Examples

The freely downloadable SsfPack Basic contains many examples in the code subdirectory. Documentation can be found in the file SsfPackV22.pdf in the doc subdirectory.

The document Illustrations of state space modelling using SsfPack contains examples of (mainly financial) times series modelling using SsfPack. These examples are taken from Chapter 11 of Analysis of Financial Time Series by Ruey S. Tsay. The document contains Ox versions of the S-Plus code in the book. The file Ox-FinancialTS.zip contains the Ox code documented in the PDF file.

Examples can also be found in State Space Methods in Ox/SsfPack by Matteo M. Pelagatti (including an introduction in SsfPack) and in A Bayesian Analysis of Unobserved Component Models Using Ox by Charles S. Bos, see also Statistical Software for State Space Methods by Jacques J.F. Commandeur, Siem Jan Koopman, and Marius Ooms.

More Examples

The following examples demonstrate the use of the SsfPack filtering, smoothing and simulation functions in Ox. The examples all consider the local linear "smooth" trend or integrated random walk model, applied to the series of nine observations: 1,9,2,5,8,4,6,7,3. More details can be found in the SsfPack documentation.

Kalman filter

ssfkf.ox

#include < oxstd.h >
#include < packages/ssfpack/ssfpack.h >

main()
{
   decl mphi = < 1,1;0,1;1,0 > ;
   decl momega = diag(< 0,0.1,1 >);
   decl msigma = < 0,0;0,0;1,.5 > ;
  
   decl mr = sqrt(momega) * rann(3, 21);
   decl md = SsfRecursion(mr, mphi, momega, msigma);
   decl myt = md[2][1:];    // 20 observations
  
   decl mkf = KalmanFil(myt, mphi, momega);
   print("mKF\' (t=10)", "%c", {"v", "K(1,1)", "K(2,1)", "F^-1"},
         mkf[][9]');

}

Output:

mKF' (t=10)

          v     K(1,1)     K(2,1)       F^-1

    -3.0347    0.76491    0.21161    0.44669

Smoothing

ssfsmo.ox

#include < oxstd.h >
#include < packages/ssfpack/ssfpack.h >

main()
{
   decl mphi = < 1,1;0,1;1,0 > ;
   decl momega = diag(< 0,0.1,1 >);
   decl msigma = < 0,0;0,0;1,.5 > ;
   
   decl mr = sqrt(momega) * rann(3, 21);
   decl md = SsfRecursion(mr, mphi, momega, msigma);
   decl myt = md[2][1:];    // 20 observations

   decl mkf = KalmanFil(myt, mphi, momega);
   print("mKF\' (t=10)", "%c", {"v", "K(1,1)", "K(2,1)", "F^-1"},
         mkf[][9]');

   decl mks = KalmanSmo(mkf, mphi, momega);
   print("Basic smoother output: mKS\' (t=10)", "%c",
         {"r","e(1,1)","e(2,1)","N","D(1,1)","D(2,2)"}, mks[][10]');

   decl msmodist = mks[0:2][0] ~ momega * mks[0:2][1:];
   print("Smoothed disturbances (t=10)", "%c",
         {"E[H.eps](1,1)","E[H.eps](2,1)","E[G.eps]"}, msmodist[][10]');

   decl msmostat = SsfRecursion(msmodist, mphi, momega);
   print("Smoothed states (t=10)",	"%c",
         {"alphahat(1,1)[t+1]","alphahat(2,1)[t+1]","y[t]"}, msmostat[][10]');
}

Output:

mKF' (t=10)

          v     K(1,1)     K(2,1)       F^-1

    -3.0347    0.76491    0.21161    0.44669

Basic smoother output: mKS' (t=10)

          r     e(1,1)     e(2,1)          N     D(1,1)     D(2,2)

   -0.64966     1.3099    -1.1358    0.60208     2.0578    0.79365

Smoothed disturbances (t=10)

E[H.eps](1,1)  E[H.eps](2,1)     E[G.eps]

      0.00000        0.13099      -1.1358

Smoothed states (t=10)

alphahat(1,1)[t+1]  alphahat(2,1)[t+1]         y[t]

           0.31542            -0.27487     -0.41452

Simulation

ssfsim.ox

#include < oxstd.h >
#include < packages/ssfpack/ssfpack.h >
							
main()
{
   decl mphi = < 1,1;0,1;1,0 > ;
   decl momega = diag(< 0,0.1,1 >);
   decl msigma = < 0,0;0,0;1,.5 > ;
   
   decl md = SsfRecursion(sqrt(momega) * rann(3, 21), mphi, momega, msigma);
   decl myt = md[2][1:];
   decl mkf = KalmanFil(myt, mphi, momega);

   decl ct = columns(myt);    // 20 observations
   decl mgamma = diag(< 0,1,0 >);
   decl mwgt = SimSmoWgt(mgamma, mkf, mphi, momega);

   print("Simulation smoother weights (t=10)", "%c",
         {"A(1,1)","A(1,2)","A(1,3)","B(1,1)"}, mwgt[][9]');

   // draw 1
   md = SimSmoDraw(mgamma, rann(1, ct), mwgt, mkf, mphi, momega);
   print("Draw 1 for slope disturbances (t=10)", "%c",
         {"H.eps(1,1)","H.eps(2,1)","G.eps"}, md[][10]');
   md = SsfRecursion(md, mphi, momega);
   print("Draw 1 for state and signal (t=10)", "%c",
         {"alpha(1,1)[t+1]","alpha(2,1)[t+1]","y[t]"}, md[][10]');

   // draw 2
   decl md2 = SimSmoDraw(mgamma, rann(1, ct), mwgt, mkf, mphi, momega);
   md2 = SsfRecursion(md2, mphi, momega);

   // draw 3
   decl md3 = SimSmoDraw(mgamma, rann(1, ct), mwgt, mkf, mphi, momega);
   md3 = SsfRecursion(md3, mphi, momega);
}

Output:

Simulation smoother weights (t=10)

     A(1,1)     A(1,2)     A(1,3)     B(1,1)

   -0.40350     1.5248  -0.014001    0.24905

Draw 1 for slope disturbances (t=10)

   H.eps(1,1)   H.eps(2,1)        G.eps

      0.00000     -0.25514      0.00000

Draw 1 for state and signal (t=10)

alpha(1,1)[t+1]   alpha(2,1)[t+1]         y[t]

      1.1744           -0.78113         1.7004