はじめに
『ベイズ推論による機械学習入門』の学習時のノートです。基本的な内容は「数式の行間を読んでみた」とそれを「RとPythonで組んでみた」になります。「数式」と「プログラム」から理解するのが目標です。
この記事は、4.4.4項の内容です。「観測モデルを平均と精度行列が未知の多次元ガウス混合分布(多変量正規混合分布)」、「事前分布をガウス・ウィシャート分布」とする混合モデルを崩壊型ギブスサンプリング(周辺化ギブスサンプリング)により推論します。
省略してある内容等ありますので、本とあわせて読んでください。初学者な自分が理解できるレベルまで落として書き下していますので、分かる人にはかなりくどくなっています。同じような立場の人のお役に立てれば幸いです。
【実装編】
www.anarchive-beta.com
www.anarchive-beta.com
【他の節一覧】
www.anarchive-beta.com
【この節の内容】
4.4.4 崩壊型ギブスサンプリング
崩壊型ギブスサンプリングを用いて、ガウス混合モデルの事後分布$p(\mathbf{X} | \mathbf{S})$から潜在変数$\mathbf{S}$をサンプルするアルゴリズムを導出する。
観測データ$\mathbf{X}$、潜在変数$\mathbf{S}$、観測モデルの平均パラメータ$\boldsymbol{\mu}$、精度行列パラメータ$\boldsymbol{\Lambda}$、混合比率パラメータ$\boldsymbol{\pi}$の同時分布$p(\mathbf{X}, \mathbf{S}, \boldsymbol{\mu},\ \boldsymbol{\Lambda}, \boldsymbol{\pi})$からパラメータを周辺化(積分消去)した
$$
p(\mathbf{X}, \mathbf{S})
= \iiint
p(\mathbf{X}, \mathbf{S}, \boldsymbol{\mu}, \boldsymbol{\Lambda}, \boldsymbol{\pi})
d\boldsymbol{\mu} d\boldsymbol{\Lambda} d\boldsymbol{\pi}
\tag{4.123}
$$
を用いる。
・潜在変数の条件付き分布の導出
$\mathbf{s}_n$をサンプルするための条件付き分布($\mathbf{s}_n$の事後分布)$p(\mathbf{s}_n | \mathbf{X}, \mathbf{S}_{\backslash n})$を求めていく。
観測データ$\mathbf{X}$に加えて、$\mathbf{s}_n$以外の潜在変数$\mathbf{S}_{\backslash n} = \{\mathbf{s}_1, \cdots, \mathbf{s}_{n-1}, \mathbf{s}_{n+1}, \cdots, \mathbf{s}_N\}$が観測された(サンプルされた)としたとき、$\mathbf{s}_n$の条件付き分布は
$$
\begin{align}
p(\mathbf{s}_n | \mathbf{X}, \mathbf{S}_{\backslash n})
&= \frac{
p(\mathbf{x}_n, \mathbf{X}_{\backslash n}, \mathbf{s}_n, \mathbf{S}_{\backslash n})
}{
p(\mathbf{x}_n, \mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n})
}
\\
&\propto
p(\mathbf{x}_n, \mathbf{X}_{\backslash n}, \mathbf{s}_n, \mathbf{S}_{\backslash n})
\\
&= p(\mathbf{x}_n | \mathbf{X}_{\backslash n}, \mathbf{s}_n, \mathbf{S}_{\backslash n})
p(\mathbf{X}_{\backslash n} | \mathbf{s}_n, \mathbf{S}_{\backslash n})
p(\mathbf{s}_n | \mathbf{S}_{\backslash n})
p(\mathbf{S}_{\backslash n})
\\
&\propto
p(\mathbf{x}_n | \mathbf{X}_{\backslash n}, \mathbf{s}_n, \mathbf{S}_{\backslash n})
p(\mathbf{s}_n | \mathbf{S}_{\backslash n})
\tag{4.124}
\end{align}
$$
で求められる。適宜$\mathbf{s}_n$に影響しない項を省いて比例関係に注目する。省略した部分については、最後
に正規化することで対応できる。また、$\mathbf{S} = \{\mathbf{s}_n, \mathbf{S}_{\backslash n}\}$、$\mathbf{X} = \{\mathbf{x}_n, \mathbf{X}_{\backslash n}\}$である。
$\mathbf{s}_n$の条件付き分布の具体的な形状を明らかにしていく。後の項は、「4.3.4：ポアソン混合モデルにおける推論：崩壊型ギブスサンプリング」で求めた
$$
p(\mathbf{s}_n | \mathbf{S}_{\backslash n})
= \prod_{k=1}^K
\eta_{k \backslash n}^{s_{n,k}}
= \mathrm{Cat}(\mathbf{s}_n | \boldsymbol{\eta}_{\backslash n})
\tag{4.74}
$$
パラメータ$\boldsymbol{\eta}_{\backslash n} = (\eta_{1 \backslash n}, \eta_{2 \backslash n}, \cdots, \eta_{K \backslash n})$を持つカテゴリ分布である。ここで
$$
\eta_{\backslash n,k}
= \frac{
\hat{\alpha}_{k \backslash n}
}{
\sum_{k'=1}^K \hat{\alpha}_{\backslash n,k'}
}
\propto
\hat{\alpha}_{k \backslash n}
\tag{4.75}
$$
である。
式(4.124)の前の項は、$\mathbf{X}_{\backslash n},\ \mathbf{S}$が与えられた下での、パラメータ$\boldsymbol{\mu},\ \boldsymbol{\Lambda}$を周辺化した$\mathbf{x}_n$の予測分布
$$
\begin{align}
p(\mathbf{x}_n | \mathbf{X}_{\backslash n}, \mathbf{s}_n, \mathbf{S}_{\backslash n})
&= \iint
p(\mathbf{x}_n, \boldsymbol{\mu}, \boldsymbol{\Lambda} | \mathbf{X}_{\backslash n}, \mathbf{s}_n, \mathbf{S}_{\backslash n})
d\boldsymbol{\mu} d\boldsymbol{\Lambda}
\\
&= \iint
p(\mathbf{x}_n | \mathbf{s}_n, \boldsymbol{\mu}, \boldsymbol{\Lambda})
p(\boldsymbol{\mu}, \boldsymbol{\Lambda} | \mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n})
d\boldsymbol{\mu} d\boldsymbol{\Lambda}
\tag{4.125}
\end{align}
$$
と解釈できる。
さらに、この式の前の項は、ガウス混合モデルの定義(4.4.1項)より
$$
p(\mathbf{x}_n | \mathbf{s}_n, \boldsymbol{\mu}, \boldsymbol{\Lambda})
= \prod_{k=1}^K
\mathcal{N}(\mathbf{x}_n | \boldsymbol{\mu}_k, \boldsymbol{\Lambda}_k^{-1})^{s_{n,k}}
\tag{4.85}
$$
である。
後の項は
$$
\begin{align}
p(\boldsymbol{\mu}, \boldsymbol{\Lambda} | \mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n})
&= \frac{
p(\mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n}, \boldsymbol{\mu}, \boldsymbol{\Lambda})
}{
p(\mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n})
}
\\
&\propto
p(\mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n}, \boldsymbol{\mu}, \boldsymbol{\Lambda})
\\
&= p(\mathbf{X}_{\backslash n} | \mathbf{S}_{\backslash n}, \boldsymbol{\mu}, \boldsymbol{\Lambda})
p(\mathbf{S}_{\backslash n})
p(\boldsymbol{\mu}, \boldsymbol{\Lambda})
\\
&\propto
p(\mathbf{X}_{\backslash n} | \mathbf{S}_{\backslash n}, \boldsymbol{\mu}, \boldsymbol{\Lambda})
p(\boldsymbol{\mu}, \boldsymbol{\Lambda})
\tag{4.126}\\
&= \left\{
\prod_{n'\neq n}
p(\mathbf{x}_{n'} | \mathbf{s}_{n'}, \boldsymbol{\mu}, \boldsymbol{\Lambda})
\right\}
\prod_{k=1}^K
p(\boldsymbol{\mu}_k, \boldsymbol{\Lambda}_k)
\\
&= \left\{
\prod_{n'\neq n} \prod_{k=1}^K
\mathcal{N}(\mathbf{x}_{n'} | \boldsymbol{\mu}, \boldsymbol{\Lambda}^{-1})^{s_{n',k}}
\right\}
\prod_{k=1}^K
\mathrm{NW}(\boldsymbol{\mu}_k, \boldsymbol{\Lambda}_k | \mathbf{m}, \beta, \nu, \mathbf{W})
\end{align}
$$
で求められる。適宜$\boldsymbol{\mu},\ \boldsymbol{\Lambda}$に影響しない項を省いている。
また、左辺の(同時)条件付き分布は
$$
p(\boldsymbol{\mu}, \boldsymbol{\Lambda} | \mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n})
= p(\boldsymbol{\mu} | \boldsymbol{\Lambda}, \mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n})
p(\boldsymbol{\Lambda} | \mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n})
$$
と分解できる。
(対数をとった)式(4.126)は、「4.4.2：ガウス混合モデルにおける推論：ギブスサンプリング」のときの$\boldsymbol{\mu},\ \boldsymbol{\Lambda}$の条件付き分布(4.96)と同じ形である。よって、式(4.126)を$\boldsymbol{\mu}$に関して整理すると式(4.97)と同じ形になるので、同様の手順で
$$
\begin{aligned}
p(\boldsymbol{\mu} | \boldsymbol{\Lambda}, \mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n})
&= \prod_{k=1}^K
p(\boldsymbol{\mu}_k | \boldsymbol{\Lambda}_k, \mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n})
\\
&= \prod_{k=1}^K
\mathcal{N}(\boldsymbol{\mu}_k | \hat{\mathbf{m}}_{k \backslash n}, (\hat{\beta}_{k \backslash n} \boldsymbol{\Lambda}_k)^{-1})
\end{aligned}
$$
$\boldsymbol{\mu}_k$の条件付き分布は、平均$\hat{\mathbf{m}}$、精度$\hat{\beta}_k \boldsymbol{\Lambda}_k$の$D$次元ガウス分布と求められる。ここで
$$
\begin{aligned}
\hat{\beta}_{k \backslash n}
&= \sum_{n'\neq n}
s_{n',k}
+ \beta
\\
\hat{\mathbf{m}}_{k \backslash n}
&= \frac{
\sum_{n'\neq n}
s_{n',k} \mathbf{x}_{n'}
+ \beta \mathbf{m}
}{
\hat{\beta}_{k \backslash n}
}
\end{aligned}
\tag{4.128.a}
$$
である。
また、式(4.126)を$\boldsymbol{\Lambda}$に関して整理すると式(4.101)と同じ形になるので、同様の手順で
$$
\begin{aligned}
p(\boldsymbol{\Lambda} | \mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n})
&= \prod_{k=1}^K
p(\boldsymbol{\Lambda}_k | \mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n})
\\
&= \prod_{k=1}^K
\mathcal{W}(\boldsymbol{\Lambda}_k | \hat{\nu}_{k \backslash n}, \hat{\mathbf{W}}_{k \backslash n})
\end{aligned}
$$
$\boldsymbol{\Lambda}_k$の条件付き分布は、パラメータ$\hat{\mathbf{W}}_k$を持つ自由度$\hat{\nu}_k$のウィシャート分布と求められる。ここで
$$
\begin{aligned}
\hat{\mathbf{W}}_{k \backslash n}^{-1}
&= \sum_{n'\neq n}
s_{n',k} \mathbf{x}_{n'} \mathbf{x}_{n'}^{\top}
+ \beta \mathbf{m} \mathbf{m}^{\top}
- \hat{\beta}_{k \backslash n} \hat{\mathbf{m}}_{k \backslash n} \hat{\mathbf{m}}_{k \backslash n}^{\top}
+ \mathbf{W}^{-1}
\\
\hat{\nu}_{k \backslash n}
&= \sum_{n'\neq n}
s_{n',k}
+ \nu
\end{aligned}
\tag{4.128.b}
$$
である。
それぞれ代入すると
$$
\begin{align}
p(\boldsymbol{\mu}, \boldsymbol{\Lambda} | \mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n})
&= p(\boldsymbol{\mu} | \boldsymbol{\Lambda}, \mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n})
p(\boldsymbol{\Lambda} | \mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n})
\\
&= \prod_{k=1}^K
\mathcal{N}(\boldsymbol{\mu}_k | \hat{\mathbf{m}}_{k \backslash n}, (\hat{\beta}_{k \backslash n} \boldsymbol{\Lambda}_k)^{-1})
\mathcal{W}(\boldsymbol{\Lambda}_k | \hat{\nu}_{k \backslash n}, \hat{\mathbf{W}}_{k \backslash n})
\tag{4.127}
\end{align}
$$
となる。
つまり、式(4.85)と式(4.127)を式(4.25)に代入して、計算を簡単にするため$s_{n,k} = 1$(それ以外の$s_{n,1}, \cdots, s_{n,k-1}, s_{n,k+1}, \cdots, s_{n,K}$は0)の場合を考えると
$$
\begin{align}
p(\mathbf{x}_n | \mathbf{X}_{\backslash n}, s_{n,k} = 1, \mathbf{S}_{\backslash n})
&= \iint
p(\mathbf{x}_n | s_{n,k} = 1, \boldsymbol{\mu}, \boldsymbol{\Lambda})
p(\boldsymbol{\mu}, \boldsymbol{\Lambda} | \mathbf{X}_{\backslash n}, \mathbf{S}_{\backslash n})
d\boldsymbol{\mu} d\boldsymbol{\Lambda}
\tag{4.125'}\\
&= \iint
\mathcal{N}(\mathbf{x}_n | \boldsymbol{\mu}_k, \boldsymbol{\Lambda}_k^{-1})
\mathcal{N}(
\boldsymbol{\mu}_k | \hat{\mathbf{m}}_{k \backslash n}, (\hat{\beta}_{k \backslash n} \boldsymbol{\Lambda}_k)^{-1}
)
\mathcal{W}(
\boldsymbol{\Lambda}_k | \hat{\nu}_{k \backslash n}, \hat{\mathbf{W}}_{k \backslash n}
)
d\boldsymbol{\mu} d\boldsymbol{\Lambda}
\end{align}
$$
となる。
この式は、「3.4.3：多次元ガウス分布の学習と予測：平均・精度が未知の場合」の予測分布(3.134)と同じ形なので、同様の手順で
$$
p(\mathbf{x}_n | \mathbf{X}_{\backslash n}, s_{n,k} = 1, \mathbf{S}_{\backslash n})
= \mathrm{St}(
\mathbf{x}_n | \hat{\boldsymbol{\mu}}^{(s)}_{k \backslash n}, \hat{\boldsymbol{\Lambda}}^{(s)}_{k \backslash n}, \hat{\nu}^{(s)}_{k \backslash n}
)
\tag{4.129}
$$
パラメータ$\hat{\boldsymbol{\mu}}^{(s)}_{k \backslash n},\ \hat{\boldsymbol{\Lambda}}^{(s)}_{k \backslash n},\ \hat{\nu}^{(s)}_{k \backslash n}$を持つ$D$次元スチューデントのt分布と求められる。ここで
$$
\begin{aligned}
\hat{\boldsymbol{\mu}}^{(s)}_{k \backslash n}
&= \hat{\mathbf{m}}_{k \backslash n}
\\
\hat{\boldsymbol{\Lambda}}^{(s)}_{k \backslash n}
&= \frac{
(1 - D + \hat{\nu}_{k \backslash n})
\hat{\beta}_{k \backslash n}
}{
1 + \hat{\beta}_{k \backslash n}
}
\hat{\mathbf{W}}_{k \backslash n}
\\
\hat{\nu}^{(s)}_{k \backslash n}
&= 1 - D + \hat{\nu}_{k \backslash n}
\end{aligned}
$$
である(右肩の$(s)$はスチューデント分布のパラメータを表している)。
よって、$x^0 = 1$であることを利用して、$k = 1, \cdots, K$の項をまとめると
$$
p(\mathbf{x}_n | \mathbf{X}_{\backslash n}, \mathbf{s}_n, \mathbf{S}_{\backslash n})
= \prod_{k=1}^K
\mathrm{St}(
\mathbf{x}_n | \hat{\boldsymbol{\mu}}^{(s)}_{k \backslash n}, \hat{\boldsymbol{\Lambda}}^{(s)}_{k \backslash n}, \hat{\nu}^{(s)}_{k \backslash n}
)^{s_{n,k}}
\tag{4.129'}
$$
となる。
したがって、式(4.129')と式(4.74)を式(4.124)に代入すると、$\mathbf{s}_n$の条件付き分布
$$
\begin{align}
p(\mathbf{s}_n | \mathbf{X}, \mathbf{S}_{\backslash n})
&\propto
p(\mathbf{x}_n | \mathbf{X}_{\backslash n}, \mathbf{s}_n, \mathbf{S}_{\backslash n})
p(\mathbf{s}_n | \mathbf{S}_{\backslash n})
\tag{4.124}\\
&= \left\{
\prod_{k=1}^K
\mathrm{St}(
\mathbf{x}_n | \hat{\boldsymbol{\mu}}^{(s)}_{k \backslash n}, \hat{\boldsymbol{\Lambda}}^{(s)}_{k \backslash n}, \hat{\nu}^{(s)}_{k \backslash n}
)^{s_{n,k}}
\right\}
\mathrm{Cat}(\mathbf{s}_n | \boldsymbol{\eta}_{\backslash n})
\end{align}
$$
が得られる。
実際に$\mathbf{s}_n$をサンプリングするには、スチューデントのt分布(4.129)とカテゴリ分布(4.74)の積が$s_{n,k} = 1$となる確率の比を表すので、$s_{n,1}$から$s_{n,K}$までそれぞれが1となる確率の比を求めて、$K$個の値の和が1となるように正規化したものが$\mathbf{s}_n$のサンプリング確率となる。
別の表現をすると、式(4.66)は$\prod_{k=1}^K \{ \mathrm{St}(\mathbf{x}_n | \hat{\boldsymbol{\mu}}^{(s)}_{k \backslash n}, \hat{\boldsymbol{\Lambda}}^{(s)}_{k \backslash n}, \hat{\nu}^{(s)}_{k \backslash n}) \eta_{n,k} \}^{s_{n,k}}$となるので、波括弧を新たなパラメータとおき正規化することで$\mathbf{s}_n$が従うカテゴリ分布となる。
参考文献
- 須山敦志『ベイズ推論による機械学習入門』(機械学習スタートアップシリーズ)杉山将監修,講談社,2017年.
おわりに
明けましておめでとうございます。哀しいかな、他にすることもなかったので元日からブログを更新することになりました。と言っても、(去年中に終わらせたかった章なので)ほとんど資料としてできていたものを確認して転載しただけですが。まぁ何はともあれ、ブログとしては幸先のいいスタートを切ることになりました。今年もよろしくお願いします。
ところで、最後の最後が本に載ってないんですが？自明なんですか？分かりません。。MLPシリーズの『ノンパラメトリックベイズ』には最後のところも載ってるようですが、パラパラと見ただけでは分かりませんでした(記号が違ったりするからね)。ここのためだけに買うのは躊躇ったので尻切れトンボです、、が一応これで4章の数式行間埋め編は終了です。
全然幸先よくないですが、今年もこんな感じで発狂しながらも楽しんでやっていきます。
【次節の内容】は、ないっ！
5.4節のトピックモデルについてはこちらで書いてます。
www.anarchive-beta.com
