Risk Analysis

Volatility Risk, GARCH Model with R and Value at Risk

  • We consider the log-return of an asset. That is, \[ r_t=\log(P_t)-\log(P_{t-1}), \] where \(P_t\) is the price of the asset at time \(t\).

  • The mean return \[ \bar{r}=\frac{1}{m}\sum_{t=1}^m r_{n-t}, \] where \(m\) is the number of observations leading upto present period.

  • If we assume that the mean return is zero, we obtain the maximum likelihood estimator of variance: \[ \sigma_n^2=\frac{1}{m}\sum_{t=1}^m r_{n-t}^2. \]

  • This formula to estimate volatility assign weght equally to each observations.

  • So more remote older return have the same influence on the estimated volatility, as returns that are more current.

  • As our aim is to estimate the present level of volatility, we would like to assign weight on current return more heavily than the older one.

  • The generic scheme can be presented as \[ \sigma_n^2=\sum_{t=1}^m \alpha_t r_{n-t}^2, ~~ \sum_{t=1}^m \alpha_t =1, \] \(\alpha_t\) is the weight of the return \(t\) days ago.

  • As our aim is to have a greater influence on current returns, then weight should decline in value for older returns, i.e., \(\alpha_t \geq \alpha_{t+1}\)

  • An addition to this representation is to assume a long-run variance.

  • The most commonly used model is the autoregressive conditional heteroskedasticity model, ARCH(m), which can be represented by: \[ \sigma_n^2=\gamma \sigma_L+ \sum_{t=1}^m \alpha_t r_{n-t}^2,~~ ~~ \gamma+ \sum_{t=1}^m \alpha_t = 1, \] where \(\sigma_L\) is the long-run variance weighted by parameter \(\gamma\).

  • We consider ARCH(1) model, which can be presented as \[ \sigma_n^2=\gamma \sigma_L+ \alpha_1 r_{n-1}^2,~~ ~~ \gamma+ \alpha_1 = 1, \]

  • We replace \(\alpha_1\) by \((1-\gamma)\) we have the model as \[ \sigma_n^2=\gamma \sigma_L+ (1-\gamma) r_{n-1}^2, \] where \(\gamma\) is weight on long-run variance.

  • Though ARCH model consider long-run variance, however it does not consider the previous volatility estimate.

  • This disadvantage of ARCH model can be overcome by the generalized autoregressive conditional heteroskedasticity, more popularly known as the GARCH model.

  • The GARCH(1,1) model can be presented as as \[ \sigma_n^2=\gamma \sigma_L+ \alpha_1 r_{n-1}^2 + \beta_1 \sigma_{n-1}^2,~~ ~~ \gamma+ \alpha_1 +\beta_1 = 1, \] where \(\alpha_1\) is the weight on previous period’s return, \(\beta_1\) is the weight on previous volatility estimate, \(\gamma\) is the weight on long-run variance.

  • We can write \(\omega=\gamma \sigma_L^2\), such that \[ \sigma_n^2=\omega+ \alpha_1 r_{n-1}^2 + \beta_1 \sigma_{n-1}^2, \] where \(\sigma_L^2=\frac{\omega}{1-\alpha_1-\beta_1}\) is the long-run variances.

  • If \(\omega=0\), \(\alpha=(1-\lambda)\) and \(\beta_1=\lambda\) then we have Exponentially Weighted Moving Average (EWMA) model, i.e., \[ \sigma_n^2=(1-\lambda)r_{n-1}^2 + \lambda \sigma_{n-1}^2 \] where \(\lambda\) is exponential decay rate. So EWMA is the special case of GARCH model.

Fitting GARCH Model

library(tseries)
VIX<-get.hist.quote(instrument = "^VIX"
                        ,start="2015-01-01"
                        ,end=Sys.Date()
                        ,quote="AdjClose"
                        ,provider = "yahoo")
## time series starts 2015-01-02
## time series ends   2016-12-16
par(mfrow=c(1,2))
plot(ts(VIX$AdjClose),ylab="VIX")

## Why are we using SPY and not ^GSPC?
SnP.500<-get.hist.quote(instrument = "SPY"
                        ,start="2015-01-01"
                        ,end=Sys.Date()
                        ,quote="AdjClose"
                        ,provider = "yahoo")
## time series starts 2015-01-02
## time series ends   2016-12-16
plot(ts(SnP.500$AdjClose),ylab="S&P 500")

## log-return
ln_rt<-as.numeric(diff(log(SnP.500)))

