Autoregressive model

The notation ${\displaystyle AR(p)}$ indicates an autoregressive model of order p. The AR(p) model is defined as

${\displaystyle X_{t}=\sum _{i=1}^{p}\varphi _{i}X_{t-i}+\varepsilon _{t}}$

where ${\displaystyle \varphi _{1},\ldots ,\varphi _{p}}$ are the parameters of the model, and ${\displaystyle \varepsilon _{t}}$ is white noise.[1][2] This can be equivalently written using the backshift operator B as

${\displaystyle X_{t}=\sum _{i=1}^{p}\varphi _{i}B^{i}X_{t}+\varepsilon _{t}}$

so that, moving the summation term to the left side and using polynomial notation, we have

${\displaystyle \phi [B]X_{t}=\varepsilon _{t}}$

An autoregressive model can thus be viewed as the output of an all-pole infinite impulse response filter whose input is white noise.

Some parameter constraints are necessary for the model to remain weak-sense stationary. For example, processes in the AR(1) model with ${\displaystyle |\varphi _{1}|\geq 1}$ are not stationary. More generally, for an AR(p) model to be weak-sense stationary, the roots of the polynomial ${\displaystyle \Phi (z):=\textstyle 1-\sum _{i=1}^{p}\varphi _{i}z^{i}}$ must lie outside the unit circle, i.e., each (complex) root ${\displaystyle z_{i}}$ must satisfy ${\displaystyle |z_{i}|>1}$ (see pages 89,92 [3]).

In an AR process, a one-time shock affects values of the evolving variable infinitely far into the future. For example, consider the AR(1) model ${\displaystyle X_{t}=\varphi _{1}X_{t-1}+\varepsilon _{t}}$. A non-zero value for ${\displaystyle \varepsilon _{t}}$ at say time t=1 affects ${\displaystyle X_{1}}$ by the amount ${\displaystyle \varepsilon _{1}}$. Then by the AR equation for ${\displaystyle X_{2}}$ in terms of ${\displaystyle X_{1}}$, this affects ${\displaystyle X_{2}}$ by the amount ${\displaystyle \varphi _{1}\varepsilon _{1}}$. Then by the AR equation for ${\displaystyle X_{3}}$ in terms of ${\displaystyle X_{2}}$, this affects ${\displaystyle X_{3}}$ by the amount ${\displaystyle \varphi _{1}^{2}\varepsilon _{1}}$. Continuing this process shows that the effect of ${\displaystyle \varepsilon _{1}}$ never ends, although if the process is stationary then the effect diminishes toward zero in the limit.

Because each shock affects X values infinitely far into the future from when they occur, any given value X_t is affected by shocks occurring infinitely far into the past. This can also be seen by rewriting the autoregression

${\displaystyle \phi (B)X_{t}=\varepsilon _{t}\,}$

(where the constant term has been suppressed by assuming that the variable has been measured as deviations from its mean) as

${\displaystyle X_{t}={\frac {1}{\phi (B)}}\varepsilon _{t}\,.}$

When the polynomial division on the right side is carried out, the polynomial in the backshift operator applied to ${\displaystyle \varepsilon _{t}}$ has an infinite order—that is, an infinite number of lagged values of ${\displaystyle \varepsilon _{t}}$ appear on the right side of the equation.

The autocorrelation function of an AR(p) process can be expressed as

${\displaystyle \rho (\tau )=\sum _{k=1}^{p}a_{k}y_{k}^{-|\tau |},}$

where ${\displaystyle y_{k}}$ are the roots of the polynomial

${\displaystyle \phi (B)=1-\sum _{k=1}^{p}\varphi _{k}B^{k}}$

where B is the backshift operator, where ${\displaystyle \phi (\cdot )}$ is the function defining the autoregression, and where ${\displaystyle \varphi _{k}}$ are the coefficients in the autoregression. The formula is valid only if all the roots have multiplicity 1.

The autocorrelation function of an AR(p) process is a sum of decaying exponentials.

• Each real root contributes a component to the autocorrelation function that decays exponentially.
• Similarly, each pair of complex conjugate roots contributes an exponentially damped oscillation.

AR(0); AR(1) with AR parameter 0.3; AR(1) with AR parameter 0.9; AR(2) with AR parameters 0.3 and 0.3; and AR(2) with AR parameters 0.9 and −0.8

The simplest AR process is AR(0), which has no dependence between the terms. Only the error/innovation/noise term contributes to the output of the process, so in the figure, AR(0) corresponds to white noise.

For an AR(1) process with a positive ${\displaystyle \varphi }$, only the previous term in the process and the noise term contribute to the output. If ${\displaystyle \varphi }$ is close to 0, then the process still looks like white noise, but as ${\displaystyle \varphi }$ approaches 1, the output gets a larger contribution from the previous term relative to the noise. This results in a "smoothing" or integration of the output, similar to a low pass filter.

