モンテカルロ法

モンテカルロ法とは、乱数を用いた数値計算法の総称である。複雑なペイオフを持つ派生証券（デリバティブ）の中には、原資産価格のモデルが単純であっても、解析解に基づく価格公式を得られないものが少なくない。そのような場合でも、モンテカルロ法を用いることで数値的な評価が可能となる。

具体的な価格評価の手順としては、まずリスク中立測度に基づいて原資産価格のシミュレーション（見本路の発生）を行う。次に、得られた各パスにおけるペイオフの平均値を算出し、それを現在価値に割り戻すことで証券価格を推定する。

原証券価格 \(S_t\) がリスク中立測度 \(Q\) の下で次の確率微分方程式（幾何ブラウン運動モデル）に従うとする。

\[ \frac{dSt}{S_t}=rdt+\sigma dW_t \]

この確率微分方程式の解は次の式で与えられる。

\[ S_t=S_0\exp\left[\left(r-\frac{\sigma^2}{2}\right)t+\sigma W_t\right] \]

ヨーロピアン派生証券のペイオフ \(f(S_T)\) は満期 \(T\) における価格 \(S_T\) のみに依存する。標準ウィーナー過程 \(W_T\) は平均 \(0\)、分散 \(T\) の正規分布に従うので、シミュレーションのために生成される標準正規乱数を \(\epsilon_t\) として \(W_t=\sqrt T \epsilon_t\) とすれば、見本路に対応するペイオフを計算することができる。生成する見本路の数を \(M\)、対応する原資産価格を \(s^{(j)}_T\) とすれば、ヨーロピアン派生証券の価格について、モンテカルロ法による推定値 \(\pi\) を以下のように決定することができる。

\[ \begin{align} \pi&=e^{-rT}\mathbf{E}^Q[f(S_T)] \\ &\approx e^{-rT}\frac{1}{M}\sum_{j=1}^Mf(s^{(j)}_T) \end{align} \]

ヨーロピアン派生証券のように、ペイオフが満期価格 \(S_T\) のみに依存する場合、期待値計算を解析的に実行できる。たとえば、ヨーロピアン・コールの価格をリスク中立測度 \(Q\) の下での期待値 \(\pi = e^{-rT}\mathbf{E}^Q[(S_T - K)^+]\) として計算すると、以下のブラック・ショールズ公式が導かれる。

\[ C(S_0, t) = S_0 \Phi(d_1) - Ke^{-r(T-t)} \Phi(d_2) \]

ここで、\(\Phi(x)\) は標準正規分布の累積分布関数であり、\(d_1, d_2\) は以下の通り定義される。

\[ \begin{aligned} d_1 &= \frac{\ln\left(\frac{S_0}{K}\right) + \left(r + \frac{\sigma^2}{2}\right)(T-t)}{\sigma\sqrt{T-t}} \\ d_2 &= d_1 - \sigma\sqrt{T-t} = \frac{\ln\left(\frac{S_0}{K}\right) + \left(r - \frac{\sigma^2}{2}\right)(T-t)}{\sigma\sqrt{T-t}} \end{aligned} \]

以下のコードでは、行使価格 \(K\) のヨーロピアン・コールの価格をモンテカルロ法によって評価し、ブラックショールズモデルに基づく解析解と比較している。

# parameters
s_0   <- 110
K     <- 100
tenor <- 1
r     <- 0.1
sigma <- 0.2
M     <- 1000 # size of sample in one simulation
N     <- 1000 # number of simulations

# Monte Carlo Simulation
set.seed(42)
drift     <- (r - sigma ^ 2 / 2) * tenor
diffusion <- sigma * sqrt(tenor) * rnorm(M * N)
s_T       <- s_0 * exp(drift + diffusion)

payoff <- pmax(s_T - K, 0)
payoff_mat <- matrix(payoff, nrow = N, ncol = M)
premium_mc <- exp(-r * tenor) * rowMeans(payoff_mat)
cat("Monte Carlo Price: ", premium_mc[1:3], "...\n")

Monte Carlo Price:  22.16895 21.16513 20.4242 ...

# Black-Scholes Analytical Solution
d1 <- (log(s_0 / K) + (r + sigma^2 / 2) * tenor) / (sigma * sqrt(tenor))
d2 <- d1 - sigma * sqrt(tenor)

premium_bs <- s_0 * pnorm(d1) - K * exp(-r * tenor) * pnorm(d2)
cat("Black-Scholes Price:", premium_bs, "\n")

Black-Scholes Price: 21.24877

# distribution of estimated prices
library(ggplot2)
theme_set(theme_bw())
ggplot() +
  geom_density(aes(premium_mc)) +
  geom_vline(xintercept = premium_bs,
             linetype = "dashed") +
  labs(subtitle = "Monte Carlo Simulation of European Call Price",
       x = "Premium", y = "Density")

モンテカルロ法の利点は、ペイオフが原資産価格の推移に依存する経路依存型オプションの評価が可能なことにある。ここではその一例として、変動行使価格型ルックバック・コールオプションの価格評価を行う。このオプションの満期におけるペイオフは、満期価格と取引期間中の最低価格の差 \(\max(S_T-\min_{0\le t\le T}S_t, 0)\) として定義される。

# parameters
tenor <- 1
ntime <- 100 # number of time spans
dt    <- tenor / ntime
time  <- seq(0, tenor, by = dt)

s_0   <- 100
r     <- 0.1
sigma <- 0.2

M     <- 1000 # size of sample in one simulation

# simulation
set.seed(42)

drift     <- (r - sigma^2 / 2) * dt
diffusion <- matrix(sigma * sqrt(dt) * rnorm(M * ntime), nrow = M)
increments <- drift + diffusion

log_paths <- cbind(rep(0, M), t(apply(increments, 1, cumsum)))
paths     <- s_0 * exp(log_paths)

s_T    <- paths[, ntime + 1]
s_min  <- apply(paths, 1, min)
payoff <- pmax(s_T - s_min, 0)

premium_lb <- exp(-r * tenor) * mean(payoff)
cat("Monte Carlo Lookback Price: ", premium_lb)

Monte Carlo Lookback Price:  18.21477

df <- data.frame(
  id = rep(seq_len(100), each = ntime + 1),
  time = time,
  path = as.numeric(t(paths[1:100, ]))
)

ggplot(df, aes(x = time, y = path, group = id, alpha = id)) +
  geom_line() +
  scale_color_continuous() +
  labs(subtitle = "Sample Paths",
       x = "Time", y = "Asset Price") +
  theme(legend.position = "none")

ggplot() +
  geom_density(aes(x = payoff)) +
  geom_vline(xintercept = premium_lb,
             linetype = "dashed") +
  labs(subtitle = "Distribution of Dayoff",
       x = "Payoff", y = "Density")