I estimated a GARCH fit to the log returns of three series (CAC 40, a french
real estate index and french T10 bond yield series) using `rugarch`. I then
manually calculated and backtested the VaR and CVaR measures. I also fitted
a DCC-GARCH(1,1) to the log returns of the 3 series using `rmgarch` and now
I would like to backtest the VaR and CVaR measures in a similar way as I did
for the univariate GARCH cases.
We'll need to specify the following functions for the CVaR before we
proceed:
#This function calculates the CVaR at a certain position gdist list
cvar <- function(p=0.05, s = "CAC", dist_params = gdist_var, pos = l, v
= df, dist = "jsu"){
ES <- abs((integrate(qdist, lower = 0, upper = p, distribution = dist,
mu = gdist_var[[s]][, 'Mu'][pos], sigma = gdist_var[[s]][, 'Sigma'][pos],
shape = gdist_var[[s]][, 'Shape'][pos], skew =
gdist_var[[s]][, 'Skew'][pos])$value)/p * v[nrow(v),s])
return(ES)
}
#This function calculates the CVaR given the arguments
cvar_df <- function(p=0.01, dist = "jsu", mu = Mu, sigma = Sigma, shape
= Shape, skew = Skew){
ES <- (integrate(qdist, lower = 0, upper = p, distribution = dist, mu
= mu, sigma = sigma, shape = shape, skew = skew)$value)/p
return(ES)
}
#This function is a vectorized form of the above
vcvar_df <- Vectorize(cvar_df)
The data can be found on dropbox under the following links (one for french
real estate index data and the other for french bonds) the CAC 40 data is
downloadable with `quantmod`:
https://www.dropbox.com/s/vy8sl88fs5opmi3/IEIF%20SIIC%20FRANCE_quote_chart.csv?dl=0https://www.dropbox.com/s/xljxk5izy6pt1ds/entre_obligations.csv?dl=0
The commented code is the following:
#loading libraries and data:
require(tidyquant)
require(reshape2)
require(astsa)
require(GGally)
require(forecast)
source("functions.R", local = T)
#
#
https://www.banque-france.fr/statistiques/taux-et-cours/les-indices-obligataires
obli_10 <- read.csv("entre_obligations.csv", sep = ";", na.strings =
"-", stringsAsFactors = F) %>%
rename(Date = 1) %>%
mutate(Date = dmy(Date)) %>%
mutate_at(vars(-Date), funs(gsub("\\,", ".", .))) %>%
mutate_at(vars(-Date), funs(as.numeric)) %>%
dplyr::select(c(1,2))
# #https://live.euronext.com/en/product/indices/QS0010980447-XPAR/quotes
indices nu de:
#
#
https://www.ieif.fr/wp-content/plugins/aa-indices/datas/histo/index.php?IndiceNu=SIICNu&IndiceNet=SIICNet&IndiceBrut=SIICBrut&Indice=Euronext%20IEIF%20SIIC%20France
reit <- read.csv("IEIF SIIC FRANCE_quote_chart.csv", stringsAsFactors =
F) %>%
dplyr::select(1,2) %>%
rename(Date = 1) %>%
mutate(Date = substr(Date, 1, 10)) %>%
mutate(Date = ymd(Date))
cac <- as.data.frame(Ad(getSymbols("^FCHI", src = "yahoo", adjust = T,
auto.assign = FALSE)))
cac <- cac %>%
mutate(Date = rownames(.)) %>%
mutate(Date = ymd(Date)) %>%
dplyr::select(Date, everything())
#Calculate the log returns
lr_df <- as.data.frame(sapply(df[2:ncol(df)], function(x) diff(log(x))))
lr_df <-cbind(df$Date[2:nrow(df)], lr_df) %>%
dplyr::rename(Date = !!names(.)[1])
#Specification of GARCH models
cac_egarch_spec <- ugarchspec(mean.model = list(armaOrder = c(3, 3),
include.mean = T, archm = F, archpow = 1),
variance.model = list(model = "eGARCH",
garchOrder = c(2, 1)),
distribution.model="jsu")
reit_egarch_spec <- ugarchspec(mean.model = list(armaOrder = c(3, 1),
include.mean = T, archm = F, archpow = 1),
variance.model = list(model = "eGARCH",
garchOrder = c(2, 1)),
distribution.model="nig")
obli_apgarch_spec <- ugarchspec(mean.model = list(armaOrder = c(2, 1),
include.mean = T, archm = F, archpow = 1),
variance.model = list(model = "apARCH",
garchOrder = c(1, 1)),
distribution.model="jsu")
#Get VaR and CVaR
cac_roll <- ugarchroll(cac_egarch_spec, lr_df[,2],n.start = 750,
refit.every = 50, refit.window = "moving",
solver = "hybrid", calculate.VaR =
TRUE, VaR.alpha = c(0.01, 0.025, 0.05), keep.coef = T,
fit.control = list(scale = 1))
reit_roll <- ugarchroll(reit_egarch_spec, lr_df[,3],n.start = 750,
refit.every = 50, refit.window = "moving",
solver = "hybrid", calculate.VaR = TRUE,
VaR.alpha = c(0.01, 0.025, 0.05), keep.coef = T,
fit.control = list(scale = 1))
obli_roll <- ugarchroll(obli_apgarch_spec, lr_df[,4],n.start = 750,
refit.every = 50, refit.window = "moving",
solver = "hybrid", calculate.VaR = TRUE,
VaR.alpha = c(0.01, 0.025, 0.05), keep.coef = T,
fit.control = list(scale = 1))
gdist_var <- list()
gdist_var[["CAC"]] <- as.data.frame(cac_roll, which = 'density')
gdist_var[["REIT"]] <- as.data.frame(reit_roll, which = 'density')
gdist_var[["OBLI_10"]] <- as.data.frame(obli_roll, which = 'density')
#VaR and CVaR calculations
p <- c(0.05, 0.025, 0.01)
l <- nrow(gdist_var[["CAC"]])
for(j in p){
for(i in 1:3){
print(paste("VaR", names(gdist_var)[i], 1-j))
print(abs(qdist(dg[[i]], p=j, mu=gdist_var[[i]]$Mu[l],
sigma=gdist_var[[i]]$Sigma[l], skew=gdist_var[[i]]$Skew[l],
shape=gdist_var[[i]]$Shape[l]))*df[nrow(df),i+1])
}
}
for(j in p){
for(i in 1:3){
print(paste("CVaR", names(gdist_var)[i], 1-j))
print(cvar(p = j, s = names(gdist_var)[i], dist_params = gdist_var,
pos = l, v = df, dist = dg[[i]]))
}
}
#VaR plots for cac only but will be done for others the same way
var_cac <- gdist_var$CAC
var_cac <- cbind.data.frame(tail(lr_df[,c("Date","CAC")],2438), var_cac)
%>%
dplyr::select(-`Shape(GIG)`, -Realized) %>%
dplyr::mutate(VaR_99 = qdist("jsu", p = 0.01, mu = Mu, sigma = Sigma,
skew = Skew, shape = Shape)) %>%
dplyr::select(-Mu, -Sigma, -Skew, -Shape)
var_cac <- melt(var_cac, id.vars = "Date")
ggplot(data = var_cac, aes(x = Date, value)) + geom_line(aes(colour =
variable)) +
ggtitle("Series with 1% 1D VaR Limit") +
theme(plot.title = element_text(hjust = 0.5))
#VaR backtesting reports using report function
report(cac_roll, type = "VaR", VaR.alpha = 0.05, conf.level = 0.95)
report(reit_roll, type = "VaR", VaR.alpha = 0.05, conf.level = 0.95)
report(obli_roll, type = "VaR", VaR.alpha = 0.05, conf.level = 0.95)
#CVaR plots for CAC only but will be done for others
cvar_cac <- gdist_var$CAC
cvar_cac <- cbind.data.frame(tail(lr_df[,c("Date","CAC")],2438),
cvar_cac) %>%
dplyr::select(-`Shape(GIG)`, -Realized) %>%
dplyr::mutate(CVaR_99 = vcvar_df(p = 0.01, dist = "jsu", mu = Mu,
sigma = Sigma, shape = Shape, skew = Skew)) %>%
dplyr::select(-Mu, -Sigma, -Skew, -Shape)
mcvar_cac <- melt(cvar_cac, id.vars = "Date")
ggplot(data = mcvar_cac, aes(x = Date, value)) + geom_line(aes(colour =
variable)) +
ggtitle("Series with 1% 1D CVaR Limit") +
theme(plot.title = element_text(hjust = 0.5))
#Bactesting CVaRby calculating nuber of times CVaR crossed
cvar_cac <- gdist_var$CAC
cvar_cac <- cbind.data.frame(tail(lr_df[,c("Date","CAC")],2438),
cvar_cac) %>%
dplyr::select(-`Shape(GIG)`, -Realized) %>%
dplyr::mutate(CVaR_99 = vcvar_df(p = 0.01, dist = "jsu", mu = Mu,
sigma = Sigma, shape = Shape, skew = Skew)) %>%
dplyr::mutate(CVaR_975 = vcvar_df(p = 0.025, dist = "jsu", mu = Mu,
sigma = Sigma, shape = Shape, skew = Skew)) %>%
dplyr::mutate(CVaR_95 = vcvar_df(p = 0.05, dist = "jsu", mu = Mu,
sigma = Sigma, shape = Shape, skew = Skew)) %>%
mutate(depasse_99 = case_when(CVaR_99 >= .[[2]] ~ 1, TRUE ~ 0)) %>%
mutate(depasse_975 = case_when(CVaR_975 >= .[[2]] ~ 1, TRUE ~ 0)) %>%
mutate(depasse_95 = case_when(CVaR_95 >= .[[2]] ~ 1, TRUE ~ 0)) %>%
mutate(sum_99 = sum(depasse_99)) %>%
mutate(sum_975 = sum(depasse_975)) %>%
mutate(sum_95 = sum(depasse_95))
#DCC GARCH of GARCH models above:
require(rmgarch)
dcc_garch <- multispec(c(cac_egarch_spec, reit_egarch_spec,
obli_apgarch_spec))
dcc_multfit <- multifit(dcc_garch, lr_df[,2:ncol(lr_df)]) #fitting many
univariate models
dcc_spec <- dccspec(uspec = dcc_garch, dccOrder = c(1,1), distribution =
"mvnorm")
dcc_fit <- dccfit(dcc_spec, lr_df[,2:ncol(lr_df)], fit.control =
list(eval.se = TRUE), fit = dcc_multfit) #fit = dcc_multfit not really
necessary but more robust
dcc_roll <- dccroll(dcc_spec, lr_df[,2:4],n.start = 750, refit.every =
50, refit.window = "moving",
solver = "solnp", calculate.VaR = TRUE, VaR.alpha
= c(0.01, 0.025, 0.05), keep.coef = T,
fit.control = list(scale = 1))
Now I want to do the backtesting and plotting steps for both 1Day VaR and
1Day CVaR measures. Ideally I would also conduct the Kupiec and
Christoffersen test just like in the function `report` of the package
`rugarch`. I am realy stumped as I tried to find an answer online but
couldn't.
Eliot TABET | Structurer/Quantitative Analyst
Les informations contenues dans ce message sont confidentielles.
Si vous receviez ce message par erreur, merci de pr?venir l'exp?diteur
et de le supprimer ainsi que ses pi?ces jointes.
N'imprimez ce mail que si c'est n?cessaire
[[alternative HTML version deleted]]