はじめに
『ベイズ推論による機械学習入門』の学習時のノートです。基本的な内容は「数式の行間を読んでみた」とそれを「Rで組んでみた」になります。「数式」と「プログラム」から理解するのが目標です。
この記事は4.3.4項の内容になります。ポアソン混合モデルにおける崩壊型ギブスサンプリングを導出し、Rで実装します(未)。
省略してある内容等ありますので、本と併せて読んでください。初学者な自分が理解できるレベルまで落として書き下していますので、分かる人にはかなりくどくなっています。同じような立場の人のお役に立てれば幸いです。
【前節の内容】
【他の節一覧】
【この節の内容】
4.3.4 崩壊型ギブスサンプリング
ポアソン混合モデルに対する崩壊型ギブスサンプリングを導出する。
同時分布からｍパラメータを周辺化する。
$$
p(\boldsymbol{X}, \boldsymbol{S})
= \int \int
p(\boldsymbol{X}, \boldsymbol{S}, \boldsymbol{\lambda}, \boldsymbol{\pi})
d\boldsymbol{\lambda} d\boldsymbol{\pi}
\tag{4.63}
$$
これを用いて、$\boldsymbol{s}_n$の事後分布は
$$
\begin{align*}
p(\boldsymbol{s}_n | \boldsymbol{X}, \boldsymbol{S}_{\backslash n})
&= \frac{
p(x_n, \boldsymbol{X}_{\backslash n}, \boldsymbol{s}_n, \boldsymbol{S}_{\backslash n})
}{
p(x_n, \boldsymbol{X}_{\backslash n}, \boldsymbol{S}_{\backslash n})
}
\\
&\propto
p(x_n, \boldsymbol{X}_{\backslash n}, \boldsymbol{s}_n, \boldsymbol{S}_{\backslash n})
\tag{4.64}\\
&= p(x_n | \boldsymbol{X}_{\backslash n}, \boldsymbol{s}_n, \boldsymbol{S}_{\backslash n})
p(\boldsymbol{X}_{\backslash n} | \boldsymbol{s}_n, \boldsymbol{S}_{\backslash n})
p(\boldsymbol{s}_n | \boldsymbol{S}_{\backslash n})
p(\boldsymbol{S}_{\backslash n})
\tag{4.65}\\
&= p(x_n | \boldsymbol{X}_{\backslash n}, \boldsymbol{s}_n, \boldsymbol{S}_{\backslash n})
p(\boldsymbol{X}_{\backslash n} | \boldsymbol{S}_{\backslash n})
p(\boldsymbol{s}_n | \boldsymbol{S}_{\backslash n})
p(\boldsymbol{S}_{\backslash n})
\\
&\propto
p(x_n | \boldsymbol{X}_{\backslash n}, \boldsymbol{s}_n, \boldsymbol{S}_{\backslash n})
p(\boldsymbol{s}_n | \boldsymbol{S}_{\backslash n})
\tag{4.66}
\end{align*}
$$
となる。
後の項は
$$
p(\boldsymbol{s}_n | \boldsymbol{S}_{\backslash n})
= \int
p(\boldsymbol{s}_n | \boldsymbol{\pi})
p(\boldsymbol{\pi} | \boldsymbol{S}_{\backslash n})
d\boldsymbol{\pi}
\tag{4.70}
$$
のような$\boldsymbol{s}_n$の予測分布であると解釈できる。
更にこの式の後ろの項は
$$
\begin{align*}
p(\boldsymbol{\pi} | \boldsymbol{S}_{\backslash n})
&= \frac{
p(\boldsymbol{\pi}, \boldsymbol{S}_{\backslash n})
}{
p(\boldsymbol{S}_{\backslash n})
}
\\
&\propto
p(\boldsymbol{S}_{\backslash n}, \boldsymbol{\pi})
\\
&= p(\boldsymbol{S}_{\backslash n} | \boldsymbol{\pi})
p(\boldsymbol{\pi})
\tag{4.71}
\end{align*}
$$
となる。
この式は、3.2.2項「カテゴリ分布の学習と予測」の事後分布(3.27)より
$$
p(\boldsymbol{\pi} | \boldsymbol{S}_{\backslash n})
= \mathrm{Dir}(\boldsymbol{\pi} | \hat{\boldsymbol{\alpha}}_{\backslash n})
\tag{4.72}
$$
となる。ただし
$$
\hat{\alpha}_{\backslash n,k}
= \sum_{n'\neq n} s_{n',k} + \alpha_k
\tag{4.73}
$$
である。
よって式(4.70)の計算を進めると、3.2.2項「カテゴリ分布の学習と予測」の予測分布(3.32)より
$$
\begin{align*}
p(\boldsymbol{s}_n | \boldsymbol{S}_{\backslash n})
&= \int
p(\boldsymbol{s}_n | \boldsymbol{\pi})
p(\boldsymbol{\pi} | \boldsymbol{S}_{\backslash n})
d\boldsymbol{\pi}
\tag{4.70}\\
&= \int
\mathrm{Cat}(\boldsymbol{s}_n | \boldsymbol{\pi})
\mathrm{Dir}(\boldsymbol{\pi} | \hat{\boldsymbol{\alpha}}_{\backslash n})
d\boldsymbol{\pi}
\\
&= \mathrm{Cat}(\boldsymbol{s}_n | \boldsymbol{\eta}_{\backslash n})
\tag{4.74}
\end{align*}
$$
となる。ただし
$$
\begin{align*}
\eta_{\backslash n,k}
&= \frac{
\hat{\alpha}_{\backslash n,k}
}{
\sum_{k'=1}^K \hat{\alpha}_{\backslash n,k'}
}
\\
&\propto \hat{\alpha}_{\backslash n,k}
\tag{4.75}
\end{align*}
$$
である。
続いて、式(4.66)の計算に戻る。式(4.66)の前の項は
$$
p(x_n | \boldsymbol{X}_{\backslash n}, \boldsymbol{s}_n, \boldsymbol{S}_{\backslash n})
= \int
p(x_n | \boldsymbol{s}_n, \boldsymbol{\lambda})
p(\boldsymbol{\lambda} | \boldsymbol{X}_{\backslash n}, \boldsymbol{S}_{\backslash n})
d\boldsymbol{\lambda}
\tag{4.76}
$$
のような$x_n$の予測分布であると解釈できる。
更にこの式の後の項は
$$
\begin{align*}
p(\boldsymbol{\lambda} | \boldsymbol{X}_{\backslash n}, \boldsymbol{S}_{\backslash n})
&= \frac{
p(\boldsymbol{X}_{\backslash n}, \boldsymbol{S}_{\backslash n}, \boldsymbol{\lambda})
}{
p(\boldsymbol{X}_{\backslash n}, \boldsymbol{S}_{\backslash n})
}
\\
&\propto
p(\boldsymbol{X}_{\backslash n}, \boldsymbol{S}_{\backslash n}, \boldsymbol{\lambda})
\\
&= p(\boldsymbol{X}_{\backslash n} | \boldsymbol{S}_{\backslash n}, \boldsymbol{\lambda})
p(\boldsymbol{\lambda})
\tag{4.77}
\end{align*}
$$
となる。
この式はポアソン分布とガンマ事前分布を用いたベイズ推論(3.2.3項「ポアソン分布の学習と予測」の事後分布(3.37))になる。
$$
\begin{align*}
\ln p(\boldsymbol{X}_{\backslash n} | \boldsymbol{S}_{\backslash n}, \boldsymbol{\lambda})
p(\boldsymbol{\lambda})
&= \ln \prod_{n'\neq n} \prod_{k=1}^K
p(x_{n'} | \boldsymbol{\lambda})^{s_{n',k}}
+ \ln \prod_{k=1}^K p(\lambda_k)
\\
&= \sum_{k=1}^K \left\{
\sum_{n'\neq n}
s_{n',k}
\ln \mathrm{Poi}(x_{n'} | \lambda_k)
+ \ln \mathrm{Gam}(\lambda_k | a, b)
\right\}
\\
&= \sum_{k=1}^K \left\{
\sum_{n'\neq n}
s_{n',k}
\ln \frac{\lambda_k^{x_{n'}}}{x_{n'}!}
\exp({-\lambda_k})
+ \ln C_G(a, b)
\lambda_k^{a-1}
\exp(- b \lambda_k)
\right\}
\\
&= \sum_{k=1}^K \left\{
\sum_{n'\neq n}
s_{n',k} (
x_{n'} \ln \lambda_k
- \lambda_k
)
+ (a - 1) \ln \lambda_k
- b \lambda_k
\right\}
+ \mathrm{const.}
\\
&= \sum_{k=1}^K \left\{
\Bigl(
\sum_{n'\neq n} s_{n',k} x_{n'} + a - 1
\Bigr)
\ln \lambda_k
- \Bigl(
\sum_{n'\neq n} s_{n',k} + b
\Bigr)
\lambda_k
\right\}
+ \mathrm{const.}
\tag{4.78}
\end{align*}
$$
この式について
$$
\begin{align*}
\hat{a}_{\backslash n,k}
&= \sum_{n'\neq n} s_{n',k} x_{n'} + a
\\
\hat{b}_{\backslash n,k}
&= \sum_{n'\neq n} s_{n',k} + b
\tag{4.80}
\end{align*}
$$
とおくと、$\lambda_k$の近似事後分布は
$$
\ln p(\boldsymbol{\lambda} | \boldsymbol{X}_{\backslash n}, \boldsymbol{S}_{\backslash n})
= \sum_{k=1}^K \left\{
\hat{a}_{\backslash n,k}
\ln \lambda_k
- \hat{b}_{\backslash n,k}
\lambda_k
\right\}
+ \mathrm{const.}
= \ln \prod_{k=1}^K \mathrm{Gam}(\lambda_k | \hat{a}_{\backslash n,k}, \hat{b}_{\backslash n,k})
\tag{4.79}
$$
パラメータ$\hat{a}_{\backslash n,k},\ \hat{b}_{\backslash n,k}$を持つガンマ分布になることが分かる。
よって式(4.76)の計算を進める。ただし計算を簡単にするため、$s_{n,k} = 1$(それ以外の1から$K$は0)の場合を考える。3.2.3項「ポアソン分布の学習と予測」の予測分布より
$$
\begin{align*}
p(x_n | \boldsymbol{X}_{\backslash n}, s_{n,k} = 1, \boldsymbol{S}_{\backslash n})
&= \int
p(x_n | \lambda_k)
p(\lambda_k | \boldsymbol{X}_{\backslash n}, \boldsymbol{S}_{\backslash n})
d\lambda_k
\tag{4.76'}\\
&= \mathrm{NB} \left(
x_n \middle| \hat{a}_{\backslash n,k}, \frac{1}{\hat{b}_{\backslash n,k} + 1}
\right)
\tag{4.81}
\end{align*}
$$
パラメータ$\hat{a}_{\backslash n,k},\ \frac{1}{\hat{b}_{\backslash n,k} + 1}$を持つ負の二項分布になる。
参考文献
- 須山敦志『ベイズ推論による機械学習入門』(機械学習スタートアップシリーズ)杉山将監修,講談社,2017年.
おわりに
…ん？式(4.66)で$\boldsymbol{s}_n$をサンプリングするんだよな？？式(4.75)と式(4.81)を組み合わせるとどうなるんだ…？という訳で実装編はまたいずれ。
【次節の内容】続く
2020.04.26：元こぶしファクトリーはまちゃん二十歳の誕生日！