I was experimenting with the most impressive forecast package. I tried to
rig up a simplest demo of auto.arima() in action. The forecast sd looks odd
to me. Here's a demo.
## Reproducibility
set.seed(101)
## Simulate a series from AR(1)
x <- arima.sim(list(order = c(1,0,0), ar = 0.8), n = 1000)
sd(x)
# This is 1.566 which seems right
## Using ar()
m1 <- ar(x)
## Use forecast::auto.arima()
library(forecast)
m2 <- auto.arima(x)
## How did they fare?
m1
m2
# Both of them got the right model.
f1 <- predict(m1, n.ahead=10)
data.frame(pointest=as.numeric(f1$pred), sd=as.numeric(f1$se))
f2 <- forecast(m2, h=10)
data.frame(pointest=cbind(as.numeric(f2$mean),
sd=(f2$upper[,2]-f2$lower[,2])/2))
Let me show you the two results.
pointest sd
1 -0.3186451 0.9628466
2 -0.2898804 1.2264023
3 -0.2671877 1.3649788
4 -0.2492853 1.4445294
5 -0.2351619 1.4918999
6 -0.2240198 1.5206374
7 -0.2152297 1.5382519
8 -0.2082952 1.5491136
9 -0.2028244 1.5558355
10 -0.1985085 1.5600043
This seems right. The forecast sd is small at first, as we're exploiting
the time series structure. But as we go deep into the future, the AR model
is useless and we're down to the unconditional sd which is ~ 1.566.
sd=(f2$upper[,2]-f2$lower[,2])/2))
pointest.V1 pointest.sd
1 -0.28090691 1.887602
2 -0.22221139 2.406792
3 -0.17578030 2.681016
4 -0.13905100 2.839174
5 -0.10999628 2.933809
6 -0.08701255 2.991505
7 -0.06883127 3.027050
8 -0.05444897 3.049081
9 -0.04307186 3.062787
10 -0.03407199 3.071333
This seems out of line. Right at n.ahead=1 the forecast sd is 1.8876 which
is > 1.566. And going beyond, the forecast sd goes up to 3.07. What am I
missing?