[R-SIG-Finance] rugarch + VineCopula for value at risk

Ole Bueker ole.bueker at outlook.com
Sun Aug 24 04:47:03 CEST 2014


Okay, I have now successfully run the simulations and have obtained + saved the mean fitted values, the mean sigmas as well as the VaR for the means of the individual 112 companies.
However, I am not quite sure how I can calculcate the VaR of a equally-weighted portfolio consisting of the 112 companies.

Here is my current approach (modified to be working with the MWE posted before):
1. Calculate the mean fitted and mean sigma for each forecasted day for all 112 companies.rvine_fitted_mean = rvine_sigma_mean = array(NA, dim=c(10,20))
for(i in 1:10){for(k in 1:20) {    rvine_fitted_mean[i,k] <- mean(rvine_fitted[i,,k])rvine_sigma_mean[i,k] <- mean(rvine_sigma[i,,k])}}remove(i,k)
2. Calculate portfolio sigma, fitted and VaR.
rvine_portfolio_VaR01 = rvine_portfolio_VaR05 = array(NA, dim=c(10))rvine_portfolio_fitted = rvine_portfolio_sigma = array(NA, dim=c(10))sigmas <- array(NA, dim=c(10,20,20))w = 1/112
for(i in 1:10){for(j in 1:20){for(l in 1:20){sigmas[i,j,l] <- w*w*rvine_sigma_mean[i,j]*rvine_sigma_mean[i,l]*cor(rvine_fitted_mean[,j], rvine_fitted_mean[,l])}}rvine_portfolio_sigma[i] <- sum(sigmas[i,,])rvine_portfolio_fitted[i] <- w*sum(rvine_fitted_mean[i,])rvine_portfolio_VaR01[i] <- rvine_portfolio_fitted[i] + rvine_portfolio_sigma[i] * qdist('ged', 0.01, mu=0, sigma=1)rvine_portfolio_VaR05[i] <- rvine_portfolio_fitted[i] + rvine_portfolio_sigma[i] * qdist('ged', 0.05, mu=0, sigma=1)}remove(w, i, j, l, k)
The rvine_portfolio_fitted values seem to be correct since the formula is quite simple (weight * returns of individual company), but for the portfolio sigma the formula is more complex and I'm not sure if this is the correct approach.
I wanted to use the formula given at http://en.wikipedia.org/wiki/Modern_portfolio_theory#Mathematical_model : sigma = \sum_i \sum_i (w * w * sigma_i * sigma_j * \rho_{ij} ), but the resulting sigma seems to be very low.

Best Regards,
Ole

