###
### simple example of estimation of a Structural Time Series Model
### via the state space apparatus. Note that:
###
### (a) this example supersedes the code discussed in Lucchetti (2011), 
### State Space Methods in gretl, Journal of Statistical Software, 041, 
### issue i11, that uses the same data, but an obsolete syntax, which 
### is no longer supported
###
### (b) this script only has instructional purposes: the StrucTiSM 
### package provides the model below and several others.
###

# open the dataset

open nile.gdt

# set up the local linear trend model with the Nile data

Nile = ksetup({nile}, {1}, {1}, {1})
Nile.diffuse = 2
Nile.obsvar = {1}

# initialize the parameters via the reduced form
# \Delta y_t = \eta_t + \Delta \epsilon_t

series dnile = diff(nile)
scalar    vd = var(dnile)
scalar  veps = -cov(dnile,dnile(-1))
scalar s_eps = sqrt(veps)
scalar s_eta = sqrt(vd - 2*veps)

# estimate the variances via MLE

scalar err = 0
mle loglik = err ? NA : Nile.llt
    matrix Nile.statevar = s_eta^2
    matrix Nile.obsvar = s_eps^2
    err = kfilter(&Nile)
    params s_eta s_eps
end mle

# extract the smoothed trend anb its estimated variance

ksmooth(&Nile)
series nile_trend = Nile.state
series nile_trend_se = sqrt(Nile.statevar)

# plot series with trend and associated 95% confidence band

list X = nile nile_trend
plot X
    options with-lines time-series
    options band=nile_trend,nile_trend_se,1.96 band-style=fill,0xdddddd
    literal set key right top
    printf "set title \"%s\"", "Nile data and trend with 95\% confidence band"
end plot --output=display