ARCH models

Warning: package 'ggplot2' was built under R version 4.3.2

We start by introducing the first parameteric volatility model, the AutoRegressive Conditional Hetereskedasticity model (ARCH). This model was first introduce by Engle (1982). When writing this, his article has 32 369 citations on Google Scholar, and Robert Engle was awarded the Nobel Memorial Prize in Economic Sciences in 2003 with Clive Granger for “for methods of analyzing economic time series with time-varying volatility (ARCH)” (see nobelprize.org).

ARCH(1)

Let $y_t$ denote the time series of interest, $\sigma_t$ the volatility and $\varepsilon_t \sim \text{iid WN}(0,1)$. Then an ARCH(1) model can be written as $$\begin{equation*} \begin{split} y_t &= \sigma_t \varepsilon_t,\\ \sigma_t^2 &= \omega + \alpha_1 y_{t-1}^2, \quad \omega>0,\, \alpha_1\ge0, \end{split} \end{equation*}$$ where $t$ is denotes the time index. As you can see, the conditional squared volatility depends on the previous squared observations. This means if we get a large return, either positive or negative, the conditional variance will be large in the next time step. With a high conditional variance the we would expect to see a large return, in either direction, also for the next time step, leading to a volatility cluster. With only one lagged variable in the model, this can “die out” quite quickly, but we will use this simple model to do some more detailed calculations before extending the model.

Moments and stationarity

Note that we assume here that $y_t$ has zero expectation, because $$\mu=\mathbb E\, y_t = \mathbb E( \underbrace{\sigma_t\, \varepsilon_t}_\text{independent}) = \mathbb E\,\sigma_t \underbrace{\mathbb E\,\varepsilon_t}_{=0} = 0$$ We thus have that the variance of $y_t$ is

$$\mathrm{Var}(y_t) = \underbrace{\mathbb E(y_t-\mu)^2}_{\mu=0} = \mathbb E\,y_t^2 = \mathbb E(\sigma_t^2\varepsilon_t^2) = \mathbb E\,\sigma_t^2 \underbrace{\mathbb E\,\varepsilon_t^2}_{=1}=\mathbb E\,\sigma_t^2$$

Let’s assume that $y_t$ is stationary, implying that $\mathrm{Var}(y_t) = \mathbb E\,\sigma_t^2 = \sigma^2$. We then have that $$\mathbb E\,\sigma_t^2 = \omega+\alpha_1\mathbb E\,y_{t-1}^2 = \omega +\alpha_1 \mathbb E\,\sigma_{t-1}^2.$$ Since $\mathbb E\,\sigma_t^2 = \mathbb E\,\sigma_{t-1}^2 = \sigma^2$, we get that $$\sigma^2 = \omega+\alpha_1\sigma^2 \quad\Leftrightarrow\quad \sigma^2 = \frac{\omega}{1-\alpha_1}.$$ That is, if $y_t$ is stationary the (unconditional) variance is $\omega/(1-\alpha_1)$, requiring $\alpha_1\in[0,1)$. This is in fact a necessary and sufficient condition for $y_t$ being stationary.

ARCH vs AR

There is a strong link between ARCH models for volatility an AR models for the level. In fact, we can write a ARCH(1) model as an AR(1) model. By adding and subtracting $\sigma_t^2$, we get $$y_t^2 = \sigma_t^2\, \varepsilon_t^2 = \sigma_t^2 + \underbrace{\sigma_t^2(\varepsilon_t^2-1)}_{v_t}=\sigma_t^2 + v_t.$$ It is easy to show that $\mathbb E\,v_t =0$ and under certain conditions ($\mathbb E\,y_t^4<\infty$) $v_t$ will be a white noise process. We can then write $$y_t^2 = \sigma_t^2+v_t = \omega+\alpha_1 y_{t-1}^2 +v_t\quad \text{or} \quad (1-\alpha_1\,B)y_t^2 = \omega+v_t.$$ Thus, we have written the squared ARCH(1) as an AR(1) time series with intercept.

