Introduction (TS 3.1, 3.3, 3.6, S6) #

Overview of main models #

We are now ready to introduce the main models with full generality, and to discuss an overall estimation, fitting, and forecasting approach. The main models we will review are:

• Autoregressive models of order $p$, denoted $AR(p)$: $$x_t=\alpha +{\varphi }_1{x}_{t-1}+{\varphi }_2{x}_{t-2}+\cdots +{\varphi }_p{x}_{t-p}+w_t,$$ where $\alpha =\mu (1-{\varphi }_1-{\varphi }_2-\cdots -{\varphi }_p)$.
• Moving averages of order $q$, denoted $MA(q)$: $$x_t=w_t+{\theta }_1{w}_{t-1}+{\theta }_2{w}_{t-2}+\cdots +{\theta }_q{w}_{t-q}.$$
• $ARMA(p,q)$ models, which are a combination of the above: $$x_t={\varphi }_1{x}_{t-1}+\cdots +{\varphi }_p{x}_{t-p}+w_t+{\theta }_1{w}_{t-1}+\cdots +{\theta }_q{w}_{t-q}$$

Variations #

• $ARIMA(p,d,q)$ models, which reduce to $ARMA(p,q)$ when differenced $d$ times.
• Multiplicative seasonal ARIMA models:
  • Pure seasonal autoregressive moving average $ARMA(P,Q{)}_s$ models
  • Multiplicative seasonal autoregressive moving average models
  • Seasonal autoregressive integrated moving average models
• Multivariate time series—Vector Auto-Regressive $VAR(p)$ models
• Extra models mentioned by the CS2 syllabus:
  • bilinear model
  • threshold autoregressive model
  • $ARCH(p)$ models

Partial Autocorrelation Function (PACF) #

Motivation #

• There is one more property we need to introduce, the Partial Autocorrelation Function (PACF).
• The ACF can help identify $q$ for a $MA(q)$ model, because it should be 0 for lags $>q$
• However, the order of an $AR(p)$ or $ARMA(p,q)$ model cannot be recognised from the ACF.
• What we need for autoregressive models is a measure of dependence at two different points in time,
  with the effect of the other points in time removed.
• We will remove the linear effect of all variables in-between, and so by increasing the lag iteratively, this partial correlation should be 0 after some lag ( $p$ for $AR(p)$)

Reminder: partial correlation #

• Let $X$, $Y$, and $Z$ be random variables.
• The partial correlation between $X$ and $Y$, given $Z$ is obtained by
  • Regressing $X$ on $Z$ to obtain $\hat{X}$;
  • Regressing $Y$ on $Z$ to obtain $\hat{Y}$,
• Then the partial correlation is $${\rho }_{XY|Z}=\text{corr}\{X-\hat{X},Y-\hat{Y}\}.$$
• In the special case where the random variables are multivariate normal (and hence dependence is linear), $${\rho }_{XY|Z}=\text{corr}\{X,Y|Z\}.$$ (Otherwise, $\hat{X}$ and $\hat{Y}$ are only linear projections, not the full story.)

The partial autocorrelation function (PACF) #

The partial autocorrelation function (PACF) of a stationary process, $x_t$, denoted ${\varphi }_{hh}$, for $h=1,2,\dots$, is $${\varphi }_{11}=\text{corr}(x_{t+1},x_t)=\rho (1),$$ and $${\varphi }_{hh}=\text{corr}(x_{t+h}-{\hat{x}}_{t+h},x_t-{\hat{x}}_t),h\ge 2.$$ Note:

• The PACF, ${\varphi }_{hh}$, is the correlation between $x_{t+h}$ and $x_t$ with the linear dependence of $\{x_{t+1},\dots ,x_{t+h-1}\}$ on each, removed.
• The PACF of an $AR(p)$ model will be 0 for $h>p$.
• The R function pacf will display a sample ${\varphi }_{hh}$ for $h>0$.

Examples of AR PACFs #

set.seed(2)
AR1 <- arima.sim(list(order = c(1, 0, 0), ar = 0.8), n = 10000)
AR2 <- arima.sim(list(order = c(2, 0, 0), ar = c(0.8, -0.1)),
  n = 10000)
AR5 <- arima.sim(list(order = c(5, 0, 0), ar = c(0.8, -0.8, 0.5,
  -0.5, 0.3)), n = 10000)
par(mfrow = c(2, 3))
acf(AR1)
acf(AR2)
acf(AR5)
pacf(AR1)
pacf(AR2)
pacf(AR5)

Examples of MA PACFs #

set.seed(2)
MA1 <- arima.sim(list(order = c(0, 0, 1), ma = 0.8), n = 10000)
MA2 <- arima.sim(list(order = c(0, 0, 2), ma = c(0.8, -0.1)),
  n = 10000)
MA5 <- arima.sim(list(order = c(0, 0, 5), ma = c(0.8, -0.8, 0.5,
  -0.5, 0.3)), n = 10000)
par(mfrow = c(2, 3))
acf(MA1)
acf(MA2)
acf(MA5)
pacf(MA1)
pacf(MA2)
pacf(MA5)

$AR(p)$ models (TS 3.1, 3.3, 3.6, S6) #

Definition #

$AR(p)$ models #

An autoregressive model of order $p$ $AR(p)$ is of the form $$x_t=\alpha +{\varphi }_1{x}_{t-1}+{\varphi }_2{x}_{t-2}+\cdots +{\varphi }_p{x}_{t-p}+w_t,$$ where $x_t$ is stationary, $w_t\sim \text{wn}(0,{\sigma }_w^2)$, ${\varphi }_1,{\varphi }_2,\dots ,{\varphi }_p$ are constants $({\varphi }_p\ne 0)$, and where $\alpha =\mu (1-{\varphi }_1-{\varphi }_2-\cdots -{\varphi }_p)$.

