ベイズ推定と最尤推定の比較

データが何かの分布に従うと仮定したら、パラメーターを求める為に何か計算をします。その時の手法として、最尤推定とベイズ推定¹があります。それぞれの手法を比較し、注意するべき事を見つけます。
分散が既知の場合の、正規分布の平均値の推定を通して検証してみます。

ベイズ推定
ベイズ推定の計算
最尤推定の計算
python による実装
予測値の精度
まとめ

ベイズ推定

ベイズ推定の手順を説明します。あるデータ$ (Y,X) $があり、 $Y=\{y_1 , \cdots , y_N \} $は, 確率分布 $ y_n =p(x_n| \theta )$に従うとします。この時、良い感じの $ \theta $ を決定すれば, $X$から$Y $が予想できる²という訳です。上手く$ p(\mu )$を指定してやって、$ p( \mu |x ) $が良く知っている確率分布になるようにしておくのが常套手段の1つです。$ p( \mu )$が指定されていると、ベイズの定理から、 $p(\mu |x ) $が求められます。
$$\begin{eqnarray}
p(\theta |x ) = \frac{ p(x| \theta ) p(\theta)}{p(x)}
\end{eqnarray}$$
この記事の中では、パラメーターの事後分布を求める事をベイズ推定と呼びます。
右辺において、$ \theta $に注目するのであれば、$ p(x) $は定数であることに注意しましょう。

ベイズ推定の計算

同じ計算をするので、多次元正規分布を考えて、データの次元を1次元にして結果を使います。分散行列$ \Sigma $は既知としましょう。 $ y_n = \mathcal{N} (x_n |\mu , \Sigma ) $の場合を考えます。また、$ \mu $ も多次元正規分布に従うとします。
$$\begin{eqnarray}
p(\mu ) \sim \mathcal{N}(\mu _0 , \Sigma _0 )
\end{eqnarray}$$
$\mu_0 , \Sigma_0 $は人間が適当に与えてやります。今回はデータを無尽蔵に増やせるので、どんな値を与えても³上手くいきます。
今回の場合だと、$p(\mu |X) $ は多次元正規分布になります。
$ p(\mu |X ) = \mathcal{N} ( \mu| \bar{\mu} , \Sigma _{\mu} ) $と置いて、$ \mu $について展開しておきます。
$$\begin{eqnarray}
\log p(\mu |X ) = -\frac{1}{2} \left[ \mu ^{T} \Sigma _{\mu} ^{-1} \mu -2 \mu ^{T} \Sigma _{\mu} ^{-1} \bar{\mu} \right] + {\rm const}
\end{eqnarray}$$
計算しやすいように$ \log $を取ったベイズの定理を考えます。
$$\begin{eqnarray}
\log p(\mu |X ) = \log p(X| \mu ) + \log p(\mu )-\log p(X)
\end{eqnarray}$$
上で定義した確率分布を代入して、$\mu $についてまとめます。
$$\begin{eqnarray}
\log p(\mu |X ) = -\frac{1}{2}\left[\mu ^{T} ( N\Sigma ^{-1} + \Sigma _{\mu _0} ^{-1} )\mu -2\mu ^{T}( \Sigma ^{-1} \sum_{n}x_n + \Sigma _{\mu _0} ^{-1} \mu _0 ) \right] + { \rm const }
\end{eqnarray}$$
上の式と比べる事で、$ \mu $が従う確率分布のパラメーターを求める事が出来ます。
$$\begin{eqnarray}
\Sigma _{\mu} ^{-1} &=& N\Sigma ^{-1} + \Sigma _{\mu _0} ^{-1}\\
\Sigma _{\mu} ^{-1} \bar{\mu} &=& \Sigma ^{-1} \sum_{n}x_n + \Sigma _{\mu _0} ^{-1} \mu _0 \\
\bar{\mu} &=& \left( N\Sigma ^{-1} + \Sigma _{\mu _0} ^{-1} \right) ^{-1} \left( \Sigma ^{-1} \sum_{n}x_n + \Sigma _{\mu _0} ^{-1} \mu _0 \right)
\end{eqnarray}$$

