1

I am trying to maximize the number N_ent through a 1x42 weighting vector (weight).

N_ent is calculated with the following function:

N_ent <- exp(-sum((((solve(pca$rotation[])) %*% t(weight))^2)*
  (pca$sdev^2)/(sum((((solve(pca$rotation[])) %*% t(weight))^2)*
  (pca$sdev^2)))*log((((solve(pca$rotation[])) %*% t(weight))^2)*
  (pca$sdev^2)/(sum((((solve(pca$rotation[])) %*% t(weight))^2)*(pca$sdev^2)))))) 

Though it looks quite complicated, the equation works fine and supplies me with N_ent = 1.0967 when equal weights of 0.0238 (1/42 = 0.0238) are used.

Further, none of the weights may be below -0.1 or above 1.

I am new to R have struggled to use both the optim() (ignoring my constraints) and constrOptim() functions, encountering the error

Error in match.arg(method) : 'arg' must be of length 1

when optim() was used and

Error in ui %*% theta : non-conformable arguments

when constrOptim() was used.

Any help on how to set up the code for such an optimization problem would be greatly appreciated.

RHertel
  • 23,412
  • 5
  • 38
  • 64
Tom.K
  • 13
  • 5
  • I use `nloptr` package and function for optimization with constraints. Please include an example dataset `pca` so I will write down the answer. – danas.zuokas Dec 29 '15 at 09:25
  • Thank you very much, I'm not sure how to upload the data but it is a 42x42 matrix and all elements are between -1 and 1 if that helps? Alternatively, perhaps I could somehow send you the data set? – Tom.K Dec 29 '15 at 10:06
  • Please read [this](http://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example) thread. – danas.zuokas Dec 29 '15 at 10:10
  • The dput command yielded output that was too large to copy from the console so I wrote it to this csv file, hope it works. https://www.dropbox.com/s/ga99tgt1q7a5rz4/pca.csv?dl=0 – Tom.K Dec 29 '15 at 10:43
  • Could you please give some detail about the problem you are solving - is it scientific or real-life problem. What is behind this entropy maximization? – danas.zuokas Dec 30 '15 at 07:25
  • Its a financial project which uses Principle Component Analysis (hence pca) in order to calculate uncorrelated returns. Consequently, entropy is a measure of portfolio diversification and the optimisation calculates the weights that maximise diversification. The methodology comes from Muecci (2010) if you're interested. http://www.master272.com/finance/GQ_docs/Meucci_10b.pdf – Tom.K Dec 30 '15 at 09:16

1 Answers1

1

Here is the solution using library nloptr.

library(nloptr)

pca <- dget('pca.csv')
#random starting point
w0 <- runif(42, -0.1, 1)
#things that do not depend on weight
rotinv <- solve(pca$rotation)
m2 <- pca$sdev^2

#function to maximize
N_ent <- function(w) {
  m1 <- (rotinv %*% w)^2
  -exp(-sum(m1 * m2 / sum(m1 * m2) * log(m1 * m2 / sum(m1 * m2))))
}

#call optimization function
optres <- nloptr(w0, N_ent, lb = rep(-0.1, 42), ub = rep(1, 42),
                 opts = list('algorithm' = 'NLOPT_LN_NEWUOA_BOUND', 'print_level' = 2, 'maxeval' = 1000, 'xtol_rel' = 0))

You can view result by optres$solution. For your particular problem I find NLOPT_LN_NEWUOA_BOUND algorithm giving best result of 42. You can view all available algorithms by nloptr.print.options(). Note that _XN_ in the names of the algorithms indicate these that do not require derivatives. In your case derivative computation is not that difficult. You can provide it and use algorithms with _XD_ in the names.

danas.zuokas
  • 4,551
  • 4
  • 29
  • 39
  • I'm trying to add an inequality constraint to this but keep getting the error: `Error: nlopt_add_inequality_mconstraint returned NLOPT_INVALID_ARGS`. Would it be possible to give me the format for that too? The constraint is A* w - B >= 0, where A=`matrix(runif(42,0.5,3), nrow=1, ncol=42)` and B is a constant. – Tom.K Jan 08 '16 at 10:25
  • You should use argument `eval_g_ineq` which takes function such that in the solution `g(w) <= 0`. So it should go something like: `g <- function(w) B - A %*% w`. – danas.zuokas Jan 08 '16 at 13:01