はじめに

　『ベイズ推論による機械学習入門』(MLSシリーズ)の独学時のノートです。各種のモデルやアルゴリズムについて「数式・プログラム・図」を用いて解説します。
　本の補助として読んでください。

　この記事では、ポアソン分布に対するベイズ推論の数式の行間を埋めます。

【前節の内容】

【他の節の内容】

【この節の内容】

3.2.3 ポアソン分布のベイズ推論の導出

　ポアソンモデル(Poisson model)に対するベイズ推論(Bayesian inference)を導出する。ポアソンモデルでは、尤度関数をポアソン分布(Poisson distribution)、事前分布をガンマ分布(Gamma distribution)とする。
　ポアソンモデルについては「3.2.3：ポアソンモデルの生成モデルの導出【緑ベイズ入門のノート】 - からっぽのしょこ」、ポアソン分布については「ポアソン分布の定義式 - からっぽのしょこ」を参照のこと。

事後分布の導出

　まずは、ポアソン分布のパラメータ $\lambda$ の事後分布(posterior distribution)を導出する。
　ガンマ分布については「準備中」を参照のこと。

　観測データ $\mathbf{X}$ が与えられた(条件とする)下でのパラメータ $\lambda$ の条件付き分布(事後分布)を求める。

$\displaystyle \begin{align} p(\lambda \mid \mathbf{X}, a, b) &= \frac{ p(\mathbf{X} \mid \lambda) p(\lambda \mid a, b) }{ p(\mathbf{X} \mid a, b) } \\ &\propto p(\mathbf{X} \mid \lambda) p(\lambda \mid a, b) \\ &= \left\{ \prod_{n=1}^N p(x_n \mid \lambda) \right\} p(\lambda \mid a, b) \\ &= \left\{ \prod_{n=1}^N \mathrm{Poi}(x_n \mid \lambda) \right\} \mathrm{Gam}(\lambda \mid a, b) \tag{3.35} \end{align}$

途中式の途中式(クリックで展開)

1: ベイズの定理 $p(y \mid x) = \frac{p(x \mid y) p(y)}{p(x)}$ より、観測変数 $\mathbf{X}$ を条件に移した式を立てる。

　事後分布は、 $\mathbf{X}, \lambda$ の結合分布と $\mathbf{X}$ の周辺分布を用いて、次のようにも求められる。

$\displaystyle \begin{aligned} p(\lambda \mid \mathbf{X}, a, b) &= \frac{ p(\mathbf{X}, \lambda \mid a, b) }{ p(\mathbf{X} \mid a, b) } \\ &= \frac{ p(\mathbf{X}, \lambda \mid a, b) }{ \int p(\mathbf{X}, \lambda \mid a, b) \mathrm{d} \lambda } \\ &= \frac{ p(\mathbf{X} \mid \lambda) p(\lambda \mid a, b) }{ \int p(\mathbf{X} \mid \lambda) p(\lambda \mid a, b) \mathrm{d} \lambda } \end{aligned}$

　1行目では、条件付き分布 $p(y \mid x) = \frac{p(x, y)}{p(x)}$ より、 $\mathbf{X}$ を条件に移している。
　2行目では、周辺化 $p(y) = \int p(x, y) \mathrm{d} x$ した $\lambda$ を明示している。
　3行目では、 $\mathbf{X}, \lambda$ の依存関係に従い項を分割している。
　生成モデル(結合分布)については「生成モデルの導出」を参照のこと。

2: $\lambda$ と無関係な項を省く。
3: 観測データ集合 $\mathbf{X}$ の生成確率を、各データ $x_n$ の生成確率の積に分解する。
4: ポアソンモデルの定義より、尤度関数をポアソン分布、事前分布をガンマ分布に置き換える。

　周辺分布(分母)は $\lambda$ に影響しないため省いて、比例関係のみに注目する。省略した項については、最後に正規化することで対応できる。

　両辺の対数をとり、指数部分の計算を分かりやすくして、 $\lambda$ に関して式を整理する。