Assume $\alpha =0$ and rewrite $$(1-{\varphi }_{1}B-{\varphi }_2{B}^2-\cdots -{\varphi }_p{B}^p)x_t=w_t\text{or}\varphi (B)x_t=w_t.$$ Note:

• Unless stated otherwise we assume $\alpha =0$. The mean of $x_t$ is then zero.
• Recall that here the regressors $x_{t-1},\dots ,x_{t-p}$ are random components (as opposed to the ‘usual’ regression).
• The autoregressive operator is defined as $$\varphi (B)=(1-{\varphi }_{1}B-{\varphi }_2{B}^2-\cdots -{\varphi }_p{B}^p),$$ the roots of which will be crucial for the analysis.

Example: AR(1) model #

An $AR(1)$ is $$\begin{array}{rcl}{x}_t& =& \varphi x_{t-1}+w_t\\ & =& \varphi (\varphi x_{t-2}+w_{t-1})+w_t\\ & ⋮& \\ & =& {\varphi }^k{x}_{t-k}+\sum _{j=0}^{k-1}{\varphi }^j{w}_{t-j}\end{array}$$ If one continues to iterate backwards, one gets (provided $|\varphi |<1$ and $\underset{t}{sup}Var(x_t)<\mathrm{\infty }$) $$x_t=\sum _{j=0}^{\mathrm{\infty }}{\varphi }^j{w}_{t-j},$$ a linear process! This is called the stationary solution of the model.

The stationary solution of the model has mean $$E[x_t]=\sum _{j=0}^{\mathrm{\infty }}{\varphi }^{j}E[w_{t-j}]=0,$$ and autocovariance function (see (3.7) in TS) $$\gamma (h)=\frac{{\sigma }_w^2{\varphi }^h}{1-{\varphi }^2},h\ge 0.$$ (remember $\gamma (h)=\gamma (-h)$) The ACF of an $AR(1)$ is $$\rho (h)=\frac{\gamma (h)}{\gamma (0)}={\varphi }^h,h\ge 0,$$ and $\rho (h)$ satisfies the recursion $$\rho (h)=\varphi \rho (h-1),h=1,2,\dots .$$

Sample paths of $AR(1)$ processes #

par(mfrow = c(2, 1))
plot(arima.sim(list(order = c(1, 0, 0), ar = 0.9), n = 100),
  ylab = "x", main = (expression(AR(1) ~ ~~phi == +0.9 ~ ~~sigma^2 ==
    1 ~ ~~rho == 0.9^h)))
plot(arima.sim(list(order = c(1, 0, 0), ar = -0.9), n = 100),
  ylab = "x", main = (expression(AR(1) ~ ~~phi == -0.9 ~ ~~sigma^2 ==
    1 ~ ~~rho == (-0.9)^h)))

Stationarity (and Causality) #

Stationarity of $AR(p)$ processes #

The time series process $x_t$ is stationary and causal if, and only if, the roots of the characteristic polynomial $$1-{\varphi }_{1}z-{\varphi }_2{z}^2-\cdots -{\varphi }_p{z}^p=\varphi (z)=0$$ (which can be complex numbers) are all greater than 1 in absolute value (outside the unit circle).

We will get a sense of why as we make our way through this module.

ACF and PACF #

Autocovariances #

Remember that for stationary $x_t=\sum _{j=1}^p{\varphi }_j{x}_{t-j}+w_t$ we have $$\gamma (k)=Cov(\sum _{j=1}^p{\varphi }_j{x}_{t-j}+w_t,x_{t-k})=\sum _{j=1}^p{\varphi }_j\gamma (k-j)\text{for }k\ge 0.$$ This is a $p$-th order difference equation with constant coefficients with solution (see TS 3.2) $$\gamma (k)=\sum _{j=1}^p{A}_j{z}_j^{-k}\text{ for all }k\ge 0,$$ where $z_1,\dots ,z_p$ are the $p$ roots of the characteristic polynomial $\varphi (z)=0$, and $A_1,\dots ,A_p$ are constants depending on the initial values.

R can calculate this for you (code is provided later).

Note:

• we expect $\gamma (k)\to 0$ as $k\to \mathrm{\infty }$ for stationary $x_t$
• this is equivalent to $|z_j|>1$ for all $j$
• this is equivalent to our iff condition for $x_t$ to be stationary and causal

Yule-Walker equations #

• We have a solution form for $\gamma (k)$, but no explicit solution yet. These need to be solved for.
• Sometimes, it is quicker to use the so-called “Yule-Walker equations”. Consider the equation for an $AR(p)$ model with $\mu =0$, multiply by $x_{t-k}$, take expectations to get $$\begin{array}{rcl}Cov(x_t,x_{t-k})& =& {\varphi }_{1}Cov(x_{t-1},x_{t-k})+{\varphi }_{2}Cov(x_{t-2},x_{t-k})+\dots \\ & & +{\varphi }_{p}Cov(x_{t-p},x_{t-k})+Cov(w_t,x_{t-k}),0\le k\le p,\end{array}$$ that is, $$\gamma (k)={\varphi }_1\gamma (k-1)+{\varphi }_2\gamma (k-2)+\dots +{\varphi }_p\gamma (k-p)+{\sigma }_w^2{1}_{\{k=0\}},0\le k\le p.$$ This is solvable thanks to the fact that $\gamma (k)=\gamma (-k)$.
• A matrix representation will be introduced later

Example #

