共役事前分布

共役事前分布の概念と代表的な例を、ベイズの定理からの導出付きで解説。ベータ-二項分布・ディリクレ-多項分布・正規-正規分布の共役性を証明し、Pythonでベータ分布の逐次事後更新を実行検証します。

ベイズ推定において、事後分布を解析的に求めるためには、尤度関数と事前分布の組み合わせが重要になります。特定の尤度関数に対して、事後分布が事前分布と同じ関数形になるような事前分布を共役事前分布 (Conjugate Prior) と呼びます。共役事前分布を用いることで、事後分布の計算が非常にシンプルになります。

本記事では、代表的な共役事前分布を天下り的に提示するのではなく、ベイズの定理から実際に事後分布を導出し、なぜ「同じ関数形」に戻ってくるのかを確認します。そのうえで、共役性が破れる注意点と、Pythonによる実行検証を行います。

代表的な共役事前分布

ディリクレ分布 (Dirichlet Distribution)

多項分布のパラメータ(各カテゴリの確率)に対する共役事前分布は、ディリクレ分布 \(\mathcal{D}(\mu|\alpha)\) です。

\[ \mathcal{D}(\mu|\alpha) = \frac{\Gamma(\sum_{j=1}^{k}\alpha_j)}{\prod_{j=1}^k\Gamma(\alpha_j)}\prod_{j=1}^k\mu_j^{\alpha_j-1} \]

ここで、\(\alpha = (\alpha_1, \dots, \alpha_k)\) はディリクレ分布のパラメータであり、\(\Gamma(\cdot)\) はガンマ関数です。ガンマ関数は階乗の概念を実数に拡張したもので、以下のように定義されます。

\[ \Gamma(x) = \int_0^\infty t^{x-1}e^{-t}dt \]

ベータ分布 (Beta Distribution)

二項分布のパラメータ(成功確率)に対する共役事前分布は、ベータ分布 \(p(\mu|a,b)\) です。

\[ p(\mu|a,b) = \frac{1}{B(a,b)}\mu^{a-1}(1-\mu)^{b-1} \]

ここで、\(a, b\) はベータ分布のパラメータであり、\(B(a,b)\) はベータ関数です。ベータ関数は以下のように定義されます。

\[ B(a,b) = \int_0^1\mu^{a-1}(1-\mu)^{b-1}d\mu \]

ガンマ分布 (Gamma Distribution)

正規分布の**分散の逆数(精度パラメータ)**に対する共役事前分布は、ガンマ分布 \(\mathcal{G}(\lambda|\kappa,\xi)\) です。

\[ \mathcal{G}(\lambda|\kappa,\xi) = \frac{\xi^\kappa}{\Gamma(\kappa)}\lambda^{\kappa-1}\exp(-\xi\lambda) \]

ここで、\(\lambda = 1/\sigma^2\) は精度パラメータ、\(\kappa\) は形状パラメータ、\(\xi\) は尺度パラメータ(またはレートパラメータの逆数)と呼ばれます。

補足: 正規分布の平均パラメータに対する共役事前分布は正規分布です。しかし、分散パラメータに対して正規分布を事前分布として用いると、事後分布が正規分布にならないため、共役ではありません。そこで、分散の逆数である精度パラメータを導入し、その共役事前分布としてガンマ分布を用いるのが一般的です。

ガンマ分布の余談

  • \(\kappa=1\) のガンマ分布は指数分布と一致します。
\[ \mathcal{G}(\lambda|1,\xi) = \frac{\xi^1}{\Gamma(1)}\lambda^{1-1}\exp(-\xi\lambda) = \xi\exp(-\xi\lambda) \]

(\(\Gamma(1)=1\) )

  • \(\xi=\frac{1}{2}\) に設定し、\(\nu=2\kappa\) で定義されるパラメータで表したものはカイ二乗分布 \(\chi^2(\lambda|\nu)\) と呼ばれます。