$\displaystyle \begin{aligned} \ln p(\lambda \mid \mathbf{X}, a, b) &= \ln \Bigl( \frac{ \left\{ \prod_{n=1}^N p(x_n \mid \lambda) \right\} p(\lambda \mid a, b) }{ p(\mathbf{X} \mid a, b) } \Bigr) \\ &= \ln \Bigl( \prod_{n=1}^N p(x_n \mid \lambda) \Bigr) + \ln p(\lambda \mid a, b) - \ln p(\mathbf{X} \mid a, b) \\ &= \sum_{n=1}^N \ln p(x_n \mid \lambda) + \ln p(\lambda \mid a, b) + \mathrm{const.} \\ &= \sum_{n=1}^N \ln \mathrm{Poi}(x_n \mid \lambda) + \ln \mathrm{Gam}(\lambda \mid a, b) + \mathrm{const.} \end{aligned}$

途中式の途中式(クリックで展開)

1: 式(3.35)に関して、対数をとった式を立てる。
2-3: 自然対数の性質 $\ln(x y) = \ln x + \ln y$ 、 $\ln \frac{x}{y} = \ln x - \ln y$ より、分数の項を展開する。

　対数の性質より、総乗 $\prod_n$ の対数をとると、対数の総和 $\sum_n$ になる。

$\displaystyle \begin{aligned} \ln p(\mathbf{X} \mid \lambda) &= \ln \prod_{n=1}^N p(x_n \mid \lambda) \\ &= \ln \Bigl( p(x_1 \mid \lambda) * p(x_2 \mid \lambda) * \cdots * p(x_N \mid \lambda) \Bigr) \\ &= \ln p(x_1 \mid \lambda) + \ln p(x_2 \mid \lambda) + \cdots + \ln p(x_N \mid \lambda) \\ &= \sum_{n=1}^N \ln p(x_n \mid \lambda) \end{aligned}$

3: $\lambda$ と無関係な項を $\mathrm{const.}$ とおく。
4: ポアソンモデルの定義より、尤度関数をポアソン分布、事前分布をガンマ分布に置き換える。

　 $\lambda$ に影響しない項を $\mathrm{const.}$ とおく。省略した項については、最後に正規化することで対応できる。

　右辺の各分布に具体的な式を代入して、式の形状を明らかにしていく。

$\displaystyle \begin{align} \ln p(\lambda \mid \mathbf{X}, a, b) &= \sum_{n=1}^N \ln \Bigl( \frac{\lambda^{x_n}}{x_n!} \exp(-\lambda) \Bigr) \\ &\quad + \ln \Bigl( \mathrm{C}_\mathrm{Gam}(a, b) \lambda^{a-1} \exp(- b \lambda) \Bigr) + \mathrm{const.} \\ &= \sum_{n=1}^N \Bigl\{ x_n \ln \lambda - \ln x_n! - \lambda \Bigr\} \\ &\quad + \ln \mathrm{C}_\mathrm{Gam}(a, b) + (a - 1) \ln \lambda - b \lambda + \mathrm{const.} \\ &= \sum_{n=1}^N x_n \ln \lambda + (a - 1) \ln \lambda - N \lambda - b \lambda + \mathrm{const.} \\ &= \left( \sum_{n=1}^N x_n + a - 1 \right) \ln \lambda - (N + b) \lambda + \mathrm{const.} \tag{3.36} \end{align}$

途中式の途中式(クリックで展開)

1: 尤度関数はポアソン分布、事前分布はガンマ分布を仮定しているので、それぞれ定義式に置き換える。

$\displaystyle \begin{align} p(x_n \mid \lambda) &= \mathrm{Poi}(x_n \mid \lambda) \tag{3.33}\\ &= \exp(-\lambda) \frac{\lambda^{x_n}}{x_n!} \\ p(\lambda \mid a, b) &= \mathrm{Gam}(\lambda \mid a, b) \tag{3.34}\\ &= \mathrm{C}_\mathrm{Gam}(a, b) \lambda^{a-1} \exp(- b \lambda) \end{align}$

　ここで、 $\mathrm{C}_\mathrm{Gam}(a, b) = \frac{b^a}{\Gamma(a)}$ は、ガンマ分布(事前分布)の正規化項である。(式変形に影響しないので簡易的に表記している。)

2: 自然対数の性質 $\ln(x y) = \ln x + \ln y$ 、 $\ln \frac{x}{y} = \ln x - \ln y$ 、 $\ln x^y = y \ln x$ 、対数と指数の関係 $\ln (\exp (x)) = x$ より、定義式の項を展開する。
3: $n$ に関する総和 $\sum_n$ の波括弧を展開する。 $n$ と無関係な項は $N$ 倍 $\sum_{n=1}^N \lambda = N \lambda$ となる。
4: $\ln \lambda, \lambda$ の項をそれぞれまとめる。

