はじめに
『ベイズ推論による機械学習入門』の学習時のノートです。「数式の行間を読んでみた」とそれを「RとPythonで組んでみた」によって、「数式」と「プログラム」から理解するのが目標です。
省略してある内容等ありますので、本とあわせて読んでください。
この記事では、「尤度関数を平均が未知の1次元ガウス分布(一変量正規分布)」、「事前分布を1次元ガウス分布」とした場合の「平均パラメータの事後分布」と「未観測値の予測分布」を導出します。
【実装編】
www.anarchive-beta.com
www.anarchive-beta.com
【前の節の内容】
www.anarchive-beta.com
【他の節一覧】
www.anarchive-beta.com
【この節の内容】
3.3.1 平均が未知の場合
尤度関数を1次元ガウス分布(Gaussian Distribution)・一変量正規分布(Normal Distribution)、事前分布を1次元ガウス分布とするモデルに対するベイズ推論を導出する。
ガウス分布については「1次元ガウス分布の定義式 - からっぽのしょこ」を参照のこと。
尤度関数の確認
ガウス分布に従うと仮定する$N$個の観測データ$\mathbf{X} = \{x_1, x_2, \cdots, x_N\}$を用いて、平均パラメータ$\mu$の事後分布と未観測データ$x_{*}$の予測分布を求めていく。
まずは、尤度関数を確認する。
$\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
\mathcal{N}(x_n | \mu, \lambda^{-1})
\right\}
\mathcal{N}(\mu | m, \lambda_{\mu}^{-1})
}{
p(\mathbf{X})
}
\\
&\propto
\left \{\prod_{n=1}^N
\mathcal{N}(x_n | \mu, \lambda^{-1})
\right\}
\mathcal{N}(\mu | m, \lambda_{\mu}^{-1})
\tag{3.49}
\end{align}
$$
となる。分母の$p(\mathbf{X})$は$\mu$に影響しないため省略して、比例関係のみに注目する。省略した部分については、最後に正規化することで対応できる。
この分布の具体的な形状を明らかにしていく。対数をとって指数部分の計算を分かりやすくして、$\mu$に関して整理すると
$$
\begin{align}
\ln p(\mu | \mathbf{X})
&= \sum_{n=1}^N
\ln \mathcal{N}(x_n | \mu, \lambda^{-1})
+ \ln \mathcal{N}(\mu | m, \lambda_{\mu}^{-1})
+ \mathrm{const.}
\\
&= \sum_{n=1}^N
\ln \left\{
\frac{1}{\sqrt{2 \pi \lambda^{-1}}}
\exp\left(
- \frac{(x_n - \mu)^2}{2 \lambda^{-1}}
\right)
\right\}
+ \ln \left\{
\frac{1}{\sqrt{2 \pi \lambda_{\mu}^{-1}}}
\exp\left(
- \frac{(\mu - m)^2}{2 \lambda_{\mu}^{-1}}
\right)
\right\}
+ \mathrm{const.}
\\
&= \sum_{n=1}^N \left\{
- \frac{1}{2}
\ln 2 \pi \lambda^{-1}
- \frac{1}{2}
(x_n - \mu)^2 \lambda
\right\}
- \frac{1}{2}
\ln 2 \pi \lambda_{\mu}^{-1}
- \frac{1}{2}
(\mu - m)^2 \lambda_{\mu}
+ \mathrm{const.}
\\
&= - \frac{1}{2} \left\{
\sum_{n=1}^N
(x_n - \mu)^2 \lambda
+ (\mu - m)^2 \lambda_{\mu}
\right\}
+ \mathrm{const.}
\\
&= - \frac{1}{2} \left\{
\sum_{n=1}^N (
x_n^2 \lambda
- 2 x_n \lambda \mu
+ \lambda \mu^2
)
+ \lambda_{\mu} \mu^2
- 2 m \lambda_{\mu} \mu
+ m^2 \lambda_{\mu}
\right\}
+ \mathrm{const.}
\\
&= - \frac{1}{2} \left(
- 2 \sum_{n=1}^N x_n
\lambda \mu
+ N \lambda \mu^2
+ \lambda_{\mu} \mu^2
- 2 m \lambda_{\mu} \mu
\right)
+ \mathrm{const.}
\\
&= - \frac{1}{2} \left\{
(N \lambda + \lambda_{\mu})
\mu^2
- 2 \left(
\sum_{n=1}^N x_n \lambda
+ m \lambda_{\mu}
\right)
\mu
\right\}
+ \mathrm{const.}
\tag{3.50}
\end{align}
$$
【途中式の途中式】
- 式(3.49)の下から2行目の式を用いている。ただし、$\mu$と関係のない$- \ln p(\mathbf{X})$を$\mathrm{const.}$とおく。
- 尤度と事前分布に、それぞれ具体的な確率分布の式を代入する。
- それぞれ対数をとる。自然対数の性質より、$\ln \sqrt{\frac{1}{x}} = \ln x^{-\frac{1}{2}} = -\frac{1}{2} \ln x$、$\ln \exp(x) = x$である。
- $-\frac{1}{2}$を括り出す。また、$\mu$と関係のない項を$\mathrm{const.}$に含める。
- 2乗に関する括弧を展開する。
- $\sum_{n=1}^N$に関する括弧を展開し、$\mu$と関係のない項を$\mathrm{const.}$に含める。
- $\mu^2$と$\mu$の項をそれぞれまとめて式を整理する。
となる。
これまでと同様に、この式が(対数をとった)1次元ガウス分布の形になっているのを確認したい。そこで
$$
\begin{align}
\hat{\lambda}_{\mu}
&= N \lambda + \lambda_{\mu}
\tag{3.53}\\
\tilde{m}
&= \lambda
\sum_{n=1}^N x_n
+ \lambda_{\mu} m
\end{align}
$$
とおき、平方完成を行う。
$$
\begin{aligned}
\ln p(\mu | \mathbf{X})
&= - \frac{1}{2} \left\{
\hat{\lambda}_{\mu}
\mu^2
- 2 \tilde{m}
\mu
\right\}
+ \mathrm{const.}
\\
&= - \frac{1}{2}
\hat{\lambda}_{\mu} \left\{
\mu^2
- 2
\frac{
\tilde{m}
}{
\hat{\lambda}_{\mu}
}
\mu
\right\}
+ \mathrm{const.}
\\
&= - \frac{1}{2}
\hat{\lambda}_{\mu} \left\{
\mu^2
- 2
\frac{
\tilde{m}
}{
\hat{\lambda}_{\mu}
}
\mu
+ \Bigl(
\frac{\tilde{m}}{\hat{\lambda}_{\mu}}
\Bigr)^2
- \Bigl(
\frac{\tilde{m}}{\hat{\lambda}_{\mu}}
\Bigr)^2
\right\}
+ \mathrm{const.}
\\
&= - \frac{1}{2}
\hat{\lambda}_{\mu} \left\{
\Bigl(
\mu - \frac{\tilde{m}}{\hat{\lambda}_{\mu}}
\Bigr)^2
- \Bigl(
\frac{\tilde{m}}{\hat{\lambda}_{\mu}}
\Bigr)^2
\right\}
+ \mathrm{const.}
\\
&= - \frac{1}{2} \Bigl(
\mu - \frac{\tilde{m}}{\hat{\lambda}_{\mu}}
\Bigr)^2
\hat{\lambda}_{\mu}
+ \frac{\tilde{m}^2}{2 \hat{\lambda}_{\mu}}
+ \mathrm{const.}
\end{aligned}
$$
【途中式の途中式】
- $\hat{\lambda}_{\mu}$を括り出す。
- 括弧内に$(\frac{\tilde{m}}{\hat{\lambda}_{\mu}})^2 - (\frac{\tilde{m}}{\hat{\lambda}_{\mu}})^2 = 0$を加える。
- 2乗の括弧でまとめる。
- 波括弧を展開する。
後の項は$\mu$と無関係なので$\mathrm{const.}$に含めて
$$
\begin{align}
\hat{m}
&= \frac{\tilde{m}}{\hat{\lambda}_{\mu}}
\\
&= \frac{
\lambda
\sum_{n=1}^N x_n
+ \lambda_{\mu} m
}{
\hat{\lambda}_{\mu}
}
\tag{3.54}
\end{align}
$$
とおき
$$
\ln p(\mu | \mathbf{X})
= - \frac{1}{2}
(\mu - \hat{m})^2
\hat{\lambda}_{\mu}
+ \mathrm{const.}
$$
さらに$\ln$を外し、$\mathrm{const.}$を正規化項に置き換える(正規化する)と
$$
p(\mu | \mathbf{X})
= \frac{1}{\sqrt{2 \pi \hat{\lambda}_{\mu}^{-1}}}
\exp \left(
- \frac{1}{2}
(\mu - \hat{m})^2
\hat{\lambda}_{\mu}
\right)
= \mathcal{N}(\mu | \hat{m}, \hat{\lambda}_{\mu}^{-1})
\tag{3.51}
$$
事後分布は式の形状から、平均$\hat{m}$、精度$\hat{\lambda}_{\mu}$の1次元ガウス分布であることが分かる。
また、式(3.53)、(3.54)が超パラメータ$\hat{m}, \hat{\lambda}_{\mu}$の計算式(更新式)である。
予測分布の計算
続いて、1次元ガウス分布に従う未観測のデータ$x_{*}$に対する予測分布を導出する。
事前分布による予測分布
事前分布(観測データによる学習を行っていない分布)$p(\mu)$を用いて、パラメータ$\mu$を周辺化することで予測分布$p(x_{*})$となる。
$$
\begin{align}
p(x_{*})
&= \int
p(x_{*} | \mu)
p(\mu)
d\mu
\\
&= \int
\mathcal{N}(x_{*} | \mu, \lambda^{-1})
\mathcal{N}(\mu | m, \lambda_{\mu}^{-1})
d\mu
\tag{3.55}
\end{align}
$$
前節で行ったようにこの計算を直接するのではなく、より簡単に計算したい。そこで、予測分布$p(x_{*})$と事前分布$p(\mu)$の関係を考える。ベイズの定理より
$$
p(\mu | x_{*})
= \frac{
p(x_{*} | \mu)
p(\mu)
}{
p(x_{*})
}
\tag{3.56}
$$
が成り立つ。この式の両辺の対数をとり
$$
\ln p(\mu | x_{*})
= \ln p(x_{*} | \mu)
+ \ln p(\mu)
- \ln p(x_{*})
$$
予測分布に関して整理すると
$$
\ln p(x_{*})
= \ln p(x_{*} | \mu)
- \ln p(\mu | x_{*})
+ \mathrm{const.}
\tag{3.57}
$$
となる。$x_{*}$に影響しない$\ln p(\mu)$を$\mathrm{const.}$とおいている。
この式から、具体的な予測分布の式を計算する。
$p(\mu | x_{*})$は、1つのデータ$x_{*}$が与えられた下での$\mu$の条件付き分布である。つまり$p(\mu | x_{*})$は、$N$個の観測データが与えられた下での事後分布$p(\mu | \mathbf{X})$と同様の手順で求められる(同様のパラメータになる)。
従って、事後分布のパラメータ(3.53)、(3.54)を用いると、$p(\mu | x_{*})$は$N = 1$より
$$
\ln p(\mu | x_{*})
= - \frac{1}{2} \Bigl\{
(\mu - m_{*})^2
(\lambda + \lambda_{\mu})
+ \ln (\lambda + \lambda_{\mu})^{-1}
+ \ln 2 \pi
\Bigr\}
= \ln \mathcal{N}(
\mu | m_{*}, (\lambda + \lambda_{\mu})^{-1}
)
\tag{3.58}
$$
平均$m_{*}$、精度$\lambda + \lambda_{\mu}$を持つ1次元ガウス分布となる。ただし
$$
m_{*}
= \frac{
\lambda x_{*} + \lambda_{\mu} m
}{
\lambda + \lambda_{\mu}
}
\tag{3.59}
$$
とおく。
尤度(3.47)と式(3.58)を式(3.57)に代入して、$x_{*}$に関して整理すると
$$
\begin{align}
\ln p(x_{*})
&= \ln p(x_{*} | \mu)
- \ln p(\mu | x_{*})
+ \mathrm{const.}
\tag{3.57}\\
&= \ln \mathcal{N}(x_{*} | \mu, \lambda^{-1})
- \ln \mathcal{N}\left(
\mu | m_{*}, (\lambda + \lambda_{\mu})^{-1}
\right)
+ \mathrm{const.}
\\
&= - \frac{1}{2} \Bigl\{
(x_{*} - \mu)^2
\lambda
+ \ln 2 \pi \lambda^{-1}
\Bigr\}
+ \frac{1}{2} \Bigl\{
(\mu - m_{*})^2
(\lambda + \lambda_{\mu})
+ \ln 2 \pi (\lambda + \lambda_{\mu})^{-1}
\Bigr\}
+ \mathrm{const.}
\\
&= - \frac{1}{2} \left\{
(x_{*} - \mu)^2
\lambda
- \Bigl(
\mu
- \frac{
\lambda x_{*} + \lambda_{\mu} m
}{
\lambda + \lambda_{\mu}
}
\Bigr)^2
(\lambda + \lambda_{\mu})
\right\}
+ \mathrm{const.}
\\
&= - \frac{1}{2} \left\{
\lambda x_{*}^2
- 2 \mu \lambda x_{*}
+ \mu^2 \lambda
- \mu^2 (\lambda + \lambda_{\mu})
+ 2 \mu (\lambda x_{*} + \lambda_{\mu} m)
- \frac{
(\lambda x_{*} + \lambda_{\mu} m)^2
}{
\lambda + \lambda_{\mu}
}
\right\}
+ \mathrm{const.}
\\
&= - \frac{1}{2} \left\{
\lambda x_{*}^2
- 2 \mu \lambda x_{*}
+ 2 \mu \lambda x_{*}
- \frac{
\lambda^2 x_{*}^2
+ 2 m \lambda \lambda_{\mu} x_{*}
}{
\lambda + \lambda_{\mu}
}
\right\}
+ \mathrm{const.}
\\
&= - \frac{1}{2} \left\{
\frac{
\lambda (\lambda + \lambda_{\mu})
- \lambda^2
}{
\lambda + \lambda_{\mu}
}
x_{*}^2
- \frac{
2 m \lambda \lambda_{\mu}
}{
\lambda + \lambda_{\mu}
}
x_{*}
\right\}
+ \mathrm{const.}
\\
&= - \frac{1}{2} \left\{
\frac{
\lambda \lambda_{\mu}
}{
\lambda + \lambda_{\mu}
}
x_{*}^2
- \frac{
2 m \lambda \lambda_{\mu}
}{
\lambda + \lambda_{\mu}
}
x_{*}
\right\}
+ \mathrm{const.}
\tag{3.60}
\end{align}
$$
となる。適宜$x_{*}$と関係のない項を$\mathrm{const.}$に含めている。
式(3.60)について
$$
\begin{align}
\lambda_{*}
&= \frac{
\lambda \lambda_{\mu}
}{
\lambda + \lambda_{\mu}
}
\tag{3.62a}\\
\tilde{\mu}_{*}
&= \frac{
2 m \lambda \lambda_{\mu}
}{
\lambda + \lambda_{\mu}
}
\end{align}
$$
とおき、平方完成を行うと
$$
\begin{aligned}
\ln p(x_{*})
&= - \frac{1}{2} \{
\lambda_{*}
x_{*}^2
- \tilde{\mu}_{*}
x_{*}
\}
+ \mathrm{const.}
\\
&= - \frac{1}{2}
\lambda_{*} \left\{
x_{*}^2
- \frac{
\tilde{\mu}_{*}
}{
\lambda_{*}
}
x_{*}
\right\}
+ \mathrm{const.}
\\
&= - \frac{1}{2}
\lambda_{*} \left\{
x_{*}^2
- \frac{
\tilde{\mu}_{*}
}{
\lambda_{*}
}
x_{*}
+ \Bigl(
\frac{
\tilde{\mu}_{*}
}{
2 \lambda_{*}
}
\Bigr)^2
- \Bigl(
\frac{
\tilde{\mu}_{*}
}{
2 \lambda_{*}
}
\Bigr)^2
\right\}
+ \mathrm{const.}
\\
&= - \frac{1}{2}
\lambda_{*} \left\{
\Bigl(
x_{*}^2
- \frac{
\tilde{\mu}_{*}
}{
2 \lambda_{*}
}
\Bigr)^2
- \Bigl(
\frac{
\tilde{\mu}_{*}
}{
2 \lambda_{*}
}
\Bigr)^2
\right\}
+ \mathrm{const.}
\\
&= - \frac{1}{2} \Bigl(
x_{*}^2
- \frac{
\tilde{\mu}_{*}
}{
2 \lambda_{*}
}
\Bigr)^2
\lambda_{*}
+ \frac{\tilde{\mu}_{*}^2}{4 \lambda_{*}}
+ \mathrm{const.}
\end{aligned}
$$
となる。後の項は$x_{*}$と無関係なので$\mathrm{const.}$に含めて
$$
\begin{align}
\mu_{*}
&= \frac{
\tilde{\mu}_{*}
}{
2 \lambda_{*}
}
\\
&= \frac{1}{2}
\frac{
\lambda + \lambda_{\mu}
}{
\lambda \lambda_{\mu}
}
\frac{
2 m \lambda \lambda_{\mu}
}{
\lambda + \lambda_{\mu}
}
\\
&= m
\tag{3.62b}
\end{align}
$$
とおき
$$
\ln p(x_{*})
= - \frac{1}{2}
(x_{*}^2 - \mu_{*})^2
\lambda_{*}
+ \mathrm{const.}
\tag{3.61}
$$
さらに$\ln$を外し、$\mathrm{const.}$を正規化項に置き換える(正規化する)と
$$
p(x_{*})
= \frac{1}{\sqrt{2 \pi \lambda_{*}^{-1}}}
\exp \left(
- \frac{1}{2}
(x_{*}^2 - \mu_{*})^2
\lambda_{*}
\right)
= \mathcal{N}(x_{*} | \mu_{*}, \lambda_{*})
\tag{3.61}
$$
予測分布は式の形状から、平均$\mu_{*}$、精度$\lambda_{*}$の1次元ガウス分布となることが分かる。
ちなみに、精度パラメータ$\lambda_{*}$の逆数である分散$\sigma_{*}^2 = \lambda_{*}^{-1}$に注目すると
$$
\begin{align}
\lambda_{*}^{-1}
&= \frac{
\lambda + \lambda_{\mu}
}{
\lambda \lambda_{\mu}
}
\\
&= \frac{
\lambda_{\mu}
}{
\lambda \lambda_{\mu}
}
+ \frac{
\lambda
}{
\lambda \lambda_{\mu}
}
\\
&= \frac{1}{\lambda}
+ \frac{1}{\lambda_{\mu}}
\\
&= \lambda^{-1}
+ \lambda_{\mu}^{-1}
\tag{3.63}
\end{align}
$$
である。
事後分布による予測分布
予測分布の計算に事前分布$p(\mu)$を用いることで、観測データによる学習を行っていない予測分布$p(x_{*})$(のパラメータ$\mu_{*}, \lambda_{*}$)を求めた。事後分布$p(\mu | \mathbf{X})$を用いると、同様の手順で観測データ$\mathbf{X}$によって学習した予測分布$p(x_{*} | \mathbf{X})$(のパラメータ)を求められる。
そこで、$\mu_{*}, \lambda_{*}$を構成する事前分布のパラメータ$m, \lambda_{\mu}$について、事後分布のパラメータ(3.53)、(3.54)に置き換えたものをそれぞれ$\hat{\mu}_{*}, \hat{\lambda}_{*}$とおくと
$$
\begin{align}
\hat{\lambda}_{*}
&= \frac{
\lambda \hat{\lambda}_{\mu}
}{
\lambda + \hat{\lambda}_{\mu}
}
\tag{3.62a}\\
&= \frac{
\lambda
(N \lambda + \lambda_{\mu})
}{
\lambda
+ N \lambda + \lambda_{\mu}
}
\\
&= \frac{
(N \lambda + \lambda_{\mu})
\lambda
}{
(N + 1) \lambda + \lambda_{\mu}
}
\end{align}
$$
また
$$
\begin{align}
\hat{\mu}_{*}
&= \hat{m}
\tag{3.62b}\\
&= \frac{
\lambda
\sum_{n=1}^N x_n
+ \lambda_{\mu} m
}{
\hat{\lambda}_{\mu}
}
\\
&= \frac{
\lambda
\sum_{n=1}^N x_n
+ \lambda_{\mu} m
}{
N \lambda + \lambda_{\mu}
}
\end{align}
$$
となり、予測分布
$$
p(x_{*} | \mathbf{X})
= \frac{1}{\sqrt{2 \pi \hat{\lambda}_{*}^{-1}}}
\exp\left(
- \frac{1}{2}
(x_{*} - \hat{\mu}_{*})^2
\hat{\lambda}_{*}
\right)
= \mathcal{N}(x_{*} | \hat{\mu}_{*}, \hat{\lambda}_{*})
$$
が得られる。
また、上の式が予測分布のパラメータ$\hat{\mu}_{*}, \hat{\lambda}_{*}$の計算式(更新式)である。
参考文献
- 須山敦志『ベイズ推論による機械学習入門』(機械学習スタートアップシリーズ)杉山将監修,講談社,2017年.
おわりに
ここまでは1日1項ペースで進んでいましたが、次とその次でかなり詰まっております…
2020/03/05：加筆修正しました。
2021/04/04：加筆修正しました。その際にRで実装編と記事を分割しました。これでも多少は読みやすくなったはず。
【次の節の内容】
www.anarchive-beta.com