## Fit GARCH(1,1)
fit.garch<-garch(ln_rt,order=c(1,1))
## 
##  ***** ESTIMATION WITH ANALYTICAL GRADIENT ***** 
## 
## 
##      I     INITIAL X(I)        D(I)
## 
##      1     7.431534e-05     1.000e+00
##      2     5.000000e-02     1.000e+00
##      3     5.000000e-02     1.000e+00
## 
##     IT   NF      F         RELDF    PRELDF    RELDX   STPPAR   D*STEP   NPRELDF
##      0    1 -2.084e+03
##      1    7 -2.084e+03  2.45e-04  1.19e-03  1.0e-04  2.5e+10  1.0e-05  1.47e+07
##      2    8 -2.085e+03  1.52e-04  1.91e-04  9.0e-05  2.0e+00  1.0e-05  6.85e+00
##      3    9 -2.085e+03  5.39e-06  4.95e-06  9.9e-05  2.0e+00  1.0e-05  6.93e+00
##      4   17 -2.095e+03  5.02e-03  8.82e-03  5.2e-01  2.0e+00  1.1e-01  6.91e+00
##      5   20 -2.107e+03  5.79e-03  6.06e-03  7.3e-01  1.8e+00  3.2e-01  3.33e-01
##      6   33 -2.108e+03  8.45e-05  1.96e-03  1.1e-05  2.3e+00  8.5e-06  6.51e-01
##      7   34 -2.109e+03  5.90e-04  4.50e-04  4.6e-06  2.0e+00  4.2e-06  3.72e-01
##      8   35 -2.109e+03  3.25e-05  4.29e-05  5.3e-06  2.0e+00  4.2e-06  5.21e-01
##      9   36 -2.109e+03  2.49e-06  2.71e-06  5.6e-06  2.0e+00  4.2e-06  4.88e-01
##     10   45 -2.116e+03  3.49e-03  5.66e-03  2.7e-01  2.0e+00  2.8e-01  4.88e-01
##     11   57 -2.118e+03  7.06e-04  2.72e-03  2.8e-06  2.8e+00  3.7e-06  3.55e-03
##     12   58 -2.118e+03  2.80e-04  2.07e-04  2.6e-06  2.0e+00  3.7e-06  3.25e-04
##     13   59 -2.118e+03  2.60e-05  3.19e-05  2.4e-06  2.0e+00  3.7e-06  9.00e-05
##     14   60 -2.118e+03  1.13e-06  1.21e-06  2.1e-06  2.0e+00  3.7e-06  4.99e-05
##     15   68 -2.119e+03  8.69e-05  4.90e-05  9.3e-03  0.0e+00  1.7e-02  4.90e-05
##     16   69 -2.119e+03  9.26e-05  1.26e-04  3.3e-02  0.0e+00  6.0e-02  1.26e-04
##     17   70 -2.119e+03  1.32e-06  1.90e-05  1.4e-02  0.0e+00  2.2e-02  1.90e-05
##     18   71 -2.119e+03  1.50e-06  1.58e-05  5.5e-03  0.0e+00  9.2e-03  1.58e-05
##     19   72 -2.119e+03  7.74e-06  7.79e-06  2.9e-03  7.1e-01  4.6e-03  9.95e-06
##     20   73 -2.119e+03  1.89e-06  2.03e-06  2.3e-03  0.0e+00  3.3e-03  2.03e-06
##     21   74 -2.119e+03  1.36e-08  1.14e-08  1.7e-04  0.0e+00  2.4e-04  1.14e-08
##     22   75 -2.119e+03 -3.25e-10  8.83e-11  5.8e-06  0.0e+00  8.2e-06  8.83e-11
## 
##  ***** RELATIVE FUNCTION CONVERGENCE *****
## 
##  FUNCTION    -2.118826e+03   RELDX        5.809e-06
##  FUNC. EVALS      75         GRAD. EVALS      22
##  PRELDF       8.825e-11      NPRELDF      8.825e-11
## 
##      I      FINAL X(I)        D(I)          G(I)
## 
##      1    8.782813e-06     1.000e+00    -4.645e+02
##      2    1.913741e-01     1.000e+00    -4.872e-03
##      3    6.973940e-01     1.000e+00    -1.982e-02
## Fitted Values
sigma_hat<-fit.garch$fitted.values[,"sigt"]

vol_hat<-sigma_hat*sqrt(252)*100

## plot GARCH fitted volatility vs. VIX

par(mfrow=c(1,1))
plot(ts(VIX$AdjClose),col="red",lwd=2,ylim=c(0,40))
lines(vol_hat,col="blue",lwd=2)
text<-c("VIX","GARCH(1,1)")
legend(350,40,text,col=c("red","blue"),lwd=c(2,2))