　適宜、 $\lambda$ に影響しない項を $\mathrm{const.}$ にまとめている。

　事後分布の式(3.36)について、次のようにおく。

$\displaystyle \begin{aligned} \hat{a} &= \sum_{n=1}^N x_n + a \\ \hat{b} &= N + b \end{aligned} \tag{3.38}$

　式(3.36)について、 $\hat{a}, \hat{b}$ で置き換える。

$\displaystyle \ln p(\lambda \mid \mathbf{X}, a, b) = (\hat{a} - 1) \ln \lambda - \hat{b} \lambda + \mathrm{const.}$

　さらに、 $\ln$ を外して $\mathrm{const.}$ を正規化項に置き換える(正規化する)と、事後分布は式の形状から、パラメータ $\hat{a}, \hat{b}$ のガンマ分布であることが分かる。

$\displaystyle \begin{align} p(\lambda \mid \mathbf{X}, a, b) &= \mathrm{C}_\mathrm{Gam}(\hat{a}, \hat{b}) \lambda^{\hat{a}-1} \exp (- \hat{b} \lambda) \\ &= \frac{\hat{b}^{\hat{a}}}{\Gamma(\hat{a})} \lambda^{\hat{a}-1} \exp (- \hat{b} \lambda) \\ &= \mathrm{Gam}(\lambda \mid \hat{a}, \hat{b}) \tag{3.37} \end{align}$

　 $\lambda$ の事後分布の式が得られた。
　ここで、 $\mathrm{C}_\mathrm{Gam}(\hat{a}, \hat{b})$ は、ガンマ分布(事後分布)の正規化項である。
　また、式(3.38)が、事後分布のパラメータ(超パラメータ) $\hat{a}, \hat{b}$ の計算式(更新式)である。

　以上で、ポアソンモデルにおける事後分布を導出した。

スポンサードリンク

予測分布の導出

　次は、ポアソン分布に従う未観測データ $x_{*}$ の予測分布(predict distribution)を導出する。
　負の二項分布については「負の二項分布の定義式 - からっぽのしょこ」を参照のこと。

事前分布による予測分布

　事前分布(観測データによる学習を行っていない $\lambda$ の分布)を用いた予測分布(事前予測分布)を求める。

$\displaystyle \begin{align} p(x_{*} \mid a, b) &= \int p(x_{*}, \lambda \mid a, b) \mathrm{d} \lambda \\ &= \int p(x_{*} \mid \lambda) p(\lambda \mid a, b) \mathrm{d} \lambda \\ &= \int \mathrm{Poi}(x_{*} \mid \lambda) \mathrm{Gam}(\lambda \mid a, b) \mathrm{d} \lambda \tag{3.39} \end{align}$

途中式の途中式(クリックで展開)

1: 未知変数 $x_{*}$ とパラメータ $\lambda$ の結合分布に対して、 $\lambda$ を周辺化した式を立てる。
2: 依存関係のある $x_{*}, \lambda$ の項を分割する。
3: ポアソンモデルの定義より、尤度関数をポアソン分布、事前分布をガンマ分布に置き換える。

　事前予測分布は、未知のデータ $x_{*}$ の生成分布(3.33)と、パラメータ $\lambda$ の事前分布(3.34)を用いた、 $x_{*}$ の周辺分布である。

　右辺の各分布について具体的な式を代入して、式の形状を明らかにしていく。

$\displaystyle \begin{align} p(x_{*} \mid a, b) &= \int \frac{\lambda^{x_{*}}}{x_{*}!} \exp(- \lambda) \mathrm{C}_{\mathrm{Gum}}(a, b) \lambda^{a-1} \exp(- b \lambda) \mathrm{d} \lambda \\ &= \frac{\mathrm{C}_{\mathrm{Gum}}(a, b)}{x_{*}!} \int \lambda^{x_{*}} \lambda^{a-1} \exp(- \lambda) \exp(- b \lambda) \mathrm{d} \lambda \\ &= \frac{\mathrm{C}_{\mathrm{Gum}}(a, b)}{x_{*}!} \int \lambda^{x_{*} + a-1} \exp \Bigl( - (1 + b) \lambda \Bigr) \mathrm{d} \lambda \tag{3.40} \end{align}$

途中式の途中式(クリックで展開)

1: 尤度関数はポアソン分布、事前分布はガンマ分布を仮定しているので、それぞれ定義式に置き換える。