For an AR(2) process, the previous two terms and the noise term contribute to the output. If both ${\displaystyle \varphi _{1}}$ and ${\displaystyle \varphi _{2}}$ are positive, the output will resemble a low pass filter, with the high frequency part of the noise decreased. If ${\displaystyle \varphi _{1}}$ is positive while ${\displaystyle \varphi _{2}}$ is negative, then the process favors changes in sign between terms of the process. The output oscillates. This can be likened to edge detection or detection of change in direction.

An AR(1) process is given by:

$${\displaystyle X_{t}=\varphi X_{t-1}+\varepsilon _{t}\,}$$

where ${\displaystyle \varepsilon _{t}}$ is a white noise process with zero mean and constant variance ${\displaystyle \sigma _{\varepsilon }^{2}}$. (Note: The subscript on ${\displaystyle \varphi _{1}}$ has been dropped.) The process is weak-sense stationary if ${\displaystyle |\varphi |<1}$ since it is obtained as the output of a stable filter whose input is white noise. (If ${\displaystyle \varphi =1}$ then the variance of ${\displaystyle X_{t}}$ depends on time lag t, so that the variance of the series diverges to infinity as t goes to infinity, and is therefore not weak sense stationary.) Assuming ${\displaystyle |\varphi |<1}$, the mean ${\displaystyle \operatorname {E} (X_{t})}$ is identical for all values of t by the very definition of weak sense stationarity. If the mean is denoted by ${\displaystyle \mu }$, it follows from

$${\displaystyle \operatorname {E} (X_{t})=\varphi \operatorname {E} (X_{t-1})+\operatorname {E} (\varepsilon _{t}),}$$

that

$${\displaystyle \mu =\varphi \mu +0,}$$

and hence

${\displaystyle \mu =0.}$

The variance is

${\displaystyle {\textrm {var}}(X_{t})=\operatorname {E} (X_{t}^{2})-\mu ^{2}={\frac {\sigma _{\varepsilon }^{2}}{1-\varphi ^{2}}},}$

where ${\displaystyle \sigma _{\varepsilon }}$ is the standard deviation of ${\displaystyle \varepsilon _{t}}$. This can be shown by noting that

${\displaystyle {\textrm {var}}(X_{t})=\varphi ^{2}{\textrm {var}}(X_{t-1})+\sigma _{\varepsilon }^{2},}$

and then by noticing that the quantity above is a stable fixed point of this relation.

The autocovariance is given by

${\displaystyle B_{n}=\operatorname {E} (X_{t+n}X_{t})-\mu ^{2}={\frac {\sigma _{\varepsilon }^{2}}{1-\varphi ^{2}}}\,\,\varphi ^{|n|}.}$

It can be seen that the autocovariance function decays with a decay time (also called time constant) of ${\displaystyle \tau =1-\varphi }$.[4]

The spectral density function is the Fourier transform of the autocovariance function. In discrete terms this will be the discrete-time Fourier transform:

${\displaystyle \Phi (\omega )={\frac {1}{\sqrt {2\pi }}}\,\sum _{n=-\infty }^{\infty }B_{n}e^{-i\omega n}={\frac {1}{\sqrt {2\pi }}}\,\left({\frac {\sigma _{\varepsilon }^{2}}{1+\varphi ^{2}-2\varphi \cos(\omega )}}\right).}$

This expression is periodic due to the discrete nature of the ${\displaystyle X_{j}}$, which is manifested as the cosine term in the denominator. If we assume that the sampling time (${\displaystyle \Delta t=1}$) is much smaller than the decay time (${\displaystyle \tau }$), then we can use a continuum approximation to ${\displaystyle B_{n}}$:

${\displaystyle B(t)\approx {\frac {\sigma _{\varepsilon }^{2}}{1-\varphi ^{2}}}\,\,\varphi ^{|t|}}$

which yields a Lorentzian profile for the spectral density:

${\displaystyle \Phi (\omega )={\frac {1}{\sqrt {2\pi }}}\,{\frac {\sigma _{\varepsilon }^{2}}{1-\varphi ^{2}}}\,{\frac {\gamma }{\pi (\gamma ^{2}+\omega ^{2})}}}$

where ${\displaystyle \gamma =1/\tau }$ is the angular frequency associated with the decay time ${\displaystyle \tau }$.

An alternative expression for ${\displaystyle X_{t}}$ can be derived by first substituting ${\displaystyle \varphi X_{t-2}+\varepsilon _{t-1}}$ for ${\displaystyle X_{t-1}}$ in the defining equation. Continuing this process N times yields

${\displaystyle X_{t}=\varphi ^{N}X_{t-N}+\sum _{k=0}^{N-1}\varphi ^{k}\varepsilon _{t-k}.}$

For N approaching infinity, ${\displaystyle \varphi ^{N}}$ will approach zero and:

${\displaystyle X_{t}=\sum _{k=0}^{\infty }\varphi ^{k}\varepsilon _{t-k}.}$

