確率微分方程式とEuler-Maruyama

🎓 レベル：発展　|　重要度：A（必須） 📎 前提：伊藤の公式、幾何ブラウン運動

要点（BLUF）

確率微分方程式（SDE） $dX_t=\mu(X_t)\,dt+\sigma(X_t)\,dB_t$ は、決定論的な傾き（ドリフト）にランダムな揺れ（拡散）を加えた発展則。多くの連続確率モデルの統一形です。
Euler-Maruyama 法： $X_{k+1}=X_k+\mu(X_k)\Delta t+\sigma(X_k)\sqrt{\Delta t}\,Z_k$ （ $Z_k\sim N(0,1)$ ）で数値的に解きます。
収束次数：強収束（経路の一致）は $\sqrt{\Delta t}$ オーダー（次数0.5）、弱収束（分布・期待値の一致）は $\Delta t$ （次数1）。

概念

常微分方程式 $dx/dt=\mu(x)$ にブラウン運動の揺れ $\sigma(x)\,dB$ を加えたのが SDE です。傾きと揺れの両方が現在の状態に依存でき、平均回帰・幾何成長・拡散などを統一的に書けます。解析解が出るのは特殊な場合だけなので、Euler-Maruyama 法で時間を刻んで数値的にたどります。確率版オイラー法ですが、揺れの項が $\Delta t$ でなく $\sqrt{\Delta t}$ で入る点が決定的に違います。

数式による定式化

SDE（伊藤型）：

dX_t = \mu(X_t)\,dt + \sigma(X_t)\,dB_t, \qquad X_0 = x_0

Euler-Maruyama スキーム（刻み $\Delta t$ ）：

X_{k+1} = X_k + \mu(X_k)\,\Delta t + \sigma(X_k)\,\Delta B_k, \qquad \Delta B_k \sim N(0,\Delta t)

強誤差（同じブラウン運動の下での経路誤差）は

\mathbb{E}\big[\,|X^{\text{EM}}_T - X_T|\,\big] = O(\sqrt{\Delta t}) \quad (\text{強収束次数 } 0.5)

例：オルンシュタイン・ウーレンベック過程（OU） $dX_t=\theta(\mu-X_t)\,dt+\sigma\,dB_t$ は平均 $\mu$ に回帰し、定常分布は $N\!\big(\mu,\ \frac{\sigma^2}{2\theta}\big)$ 。

直観

要するに「毎ステップ、傾きで少し進め、サイコロで少し揺らす」。決定論的オイラー法に乱数を足すだけですが、揺れが $\sqrt{\Delta t}$ で入るのがミソ（ブラウン増分の大きさ）。だから刻みを細かくしても揺れは $\sqrt{\Delta t}$ でしか縮まず、経路を正確に追う（強収束）のは遅い（次数0.5）。一方、期待値や分布だけ合わせる（弱収束）なら次数1で速い。「個々の経路を当てる」のは「平均を当てる」より難しい、というのが収束次数の差です。

具体例

OU 過程が定常分布 $N(\mu,\sigma^2/2\theta)$ に達すること、Euler-Maruyama 法の強誤差が刻みを1/4にするたびに半分（次数0.5）になることを確認します。

import numpy as np
rng = np.random.default_rng(3)
theta, mu_ou, sigma = 1.0, 2.0, 0.5            # OU：平均回帰
T, n_steps, n_paths = 10.0, 2000, 100000
dt = T/n_steps
X = np.zeros(n_paths)
for _ in range(n_steps):
    X = X + theta*(mu_ou - X)*dt + sigma*rng.normal(0, np.sqrt(dt), n_paths)
print(f"OU 終端: 平均={X.mean():.4f}(定常{mu_ou}) 分散={X.var():.4f}"
      f"(定常 sigma^2/(2theta)={sigma**2/(2*theta):.4f})")

S0, muG, sigG, T2 = 1.0, 0.05, 0.2, 1.0        # GBM で強誤差（厳密解と比較）
def em_strong_error(ns):
    d = T2/ns
    r = np.random.default_rng(7)               # 厳密解とEMで同じBM
    npath = 20000
    dBk = r.normal(0, np.sqrt(d), size=(npath, ns))
    exact = S0*np.exp((muG - 0.5*sigG**2)*T2 + sigG*dBk.sum(axis=1))
    S = np.full(npath, S0)
    for k in range(ns):
        S = S + muG*S*d + sigG*S*dBk[:, k]
    return np.abs(S - exact).mean()
prev = None
for ns in [16, 64, 256, 1024]:
    e = em_strong_error(ns)
    print(f"EM 強誤差(n={ns:4d})={e:.5f}" + ("" if prev is None else f"  前との比={prev/e:.2f}"))
    prev = e
# OU 終端: 平均=1.9997(定常2.0) 分散=0.1262(定常 sigma^2/(2theta)=0.1250)
# EM 強誤差(n=  16)=0.00586
# EM 強誤差(n=  64)=0.00299  前との比=1.96
# EM 強誤差(n= 256)=0.00148  前との比=2.02
# EM 強誤差(n=1024)=0.00074  前との比=2.01

OU は平均2・分散0.125へ収束（平均回帰の定常分布）。EM の強誤差は刻み数を4倍（ $\Delta t$ を1/4）にするたびに約1/2になり、 $\sqrt{\Delta t}$ の収束次数0.5が確認できます。

他過程との関係

GBM（幾何ブラウン運動）は $\mu(x)=\mu x,\ \sigma(x)=\sigma x$ の SDE。OU 過程は AR(1)（定常性と独立増分）の連続時間版で、定常・平均回帰する点が共通です。
SDE の解は拡散過程（連続状態のマルコフ過程）で、その分布の発展は連続時間マルコフ連鎖と生成行列の生成行列を連続化したフォッカー＝プランク方程式に従います。

数式の直観的意味

強収束次数が0.5に留まるのは、ブラウン増分の大きさ $\sqrt{\Delta t}$ と、Euler-Maruyama が拡散項を1次近似しかしないことの帰結です。より高次（ミルシュタイン法）にするには伊藤の公式由来の補正項 $\frac12\sigma\sigma'((\Delta B)^2-\Delta t)$ を足す必要があり、ここでも伊藤の公式の二次変分が顔を出します。弱収束が次数1なのは、分布だけ見れば増分の細部が平均化されて効かなくなるためです。

⚠️ よくある誤解

拡散項は $\sigma\Delta t$ ではなく $\sigma\sqrt{\Delta t}\,Z$ 。決定論的オイラーの感覚で $\Delta t$ を掛けると揺れが消えてしまいます。
強収束と弱収束を混同しない。経路を当てたい（パスごとの誤差）なら次数0.5、期待値だけ当てたい（オプション価格など）なら次数1。目的で必要な刻みが変わります。
EM は安定とは限らない。剛性の強い SDE では陰的スキームが必要。刻みが粗いと数値的に発散することがあります。

対応シミュレーション

本文の OU・GBM に加え、ミルシュタイン法との収束次数比較（0.5 対 1.0）を stochastic-processes-study/simulations/ に置きます。