Library: | kalman |
See also: | gkalarray gkalfilter gkallag gkalresiduals |
Quantlet: | gkalsmoother | |
Description: | calculates a smoothed time series for a state space model (uni- or multivariate) using the Kalman smoother. The quantlet gkalsmoother needs a pre-run of gkalfilter. |
Usage: | gkalsmoothOut = gkalsmoother(Y,Ta,Ra,gkalfilOut {,ca}) | |
Input: | ||
Y | N_max x TIME matrix, observed time series; N_max is the maximal number of variables observed at any time t | |
Ta | K x K x TIME array, observations for T_t | |
Ra | K x K x TIME array, covariance matrices R_t | |
gkalfilOut | K x (K+1) x (TIME+1) array, output of gkalfilter | |
ca | optional K x 1 x TIME array, state constants c_t; default ca = 0*Ta[,1,] | |
Output: | ||
gkalsmoothOut | K x (K+1) x (TIME+1) array, smoothed state space series a_(t|T)~P_(t|T). First entry is for t=0. |
State equation
alpha_t = c_t + T_t*alpha_t-1 + e^s_t
Measurement equation
y_t = d_t + Z_t*alpha_t + e^m_t
with
alpha_0 ~ (mu,Sig), e^s_t ~ (0,R_t), e^m_t ~ (0,H_t)
All parameters are known.
We denote with as_(t|T) the filtered series of the unobservable state vector alpha_t and with Ps_(t|T) the corresponding covariance matrix. The output of the procedure is an array with as_(t|T) and Ps_(t|T) for all t.
To run the smoother, it is necessary to run gkalfilter first. This produces the array gkalfilOut.
library("kalman") library("xplore") ; loads the quantlets from xplore library library("plot") ; loads the quantlets from plot library Y = read("houseprice.dat") Aa = read("houseAa.dat") dimAa = read("housedimAa.dat") Aa = reshape(Aa',dimAa) Pa = read("housePa.dat") dimPa = read("housedimPa.dat") Pa = reshape(Pa',dimPa) gkalfilOut = Aa~Pa T =(#(1.4,1)'|#(-0.4,0)')~0*matrix(2,3)|(0*matrix(3,2)~unit(3)) Ta = gkalarray(Y,T,0,0) Ra = gkalarray(Y,diag(#(0.4,0,0,0,0)),0,0) gkalsmoothOut = gkalsmoother(Y,Ta,Ra,gkalfilOut) Time =(1:cols(Y)) T = rows(Time) state = reshape(gkalsmoothOut[1,1,2:(T+1)],#(1,T,1))' state = Time~state state = setmask(state,"line","red","medium") Y1 = Time~(Y[1,]') setmaskp(Y1,0,8,5) Y2 = Time~(Y[2,]') setmaskp(Y2,0,8,5) Y3 = Time~(Y[3,]') setmaskp(Y3,0,8,5) disp = createdisplay(1,1) show(disp,1,1,state,Y1,Y2,Y3) setgopt(disp,1,1,"title","Kalman smoother","xlabel","Time") setgopt(disp,1,1,"ylabel","y,a","border",0)
The first component of the smoothed state vector is shown as a red line; the observations of Y are shown as black dots.