2012-09-29 24 views
5

Aquí está el código (lo siento si es tan largo, pero fue el primer ejemplo que tuve); Estoy usando el ejemplo de CVaR CreditMetrics paquete por A. Wittmann y DEoptim solucionador de optimizar:Cómo establecer la suma de los parámetros en 1 en la optimización restringida

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)] 

Ahora escribo mi función ...

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

... y quiero para optimizar que:

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

¿me puede decir lo que tengo que insertar en # ... hacer sum(w) = 1 durante la optimización?

A continuación os muestro los resultados de optimización de acuerdo con los consejos de flodel:

# 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 

... y ...

Como se puede ver, el primer truco funciona mejor en que encuentra un resultado el cual es más pequeño que el segundo. Desafortunadamente, parece que lleva más tiempo.

# 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 

¿Algún comentario?

¿Debo aumentar el número de iteraciones debido a una función "demasiado estocástica" para ser optimizada?

+0

En el código para el "segundo truco" anterior, olvidaste agregar 'w <- w/sum (w)' en el cuerpo de 'fun'.¿Puedes actualizar tu código y resultados? – flodel

+0

Gracias, acabo de actualizar. ¿Cuál sería el mejor método según estos resultados? –

+0

Gracias. Tenga en cuenta que nunca sugerí lo que todavía está etiquetado como "segundo truco": está mal y debe eliminarlo. Lo que sugerí es lo que actualmente se etiqueta como "primer truco" y "tercer truco". Sugerí otro método: para usar 'w <- c (w, 1-sum (w))' en el cuerpo de 'fun' y afuera, ¿has probado eso también? Creo que este último método podría ser un poco más robusto y rápido. – flodel

Respuesta

6

Probar:

w <- w/sum(w) 

y si DEoptim le da una solución óptima w* tal que sum(w*) != 1 continuación w*/sum(w*) debería ser la solución óptima.

Otro enfoque es resolver todas sus variables excepto una. Sabemos que el valor de la última variable debe ser 1 - sum(w) por lo que en el cuerpo de la función, tener:

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

y hacer lo mismo a la solución óptima devuelto por DEoptim: w* <- c(w*, 1-sum(w*))

Ambas soluciones requieren que se vuelva a formular su problema en una optimización no restringida (sin contar para límites variables) para que se pueda usar DEoptim; lo que te obliga a hacer un poco de trabajo extra fuera del DEoptim para recuperar la solución al problema original.

En respuesta a su comentario, si quiere DEoptim para darle la respuesta correcta de inmediato (es decir, sin necesidad de una transformación posterior), también puede intentar incluir un costo de penalización para su función objetivo: por ejemplo agregue B * abs(sum(w)-1) donde B es un número arbitrario grande por lo que sum(w) se forzará a 1.

+0

Gracias, flodel, esa era la solución a la que estaba acostumbrado. Pero me gustaría que la salida de opitmización ya sume a 1. Hace pocos meses recuerdo de un ejemplo de función que tenía esta restricción incorporada en la función en sí, pero no puedo recuperarla :( –

+0

Gracias por su sugerencia, flode Desafortunadamente no soy experto en optimización limitada, entonces voy a mostrarle lo que obtuve utilizando ambas sugerencias. El «gran truco B» fue más lento en la convergencia: de hecho, trató de minimizar 'B' luego siguió la primera parte de mi ecuación, la normalización a 1 de los parámetros permitió una convergencia más rápida pero, después de la normalización, obtuve una solución subóptima. Voy a editar mi pregunta para mostrarle mis códigos ... –

0

Creo que debería agregar una penalización por cualquier desviación de uno. Agregue a su problema de minimización el término +(sum(weights) - 1)^2 * 1e10. ¡Debería ver que esta gran penalización forzará los pesos a sumar a 1!

Cuestiones relacionadas