The $\varepsilon_t$ is often referred to as the innovations of the ARCH process, similar to innovation residuals in the fable environment. We can write a ARCH(1) only in terms of the parameters and the infinite past innovations. $$\sigma_t^2 = \omega + \alpha_1 y_{t-1}^2= \omega+\alpha_1\sigma_{t-1}^2\varepsilon_{t-1}^2=\omega+\alpha_1(\omega+\alpha_1 \sigma_{t-2}^2\varepsilon_{t-2}^2)\varepsilon_{t-1}^2.$$ If we continue this iteration backwards, we end up with $$\sigma_t^2=\cdots = \omega+\omega\sum_{j=1}^\infty\alpha_1^j\varepsilon_{t-j}^2.$$ This representation is actually the MA$(\infty)$ representation corresponding to the AR(1). One can take expectations on both sides of this equation and derive the formula for the variance of a stationary ARCH(1) process.

ARCH(r)

Similar to the relationship between an AR(1) and an AR(p), we can set up a ARCH(r). We use $r$ to denote the order of the ARCH model to distinguish it from the p and q of ARIMA, but in the literature it is quite common to use p to denote the order of an ARCH model.

An ARCH(r) model can be set up as $$y_t = \sigma_t\,\varepsilon_t,\quad \sigma_t^2 = \omega + \alpha_1\,y_{t-1}^2+\alpha_2\,y_{t-2}^2+\cdots \alpha_r\,y_{t-r}^2=\omega + \sum_{j=1}^r\alpha_jy_{t-j}^2,$$ where $\omega>0,\,\alpha_j\ge0, j=1,\ldots, r.$

Including more lags, will make the ARCH model’s volatility more persistent, because a large return will stay in the model for some time. We could do similar calculations as above to show that the expectation of $y_t$ is zero and if $y_t$ is stationary the variance is given by $$\sigma^2= \frac{\omega}{1-\sum_{j=1}^r\alpha_j},\quad \sum_{j=1}^r\alpha_j<1$$ You can also find an AR(r) representation of $y_t^2$ using the same approach and assumptions as above.

Simulation

Let us simulate an ARCH(1) model by basic R code, assuming Gaussian innovations.

library(fpp3)
library(tidyverse)
# Setting seed for reproduciblity
set.seed(12345)
nT = 1000
omg = 1.4
a1 = 0.95
# Initiating y and sigma
y <- sig <- rep(sqrt(omg/(1-a1)),
                nT)
# Simulation: 
for(t in 2:nT){
  sig[t] <- sqrt(omg+a1*y[t-1]^2)
  y[t]   <- sig[t]*rnorm(1)
}
# tsibble object: 
arch <- tsibble(
  t = 1:nT,
  y = y,
  sig = sig,
  index = t
) 
# Plotting: 
arch %>% pivot_longer(-t) %>% 
  mutate(name = factor(name, levels = c("y","sig"))) %>% 
  ggplot(aes(x=t,y=value)) + 
  geom_line()+ 
  facet_wrap(~name, scales = "free_y", nrow=2, 
             strip.position = "left")+
  theme(axis.title =element_blank())

In the top panel you see the time series (think of the Microsoft returns from the previous chapter). When $\sigma_t$ is high, you see high variation in $y_t$. When $\sigma_t$ is high, it stays high for some time, leading to volatility clusters.

Let’s look at the autocorrelation of the simulated time series and its squared values. We also do a qq-plot and a distribution plot.

arch %>% ACF(y) %>% autoplot() + labs(y = "ACF of y")
arch %>% ACF(y^2) %>% autoplot() + labs(y = "ACF of squared y")
arch %>%  ggplot(aes(sample = y)) + geom_qq() +geom_qq_line()
arch %>%  ggplot(aes(x = y)) + stat_density(fill = "blue", alpha = .2)

The patterns are very similar to the pattern we saw for Microsoft returns - as is the main goal of these models!

References

• Engle, R. F. (1982). Autoregressive conditional heteroscedasticity with estimates of the variance of United Kingdom inflation. Econometrica: Journal of the econometric society, 987-1007.