• For an $AR(3)$ model we have $$\begin{array}{rcl}\gamma (3)& =& {\varphi }_1\gamma (2)+{\varphi }_2\gamma (1)+{\varphi }_3\gamma (0),\\ \gamma (2)& =& {\varphi }_1\gamma (1)+{\varphi }_2\gamma (0)+{\varphi }_3\gamma (1),\\ \gamma (1)& =& {\varphi }_1\gamma (0)+{\varphi }_2\gamma (1)+{\varphi }_3\gamma (2),\\ \gamma (0)& =& {\varphi }_1\gamma (1)+{\varphi }_2\gamma (2)+{\varphi }_3\gamma (3)+{\sigma }_w^2.\end{array}$$
• The second and third equations can be solved linearly to obtain expressions for $\gamma (1)$ and $\gamma (2)$ as a constant times $\gamma (0)$, which yields explicitly $\rho (1)$ and $\rho (2)$. ( Remember $\rho (k)=\gamma (k)/\gamma (0)$, $\rho (0)=1$ ) $$\begin{array}{rcl}\rho (2)& =& {\varphi }_1\rho (1)+{\varphi }_2+{\varphi }_3\rho (1),\\ \rho (1)& =& {\varphi }_1+{\varphi }_2\rho (1)+{\varphi }_3\rho (2),+{\sigma }_w^2.\end{array}$$
• For numerical values of $\gamma (1)$ and $\gamma (2)$, and indeed the others, too, the system needs to be solved for.
  The matrix representation will help.

ACF #

As illustrated above $$\rho (h)=\sum _{j=1}^p{\varphi }_j\rho (k-j),$$ because $\rho (h)=\gamma (h)/\gamma (0).$

Note:

• For stationary and causal AR, $|z_i|>1$, $i=1,\dots ,r$ ( $r\le p$ distinct roots)
• If all the roots are real, then $\rho (h)$ dampens exponentially fast to zero as $h\to \mathrm{\infty }$.
• If some of the roots are complex, then they will be in conjugate pairs and $\rho (h)$ will dampen, in a sinusoidal fashion, exponentially fast to zero as $h\to \mathrm{\infty }$. In this case, the time series may appear to be cyclic in nature.
• This property flows on to ARMA models.

PACF #

When $h>p$ the regression of $x_{t+h}$ on $\{x_{t+1},\dots ,x_{t+h-1}\}$ is $${\hat{x}}_{t+h}=\sum _{j=1}^p{\varphi }_j{x}_{t+h-j},h>p.$$ This is not a typo, it is a really nice result! (see TS for a proof).

[why? because this means that ${\hat{x}}_{t+h}=x_{t+h}-w_{t+h}$]

This means that $x_t-{\hat{x}}_t$, which will depend only on $\{w_{t+h-1},w_{t+h-2},\dots \}$, has no overlap with $x_{t+h}-{\hat{x}}_{t+h}=w_{t+h}$

Hence the PACF $${\varphi }_{hh}=\text{corr}(x_{t+h}-{\hat{x}}_{t+h},x_t-{\hat{x}}_t)=\text{corr}(w_{t+h},x_t-{\hat{x}}_t)=0,h>p.$$

In summary, we know that $${\varphi }_{hh}=0,\text{ for all }h>p$$ “by design” (we wanted a measure that would die down beyond $p$ for diagnostic reasons).

“In-between”, we have that $${\varphi }_{11},\dots ,{\varphi }_{pp}$$ are not necessarily 0.

Furthermore, it can be shown that $${\varphi }_{pp}={\varphi }_p,$$ a really nice feature.

Explosive AR models and causality #

Remember the random walk (an $AR(1)$ with $\varphi =1$) $$x_t=x_{t-1}+w_t$$ is not stationary. This was because $x_{t-1}$ includes all past $w_t$’s, leading to an exploding variance.

Examination of the $AR(1)$ representation $$x_t=\sum _{j=0}^{\mathrm{\infty }}{\varphi }^j{w}_{t-j}$$ suggests that the key to contain explosion is with $\varphi$. If $|{\varphi }^j|$ increases without bound as $j\to \mathrm{\infty }$ the expected value of this quantity won’t converge. Explosive processes quickly become large in magnitude, leading to unstationarity.

Unfortunately not all AR models are stationary. For the random walk, $\varphi =1$ and it is not stationary. What are the values of $\varphi$ so that this does not happen?

Assume $|\varphi |>1$. Write (by iterating forward $k$ steps) $$\begin{array}{rcl}{x}_t& =& {\varphi }^{-1}{x}_{t+1}-{\varphi }^{-1}{w}_{t+1}\\ & =& {\varphi }^{-1}({\varphi }^{-1}{x}_{t+2}-{\varphi }^{-1}{w}_{t+2})-{\varphi }^{-1}{w}_{t+1}\\ & ⋮& \\ & =& {\varphi }^{-k}{x}_{t+k}-\sum _{j=1}^{k-1}{\varphi }^{-j}{w}_{t+j}\end{array}$$

Because $|\varphi {|}^{-1}<1$, this result suggests the future dependent $AR(1)$ model $$x_t=-\sum _{j=1}^{\mathrm{\infty }}{\varphi }^{-j}{w}_{t+j}.$$ This is of the $AR(1)$ form $x_t=\varphi x_{t-1}+w_t$, but it is useless because it requires us to know the future to be able to predict the future!

When a process does not depend on the future— such as $AR(1)$ when $|\varphi |<1$— we will say the process is causal.

The model above with $|\varphi |>1$ is stationary, but it is also future dependent, and hence is not causal.

Here is the lesson for $p>1$:

• stationary and causal are not equivalent conditions
• depending on the parameters of your $AR(p)$ model, you may have a future dependent model without knowing it: when $p>1$ it is not obvious by just looking at the parameters ${\varphi }_j$
• that’s why the condition above stated “stationary and causal”
• this is further discussed in the ARMA section, where the argument above is further generalised to $p>1$

$MA(q)$ models (TS 3.1, 3.3, 3.6, S6) #

