8a: GARCH model

1. Stylized Facts of Financial Return

1a. Financial (log) Return

\[\begin{align*} Y_t &= \mbox{ Stock Price (observation) } \\ \\ X_t &= \ln(Y_t) - \ln(Y_{t-1}) \hspace{5mm} : \mbox{ log-return } \end{align*}\]

1b. Ex: Apple

  library(quantmod)

## Loading required package: xts

## Loading required package: zoo

## 
## Attaching package: 'zoo'

## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric

## Loading required package: TTR

## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo

  source('https://nmimoto.github.io/R/TS-00.txt')

  getSymbols("AAPL")    #- download from Yahoo!

## [1] "AAPL"

  Apple = Ad(AAPL)['2017::2019']

  plot(Apple)

  plot(log(Apple))

  plot(diff(log(Apple)))

  hist(diff(log(Apple)), 20)

1c. Stylized Facts

Not normal - Heavy Tailed unconditional and conditional distribution

Uncorrelated

Squares are correlated

Clustering

Asymmetry

  Randomness.tests(as.numeric(diff(log(Apple))))

##   B-L test H0: the series is uncorrelated
##   M-L test H0: the square of the series is uncorrelated
##   J-B test H0: the series came from Normal distribution
##   SD         : Standard Deviation of the series

##       BL15 BL20  BL25 ML15 ML20 JB    SD
## [1,] 0.189 0.15 0.226    0    0  0 0.016

1d. Ex: SP500

We need a model to capture these characteristics.

  library(quantmod)
  getSymbols("^GSPC")    #- SP500 download from Yahoo!

## [1] "GSPC"

  SP500 = Ad(GSPC)['2000::']

  plot(SP500, type="l")

  plot(log(SP500), type="l")

  plot(diff(log(SP500)), type="l")

  Randomness.tests(as.numeric(diff(log(SP500))))

##   B-L test H0: the series is uncorrelated
##   M-L test H0: the square of the series is uncorrelated
##   J-B test H0: the series came from Normal distribution
##   SD         : Standard Deviation of the series

##      BL15 BL20 BL25 ML15 ML20 JB    SD
## [1,]    0    0    0    0    0  0 0.012

2. ARCH model

Engle (1985) AutoRegressive Conditionally Heteroscedastic Model

Won Nobel Prize in Economics

\[\begin{align*} Y_t &= \sigma_t \epsilon_t \hspace{10mm} \epsilon_t \sim_{iid} N(0,1) \\\\ \sigma_t^2 &= \omega+ \alpha Y_{t-1}^2 \end{align*}\]

(unconditional) Mean \[ E(Y_t) = \sigma_t E(\epsilon_t) = 0 \]

(unconditional) Variance \[ V(Y_t) = \frac{\omega}{1-\alpha} \hspace{5mm} 0<\alpha<1 \]

2a. Cond’l Mean and Var of ARCH

Conditional mean \[ E[Y_t \Big| Y_{t-1}, \epsilon_{t-1}, \ldots] \hspace{3mm} = \hspace{3mm} \sigma_t E[\epsilon_t] = 0 \]

Conditional variance \[ V[Y_t \Big| Y_{t-1}, \epsilon_{t-1}, \ldots] \hspace{3mm} = \hspace{3mm} \sigma_t^2 V[\epsilon_t] = \sigma_t^2 \]

2b. ARCH is uncorrelated

\(Y_t\) is uncorrelated, so it will pass the Ljung-Box test.

But \(Y_t^2\) is correlated, so it will NOT pass the McLeod-Li test.

2c. Cond’l Mean and Var of ARMA

Suppose \(Y_t\) is AR(1) series observation

(unconditional) Mean \[ E( Y_t ) = 0 \]

(unconditional) Variance \[ V( Y_t ) = \gamma(0) = (1+\phi_1^2) \sigma^2 \]

2d. Conditional Mean

Conditonal Mean: \[ E\Big(Y_t \hspace{2mm} \Big| \hspace{2mm} \mbox{all variables realized by yesterday. }\Big) \]

Conditonal mean of AR(1): \(Y_t = \phi Y_{t-1} + \epsilon_t\) \[\begin{align*} E\Big(Y_t \hspace{2mm} \Big| \hspace{2mm} Y_{t-1}, Y_{t-2}, \ldots, \epsilon_{t-1}, \epsilon_{t-2}, \ldots\Big) &= E \Big(\phi Y_{t-1} + \epsilon_t \hspace{2mm} \Big| \hspace{2mm} Y_{t-1}, Y_{t-2}, \ldots, \epsilon_{t-1}, \epsilon_{t-2}, \ldots\Big)\\ \\ &= \phi Y_{t-1} + E(\epsilon_t) \hspace{2mm} = \hspace{2mm} \phi Y_{t-1}. \end{align*}\]

