Introduction

In this brief tutorial we will go over the ARCH model for measuring volatility in time-series data. The ARCH model has become an exteremly significant econometric model because it is able to capture the features of real-world volatility.

1. The ARCH Model

ARCH stands for auto-regressive conditional heteroskedasticity. To understand the ARCH model, let us first start with a discussion of the concepts of conditional and unconditional means and variances of the error terms in time-series data. First, consider the auto-regressive AR(1) model

\[ \begin{aligned} & y_{t} = \phi + e_{t} &&&& (1.1a) \\\\ & e_{t} = \rho_{t-1} + v_{t} &&&& (1.1b) \\\\ & v_{t} \sim N(0, \sigma_{v}^{2}) &&&& (1.1c) \end{aligned} \]

For convenience of clarifying the concept, first perform some successive substitution to obtain \(e_{t}\) as the sum of an infinite series of the error term \(v_{t}\). To do so, note that if \(e_{t} = \rho e_{t-1} + v_{t}\), then \(e_{t-1} = \rho e_{t-2} + v_{t-1}\) and \(e_{t-2} = \rho e_{t-3} + v_{t - 2}\) and so on. Hence, \(e_{t} = v_{t} + \rho v_{t-1} + \rho^{2} v_{t-2} + \ldots + \rho^{t} e_{0}\) is negligible.

The unconditional mean of the error term is

\[ \begin{aligned} \mathbb{E}(e_{t}) = \mathbb{E}[v_{t} +\rho v_{t-1} + \rho^{2} v_{t-2} + \ldots] = 0 \end{aligned} \]

because \(\mathbb{E}[v_{t-j}] = 0\) for all \(j\), where as the conditional mean for the error is

\[ \begin{aligned} \mathbb{E}[e_{t} \vert I_{t-1}] = \mathbb{E}[\rho e_{t-1} \vert I_{t-1}] + \mathbb{E}[v_{t}] = \rho e_{t-1} \end{aligned} \] because the information set at time \(t-1\), \(I_{t-1}\) includes knowing \(\rho e_{t-1}\). Simply put, unconditional describes a situation in which there is no information, whereas conditional describes a situation in which there is information up to a certain point in time.

The unconditional variance of the error term is

\[ \begin{aligned} \mathbb{E}[e_{t} - 0]^{2} & = \mathbb{E}[v_{t} + \rho v_{t-1} + \rho^{2} v_{t-2} + \ldots]^{2} \\\\ & = \mathbb{E}[v^{2}_{t} + \rho^{2} v^{2}_{t-1} + \rho^{4} v_{t-2}^{4} + ...] \\\\ & = \sigma^{2}_{v}[1 + \rho^{2} + \rho^{4} + \ldots] = \frac{\sigma_{v}^{2}}{1 - \rho^{2}} \end{aligned} \]

because \(\mathbb{E}[v_{t-j}v_{t-i}] = \sigma^{2}_{v}\) when \(i = j\); \(\mathbb{E}[v_{t-j}v_{t-i}] = 0\) when \(i \neq j\) and the sum of the geometric series \([1 + \rho^{2} + \rho^{4} + \ldots]\) is \(1/(1 - \rho^{2})\). The conditional variance for the error term is

\[ \begin{aligned} \mathbb{E}\big[(e_{t} - \rho e_{t-1})^{2} \vert I_{t-1}\big] = \mathbb{E}\big[v_{t}^{2} \vert I_{t-1}\big] = \sigma_{v}^{2} \end{aligned} \]

Notice that for the ARCH model, the conditional mean of the error term varies over time, whereas the conditional variance does not. Now, suppose instead that rather than a conditional mean that changes over time, we have a conditional variance that does change over time. To introduce modification, let us consider a variant of models (1.1a), (1.1b) and (1.1c)

\[ \begin{aligned} & y_{t} = \phi + e_{t} &&&& (1.2a) \\\\ & e_{t} \vert I_{t-1} \sim N(0, h_{t}) &&&& (1.2b) \\\\ & h_{t} = \alpha_{0} + \alpha_{1}e^{2}_{t-1} &&&& (1.2c) \end{aligned} \]

Equations (1.2b) and (1.2c) describe the autoregressive conditional heteroskedastic (ARCH) class of models. Equation (1.2b) states that the error term is conditionally normal \(e_{t} \vert I_{t-1} \sim N(0, h_{t})\), where \(I_{t-1}\) represents the information available at time \(t-1\) with mean \(0\) and time-varying variance, denoted as \(h_{t}\).

