# Cut-and-Paste code below into window above and Run
#
# BR20 Model
#
# Measurement Matrix (Growth) (LU- Q - CO2) (CO2 - E - N)
#
# EN.ATM.CO2E.KT NY.GDP.MKTP.KD EG.USE.COMM.KT.OE SL.TLF.TOTL.IN
#[1,] 0.409 0.409 0.411 0.4135
#[2,] -0.245 -0.380 -0.159 -0.0514
#[3,] 0.711 0.206 -0.436 -0.2619
# SP.POP.TOTL SL.UEM.TOTL.ZS
#[1,] 0.41276 0.393
#[2,] -0.00498 0.876
#[3,] -0.39564 0.191
#
# Fraction of Variance
#[1] 0.971 0.992 0.997 0.999 1.000 1.000
#
stability <- function(obj,...) UseMethod("stability")
stability.SS <- function (obj,...)
{
roots(obj, fuzz = 1e-04, randomize = FALSE)
}
SSmodel <- function(obj1,obj2) UseMethod("SSmodel")
SSmodel.TSmodel <- function (obj1,obj2) return(l(obj1,obj2))
SSinnov <- function(F. = NULL, G = NULL, H = NULL, K = NULL, z0 = NULL, input.names = NULL, output.names = NULL,...) UseMethod("SSinnov")
SSinnov.matrix <- function(F. = NULL, G = NULL, H = NULL, K = NULL, z0 = NULL, input.names = NULL, output.names = NULL,...) {
require(matlab)
if (is.null(F.) )
stop("System Matrix F must be specified as a matrix.")
m <- ncol(G)
n <- nrow(F.)
p <- nrow(H)
if (n != ncol(F.))
stop("Model F matrix must be square.")
if (is.null(H)) H <- eye(n-1,n)
if (n != ncol(H))
stop("Model H matrix have second dimension consistent with matrix F.")
if (is.null(K)) K <- F.[,1:(n-1),drop=FALSE]
model <- list(F = F., G = G, H = H, K = K,
z0 = z0, P0 = NULL, rootP0 = NULL, constants = NULL,
description = "SS innovation model")
class(model) <- c("innov", "SS", "TSmodel")
return(model)
}
SSnonInnov <- function(F. = NULL, G = NULL, H = NULL, Q = NULL,
R = NULL, z0 = NULL, input.names = NULL, output.names = NULL,...) UseMethod("SSnonInnov")
SSnonInnov.matrix <- function(F. = NULL, G = NULL, H = NULL, R = NULL, Q = NULL, z0 = NULL, input.names = NULL,
output.names = NULL,...) {
require(matlab)
if (is.null(F.) )
stop("System Matrix F must be specified.")
m <- ncol(G)
n <- nrow(F.)
p <- nrow(H)
if (n != ncol(F.))
stop("Model F matrix must be square.")
if (is.null(H)) H <- eye(n-1,n) # Assumes constant in F
if (n != ncol(H))
stop("Model H matrix have second dimension consistent with matrix F.")
if (is.null(Q)) Q <- eye(n,n-1)
if (is.null(R)) R <- eye(n-1,n-1)
if (nrow(Q)!=n || ncol(Q)!=(n-1)) stop("Model Q matrix dimensions inconsistent with matrix F.")
if (nrow(R)!=(n-1) || ncol(R)!=(n-1)) stop("Model R matrix dimensions inconsistent with matrix F.")
model <- list(F = F., G = G, H = H, R=R, Q=Q,
z0 = z0, P0 = NULL, rootP0 = NULL, constants = NULL,
description = "SS nonInnov Model")
class(model) <- c("nonInnov", "SS", "TSmodel")
return(model)
}
merge.forecast <- function (fx,n=1) {
x <- splice(fx$pred,fx$forecast[[n]])
colnames(x) <- seriesNames(fx$data$output)
return(x)
}
AIC <- function(model) {informationTestsCalculations(model)[3]}
require(dse)
require(matlab)
f <- matrix( c(1.005057242, -0.04445958, -0.09497071, 0.15544238,
-0.007611132, 0.83908479, -0.07477594, -0.01690898,
0.016053556, -0.12791781, 1.09030592, 0.01232127,
0.00000000, 0.0000000, 0.0000000, 1.0000000000
),byrow=TRUE,nrow=4,ncol=4)
#
# to Stabilize, Uncomment following line
# f[1,1] <- 0.899866196; f[2,2] <- 0.75126471; f[3,3] <- 0.97619260;
#
h <- eye(3,4)
k <- (f[,1:3,drop=FALSE])
#BR20 <- SSinnov(F=f,H=h,K=k,z0=c( 0.15544238, -0.01690898, 0.01232127, 1.0000000000),
BR20 <- SSinnov(F=f,z0=c( 0.15544238, -0.01690898, 0.01232127, 1.0000000000), output.names=c("BR1","BR2","BR3"))
print(BR20)
is.SS(BR20)
stability(BR20)
# tfplot(simulate(BR20,sampleT=150))
BR20.data <- simulate(BR20,sampleT=50,noise=matrix(0,50,3),start=1950)
BR20.f <- forecast(BR20,BR20.data,horizon=50)
tfplot(BR20.f)
AIC(SSmodel(BR20,BR20.data))
#BR20x <- SSnonInnov(F=f,H=h,Q=eye(4,3),R=eye(3,3),z0=c( 0.15544238, -0.01690898, #0.01232127, 1.0000000000),
BR20x <- SSnonInnov(F=f,z0=c( 0.15544238, -0.01690898, 0.01232127, 1.0000000000), output.names=c("BR1","BR2","BR3"))
shockDecomposition(BR20x)
#
# LCI Parameter UCI P>=T[1] P< T[1] Std. Dev. Bias Bias-z
# [1,] 0.992797 1.005057 1.017774 0.43 0.57 0.009276 0.001966 0.2119
# [2,] -0.015245 -0.007611 0.001641 0.62 0.38 0.006725 1.010342 150.2425
# [3,] 0.008530 0.016054 0.021060 0.18 0.82 0.004894 0.993139 202.9282
# [4,] -0.115408 -0.044460 0.036786 0.72 0.28 0.056672 1.019969 17.9977
# [5,] 0.776091 0.839085 0.904875 0.24 0.76 0.056028 0.203818 3.6378
# [6,] -0.171281 -0.127918 -0.064199 0.85 0.15 0.038658 1.095055 28.3268
# [7,] -0.270216 -0.094971 0.071573 0.31 0.69 0.130951 1.151806 8.7957
# [8,] -0.205173 -0.074776 0.121551 0.82 0.18 0.135977 0.979619 7.2043
# [9,] 1.001884 1.090306 1.198789 0.04 0.96 0.073040 0.054890 0.7515
#[10,] 0.145919 0.155442 0.164191 0.72 0.28 0.007699 0.846993 110.0123
#[11,] -0.030395 -0.016909 -0.002869 0.57 0.43 0.012043 1.019858 84.6861
#[12,] 0.001821 0.012321 0.021176 0.23 0.77 0.007321 0.998233 136.3526
#