# ar() time series function in R, manually checking the residuals/predicted values

I am using the ar() function to fit an AR model to some data, and this object will return the in sample residuals. I also know the syntax for how to get the corresponding predicted values, but I want to compute these predicted values manually (just to check my own understanding) for a simple AR(1) example.

My problem is that my manually computed residuals (based on my manually computed predictions) do not match the in sample residuals stored in the ar object (well the 1st residual does match, but not the rest).

From the documentation, I see that

x[t]  = m + a[1]*(x[t-1] - m)


where m is the sample mean of the series. Here is an example of what I am doing manually.

# Create some true AR(1) data
set.seed(123)
x = w = rnorm(30) + 3 ; for (t in 2:30) x[t] = .60*x[t-1] + w[t]
# Fit an ar model
x.model = ar(x) # coefficient is .49, mean value of x is 6.98
# Manually create predictions
x.MyPred = rep(0,30) ; x.MyPred[1] = x[1]
for (t in 2:30) x.MyPred[t] = 6.984234 + .4988327*(x.MyPred[t-1] - 6.984234)
MyResid = x - x.MyPred
cbind(MyResid, x.model\$res) # does not match


And interestingly, the first residual (at observation 2) does match, but the rest do not. Thanks in advance.

Update

The answer given below is basically highlighting the difference between so called static forecasting and dynamic forecasting, here are a few more details.

The two possible choices to make the fitted values are

# Method 1:
for (t in 2:30) x.MyPred[t] = 6.984234 + .4988327*(x.MyPred[t-1] - 6.984234)

# Method 2:
for (t in 2:30) x.MyPred[t] = 6.984234 + .4988327*(x[t-1] - 6.984234)


Method 1 is taking the forecasted value from the previous step and using this to plug in for the lagged values on the RHS. This could feasibly be used to create out of sample predictions forever. Also, depending on what time point you start to make predictions, the resulting predictions for a given time period can be different.

Method 2 is taking the actual known historical values to plug in for the lagged values on the RHS. It will never be able to forecast more than 1 step out of sample. Also, the forecast values will always be the same no matter where you start.

Both of these will give the same 1 step ahead forecast value. Method 2 is what produces the residuals returned in the ar() object.

The Eviews software documentation has a good discussion of this. This R documentation was problematic in not specifying which of the x values are fitted values versus known historical values. A better expression to show in the documentation would be

\hat{x}[t]  = m + a[1]*(x[t-1] - m)


and, dare I say, possibly a few sentences on this very topic. But good documentation is hardly something to expect in R.

The reason is that you are using x.MyPred[t-1] when you should use x[t-1]. What you are doing is creating a series of predictions over the entire sample period, with only the initial value given. The residuals are one-step ahead predictions, meaning that the lagged value is taken as given.