The RATS Software Forum

This is an updated version of the SeasonalDLM procedure, which creates one of two state space representations for a seasonal. The main change has been to switch to the use of named options (A, C, and F) rather than the unnamed parameters, which will make it easier to use. (The old syntax still works, however). Also, it now defines the "F" matrix (mapping from fundamental shocks to states), rather than the SW matrix. This also will make it easier to use.

Code: Select all

*
* @SeasonalDLM(options)
*
* Creates the "A", "C" and "F" matrices for the seasonal component of a DLM.
*
* Options:
*   TYPE=[FOURIER]/ADDITIVE
*   SPAN=seasonal span [calendar seasonal]
*   A (output)    the state transition matrix
*   C (output)    the coefficient matrix
*   F (output)    the loadings from the state shocks to the states
*
* TYPE=FOURIER models the seasonal states as cosine and sine waves at the harmonic
* frequencies of the seasonal (2 x pi x j/S), j=1,...,S/2, with no sine component
* when j=S/2 for S even. The "C" vector adds up the cosine components. The
* amplitudes are independent random walks.
*
* TYPE=ADDITIVE models the seasonal states as additive components for each month which
* approximately add up to zero over the course of one year. The "C" vector extracts
* just the current month.
*
* See Durbin and Koopman, "Time Series Analysis by State Space Methods", Oxford
* University Press 2001, pp 40-42; West and Harrison, "Bayesian Forecasting and
* Dynamic Models" 2nd ed, Springer 1997, chapter 8 for more information.
*
* Parameters (still supported, in order):
*   A (output)    the (span-1) x (span-1) state transition matrix
*   C (output)    the (span-1) x 1 coefficient matrix
*   SW (output)   the (span-1) x (span-1) symmetric state disturbance matrix
*                  (shape only - you'll have to multiply by the component variance)
*
* Revision Schedule:
*   02/2007 Written by Tom Doan, Estima
*   07/2009 Revised to use options rather than parameters and to generate "F"
*           (loadings from shocks to states) rather than SW.
*
procedure SeasonalDLM ap cp swp
type rect *ap *cp
type symm *swp
*
option integer span
option choice  type  1  fourier additive
option rect    *a
option rect    *c
option rect    *f
*
local integer spanl i j k
local real    freq
local rect    al cl fl
*
if .not.%defined(span)
   inquire(seasonal) spanl
else
   compute spanl=span

if type==1 {
   *
   * Fourier model.
   *
   compute al=%zeros(spanl-1,spanl-1)
   *
   * Fill in the matched pairs
   *
   do j=1,(spanl-1)/2
      compute k=j*2-1,freq=2*%pi*j/spanl
      compute al(k,k)=cos(freq),al(k,k+1)=sin(freq),al(k+1,k)=-sin(freq),al(k+1,k+1)=cos(freq)
   end do j
   *
   * Now do the cosine only (frequency pi)
   *
   if %clock(spanl,2)==2
      compute al(spanl-1,spanl-1)=-1.0
   *
   * Generate C and F
   *
   dim cl(spanl-1,1)
   ewise cl(i,j)=%clock(i,2)==1
   compute fl=%identity(spanl-1)
}
else {
   dim al(spanl-1,spanl-1)
   ewise al(i,j)=%if(i==1,-1,i==j+1)
   compute cl=%unitv(spanl-1,1)
   compute fl=%unitv(spanl-1,1)
}
*
* Save whatever was requested
*
if %defined(a)
   compute a=al
if %defined(ap)
   compute ap=al
if %defined(c)
   compute c=cl
if %defined(cp)
   compute cp=cl
if %defined(f)
   compute f=fl
if %defined(swp)
   compute swp=%outerxx(fl)
end