Conditonal variance \[\begin{align*} V\Big(Y_t \hspace{2mm} \Big| \hspace{2mm} Y_{t-1}, Y_{t-2}, \ldots, \epsilon_{t-1}, \epsilon_{t-2}, \ldots\Big) &= V\Big(\phi Y_{t-1} + \epsilon_t \hspace{2mm} \Big| \hspace{2mm} Y_{t-1}, Y_{t-2}, \ldots, \epsilon_{t-1}, \epsilon_{t-2}, \ldots\Big)\\ \\ &= V(\epsilon_t) \hspace{2mm} = \hspace{2mm} \sigma^2 \end{align*}\]

    n=200
    Y = arima.sim(n = n, list(ar = c(0.8) ))

    v.AR1 = (1+.8^2)*1

    plot(Y, type="o", main="AR(1)", xlim=c(0,n*1.1))
    abline(h=0)
    abline(h=c(1,-1)*1.96*v.AR1, col="red")

    lines(n+1, .8*Y[n], type="p", col="blue")
    lines(c(n+1,n+1), .8*Y[n]+c(1.96,-1.96), type="p", col="red")

Suppose \(Y_t\) is MA(1) series observation

(unconditional) Mean \[ E( Y_t ) = 0 \]

(unconditional) Variance \[ V( Y_t ) = \gamma(0) = (1+\theta_1^2) \sigma^2 \]

Conditonal mean of MA(1): \(Y_t = \epsilon_t + \theta_1 \epsilon_{t-1}\) \[\begin{align*} E\Big(Y_t \hspace{2mm} \Big| \hspace{2mm} Y_{t-1}, Y_{t-2}, \ldots, \epsilon_{t-1}, \epsilon_{t-2}, \ldots\Big) &= E \Big( \epsilon_t + \theta_1 \epsilon_{t-1} \hspace{2mm} \Big| \hspace{2mm} Y_{t-1}, Y_{t-2}, \ldots, \epsilon_{t-1}, \epsilon_{t-2}, \ldots\Big) \\\\ &= E ( \epsilon_t )+ \theta_1 \epsilon_{t-1} \hspace{2mm} = \hspace{2mm} \theta_1 \epsilon_{t-1} \end{align*}\]

Note that \(\epsilon_{t-1}\) is not observable.

Conditonal variance of MA(1): \(Y_t = \epsilon_t + \theta_1 \epsilon_{t-1}\) \[\begin{align*} Var\Big(Y_t \hspace{2mm} \Big| \hspace{2mm} Y_{t-1}, Y_{t-2}, \ldots, \epsilon_{t-1}, \epsilon_{t-2}, \ldots\Big) \ &= Var\Big( \epsilon_t + \theta_1 \epsilon_{t-1} \hspace{2mm} \Big| \hspace{2mm} Y_{t-1}, Y_{t-2}, \ldots, \epsilon_{t-1}, \epsilon_{t-2}, \ldots\Big)\\ &= Var( \epsilon_t ) \hspace{2mm} = \hspace{2mm} \sigma^2. \end{align*}\]ilon_t ) = ^2. \end{align*}

AR(1) \[\begin{align*} \mbox{ Uncond'l } \hspace{20mm} \mbox{ cond'l } \\ E(Y_t) &= 0, \hspace{28mm} E(Y_t|\omega_{t-1}) = \phi_1 Y_{t-1}, \\ Var(Y_t) &= (1+\phi_1^2)\sigma^2 \hspace{10mm} Var(Y_t|\omega_{t-1}) = \sigma^2 \end{align*}\]

MA(1) \[\begin{align*} \mbox{ Uncond'l } \hspace{20mm} \mbox{ cond'l } \\ E(Y_t) &= 0, \hspace{28mm} E(Y_t|\omega_{t-1}) = \theta_1 \epsilon_{t-1}, \\ Var(Y_t) &= (1+\theta_1^2)\sigma^2 \hspace{10mm} Var(Y_t|\omega_{t-1}) = \sigma^2 \end{align*}\]

For ARMA(p,q) model, conditional mean changes, but conditional variance is constant.

  x  = arima.sim(n = 50, list(ar = c(0.8) ))
  plot(x, type="o", main="AR(1)")
  abline(h=0)

  x  = arima.sim(n = 50, list(ma = c(0.5) ))
  plot(x, type="o", main="MA(1)")
  abline(h=0)

