[R] Problem with arima function with xreg argument to forecast using predict - How to write the forecast function?

Richard A. Bilonick rab at nauticom.net
Thu Apr 17 21:44:27 CEST 2003


Using the LakeHuron data, fit a simple AR(2) model:

 > data(LakeHuron) 
 > ar.lh <- arima(LakeHuron, order = c(2,0,0))
 > ar.lh

Call:
arima(x = LakeHuron, order = c(2, 0, 0))

Coefficients:
         ar1      ar2  intercept
      1.0436  -0.2495   579.0473
s.e.  0.0983   0.1008     0.3319

sigma^2 estimated as 0.4788:  log likelihood = -103.63,  aic = 215.27

Make a 1-step ahead forecast:

 > predict(ar.lh,1)[[1]]
Time Series:
Start = 1973
End = 1973
Frequency = 1
[1] 579.7896

Compute the forecast manually:

 > sum(ar.lh$coef*c(c(579.96,579.89)-ar.lh$coef[3],1))
[1] 579.7896

This just says that the forecast for the next period (after the end of 
the data) is 579.0473 + 1.0436*(579.96 - 579.0473) - 0.2495*(579.89 - 
579.0473). In other words: the forecast is the intercept plus the AR 
coefficients times the  (previous ts values minus the intercepts).

Now add an exogenous variable (in this case, the (year - 1920):

 > ar.lh <- arima(LakeHuron, order = c(2,0,0), xreg = time(LakeHuron)-1920)
 > ar.lh

Call:
arima(x = LakeHuron, order = c(2, 0, 0), xreg = time(LakeHuron) - 1920)

Coefficients:
         ar1      ar2  intercept  time(LakeHuron) - 1920
      1.0048  -0.2913   579.0993                 -0.0216
s.e.  0.0976   0.1004     0.2370                  0.0081

sigma^2 estimated as 0.4566:  log likelihood = -101.2,  aic = 212.4

The prediction is:

 > predict(ar.lh,1,newxreg=53)[[1]]
Time Series:
Start = 1973
End = 1973
Frequency = 1
[1] 579.3972


Now try to manually forecast where the next time period is 53:

 > sum(ar.lh$coef*c(c(579.96,579.89)-ar.lh$coef[3],1,53))
[1] 578.5907

What am I doing wrong? I've tried this with numerous examples and 
whenever there is an exogenous variable I cannot get the manual forecast 
to agree with predict. Is it not correct to just add (-0.0216 times 53) 
to the rest? I need to know how to write the model correctly. Obviously 
there is something I am overlooking. R's arima function and predict 
function work correctly - at least they agree with SAS for example so 
I'm not doing something right.

I would really appreciate some insight here.

Rick B.



More information about the R-help mailing list