[R] how to predict dynamic model in R

[Gabor Grothendieck]

Best thing is to test it out on simulated data for which you
already know the answer.

library(dyn)
set.seed(123)
DF <- data.frame(x = rnorm(10), z = rnorm(10))
DF$Y <- 0
for(i in 2:10) {
  DF$Y[i] <- with(DF, 2*Y[i-1] + 3*z[i] +4* x[i] + 5*x[i-1] + rnorm(1))
}
DF.zoo <- do.call(merge, lapply(DF, zoo))
L <- function(x, k = 1) lag(x, -k)
fit <- dyn$lm(Y ~ L(Y) + z + L(x, 0:1), DF.zoo[-10, ])
fit
pred <- predict(fit, DF.zoo)
cbind(pred, DF.zoo$Y)
plot(cbind(pred, act = DF.zoo$Y), type = c("l", "p"), screen = 1)




On Thu, Jul 23, 2009 at 2:48 PM, Hongwei Dong<pdxdong at gmail.com> wrote:
> Hi, Gabor, got it. Thanks a lot.
> I have one more question about how the "predict" function works here,
> especially for the lag(Y,-1) part.
> In my model, I assume I know predictors x and z in the next two years, and
> use them to predict Y. For each forecast step at time t, the lag(Y,-1) in
> the model should be the estimated result from the last forecast step. For
> example, the prediction of the Y in 2007 should be based on the estimated Y
> for 2006, not the actual number already in my dataset. Before I report my
> results, I want to make sure the "predict" function works in this way.
> Thanks.
> Harry
>
>
> On Thu, Jul 23, 2009 at 10:28 AM, Gabor Grothendieck
> <ggrothendieck at gmail.com> wrote:
>>
>> Please provide your code in a reproducible form (as requested previously).
>>
>> Here we fit with first 9 points and then add a point for prediction.  (Of
>> course your model can only predict the current value of Y so you may
>> have to rethink your model even aside from the implementation if you
>> really want to predict future values of Y.)
>>
>> > library(dyn)
>> > set.seed(123)
>> > tt <- ts(cbind(Y = 1:10, x = rnorm(10), z = rnorm(10)))
>> > L <- function(x, k = 1) lag(x, -k)
>> > tt.zoo <- as.zoo(tt)
>> > fit <- dyn$lm(Y ~ L(Y) + L(x, 0:2) + z, tt.zoo[-10, ])
>> > predict(fit, tt.zoo)
>>  1  2  3  4  5  6  7  8  9 10 11 12
>> NA NA  3  4  5  6  7  8  9 10 NA NA
>>
>>
>>
>>
>> On Thu, Jul 23, 2009 at 1:01 PM, Hongwei Dong<pdxdong at gmail.com> wrote:
>> > Hi, Gabor and Other R users,
>> >   I'm re-posting my script and the results I got.
>> >
>> >   here is the dynamic model I used to estimate in-sample model
>> > (1996-2006)
>> > and it works:
>> >
>> >   fit<-dyn$lm(Y~lag(Y,-1)+z+x+lag(x,-1)+lag(x,-2)+lag(x,-3)+lag(x,-4))
>> >   Then I used this model to do out sample forecast with the following
>> > scripts, which do not work:
>> >
>> >  z<-ts(Z[41:52],start=2006,frequency=4)
>> >  x<-ts(X[41:52],start=2006,frequency=4)
>> > newdata<-data.frame(cbind(z,x))
>> > newdata<-ts(newdata,start=2006,frequency=4)
>> > pred<-predict(fit,newdata)
>> > Here is the results I got from R:
>> >              Qtr1     Qtr2     Qtr3     Qtr4
>> > 2006       NA       NA       NA       NA
>> > 2007 3083.362       NA       NA       NA
>> > 2008       NA       NA       NA       NA
>> > 2009       NA       NA       NA       NA
>> >
>> > I got only one prediction for the first quarter in 2007. Intuitively,
>> > there
>> > might be two problems: the definition of the newdata and how to define Y
>> > in
>> > newdata. But I just can't figure this out. It will greatly appreciated
>> > if
>> > someone can give me some help. Thanks.
>> > Harry
>> >
>> >
>> > On Thu, Jul 23, 2009 at 5:15 AM, Gabor Grothendieck
>> > <ggrothendieck at gmail.com> wrote:
>> >>
>> >> You have to use consistent classes.  You can't start out using
>> >> ts class and then suddenly switch to zoo class as if you had been
>> >> using zoo all along.   Either use ts everywhere or zoo everywhere.
>> >>
>> >> Also in the future please post reproducible examples as requested
>> >> at the bottom of every message to r-help.  That means include
>> >> a minimal amount of data so we can get exactly what you
>> >> are getting.
>> >>
>> >> On Thu, Jul 23, 2009 at 4:48 AM, Hongwei Dong<pdxdong at gmail.com> wrote:
>> >> > Thanks, Gabor. This is really helpful.
>> >> > When the regressive part, lag(Y,-1), is not included, my sytax works
>> >> > well.
>> >> > However, when I include the lag(Y) in the model, the prediction part
>> >> > does
>> >> > not work. Here is my sytax for in-sample estimation and it works
>> >> > well:
>> >> > fit<-dyn$lm(Y~lag(Y,-1)+x+lag(x,-1)+lag(x,-2)+lag(x,-3)+lag(x,-4)+z).
>> >> > Then I use this model to do out of sample prediction:
>> >> > x<-ts(X[41:52],start=2006,frequency=4)
>> >> > z<-ts(Z[41:52],start=2006,frequency=4)
>> >> > newdata<-data.frame(cbind(x,z))
>> >> > newdata<-zooreg(newdata)
>> >> > pred<-predict(fit,newdata)
>> >> > With these, I got weird result, a prediction for each year from year
>> >> > 1
>> >> > to
>> >> > the first quarter of year 2007 (all "NA"). What I expect is a
>> >> > prediction
>> >> > for
>> >> > the 8 quarters from 2007 to 2008. Intuitively, I know there must be
>> >> > something wrong with my newdata definition. But I can't figure it
>> >> > out.
>> >> > I'll
>> >> > appreciate it very much if you can give some suggestions to modify my
>> >> > syntax. Thanks.
>> >> >
>> >> > Harry
>> >> >
>> >> >
>> >> >
>> >> >
>> >> >
>> >> >
>> >> >
>> >> >
>> >> >
>> >> >
>> >> >
>> >> > On Wed, Jul 22, 2009 at 10:53 PM, Gabor Grothendieck
>> >> > <ggrothendieck at gmail.com> wrote:
>> >> >>
>> >> >> Here is an example closer to yours.
>> >> >>
>> >> >> > library(dyn)
>> >> >> > set.seed(123)
>> >> >> > x <- zooreg(rnorm(10))
>> >> >> > y <- zooreg(rnorm(10))
>> >> >> > L <- function(x, k = 1) lag(x, k = -k)
>> >> >> > mod <- dyn$lm(y ~ L(y) + L(x, 0:2))
>> >> >> > mod
>> >> >>
>> >> >> Call:
>> >> >> lm(formula = dyn(y ~ L(y) + L(x, 0:2)))
>> >> >>
>> >> >> Coefficients:
>> >> >> (Intercept)         L(y)   L(x, 0:2)1   L(x, 0:2)2   L(x, 0:2)3
>> >> >>    0.06355     -0.74540      0.63649      0.44957     -0.41360
>> >> >>
>> >> >> > newdata <- cbind(x = c(coredata(x), rnorm(1)), y = c(coredata(y),
>> >> >> > rnorm(1)))
>> >> >> > newdata <- zooreg(newdata)
>> >> >> > predict(mod, newdata)
>> >> >>         1          2          3          4          5          6
>> >> >>  7
>> >> >>        NA         NA  0.9157808  0.6056333 -0.5496422  1.5984615
>> >> >> -0.2574875
>> >> >>         8          9         10         11         12         13
>> >> >> -1.6148859  0.3329285 -0.5284646 -0.1799693         NA         NA
>> >> >>
>> >> >>
>> >> >> On Thu, Jul 23, 2009 at 1:04 AM, Gabor
>> >> >> Grothendieck<ggrothendieck at gmail.com> wrote:
>> >> >> > Use dyn.predict like this:
>> >> >> >
>> >> >> >> library(dyn)
>> >> >> >> x <- y <- zoo(1:5)
>> >> >> >> mod <- dyn$lm(y ~ lag(x, -1))
>> >> >> >> predict(mod, list(x = zoo(6:10, 6:10)))
>> >> >> >  7  8  9 10
>> >> >> >  7  8  9 10
>> >> >> >
>> >> >> >
>> >> >> > On Thu, Jul 23, 2009 at 12:54 AM, Hongwei Dong<pdxdong at gmail.com>
>> >> >> > wrote:
>> >> >> >> I have a dynamic time series model like this:
>> >> >> >> dyn$lm( y ~ lag(y,-1) + x + lag(x,-1)+lag(x,-2) )
>> >> >> >>
>> >> >> >> I need to do an out of sample forecast with this model. Is there
>> >> >> >> any
>> >> >> >> way I
>> >> >> >> can do this with R?
>> >> >> >> It would be greatly appreciated if some one can give me an
>> >> >> >> example.
>> >> >> >> Thanks.
>> >> >> >>
>> >> >> >>
>> >> >> >> Harry
>> >> >> >>
>> >> >> >>        [[alternative HTML version deleted]]
>> >> >> >>
>> >> >> >> ______________________________________________
>> >> >> >> 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.
>> >> >> >>
>> >> >> >
>> >> >
>> >> >
>> >
>> >
>
>