Forecasting error in R when passing around argumen

Question

When trying to compose a function from smaller ones using Rob Hyndman's forecast library, like so:

> library('forecast')
> arf <- function(data, ...) forecast(ar(data, order.max=1, method="ols"), ...)

I get an error when trying to plug in some data:

> arf(ts(1:100, start=c(2000,1), frequency=4))
Error in ts(x, frequency = 1, start = 1) : object is not a matrix

However, using the body of arf directly works perfectly:

> forecast(ar(ts(1:100, start=c(2000,1), frequency=4), order.max=1,method="ols"))
        Point Forecast Lo 80 Hi 80 Lo 95 Hi 95
2025 Q1            101   101   101   101   101
2025 Q2            102   102   102   102   102
2025 Q3            103   103   103   103   103
2025 Q4            104   104   104   104   104
2026 Q1            105   105   105   105   105
2026 Q2            106   106   106   106   106
2026 Q3            107   107   107   107   107
2026 Q4            108   108   108   108   108
2027 Q1            109   109   109   109   109
2027 Q2            110   110   110   110   110

Why does arf not work as it should?

Answer 1

the problem is a bug in the S3 method predict for ar class. predict.ar try to evaluate the newdata parameter using ar object. Showing the first lines of getS3method('predict','ar')

function (object, newdata, n.ahead = 1L, se.fit = TRUE, ...) 
{
    if (n.ahead < 1L) 
        stop("'n.ahead' must be at least 1")
    if (missing(newdata)) {
        newdata <- eval.parent(parse(text = object$series))
        if (!is.null(nas <- object$call$na.action)) 
            newdata <- eval.parent(call(nas, newdata))
    }

  .....
}

the relevant/bugged line is:

 newdata <- eval.parent(parse(text = object$series))

But object$series don't have the right expression/character. since it is hidden by the new level of wrapper function arf. Here a workaround, is to set the right expression for this term:

arf <- function(dat, ...) 
  {
    object <- ar(dat, order.max=1, method="ols")
    object$series <- as.character(as.expression(as.list(match.call())$dat))
    forecast(object,...)
}
arf( ts(1:100, start=c(2000,1), frequency=4)

Note that ; this solution works also with:

aa <- ts(1:100, start=c(2000,1), frequency=4)
arf(aa)

Answer 2

This is a problem (not really a bug) in forecast.ar(). All forecast.xxx() functions try to store the data used to estimate the time series model as it is required for plots and accuracy calculations. However, ar() does not return the data, so forecast.ar() attempts to find the data in the calling environment, or in the parent environment. When you call forecast(ar(...)) directly, the function manages to find the data, but arf() places the call to ar() one level deeper making it that much harder for forecast to figure out what data was being used.

I can modify the function to make it look harder for the data (i.e., look in the grandparent environment also), but the construction will still fail because predict.ar() (part of the stats package) will then cause a similar error. It would be much better if ar() returned the data, but ar() is part of the stats package and I have no control over it.

There are several possible solutions.

1. You could replace ar with Arima:

   arf <- function(dat, ...) forecast(Arima(dat, order=c(1,0,0)), ...)
   

   That should return the same model if the data are stationary (although the parameter estimates will be slightly different). It won't return the same answer for the example in your question because the time series is non-stationary.

2. You could use auto.arima() instead if you were willing to use more general ARIMA models than AR(1).

   arf <- function(dat, ...) forecast(auto.arima(dat, ...)
   

3. (Based on a suggestion from @agstudy). A workaround is to ensure the data is stored within the ar object:

   arf <- function(dat, ...) 
   {
     object <- ar(dat, order.max=1, method="ols")
     object$x <- dat
     forecast(object,...)
   }