[R] Setting up a State Space Model in dlm

Gavin Simpson gavin.simpson at ucl.ac.uk
Fri Jun 10 09:22:52 CEST 2011


On Tue, 2011-06-07 at 17:24 +0100, Michael Ash wrote:
> This question pertains to setting up a model in the package "dlm"
> (dynamic linear models,
> http://cran.r-project.org/web/packages/dlm/index.html

The author of the dlm package has just published a paper on state space
models in R including details on setting up dlm:

http://www.jstatsoft.org/v41/i04

That might help with your question - I haven't seen a reply on list, but
am unable to help answer it either.

HTH

G

> I have read both the vignette and "An R Package for Dynamic Linear
> Models" (http://www.jstatsoft.org/v36/i12/paper), both of which are
> very helpful.  There is also some discussion at
> https://stat.ethz.ch/pipermail/r-help/2009-May/198463.html
> 
> I have what I think is a relatively straightforward state-space model
> but am unable to translate it into the terms of dlm.   It would be
> very helpful to get a basic dlm setup for the problem and I would
> guess that I can then modify it with more lags, etc., etc.
> 
> The main equation is
> pi[t] = a * pi[t-1] + b*(U[t] - UN[t]) + e[t]
> 
> (see http://chart.apis.google.com/chart?cht=tx&chl=%5Cpi_t=a%5Cpi_{t-1}%2bb%28U_t-U^N_{t}%29%2Be_t
> for a pretty version)
> 
> with pi and U observed, a and b fixed coefficients, and e a
> well-behaved error term (gaussian, say, variance unknown).
> The object of interest is the unobserved and time-varying component UN
> which evolves according to
> 
> UN[t] = UN[t-1] + w[t]
> 
> (see http://chart.apis.google.com/chart?cht=tx&chl=U%5EN_%7Bt%7D%20=%20U%5EN_%7Bt-1%7D%20%2B%20%5Cepsilon_t
> for a pretty version)
> that is, a random walk with well-behaved error term with var(w) known
> (or assumed).
> 
> I'm interested in the estimates of a and b and also in estimating the
> time series of UN.
> 
> Note that the term b*(U[t] - UN[t]) makes this a nonlinear model.
> 
> Below is code that does not work as expected.  I see the model as
> having four parameters, a, b, var(e), and UN.  (Or do I have a
> parameter UN[t] for every period?)
> 
> I do not fully understand the dlm syntax.  Is FF specified properly?
> What should X look like?  How does m0 relate to parm()?
> 
> 
> I would be grateful if someone would be willing to glance at the code.
>  Thanks. Michael
> 
> library(quantmod)
> library(dlm)
> ## Get and organize the data
> getSymbols("UNRATE",src="FRED")  ## Unemployment rate
> getSymbols("GDPDEF",src="FRED")  ## Quarterly GDP Implicit Price Deflator
> u <- aggregate(UNRATE,as.yearqtr,mean)
> gdpdef <- aggregate(GDPDEF,as.yearqtr,mean)
> pi <- diff(log(gdpdef))*400
> pilag <- lag(pi,-1)
> tvnairu <- cbind(pi,pilag,u)
> tvnairu.df <- subset(data.frame(tvnairu), !is.na(pi) & !is.na(u) &
> !is.na(pilag))
> 
> 
> ## First attempt
> buildNAIRU <- function(x) {
>   modNAIRU <- dlm(FF=t(matrix(c(1,1,1,0))),
>                   GG=diag(4),
>                   W=matrix(c(0,0,0,0,  0,0,0,0,  0,0,0.04,0,  0,0,0,0),4,4),
>                   V=exp(x[4]), m0=rep(0,4), C0=diag(1e07,4),
>                   JFF = t(matrix(c(1,1,0,0))),
>                   X=cbind( tvnairu.df$pilag, tvnairu.df$u))
>   return(modNAIRU)
> }
> 
> (fitNAIRU <- dlmMLE(tvnairu.df$pi, parm=c(0,0,0,0) , build=buildNAIRU,
> hessian=TRUE, control=list(maxit=500)))
> (dlmNAIRU <- buildNAIRU(fitNAIRU$par))
> 
> 
> ## Second attempt
> buildNAIRU <- function(x) {
>   modNAIRU <- dlm(FF=t(matrix(c(1,1,0,0))),
>                   GG=diag(4),
>                   W=matrix(c(0,0,0,0, 0,0,0,0, 0,0,0.04,0, 0,0,0,0 ),4,4),
>                   V=exp(x[4]), m0=c(0.7,-0.3,5,1), C0=diag(100,4),
>                   JFF = t(matrix(c(1,1,0,0))),
>                   X=cbind( tvnairu.df$pilag, tvnairu.df$u - x[3] ))
>   return(modNAIRU)
> }
> 
> (fitNAIRU <- dlmMLE(tvnairu.df$pi, parm=c(0,0,0,0) , build=buildNAIRU,
> hessian=TRUE, control=list(maxit=500)))
> (dlmNAIRU <- buildNAIRU(fitNAIRU$par))
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

-- 
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%
 Dr. Gavin Simpson             [t] +44 (0)20 7679 0522
 ECRC, UCL Geography,          [f] +44 (0)20 7679 0565
 Pearson Building,             [e] gavin.simpsonATNOSPAMucl.ac.uk
 Gower Street, London          [w] http://www.ucl.ac.uk/~ucfagls/
 UK. WC1E 6BT.                 [w] http://www.freshwaters.org.uk
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%



More information about the R-help mailing list