事前・尤度・事後・周辺尤度

🎓 レベル：基礎　|　重要度：A（必須）

📎 前提：頻度論とベイズの違い　|　数理：事前分布・事後分布・ベイズ更新（統計）・最尤法・モーメント法（推定量の作り方と最尤推定量の漸近論）（統計）

要点（BLUF）

ベイズの定理 $P(\theta\mid D)=\dfrac{P(D\mid\theta)\,P(\theta)}{P(D)}$ は、尤度・事前・事後・周辺尤度の4部品でできています。
尤度 $P(D\mid\theta)$ は「データ生成モデルを $\theta$ の関数として見たもの」。 $\theta$ について積分しても1にならず、確率分布ではありません。
周辺尤度 $P(D)$ は事後を1に正規化する定数であると同時に、「そのモデルがデータをどれだけ説明できるか」を表し、モデル比較（第8章）の主役です。計算が難しいことが、MCMC（第4章）が要る理由でもあります。

1. 尤度 P(D|θ)：θ の関数として見る

データ生成の確率モデルを書き、そこに観測データを代入して $\theta$ の関数とみなしたものが尤度です。コインを $n$ 回投げて $k$ 回表なら、

P(D\mid\theta) = \theta^{k}(1-\theta)^{n-k}

ここがつまずきやすい点： $\theta$ を固定すれば「データの確率分布」ですが、データを固定して $\theta$ を動かすと、もはや確率分布ではありません（ $\theta$ について積分しても1になりません）。だから尤度は「高さの比」だけに意味があり、絶対値そのものには意味がありません。

2. 事前 P(θ)：データを見る前の信念

$\theta$ に対する観測前の信念を確率分布で表したものが事前分布です。これは確率分布（積分して1）です。「どんな値もありうる」なら一様 $\mathrm{Beta}(1,1)$ 、「0.5 付近が濃いだろう」なら $\mathrm{Beta}(2,2)$ のように、信念の強さを形で表します（事前の選び方はベイズ更新と逐次推論・無情報事前と弱情報事前で深掘り）。

3. 事後 ∝ 尤度 × 事前

更新後の信念が事後分布です。形は「尤度 × 事前」で決まります。

P(\theta\mid D)\;\propto\;P(D\mid\theta)\,P(\theta)

flowchart LR
  PR["事前 P(θ)"] --> MUL["かけ算"]
  LK["尤度 P(D|θ)"] --> MUL
  MUL --> U["尤度×事前（未正規化）"]
  U --> DIV["÷ 周辺尤度 P(D)"]
  DIV --> PO["事後 P(θ|D)"]

4. 周辺尤度 P(D)：正規化定数にして、モデルの説明力

「 $\propto$ 」を「 $=$ 」に戻す正規化定数が周辺尤度です。 $\theta$ を積分で消し去ります。

P(D) = \int P(D\mid\theta)\,P(\theta)\,d\theta

これは「 $\theta$ が何であれ、このモデルからこのデータが出る平均的な確率」＝モデルのデータ説明力です。事後分布の形を見るだけなら無視できますが、複数モデルを比べるとき（第8章のベイズファクター）に主役になります。一般には高次元の積分で解析的に解けない——この困難を回避・近似するのが共役（第2章）や MCMC（第4章）です。

5. コードで「事前 × 尤度 → 事後」を組み立てる

数値的に事後を作り、周辺尤度も求めてみます。共役の公式を使わず、定義どおり積で作って正規化する素朴な実装です。

import numpy as np
from scipy import stats
from scipy.integrate import trapezoid  # numpy 2.0+ で np.trapz は廃止のため
import matplotlib.pyplot as plt
import japanize_matplotlib  # 日本語ラベル用

# コイン：n 回中 k 回表
n, k = 20, 14
theta = np.linspace(0, 1, 500)

# --- 尤度 P(D|θ) = θ^k (1-θ)^(n-k)：θ の関数（確率分布ではない） ---
likelihood = theta**k * (1 - theta)**(n - k)

# --- 事前 P(θ)：弱情報事前 Beta(2,2)（0.5 付近を少しだけ信じる） ---
prior = stats.beta(2, 2).pdf(theta)

# --- 事後 ∝ 尤度 × 事前。台形則で積分し正規化 ---
unnorm = likelihood * prior
evidence = trapezoid(unnorm, theta)     # 周辺尤度 P(D) = ∫ 尤度×事前 dθ
posterior = unnorm / evidence

print(f"周辺尤度 P(D) = {evidence:.4e}")
print(f"事後の最頻値(MAP) ≈ {theta[np.argmax(posterior)]:.3f}")

# 3つを重ねて描く（尤度は形を見たいので事後の高さにスケール）
plt.figure(figsize=(7, 4))
plt.plot(theta, prior, ls="--", label="事前 Beta(2,2)")
plt.plot(theta, likelihood / likelihood.max() * posterior.max(),
         ls=":", label="尤度（形のみ・規格化）")
plt.plot(theta, posterior, lw=2, label="事後 P(θ|D)")
plt.xlabel("θ（表の確率）"); plt.ylabel("密度")
plt.title("事前 × 尤度 → 事後")
plt.legend(); plt.tight_layout(); plt.show()

出力はおよそ次のようになります。

周辺尤度 P(D) = 1.5296e-06
事後の最頻値(MAP) ≈ 0.681

出力の意味：事後分布（実線）は、事前（破線）と尤度（点線）のあいだに位置します。尤度の最頻 $0.70$ より少し 0.5 寄りに引かれているのは、弱情報事前 $\mathrm{Beta}(2,2)$ が「0.5 付近だろう」と軽く主張しているからです。周辺尤度 $P(D)\approx 1.5\times10^{-6}$ は、この $\theta$ の積分でデータが出る平均確率——それ自体は小さな数ですが、別のモデル（別の事前や別の生成過程）と比べたときの比に意味が出ます（第8章）。

⚠️ よくある誤解

「尤度は確率分布」ではない： $P(D\mid\theta)$ を $\theta$ で積分しても1になりません。確率分布なのは事前・事後・データ側 $P(D\mid\theta)$ （ $D$ について）です。
「周辺尤度は無視してよい」わけではない：事後の形を見るだけなら定数として無視できますが、モデル比較では $P(D)$ そのものが評価指標になります。
「事後は必ず事前と尤度の間にある」とは限らない：1次元の素直な例では間に来やすいですが、多次元や非共役では直観に反する形になることもあります。だからこそ数値計算（MCMC）で確かめます。