$\displaystyle \begin{align} p(x_{*} \mid \lambda) &= \mathrm{Poi}(x_{*} \mid \lambda) \tag{3.33}\\ &= \exp(-\lambda) \frac{\lambda^{x_{*}}}{x_{*}!} \\ p(\lambda \mid a, b) &= \mathrm{Gam}(\lambda \mid a, b) \tag{3.34}\\ &= \mathrm{C}_\mathrm{Gam}(a, b) \lambda^{a-1} \exp(- b \lambda) \end{align}$

2: $\lambda$ と無関係な項を $\int \mathrm{d} \lambda$ の外に出す。
3: 指数の性質 $e^{x+y} = e^x e^y$ より、 $\lambda, \exp(\lambda)$ の項をそれぞれまとめる。

　 $\lambda$ に関する積分の項に注目すると、パラメータ $x_{*} + a, 1 + b$ の正規化項のないガンマ分布の形をしている。
　この積分の項は、ガンマ分布の正規化項の逆数

$\displaystyle \int \lambda^{x_{*} + a-1} \exp \Bigl( - (1 + b) \lambda \Bigr) \mathrm{d} \lambda = \frac{1}{\mathrm{C}_{\mathrm{Gum}}(x_{*} + a, 1 + b)} \tag{3.41}$

に変形できる(そもそも確率分布の関数(正規化項以外)の項を積分した逆数が正規化項である)。
　予測分布の式(3.40)について、正規化項の逆数の式(3.41)で置き換える。

$\displaystyle p(x_{*} \mid a, b) = \frac{ \mathrm{C}_{\mathrm{Gum}}(a, b) }{ x_{*}! \mathrm{C}_{\mathrm{Gum}}(x_{*} + a, 1 + b) } \tag{3.42}$

　右辺の各正規化項について具体的な式を代入して、式の形状を明らかにしていく。

$\displaystyle \begin{align} p(x_{*} \mid a, b) &= \frac{1}{x_{*}!} \frac{b^a}{\Gamma(a)} \frac{ \Gamma(x_{*} + a) }{ (1 + b)^{x_{*} + a} } \\ &= \frac{ \Gamma(x_{*} + a) }{ x_{*}! \Gamma(a) } \frac{b^a}{(1 + b)^a} \frac{1}{(1 + b)^{x_{*}}} \\ &= \frac{ \Gamma(x_{*} + a) }{ x_{*}! \Gamma(a) } \Bigl( \frac{b}{1 + b} \Bigr)^{a} \Bigl( \frac{1}{1 + b} \Bigr)^{x_{*}} \\ &= \frac{ \Gamma(x_{*} + a) }{ x_{*}! \Gamma(a) } \Bigl( 1 - \frac{1}{1 + b} \Bigr)^{a} \Bigl( \frac{1}{1 + b} \Bigr)^{x_{*}} \tag{1} \end{align}$

途中式の途中式(クリックで展開)

1: ガンマ分布の正規化項(2.57)より、それぞれ置き換える。
2: 指数の性質 $x^{a + b} = x^a x^b$ より、 $1 + b$ の項を分割する。
3: 指数の性質 $(\frac{x}{y})^a = \frac{x^a}{y^a}$ より、 $a, x_{*}$ 指数の項をそれぞれまとめる。
4: $\frac{b}{1 + b} = \frac{1 + b}{1 + b} - \frac{1}{1 + b}$ より、分数の項の形を揃える。

　(次の式の) $p$ を成功確率とすると、 $q = 1-p$ は失敗確率と言える。

　予測分布の式(1)について、次のようにおく。

$\displaystyle \begin{aligned} r &= a \\ q &=\frac{1}{1 + b} \\ p &= 1 - q \\ &= \frac{b}{1 + b} \end{aligned} \tag{3.44}$

　式(1)について、 $r, p$ で置き換えると、予測分布は式の形状から、パラメータ $r, p$ の負の二項分布であることが分かる。

$\displaystyle \begin{align} p(x_{*} \mid a, b) &= \frac{ \Gamma(x_{*} + r) }{ x_{*}! \Gamma(r) } (1 - q)^r q^{x_{*}} \\ &= \frac{ \Gamma(x_{*} + r) }{ x_{*}! \Gamma(r) } p^r (1 - p)^{x_{*}} \\ &= \mathrm{NB}(x_{*} \mid r, p) \tag{3.43} \end{align}$

