モンテカルロシミュレーションでVasicekモデルの債券価格

Vasicekモデルでパスの発生方法を比較の続き、モンテカルロシミュレーションから債券価格を求めてみる

パラメータなどは前回と同じで、それぞれ1000回シミュレーション

nsim <- 1000
p1 <- replicate(nsim, Vasicek_euler(kappa, mu, sigma, r0, dt, n))
p2 <- replicate(nsim, Vasicek_exact(kappa, mu, sigma, r0, dt, n))
p3 <- replicate(nsim, sde.sim(X0 = r0, N = n, delta = dt, theta = c(kappa * mu, kappa, sigma), model = "OU")[-(n+1)])
p4 <- replicate(nsim, sde.sim(X0 = r0, N = n, delta = dt, method = "cdist", rcdist = rcOU, theta = c(kappa * mu, kappa, sigma))[-(n+1)])
p5 <- replicate(nsim, sde.sim(X0 = r0, N = n, delta = dt, drift = expression(0.2 * (0.015 - x)), sigma = expression(0.003))[-(n+1)])
p6 <- replicate(nsim, simulate(model, true.param = list(kappa = kappa, mu = mu, sigma = sigma), sampling = sampling)@data@original.data[-(n+1)])

replicateはシミュレーションしたパスをmatrixの形で返してくる。ここから、債券価格を求める関数を作る

MCZCBond <- function(r, dt) {
  rsum <- colSums(r)
  sbond <- exp(-rsum * dt)
  mean(sbond)
}

それぞれの方法について債券価格を計算する

> MCZCBond(p1, dt)
[1] 0.9162368
> MCZCBond(p2, dt)
[1] 0.9140108
> MCZCBond(p3, dt)
[1] 0.9145012
> MCZCBond(p4, dt)
[1] 0.9144614
> MCZCBond(p5, dt)
[1] 0.9147488
> MCZCBond(p6, dt)
[1] 0.9155063

明示的に求まるVasicekモデルの債券価格も計算する

> t <- 0
> T <- dt * n
> VasicekZCBond(kappa, mu, sigma, r0, t, T)$price
[1] 0.9148043

今回はsdeパッケージを使ったドリフト項とディフュージョン項を指定する方法が最も差が小さいという結果に

何回かシミュレーションをやってみると、必ずどの方法が明示的に求まる債券価格に近いということはなく、パスの発生方法の違いによって生じる誤差よりもシミュレーションした時に生じる誤差の方が大きいと言えそうです。モンテカルロシミュレーションに限らず、Vasicekモデルは基本を確認するには大変役立ちます