周辺尤度

ベイズ推定の周辺尤度による表現

復習すると、ベイズ予測分布は

\begin{align*}
p(X|\{x_n\})
 &= \frac{p(X,\{x_n\}_{n=1}^N)}{p(\{x_n\}_{n=1}^N)} \\
 &\approx \frac{1}{p(\{x_n\}_{n=1}^N)} \int p(X,\{x_n\}_{n=1}^N,\theta) d\theta \\
 &= \frac{1}{p(\{x_n\}_{n=1}^N)} \int p(X,\{x_n\}_{n=1}^N|\theta) p(\theta) d\theta \\
 &= \frac{1}{p(\{x_n\}_{n=1}^N)} \int p(X|\{x_n\}_{n=1}^N,\theta) p(\{x_n\}_{n=1}^N|\theta) p(\theta) d\theta \\
 &= \int p(X|\theta) \frac{p(\{x_n\}_{n=1}^N|\theta) p(\theta)}{p(\{x_n\}_{n=1}^N)} d\theta \\
 &= \int p(X|\theta) p(\theta|\{x_n\}_{n=1}^N) d\theta
\end{align*}

でした。事後分布は

\begin{align}
p(\theta|\{x_n\}_{n=1}^N) = \frac{p(\{x_n\}_{n=1}^N|\theta)p(\theta)}{p(\{x_n\}_{n=1}^N)}
\end{align}

で、右辺分母の

\begin{align}
p(\{x_n\}_{n=1}^N)
 &=\int p(\{x_n\}_{n=1}^N,\theta) d\theta \\
 &=\int p(\{x_n\}_{n=1}^N|\theta) p(\theta) d\theta
\end{align}

は正規化定数なので、事後分布の関数形には関わらない、としてそれ以上の説明をしていませんでした。この正規化定数は**周辺尤度（Marginal Likelihood）**と言います。

ここまでベイズ推定の理解のために、観測値$\{x_n\}_{n=1}^N$と予測したい確率変数$X$を区別して記載していましたが、観測値も確率変数も$x_n$で表記し、観測データを$\mathcal{D}_N=\{x_n\}_{n=1}^N$、予測したい確率変数を$x_{N+1}$とするとベイズ予測分布は

\begin{align}
p(x_{N+1}|\mathcal{D}_N)
 = \frac{p(x_{N+1},\mathcal{D}_N)}{p(\mathcal{D}_N)}
 = \frac{p(\mathcal{D}_{N+1})}{p(\mathcal{D}_N)}
\end{align}

と周辺尤度で書けてしまいます！ 実際には$x_{N+1}$は未確定な確率変数で、$\mathcal{D}_N=\{x_n\}_{n=1}^N$だけが確定した観測値なので、分子の $p(\mathcal{D}_{N+1}) = p(x_{N+1},\mathcal{D}_N)$ は周辺尤度ではなく$x_{N+1}$についての関数で、**「$x_{N+1}$という値に対して『$x_{N+1}$も含めた観測値$\mathcal{D}_{N+1}=\{x_n\}_{n=1}^{N+1}$に対する周辺尤度値』を取る関数」**です。

すなわち周辺尤度 $p(\mathcal{D}_N)$ は、単なる正規化定数ではなく、予測分布の情報を持っている「ベイズ推定の本質」と言っても過言ではないとても重要な量です。

階層ベイズと周辺尤度

$\psi$をモデルを指すインデックス（あるいはモデルのハイパーパラメータ）として、モデルパラメータ$\theta$の事前分布を$p(\theta|\psi)$とします。ここで$\psi$が違えば全然違うモデルかもしれないし、$\theta$の定義域も違うかもしれないのですが、同じ$\theta$という表記を用いてしまいます。($\theta$について積分消去するので問題ない)
$\psi$も含めた階層ベイズによる予測分布は