最尤推定の計算

正規分布の最尤推定については、以下の記事で解説しています。

Fisher情報量(Fisher情報行列)

Fisher 情報量(Fisher 情報行列)の定義と、その役割について解説します。統計量の推定をした時、推定量の分散の大きさについての情報を与えてくれるのがFisher 情報量(Fisher 情報行列)です。

結果的に、正規分布の平均値の最尤推定量は$ \bar{\mu} =\frac{1}{N} \sum x_n $です。

python による実装

pythonで適当な正規分布を生成し、最尤推定、ベイズ推定それぞれで平均値を推定します。実験に使うコードを載せます。

import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns

mu = 10
sig =5
Y = np.random.normal(loc=mu, scale=sig, size=100)
ys=[]
mu_ml=[]
fis=[]
mu_by=[]
lam =[]
for y in Y:
    ys.append(y)
    mu_ml.append(np.sum(ys)/len(ys) ) #最尤推定
    fis=np.append(fis, sig**2 /len(ys) )
　　lam_by = len(ys)/sig+ 1/sig_zero
    lam = np.append(lam,lam_by)
    mu_by=np.append( mu_by,(np.sum(ys)/sig +mu_zero/sig_zero )/lam_by ) #ベイズ推定

plt.plot(mu_ml,label="Likelyhood")
plt.plot(mu_by, label="Byes")
plt.plot(10*np.ones(len(mu_ml)),c="r",lw=0.5)
plt.legend()
plt.ylim(5,15)
plt.title("Byes vs Lilelyhood")

縦軸に推定された平均値、横軸に更新の回数を取ってグラフを描きました。
結果を見ると、ベイズ推定と最尤推定は殆ど同じ動きをしながらパラメーターを推定しています。

予測値の精度

予測の精度はどのくらいあるのでしょうか。
ベイズ推定なら$ \mu $が従う正規分布の分散、最尤推定ならフィッシャー情報量の逆数が予測の精度の目安になります。フィッシャー情報量を知らない方は以下の記事をどうぞ。

Fisher情報量(Fisher情報行列)

正規分布の場合、平均値の最尤推定量の⁴分散は以下の式で下から抑えられています。
$$\begin{eqnarray}
{ \rm var }[\mu] \geq \sigma ^2 /N
\end{eqnarray} $$
この2つの目安を予測値の信頼度として、予測値のグラフに描き入れてみます。

最尤推定の方が信頼度が小さく、ベイズ推定は最尤推定より精度が高い、という結果になりました。
ベイズ推定の場合は、$ \mu $を計算しやすいように色々仮定を置いていました。その仮定の下では自信を持って予測しているという事です。
このことは、上手く現実を再現しているとは限らない事に注意しましょう。実際、繰り返しが20回未満では全く違う平均値を予測していますが、ベイズ推定の予測精度は高いことになっています。⁵ついでに書いておくと、事後分布の分散は、既知の分散$\times $ データ数と、人が設定した分散で決まるので、データ数が少ない場合は非常に恣意性が高くなります。
一方、最尤推定の誤差は、Fisher 情報量 (既知の分散/データ数)で下から抑えられています。更新をいくら続けても誤差が0になることはありません。その為、ベイズ推定より予測の精度が小さいという結果になっています。
今回の場合は、不偏推定量と最尤推定量が一致しているので、最尤推定を使っておけば良いということになります。
一般の場合では、最尤推定量が不偏推定量になるとは限らないのと、最尤推定の計算が出来るとは限らりません。ベイズ推定は、計算のしやすい仮定を置いたり、近似手法を使って計算を進める事が出来たりするので、難しい分布が出てくる問題にこそ、ベイズ推定が使えるかもしれません。