# Cut-and-Paste code below into window above and Run
#
# W_LM (Late Modern) Model
#
# Measurement Matrix (Growth) (L-T-X-XREAL), (T+L-XREAL-X)
#
# Q N XREAL X L T
#[1,] 0.4343 0.433 0.418 0.424 0.357 0.375
#[2,] 0.0827 0.122 -0.291 -0.244 0.815 -0.413
#[3,] -0.1689 -0.162 -0.386 -0.226 0.254 0.825
#
# Fraction of Variance
#
periodsOutput <-
function (obj,...) {
return(dim(obj$data$output)[1])
}
attractorPath <- function(obj,data,series=1) UseMethod("attractorPath")
#
# Plot the Attractor Path for a TSestModel
# series = output series to plot
#
attractorPath.SS <- function(obj,data,series=1) {
obj1 <- l(obj,data)
n <- periodsOutput(obj1)
m <- nseriesOutput(obj1)
tstart=start(obj1)
obj2 <- simulate(obj1,noise=matrix(0,n,m),start=tstart)
tfOnePlot(tbind(outputData(obj1,series=series),outputData(obj2,series=series)))
}
phaseSpace <- function(obj,data,n=500) UseMethod("phaseSpace")
phaseSpace.SS <-
function(obj,data,n=500) {
require(scatterplot3d)
obj1 <- l(obj,data)
m <- dim(outputData(data))[2]
model <- simulate(obj1,sampleT=n,noise=matrix(0,n,m))
data <- outputData(obj1)
m <- dim(data)[2]
if (m > 3) m <- 3
if (m == 3) scatterplot3d(data[,1:m], type="l") else if (m == 2) plot(data[,1],data[,2],xlab=seriesNames(data)[1],ylab=seriesNames(data)[2], type="l") else stop("no phase space")
}
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]}
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.00713465, 0.1136613, -0.3071767, 0.15009523,
-0.03424590, 0.9703990, 0.3867074, -0.02887742,
-0.03421009, 0.2132679, 0.3331558, -0.00712312,
0.00000000, 0.0000000, 0.0000000, 1.0000000000
),byrow=TRUE,nrow=4,ncol=4)
#
# To stabilize the model, uncomment the next line:
# f[1,1] <- 0.92932816; f[2,2] <- 0.8954306; f[3,3] <- 0.3074178
#
h <- eye(3,4)
k <- (f[,1:3,drop=FALSE])
W_LM <- SS(F=f,H=h,K=k,z0=c(0.15009523, -0.02887742, -0.0071231, 1.0000000000),
output.names=c("W1","W2","W3"))
print(W_LM)
is.SS(W_LM)
stability(W_LM)
#W_LM.data <- (simulate(W_LM,sampleT=150,start=1950))
W_LM.data <- simulate(W_LM,sampleT=150,noise=matrix(0,20,3),start=1950)
W_LM.f <- forecast(W_LM,W_LM.data,horizon=150)
tfplot(W_LM.f)
W_LM.fx <- merge.forecast(W_LM.f)
W_LMx <- SS(F=f,H=h,Q=eye(4,3),R=eye(3,3),z0=c(0.15009523, -0.02887742, -0.0071231, 1.0000000000),
output.names=c("W1","W2","W3"))
attractorPath(W_LM,W_LM.data)
phaseSpace(W_LM,W_LM.data)
shockDecomposition(W_LMx)