It is seen that ${\displaystyle X_{t}}$ is white noise convolved with the ${\displaystyle \varphi ^{k}}$ kernel plus the constant mean. If the white noise ${\displaystyle \varepsilon _{t}}$ is a Gaussian process then ${\displaystyle X_{t}}$ is also a Gaussian process. In other cases, the central limit theorem indicates that ${\displaystyle X_{t}}$ will be approximately normally distributed when ${\displaystyle \varphi }$ is close to one.

For ${\displaystyle \varepsilon _{t}=0}$, the process ${\displaystyle X_{t}=\varphi X_{t-1}}$ will be a geometric progression (exponential growth or decay). In this case, the solution can be found analytically: ${\displaystyle X_{t}=a\varphi ^{t}}$ whereby ${\displaystyle a}$ is an unknown constant (initial condition).

The partial autocorrelation of an AR(p) process equals zero at lags larger than p, so the appropriate maximum lag p is the one after which the partial autocorrelations are all zero.

There are many ways to estimate the coefficients, such as the ordinary least squares procedure or method of moments (through Yule–Walker equations).

The AR(p) model is given by the equation

${\displaystyle X_{t}=\sum _{i=1}^{p}\varphi _{i}X_{t-i}+\varepsilon _{t}.\,}$

It is based on parameters ${\displaystyle \varphi _{i}}$ where i = 1, ..., p. There is a direct correspondence between these parameters and the covariance function of the process, and this correspondence can be inverted to determine the parameters from the autocorrelation function (which is itself obtained from the covariances). This is done using the Yule–Walker equations.

The power spectral density (PSD) of an AR(p) process with noise variance ${\displaystyle \mathrm {Var} (Z_{t})=\sigma _{Z}^{2}}$ is[8]

${\displaystyle S(f)={\frac {\sigma _{Z}^{2}}{|1-\sum _{k=1}^{p}\varphi _{k}e^{-i2\pi fk}|^{2}}}.}$

• R, the stats package includes an ar function.[14]
• R, the package astsa includes an sarima function to fit various models including AR models.[15]
• MATLAB's Econometrics Toolbox[16] and System Identification Toolbox[17] includes autoregressive models[18]
• Matlab and Octave: the TSA toolbox contains several estimation functions for uni-variate, multivariate and adaptive autoregressive models.[19]
• PyMC3: the Bayesian statistics and probabilistic programming framework supports autoregressive modes with p lags.
• bayesloop supports parameter inference and model selection for the AR-1 process with time-varying parameters.[20]
• Python: implementation in statsmodels.[21]

The impulse response of a system is the change in an evolving variable in response to a change in the value of a shock term k periods earlier, as a function of k. Since the AR model is a special case of the vector autoregressive model, the computation of the impulse response in vector autoregression#impulse response applies here.

Once the parameters of the autoregression

${\displaystyle X_{t}=\sum _{i=1}^{p}\varphi _{i}X_{t-i}+\varepsilon _{t}\,}$

have been estimated, the autoregression can be used to forecast an arbitrary number of periods into the future. First use t to refer to the first period for which data is not yet available; substitute the known preceding values X_t-i for i=1, ..., p into the autoregressive equation while setting the error term ${\displaystyle \varepsilon _{t}}$ equal to zero (because we forecast X_t to equal its expected value, and the expected value of the unobserved error term is zero). The output of the autoregressive equation is the forecast for the first unobserved period. Next, use t to refer to the next period for which data is not yet available; again the autoregressive equation is used to make the forecast, with one difference: the value of X one period prior to the one now being forecast is not known, so its expected value—the predicted value arising from the previous forecasting step—is used instead. Then for future periods the same procedure is used, each time using one more forecast value on the right side of the predictive equation until, after p predictions, all p right-side values are predicted values from preceding steps.

There are four sources of uncertainty regarding predictions obtained in this manner: (1) uncertainty as to whether the autoregressive model is the correct model; (2) uncertainty about the accuracy of the forecasted values that are used as lagged values in the right side of the autoregressive equation; (3) uncertainty about the true values of the autoregressive coefficients; and (4) uncertainty about the value of the error term ${\displaystyle \varepsilon _{t}\,}$ for the period being predicted. Each of the last three can be quantified and combined to give a confidence interval for the n-step-ahead predictions; the confidence interval will become wider as n increases because of the use of an increasing number of estimated values for the right-side variables.

• Moving average model
• Linear difference equation
• Predictive analytics
• Linear predictive coding
• Resonance
• Levinson recursion
• Ornstein–Uhlenbeck process
• Infinite impulse response

• Mills, Terence C. (1990). Time Series Techniques for Economists. Cambridge University Press. ISBN 9780521343398.
• Percival, Donald B.; Walden, Andrew T. (1993). Spectral Analysis for Physical Applications. Cambridge University Press. Bibcode:1993sapa.book.....P.
• Pandit, Sudhakar M.; Wu, Shien-Ming (1983). Time Series and System Analysis with Applications. John Wiley & Sons.

• AutoRegression Analysis (AR) by Paul Bourke
• Econometrics lecture (topic: Autoregressive models) on YouTube by Mark Thoma