\begin{align}
p(x_{N+1}|\mathcal{D}_N)
 &= \frac{p(\mathcal{D}_{N+1})}{p(\mathcal{D}_N)} \\
 &\approx \frac{1}{p(\mathcal{D}_N)} \int \int p(\mathcal{D}_{N+1},\theta,\psi)  d\theta d\psi \\
 &= \frac{1}{p(\mathcal{D}_N)} \int \int p(\mathcal{D}_{N+1}|\theta) p(\theta|\psi) p(\psi)  d\theta d\psi \\
 &= \frac{1}{p(\mathcal{D}_N)} \int \left( \int p(\mathcal{D}_{N+1}|\theta) p(\theta|\psi) d\theta \right) p(\psi) d\psi \\
 &= \frac{1}{p(\mathcal{D}_N)} \int p(\mathcal{D}_{N+1}|\psi) p(\psi) d\psi \\
 &= \frac{1}{p(\mathcal{D}_N)} \int p(x_{N+1}|\mathcal{D}_N,\psi) p(\mathcal{D}_N|\psi) p(\psi) d\psi \\
 &= \int p(x_{N+1}|\mathcal{D}_N,\psi) \frac{p(\mathcal{D}_N|\psi) p(\psi)}{p(\mathcal{D}_N)} d\psi \\
 &= \int p(x_{N+1}|\mathcal{D}_N,\psi) p(\psi|\mathcal{D}_N) d\psi
\tag{1}
\end{align}

で、

\begin{align}
p(\mathcal{D}_N|\psi) &= \int p(\mathcal{D}_N|\theta) p(\theta|\psi) d\theta
\end{align}

は$\psi$が指すモデルの、モデルパラメータ$\theta$について周辺化した周辺尤度であり、

\begin{align}
p(\psi|\mathcal{D}_N) &= \frac{p(\mathcal{D}_N|\psi) p(\psi)}{p(\mathcal{D}_N)}
\tag{2}
\end{align}

は$\psi$の事後分布です。$\theta$の事後分布では分子は尤度と事前分布の積でしたが、$\psi$の事後分布の分子は周辺尤度と事前分布の積になっています。分母は正規化定数（$\psi$も含めて周辺化した周辺尤度）で、関数形は分子で決まります。$p(x_{N+1}|\mathcal{D}_N,\psi)$は$\psi$の下でのベイズ予測分布で

\begin{align}
p(x_{N+1}|\mathcal{D}_N,\psi)
 &= \int p(x_{N+1},\theta|\mathcal{D}_N,\psi) d\theta \\
 &= \int p(x_{N+1}|\theta) p(\theta|\mathcal{D}_N,\psi) d\theta
\end{align}

です。(1)に代入すると

\begin{align}
p(x_{N+1}|\mathcal{D}_N)
 &= \int \left( \int p(x_{N+1}|\theta) p(\theta|\mathcal{D}_N,\psi) d\theta \right) p(\psi|\mathcal{D}_N) d\psi
\end{align}

と書けます。カッコ内は$\psi$の下でのベイズ予測分布で「モデル$p(x_{N+1}|\theta)$の、$\psi$の下での$\theta$の事後分布$p(\theta|\mathcal{D}_N,\psi)$による期待値」であり、**階層ベイズ全体の予測分布は「$\psi$の下でのベイズ予測分布の、$\psi$の事後分布$p(\psi|\mathcal{D}_N)$による期待値」**です。

モデルパラメータ$\theta$についてもハイパーパラメータ$\psi$についても、統一的に扱うことができていて美しいです。

経験ベイズ法は「モデル選択」、階層ベイズは「アンサンブル学習」

$\psi$の事後分布(2)による期待値計算において、$\psi$の事前分布が無情報であれば、モデル$\psi$に対する重みは周辺尤度$p(\mathcal{D}_N|\psi)$になります。このことから、周辺尤度はモデル（と$\theta$の事前分布）の良さを量るベイズ的な指標と言え、証拠（エビデンス）(Evidence)とも呼ばれます。
もし$\psi$を一つ選ぶとするなら$p(\mathcal{D}_N|\psi)$が最大の一点を選ぶことがリーズナブルでしょう。周辺尤度を$\psi$について最大化することを**経験ベイズ法（第二種最尤推定）**と呼びます。ここでは詳しい説明を省きますが、モデル選択規準の一つであるBayesian Information Criterion（BIC）も、事後分布をモード（MAP推定値）周りで正規分布で近似（ラプラス近似）して算出した対数周辺尤度です。最尤推定ではモデル選択に最尤推定値の尤度を用いることはできず、情報量規準やCross-Validationを用いる必要がありますが、ベイズ推定では周辺尤度はそのままモデル（ハイパーパラメータ）選択に用いることができます。（最適性を言うにはいろいろ条件が要ります。）
階層ベイズはモデルを選択するのではなくさらに周辺化します。**階層ベイズは、各モデルの予測分布の周辺尤度（×事前分布）による重み付き和の「アンサンブル学習」**と解釈することもできます。

次回は「ベイズ推定７：EMアルゴリズム」です。