[R-SIG-Finance] Quadratic programming solve.QP's Lagrangians

[Bob Jansen]

Hi,

Careful this does rm(list=ls()) and has been cross posted to Quant.SE:

http://quant.stackexchange.com/q/17997/848

On Wed, May 27, 2015 at 2:41 PM, Ton Jean-Claude <jean-claude.ton at epfl.ch>
wrote:

> Hi I am currently playing a little with Markowitz's portfolio optimization
> problem and from my computations I get that if we start from the following
> problem :
>
> $$    \underset{w}{\arg \max} \quad w'r - \frac{\gamma}{2} w' Q w \quad
> \text{s.t}  \sum_i w_i=1, \quad w_i \geq 0\quad
> $$
>
> Then $\gamma$ should be of the form :
> $$\gamma = \frac{w'(r-\beta)-\alpha}{w' \Sigma w},$$
>
> where $\alpha$ is the Lagrangian multiplier due to the sum of weights
> being equal to one and $\beta$ the Lagrangian multipliers due to the no
> short selling condition.
> Using `solve.QP` from the package `quadprog` I have coded the
> mean-variance formulation. And the function `ComputeRho` which you may find
> attached below computes the $\gamma$ (Sorry for the bad notations)
>
>
> My problem is the following
> depending on the r.a (risk aversion that is $\gamma$ in my case),
> `ComputeRho` gets me the correct answer or not.
>
>     rm(list=ls())
>         ###################################### begin : packages
> ####################################
>         require(tseries)
>         require(quantmod)
>         require(quadprog)
>         require(MASS)
>         require(PerformanceAnalytics)
>         require(corpcor)
>         require(quantmod)
>         require(parallel)
>         require(ggplot2)
>         #require(mvtnorm)
>         ###################################### end : packages
> ######################################
>         LogRet <- function(x,lag=5){
>                 #require(quantmod)
>                 #return(na.omit(apply(x,2,function(x) Delt(x,k=lag))) )
>                 require(zoo)
>                 return(na.omit(diff(log(x),lag=lag)))
>         }
>
> ##########################################################################
>         # Input :        sol = solution from solve.qp
>         #                lr = logreturns
>         # Output :      rho
>         ComputeRho<-function(sol,lr){
>                 lr<-as.matrix(lr)
>                 alpha<-sol$lagrange_mults[1] # Lagrangian due to
> normalizing constraint
>                 beta<-sol$lagrange_mults[2:length(sol$lagrange_mults)] #
> Lagrangians due to no short sell constraint
>
>                 weights<-sol$weights
>                 cov.mat<-cov(lr)
>                 # unconstrained form
>                 uncons
> <-(t(weights)%*%colMeans(lr))/(t(weights)%*%cov.mat%*%weights)
>                 # due to Lagrangian multipliers
>                 cons1 <-
> (t(weights)%*%beta)/(t(weights)%*%cov.mat%*%weights)
>                 cons2 <- alpha/(t(weights)%*%cov.mat%*%weights)
>                 cons<-cons1+cons2
>
>                 rho<-uncons+cons # lagrangians change sign...
>                 rho
>         }
>
>
>         ############################################################
>         # Input :       data = usually logreturns
>         #               r.aversion = risk aversion coefficient usually
> between 1 to 10
>         # Output :      er = expected return
>         #               sd = standard deviation
>         #               weights = weights
>         #               Lagrangians = Lagrange multipliers
>
>         MeanVariance<-function(data,r.aversion=10,method="MeanVariance")
>         {
>                 # compute probability beliefs mean and risk ( here sample
> covariance)
>                 er<-colMeans(data)
>                 cov.mat<-cov(data)
>                 # get length of data
>                 n_asset <- length(er)
>
>                 # quad problem
>                 Dmat <- r.aversion*as.matrix(cov.mat)
>                 dvec <- er
>                 meq <-  1
>
>                 bvec <- c(1, rep(0,n_asset))
>                 Amat <- matrix( c(rep(1,n_asset), diag(nrow=n_asset)),
> nrow=n_asset)
>
>                 sol<-solve.QP(Dmat, dvec,Amat, bvec=bvec, meq=meq)
>
>                 # fetch returns from solve.QP
>                 weights<-sol$solution
>                 exp.ret <- t(er)%*%weights
>                 std.dev <- sqrt(weights %*% cov.mat %*% weights)
>
>                 # return results
>                 ret <- list(er = as.vector(exp.ret),
>                             sd = as.vector(std.dev),
>                             weights = weights,
>                             lagrange_mults=sol$Lagrangian,
>                             weights2=sol$unconstrained.solution)
>
>         }
>         # load asset prices
>         data <- read.table("~/data.csv", header=TRUE, quote="\"")
>         # risk aversion
>         r.a<- 17 # for values below 18 r.a is different from rho
>         lr<-LogRet(as.matrix(data))
>
>         sol<-MeanVariance(lr,r.aversion = r.a)
>         rho<-ComputeRho(sol,lr)
>         rho
>
>         r.a2<- 18 # for values below 18 r.a is different from rho
>         sol2<-MeanVariance(lr,r.aversion = r.a2)
>         rho2<-ComputeRho(sol2,lr)
>         rho2
>
>
> I was able to cern the problem to the sign of the Lagranians (here `cons1`
> and `cons2` in `ComputeRho` ) which switch signs depending on the `r.a.` I
> choose.
>
> Does someone knows how the Lagrangians are computed explicitly ? Since It
> seems that this is done in FORTRAN.
>
> Best,
> Jc
>
> here is the data set
> https://www.dropbox.com/s/ynlnmdgfm9w4rq0/data.csv?dl=0
>
>         [[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.
>



-- 
Met vriendelijke groeten,

Bob Jansen

	[[alternative HTML version deleted]]