The name of the ARCH model implies that the model works with time-varying variances (i.e. heterskedasticity) that depend on lagged effects (i.e. autocorrelation). Our example is an ARCH(1) model since the time-varying variance \(h_{t}\) is a function of a constant term \(\alpha_{0}\) plus one lagged term, which is the square of the error in the previous period \(\alpha_{1}e^{2}_{t-1}\). The coefficients \(\alpha_{0}\) and \(\alpha_{1}\) must be positive to ensure a positive variance. The coefficient \(\alpha_{1}\) must be \(\lt1\). Otherwise, \(h_{t}\) will continue to increase over time, eventually exploding (kabooooom!). Conditional normality means that the normal distribution is a function of known information at time \(t-1\). That is, when \(t=2\), \(e_{2} \vert I_{1} \sim N(0, \alpha_{0} + \alpha_{1}e^{2}_{1})\) and when \(t=3\), then \(e_{3} \vert I_{2} \sim N(0, \alpha_{0} + \alpha_{1}e_{2}^{2})\), and so on. It is also important to note that while the conditional distribution of the error term \(e_{t}\) is assumed to be normal, the unconditional distribution of the error \(e_{t}\) will not. However, this is an inconsequential consideration given thatn many real-world data ppear to be drawn from non-normal distributions.

To find the mean and variance of the unconditional distribution of \(e_{t}\), we note that, conditional on \(e^{2}_{t-1}\), the standardized errors are standard normal. In other words,

\[ \begin{aligned} \Bigg(\frac{e_{t}}{\sqrt{h_{t}}} \Bigg\vert I_{t-1}\Bigg) = z \sim N(0,1) \end{aligned} \]

Since this distribution does not depend on \(e_{t-1}^{2}\), it follows that the unconditional distribution of \(z = (e_{t}/\sqrt{h_{t}})\) is also \(N(0,1)\), and that \(z\) and \(e^{2}_{t-1}\) are independent. Thus, we can write

\[ \begin{equation} \mathbb{E}(e_{t}) = \mathbb{E}(z_{t})\mathbb{E} \Bigg(\sqrt{\alpha_{0} + \alpha_{1}e_{t-1}^{2}}\Bigg) \end{equation} \]

and

\[ \begin{aligned} \mathbb{E}(e_{t}^{2}) = \mathbb{E}(z_{t}^{2}) \mathbb{E}(\alpha_{0} + \alpha_{1}e^{2}_{t-1}) = \alpha_{0} + \alpha_{1} \mathbb{E}(e^{2}_{t-1}) \end{aligned} \]

From the first of these equations, we get \(\mathbb{E}(e_{t}) = 0\) because \(\mathbb{E}(z_{t}) = 0\). From the second equations, we get \(\mathbf{var}(e_{t}^{2}) = \mathbb{E}(e^{2}_{t}) = \alpha_{0}/(1 - \alpha_{1})\) because \(\mathbb{E}(z^{2}_{t}) = 1\) and \(\mathbb{E}(e^{2}_{t}) = \mathbb{E}\big(e^{2}_{t-1}\big)\).

1.2 Time-Varying Volatility

The ARCH model is a popular model because its variance specifications can capture commonly observed features of the time series of financial variables. It is useful for modeling volatility and especially changes in volatility over time.

Let us examine a time-series data set on returns on four major indices to examine their respective volatilities.

library(foreign)
library(ggplot2)
library(grid)
library(gridExtra)
theme_set(theme_bw())

returns <- read.dta('/Users/czar.yobero/R/Data Sets/returns.dta')

# conver to time-series objects
returns$nasdaq <- ts(returns$nasdaq, start = c(1988, 1), frequency = 12)
returns$allords <- ts(returns$allords, start = c(1988, 1), frequency = 12)
returns$ftse <- ts(returns$ftse, start = c(1988, 1), frequency = 12)
returns$nikkei <- ts(returns$nikkei, start = c(1988, 1), frequency = 12)
returns$date <- seq(as.Date('1988-01-01'), by = 'month', length.out = length(returns$nasdaq))

