# different coefficients in ARX(1) and ln. regression with lagged y

According to the theory ARX(1) is the same as a linear regression with a lagged y. For some reason the coefficients can be very different though. Here is a simple dataset for demonstration purposes:

y=c(3,4,6,3,3,5,7,8,5,7,3,2,5,7,4,3,5,6,8,4,3,5,6,3,2,5,6,7,8) y.lag1=c(3,3,4,6,3,3,5,7,8,5,7,3,2,5,7,4,3,5,6,8,4,3,5,6,3,2,5,6,7) x1=c(7,4,9,6,4,5,3,5,7,9,7,8,5,6,8,4,3,3,6,3,1,4,5,3,4,5,6,7,4) x2=c(6,3,5,6,7,4,6,8,9,6,0,5,7,4,6,7,8,5,4,5,6,4,1,2,4,7,8,0,8)

Fitting an ARX(1) model gives me the following coefficients:

> arx.model=arima(y, xreg=data.frame(x1, x2), order=c(1,0,0), method="ML")

Coefficients: ar1 intercept x1 x2 0.331 3.6976 0.1967 0.0430

Fitting a linear Regression with a lagged y variable gives me the following coefficients:

> lm.model=lm(y~y.lag1+x1+x2)

Coefficients: (Intercept) y.lag1 x1 x2 2.44489 0.31319 0.09048 0.10077

I'm wondering why there is such a big difference in the coefficients. I assume that the different estimation procedures (OLS in the linear regression and ML in ARX) do not cause such a big difference. Does anyone know what causes the difference?

Edit:

This is how I understand the concept of a lagged y variable:

Instead of adding the number 3 to the first missing observation in y_(t-1) I now ignore the first row (as suggested in the comments). That means, I use the sample size T-1. Using this logic I refitted the model:

disign.matrix=data.frame(y.lag1=y[1:28], x1=x1[2:29], x2=x2[2:29]) lm.model=lm(y[2:29]~y.lag1+x1+x2, data=disign.matrix)

and get the following coefficients:

(Intercept) y.lag1 x1 x2

2.4645 0.2722 0.1265 0.1106

They are still quite different from the ARX coefficients. Any further suggestions?

• At this point I think the question is really becoming off-topic as it's clearly about R not the general statistics. MLE and OLS should give the same parameter estimates in ARX model, the only reason they're not is something to do with R code and its packages Commented Jun 1, 2017 at 19:00

The ARX(1) model that you are considering, that is, $$y_t = \theta_0 + \phi_1 y_t + \beta_1 x_t + w_t \tag{1}$$ is not the same as a regression model $$y_t = \beta_0 + \beta_1 x_t + e_t \tag{2a}$$ with errors $e_t$ following an autoregressive model $$e_t = \phi_1 e_t + w_t. \tag{2b}$$ An important difference is that for model (2), $y_t$ is only influenced only by the current value of $x_t$ whereas for model (1) $y_t$ is influenced by the current and all past values of $x_t$.

The arima function in R with covariates contained in the xreg argument fits models such as (2).

To fit ARIMAX models such as the above ARX(1) model use the arimax function in the TSA package.

In addition, method="ML" (maximum likelihood) is not the same as OLS equivalent to the conditional sum of square method for models with no MA part (method="CSS") as the latter method ignores information contained in the first initial observations. Methods based on the exact likelihood computes the probability of the first observations based on the stationary distribution of the process, via the Kalman filter applied to the state space representation of the model. So this is another reason for the difference you see.

If you line up the series right you certainly should get the same results in MLE and OLS estimations in terms of the coefficients of ARX(1) model.

I think in your case there's something funky that R arima function does with pre-sampling. It's the very first observation of the dependent variable. In MATLAB you either supply that observation with Y0 parameter or MATLAB will back forecast it. I'll show you the code where I'm providing this value below. If I let MATLAB backforecast the first observation then OLS and ARIMA estimates will not match.

Here's how it works in MATLAB:

y=[3,4,6,3,3,5,7,8,5,7,3,2,5,7,4,3,5,6,8,4,3,5,6,3,2,5,6,7,8]';
x1=[7,4,9,6,4,5,3,5,7,9,7,8,5,6,8,4,3,3,6,3,1,4,5,3,4,5,6,7,4]';
x2=[6,3,5,6,7,4,6,8,9,6,0,5,7,4,6,7,8,5,4,5,6,4,1,2,4,7,8,0,8]';
X=[x1 x2];

model = arima(1,0,0);
fit = estimate(model,y(2:end),'X',X(2:end,:),'Y0',y(1))


With output

                              Standard          t
Parameter       Value          Error       Statistic
-----------   -----------   ------------   -----------
Constant        2.46451       1.69782        1.45157
AR{1}       0.272167      0.233228        1.16695
Beta1       0.126459      0.189871       0.666028
Beta2       0.110629      0.136319       0.811544
Variance        2.91062       1.76726        1.64696


Then we do OLS:

fitlm([y(1:end-1) X(2:end,:)],y(2:end))


and output:

Estimated Coefficients:
Estimate      SE        tStat     pValue
________    _______    _______    _______

(Intercept)     2.4645      1.5453     1.5949    0.12383
x1             0.27217     0.20217     1.3462     0.1908
x2             0.12646     0.18162    0.69627    0.49295
x3             0.11063      0.1483      0.746    0.46291

• In addition to this reply, you must remove the first row of the data! see dspace.brunel.ac.uk/handle/2438/14135 Commented May 31, 2017 at 13:48
• hang on! lag looks okay, except the first one that is not important! Commented May 31, 2017 at 14:05
• @Aksakal and TPArrow: Thanks for the quick response. Please check my edit above.
– shb
Commented Jun 1, 2017 at 17:28