This R Markdown document is intended to show my R skills and interest for financial econometrics in an applied and tangible manner.

Value-at-risk (VaR) with all its associated problems remains the workhorse model in risk management. One key problem of VaR is that it does not properly account for volatility clustering, which means that VaR limits are breached in serial dependence across time. As a result, risk is underestimated during a crisis.

A powerful approach to solve this problem is to combine VaR with GARCH models, which take conditional volatility into account. In order to illustrate this method, I apply a GARCH(1,1) with a normal distribution to the Swiss equity market index SMI.

1. Getting Data from Yahoo and Plots

##Initialisation
#Load Packages
library(fImport)
library(fPortfolio)
library(ggplot2)

#Inputs
from = "1995-11-20"
to = "2015-12-17"
symbol = "^SSMI"

#Get Data from Yahoo
TS <- yahooSeries(symbol, from = from, to = to)
SMI <- TS[,ncol(TS)]
SMI <- returns(SMI, method = "continuous")

#Plot SMI Returns
seriesPlot(TS[,4])

histPlot(SMI, main = "SMI Returns")



2. Some Math

In standard textbooks VaR is mathematically defined as

\[VaR(X) = \mu + \sigma N^{-1}(X)\]

where \(\mu\) is the mean of daily returns, \(\sigma\) is the standard deviation and \(N^{-1}(.)\) is the inverse of the cumulative normal distribution at a given confidence level \(X\). Unfortunately, it is not clear how the standard deviation \(\sigma\) is defined. In the case of the VaR formula, it is the unconditional volatility (some sort of long term average) and usually the sample standard deviation.

Since the mean is approximately zero, I will use the formula below to calculate the volatility:

\[\sigma_n^2 = \frac{1}{n} \sum_{i=1}^n R_i^2\]

\(n\) is the number of observations and I will use both, the full sample and a 150-day moving average. This formula ignores “volatility clustering”. By doing so, risk is underestimated/overestimated during periods with high/low volatility periods.

In order to incorporate the conditional volatility into the VaR formula, I use a GARCH(1,1), which is defined as \[\sigma_t^2 = \omega + \alpha R_{t-1}^2 + \beta \sigma_{t-1}^2\]



3. Model Estimation

The SMI returns data has 5078 observations. I use the first 3078 observations to make the initial estimation for the GARCH model. The remaining 2000 observations are used for validation and testing.

library(rugarch)
library(zoo)

SMIdf <- as.data.frame(SMI)

#Historical Volatility (Moving Average)
SMIdf <- as.data.frame(SMI)
hv <- rollapply(SMIdf, 150, sd)
l <- length(hv)
hv <- hv[(l-2000+1):l]


#GARCH
#GARCH Spec - (Change Distribution here)
gspec11 <- ugarchspec(variance.model = list(model = "sGARCH", 
                      garchOrder = c(1, 1)),
                      mean.model=list(armaOrder=c(0,0), 
                      include.mean = FALSE), 
                      distribution="norm")

#Rolling Estimation
roll11 <- ugarchroll(gspec11, SMIdf, n.start=3048,
                     refit.every = 25, refit.window = "moving",
                     VaR.alpha = c(0.025, 0.05))



4. Results

Returns <- SMI[3048:5047]  #Returns for Validation


#VaR 5%
VaRStatic <- sd(Returns) * qnorm(0.05)   #Static, unconditional
VaRMA <- hv * qnorm(0.05)                       #Moving Average
VaRGARCH <- roll11@forecast$VaR[,2]             #GARCH Extraxtion

#xaxis
xaxis <- rownames(SMI)
xaxis <- xaxis[3048:5047]
xaxis <- c(xaxis[1], xaxis[500], xaxis[1000], xaxis[1500], xaxis[2000])

#VaR Plot
plot(Returns, type = "l", pch = 16, cex = 0.8,  col = gray(0.2, 0.5),
     ylab = "Returns", main = "95% VaR Forecasting", xaxt = "n")
axis(1, at=c(1, 500, 1000, 1500, 2000), labels=xaxis)
lines(VaRGARCH, col = 1)
lines(VaRMA, col = 4)
abline(h=VaRStatic, col = 2)
legend('topright', c("GARCH(1,1)", "MA 150 days", "Static VaR") , 
       lty=1, col=c(1,4,2), bty='n', cex=.75)

As we can see from the plot, the VaR-GARCH (black line) combination is way more realistic and lowers the VAR limit when volatility clustering occurs, whereas for the static VaR (red line) we observe serial limit breaches. 

#Volatility Plot
plot(abs(Returns), type = "l", col = grey(0.4, 0.5), 
     ylab = "Absolute Returns", main = "Volatility Forecasting", xaxt = "n")
axis(1, at=c(1, 500, 1000, 1500, 2000), labels=xaxis)
Sigma11 <- roll11@forecast$density$Sigma
lines(Sigma11, col = 1)
lines(hv, col = 4)
abline(h=sd(Returns), col = 2)
legend('topright', c("GARCH(1,1)", "MA 150 days", "Unconditional") , 
       lty=1, col=c(1,4, 2), bty='n', cex=.75)