nasdaq.plot <- ggplot(returns, aes(y = nasdaq, x = date)) + geom_line(col = 'blue') + labs(ylab = 'return', xlab = 'Time', title = 'Nasdaq Returns')
allords.plot <- ggplot(returns, aes(y = allords, x = date)) + geom_line(col = 'blue') + labs(ylab = 'return', xlab = 'Time', title = 'Australia All Ordinaries Returns')
ftse.plot <- ggplot(returns, aes(y = ftse, x = date)) + geom_line(col = 'blue') + labs(ylab = 'return', xlab = 'Time', title = 'FTSE Returns')
nikkei.plot <- ggplot(returns, aes(y = nikkei, x = date)) + geom_line(col = 'blue') + labs(ylab = 'return', xlab = 'Time', title = 'Nikkei Returns')

grid.arrange(nasdaq.plot, allords.plot, ftse.plot, nikkei.plot, ncol = 2 )

As you can see, the values of these series change rapidly from period to period in an seemingly unpredictable manner, which we can say is volatile. Furthermore, there are periods when large changes are followed by further large changes and periods in which small changes are followed by further small changes, another evidence of volatility.

1.3 Testing, Estimating, and Forecasting

1.3.1 Testing for ARCH Effects

A Langrange multiplier (LM) test is often used to test for the presence of ARCH effects. To perform the LM test, we first estimate the mean equation, which can be a regression of the varriable on a constant like equation (1.2) or may include other variables as well. Then, we save the estimated residuals \(\hat{e}_{t}\) and obtain their squares \(\hat{e}_{t}^{2}\). To test for the first-order ARCH, we regres \(\hat{e}_{t}^{2}\) on the squared residuals lagged \(\hat{e}_{t-1}^{2}\) as follows

\[ \begin{aligned} \hat{e}_{t}^{2} = \gamma_{0} + \gamma_{1}\hat{e}_{t-1}^{2} + v_{t} &&&& (1.3) \end{aligned} \]

where \(v_{t}\) is a random error term. the null and alternative hypotheses are

\[ \begin{aligned} H_{0} : \gamma_{1} = 0 \qquad H_{1} : \gamma_{1} \neq 0 \end{aligned} \]

If there are no ARCH effects, then \(\gamma_{1} = 0\) and the fit of (1.3) will be poor, and the equation \(R^{2}\) will be low. If there exists ARCH effects, we expect the magnitude of \(\hat{e}^{2}_{t}\) to depend on its lagged values, and the \(R^{2}\) will be relatively high. The \(LM\) test test statistic is \((T - q)R^{2}\) where \(T\) is the sample size, \(q\) is the number of \(\hat{e}^{2}_{t-j}\) terms of the right-hand side of (1.3), and \(R^{2}\) is distributed as \(\chi^{2}_{(q)}\), where \(q\) is the order of lag, and \(T - q\) is the number of complete obervations. In thise case, \(q = 1\). If \((T - q)R^{2} \geq \chi^{2}_{(1 - \alpha, q)}\), then we reject the null hypothesis that \(\gamma_{1} = 0\) and conclude that ARCH effects are present.

1.2.1 An Example

For our example, let us examine the stock returns of a fictitous company called “Bring Your Day Lighting” (BYD). Let’s first load the data set and plot the returns against time.

library(kableExtra)
library(knitr)
byd <- read.dta('/Users/czar.yobero/R/Data Sets/byd.dta')
byd$r <- ts(byd$r)
byd$date <- seq.Date(as.Date('2010-01-01'), by = 'day', length.out = length(byd$r))
ggplot(byd, aes(y = r, x = date )) + geom_line(col = 'red') +
  labs(title = 'BYD Daily Stock Returns', ylab = 'return')

To perform the test for ARCH effects, we must

  1. Estimate a mean equation that, in this example, is \(r_{t} = \beta_{0} + e_{t}\)
  2. Retrieve the estimated residuals \(\hat{e}_{t}\)
  3. Estimate (1.3).
library(dynlm)
# Step 1: Estimate mean equation r = beta + error
byd.mean <- dynlm(r ~ 1, data = byd)

# Step 2: Retrieve the residuals from the former model and square them
ehatsq <- ts(resid(byd.mean)^2)

# Step 3: regress squared residuals on one-lagged squared residuals
byd.arch <- dynlm(ehatsq ~ L(ehatsq), data = ehatsq)

