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]]
_______________________________________________
[email protected] 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.
