*
*  Illustration from sectin 2.7.1
*  Filtering and smoothing with missing values
*
open data nile.dat
calendar 1871
data(format=free,org=columns,skips=1) 1871:1 1970:1 nile
*
* Create the series with missing values for two entry ranges
*
set withmiss = %if(t>=21.and.t<=40.or.t>=61.and.t<=80,%na,nile)
*
*  Get filtered estimates
*
dlm(y=withmiss,a=1.0,c=1.0,sx0=1.e+7,sv=15099.0,sw=1469.1,vhat=vhat,svhat=fhat) / xstates vstates
set a = xstates(t)(1)
set p = vstates(t)(1,1)
*
*  Get smoothed estimates
*
dlm(type=smooth,y=withmiss,a=1.0,c=1.0,sx0=1.e+7,sv=15099.0,sw=1469.1,vhat=vhat,svhat=fhat) / xstates vstates
set as = xstates(t)(1)
set ps = vstates(t)(1,1)
*
*
spgraph(footer="Figure 2.5. Filtering and smoothing output when observations are missing",vfields=2,hfields=2)
graph(header="Filtered state (extrapolation)") 2
# a 1871:1 *
# withmiss 1871:1 *
graph(header="Smoothed state (interpolation)") 2
# as 1871:1 *
# withmiss 1871:1 *
graph(header="Filtered state variance")
# p 1871:1 *
graph(header="Smoothed state variance")
# ps 1871:1 *
spgraph(done)