> From: ole.bueker at outlook.com
> To: alexios at 4dscape.com; r-sig-finance at r-project.org
> Date: Fri, 22 Aug 2014 06:34:10 +0200
> Subject: Re: [R-SIG-Finance] rugarch + VineCopula for value at risk
> 
> Hey Alexios,thanks a lot for the corrected version, seems to be working faster than before!'returns_crisis' are indeed the same as 'returns', forgot to replace the name in that line - should be fixed now in the dropbox file.I'll remember to post the required packages next time! The qged and pged functions are from the fGarch package I believe.Good to see that my approach seems to be correct - although in practice I would increase the moving window to 500 days (which isn't an option for me unfortunately). I don't think there is anything special with the VineCopula simulation method, not sure if the authors/maintainers of the package follow this list. Guess I can't fix that some of the values are positive, however the mean VaR should be negative - therefore I will increase the number of simulations and then just take the average for all the simulations to obtain 25 1-day ahead VaR forecasts for all 112 companies. (so VaR01 is a 25 x 112 matrix of the mean VaR forecasts). Thank!
>  s again for the help!Best Regards,Ole
> > -----Original Message-----
> > From: alexios ghalanos [mailto:alexios at 4dscape.com] 
> > Sent: Friday, August 22, 2014 11:40 AM
> > To: Ole Bueker; r-sig-finance at r-project.org
> > Subject: Re: [R-SIG-Finance] rugarch + VineCopula for value at risk
> > 
> > Ole,
> > 
> > See a 'corrected' extract of the relevant code section below. You did not define 'returns_crisis' but I assume this is the same as 'returns'.
> > You also did not tell us the required packages (please do so in future).
> > 
> > There is absolutely nothing wrong (at least that I can see with the approach and results).
> > 1. You are filtering the returns by an AR-GARCH model.
> > 2. You estimate the vine copula on the uniform margins based on the standardized residuals transformed by the conditional distribution from (1).
> > 3. You simulate from this unconditional distribution 100x20 uniform variates which you transform back to standardized values using the shape parameter from (1).
> > 3. You then simulate the T+1 conditional distribution (with 100 samples for each asset) using the GARCH model from (1).
> > 4. You obtain the T+1 distribution of the quantiles for each asset.
> > 
> > Are you very surprised that a particular point in time distribution of the 5% quantile of certain assets is not entirely in the negative domain? Think about what it is you are doing, the assumptions, high degree of uncertainty (you are using 250 points to estimate a GARCH
> > model) etc The surprise would be if the MEAN of those value was positive.
> > 
> > Am not too familiar with Vine Copulas and the particular package (or if there is something particular about the way the simulation is done) so I may have missed something...
> > 
> > -Alexios
> > 
> > Note also the following:
> > a. For T+1 there is not uncertainty about sigma (see http://unstarched.net/r-examples/rugarch/a-note-on-simulation-in-the-rugarch-package/
> > if you are confused about this point).
> > b. All distributions in rugarch (and their functions qdist, rdist, pdist, ddist) are based on the standardized (0,1) representation (all have the location-scale invariance property).
> > 
> > 
> > ################################################################################
> > for(i in 1:10)
> > {
> > print(i)
> > windows <- window(returns_crisis[,1:20], start=times[376-250-24+i,1],
> > end=times[376-25+i,1])
> > fit <- lapply(windows, ugarchfit, spec=model, solver="hybrid") print("rugarch fitting done") residuals <- sapply(fit, residuals, standardize=TRUE) shape <- sapply(fit, coef) shape <- shape[5,] # # I don't know where you found pged, but the correct version is # pdist("ged",....) UniformResiduals <- sapply(1:ncol(residuals), function(i){pdist("ged", residuals[,i], mu=0, sigma=1, shape = shape[i])}) rvine <- RVineStructureSelect(UniformResiduals, indeptest=TRUE,
> > familyset=familyset)
> > print(paste(i,"RVine fitting done"))
> > set.seed(100)
> > sim <- RVineSim(100, rvine)
> > print(paste(i,"RVine simulation done"))
> > for(k in 1:20)
> > {
> > # I don't know where you found qged, but the correct version is # qdist("ged",....)
> > residuals2 <- qdist("ged", mu=0, sigma=1, sim[,k], shape = shape[k]) # Comment: custom.dist$name must be a character (anything but not NA or # NULL).
> > # Also use startMethod="sample" to use the last values of the estimated # object for the simulation rather than the unconditional values.
> > rvine_sim <- ugarchsim(fit[[k]], n.sim=1, m.sim=100, custom.dist = list(name="sample", distfit=matrix(residuals2, ncol=100)),
> > startMethod="sample")
> > # cleaner to calculate at each step
> > rvine_fitted[i,,k] <- as.numeric(fitted(rvine_sim)) rvine_sigma[i,,k] <- as.numeric(sigma(rvine_sim)) VaR01[i,,k] <- rvine_fitted[i,,k] + rvine_sigma[i,,k] * qdist('ged', 0.01, mu=0, sigma=1, shape = shape[k]) VaR05[i,,k] <- rvine_fitted[i,,k] + rvine_sigma[i,,k] * qdist('ged', 0.05, mu=0, sigma=1, shape = shape[k]) } print(paste(i, "rugarch forecast done")) } ################################################################################
> > 
> > On 22/08/2014 00:44, Ole Bueker wrote:
> > > Sorry for sending another message, but it seems like there is still an 
> > > error in the formatting of the code..
> > > I'll share the R code using dropbox, in case the code is still not 
> > > readable in this
> > > email: https://www.dropbox.com/s/rc9kpanule52zu8/minimum_working.R
> > > 
> > > #Load data and define variables
> > > returns <- read.zoo("E:/Dropbox/my own/Programming/R/returns.csv", 
> > > header=TRUE, sep=",", format="%d-%m-%y")
> > > model<-ugarchspec(variance.model=list(model="sGARCH",garchOrder=c(1,1)
> > > ),mean.model=list(armaOrder=c(1,0),include.mean=FALSE),distribution.mo
> > > del="ged") times <- as.data.frame(time(returns_crisis))
> > > windows <- matrix(0, 20, 250)
> > > familyset <- c(1:5, 7, 10, 13, 14, 17, 20) sim <- array(0, dim = 
> > > c(100, 20))
> > > residuals2 <- array(0, dim = c(100, 20)) rvine_fitted <- array(0, dim 
> > > = c(10,100,20)) rvine_sigma <- array(0, dim = c(10,100,20))
> > > VaR01 = VaR05 = array(0, dim = c(10,100,20))
> > > 
> > > #Main calculation
> > > for(i in 1:10)
> > > {
> > >   print(i)
> > >   windows <- window(returns_crisis[,1:20], start=times[376-250-24+i,1],
> > > end=times[376-25+i,1])     # Step 1
> > >   fit <- lapply(windows, ugarchfit, spec=model, solver="hybrid")
> > >   print("rugarch fitting done")
> > > 
> > >   residuals <- sapply(fit, residuals, standardize=TRUE)  # Step 2     
> > >   shape <- sapply(fit, coef)
> > >   shape <- shape[5,]
> > > 
> > >   UniformResiduals <- pged(residuals, nu = shape)     # Step 3
> > >   if(any(UniformResiduals > 0.99999))
> > >   {
> > >     ix = which(UniformResiduals > 0.99999)
> > >    UniformResiduals [ix] = 0.99999
> > >   }
> > >   if(any(UniformResiduals < .Machine$double.eps))
> > >   {
> > >     ix = which(UniformResiduals < (1.5*.Machine$double.eps))
> > >     UniformResiduals [ix] = .Machine$double.eps
> > >   }
> > >   rvine <- RVineStructureSelect(UniformResiduals, indeptest=TRUE,
> > > familyset=familyset)   # Step 4
> > >   print(paste(i,"RVine fitting done"))
> > >   for(j in 1:100)
> > >   {
> > >   sim[j,] <- RVineSim(1, rvine)      # Step 5
> > >   }
> > >   print(paste(i,"RVine simulation done"))
> > >   for(k in 1:20)
> > >   {
> > >   residuals2[,] <- qged(sim[,], nu = shape[k])         # Step 6
> > > 
> > >   residuals_temp <- residuals2[,k]
> > >   rvine_sim <- ugarchsim(fit[[k]], n.sim=1, m.sim=100, custom.dist = 
> > > list(name=NA, distfit=residuals_temp))  # Step 7
> > > 
> > >   rvine_fitted[i,,k] <- fitted(rvine_sim)  # Step 8
> > >   rvine_sigma[i,,k] <- sigma(rvine_sim)  
> > >   if (i==10)
> > >     {
> > >     for(j in 1:100)
> > >     {
> > >     VaR01[,j,k] <- rvine_fitted[,j,k] + rvine_sigma[,j,k] * 
> > > qdist('ged', 0.01, mu=0, sigma=1, shape = shape[k])  # Step 9
> > >     VaR05[,j,k] <- rvine_fitted[,j,k] + rvine_sigma[,j,k] * 
> > > qdist('ged', 0.05, mu=0, sigma=1, shape = shape[k])
> > >     }
> > >     }else {}
> > >   }
> > >   print(paste(i, "rugarch forecast done")) }
> > > 
> > > 
> > > #Cleanup
> > > remove(i, j, k, ix, familyset, model, sim, rvine_sim, residuals,
> > > residuals_temp)
> > > 
> > > 
> > > ----------------------------------------------------------------------
> > > --
> > > 
> > >> From: ole.bueker at outlook.com
> > >> To: alexios at 4dscape.com; r-sig-finance at r-project.org
> > >> Date: Fri, 22 Aug 2014 01:10:55 +0200
> > >> Subject: Re: [R-SIG-Finance] rugarch + VineCopula for value at risk
> > >>
> > >>
> > >>
> > >>
> > >> I have no idea what happened to the formatting, looked fine to me 
> > >> when
> > > I typed it in Outlook.
> > >> Here's the revised version, with amended code so that it is minimally
> > > reproducible (for this purpose, I will only use the first 25 
> > > companies,
> > > 100 simulations for 10 trading days).
> > >> The revised code should take around 5 minutes,
> > >>
> > >>
> > >> Again, heres the code summary:
> > >> 1. Fit GARCH models to each series.
> > >> 2. Extract standardized returns.
> > >> 3. Transform standardized returns to uniform marginals using the
> > > parametric IFM method by Joe.
> > >> 4. Fit vine copulas.
> > >> 5. Generate 100 1-day ahead forecasts from the vine copulas.
> > >> 6. Reverse transform the simulated values.
> > >> 7. Use these transformed forecasts in ugarchsim (using custom.dist) 
> > >> 8. Extract forecasted mu and sigma.
> > >> 9. Calculate 95% and 99% VaR.
> > >>
> > >> The full code:
> > >> #Load data and define variablesreturns <- read.zoo("E:/Dropbox/my
> > > own/Programming/R/returns.csv", header=TRUE, sep=",", 
> > > format="%d-%m-%y")
> > >>
> > > model<-ugarchspec(variance.model=list(model="sGARCH",garchOrder=c(1,1)
> > > ),mean.model=list(armaOrder=c(1,0),include.mean=FALSE),distribution.mo
> > > del="ged")times
> > > <- as.data.frame(time(returns_crisis))windows <- matrix(0, 20, 
> > > 250)familyset <- c(1:5, 7, 10, 13, 14, 17, 20)sim <- array(0, dim = 
> > > c(100, 20))residuals2 <- array(0, dim = c(100, 20))rvine_fitted <- 
> > > array(0, dim = c(10,100,20))rvine_sigma <- array(0, dim =
> > > c(10,100,20))VaR01 = VaR05 = array(0, dim = c(10,100,20))
> > >> #Main calculation
> > >> for(i in 1:10){ print(i) windows <- window(returns_crisis[,1:20],
> > > start=times[376-250-24+i,1], end=times[376-25+i,1]) # Step 1 fit <- 
> > > lapply(windows, ugarchfit, spec=model, solver="hybrid") print("rugarch 
> > > fitting done")
> > >> residuals <- sapply(fit, residuals, standardize=TRUE) # Step 2 shape
> > > <- sapply(fit, coef) shape <- shape[5,]
> > >> UniformResiduals <- pged(residuals, nu = shape) # Step 3
> > > if(any(UniformResiduals > 0.99999)) { ix = which(UniformResiduals >
> > > 0.99999) UniformResiduals [ix] = 0.99999 } if(any(UniformResiduals <
> > > .Machine$double.eps)) { ix = which(UniformResiduals <
> > > (1.5*.Machine$double.eps)) UniformResiduals [ix] = .Machine$double.eps 
> > > } rvine <- RVineStructureSelect(UniformResiduals, indeptest=TRUE,
> > > familyset=familyset) # Step 4 print(paste(i,"RVine fitting done")) 
> > > for(j in 1:100) { sim[j,] <- RVineSim(1, rvine) # Step 5 } 
> > > print(paste(i,"RVine simulation done")) for(k in 1:20) { residuals2[,]
> > > <- qged(sim[,], nu = shape[k]) # Step 6
> > >> residuals_temp <- residuals2[,k] rvine_sim <- ugarchsim(fit[[k]],
> > > n.sim=1, m.sim=100, custom.dist = list(name=NA, 
> > > distfit=residuals_temp)) # Step 7
> > >> rvine_fitted[i,,k] <- fitted(rvine_sim) # Step 8 rvine_sigma[i,,k] <-
> > > sigma(rvine_sim) if (i==10) { for(j in 1:100) { VaR01[,j,k] <- 
> > > rvine_fitted[,j,k] + rvine_sigma[,j,k] * qdist('ged', 0.01, mu=0, 
> > > sigma=1, shape = shape[k]) # Step 9 VaR05[,j,k] <- rvine_fitted[,j,k] 
> > > + rvine_sigma[,j,k] * qdist('ged', 0.05, mu=0, sigma=1, shape = 
> > > shape[k]) } }else {} } print(paste(i, "rugarch forecast done"))}
> > >>
> > >> #Cleanupremove(i, j, k, ix, familyset, model, sim, rvine_sim,
> > > residuals, residuals_temp)
> > >>
> > >> Hope this time everything is formatted correctly!
> > >> After running the code, I export VaR01 and VaR05 to Excel, and notice
> > > that there's usually quite a few positive values.
> > >> This seems to happen because sometimes the sigma is too low compare 
> > >> to
> > > the (positive) forecast, and therefore the VaR stays in the positive.
> > >> My guess is that either there's a mistake in my code at step 3 (maybe
> > > the wrong shape is used?).Or more likely there's a mistake at step 6/7 
> > > (wrong transformation or transformed values are not inserted correctly 
> > > into ugarchsim).
> > >> > Date: Thu, 21 Aug 2014 18:21:49 +0100
> > >> > From: alexios at 4dscape.com
> > >> > To: ole.bueker at outlook.com; r-sig-finance at r-project.org
> > >> > Subject: Re: [R-SIG-Finance] rugarch + VineCopula for value at risk
> > >> >
> > >> > Ole,
> > >> >
> > >> > The email you've sent is really badly formatted and all over the place.
> > >> >
> > >> > 1. Try to send a well structured text-only email (with no special 
> > >> > characters).
> > >> > 2. Make an effort to send a minimally reproducible example.
> > >> > Telling us that your problem took 8 hours to run is not the way to go.
> > >> > Take a small subsample of your data which generates the behavior 
> > >> > you want to highlight (try to find one that does), and use that so 
> > >> > that we can investigate with minimal fuss.
> > >> >
> > >> > -Alexios
> > >> >
> > >> > On 21/08/2014 17:50, Ole Bueker wrote:
> > >> > > Anyway, here s my code so far:# Load Data and define variables 
> > >> > > returns <- read.zoo("E:/Dropbox/my 
> > >> > > own/Programming/R/returns.csv",
> > > header=TRUE, sep=",",
> > > format="%d-%m-%y")model<-ugarchspec(variance.model=list(model="sGARCH"
> > > ,garchOrder=c(1,1)),mean.model=list(armaOrder=c(1,0),include.mean=FALS
> > > E),distribution.model="ged")times
> > > <- as.data.frame(time(returns))windows <- matrix(0, 112, 250)familyset
> > > <- c(1:5, 7, 10, 13, 14, 17, 20) # The vine copulas to be testedsim <- 
> > > array(0, dim = c(1000, 112))residuals2 <- array(0, dim = c(1000, 
> > > 112))rvine_fitted <- array(0, dim = c(25,1000,112))rvine_sigma <- 
> > > array(0, dim = c(25,1000,112))VaR01 = VaR05 = array(0, dim = 
> > > c(25,1000,112))
> > >> > >
> > >> > > #Main calculation
> > >> > > for(i in 1:25){ print(i) windows <- window(returns_crisis,
> > > start=times[376-250-24+i,1], end=times[376-25+i,1]) #Define the moving 
> > > window fit <- lapply(windows, ugarchfit, spec=model, solver="hybrid") 
> > > #Fit the garch models print("rugarch fitting done") residuals <- 
> > > sapply(fit, residuals, standardize=TRUE) #Extract residuals & shape 
> > > parameters shape <- sapply(fit, coef) shape <- shape[5,] 
> > > UniformResiduals <- pged(residuals, nu = shape) #Transform residuals 
> > > into uniform marginals if(any(UniformResiduals > 0.99999)) { ix = 
> > > which(UniformResiduals > 0.99999) UniformResiduals [ix] = 0.99999 } 
> > > if(any(UniformResiduals < .Machine$double.eps)) { ix = 
> > > which(UniformResiduals < (1.5*.Machine$double.eps)) UniformResiduals 
> > > [ix] =
> > >> > .Machine$double.eps } rvine <-
> > > RVineStructureSelect(UniformResiduals, indeptest=TRUE,
> > > familyset=familyset) #Fit the Vine copulas print(paste(i,"RVine 
> > > fitting
> > > done")) for(j in 1:1000) #Simulate 1000 1-day ahead using VineCopula { 
> > > sim[j,] <- RVineSim(1, rvine) # 1000 x 112 matrix of forecasts 
> > > }print(paste(i,"RVine simulation done")) for(k in 1:112) #Next:
> > > ugarchsimulation for all 112 companies { residuals2[,] <- qged(sim[,], 
> > > nu = shape[k])
> > >> > # 1000 x 112 matrix of standardized residuals residuals_temp <-
> > > residuals2[,k] # 1000 x 1 vector of standardized residuals for 
> > > individual company rvine_sim <- ugarchsim(fit[[k]], n.sim=1, 
> > > m.sim=1000, custom.dist = list(name=NA, distfit=residuals_temp)) #1000 
> > > simulations using the standardized residuals from Vine copula models 
> > > for ugarchfit rvine_fitted[i,,k] <- fitted(rvine_sim) #Extract 
> > > forecasted values - 25 x 1000 x 112 rvine_sigma[i,,k] <- 
> > > sigma(rvine_sim) #Extract forecasted sigmas - 25 x 1000 x 112 for(j in 
> > > 1:1000)
> > >> > #Next: Value at risk { VaR01[,j,k] <- rvine_fitted[,j,k] +
> > > rvine_sigma[,j,k] * qdist('ged', 0.01, mu=0, sigma=1, shape = 
> > > shape[k]) #Value at risk for 99% quantile VaR05[,j,k] <- 
> > > rvine_fitted[,j,k] + rvine_sigma[,j,k] * qdist('ged', 0.05, mu=0, 
> > > sigma=1, shape = shape[k]) #Value at risk for 95% quantile } 
> > > }}remove(i, j, k) #Cleanupremove(windows, fit, residuals, shape, 
> > > residuals2, residuals_temp, rvine, sim, rvine_sim) #Cleanup Hope I 
> > > didn t make any mistakes in my approach, but it seems like this is the 
> > >  standard  copula
> > > + rugarch approach   if anyone is familiar with this, I am open to
> > > suggestions on how to speed up the si
> > >> > mulations.
> > >> > > So far, so good   the problem I am facing now:
> > >> > > Some (only a few) of my value at risk values are positive..I have
> > > manually checked and it seems like the fitted value is much larger 
> > > than the sigma, so Value at Risk is positive   which doesn t really 
> > > make any economic sense to me.
> > >> > > Here s a dropbox link to the returns.csv, in case anyone is
> > > interested in running my code:
> > > https://www.dropbox.com/s/69i5959f3h4kweb/returns.csv
> > >> > >
> > >> > > Best Regards,
> > >> > > Ole
> > >> > > [[alternative HTML version deleted]]
> > >>
> > >>
> > >> [[alternative HTML version deleted]]
> > >>
> > >> _______________________________________________
> > >> R-SIG-Finance at r-project.org mailing list 
> > >> https://stat.ethz.ch/mailman/listinfo/r-sig-finance
> > >> -- Subscriber-posting only. If you want to post, subscribe first.
> > >> -- Also note that this is not the r-help list where general R
> > > questions should go.
> > 
>  		 	   		   		 	   		  
> 	[[alternative HTML version deleted]]
> 
> _______________________________________________
> R-SIG-Finance at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-finance
> -- Subscriber-posting only. If you want to post, subscribe first.
> -- Also note that this is not the r-help list where general R questions should go.
 		 	   		  
	[[alternative HTML version deleted]]



More information about the R-SIG-Finance mailing list