はじめに
『ベイズ推論による機械学習入門』の学習時のノートです。基本的な内容は「数式の行間を読んでみた」とそれを「RとPythonで組んでみた」になります。「数式」と「プログラム」から理解するのが目標です。
この記事は、3.2.1項の内容です。「尤度関数をベルヌーイ分布」、「事前分布をベータ分布」とした場合の「パラメータの事後分布」と「未観測値の予測分布」を導出します。
省略してる内容等ありますので、本とあわせて読んでください。初学者な自分が理解できるレベルまで落として書き下していますので、分かる人にはかなりくどくなっています。同じような立場の人のお役に立てば幸いです。
【実装編】
www.anarchive-beta.com
www.anarchive-beta.com
【他の節の内容】
www.anarchive-beta.com
【この節の内容】
3.2.1 ベルヌーイ分布の学習と予測
ベルヌーイ分布に従うと仮定する$N$個の観測データ$\mathbf{X} = \{x_1, x_2, \cdots, x_N\}$を用いて、パラメータ$\mu$の事後分布と未観測のデータ$x_{*}$の予測分布を求めていく。
各データ$x_n$は、0か1の2値をとる。また$\mu$は、$x_n = 1$となる確率を表すパラメータで、$0 \leq \mu \leq 1$である。
$\mathbf{X}$の各データが独立に発生しているとの仮定の下で、尤度$p(\mathbf{X} | \mu)$は
$$
\begin{aligned}
p(\mathbf{X} | \mu)
&= p(x_1, x_2, \cdots, x_N | \mu)
\\
&= p(x_1 | \mu) p(x_2 | \mu) \cdots p(x_N | \mu)
\\
&= \prod_{n=1}^N p(x_n | \mu)
\end{aligned}
$$
と分解できる。
・事後分布の計算
まずは、パラメータ$\mu$の事後分布$p(\mu | \mathbf{X})$を導出する。
$\mathbf{X}$が与えられた下での$\mu$の事後分布は、ベイズの定理より
$$
\begin{align}
p(\mu | \mathbf{X})
&= \frac{
p(\mathbf{X} | \mu)
p(\mu)
}{
p(\mathbf{X})
}
\\
&= \frac{
\left\{
\prod_{n=1}^N p(x_n | \mu)
\right\}
p(\mu)
}{
p(\mathbf{X})
}
\\
&= \frac{
\left\{
\prod_{n=1}^N \mathrm{Bern}(x_n | \mu)
\right\}
\mathrm{Beta}(\mu | a, b)
}{
p(\mathbf{X})
}
\\
&\propto
\left\{
\prod_{n=1}^N \mathrm{Bern}(x_n | \mu)
\right\}
\mathrm{Beta}(\mu | a, b)
\tag{3.12}
\end{align}
$$
となる。分母の$p(\mathbf{X})$は$\mu$に影響しないため省略して、比例関係のみに注目する。省略した部分については、最後に正規化することで対応できる。
この分布の具体的な形状を明らかにしていく。対数をとって指数部分の計算を分かりやすくして、$\mu$に関して整理すると
$$
\begin{align}
\ln p(\mu | \mathbf{X})
&= \ln \prod_{n=1}^N
\mathrm{Bern}(x_n | \mu)
+ \ln \mathrm{Beta}(\mu | a, b)
- \ln p(\mathbf{X})
\\
&= \sum_{n=1}^N
\ln \mathrm{Bern}(x_n | \mu)
+ \ln \mathrm{Beta}(\mu | a, b)
+ \mathrm{const.}
\\
&= \sum_{n=1}^N \ln \left\{
\mu^{x_n}
(1 - \mu)^{1-x_n}
\right\}
+ \ln \left\{
C_B(a, b)
\mu^{a-1}
(1 - \mu)^{b-1}
\right\}
+ \mathrm{const.}
\\
&= \sum_{n=1}^N \left\{
x_n \ln \mu
+ (1 - x_n) \ln(1 - \mu)
\right\}
+ \ln C_B(a, b)
+ (a - 1) \ln \mu
+ (b - 1) \ln (1 - \mu)
+ \mathrm{const.}
\\
&= \sum_{n=1}^N
x_n \ln \mu
+ \left(N - \sum_{n=1}^N x_n\right) \ln(1 - \mu)
+ (a - 1) \ln \mu
+ (b - 1) \ln (1 - \mu)
+ \mathrm{const.}
\\
&= \left(
\sum_{n=1}^N x_n + a - 1
\right)
\ln \mu
+ \left(
N - \sum_{n=1}^N x_n + b - 1
\right)
\ln (1 - \mu)
+ \mathrm{const.}
\tag{3.13}
\end{align}
$$
【途中式の途中式】
- 式(3.12)の下から2行目の関係を用いている。
- $\mu$に影響しない$- \ln p(\mathbf{X})$を$\mathrm{const.}$とおく。
総乗部分の対数をとると、自然対数の性質$\ln x y = \ln x + \ln y$より
$$
\begin{aligned}
\ln p(\mathbf{X} | \mu)
&= \ln \prod_{n=1}^N p(x_n | \mu)
\\
&= \ln \Bigl\{
p(x_1 | \mu) * p(x_2 | \mu) * \cdots * p(x_N | \mu)
\Bigr\}
\\
&= \ln p(x_1 | \mu) + \ln p(x_2 | \mu) + \cdots + \ln p(x_N | \mu)
\\
&= \sum_{n=1}^N \ln p(x_n | \mu)
\end{aligned}
$$
総和になる。
- 尤度と事前分布に、それぞれ具体的な確率分布の式を代入する。
- それぞれ対数をとる。$\ln x^a = a \ln x$である。
- $\sum_{n=1}^N$に関する括弧を展開し、また$\mu$に影響しない$\ln C_B(a, b)$を$\mathrm{const.}$にまとめる。
- $\mu$と$1 - \mu$の項をそれぞれまとめて式を整理する。
となる。
式(3.13)について
$$
\begin{aligned}
\hat{a}
&= \sum_{n=1}^N x_n + a
\\
\hat{b}
&= N - \sum_{n=1}^N x_n + b
\end{aligned}
\tag{3.15}
$$
とおき
$$
\ln p(\mu | \mathbf{X})
= (\hat{a} - 1)
\ln \mu
+ (\hat{b} - 1)
\ln (1 - \mu)
+ \mathrm{const.}
$$
さらに$\ln$を外して、$\mathrm{const.}$を正規化項に置き換えると(正規化する)と
$$
p(\mu | \mathbf{X})
= C_B(\hat{a}, \hat{b})
\mu^{\hat{a} - 1}
(1 - \mu)^{\hat{b} - 1}
= \mathrm{Beta}(\mu | \hat{a}, \hat{b})
\tag{3.14}
$$
事後分布は式の形状から、パラメータ$\hat{a},\ \hat{b}$を持つベータ分布であることが分かる。
また、式(3.15)が超パラメータ$\hat{a},\ \hat{b}$の計算式(更新式)である。
・予測分布の計算
続いて、ベルヌーイ分布に従う未観測のデータ$x_{*}$に対する予測分布を導出する。
事前分布(観測データによる学習を行っていない分布)$p(\mu)$を用いて、パラメータ$\mu$を周辺化することで予測分布$p(x_{*})$となる。
$$
\begin{align}
p(x_{*})
&= \int
p(x_{*} | \mu)
p(\mu)
d\mu
\\
&= \int
\mathrm{Bern}(x_{*} | \mu)
\mathrm{Beta}(\mu | a, b)
d\mu
\\
&= \int
\mu^{x_{*}}
(1 - \mu)^{1-x_{*}}
C_B(a, b)
\mu^{a-1}
(1 - \mu)^{b-1}
d\mu
\\
&= C_B(a, b)
\int
\mu^{x_{*}+a-1}
(1 - \mu)^{1-x_{*}+b-1}
d\mu
\tag{3.16}
\end{align}
$$
$\mu$と$1 - \mu$の項をまとめて積分部分に注目すると、パラメータ$x_{*} + a,\ 1 - x_{*} + b$を持つ正規化項のないベータ分布の形をしている。
これはベータ分布の正規化項の逆数に変形できる(そもそもこの部分の逆数が正規化項である。詳しくは「1.2.3：ベータ分布【『トピックモデル』の勉強ノート】 - からっぽのしょこ」を参照のこと)。
$$
\int
\mu^{x_{*}+a-1}
(1 - \mu)^{1-x_{*}+b-1}
d\mu
= \frac{1}{C_B(x_{*} + a, 1 - x_{*} + b)}
\tag{3.17}
$$
これを式(3.16)に代入すると、予測分布は
$$
p(x_{*})
= \frac{C_B(a, b)}{C_B(x_{*} + a, x_{*} + b)}
$$
となる。
さらにベータ分布の正規化項(2.42)を用いると
$$
\begin{align}
p(x_{*})
&= C_B(a, b)
\frac{1}{C_B(x_{*} + a, x_{*} + b)}
\\
&= \frac{
\Gamma(a + b)
}{
\Gamma(a)
\Gamma(b)
}
\frac{
\Gamma(x_{*} + a)
\Gamma(1 - x_{*} + b)
}{
\Gamma(x_{*} + a + 1 - x_{*} + b)
}
\\
&= \frac{
\Gamma(a + b)
}{
\Gamma(a)
\Gamma(b)
}
\frac{
\Gamma(x_{*} + a)
\Gamma(1 - x_{*} + b)
}{
\Gamma(a + b + 1)
}
\\
&= \frac{
\Gamma(a + b)
}{
\Gamma(a)
\Gamma(b)
}
\frac{
\Gamma(x_{*} + a)
\Gamma(1 - x_{*} + b)
}{
(a + b) \Gamma(a + b)
}
\\
&= \frac{
\Gamma(x_{*} + a)
\Gamma(1 - x_{*} + b)
}{
\Gamma(a)
\Gamma(b)
(a + b)
}
\tag{3.18}
\end{align}
$$
となる。ガンマ関数の性質$\Gamma(x + 1) = x \Gamma(x)$を用いている。
$x_{*}$は0か1しかとらないため、場合分けしてこの式を整理する。$x_{*} = 1$のとき
$$
\begin{align}
p(x_{*} = 1)
&= \frac{
\Gamma(1 + a)
\Gamma(1 - 1 + b)
}{
\Gamma(a)
\Gamma(b)
(a + b)
}
\tag{3.18'}\\
&= \frac{
a \Gamma(a)
\Gamma(b)
}{
\Gamma(a)
\Gamma(b)
(a + b)
}
\\
&= \frac{a}{a + b}
\tag{3.19}
\end{align}
$$
となる。また$x_{*} = 0$のとき
$$
\begin{align}
p(x_{*} = 0)
&= \frac{
\Gamma(0 + a)
\Gamma(1 - 0 + b)
}{
\Gamma(a)
\Gamma(b)
(a + b)
}
\tag{3.18'}\\
&= \frac{
\Gamma(a)
b \Gamma(b)
}{
\Gamma(a)
\Gamma(b)
(a + b)
}
\\
&= \frac{b}{a + b}
\tag{3.20}
\end{align}
$$
となる。これは
$$
\frac{b}{a + b}
= 1 - \frac{a}{a + b}
$$
と置き換えられる。
従って、$x^0 = 1$であることを利用して式(3.19)と式(3.20)をまとめると、式(3.18)は
$$
\begin{aligned}
p(x_{*})
&= \Bigl( \frac{a}{a + b} \Bigr)^{x_{*}}
\Bigl( \frac{b}{a + b} \Bigr)^{1-x_{*}}
\\
&= \Bigl( \frac{a}{a + b} \Bigr)^{x_{*}}
\Bigl( 1 - \frac{a}{a + b} \Bigr)^{1-x_{*}}
\end{aligned}
$$
と書き換えられる。($x_{*} = 0$のとき前の項が1となり、$x_{*} = 1$のとき後ろの項が1となるので、それぞれ式(3.20)と式(3.19)が成り立つ。)
この式について
$$
\mu_{*}
= \frac{a}{a + b}
$$
とおくと
$$
p(x_{*})
= \mu_{*}^{x_{*}}
(1 - \mu_{*})^{1-x_{*}}
= \mathrm{Bern}(x_{*} | \mu_{*})
\tag{3.21}
$$
予測分布は式の形状から、パラメータ$\mu_{*}$を持つベルヌーイ分布であることが分かる。ちなみに、$\frac{a}{a + b}$はベータ分布$\mathrm{Beta}(\mu | a, b)$の期待値$\mathbb{E}[\mu]$である。
予測分布の計算に事前分布$p(\mu)$を用いることで、観測データによる学習を行っていない予測分布$p(x_{*})$(のパラメータ$\mu_{*}$)を求めた。事後分布$p(\mu | \mathbf{X})$を用いると、同様の手順で観測データ$\mathbf{X}$によって学習した予測分布$p(x_{*} | \mathbf{X})$(のパラメータ)を求められる。
そこで、$\mu_{*}$を構成する事前分布のパラメータ$a,\ b$について、事後分布のパラメータ(3.15)に置き換えたものを$\hat{\mu}_{*}$とおくと
$$
\begin{aligned}
\hat{\mu}_{*}
&= \frac{\hat{a}}{\hat{a} + \hat{b}}
\\
&= \frac{
\sum_{n=1}^N x_n + a
}{
\sum_{n=1}^N x_n + a
+ N - \sum_{n=1}^N x_n + b
}
\\
&= \frac{\sum_{n=1}^N x_n + a}{N + a + b}
\end{aligned}
$$
となり、予測分布
$$
p(x_{*} | \mathbf{X})
= \hat{\mu}_{*}^{x_{*}}
(1 - \hat{\mu}_{*})^{1-x_{*}}
= \mathrm{Bern}(x_{*} | \hat{\mu}_{*})
\tag{3.22}
$$
が得られる。
また上の式が、予測分布のパラメータ$\hat{\mu}_{*}$の計算式(更新式)である。
参考文献
- 須山敦志『ベイズ推論による機械学習入門』(機械学習スタートアップシリーズ)杉山将監修,講談社,2017年.
おわりに
ずっと辞書的に使っていた須山ベイズ本の精読を始めました。
実装編とで口調が変わることに深い意味はありません。正直ですます調の方が書きやすいです。
暫く続きます。宜しくお願いします。
2020年2月19日：モーニング娘。'20森戸知沙希ちゃん二十歳の誕生日！！！
Enjoy!
【次節の内容】
www.anarchive-beta.com
2020/03/04：加筆修正しました。
2020/02/19：加筆修正して、Rで実装編とは記事を分割しました。ブログ投稿から丁度一年後に更新することになるとは。
