LakeHuron
plot.ts ( LakeHuron )

  # load the dlm package
library(dlm)

  # create a locally constant model, the simplest dlm
mod1 <- dlmModPoly ( order=1 )
mod2 <- dlmModPoly ( order=2 )
  # show what it is
mod1
names(mod1)
unlist(mod1)
FF(mod1)
mod1$FF
GG(mod1)
mod1$GG

  # filter
filt1 <- dlmFilter ( LakeHuron, mod1 )
  # show what it is
names(filt1)
filt1
filt1$y
filt1$mod
filt1$m
length(filt1$m)
cbind(LakeHuron,filt1$m[2:99])
dlmSvd2var ( filt1$U.C, filt1$D.C )

filt1$a
dlmSvd2var ( filt1$U.R, filt1$D.R )

  # how well does the model fit?
    # graphically
plot ( LakeHuron, type="p", xlab="Year", ylab="Depth",
       main="Lake Huron" )
lines ( 1875:1972, filt1$m[-1] )
var <- dlmSvd2var(filt1$U.C, filt1$D.C)
sd <- sqrt(unlist(var))
lines ( 1875:1972, filt1$m[-1] + 2*sd[-1], lty=3, col="red" )
lines ( 1875:1972, filt1$m[-1] - 2*sd[-1], lty=3, col="red" )

points ( 1875:1972, filt1$a, lty=2, col="green" ) # forecast

    # numerically
dlmLL ( LakeHuron, mod1 ) # negative loglikelihood


  # do V and W matter?
mod2 <- dlmModPoly ( order=1, dV=10, dW=1 )
filt2 <- dlmFilter ( LakeHuron, mod2 )
mod3 <- dlmModPoly ( order=1, dV=1, dW=10 )
filt3 <- dlmFilter ( LakeHuron, mod3 )
plot ( LakeHuron, type="p", xlab="Year", ylab="Depth",
       main="Lake Huron" )
lines ( 1875:1972, filt1$m[-1], col="red" )
lines ( 1875:1972, filt2$m[-1], col="blue" )
lines ( 1875:1972, filt3$m[-1], col="green" )

dlmLL ( LakeHuron, mod1 )
dlmLL ( LakeHuron, mod2 )
dlmLL ( LakeHuron, mod3 )

  # find V and W by maximum likelihood
buildit <- function ( par ) {
  mod <- dlmModPoly ( order=1, dV=exp(par[1]), dW=exp(par[2]) )
  return ( mod )
}
dlmMLE ( LakeHuron, c(1,1), buildit )
  # try different starting values
mle1 <- dlmMLE ( LakeHuron, c(1,1), buildit )
exp(mle1$par)
modmle <- dlmModPoly ( order=1, dV=exp(mle1$par[1]),
                       dW=exp(mle1$par[2]) )
dlmLL ( LakeHuron, modmle )
filtmle <- dlmFilter ( LakeHuron, modmle )
lines ( 1875:1972, filtmle$m[-1], col="yellow" )

# still to do:
#   - smoothing
#   - locally linear model
#   - Bayesian treatment of V and W
#   - more data sets
