One-Factor Gaussian Copulaモデルでデフォルト相関を計算してみた

R

CDOのプライシングに使われる（ていた？）One-Factor Gaussian Copulaモデルを使ってデフォルト相関を計算。
具体的にはThe One-Factor Gaussian Copula Applied To CDOs: Just Say NO (Or, If You See A Correlation Smile, She Is Laughing At Your “Results”)のFigure 4.を再現してみた。

このモデルを世に知らしめたのは「On Default Correlation: A copula function approach 」by David X. Li of the RiskMetrics Groupの論文っぽい。

んで、今回やる問題の設定としては「デフォルト確率が等しい２つの企業があって、その企業価値の相関とデフォルト確率に対して、デフォルト確率の相関はどう振舞うのか？」ということ。
文章だと非常に分かり辛いんで上の論文のFigure 2.を見てもらった方がはやいかも。
ちなみにこの著者は「こんなクソみたいなモデリングしてんじゃねー！」って立ち位置の模様。僕もそう思う。

モデル自体が正規分布なんでちゃんと計算できるんだろうけど、モンテカルロが好きなのでモンテカルロで計算。
マシンはwindowsXP32bit,Core2のメモリ2GB。
まずは1000サンプルで計算

結構ギザギザしてる。

次に1000サンプルで計算（5分くらい）

ちょっとはマシ？これ以上Nをデカくすると計算が終わらん・・・

以下、ソースコード。Nの値は適当に。

N   <- 10000
#デフォルト相関計算する関数
CorrelationOfDefault <- function(probability.default,cor.copula){
  #シングルファクターモデル用の共通ファクター
  common.factor <-  matrix(rep(rnorm(N),2),N)
  #構造型アプローチなので基準化された資産価値出す
  value.asset   <- sqrt(cor.copula) * common.factor + sqrt(1 - cor.copula) * matrix(rnorm(2*N),N)
  #倒産判定
  index.default <- ifelse(pnorm(value.asset) < probability.default,1,0)
  return(cor(index.default)[1,2])
}
#デフォルト確率と資産相関（アセット相関）
probability.default <- seq(0.05,0.95,0.01)
cor.copula          <- seq(0.05,0.95,0.01)
#outerで回すためのラッパー関数
wrapper <- function(x, y, my.fun, ...) {
  sapply(seq(along = x), FUN = function(i)my.fun(x[i], y[i], ...))
}
#各デフォルト確率・資産相関に対してデフォルト相関を計算
cor.default <- outer(probability.default, cor.copula,FUN=wrapper,my.fun=CorrelationOfDefault)
#描画
persp(probability.default, cor.copula, cor.default, theta = -60, phi = 30,ltheta = -30,
  expand = 0.5, col = rainbow(length(probability.default)),shade = 1, ticktype = "detailed",)

outer関数をそのまま使うとベクトル丸ごと渡されちゃって、上のコードでいうifelseの箇所がうまく計算できないので、ラッパー関数を噛ませている。
ラッパー関数自体は7.19 ある関数について outer() を実行すると挙動が変なのですがなぜでしょうか? (Why does outer() behave strangely with my function?)からもらってきてた。
ありがたい。