> from Ingmar Visser:
>I would like to optimize a (log-)likelihood function subject to a number of
>linear constraints between parameters. These constraints are equality
>constraints of the form A%*%theta=c, ie (1,1) %*% 0.8,0.2)^t = 1 meaning
>that these parameters should sum to one. Moreover, there are bounds on the
>individual parameters, in most cases that I am considering parameters are
>bound between zero and one because they are probabilities.
>My problems/questions are the following:
>1) constrOptim does not work in this case because it only fits inequality
>constraints, ie A%*%theta > = c
---
I was struggling with the same problem a few weeks ago in the portfolio optimization
context. You can impose equality constraints by using inequality constraints >= and <=
simultaneously. See the example bellow.
>As a result when providing starting values constrOptim exits with an error
>saying that the initial point is not feasible, which it isn't because it is
>not in the interior of the constraint space.
In constrOptim the feasible region is defined by ui%*%theta-ci >=0. My first attempt
to obtain feasible starting values was based on solving for theta from ui%*%theta =
ci. Some of the items in the right hand side of the feasible region are, however, very
often very small negative numbers, and hence, theta is not feasible. Next, I started
to increase ci by epsilon ("a slack variable") and checked if the result was feasible.
The code is in the example bellow. If ui is not a square matrix, theta is obtained by
simulation. It is helpfull to know the upper and lower limits of the theta vector. In
my case they are often [0,1]. Usually only 2-3 simulations are required.
In the example, the weights (parameters) have limits [0,1] such that their sum is
limited to unity, as in your case. See, how Amat and b0 are defined.
V <- matrix(c(0.12,0,0,0,0.21,0,0,0,0.28),3,3) # variances
Cor <- matrix(c(1,0.25,0.25,0.25,1,0.45,0.05,0.45,1),3,3) # correlations
sigma <- t(V)%*%Cor%*%V # covariance matrix
ER <- c(0.07,0.12,0.18) # expected returns
Dmat <- sigma # adopted from solve.QP
dvec <- rep(0,3) # " " "
k <- 3 # number of assets
reslow <- rep(0,k) # lower limits for
portfolio weights
reshigh <- rep(1,k) # upper limits for
portfolio weights
pm <- 0.10 # target return
rfree <- 0.05 # risk-free return
risk.aversion <- 2.5 # risk aversion parameter
####### RISKLESS = FALSE; SHORTS = TRUE ; CONSTRAINTS = TRUE
########################################
a1 <- rep(1, k)
a2 <- as.vector(ER)+rfree
a3 <- matrix(0,k,k)
diag(a3) <- 1
Amat <- t(rbind(a1, a2, a3, -a3))
b0 <- c(1,pm,reslow, -reshigh)
objectf.mean <- function(x) # object
function
{
return(sum(t(x)%*%ER-1/2*risk.aversion*t(x)%*%sigma%*%x)-(t(x)%*%ER-pm))
}
# Getting starting values for constrOptim
ui <- t(Amat)
ci <- b0
if (dim(ui)[1] == dim(ui)[2])
{
test <- F
cieps <- ci
while(test==F)
{
theta <- qr.solve(ui,cieps)
foo <- (ui%*%theta-ci) # check if the initial values are in
the feasible area
test <- all(foo > 0)
if(test==T) initial <- theta # initial values for portfolio weights
cieps <- cieps+0.1
}
}
if (dim(ui)[1] != dim(ui)[2]) # if Amat is not square, simulate the
starting values
{
test <- F
i <- 1
while(test==F)
{
theta <- runif(k, min = 0, max = 1)
foo <- (ui%*%theta-ci)
test <- all(foo > 0) # and check that theta is feasible
if (test == T) initial <- theta
i <- i+1
}
}
initial
res <- constrOptim(initial, objectf.mean, NULL, ui=t(Amat), ci=b0, control =
list(fnscale=-1))
res$par # portfolio weights (parameters)
y <- t(res$par)%*%ER # return on the portfolio
y
I hope this helps.
Hannu Kahra
Progetti Speciali
Monte Paschi Asset Management SGR S.p.A.
Via San Vittore, 37 - 20123 Milano, Italia
Tel.: +39 02 43828 754
Mobile: +39 333 876 1558
Fax: +39 02 43828 247
E-mail: [EMAIL PROTECTED]
Web: www.mpsam.it
"Le informazioni trasmesse sono da intendersi per esclusivo uso del destinatario e
possono contenere informazioni e materiale confidenziale e privilegiato. Qualsiasi
correzione, inoltro e divulgazione in qualsiasi forma e modo anche a tenore generale �
del tutto proibita. Se avete ricevuto per errore il presente messaggio, cortesemente
contattate subito il mittente e cancellate da qualsiasi supporto il messaggio e gli
allegati a voi giunti. Tutti gli usi illegali saranno perseguiti penalmente e
civilmente da Monte Paschi Asset Management SGR S.p.A."
"The information transmitted are intended only for the person or entity to wich it is
addressed and may contain confidential and/or privileged material. Any review,
retransmission, dissemination, taking of any action in reliance upon, or other general
use are strictly prohibited. If you received this in error, please contact immediately
the sender and delete the material from any computer. All the illegal uses will be
persecuted by Monte Paschi Asset Management SGR S.p.A."
______________________________________________
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