summary(byd.arch)
## 
## Time series regression with "ts" data:
## Start = 2, End = 500
## 
## Call:
## dynlm(formula = ehatsq ~ L(ehatsq), data = ehatsq)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -10.8023  -0.9496  -0.7055   0.3204  31.3468 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.90826    0.12440   7.301 1.14e-12 ***
## L(ehatsq)    0.35307    0.04198   8.410 4.39e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.45 on 497 degrees of freedom
## Multiple R-squared:  0.1246, Adjusted R-squared:  0.1228 
## F-statistic: 70.72 on 1 and 497 DF,  p-value: 4.387e-16

We can check for ARCH effects by using the \(\text{ArchTest}()\) function from the \(\text{FinTS}\) package. We will use a significance level of \(\alpha=0.05\) for our null hypothesis test.

library(FinTS)
byd.archTest <- ArchTest(byd$r, lags = 1, demean = TRUE)
byd.archTest
## 
##  ARCH LM-test; Null hypothesis: no ARCH effects
## 
## data:  byd$r
## Chi-squared = 62.16, df = 1, p-value = 3.22e-15

Because the \(p\)-value is \(\lt\) \(0.05\), we reject the null hypothesis and conclude the presence of ARCH(1) effects.

1.2.2 Estimating ARCH Models

We can estimate an ARCH(1) model using the \(\text{garchFit}()\) function from the \(\text{fGarch}\) package in R. Specifically, we need to estimate the variance given in equation (1.2c). ARCH models are estiamted using the maximum likelihood method.

library(fGarch)
arch.fit <- garchFit(~garch(1,0), data = byd$r, trace = F)
summary(arch.fit)
## 
## Title:
##  GARCH Modelling 
## 
## Call:
##  garchFit(formula = ~garch(1, 0), data = byd$r, trace = F) 
## 
## Mean and Variance Equation:
##  data ~ garch(1, 0)
## <environment: 0x7f857fa1bed0>
##  [data = byd$r]
## 
## Conditional Distribution:
##  norm 
## 
## Coefficient(s):
##      mu    omega   alpha1  
## 1.06394  0.64214  0.56935  
## 
## Std. Errors:
##  based on Hessian 
## 
## Error Analysis:
##         Estimate  Std. Error  t value Pr(>|t|)    
## mu       1.06394     0.03992   26.649  < 2e-16 ***
## omega    0.64214     0.06482    9.907  < 2e-16 ***
## alpha1   0.56935     0.09131    6.235 4.52e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Log Likelihood:
##  -740.7932    normalized:  -1.481586 
## 
## Description:
##  Thu Jan  4 13:52:38 2018 by user:  
## 
## 
## Standardised Residuals Tests:
##                                 Statistic p-Value    
##  Jarque-Bera Test   R    Chi^2  0.2026129 0.9036561  
##  Shapiro-Wilk Test  R    W      0.9985749 0.9614168  
##  Ljung-Box Test     R    Q(10)  12.09516  0.2787384  
##  Ljung-Box Test     R    Q(15)  17.83393  0.2714968  
##  Ljung-Box Test     R    Q(20)  20.77658  0.410386   
##  Ljung-Box Test     R^2  Q(10)  16.73949  0.08033087 
##  Ljung-Box Test     R^2  Q(15)  21.76655  0.1140729  
##  Ljung-Box Test     R^2  Q(20)  40.54622  0.004258137
##  LM Arch Test       R    TR^2   17.4092   0.1348419  
## 
## Information Criterion Statistics:
##      AIC      BIC      SIC     HQIC 
## 2.975173 3.000460 2.975101 2.985096

Based on the output, the estimated mean of the series is

\[ \begin{aligned} \hat{r}_{t} = \hat{\beta}_{0} = 1.064 \end{aligned} \]

and the estimated variance is

\[ \begin{aligned} \hat{h}_{t} = \hat{\alpha}_{0} + \hat{\alpha}_{1}\hat{e}^{2}_{t-1} = 0.642 + 0.569 \end{aligned} \]

Plotting the conditional variance looks something like

byd$ht <- arch.fit@h.t
ggplot(byd, aes(y = ht, x = date)) + geom_line(col = '#ff9933') + ylab('Conditional Variance') + xlab('Date')

Looking at the graph, you can see the periods in which volatility was high.

Exercise

Now, that you have a rudimentary understanding of the ARCH(1) model, as an exercise, try modeling the volatility of bitcoin or other cryptocurrencies and see if you can spot anything interesting.