\[ \chi^2(\lambda|\nu) = \mathcal{G}\left(\lambda \middle| \frac{\nu}{2},\frac{1}{2}\right) \]

正規-ガンマ分布 (Normal-Gamma Distribution)

正規分布の平均 \(\mu\) と精度 \(\lambda\) の両方を同時に推定する場合、その同時事前分布として正規-ガンマ分布 \(\mathcal{NG}(\mu,\lambda|\psi,\beta,\kappa,\xi)\) を用いることができます。これは、平均 \(\mu\) が正規分布に従い、精度 \(\lambda\) がガンマ分布に従うという構造を持ちます。

\[ \mathcal{NG}(\mu,\lambda|\psi,\beta,\kappa,\xi) = \mathcal{N}(\mu|\psi,(\beta \lambda)^{-1}) \mathcal{G}(\lambda|\kappa,\xi) \]

なぜ共役性が成り立つのか:ベータ-二項分布の導出

「事後分布が事前分布と同じ関数形になる」という主張は、天下りに信じるものではなく、ベイズの定理を実際に計算すれば確認できます。

コインを \(n\) 回投げて表が \(k\) 回出たとします(二項分布に従う尤度)。表が出る確率を \(\mu\) とすると、尤度関数は

\[ p(D|\mu) = \binom{n}{k}\mu^k(1-\mu)^{n-k} \]

です。事前分布としてベータ分布 \(p(\mu) = \frac{1}{B(a,b)}\mu^{a-1}(1-\mu)^{b-1}\) を置くと、ベイズの定理より事後分布は

\[ p(\mu|D) = \frac{p(D|\mu)p(\mu)}{p(D)} \propto p(D|\mu)p(\mu) \]

となります(\(p(D)\) は \(\mu\) に依存しない正規化定数なので、比例関係の議論では無視できます)。尤度と事前分布を代入すると、

\[ p(\mu|D) \propto \binom{n}{k}\mu^k(1-\mu)^{n-k}\cdot\frac{1}{B(a,b)}\mu^{a-1}(1-\mu)^{b-1} \]

ここで \(\binom{n}{k}\) と \(\frac{1}{B(a,b)}\) はどちらも \(\mu\) に依存しない定数なので比例関係にまとめて吸収でき、\(\mu^k\) と \(\mu^{a-1}\) 、\((1-\mu)^{n-k}\) と \((1-\mu)^{b-1}\) をそれぞれ指数法則でまとめると、

\[ p(\mu|D) \propto \mu^{k+a-1}(1-\mu)^{n-k+b-1} \]

が得られます。この式は、\(\mu^{\alpha-1}(1-\mu)^{\beta-1}\) というベータ分布の核(カーネル)そのものの形をしています。確率分布は \([0,1]\) 上で積分すると1になる必要がありますが、ベータ関数の定義 \(B(\alpha,\beta) = \int_0^1 \mu^{\alpha-1}(1-\mu)^{\beta-1}d\mu\) から、この核を正規化する定数は必ず \(B(\alpha,\beta)\) になることが保証されています。したがって、正規化定数を計算し直すまでもなく、

\[ p(\mu|D) = \frac{1}{B(a+k,\, b+n-k)}\mu^{a+k-1}(1-\mu)^{b+n-k-1} = \text{Beta}(\mu \mid a+k,\, b+n-k) \]

であることが分かります。事前分布が \(\text{Beta}(a,b)\) 、事後分布が \(\text{Beta}(a+k, b+n-k)\) となり、確かに同じ関数形(ベータ分布)に戻ってきていることが導出できました。\(a, b\) はそれぞれ「観測前から仮定していた成功・失敗の疑似カウント(pseudo-count)」と解釈でき、観測データの成功数 \(k\) ・失敗数 \(n-k\) がそこに単純加算される形になっています。

ディリクレ-多項分布への一般化

