Post-Sample predictive test in State Space Model

[Jennylai]

Dear Tom,

May I ask you another question about State space model? In Harvey(1989)'s book, on page 270 there is a test called " post-sample predictive testing and model evaluation". May I know if RATS has a built-in test for this "post-sample predictive test"

Thank you very much for you kind help all the time!

[Jennylai]

Dear Tom,

I am now doing the predictive testing for State Space Model. May I ask for your advice and help? The thing is this: I chop my sample into two subsamples: sub1 for insample estimation using DLM, and sub2 for out-of-sample predictive testing (which is decribed by Harvey(1989) on P271). What I want to do is to compare the out-of-sample one-step-ahead prediction error (vps1) and the insample one-step-ahead prediction error (vt1) in a formula like this:

vps1=vps1/sqrt(svps1), vt1=vt1/sqrt(gt(1,1)) (standardize those two groups of prediction errors so that they will have the distribution of NID(0, sigmastar)).

epslon(l) = ([sum from j=1 to l (vps1(j)**2)]/l)/[sume from i=d+1 to T(vt1(i)**2)/(T-d)]
and epslon(l) has an F(l,T-d) distribution under the null.

I take reference from durkp162.r and the programme is as follows:

"
cal(m) 1999
open data FDI_EXPV.xls
data(format=xls,org=columns) 1999:01 2009:03 FDI AR1001MB AR1001MP AR1001MS NDF1M

...... (the middle part are the same as in durkp162.r, except that I use 'smpl' to get my sub1 period)
*
* We will 'smpl' the orginal series FDI to make way for insample estimation and out-of-sample
* predictive test. Here we use time period 1999:01 to 2008:01 for estimation.
*
smpl 1999:01 2008:01
*
dlm(method=bfgs,start=(sw=exp(2*phil)*swl~\exp(2*phis)*sws),$
y=FDI,a=a,c=c,exact,scale,sv=1.0,sw=sw,type=filter,vhat=vt,svhat=gt) / xfilter
compute sigmastar=%variance

set gt1 = gt(t)(1,1)*sigmastar
set vt1 = vt(t)(1)
set ehat = vt1/sqrt(gt1)

* Now we do the post sample predictive test for this model 2 + seasonality.
*(Here I basically want to use DLM and the estimated parameters I obtained from the above "DLM" command line to do the algorithm to obtain the X(T+j/T), Y(T+j/T), and one-step ahead prediction error during this prediction period. T+j/T reads as T+j conditional on time T)

smpl 2008:02 2009:03

dlm(y=FDI,a=a,c=c,exact,sv=sigmastar,sw=sw,x0=xsmooth(2008:01),sx0=psmooth(2008:01),$
type=filter,vhat=vps,svhat=svps) /

*
* calculate the out-of-sample one-step-ahead prediction error
*
set svps1 = svps(t)(1,1)
set vps1 = vt(t)(1)