　 $x_{*}$ の事前予測分布の式が得られた。(成功確率 $p$ と失敗確率 $q$ の表記が本と異なる点に注意。)
　また、式(3.44)が、予測分布のパラメータ $r, p$ の計算式である。

事後分布による予測分布

　予測分布の計算に事前分布 $p(\lambda \mid a, b)$ を用いて、観測データ $\mathbf{X}$ による学習を行っていない予測分布(事前予測分布) $p(x_{*} \mid a, b)$ (のパラメータ $r, p$ )を求めた。事後分布 $p(\lambda \mid \mathbf{X}, a, b)$ を用いると、観測データ $\mathbf{X}$ によって学習した予測分布(事後予測分布) $p(x_{*} \mid \mathbf{X}, a, b)$ (のパラメータ $\hat{r}, \hat{p}$ )を求められる。

$\displaystyle \begin{align} p(x_{*} \mid \mathbf{X}, a, b) &= \int p(x_{*}, \lambda \mid \mathbf{X}, a, b) \mathrm{d} \lambda \\ &= \int p(x_{*} \mid \lambda) p(\lambda \mid \mathbf{X}, a, b) \mathrm{d} \lambda \\ &= \int \mathrm{Poi}(x_{*} \mid \lambda) \mathrm{Gam}(\lambda \mid \hat{a}, \hat{b}) \mathrm{d} \lambda \tag{3.39'} \end{align}$

途中式の途中式(クリックで展開)

1: 観測変数 $\mathbf{X}$ を条件として、未知変数 $x_{*}$ とパラメータ $\lambda$ の結合分布に対して、 $\lambda$ を周辺化した式を立てる。
2: 依存関係のある $x_{*}, \lambda$ の項を分割する。
3: ポアソンモデルの定義より、尤度関数をポアソン分布、事後分布をガンマ分布に置き換える。

　事後予測分布は、未知のデータ $x_{*}$ の生成分布(3.33)と、パラメータ $\lambda$ の事後分布(3.37)を用いた、 $x_{*}$ の周辺分布である。

　事後分布は事前分布と同じくガンマ分布なので、事前予測分布の式(3.43)と、同様の手順で事後予測分布の式も求められる。
　そこで、事前予測分布のパラメータ $r, p$ の式(3.44)を構成する事前分布のパラメータ $a, b$ について、事後分布のパラメータ $\hat{a}, \hat{b}$ の式(3.38)に置き換えたものを事後予測分布のパラメータ $\hat{r}, \hat{p}$ とおく。

$\displaystyle \begin{aligned} \hat{r} &= \hat{a} \\ &= \sum_{n=1}^N x_n + a \\ \hat{q} &= \frac{1}{1 + \hat{b}} \\ &= \frac{1}{N + 1 + b} \\ \hat{p} &= 1 - \hat{q} \\ &= \frac{\hat{b}}{1 + \hat{b}} \\ &= \frac{N + b}{N + 1 + b} \end{aligned} \tag{3.44'}$

　予測分布の式(3.43)についても置き換える(同様の手順で導出する)と、パラメータ $\hat{r}, \hat{p}$ の負の二項分布となる。

$\displaystyle \begin{align} p(x_{*} \mid \mathbf{X}, a, b) &= \frac{ \Gamma(x_{*} + \hat{r}) }{ x_{*}! \Gamma(\hat{r}) } \hat{p}^{\hat{r}} (1 - \hat{p})^{x_{*}} \\ &= \mathrm{NB}(x_{*} \mid \hat{r}, \hat{p}) \tag{3.43'} \end{align}$

　 $x_{*}$ の事後予測分布の式が得られた。(成功確率 $\hat{p}$ と失敗確率 $\hat{q}$ の表記が本と異なる点に注意。)
　また、式(3.44')が、予測分布のパラメータ $\hat{r}, \hat{p}$ の計算式(更新式)である。

　ちなみに、予測分布の期待値は

$\displaystyle \begin{align} \mathbb{E}_{\mathrm{NB}(x_{*} \mid \hat{r}, \hat{p})}[x] &= \frac{\hat{r} (1 - \hat{p})}{\hat{p}} \\ &= \hat{a} \frac{1}{1 + \hat{b}} \frac{1 + \hat{b}}{\hat{b}} \\ &= \frac{\hat{a}}{\hat{b}} \end{align}$