I have been looking at CVaR function in the PerformanceAnalytics package as an option to use in portfolio optimization. However, I cannot figure out the reason behind discrepancy that the function is giving for my test portfolio.

Here is a reproducable example:


a = rnorm(n = 10, mean = 0, sd = 0.02)
b = rnorm(n = 10, mean = 0, sd = 0.03)

wght = c(0.4, 0.6)  ## assign weights

dt = cbind(a, b)

startDate = as.Date("2016-01-01") 
h = seq(1:10)
myDates = startDate + h

dt = as.data.frame(dt)  ## dataframe of single assets
rownames(dt) = myDates

tmp = sweep(dt, 2, wght, `*`)

P = as.data.frame(rowSums(tmp))  ## resultant portfolio
colnames(P) = "Portfolio"

My issue arises here:

MES1 = ES(dt, p = 0.95, method = "modified", portfolio_method = "component",
       weights = wght)
MES2 = ES(P, p = 0.95, method = "modified", portfolio_method = "component",
      weights = 1)

I get the same results of CVaR for the portfolio in both cases. However, when I do not mention portfolio_method = "component" and weights = 1 should I not get the same result?

ES(P, p = 0.95, method = "modified")

Here the function returns completely different result. I am currently inspecting the code on the authors' page, but would also like to know the views of those who went through this.

Some background information regarding the methodology and implementation can be found here and here.


It looks like there's a bit of a problem with the code when the option portfolio_method = "single" is used. In the code below, your three calculations with ES are repeated using method="gaussian" and then compared with expected shortfall calcuated directly from the expression using mu and sigma. I've done the direct calculation with dt using the cov() function and with P using the sd() function. The resuls of all calcuations should be identical. Similar to your case, all agree except the ES(P, portfolio_method = "single") calculation.

Looking into the ES code, the different portfolio_method options use different code to calculate sigma and higher moments. In particular, for gaussian calculations, ES calculates sigma using cov() for portfolio_method = "component" and a general higher moment calculation function centeredmoment() for portfolio_method = "single". The problem is that for n = number of rows in the return vector, cov() and sd() normalize the moment calculations using 1/(n-1) while the calculations using centeredmoment() use 1/n. I've demonstrated this in the code below where the calculation for ES_P_modied_sigma rescales sigma to use the 1/n normalization and gets the same result as ES(P, portfolio_method ="single") for the Gausian distributions.

The discrepancy is more apparent in your test calculation since you're using a small number of rows. As the number of rows increases, the descrepancy becomes much smaller. Nevertheless, it's something that should be corrected.

#  repeat ES calculations using Gaussian distributions
  MES1g = ES(dt, p = 0.95, method = "gaussian", portfolio_method = "component",
            weights = wght)
  MES2g = ES(P, p = 0.95, method = "gaussian", portfolio_method = "component",
            weights = 1)
  MES3g <- ES(P, p = 0.95, method = "gaussian")
#  calculate Gaussian shortfalls directly using dt and wght and then P 
  prob <- .05
  mean_dt <-  apply(dt,2,mean)%*%wght
  sigma_dt <-  sqrt(t(wght)%*%cov(dt)%*%wght)
  ES_dt <- -mean_dt + dnorm(qnorm(prob))*sigma_dt/prob

  mean_P <- mean(P[,1])
  sigma_P <- sd(P[,1])
  ES_P <- -mean_P + dnorm(qnorm(prob))*sigma_P/prob 
# redo calculation with P using sigma scaled by sqrt((nrow(P) - 1)/nrow(P))
  ES_P_modied_sigma <- -mean_P + dnorm(qnorm(prob))*sigma_P/prob*(sqrt((nrow(P)-1)/nrow(P))) 
# compare all results
# note that result ES_P_modified_sigma reproduces MES3g showing problem with normalization for sigma
  c(MES1g$ES, MES2g$ES, MES3g)
  c(ES_dt, ES_P, ES_P_modied_sigma)
  • $\begingroup$ Thank you WaltS. I cannot seem to find the calculation function centeredmoment() and its code in the repo, can you please give me a link? $\endgroup$
    – AK88
    Jun 16 '17 at 5:01
  • 1
    $\begingroup$ With PerformanceAnalytics loaded into your R session, type centeredmoment with no trailing parentheses at the command line. The code for the function should appear. centeredmoment is also described in the help docs for PerformananceAnalytics. $\endgroup$
    – WaltS
    Jun 16 '17 at 11:23

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.