# Cut-and-Paste Code below into Window Above and Run
#
# W_Malthus (0-2000)
#
require(dse)
require(matlab)
# Measurement Matrix (Growth) (Q-N)
# Q N
#[1,] 0.707 0.707
#[2,] 0.707 -0.707
#
# Fraction of Variance
#[1] 0.96 1.00
#
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]}
f <- matrix( c( 1.023336956, 0.008819228, 0.0069498915,
0.006613558, 1.012095445, 0.0009444344,
0.0000000, 0.0000000, 1.0000000
),byrow=TRUE,nrow=3,ncol=3)
#
# Create a Random Walk (uncomment next line)
# f[1,1] <- f[2,2] <- 1
# Stabilize (uncomment next line)
# f[1,1] <- f[2,2] <- .99
# Classical Malthus Model
# f[2,1] <- 0; f[2,2] <- 1;
#
h <- eye(2,3)
k <- f[1:3,1:2,drop=FALSE]
Malthus <- SS(F=f,H=h,K=k,z0=c( 0.0069498915, 0.0009444344, 1.0000000),
output.names=c("Growth","Q-N"))
stability(Malthus)
#shockDecomposition(toSSChol(Malthus))
#tfplot(simulate(Malthus,sampleT=2000,noise=matrix(0,2000,2),start=1))
Malthus.data <- simulate(Malthus,sampleT=2000,start=1)
m <- l(Malthus,Malthus.data)
#tfplot(m)
Malthus.f <- forecast(m,horizon=500)
tfplot(Malthus.f)
AIC(m)
MalthusX<- SS(F=f,H=h,Q=eye(3,2),R=eye(2,2),z0=c( 0.0069498915, 0.0009444344, 1.0000000),
output.names=c("Growth","Q-N"))
shockDecomposition(MalthusX)
#
# LCI Parameter UCI P>=T[1] P< T[1] Std. Dev. Bias Bias-z
# [1,] 1.0197290 1.0233370 1.028431 0.34 0.66 0.0033912 0.00101 0.2979
# [2,] 0.0054137 0.0066136 0.008318 0.33 0.67 0.0010994 1.01704 925.0700
# [3,] -0.0037078 0.0088192 0.025156 0.32 0.68 0.0106613 1.01916 95.5944
# [4,] 1.0079072 1.0120954 1.016778 0.33 0.67 0.0032656 0.01255 3.8423
# [5,] 0.0062602 0.0069499 0.007678 0.42 0.58 0.0006533 1.01647 1555.8314
# [6,] 0.0007261 0.0009444 0.001218 0.34 0.66 0.0001961 1.02243 5214.0136
#