マルコフ連鎖と定常分布｜シミュレーション

🎓 レベル：標準　|　重要度：A（必須）

📎 関連：メトロポリス・ヘイスティングス法　|　収束診断と実装の注意

要点（BLUF）

マルコフ連鎖：次の状態が現在の状態だけで決まり、過去の経路によらない（マルコフ性）確率過程。遷移は遷移行列 $P$ で表す。
定常分布 $\pi$ ：遷移を繰り返しても変わらない分布 $\pi P = \pi$ 。既約・非周期な連鎖は、初期状態によらず唯一の定常分布へ収束します。
3状態連鎖の定常分布を、固有ベクトル・べき乗法・長時間頻度の3通りで計算し、すべて $[0.457, 0.282, 0.261]$ に一致することを確かめます。

1. マルコフ連鎖

状態列 $X_0, X_1, X_2, \dots$ がマルコフ性を持つとは、

P(X_{t+1} = j \mid X_t = i, X_{t-1}, \dots, X_0) = P(X_{t+1} = j \mid X_t = i) = P_{ij}

「次は今だけで決まり、過去の履歴は要らない」こと。遷移確率を並べた行列 $P$ （各行の和が1）が連鎖を定義します。分布 $\mu_t$ （行ベクトル）の1ステップ後は $\mu_{t+1} = \mu_t P$ です。

2. 定常分布

分布 $\pi$ が定常分布とは、遷移しても変わらないこと：

\pi P = \pi, \qquad \sum_i \pi_i = 1,\ \pi_i \ge 0

$\pi$ は $P^\top$ の固有値1に対する左固有ベクトル（を正規化したもの）です。連鎖が

既約（irreducible）：どの状態からどの状態へも有限ステップで到達可能
非周期（aperiodic）：周期的なループに嵌らない

を満たせば、定常分布は唯一存在し、任意の初期分布が $t\to\infty$ で $\pi$ に収束します（エルゴード定理）。MCMC はこの事実を逆用し、「サンプルしたい分布 $\pi$ を定常分布に持つ連鎖」を設計します（メトロポリス・ヘイスティングス法）。

flowchart LR
    A["初期分布 μ0"] -->|"xP"| B["μ1"]
    B -->|"xP"| C["μ2"]
    C -->|"… 収束 …"| D["定常分布 π （πP=π）"]

3. 3通りで定常分布を求める

3状態の遷移行列で、定常分布を (1) 固有ベクトル、(2) べき乗法（ $\mu_0 P^{100}$ ）、(3) 連鎖を実際に歩かせた長時間頻度——の3通りで求め、一致を確かめます。

import numpy as np

P = np.array([[0.7, 0.2, 0.1],
              [0.3, 0.4, 0.3],
              [0.2, 0.3, 0.5]])

# (1) 固有ベクトル（P^T の固有値1）
vals, vecs = np.linalg.eig(P.T)
idx = np.argmin(np.abs(vals - 1))
pi = np.real(vecs[:, idx]); pi = pi / pi.sum()
print(f"定常分布（固有ベクトル）= {np.round(pi, 4)}")

# (2) べき乗法：初期分布に P を100回かける
v = np.array([1.0, 0.0, 0.0])
for _ in range(100):
    v = v @ P
print(f"100ステップ後の分布     = {np.round(v, 4)}")

# (3) 連鎖を歩かせた長時間頻度
rng = np.random.default_rng(40)
state = 0; counts = np.zeros(3); N = 2_000_000
for _ in range(N):
    state = rng.choice(3, p=P[state]); counts[state] += 1
print(f"長時間頻度（シミュ）    = {np.round(counts/N, 4)}")

出力：

定常分布（固有ベクトル）= [0.4565 0.2826 0.2609]
100ステップ後の分布     = [0.4565 0.2826 0.2609]
長時間頻度（シミュ）    = [0.4571 0.2822 0.2607]

出力の意味：3通りの方法がすべて $[0.457, 0.282, 0.261]$ に一致。とくに重要なのが (3)——状態0からスタートしても、連鎖を200万歩歩かせると、各状態を訪れた頻度が定常分布に一致します。これが MCMC の核心原理：「連鎖を長く歩かせれば、訪問頻度が定常分布のサンプルになる」。(2) のべき乗法は初期状態(1,0,0)が100歩で完全に定常へ収束していることを示し、この連鎖の収束が速いことも分かります。

4. MCMC への橋渡し

ここでは遷移行列 $P$ が与えられ、そこから定常分布 $\pi$ を求めました。MCMC はこの問題を逆に解きます——「サンプルしたい $\pi$ が定常分布になるように $P$ を設計する」。設計の指針が詳細釣り合い（メトロポリス・ヘイスティングス法）で、これを満たせば $\pi P = \pi$ が自動的に成り立ちます。連続状態空間でも考え方は同じで、遷移行列が遷移核（カーネル）に変わるだけです。

数式の直観的意味

定常分布は「確率の流れが釣り合った状態」です。 $\pi P = \pi$ は「各状態に流れ込む確率と流れ出る確率が等しく、分布が動かない」こと。既約性は「全状態が繋がっている（孤立した島がない）」、非周期性は「決まったリズムで循環しない」を保証し、この2つがあれば連鎖は初期状態を忘れて唯一の平衡へ落ち着く。長時間頻度が $\pi$ に一致するのはエルゴード性——時間平均（1本の長い軌道での頻度）が空間平均（定常分布）に等しい、という性質で、これがあるから「1本の連鎖を長く回す」MCMC が正当化されます。収束の速さは $P$ の第2固有値の大きさで決まり、これが収束診断（収束診断と実装の注意）の理論的背景です。

⚠️ よくある誤解・落とし穴

「定常分布はいつも一意」ではない：既約でないと複数の定常分布があり得ます（孤立した部分集合ごとに）。
「すぐ収束する」ではない：収束の速さは連鎖次第（第2固有値）。遅い連鎖はバーンインを長く取る必要があります。
「マルコフ性は記憶ゼロ」ではない：「1ステップ前までは使える」。状態に過去の要約を含めれば高次の依存も表せます。
「周期的でも頻度は定常分布」ではない：周期連鎖は分布が振動して収束しません。非周期性が必須。
「定常分布＝初期分布」ではない：初期分布は任意。収束した先が定常分布で、初期の影響（過渡）はバーンインで捨てます。

対応シミュレーション参照

本文の3状態連鎖の定常分布（固有ベクトル／べき乗法／長時間頻度 default_rng(40)）。