同じ議論はカテゴリ数 \(K\) の多項分布にそのまま拡張できます。各カテゴリの観測回数を \(k_1,\dots,k_K\) (\(\sum_j k_j = n\) )、事前分布を \(\mathcal{D}(\mu|\alpha)\) とすると、尤度と事前分布の積は

\[ p(\mu|D) \propto \left(\prod_{j=1}^K \mu_j^{k_j}\right)\left(\prod_{j=1}^K \mu_j^{\alpha_j-1}\right) = \prod_{j=1}^K \mu_j^{\alpha_j+k_j-1} \]

となり、これはディリクレ分布の核 \(\prod_j \mu_j^{\alpha_j'-1}\) の形そのものです。よって事後分布は \(\mathcal{D}(\mu \mid \alpha_1+k_1, \dots, \alpha_K+k_K)\) となります。ベータ-二項分布の場合(\(K=2\) )の特殊ケースになっていることが確認できます。

正規分布の平均に対する共役性の導出(分散既知の場合)

離散分布だけでなく、連続分布でも同様の議論が成り立つことを確認します。分散 \(\sigma^2\) が既知で、平均 \(\mu\) のみが未知の正規分布から \(n\) 個のデータ \(x_1,\dots,x_n\) が観測されたとします。事前分布として \(p(\mu) = \mathcal{N}(\mu \mid \mu_0, \tau_0^2)\) を置くと、事後分布は

\[ p(\mu \mid D) \propto \left[\prod_{i=1}^n \exp\left(-\frac{(x_i-\mu)^2}{2\sigma^2}\right)\right]\exp\left(-\frac{(\mu-\mu_0)^2}{2\tau_0^2}\right) \]

指数部分を \(\mu\) について整理します(平方完成)。\(\mu\) の2次の項の係数と1次の項の係数をそれぞれ集めると、

\[ -\frac{1}{2}\left[\left(\frac{n}{\sigma^2}+\frac{1}{\tau_0^2}\right)\mu^2 - 2\left(\frac{n\bar{x}}{\sigma^2}+\frac{\mu_0}{\tau_0^2}\right)\mu\right] + \text{const} \]

(\(\bar{x} = \frac{1}{n}\sum_i x_i\) は標本平均、\(\mu\) を含まない項は正規化定数に吸収されます)。この形はガウス分布の指数部 \(-\frac{1}{2\tau_n^2}(\mu-\mu_n)^2\) を展開したものと係数比較でき、

\[ \frac{1}{\tau_n^2} = \frac{1}{\tau_0^2} + \frac{n}{\sigma^2}, \qquad \mu_n = \tau_n^2\left(\frac{\mu_0}{\tau_0^2} + \frac{n\bar{x}}{\sigma^2}\right) \]

という更新式が導かれます。事後分布は \(\mathcal{N}(\mu \mid \mu_n, \tau_n^2)\) となり、正規分布の平均に対する共役事前分布が正規分布であることが確認できます。精度(分散の逆数)どうしが単純加算される点(\(1/\tau_n^2 = 1/\tau_0^2 + n/\sigma^2\) )は、精度がデータ数に応じて線形に蓄積していく直感とも整合します。

一般的な確率分布と共役事前分布の対応

尤度関数(データ分布)パラメータ共役事前分布
ベルヌーイ分布成功確率 \(\mu\)ベータ分布
二項分布成功確率 \(\mu\)ベータ分布
カテゴリカル分布各カテゴリの確率 \(\mu\)ディリクレ分布
多項分布各カテゴリの確率 \(\mu\)ディリクレ分布
正規分布平均 \(\mu\)正規分布
正規分布精度 \(\lambda\)ガンマ分布
正規分布平均 \(\mu\) , 精度 \(\lambda\)正規-ガンマ分布

共役事前分布

エッジケース・注意点

  • 一様事前分布との関係: \(\text{Beta}(1,1)\) は \([0,1]\) 上の一様分布(Bayes–Laplace の事前分布)に一致します。一方 \(\text{Beta}(0,0)\) (Haldane 事前分布)は積分が発散する非正則(improper)事前分布であり、単独では確率分布として成立しません。ただし観測データに成功・失敗が最低1回ずつ含まれれば、事後分布 \(\text{Beta}(k, n-k)\) は正則になります。
  • 疑似カウントとしての解釈: 上記の導出からも分かる通り、事前分布のパラメータ \(a, b\) は「観測前に仮定していた成功 \(a-1\) 回・失敗 \(b-1\) 回」という疑似的な観測データとみなせます。事前分布が強い(\(a, b\) が大きい)ほど、実データの影響が相対的に小さくなります(後述のPython実行例で数値確認します)。
  • データ数 \(n \to \infty\) での漸近的な挙動: 共役事前分布を使うと事後分布は解析的に更新できますが、\(n\) が大きくなるにつれて事後分布は尤度(=最尤推定量周辺)に支配され、事前分布の影響は消えていきます。これはBernstein–von Mises の定理として知られる一般的な漸近正規性の結果で、正則条件のもとでは事前分布の選択によらず事後分布が真のパラメータ周りの正規分布に漸近的に収束することを保証します。
  • 共役性が使えない場合: 共役事前分布が存在するのは、尤度が指数型分布族に属する場合に限られます。例えばロジスティック回帰の尤度(シグモイド関数を含む)には解析的な共役事前分布が存在せず、Pólya-Gamma 補助変数によるデータ拡大法(Polson, Scott & Windle, 2013)のような特殊な工夫か、後述する変分ベイズ法変分ベイズ法 )やMCMCマルコフ連鎖モンテカルロ法 )による近似推論が必要になります。

