Models have a compact R-code in dse and can be run easily in Snippets. Note that the RUL20 model is only unstable in the feedback components $F[3,3] = 1.03$. If you stabilize the third feedback component (e.g.,$F[3,3] = 0.98$) the entire system is stable.
===========
merge.forecast <-
function (fx,n=1) {
#
# Merges a forecast with the output data
x <- splice(fx$pred,fx$forecast[[n]])
colnames(x) <- seriesNames(fx$data$output)
return(x)
}
#
#
# RUL20 Russia (1950-2000)
#
# Measurement Matrix
# EN.ATM.CO2E.KT EG.USE.COMM.KT.OE NY.GDP.MKTP.KD SL.TLF.TOTL.IN SP.POP.TOTL
#[1,] 0.28209 0.4060 0.3336 0.4343 0.3721
#[2,] -0.43993 -0.2199 -0.1359 0.1328 0.2266
#[3,] 0.06433 -0.1446 0.5944 0.0494 -0.3944
# SL.UEM.TOTL.ZS HDI EF KOF
#[1,] 0.1893 0.43400 0.2348 0.1948
#[2,] 0.4971 0.03409 -0.3950 0.5161
#[3,] -0.3374 0.22538 -0.4878 0.2470
#
#Fraction of Variance
#[1] 0.5489 0.8484 0.9677 0.9814 0.9921 0.9975 0.9988 0.9995 1.0000
require(dse)
require(matlab)
AIC <- function(model) {informationTestsCalculations(model)[3]}
f <- matrix( c( 0.976658001, 0.02127754, 0.2249071, 0.161990376,
-0.008392259, 0.97543675, -0.2319674, 0.001666581,
-0.005587383, 0.09718186, 1.0312754, 0.035854885,
0.000000000, 0.00000000, 0.0000000, 1.000000000
),byrow=TRUE,nrow=4,ncol=4)
h <- eye(3,4)
k <- f[1:4,1:3,drop=FALSE]
RUL20 <- SS(F=f,H=h,K=k,
z0=c(0.161990376, 0.001666581, 0.035854885, 1.00000000),
output.names=c("RU1","RU2","RU3"))
stability(RUL20)
#tfplot(simulate(RUL20,sampleT=20))
shockDecomposition(toSSChol(RUL20))
RUL20.data <- simulate(RUL20,sampleT=50,start=1950,noise=matrix(0,50,3))
m <- l(RUL20,RUL20.data)
#tfplot(m)
AIC(m)
tfplot(RUL20.f <- forecast(m,horizon=50))
RUL20.fx <- merge.forecast(RUL20.f)
# Models have a compact R-code in dse and can be run easily in Snippets.
# The RU20 must be run first to provide input for #ETH20.
#
# W20 ETH (Ethiopia) Russia Input
#
# Measurement Matrix
# N Q
#[1,] 0.707 0.707
#[2,] -0.707 0.707
#
#Fraction of Variance
#[1] 0.958 1.000
#
require(dse)
require(matlab)
AIC <- function(model) {informationTestsCalculations(model)[3]}
f <- matrix( c( 1.04647269, -0.1059723, 0.10155553,
0.01918655, 0.9306230, 0.01335861,
0.000000000, 0.00000000, 1.00000000
),byrow=TRUE,nrow=3,ncol=3)
g <- matrix(c(-0.003819938, -0.01406381, 0.009235749,
-0.010773904, -0.01627134, 0.017135870,
0.000000000 , 0.00000000, 0.000000000
),byrow=TRUE,nrow=3,ncol=3)
h <- eye(2,3)
k <- f[1:3,1:2,drop=FALSE]
ETH20 <- SS(F=f,H=h,K=k,z0=c(0.11991149, 0.04457381, 1.00000000),
output.names=c("Growth","(Q-N)"))
stability(ETH20)
shockDecomposition(toSSChol(ETH20))
ETH20.data <- simulate(ETH20,sampleT=50,start=1950,noise=matrix(0,50,2))
m <- l(ETH20,ETH20.data)
#tfplot(m)
AIC(m)
tfplot(forecast(m,horizon=50))
ETH20x <- SS(F=f,H=h,K=k,G=g,z0=c(0.11991149, 0.04457381, 1.00000000),
output.names=c("Growth","(Q-N)"))
ETH20x
data <- TSdata(output=outputData(ETH20.data),input=RUL20.fx)
m <- l(ETH20x,data)
#tfplot(m)
AIC(m)
shockDecomposition(m)
tfplot(forecast(m,horizon=50,conditioning.inputs=RUL20.fx))