How to create a random walk model using {forecast} R package

Question

I have a good understanding of ARIMA models but I've always found significant spikes in ACFs and PACFs that gave me the appropriate AR and MA parameters.

Now I'm dealing with a series that is more like an ARIMA(0,1,0) which I think is the same as a random walk? This random walk concept is a little new to me but I sort of understand it. My question now is how do I even create a model for this series?

I've seen some methods in the {forecast} package that might be the ones I'm looking for but I want to understand how they are different.

The functions that I am confused about are naive() and rwf(). It seems like both try to address the same random walk problem. But also, what happens if I just build an Arima() model with parameters 0,1,0? how is that different?

Answer 1

The forecasts from a random walk are flat and equal to the last observation. Adding a drift term, a trend pattern can be captured. This answer shows that a constant in a random walk has the effect of a deterministic linear trend. Some illustrations and related comments are given in this post and this post.

---

If you choose taking logarithms on the data, the original scale of the data needs to be recovered. I think that naive and rwf undo the logarithms by taking the exponential function on the forecasts. While this is the common approach, it is more accurate to rescale the forecasts taking into account that $E(\exp(\epsilon_t)) = \exp(\sigma_\epsilon^2/2)$: $$ \exp\{E_T(y_{t+k}) + 0.5 \times \hat{\sigma}^2_{e_T(k)}\} \,, $$ i.e., the exponential of the forecasts times half the variance of the forecast errors. This rescaling transformation is suggested in Novales (1993, Sec. 13.9) Econometría, Mc Graw-Hill. This comes from Jensen's inequality, discussed for example here.

Here is a small simulation exercise for a random walk fitted in logs:

iter <- 10000 # number of iterations
# storage matrices for accuracy measures
m1 <- m2 <- matrix(nrow = iter, ncol = 3)
colnames(m1) <- colnames(m2) <- c("ME", "RMSE", "MAE")
set.seed(123)
for (i in seq_len(iter))
{
  # generate 210 observations from a random walk
  y_210 <- cumsum(rnorm(210, mean = 10))
  y_200 <- y_210[1:200]
  # fit a random walk,
  # using arima for convenience to get forecast standard errors
  fit <- arima(log(y_200), order = c(0,1,0))
  pred <- predict(fit, n.ahead=10)
  fcast <- exp(pred$pred)
  # get accuracy measures for the standard forecasts
  m1[i,] <- c(
mean(y_210[201:210] - fcast), # Mean Error
sqrt(mean((y_210[201:210] - fcast)^2)), # Root Mean Squared Error
mean(abs(y_210[201:210] - fcast))) # Mean Absolute Error
  # rescale (undo logs) adjusting for the variance
  fcast.adjusted <- exp(pred$pred + 0.5 * pred$se^2)
  # get accuracy measures for the adjusted forecasts
  m2[i,] <- c(
    mean(y_210[201:210] - fcast.adjusted), # Mean Error
    sqrt(mean((y_210[201:210] - fcast.adjusted)^2)), # Root Mean Squared Error
    mean(abs(y_210[201:210] - fcast.adjusted))) # Mean Absolute Error
}
cbind("fcast" = colMeans(m1), "fcast.adjusted" = colMeans(m2))
#         fcast fcast.adjusted
# ME   55.03849       27.52379
# RMSE 62.09735       31.03224
# MAE  55.03849       27.52379

In this exercise we see that the forecasts based on the transformation given above have smaller accuracy measures than those rescaled taking the exponential. Further insight would be required in order to assess to what extent the improvement is actually relevant.

For a random walk with drift the same should be observed, but I didn't get expected results in the simulation; maybe the forecast errors that I used did not account for the presence of the drift term.