Cross Validated is a question and answer site for people interested in statistics, machine learning, data analysis, data mining, and data visualization. Join them; it only takes a minute:

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

I'm trying to understand how fitted values are calculated for ARMA(p,q) models. I've already found a question on here concerning fitted values of ARMA processes but haven't been able to make sense of it.

If I have a ARMA(1,1) model, i. e.

$$X_t = \alpha_1X_{t-1}+\epsilon_t - \beta_1 \epsilon_{t-1}$$

and am given a (stationary) time series I can estimate the parameters. How would I calculate the fitted values using those estimates. For an AR(1) model the fitted values are given by

$$\hat{X_t} = \hat{\alpha_1}X_{t-1} . $$

Since the innovations in an ARMA model are unobservable how would I use the estimate of the MA parameter? Would I just ignore the MA part and calculate the fitted values of the AR part?


share|improve this question
up vote 8 down vote accepted

To answer your questions, you basically need to know how the residuals i.e. $e_t$ are calculated in an arma model. Because then $\hat{X_{t}}=X_{t}-e_{t}$. Let's first generate a fake data ($X_t$) from arima(.5,.6) and fit the arma model (without mean):

> library(stats)
> library(forecast)
   > n=1000
    > ts_AR <- arima.sim(n = n, list(ar = 0.5,ma=0.6))
    > f=arima(ts_AR,order=c(1,0,1),include.mean=FALSE)
    > summary(f)
    Series: ts_AR 
    ARIMA(1,0,1) with zero mean     

             ar1     ma1
          0.4879  0.5595
    s.e.  0.0335  0.0317

    sigma^2 estimated as 1.014:  log likelihood=-1426.7
    AIC=2859.4   AICc=2859.42   BIC=2874.12

    Training set error measures:
                         ME    RMSE       MAE      MPE     MAPE      MASE
    Training set 0.02102758 1.00722 0.8057205 40.05802 160.1078 0.6313145

Now I create the residuals as follows: $e_1=0$ (since there is no residual at 1) and for $t=2,...,n$ we have: $e_t=X_t-Ar*X_{t-1}-Ma*e_{t-1}$, where $Ar$ and $Ma$ are the estimated auto-regressive and moving average part in above fitted model. Here is the code:

>  e = rep(1,n)
>       e[1] = 0 ##since there is no residual at 1, e1 = 0
>       for (t in (2 : n)){
+         e[t] = ts_AR[t]-coef(f)[1]*ts_AR[t-1]-coef(f)[2]*e[t-1]
+       }

Once you find the residuals $e_{t}$, the fitted values are just $\hat{X_{t}}=X_{t}-e_{t}$. So in the following, I compared the first 10 fitted values obtained from R and the ones I can calculate from $e_{t}$ I created above (i.e. manually).

> cbind(fitted.from.package=fitted(f)[1:10],fitted.calculated.manually=ts_AR[1:10]-e[1:10])
      fitted.from.package fitted.calculated.manually
 [1,]          -0.4193068                 -1.1653515
 [2,]          -0.8395447                 -0.5685977
 [3,]          -0.4386956                 -0.6051324
 [4,]           0.3594109                  0.4403898
 [5,]           2.9358336                  2.9013738
 [6,]           1.3489537                  1.3682191
 [7,]           0.5329436                  0.5219576
 [8,]           1.0221220                  1.0283511
 [9,]           0.6083310                  0.6048668
[10,]          -0.5371484                 -0.5352324

As you see there are close but not exactly the same. The reason is that when I created the residuals I set $e_{1}=0$. There are other choices though. For example based on the help file to arima, the residuals and their variance found by a Kalman filter and therefore their calculation of $e_t$ will be slightly different from me. But as time goes on they are converging.
Now for the Ar(1) model. I fitted the model (without mean) and directly show you how to calculate the fitted values using the coefficients. This time I didn't calculate the residuals. Note that I reported the first 10 fitted values removing the first one (as again it would be different depending on how you define it). As you can see, they are completely the same.

> f=arima(ts_AR,order=c(1,0,0),include.mean=FALSE)
> cbind(fitted.from.package=fitted(f)[2:10],fitted.calculated.manually=coef(f)*ts_AR[1:9])
      fitted.from.package fitted.calculated.manually
 [1,]          -0.8356307                 -0.8356307
 [2,]          -0.6320580                 -0.6320580
 [3,]           0.0696877                  0.0696877
 [4,]           2.1549019                  2.1549019
 [5,]           2.0480074                  2.0480074
 [6,]           0.8814094                  0.8814094
 [7,]           0.9039184                  0.9039184
 [8,]           0.8079823                  0.8079823
 [9,]          -0.1347165                 -0.1347165
share|improve this answer
In the help file to arima they say: "(...) the innovations and their variance found by a Kalman filter." So the function apparently somehow uses the Kalman filter for the initial values. – DatamineR May 16 '14 at 20:49
@RStudent, thanks I updated the answer. – Stat May 17 '14 at 14:59

Your Answer


By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.