2e. Heteroschedasticity

Don’t confuse the conditional heteroscedasticity with (unconditonal) heteroscedasticity:

  D = read.csv("https://nmimoto.github.io/datasets/gas.csv")
  D1 = ts(D[,2], start=c(1956, 1), freq=12)
  plot(D1, type='o')

3. GARCH model

GARCH(1,1) model

\[\begin{align*} Y_t &= \sigma_t \epsilon_t \hspace{10mm} \epsilon_t \sim_{iid} N(0,1) \\\\ \sigma_t^2 &= \omega + \alpha Y_{t-1}^2 + \beta \sigma^2_{t-1} \end{align*}\]

3a. Conditional Mean and Var of GARCH

Conditional Mean: \(0\)

Conditional Variance : \(\sigma_t^2\)

3b. Ex: Daily SPY

Daily Price of SP500 ETF (SPY) from Jan 02 2000 to Dec 31 2014

  library(quantmod)
  source('https://nmimoto.github.io/R/TS-00.txt')

  getSymbols("SPY")    #- SP500 download from Yahoo!

## [1] "SPY"

  SPY = Ad(SPY)['2010::']

  is.ts(SPY)      # not ts object

## [1] FALSE

  is.xts(SPY)     # its xts object

## [1] TRUE

  plot( SPY )

  plot( log(SPY) )

  plot( diff( log(SPY) ) )

  Randomness.tests( diff(log(SPY)) )

##   B-L test H0: the series is uncorrelated
##   M-L test H0: the square of the series is uncorrelated
##   J-B test H0: the series came from Normal distribution
##   SD         : Standard Deviation of the series

##      BL15 BL20 BL25 ML15 ML20 JB    SD
## [1,]    0    0    0    0    0  0 0.011

  # Fit GARCH model
  library(fGarch)

## NOTE: Packages 'fBasics', 'timeDate', and 'timeSeries' are no longer
## attached to the search() path when 'fGarch' is attached.
## 
## If needed attach them yourself in your R script by e.g.,
##         require("timeSeries")

## 
## Attaching package: 'fGarch'

## The following object is masked from 'package:TTR':
## 
##     volatility

  Y = diff( log(SPY) )[-1]     # remove the first diff for NA

    Fit01 =  garchFit(~ garch(1,1), data=Y, cond.dist="norm", include.mean = FALSE, trace = FALSE)
    Fit01

## 
## Title:
##  GARCH Modelling 
## 
## Call:
##  garchFit(formula = ~garch(1, 1), data = Y, cond.dist = "norm", 
##     include.mean = FALSE, trace = FALSE) 
## 
## Mean and Variance Equation:
##  data ~ garch(1, 1)
## <environment: 0x10e015468>
##  [data = Y]
## 
## Conditional Distribution:
##  norm 
## 
## Coefficient(s):
##      omega      alpha1       beta1  
## 3.7798e-06  1.5098e-01  8.1387e-01  
## 
## Std. Errors:
##  based on Hessian 
## 
## Error Analysis:
##         Estimate  Std. Error  t value Pr(>|t|)    
## omega  3.780e-06   4.755e-07    7.949 1.78e-15 ***
## alpha1 1.510e-01   1.290e-02   11.703  < 2e-16 ***
## beta1  8.139e-01   1.396e-02   58.291  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Log Likelihood:
##  13917.88    normalized:  3.316149 
## 
## Description:
##  Sat Sep 12 23:10:40 2026 by user:

    ##  Fit01@fit$par                  # estimated parameters
    ##  Fit01@residuals                # this is not GARCH residuals! This is same as Y.
    ##  Fit01@sigma.t                  # estimated sig_t
    ##  Fit01@residuals/Fit1@sigma.t   # this is the (standardized) GARCH residuals

    Fit01@fit$ics                   # AIC and BIC are here

##       AIC       BIC       SIC      HQIC 
## -6.630869 -6.626336 -6.630870 -6.629266

    plot(Y)

    plot(Fit01@sigma.t, type="l")   # estimated Sigma_t

    sigma.upper = xts(1.96*Fit01@sigma.t, order.by=index(Y))
    sigma.lower = xts(-1.96*Fit01@sigma.t, order.by=index(Y))

    plot( cbind(Y, sigma.upper, sigma.lower),
        col=c("black", "red", "red"), lwd=c(2,1,1) )

    plot(  as.numeric(          Y["2018::"]), type="h", lwd=2)
    lines( as.numeric(sigma.upper["2018::"]), col="red")
    lines( as.numeric(sigma.lower["2018::"]), col="red")