2012-09-29 7 views
5

コードは次のとおりです(長すぎると申し訳ありませんが、最初の例でした)。 ...今、私は私の関数を記述制約のある最適化でパラメータの合計を1に設定する方法

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

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

を...と私は最適化したい:私は最適化するために、A.ヴィットマンによってCreditMetricsパッケージとDEoptimソルバーからCVARの例を使用していますそれ:

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

私は最適化中sum(w) = 1を作るために# ...に挿入する必要が何をすべきかを教えてもらえますか?あなたが見ることができるように、最初のトリックは、彼がどの結果を見つけたことで、より良い作品

# 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 

...と...

I以下

はあなたの最適化のflodelのヒントに従った結果を示しています2番目のものよりも小さい。残念ながら、彼は長くかかるようです。

# 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 

最適化機能が「あまりにも確率的」であるため、反復回数を増やす必要がありますか?

+0

上記の「第2のトリック」のコードでは、 'fun'の本体に' w < - w/sum(w) 'を追加するのを忘れてしまいました。コードと結果を更新していただけますか? – flodel

+0

ありがとう、私はちょうど更新しました。これらの結果によると、最良の方法は何でしょうか? –

+0

ありがとうございます。私はまだ「第2のトリック」とラベル付けされているものは提案していないことに注意してください。間違っていて、削除する必要があります。私が提案したのは、現在「第1のトリック」と「第3のトリック」とラベル付けされているものです。私はさらに別の方法を提案しました: 'fun'と体外で 'w < - c(w、1-sum(w))'を使用するには、私はこの最後の方法がもう少し堅牢で速くなると思う。 – flodel

答えて

6

試してみてください。

w <- w/sum(w) 

DEoptimはあなたに最適なソリューションを提供している場合 w* sum(w*) != 1は、その後 w*/sum(w*)あなたの最適なソリューションでなければならないことなど。

もう1つの方法は、1つではなくすべての変数を解決することです。

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

DEoptimによって返された最適解に同じことを実行します:w* <- c(w*, 1-sum(w*))

どちらのソリューションは、あなたのことを必要と我々は最後の変数の値は、関数の本体でとても1 - sum(w)ていなければなりません知っています問題を無制約(変数の境界を考慮しない)に再定式化するので、DEoptimを使用できます。この問題を解決するために、DEoptim以外の作業を余儀なくされることになります。

DEoptimがあなたに正しい答えをすぐに(換算の必要なしに)したい場合は、目的関数にペナルティコストを含めるようにすることもできます。 B * abs(sum(w)-1)を追加してください。Bは任意の大きい数値ですので、sum(w)は、1に強制されます。

+0

ありがとう、私が慣れていた解決策だったフラデル。しかし、私はopitmization ouputがすでに1になっていることを望んでいます。数ヶ月前に、この制約を関数自体に埋め込んだ関数の例を思い出しましたが、回復することはできません:( –

+0

あなたの提案をいただきありがとうございます。残念ながら私は制約の厳しい最適化に熟練していないので、私はあなたの提案を使って得たものを見せてくれるでしょう "大きな" Bトリック»は収束の方が遅いです:実際には、それは私の方程式の最初の部分に進みました; 1つのパラメータの正規化はより速い収束を可能にしましたが、正規化した後、私は下位最適解を明らかにしました。 –

0

私はあなたからの逸脱に対してペナルティを加えるべきだと思います。 最小化問題に「+(sum(weights) - 1)^2 * 1e10」という語句を追加してください。この巨額のペナルティは、合計が1になるよう強制されます。

関連する問題