Definition #

$MA(q)$ models #

The moving average model of order $q$ $MA(q)$ is of the form $$x_t=w_t+{\theta }_1{w}_{t-1}+{\theta }_2{w}_{t-2}+\cdots +{\theta }_q{w}_{t-q}=\theta (B)w_t,$$ where $w_t\sim \text{wn}(0,{\sigma }_w^2)$ and ${\theta }_1,{\theta }_2,\dots ,{\theta }_q$ $({\theta }_q\ne 0)$ are parameters.

Note:

• The AR combines linearly the $x_t$’s, whereas the MA combines linearly the $w_t$’s.
• The moving average operator is defined as $$\theta (B)=1+{\theta }_{1}B+{\theta }_2{B}^2+\cdots +{\theta }_q{B}^q.$$
• Unlike the AR, the MA is stationary for any value of the parameters $\theta$’s

Example: $MA(1)$ model #

An $MA(1)$ is $$x_t=w_t+\theta w_{t-1}.$$ Hence, $$E[x_t]=0,$$ and $$\gamma (h)=\{\begin{array}{lc}(1+{\theta }^2){\sigma }_w^2& h=0,\\ \theta {\sigma }_w^2& h=1,\\ 0& h>1,\end{array}$$ and the ACF is $$\rho (h)=\{\begin{array}{lc}(1& h=0)\\ \frac{\theta }{(1+{\theta }^2)}& h=1,\\ 0& h>1.\end{array}$$ Furthermore $|\rho (1)|\le 1/2$ for all values of $\theta$.

Sample paths of $MA(1)$ processes #

par(mfrow = c(2, 1))
plot(arima.sim(list(order = c(0, 0, 1), ma = 0.5), n = 100),
  ylab = "x", main = (expression(MA(1) ~ ~~theta == +0.5)))
plot(arima.sim(list(order = c(0, 0, 1), ma = -0.5), n = 100),
  ylab = "x", main = (expression(MA(1) ~ ~~theta == -0.5)))

Non-uniqueness and invertibility #

• Note that (multiply the RHS by ${\theta }^2$ to get the original formula): $$\rho (h)=\frac{\theta }{(1+{\theta }^2)}=\frac{1/\theta }{(1+1/{\theta }^2)}.$$
• Furthermore, the pair $(\theta =5,{\sigma }_w^2=1)$ yield the same autocovariance function as the pair $(\theta =1/5,{\sigma }_w^2=25)$, namely, $$\gamma (h)=\{\begin{array}{lc}26& h=0,\\ 5& h=1,\\ 0& h>1.\end{array}$$

• Thus, the $MA(1)$ processes $$x_t=w_t+\frac{1}{5}{w}_{t-1},w_t\sim \text{N}(0,25)$$ and $$y_t=v_t+5v_{t-1},v_t\sim \text{N}(0,1)$$ are the same because of normality (same first two moments, and the $\gamma$ fully determine their dependence structure).
• We can only observe $x_t$ or $y_t$, but not $w_t$ or $v_t$ so we cannot distinguish between them.
• We encountered this phenomenon already in Tutorial Exercise TS5.

For convenience, we will systematically choose the version that is invertible, which means a process that has an infinite AR representation, as explained below.

Example:

• Consider the following inversion of the roles of $x_t$ and $w_t$ in the specific case of an MA(1) model: $$\begin{array}{rcl}{w}_t& =& -\theta w_{t-1}+x_t\\ & =& \sum _{j=0}^{\mathrm{\infty }}(-\theta {)}^j{x}_{t-j}\text{ if }|\theta |<1,\end{array}$$ which is an infinite AR representation of the model.
• Since we need $|\theta |<1$ for this to work, we will choose the version with $(\theta =1/5,{\sigma }_2^2=25)$.

How can we generalise this to $q>1$?

• As in the AR case, the polynomial $\theta (z)$ is key. The inversion in general is $$x_t=\theta (B)w_t⟺\pi (B)x_t=w_t\text{where }\pi (B)={\theta }^{-1}(B)$$
• Just as we required $|\theta |<1$ for , we will

Example: in the $MA(1)$ case, $$\theta (z)=1+\theta z⟺\pi (z)={\theta }^{-1}(z)=\frac{1}{1+\theta z}=\sum _{j=0}^{\mathrm{\infty }}(-\theta {)}^j{z}^j\text{ if }|\theta z|<1.$$ Consequently, $$\pi (B)=\sum _{j=0}^{\mathrm{\infty }}(-\theta {)}^j{B}^j.$$

ACF and PACF #

Autocovariance #

We have $$\begin{array}{rcl}\gamma (h)& =& Cov(x_{t+h},x_t)\\ & =& Cov(\sum _{j=0}^q{\theta }_j{w}_{t+h-j},\sum _{k=0}^q{\theta }_k{w}_{t-k})\\ & =& \{\begin{array}{ll}{\sigma }_w^2\sum _{j=0}^{q-h}{\theta }_j{\theta }_{j+h},& 0\le h\le q,\\ 0& h>q.\end{array}\\ & =& \gamma (-h)\end{array}$$ Note

• $\gamma (q)\ne 0$ because ${\theta }_q\ne 0$
• the cutting off of $\gamma (h)$ after $q$ lags is the signature of the $MA(q)$ model.

ACF #

The ACF is then $$\rho (h)=\{\begin{array}{ll}\frac{\sum _{j=0}^{q-h}{\theta }_j{\theta }_{j+h}}{1+{\theta }_1^2+\cdots +{\theta }_q^2},& 1\le h\le q,\\ 0& h>q.\end{array}$$

PACF #

An invertible $MA(q)$ can be written as $$x_t=-\sum _{j=1}^{\mathrm{\infty }}{\pi }_j{x}_{t-j}+w_t.$$ No finite representation exists, and hence the PACF will never cut off (as opposed to the case of $AR(p)$ ).

Example: For an $MA(1)$, $$\begin{array}{rcl}{\varphi }_{11}& =& \rho (1),\\ {\varphi }_{22}& =& -\frac{{\theta }^2}{1+{\theta }^2+{\theta }^4},\\ {\varphi }_{hh}& =& -\frac{(-\theta {)}^h(1-{\theta }^2)}{1-{\theta }^{2(h+1)}},h\ge 1.\end{array}$$ This is derived in the textbook.

$ARMA(p,q)$ models (TS 3.1, 3.3, 3.6, S6) #

Definition #

A time series $\{x_t;t=0,\pm 1,\pm 2,\dots \}$ is $ARMA(p,q)$ if it is stationary and $$x_t={\varphi }_1{x}_{t-1}+\cdots +{\varphi }_p{x}_{t-p}+w_t+{\theta }_1{w}_{t-1}+\cdots +{\theta }_q{w}_{t-q},$$ with ${\varphi }_p\ne 0$, ${\theta }_q\ne 0$, and ${\sigma }_w^2>0$. The parameters $p$ and $q$ are called the autoregressive and the moving average orders, respectively. If $x_t$ has a nonzero $\mu$, we set $\alpha =\mu (1-{\varphi }_1-\cdots -{\varphi }_p)$ and write the model as $$x_t=\alpha +{\varphi }_1{x}_{t-1}+\cdots +{\varphi }_p{x}_{t-p}+w_t+{\theta }_1{w}_{t-1}+\cdots +{\theta }_q{w}_{t-q},$$ where $w_t\sim \text{wn}(0,{\sigma }_w^2)$.

This can be rewritten very concisely as $$\varphi (B)x_t=\theta (B)w_t.$$

Parameter redundancy #

Multiplying both sides by a third operator $\eta (B)$ leads to $$\eta (B)\varphi (B)x_t=\eta (B)\theta (B)w_t,$$ which doesn’t change the dynamics. One could then have redundant parameters. For example:

• Consider the white noise process $x_t=w_t$, which is “$ARMA(0,0)$”
• If we multiply both sides by $\eta (B)=1-0.5B$ then the model becomes $$x_t=0.5x_{t-1}-0.5w_{t-1}+w_t,$$ which looks like $ARMA(1,1)$, even though it is still white noise.
• We have hidden that fact through over-parametrisation, leading to parameter redundancy.
• Fitting results are likely to suggest that all parameters (including
  the unnecessary ones) are significant.

The AR and MA polynomials #

The AR and MA polynomials are defined as $$\varphi (z)=1-{\varphi }_{1}z-\cdots -{\varphi }_p{z}^p,{\varphi }_p\ne 0,$$ and $$\theta (z)=1+{\theta }_{1}z+\cdots +{\theta }_q{z}^q,{\theta }_q\ne 0,$$ respectively, where $z$ may be a complex number.

These will be used extensively to study the properties of ARMA processes.

Properties #

Problems to keep in mind #

In the previous sections, we identified the following potential “problems”, or potential “issues to keep in mind”:

• parameter redundant models,
• stationary AR models that depend on the future, and
• MA models that are not unique.

We now discuss how to avoid those issues.

Parameter redundant models #

We will require:

• $\varphi (z)$ and $\theta (z)$ cannot have common factors

This will ensure that when referring to an $ARMA(p,q)$ model, there can’t be a reduced form of it.

Future dependence and causality #

[What does causal mean?] An $ARMA(p,q)$ model is said to be causal, if the time series $\{x_t;t=0,\pm 1,\pm 2,\dots \}$ can be written as a one-sided linear process: $$x_t=\sum _{j=0}^{\mathrm{\infty }}{\psi }_j{w}_{t-j}=\psi (B)w_t,$$ where $\psi (B)=\sum _{j=0}^{\mathrm{\infty }}{\psi }_j{B}^j$, and $\sum _{j=0}^{\mathrm{\infty }}|{\psi }_j|<\mathrm{\infty }$; we set ${\psi }_0=1$.

Note:

• key here is that the sum starts at $j=0$ (so that it depends only on the past);
• the $\psi$ parameters are new ones, which are defined in the equation above (“such that” $x_t$ can be written as such a sum).

[When is it causal?] An $ARMA(p,q)$ model is causal if and only if $\varphi (z)\ne 0$ for $|z|\le 1$.

To see this, note that the coefficients of the linear process given above can be determined by solving $$\psi (z)=\sum _{j=0}^{\mathrm{\infty }}{\psi }_j{z}^j=\frac{\theta (z)}{\varphi (z)},|z|\le 1.$$

Equivalently, an ARMA process is causal only when the roots of $\varphi (z)$
lie outside the unit circle, that is, $\varphi (z)=0$ only when $|z|>1$.

Invertibility #

[What does invertible mean?] An $ARMA(p,q)$ model is said to be invertible, if the time series $\{x_t;t=0,\pm 1,\pm 2,\dots \}$ can be written as $$\pi (B)x_t=\sum _{j=0}^{\mathrm{\infty }}{\pi }_j{x}_{t-j}=w_t,$$ where $\pi (B)=\sum _{j=0}^{\mathrm{\infty }}{\pi }_j{B}^j$, and $\sum _{j=0}^{\mathrm{\infty }}|{\pi }_j|<\mathrm{\infty }$; we set ${\pi }_0=1$.

Note:

• key here is that we express the model as “$w_t$ is a function of $x_t$’s” even though our focus (what we model and observe) is $x_t$. That’s the “invertible” idea;
• the $\pi$ parameters are new ones, which are defined in the equation above (“such that” $w_t$ can be written as such a sum of $x_t$’s).

[When is it invertible?] An $ARMA(p,q)$ model is invertible if and only if $\theta (z)\ne 0$ for $|z|\le 1$. The coefficients ${\pi }_j$ of $\pi (B)$ given above can be determined by solving $$\pi (z)=\sum _{j=0}^{\mathrm{\infty }}{\pi }_j{z}^j=\frac{\varphi (z)}{\theta (z)},|z|\le 1.$$ Equivalently, an ARMA process is invertible only when the roots of
$\theta (z)$ lie outside the unit circle; that is, $\theta (z)=0$ only when $|z|>1$.

Example: parameter redundancy, causality, invertability #

Consider the process: $$x_t=0.4x_{t-1}+0.45x_{t-2}+w_t+w_{t-1}+0.25w_{t-2}.$$ The first step here is to express that with the help of the backshift operator: $$(1-0.4B-0.45B^2)x_t=(1+B+0.25B^2)w_t$$ This looks like an $ARMA(2,2)$ process, but there’s a trick! Write the AR and MA polynomials, and try to factorise them: $$\begin{array}{rcl}\varphi (B)& =& 1-0.4B-0.45B^2=(1+0.5B)(1-0.9B)\\ \theta (B)& =& 1+B+0.25B^2=(1+0.5B{)}^2\end{array}$$ There is a common factor, leading to parameter redundancy!

Factorise the common factor $(1+0.5B)$ out, and get $$\begin{array}{rcl}\varphi (B)& =& 1-0.9B\\ \theta (B)& =& 1+0.5B\end{array}$$ so that our model is in fact $$x_t=0.9x_{t-1}+0.5w_{t-1}+w_t,$$ an $ARMA(1,1)$ model.
[Parameter redundancy: checked!]

We next check whether the process is causal. We need the root of $$\varphi (z)=1-0.9z=0$$ to be outside the unit circle, which it is as the solution is $z=10/9>1$. [Causality: checked!]

We next check wether the model is invertible. We need the root of $$\theta (z)=1+0.5z=0$$ to be outside the unit circle, which it is as the solution is $z=-2$. [Invertibility: checked!]

Now, we would like to find the linear representation of the process, that is, get the $\psi$ weights. Because $$\begin{array}{rcl}\varphi (z)\psi (z)& =& \theta (z)\\ (1-0.9z)(1+{\psi }_{1}z+{\psi }_2{z}^2+\cdots +{\psi }_j{z}^j+\cdots )& =& 1+0.5z\\ (\text{regrouping coefficients of powers of }z)& & \\ 1+({\psi }_1-0.9)z+({\psi }_2-0.9{\psi }_1)z^2+& & \\ \cdots +({\psi }_j-0.9{\psi }_{j-1})z^j+\cdots & =& 1+0.5z\end{array}$$ We compare coefficients of the powers of $z$, and note that all coefficients of $z^j$ are 0 for $j>1$ on the RHS.

We obtain then $$\begin{array}{rcl}{\psi }_1-0.9=0.5& ⟹& {\psi }_1=1.4,\\ {\psi }_j-0.9{\psi }_{j-1}=0& ⟹& \frac{{\psi }_j}{{\psi }_{j-1}}=0.9,j>1.\end{array}$$ and thus $${\psi }_j=1.4(0.9{)}^{j-1}\text{ for }j\ge 1,$$ and hence $$x_t=w_t+1.4\sum _{j=1}^{\mathrm{\infty }}{0.9}^{j-1}{w}_{t-j}.$$ In R, this is much quicker! just use ( $x_t=0.9x_{t-1}+0.5w_{t-1}+w_t$ )

format(ARMAtoMA(ar = 0.9, ma = 0.5, 10), digits = 2)  # first 10 psi-weights
##  [1] "1.40" "1.26" "1.13" "1.02" "0.92" "0.83" "0.74" "0.67" "0.60"
## [10] "0.54"

[Linear representation: checked!]

(no, it’s not over yet!)

Next, we want to determine the invertible representation of the model. Because $$\begin{array}{rcl}\theta (z)\pi (z)& =& \varphi (z)\\ (1+0.5z)(1+{\pi }_{1}z+{\pi }_2{z}^2+{\pi }_3{z}^3+\cdots )& =& 1-0.9z\\ (\text{regrouping coefficients of powers of }z)& & \\ 1+({\pi }_1+0.5)z+({\pi }_2+0.5{\pi }_1)z^2+& & \\ \cdots +({\pi }_j+0.5{\pi }_{j-1})z^j+\cdots & =& 1-0.9z\end{array}$$ We compare coefficients of the powers of $z$.

We get $${\pi }_j=(-1{)}^{j}1.4(0.5{)}^{j-1}\text{ for }j\ge 1$$ and then $$w_t=\sum _{j=0}^{\mathrm{\infty }}{\pi }_j{x}_{t-j}=x_t+\sum _{j=1}^{\mathrm{\infty }}{\pi }_j{x}_{t-j}$$ so that $$x_t=-\sum _{j=1}^{\mathrm{\infty }}{\pi }_j{x}_{t-j}+w_t.$$

Again, this is much quicker in R: ( $w_t=-0.5w_{t-1}-0.9x_{t-1}+x_t$ )

format(ARMAtoMA(ar = -0.5, ma = -0.9, 10), digits = 1)  # first 10 pi-weights
##  [1] "-1.400" " 0.700" "-0.350" " 0.175" "-0.087" " 0.044" "-0.022"
##  [8] " 0.011" "-0.005" " 0.003"

[Invertible representation: checked!]

Stationarity, Causality and Invertibility #

Wrapping it up #

First, it is helpful to rewrite the ARMA representation as $$x_t=\frac{\theta (B)}{\varphi (B)}{w}_t=\psi (B)w_t$$ To summarise:

• first, we require $\theta (B)$ and $\varphi (B)$ to not have common factors. If they do, these will obviously cancel out in the ratio $\theta (B)/\varphi (B)$ and we, really, are only dealing with the “reduced” model (without the redundant parameters).
• pure $AR(p)$ models will be stationary as long as $\psi (B)$ is well behaved (say, finite), which will happen as long as the roots of $\varphi (B)$ are all greater than one in modulus $(|z_j|>1)$; this is because $\varphi (B)$ appears in the denominator. Since this is not impacted by $\theta (B)$, which is in the numerator, this result also holds in the case of $ARMA(p,q)$
  models.

• pure $MA(q)$ models—where $\varphi (B)=1$—are always stationary (under some mild conditions on the coefficients which we ignore here), because by definition they include a finite number of $w$’s (all covariances are finite). They are also causal by definition.
• now, even though $MA(q)$ are always causal, establishing causality for $AR(p)$ and $ARMA(p,q)$ models—where $p>0$—is a little trickier. We not only need to to be well behaved, but we also need it to depend on past values only (an example of stationary but non-causal model is provided on slide 26 of Module 9 above; see also next subsection). It turns out that you only need the roots to not be on the unit circle for the process to be stationary. For it to be causal, the additional requirement is that it needs to be outside the unit circle. In other words, while you could have roots inside the unit circle and still achieve stationarity, the process would not be causal. This means that causality will imply stationarity, but not the other way around.

• Again, stationarity does not require the roots of $\theta (B)$ to be greater than one in modulus. This is required for invertibility, which aims to flip things around (“invert” the process): $$w_t=\frac{\varphi (B)}{\theta (B)}{x}_t.$$ It becomes clear why, now, it is the roots of $\theta (B)$ that need to be well behaved (as it is now $\theta (B)$ that is in the denominator).

Revisiting the future-dependent example #

Let us revisit that example: We have $$x_t=\varphi x_{t-1}+w_t.$$ In this case the characteristic equation is $$\varphi (z)=1-\varphi z.$$ This has root $$z_0=1/\varphi .$$

We distinguish three cases:

1. $\varphi <1$ (say, $\varphi =0.5$ so that $z_0=1/\varphi =2$): the root is not on the unit circle, and also outside the unit circle, so that the process—which is $AR(1)$—is stationary and causal.
2. $\varphi =1$ $(z_0=1)$: the root is on the unit circle, and hence the process is not stationary. This is the random walk case. Note it is not causal either, because a causal process is a linear process that depends on the past only, and while the random walk depends on the present and past only, it does not satisfy the requirement (for it to be a linear process) that the sum of the absolute value of weights is finite (the sum of an infinite number of 1’s is infinite) - see Module 7.
3. $\varphi >1$ (say, $\varphi =3$ so that $z_0=1/3<1$): the root is not on the unit circle, and hence the process is stationary. However, the root is inside the unit circle, which implies it is not causal. The process will depend on the future as demonstrated earlier.

ACF and PACF #

Autocovariance and ACF #

For a causal $ARMA(p,q)$ model $$\varphi (B)x_t=\theta (B)w_t$$ we use the linear representation $$x_t=\sum _{j=0}^{\mathrm{\infty }}{\psi }_j{w}_{t-j},$$ from which it follows immediately that $E[x_t]=0$, and the autocovariance function of $x_t$ is $$\gamma (h)=Cov(x_{t+h},x_t)={\sigma }_w^2\sum _{j=0}^{\mathrm{\infty }}{\psi }_j{\psi }_{j+h},h\ge 0.$$ This approach requires solving for the $\psi$’s.

Alternatively, it is possible to write a general homogeneous equation for the ACF of a causal ARMA process to solve for the $\gamma$’s directly
(The proof is outside scope but available in TS): $$\gamma (h)-{\varphi }_1\gamma (h-1)-\cdots -{\varphi }_p\gamma (h-p)=0,h\ge max(p,q+1),$$ with initial conditions $$\gamma (h)-\sum _{j=1}^p{\varphi }_j\gamma (h-j)={\sigma }_w^2\sum _{j=h}^q{\theta }_j{\psi }_{j-h},0\le h<max(p,q+1).$$ Finally, the ACF is $$\rho (h)=\frac{\gamma (h)}{\gamma (0)}.$$ In general, the ACF cannot distinguish between AR and ARMA, which is why PACF is useful (in presence of pure AR, it will cut off).

Example: ACF of $ARMA(1,1)$ #

Consider the $ARMA(1,1)$ process $$x_t=\varphi x_{t-1}+\theta w_{t-1}+w_t,\text{ where }|\varphi |<1.$$ The autocovariance then satisfies $$\gamma (h)-\varphi \gamma (h-1)=0,h=2,3,\dots ,$$ which has general solution $$\gamma (h)=c{\varphi }^h,h=1,2,\dots .$$ Initial conditions are $$\begin{array}{rcl}\gamma (0)& =& \varphi \gamma (1)+{\sigma }_w^2({\theta }_0{\psi }_0+{\theta }_1{\psi }_1)=\varphi \gamma (1)+{\sigma }_w^2[1+\theta \varphi +{\theta }^2]\\ \gamma (1)& =& \varphi \gamma (0)+{\sigma }_w^2\theta .\end{array}$$ Note that ${\psi }_1=\theta +\varphi$ for an $ARMA(1,1)$ model
(see Example 3.12 in the textbook).

Solving for $\gamma (0)$ and $\gamma (1)$ yields $$\begin{array}{rcl}\gamma (0)& =& {\sigma }_w^2\frac{1+2\theta \varphi +{\theta }^2}{1-{\varphi }^2}\\ \gamma (1)& =& {\sigma }_w^2\frac{(1+\theta \varphi )(\varphi +\theta )}{1-{\varphi }^2}.\end{array}$$ Since $\gamma (1)=c\varphi$, we get $c=\gamma (1)/\varphi$ and $$\gamma (h)=\frac{\gamma (1)}{\varphi }{\varphi }^h={\sigma }_w^2\frac{(1+\theta \varphi )(\varphi +\theta )}{1-{\varphi }^2}{\varphi }^{h-1},h\ge 1,$$ from which we obtain the ACF $$\rho (h)=\frac{\gamma (h)}{\gamma (0)}=\frac{\gamma (1)}{\gamma (0)}{\varphi }^{h-1}=\frac{(1+\theta \varphi )(\varphi +\theta )}{1+2\theta \varphi +{\theta }^2}{\varphi }^{h-1},h\ge 1.$$ This has the same pattern as an $AR(1)$.

$ARIMA(p,d,q)$ models (TS 3.1, 3.3, 3.6, S6) #

Definition #

A process $x_t$ is said to be $ARIMA(p,d,q)$ if $${\mathrm{\nabla }}^d{x}_t=(1-B{)}^d{x}_t$$ is $ARMA(p,q)$. In general, we will write the model as $$\varphi (B)(1-B{)}^d{x}_t=\theta (B)w_t.$$ If $E[{\mathrm{\nabla }}^d{x}_t]=\mu$, we write the model as $$\varphi (B)(1-B{)}^d{x}_t=\delta +\theta (B)w_t,$$ where $$\delta =\mu (1-{\varphi }_1-\cdots -{\varphi }_p).$$ The integrated ARMA, or ARIMA, is a broadening of the class of
ARMA models to include differencing.

Remarks #

• Because of nonstationarity, care must be taken when deriving forecasts. It is best to use so-called state-space models for handling nonstationary models (but these are outside the scope of this course).
• Since $$y_t={\mathrm{\nabla }}^d{x}_t$$ is ARMA, it suffices to discuss how to fit and forecast ARMA models. For instance, if $d=1$, given forecast $y_{n+m}^n$ for $m=1,2,\dots$, we have $$y_{n+m}^n=x_{n+m}^n-x_{n+m-1}^n\text{ so that }{x}_{n+m}^n=y_{n+m}^n+x_{n+m-1}^n$$ with initial condition $x_{n+1}^n=y_{n+1}^n+x_n$ (noting $x_n^n=x_n$).
• Derivation of prediction errors is slightly more involved (but also outside the scope of this course).

The IMA(1,1) model and exponential smoothing #

IMA(1,1) and EWMA #

• The $ARIMA(0,1,1)$, or $IMA(1,1)$ model is of interest because many economic time series can be successfully modelled this way.
• The model leads to a frequently used— and abused! —forecasting method called exponentially weighted moving averages (EWMA) (“exponential smoothing” in S6).
• One should not use this method unless there is strong statistical evidence that this is the right model to use!
  (as per the following sections)

EWMA example #

set.seed(666)
x <- arima.sim(list(order = c(0, 1, 1), ma = -0.8), n = 100)  #simulate IMA(1,1)
x.ima <- HoltWinters(x, beta = FALSE, gamma = FALSE)  # fit EWMA
plot(x.ima)

Multiplicative seasonal ARIMA models (TS 3.9) #

Introduction #

• We observed seasonality in many time series so far. This is manifested when the dependence on the past tends to occur most strongly at multiples of some underlying seasonal lag $s$. For example:
  • monthly economic data: strong yearly component occurring at lags that are multiples of $s=12$;
  • data taken quarterly: will exhibit the yearly repetitive period at $s=4$ quarters.
• In this section, we introduce several modifications made to the ARIMA model to account for seasonal and nonstationary behavior.
• The reference for this section is TS 3.9.
• Note that seasonality can be due to:
  • Deterministic reasons: include in trend
  • Stochastic dependence at seasonal lags $s$: use SARIMA models

Pure seasonal $ARMA(P,Q{)}_s$ models #

The pure seasonal autoregressive moving average model $ARMA(P,Q{)}_s$ takes the form $${\Phi }_P(B^s)x_t={\Theta }_Q(B^s)w_t,$$ where $${\Phi }_P(B^s)=1-{\Phi }_1{B}^s-{\Phi }_2{B}^{2s}-\cdots -{\Phi }_P{B}^{Ps}$$ is the seasonal autoregressive operator of order $P$ with seasonal period $s$, and where $${\Theta }_Q(B^s)=1+{\Theta }_1{B}^s+{\Theta }_2{B}^{2s}+\cdots +{\Theta }_Q{B}^{Ps}$$ is the seasonal moving average operator of order $Q$ with seasonal period $s$.

Here, inter-temporal correlations are exclusively seasonal (there are no non-seasonal correlations).

Note: (analogous to nonseasonal ARMA)

• It will be causal only when the roots of ${\Phi }_P(z^s)$ lie outside the unit circle
• It will be invertible only when the roots of ${\Theta }_Q(z^s)$ lie outside the unit circle

Example: first order seasonal $(s=12)$ AR model #

A first-order seasonal autoregressive series that might run over months could be written as $$(1-\Phi B^{12})x_t=w_t\text{or}{x}_t=\Phi x_{t-12}+w_t.$$` Note:

• $x_t$ is expressed in terms of past lags at multiple of the (yearly) seasonal period $s=12$ months
• very similar to the unit lag model $AR(1)$ that we know
• causal if $|\Phi |<1$
  Simulated example (with $\Phi =0.9$):

set.seed(666)
phi <- c(rep(0, 11), 0.9)
sAR <- arima.sim(list(order = c(12, 0, 0), ar = phi), n = 37)
sAR <- ts(sAR, freq = 12)