Solved – Fit an ARMAX model in R

arimaarmaxr

I would like to fit an ARMAX model in R of the form that is mostly used in literature:

$$y_t = \beta_1 x_t+\cdots+ \beta_{k} x_{t_{k-1}}+ \phi_1 y_{t-1}+\cdots+\phi_p y_{t-p}+ \theta_1 z_{t-1}+\cdots \theta_{q} z_{t-q} + z_t$$
where $x$ is exogenous input and $z$ is white noise.
Using the backshift operator $B$ with $B^k(y_t)=y_{t-k}$ the equation becomes

\begin{align}\tag{1}\label{a}
\phi(B)y_t = \beta(B)x_t + \theta(B)z_t
\end{align}

I know that R's built-in arima function includes exogenous input as follows:

$$ y_t = \beta x_t + n_t$$
$$ \varphi(B) n_t = \theta(B) z_t $$

But this is not want I want.
I did some research and found out that there are (at least) three possible functions that fit ARMA models with exogenous variables:

1) stats:::arima (built-in)

2) forecast:::Arima

3) TSA:::arima/arimax (2 identical functions with different names)

As @Rob Hyndman writes in 2010 (https://robjhyndman.com/hyndsight/arimax/),
functions 1) and 2) build a model of the form I don't want, whereas function 3) builds a transfer function model of the form

$$ y_t = \frac{\beta(B)}{\nu(B)}x_t + \frac{\theta(B)}{\phi(B)}z_t$$

An ARMAX model is a special case of the above transfer function model with $\nu=\phi$. However, (as far as I know) the TSA:::arimax function does not provide a possibility to respect this constraint. Therefore, it is not possible to fit a "true" ARMAX model.

My questions are:

1) Is what I write correct?

2) Since most of what I have read about this is 6-8 years old, is there a way to fit an ARMAX model in the form of (\ref{a}) in R?

Best Answer

As far as my research on this topic has taken me, I agree that that the arima/Arima functions from the stats and forecast packages, does fit transfer functions as you mention, but instead a linear model with ARMA errors.

I don't see the possibility to tell the TSA::arimax function that the $\nu(B)$ should be equal to $\phi(B)$. But you can force the order of all the individual polynomial in $B$ to be what you want, including $0$. But that does not really give you your ARMAX model.

The last thing I can suggest is to take a look at the MARIMA package. It should somehow be able to fit an ARMAX model, but I am not 100 % sure about the procedure.

Related Question