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=0
https://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 VaR and 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.
