ベルヌーイ分布のパラメータを最尤推定する

Posted: , Category: 機械学習 , 統計学

機械学習のパラメータ推定で頻繁に登場する、ベルヌーイ分布(Bernoulli distribution)です。

よくコイン投げの思考に例えられるこの分布ですが、基本的な確率分布としてよく機械学習や統計学のテキストで登場します。一方、データが得られた時に、尤もらしいパラメータを推定する最尤推定法も、非常によく登場しますよね。

今回は、機械学習や統計分野で最も基本的なベルヌーイ 分布と最尤推定という具体例を用いて、データが得られた際にモデルをどのように学習するかについて、学んでいきます。

また実装はほとんどないですが、最尤推定が実際にどのようなことをしているのか理解しやすいように、Pythonを用いてビジュアライズしているので、ぜひ軽く見通して、データからモデルを学習するとはどのようなことかを理解してみてください。

本記事の内容
  • ベルヌーイ 分布の最尤推定する方法
  • ベルヌーイ 分布の例を通して、尤度関数と対数尤度関数の形状を確かめる

問題設定

今回はベルヌーイ 分布を用いて最尤推定をやっていくので、コイン投げの試行で問題設定をしてみます。

ベルヌーイ 分布の最尤推定における問題設定

表(1)と裏(0)がどの割合で出るかわからない、謎のコインがあるとして、この表が出る確率$\mu$を知りたい。

今、このコインを$N$回振ったことで、得られた$mathcal{D} = \{ x_1, x_2, \dots, x_N\}$のデータを用いて、このパラメータ$\mu$を最尤推定することを考える。

ベルヌーイ 分布のパラメータを最尤推定する

では、早速ベルヌーイ 分布のパラメータを最尤推定していきましょう。

ベルヌーイ 分布の確率密度関数は、パラメータ$\mu$を用いてこのように表されていました。

ベルヌーイ分布の定義
\begin{equation}
p(x | \mu) = μ^x(1-μ)^{1-x}
\end{equation}

今、手元にある観測データ$mathcal{D} = \{ x_1, x_2, \dots, x_N\}$が、(1)のベルヌーイ 分布に従って独立に$N$回の試行によって得られたと仮定すると、次に尤度関数を考えることができます。

\begin{equation}
\mathcal{L}(\mu) = \prod_{n=1}^{N} μ^{x_n}(1-μ)^{1-{x_n}}
\end{equation}

ここで、尤度関数はLikelihodd Functionと呼ぶため(1)式は先頭のLをとって、$\mathcal{L}(\mu)$と書きました。

最尤推定の考え方では、(2)式を最大化するパラメータ$\mu$が最尤推定値となります。

ここで、(2)式は積形をしており、扱いにくいので、(2)式に対数を取った関数の最大値問題を考えます。

よく尤度関数の最大化するパラメータの問題を考える際に、対数を取った関数を考えますが、これは正の値をとる関数に、対数を取ってもその大小関係は変わらないという、自然対数に関する性質があるためです。

(2)式に対数を取った対数尤度関数を考えてみます。

\begin{equation}
\begin{split}
ln \mathcal{L}(\mu) = \sum_{n=1}^{N}
\{x_n ln μ + (1-{x_n}) ln (1-μ)\}
\end{split}
\end{equation}

対数尤度関数はこのようになりました。

この(3)式を$\mu$で偏微分した導関数を0とすることで、最尤推定値$\mu_{ML}$を得ることができます。

\begin{equation}
\begin{split}
\frac{\partial  ln \mathcal{L}(\mu)}{\partial \mu}  &= 

\sum_{n=1}^{N} 
\{ 
x_n \frac{1}{\mu} - (1-x_n) \frac{1}{(1 -\mu)} 
\}
= 0 \\
\sum_{n=1}^{N} x_n
\biggl \{
\frac{1}{\mu(1 - \mu)} 
\biggr  \}
&= 
\sum_{n=1}^{N}
\biggl \{
\frac{1}{(1 - \mu)} 
\biggr  \}  = \frac{N}{1-\mu} \\
\therefore \quad \mu & = \frac{1}{N}  \sum_{n=1}^{N}x_n

\end{split}
\end{equation}

となり、最尤推定値

\begin{equation}
\mu_{ML} = \frac{1}{N}  \sum_{n=1}^{N}x_n
\end{equation}

を得ることができました。$x_n$は表が出た時に1になるので、表に出た回数を$N$回の試行のうち、表が出た回数を$m$回とすると、(5)式は次のように書き換えることができる

\begin{equation}
\mu_{ML} = \frac{m}{N} 
\end{equation}

つまり、N回投げた時に表が出る確率を$\mu$の推定値に使うのが、最尤推定法の結論ということになった。

補足: 尤度関数を対数尤度関数に置き換えてよい理由

今回、最尤推定するにあたって、(2)式の尤度関数(Likelihood関数)を(3)式の対数尤度関数にした上で、パラメータ$\mu$で偏微分することで最適なパラメータ$\mu_{ML}$を見つけました。

しかし、実際に尤度関数を対数尤度関数に置き換えても大丈夫なのかについては、いまいち腑に落ちない人もいるかと思います。そういった人向けに、尤度関数と対数尤度関数を可視化してみます。

問題設定

簡単のため、問題を具体化します。今回、10回コインを投げて、5回表が出た場合だと、この場合、(6)式ですでに解析的に解は求められているので、最尤推定値$\mu_{ML} = 0.5$になるはずです。

それをグラフで確認してみましょう。

まず(2)の尤度関数をPythonでグラフ化するとこのようになります。

グラフ化には下記のコードを実装します。

import numpy as np
import matplotlib.pyplot as plt

plt.figure(figsize=(10, 10), dpi=50)
plt.xlim([0, 1.0])
plt.ylim([0, 1.0e-3])

mu = np.linspace(0, 1.0, 1000)

likelihood = (mu ** 5) * ((1 - mu) ** 5)
ln_likelihood = np.log(likelihood)

plt.plot(mu, likelihood, label="$\mathcal{L}$", color="red")
plt.vlines(0.5, 0, 5.0e-3, linestyle='dashed', linewidth=1, color="black")
plt.legend()
plt.show()

この場合の尤度関数(2)の形状はこのようになりました。偏微分するまでもなく、グラフが単峰でできており、$\mu = 0.5$が最尤推定値であることがわかります。

では、今度は対数尤度関数にすると、グラフの形状はどのように変わるのでしょうか。図示して確認してみます。

plt.figure(figsize=(10, 10), dpi=50)
plt.xlim([0, 1.0])
plt.ylim([-8, -6])

mu = np.linspace(0, 1, 1000)

likelihood = (mu ** 5) * ((1 - mu) ** 5)
ln_likelihood = np.log(likelihood)

plt.vlines(0.5, -8, -6 , linestyle='dashed', linewidth=1, color="black")
plt.plot(mu, ln_likelihood, label="log$\mathcal{L}$", color="blue")

plt.legend()
plt.show()

対数尤度関数を図示した結果がこちらです。こちらでも、$\mu = 0.5$が最尤推定地になっていることがわかりますね。

確かに尤度関数と対数尤度関数だと、y方向の値がかなり変わっており、対数尤度関数だと負の値になっていますが、関数の大まかな形は変わらず、最大、最小の場所も変化しないことがわかり、尤度関数の代わりに対数尤度関数を利用して、最大値を得るパラメータを見つけても良いことがわかりました。

【広告】
統計学的にあなたの悩みを解決します。
仕事やプライベートでお悩みの方は、ベテラン占い師 蓮若菜にご相談ください。

機械学習と情報技術