CODEWITHSUNDEEP
roptimizationconstraints
Lisa Ann
Lisa Ann

Reputation: 3485

How to set parameters' sum to 1 in constrained optimization

Here's the code (I'm sorry if it's so long, but it was the first example I had); I'm using the CVaR example from CreditMetrics package by A. Wittmann and DEoptim solver to optimize:

library(CreditMetrics)
library(DEoptim)

N <- 3
n <- 100000
r <- 0.003
ead <- rep(1/N,N)
rc <- c("AAA", "AA", "A", "BBB", "BB", "B", "CCC", "D")
lgd <- 0.99
rating <- c("BBB", "AA", "B")   
firmnames <- c("firm 1", "firm 2", "firm 3")
alpha <- 0.99

# correlation matrix
rho <- matrix(c(  1, 0.4, 0.6,
                  0.4,   1, 0.5,
                  0.6, 0.5,   1), 3, 3, dimnames = list(firmnames, firmnames),
              byrow = TRUE)

# one year empirical migration matrix from standard&poors website
rc <- c("AAA", "AA", "A", "BBB", "BB", "B", "CCC", "D")
M <- matrix(c(90.81,  8.33,  0.68,  0.06,  0.08,  0.02,  0.01,   0.01,
              0.70, 90.65,  7.79,  0.64,  0.06,  0.13,  0.02,   0.01,
              0.09,  2.27, 91.05,  5.52,  0.74,  0.26,  0.01,   0.06,
              0.02,  0.33,  5.95, 85.93,  5.30,  1.17,  1.12,   0.18,
              0.03,  0.14,  0.67,  7.73, 80.53,  8.84,  1.00,   1.06,
              0.01,  0.11,  0.24,  0.43,  6.48, 83.46,  4.07,   5.20,
              0.21,     0,  0.22,  1.30,  2.38, 11.24, 64.86,  19.79,
              0,     0,     0,     0,     0,     0,     0, 100
)/100, 8, 8, dimnames = list(rc, rc), byrow = TRUE)

cm.CVaR(M, lgd, ead, N, n, r, rho, alpha, rating)

y <- cm.cs(M, lgd)[which(names(cm.cs(M, lgd)) == rating)]

Now I write my function...

fun <- function(w) {
  # ... 
  - (t(w) %*% y - r) / cm.CVaR(M, lgd, ead = w, N, n, r, 
                           rho, alpha, rating)
}

...and I want to optimize it:

DEoptim(fn = fun, lower = rep(0, N), upper = rep(1, N), 
        control = DEoptim.control())

Can you tell me what do I have to insert in # ... to make sum(w) = 1 during optimization?

Below I show you optimization results according to flodel's tips:

# The first trick is to include B as large number to force the algorithm to put sum(w) = 1

fun <- function(w) {
  - (t(w) %*% y - r) / cm.CVaR(M, lgd, ead = w, N, n, r, rho, alpha, rating) + 
    abs(10000 * (sum(w) - 1))
}

DEoptim(fn = fun, lower = rep(0, N), upper = rep(1, N), 
        control = DEoptim.control())

$optim$bestval
[1] -0.05326055

$optim$bestmem
par1        par2        par3 
0.005046258 0.000201286 0.994752456

parsB <- c(0.005046258, 0.000201286, 0.994752456)

> fun(parsB)
            [,1]
[1,] -0.05326089

...and...

As you can see, the first trick works better in that he finds a results which is smaller than the second one. Unfortunately it seems he takes longer.

# The second trick needs you use w <- w / sum(w) in the function itself

fun <- function(w) {
  w <- w / sum(w)
  - (t(w) %*% y - r) / cm.CVaR(M, lgd, ead = w, N, n, r, rho, alpha, rating) #+ 
    #abs(10000 * (sum(w) - 1))
}

DEoptim(fn = fun, lower = rep(0, N), upper = rep(1, N), 
        control = DEoptim.control())

$optim$bestval
[1] -0.0532794

$optim$bestmem
par1         par2         par3 
1.306302e-15 2.586823e-15 9.307001e-01

parsC <- c(1.306302e-15, 2.586823e-15, 9.307001e-01)
parC <- parsC / sum(parsC)

> fun(parC)
           [,1]
[1,] -0.0532794

Any comment?

Should I increase the number of iterations because of a "too-stochastic" to-be-optimized-function?

Upvotes: 6

Views: 7889

Answers (3)

Karan Samlal
Karan Samlal

Reputation: 1

With the trick you applied:

fun <- function(w) {
  w <- w / sum(w)
  - (t(w) %*% y - r) / cm.CVaR(M, lgd, ead = w, N, n, r, rho, alpha, rating) #+ 
    #abs(10000 * (sum(w) - 1))
}

Why would you not use optim in this case? I think it will be much faster.

Upvotes: 0

Marco
Marco

Reputation: 21

I think you should add a penalty for any deviation from one. Add to your minimizing problem the term +(sum(weights) - 1)^2 * 1e10. You should see that this huge penalty will force the weights to sum to 1!

Upvotes: 2

flodel
flodel

Reputation: 89057

Try:

w <- w / sum(w)

and if DEoptim gives you an optimal solution w* such that sum(w*) != 1 then w*/sum(w*) should be your optimal solution.

Another approach is to solve over all your variables but one. We know the value of the last variable must be 1 - sum(w) so in the body of the function, have:

w <- c(w, 1-sum(w))

and do the same to the optimal solution returned by DEoptim: w* <- c(w*, 1-sum(w*))

Both solutions require that you re-formulate your problem into an unconstrained (not counting for variable bounds) optimization so DEoptim can be used; which forces you to do a little extra work outside of DEoptim to recover the solution to the original problem.

In reply to your comment, if you want DEoptim to give you the correct answer right away (i.e. without the need for a post-transformation), you could also try to include a penalty cost to your objective function: for example add B * abs(sum(w)-1) where B is some arbitrary large number so sum(w) will be forced to 1.

Upvotes: 9

Related Questions