Python実行検証

ベータ-二項分布の共役性を、実際にコインを投げるシミュレーションで確認します。真の成功確率 \(p=0.7\) のコインを50回投げ、事前分布 \(\text{Beta}(2,2)\) から事後分布を逐次更新します。

import numpy as np
from scipy import stats

rng = np.random.default_rng(42)

true_p = 0.7
n_total = 50
flips = rng.binomial(1, true_p, size=n_total)

a0, b0 = 2, 2  # 弱情報事前分布 Beta(2,2)
checkpoints = [0, 5, 10, 20, 50]
for n in checkpoints:
    k = flips[:n].sum()
    a, b = a0 + k, b0 + (n - k)
    mean = a / (a + b)
    mode = (a - 1) / (a + b - 2) if a > 1 and b > 1 else float("nan")
    lo, hi = stats.beta.ppf([0.025, 0.975], a, b)
    print(f"n={n:3d}  successes={k:3d}  posterior=Beta({a},{b})  "
          f"mean={mean:.4f}  MAP={mode:.4f}  95%CI=[{lo:.4f}, {hi:.4f}]")

実行結果は以下の通りです。

n=  0  successes=  0  posterior=Beta(2,2)   mean=0.5000  MAP=0.5000  95%CI=[0.0943, 0.9057]
n=  5  successes=  3  posterior=Beta(5,4)   mean=0.5556  MAP=0.5714  95%CI=[0.2449, 0.8430]
n= 10  successes=  5  posterior=Beta(7,7)   mean=0.5000  MAP=0.5000  95%CI=[0.2513, 0.7487]
n= 20  successes= 12  posterior=Beta(14,10) mean=0.5833  MAP=0.5909  95%CI=[0.3854, 0.7681]
n= 50  successes= 32  posterior=Beta(34,20) mean=0.6296  MAP=0.6346  95%CI=[0.4980, 0.7521]

観測数 \(n\) が増えるにつれて、事後平均が真の値 \(p=0.7\) に近づき(\(n=50\) で \(0.6296\) )、95%信用区間の幅も \([0.0943, 0.9057]\) (幅 \(0.8114\) )から \([0.4980, 0.7521]\) (幅 \(0.2541\) )へと単調に縮小していることが確認できます。

次に、事前分布の「強さ」が事後分布に与える影響を確認します。\(n=8\) 回の観測で成功 \(k=6\) 回(標本比率 \(0.75\) )というデータに対し、事前分布の強さを変えて事後平均を比較します。

