状态空间模型计算和例子

1965
1970
1975
1980
Series 2
0
5
10
1955
1960
1965
1970
1975
1980
#library(dse) my2=TSdata(output=newdat[,3,drop = F]) my2=tframed(my2, list(start=c(1953,1), frequency=4))
library(dse) my=TSdata(input= newdat[,3,drop = F],output=newdat[, 2,drop = F]) my=tframed(my, list(start=c(1953,1), frequency=4)) seriesNamesInput(my)="CPI" seriesNamesOutput(my)="BILLS"
seriesNamesInput(my)="CPI" seriesNamesOutput(my)="BILLS“ bb=estBlackBox(my,g=3) #lag=?!!! tfplot(bb) bb$model #Giving F, G, H, K bb$estimates #giving $cov, $like, $pred attributes(bb)#model #看结果
Step 4: Kalman Filter Filtering, Smoothing, and Forecasting
• FKF: using R function fkf • Shumway using functions Kfilter0, Ksmooth0, Kfilter1, Ksmooth1, Kfilter2, Ksmooth2 • dlm: using function dlmFilter, dlmSmooth, dlmFilterdlmForecast • dse: using functions such as forecast, featherForecasts, horizonForecasts, etc Note: you have to pay attention to the meanings of the output!
Step 2: Log Likelihood Function
• FKF: write function such as objective to create log likelihood function • Shumway: write function via Kfilter0, Kfilter1, Kfilter2 to create log likelihood function • dlm: write function such as buildFun
状态空间模型及例子
吴喜之
Models and Computation Steps for State-Space Model
Sources (for four packages): Package: FKF, related article (not for R, but having formulas): /p/dgr/kubcen/1998141.html (file3826.pdf) Package: dlm dlm.pdf (in R project or your R directory) Package: Shumway: /stoffer/tsa2/chap6.htm Package: dse In your R directory : dse-guide.pdf and article: http://www.bank-banque-canada.ca/pgilbert/gil93.pdf
ans <- fkf(a0 = c(0, 0), P0 = diag(c(10, 10)), dt = rbind(0, 0), ct = matrix(0), Tt = matrix(c(ar1, 1, ar2, 0), ncol = 2), Zt = cbind(1, 0), HHt = matrix(c(sigma^2, 0, 0, 0), ncol = 2), GGt = matrix(0), yt = rbind(a))
#library(dse) my1=TSdata(output=newdat[,c(3,2),drop = F]) my1=tframed(my1, list(start=c(1953,1), frequency=4)) bb1=estBlackBox(my1,g=3) tfplot(bb1) attributes(bb1)#model
Series 1
One step ahead predictions (dotted) and actual data (solid)
12
bb1$model #Giving F, H, K bb1$estimates #giving $cov, $like, $pred #看结果
0
4
8
1955
1960
Filtering, Smoothing, Forecast
• s<t: forecasting or prediction • s=t: filtering • s>t: smoothing
Kalman Filter (Forecast ): with
Prediction equations
with
0 2
4
6
bb2$model #Giving F, H, K (one dimention) bb2$estimates #giving $cov, $like, $pred
8
10
12
14
1955
1960
1965
1970
1975
1980
bb4=estVARXls(my1,g=3)#two variables (TWO output) bb4$model #Giving A(L), B(L), C(L) bb4$estimates #giving $cov, $like, $pred
Model (Shumway notation)
Transition equation and observation equation
• • • • • • •
xt : p×1 state vector N(m, S) yt : q×1 observation/measurement vector wt : p×1 system noise iid N(0,Q) vt : q×1 observation noise iid N(0,R) ut : r×1 fixed input At : q×p observation/measurement matrix F : p×p ; G : q×r, U: p×r, Q: p×p, R: q×q
updating equations where
The MSE of
is Kalman gain
Kalman Smoother: with initial
The Lag-one Covariance Smoother
With initial condition
for t=n,n-1,…,2
Three “levels” of Shumway’s models
One step ahead predictions (dotted) and actual data (solid)
bb2=estBlackBox(my2,g=3) tfplot(bb2) attributes(bb2)#model
Series 1
bb3=estVARXls(my,g=1)#two variables (one input one output) bb3$model #Giving A(L), B(L), C(L) bb3$estimates #giving $cov, $like, $pred
– dlmModARMA , dlmModPoly, dlmModReg , dlmModSeas , dlmModTrig
• dse: using functions such as ARMA, SS (can transfer to each other) to build an empty model (but)
Example: FKF
Example: FKF


#Shumway Example 6.12 newdat = read.table("f:/xzwu/2009/berkeley/ts/shumway/mydata/Newbold1.txt", comment.char="#") y = newdat[1:50,2] # quarterly inflation (column 2) z = newdat[1:50,3] # quarterly interest rates (column 3) plot(ts(y,start=c(1953,1),frequency=1)) points(ts(y,start=c(1953,1),frequency=1),type='p',pch=15) lines(ts(z,start=c(1953,1),frequency=1),lty=2) points(ts(z,start=c(1953,1),frequency=1),type='p',pch=18)
One step ahead predictions (dotted) and actual data (solid) bb=estBlackBox(my,g=3) #lag=?!!! tfplot(bb)
15 Series 1 0 5 10