(Exactly below this line, RATS returns an error message: ## MAT15. Subscripts Too Large or Non-Positive)

As a result, may I ask if you would give me some advice or idea of how to fix this problem?

Thank you very much for your kind help all the time!

[Jennylai]

Dear Tom,

I've solved the error that RATS reports in my above post by using :
"
smpl 2008:02 2009:03

dlm(y=FDI,a=a,c=c,sv=sigmastar,sw=sw*sigmastar,x0=xsmooth(2008:01),sx0=psmooth(2008:01),$
type=filter,vhat=vps,svhat=svps) 2008:02 2009:03

set svps1 = svps(t)(1,1)
set vps1 = vps(t)(1)
set eohat = vps1/sqrt(svps1)
"
. Thank you all the same.
And now my biggest concern is: Can I use this DLM to produce the out-of-sample one step ahead prediction error V(T+j/T) ?
Is it right?

Thank you! Waiting for your reply.

[TomDoan]

That looks like you have that correct. You run the Kalman filter starting at the final state and variance from the estimation sample. An alternative which you might find easier is to use

Code: Select all

dlm(method=bfgs,start=(sw=exp(2*phil)*swl~\exp(2*phis)*sws),$
y=FDI,a=a,c=c,exact,scale,sv=1.0,sw=sw,type=filter,vhat=vt,svhat=gt) 1999:01 2008:01 xfilter
compute sigmastar=%variance

set gt1 = gt(t)(1,1)*sigmastar
set vt1 = vt(t)(1)
set ehat = vt1/sqrt(gt1)

dlm(method=solve,start=(sw=exp(2*phil)*swl~\exp(2*phis)*sws),$
y=FDI,a=a,c=c,exact,scale,sv=1.0,sw=sw,type=filter,vhat=vt,svhat=gt) 1999:01 2009:02 xfilter

That way you just use the same DLM instruction - you estimate (METHOD=BFGS) on the initial sample, and just filter (METHOD=SOLVE) over the full sample.

[Jennylai]

Dear Tom,

Thank you very much for your kind reply! It helps a lot!

[TomDoan]

This is an example from Harvey which does the post-sample predictive test (among other things). As mentioned earlier, you do that by estimating through the initial range, then doing the same DLM instruction through the full range with METHOD=SOLVE rather than one of the estimation options. The test is constructed by comparing the sums of squared recursive residuals over the estimation period and predictive period.

This requires the updated versions of LocalDLM and SeasonalDLM which you can get from the RATS Procedures forum. The data file is attached.

Code: Select all

*
* Example from pp 82-84
*
open data seatbelt.prn
calendar(m) 1969
data(format=prn,org=columns) 1969:01 1984:12 drivers front rear kms petrol rate law
*
@localdlm(type=trend,a=at,c=ct,f=ft,shocks=both)
@seasonaldlm(type=additive,a=as,c=cs,f=fs)
*
compute a=at~\as,c=ct~~cs,f=ft~\fs
*
* The success of the estimation is quite sensitive to the initial guess
* values. The following is a bit cruder than the guess values computed
* using the reference in the book, but should work in practice
*
* (a) get rid of the seasonal by regression. (This corresponds to
*     a zero variance for that component)
*
filter(remove=seasonal) drivers / deseas
*
* (b) do a short span filter on the deseasonalized data to get an
*     estimate of the local trend
*
filter(type=centered,width=5) deseas / trend
*
* (c) compute the approximate irregular component and its variance
*
set irreg = deseas-trend
stats(noprint) irreg
compute vi=%variance
*
* (d) compute the first difference of the local level component and
*     its variance
*
set difft = trend-trend{1}
stats(noprint) difft
compute vl=%variance
*
* Because we're concentrating out the irregular variance, the <<phil>>
* parameter needs to be initialized based upon the ratio. We'll make
* <<phis>> small.
*
nonlin sigsqeps sigsqeta sigsqzeta sigsqomega
compute sigsqeps=vi,sigsqeta=vl,sigsqzeta=vl*.1,sigsqomega=vl*.1
*
function %%GlueSWBSM
type symm %%GlueSWBSM
compute %%GlueSWBSM=%diag(||sigsqeta,sigsqzeta,sigsqomega||)
end
*
nonlin sigsqeps sigsqeta sigsqzeta=0.00 sigsqomega=0.00
dlm(method=bfgs,start=(sw=%%GlueSWBSM()),type=filter,$
   y=drivers,a=a,c=c,f=f,exact,var=known,sv=sigsqeps,sw=sw,$
   vhat=vhat,svhat=svhat) * 1981:12 xfilter
*
* Compute recursive residuals
*
set resids * 1981:12 = vhat(t)(1)/sqrt(svhat(t)(1,1))
@regcorrs(number=15,dfc=1,qstats) resids
*
report(action=define)
report(atrow=1,atcol=1,align=center) "J" "F" "M" "A" "M" "J" "J" "A" "S" "O" "N" "D"
*
* The seasonals are in "reversed" order. The first seasonal state for 1981:12
* is for December, the second for November, etc. January is - the sum of the
* other months.
*
do i=2,12
   report(atrow=2,atcol=i) exp(xfilter(1981:12)(15-i))
   report(atrow=2,atcol=1) exp(-%sum(%xsubvec(xfilter(1981:12),3,13)))
end do i
report(action=format,picture="*.###")
report(action=show)
*
* Post sample predictive test. Run the filter through 1982:12 using the
* coefficients estimated through 1981:12. (Same instruction with
* METHOD=SOLVE replacing METHOD=BFGS)
*
dlm(method=solve,start=(sw=%%GlueSWBSM()),type=filter,$
   y=drivers,a=a,c=c,f=f,exact,var=known,sv=sigsqeps,sw=sw,$
   vhat=vhat,svhat=svhat) * 1982:12 xfilter
*
set resids * 1982:12 = vhat(t)(1)/sqrt(svhat(t)(1,1))
*
* Residuals aren't defined until 1970:2 because of the 13 diffuse states
*
sstats 1970:2 1981:12 resids^2>>base
sstats 1982:1 1982:12 resids^2>>forecast
cdf(title="Post-sample predictive test") ftest (forecast/12)/(base/143) 12 143
*
* Re-estimate through 1982:12
*
dlm(method=bfgs,start=(sw=%%GlueSWBSM()),type=filter,$
   y=drivers,a=a,c=c,f=f,exact,var=known,sv=sigsqeps,sw=sw,$
   vhat=vhat,svhat=svhat) * 1982:12 xfilter
*
* Do diagnostics again
*
set resids * 1982:12 = vhat(t)(1)/sqrt(svhat(t)(1,1))
@regcorrs(number=15,dfc=1,qstats) resids
*
* Estimate the "airline" model - doesn't really work
*
boxjenk(diffs=1,sdiffs=1,ma=1,sma=1,maxl) drivers * 1982:12

[Jennylai]

Dear Tom,

Thank you very much for your help all the time! It helped me a lot!