for (a0s, b0s, label) in [
    (1, 1, "Uniform Beta(1,1)"),
    (2, 2, "Weak Beta(2,2)"),
    (20, 20, "Strong Beta(20,20)"),
]:
    a, b = a0s + 6, b0s + (8 - 6)
    mean = a / (a + b)
    print(f"{label:22s} -> posterior=Beta({a},{b})  mean={mean:.4f}")
Uniform Beta(1,1)      -> posterior=Beta(7,3)   mean=0.7000
Weak Beta(2,2)         -> posterior=Beta(8,4)   mean=0.6667
Strong Beta(20,20)     -> posterior=Beta(26,22)  mean=0.5417

事前分布が強くなる(\(a_0=b_0=20\) )ほど、同じデータに対して事後平均が事前分布の中心 \(0.5\) に強く「引き寄せられる」ことが数値的に確認できます。これは疑似カウント解釈と整合する結果です(\(\text{Beta}(20,20)\) は「20回ずつの成功・失敗を既に観測した」のと同等の情報量を持つため、8回・6成功程度の新規データでは大きく動かされません)。

最後に、正規-正規分布の共役性についても、平均3・分散4の正規分布から30個のデータを生成し、事前分布 \(\mathcal{N}(0, 9)\) から逐次更新した結果を示します。

rng_normal = np.random.default_rng(42)
true_mean, true_var = 3.0, 4.0
n_total = 30
data = rng_normal.normal(true_mean, np.sqrt(true_var), size=n_total)

prior_mean, prior_var = 0.0, 9.0
for n in [0, 5, 10, 30]:
    xs = data[:n]
    if n == 0:
        post_mean, post_var = prior_mean, prior_var
    else:
        prec0 = 1 / prior_var
        precl = n / true_var
        post_var = 1 / (prec0 + precl)
        post_mean = post_var * (prior_mean * prec0 + xs.sum() / true_var)
    print(f"n={n:2d}  posterior mean={post_mean:.4f}  posterior var={post_var:.4f}")
n= 0  posterior mean=0.0000  posterior var=9.0000
n= 5  posterior mean=2.3895  posterior var=0.7347
n=10  posterior mean=2.2298  posterior var=0.3830
n=30  posterior mean=2.9893  posterior var=0.1314

事後分散が \(9.0 \to 0.1314\) と単調に減少し、精度(分散の逆数)が線形加算されるという導出結果 \(1/\tau_n^2 = 1/\tau_0^2 + n/\sigma^2\) が数値的にも裏付けられています。事後平均は真の平均 \(3.0\) の近傍でサンプルのばらつきに応じて上下しますが(\(n=30\) で \(2.9893\) )、これは有限サンプルにおける正規分布特有の揺らぎであり、事後分散が確実に縮小していく点こそが共役更新の理論的に保証された性質です。

上記のベータ分布の逐次更新を可視化すると、以下のようになります。

Beta事後分布の逐次更新(真の成功確率p=0.7、事前分布Beta(2,2))。観測数が増えるにつれて分布が真の値の周りに集中していく

発展:共役性の限界と一般化ベイズ推論

共役事前分布は計算上非常に便利ですが、その適用範囲は尤度が指数型分布族に属する場合に限られます。近年は、尤度そのものの計算が困難(intractable)なモデルに対しても妥当な事後分布を構成する一般化ベイズ推論 (Generalized Bayesian Inference) の研究が進んでいます。例えば Matsubara ら (2024) は、正規化定数の計算が困難な離散尤度モデルに対し、離散フィッシャーダイバージェンスに基づく一般化事後分布を提案し、事後一致性・漸近正規性を保証する条件を示しています。これは、共役事前分布が使えない現実的なモデル(グラフィカルモデル、空間統計モデルなど)に対しても、ベイズ的な不確実性定量化を可能にする方向性の一つです。

参考