#
# WL19D World Model (1872-1908)
#
AIC <- function(model) {informationTestsCalculations(model)[3]}
require(dse)
require(matlab)
#
# Measurement Matrix (Growth-T), (T), (XREAL-N)
#
# Q N XREAL X L T
#[1,] 0.4255 0.420 0.4192 0.428 0.422 -3.25e-01
#[2,] 0.1844 0.255 -0.0705 0.101 0.231 9.12e-01
#[3,] -0.0591 -0.376 0.8067 0.162 -0.342 2.48e-01
#
#[4,] 0.7518 0.147 -0.0360 -0.331 -0.549 -2.04e-02
#[5,] 0.0753 -0.125 -0.3862 0.819 -0.399 1.39e-16
# Fraction of Variance
#[1] 0.904 0.989 1.000 1.000 1.000 1.000
#
f <- matrix( c(0.95938093, 0.1108599, 0.3349007, 0.2038538953,
0.05563723, 0.6056210, -0.7238863, 0.0000852837,
-0.01263019, -0.1333194, 0.8623076, -0.0135689120,
0.00000000, 0.0000000, 0.0000000, 1.0000000000
),byrow=TRUE,nrow=4,ncol=4)
#
# There are two ways to stabilize the system:
# (1) Reduce Growth Rates (uncomment next line)
# f[1,1] <- 0.93; f[2,2] <- 0.59; f[3,3] <- 0.84
# (2) Moderate Global Temperature Feedback (uncomment next line)
# f[2,3] <- -.1
#
h <- eye(3,4)
k <- (f[,1:3,drop=FALSE])
WL19D <- SS(F=f,H=h,K=k,z0=c(0.2038538953, 0.0000852837, -0.0135689120, 1.0000000000),
output.names=c("W1","W2","W3"))
print(WL19D)
is.SS(WL19D)
stability(SS(F=f[1:3,1:3,drop=FALSE],H=eye(3),Q=eye(3),R=eye(3)))
# tfplot(simulate(WL19D,sampleT=100))
WL19D.data <- simulate(WL19D,sampleT=20,noise=matrix(0,20,3),start=1872)
WL19D.f <- forecast(l(WL19D,WL19D.data),horizon=50)
tfplot(WL19D.f)