Interpretation of TSA::arimax output model is presented in R I have a relatively simple problem, but yet taking some time to solve it. I am using the arimax() function from the TSA package. (Note: not arima() from the stats package.) This is the model: 
out <- arimax(sub_s_t_series, order=c(2,0,1), xreg=sub_r_t_series, method=c("ML"))

and these are my coefficients:
Call
arimax(x = sub_s_t_series, order = c(2, 0, 1), xreg = sub_r_t_series, method = c("ML"))

Coefficients:
         ar1      ar2      ma1  intercept     xreg
      1.4825  -0.6613  -0.8516  52745.107  -1.0132
s.e.  0.0295   0.0294   0.0064     40.828   0.0012

sigma^2 estimated as 0.08929:  log likelihood = -105.98,  aic = 221.97

All I am trying to do is to interpret the results. According to my understanding and the help given in the TSA package, the above ARIMAX(2,0,1) model is represented as follows:
$$
{\rm sub\_s\_t\_series\_hat[k]} = {\rm intercept} + xreg\times {\rm sub\_r\_t\_series[k]} + 
\frac{a_{t[k]}+ma1*a_{t[k-1]}}{a_{t[k]}-ar1*a_{t[k-1]}-ar2*a_{t[k-2]}}  \tag{1} 
$$
where $a_t$ are the residuals. When I use e_t = fitted(out)-sub_s_t_series_hat to measure the error / residuals myself, e_t matches exactly to the values obtained by out[["residuals"]].
But when I use (1) as follows: e_t_hat = sub_s_t_series_hat - sub_s_t_series, 
e_t_hat does not match with out[["residuals"]], in fact the results deviate by a magnitude of almost 4.
My questions is: did an ARIMAX(2,0,1) fit would result in (1) or am I missing something?
 A: This is an old post, but since there is no accepted answer, I still want to provide an explanation for future readers. You have the following R output:
Coefficients:
         ar1      ar2      ma1  intercept     xreg
      1.4825  -0.6613  -0.8516  52745.107  -1.0132
s.e.  0.0295   0.0294   0.0064     40.828   0.0012

sigma^2 estimated as 0.08929:  log likelihood = -105.98,  aic = 221.97

(1) y = intercept + xreg * x series name + n_t(error term)
(2) n_t(error term) = ar1 * n(t-1) + ar2 * n(t-2) + white noise epsilon + ma1 * epsilon(t-1)
(3) white noise epsilon ~ NID(0,0.08929)
To explain it further:
(1) This formula focuses on xreg coefficient. The R outputs gives you a coefficient for xreg -- the effect of the exogenous variable on Y at it's current value.
(2) This formula focuses on unpacking the error term left unexplained by formula (1). Formula 1 and 2 together is also called regression with ARIMA error (ARMA error in this case). This error term itself is an ARIMA (2,0,1) process, according to your output. This is what ar1, ar2, ma2 in the R output corresponds to.
(3) This formula is describing what is left after formula 1 and 2. sigma^2 describes the variance of the white noise series -- epsilon in most Econometric textbooks. 
