変分ベイズ法の解説

# 変分ベイズ法の解説
この記事は[古川研究室 Workout_calendar](https://qiita.com/flab5420/private/aee38a24be700d146817) 6日目の記事です。
本記事は古川研究室の学生が学習の一環として書いたものです。

本記事では、隠れ変数(潜在変数)が出てくるモデルにおけるパラメータの分布推定法としての変分ベイズ法を解説します.
本記事で主に以下の2点について解説していきます.
* 変分ベイズの問題の枠組みについて
* 変分ベイズ法のアルゴリズムの導出



## はじめに
変分ベイズ法とは潜在変数を伴うモデルにおける解法の１つです.通常のベイズ学習では、データ集合$D$が与えられた時のデータ尤度$p(D|\theta)$を最大にするようなパラメータ$\theta$を推定することを考えます.その際に、混合ガウスモデルや潜在ディリクレ配分法などの隠れ変数(潜在変数)を仮定するようなモデルの場合は、積分計算が解析的に計算することができません.そのような場合に、潜在変数やモデルのパラメータをベイズ推定する解法が変分ベイズ法です.




変分ベイズ法の問題の枠組み
---
変分ベイズ法では、推定したいパラメータの事後分布$p(\boldsymbol{\theta}|D)$を直接推定することは諦めて別の確率分布$q(\boldsymbol{\theta})$で近似できると仮定します.出来るだけ近似した確率分布$q(\boldsymbol{\theta})$は出来るだけ$p(\boldsymbol{\theta}|D)$を表現できるような確率分布であることが望まれます.そこで変分ベイズ法では以下の式を最小化することで確率分布$q(\boldsymbol{\theta})$を推定します。

```math
 \hat{q}(\boldsymbol{\theta})=\mathop{\rm arg~min}\limits_{{q}(\boldsymbol{\theta})} D_{KL}(q(\boldsymbol{\theta})||p(\boldsymbol{\theta}|D))
```

また$D_{KL}(q(\boldsymbol{\theta})||p(\boldsymbol{\theta}|D))$はKL情報量という確率分布間の類似度を測る方法の1つです.
$D_{KL}(q(\boldsymbol{\theta})||p(\boldsymbol{\theta}|D))$の式をみてみると

```math
D_{KL}(q(\boldsymbol{\theta})||p(\boldsymbol{\theta}|D))=\int{q(\boldsymbol{\theta})\log\frac{q(\boldsymbol{\theta})}{p(\boldsymbol{\theta}|D)} d\boldsymbol{\theta}}\tag{1}
```

となります.

ここで式(1)をよくみてみると$D_{KL}(q(\boldsymbol{\theta})||p(\boldsymbol{\theta}|D))$式中には$p(\boldsymbol{\theta}|D)$が入っているためこのままでは直接計算できません.そのため変分ベイズ法では以下のように式変形を行い$p(\boldsymbol{\theta}|D)$を消し、$D_{KL}(q(\boldsymbol{\theta})||p(\boldsymbol{\theta}|D))$の下限を計算できるようにしています。
まず、ベイズの定理$p(\theta|D)=\frac{p(D,\theta)}{p(D)}$より$D_{KL}(q(\boldsymbol{\theta})||p(\boldsymbol{\theta}|D))$は以下のように変形できます.

```math
\begin{align}
D_{KL}(q(\theta)||p(\theta|D))&= \int  q(\boldsymbol{\theta})\log\frac{q(\boldsymbol{\theta})p(D)}{p(D,\boldsymbol{\theta})}d\boldsymbol{\theta}\\
&= -\int q(\boldsymbol{\theta})\log \frac{p(D,\boldsymbol{\theta})}{q(\boldsymbol{\theta})} d\boldsymbol{\theta}+\log p(D)
\end{align}

```

$F[q(\boldsymbol{\theta})]=\int{q(\boldsymbol{\theta})\log \frac{p(D,\boldsymbol{\theta})}{q(\boldsymbol{\theta)}} d \boldsymbol{\theta}}$とすると以下のような関係式が出てきます.

```math
\log p(D)=F[q(\theta)]+D_{KL}(q(\theta)||p(\theta|D))\tag{2}

```

この式は、周辺対数尤度$\log p(D)=\log\left(\int{ p(D,\theta)}d \boldsymbol{\theta}\right)$が$F[q(\theta)]$と$D_{KL}(q(\theta)||p(\theta|D))$の足し合わせで表現できることを表しています.ここで重要なのは周辺対数尤度$\log p(D)$は$q(\boldsymbol{\theta})$に依存しないということです.$q(\boldsymbol{\theta})$を変化させると周辺対数尤度$\log p(D)$は変わらず、$F[q(\theta)]$と$D_{KL}(q(\theta)||p(\theta|D))$の比が変化します.
つまり、式(2)から見ると、$D_{KL}(q(\theta)||p(\theta|D))$を最小化するような$q(\boldsymbol{\theta})$を推定するのと$F[q(\boldsymbol{\theta})]$を最大化するような$q(\boldsymbol{\theta})$を推定する問題は等価になります.

よって変分ベイズ法は以下の式(3)のようにも定式化ができます.

```math
 \hat{q}(\boldsymbol{\theta})=\mathop{\rm arg~max}\limits_{{q}(\boldsymbol{\theta})} F[q(\boldsymbol{\theta})]\tag{3}
```

ちなみに$F[q(\boldsymbol{\theta})]$は変分下限といい、周辺対数尤度$\log p(D)$をイェンセンの不等式で評価したものになります.変分ベイズ法の導出は式(3)の変分下限最大化の問題として定式化したものから行なっていきたいと思います.

導出
---
潜在変数を$\boldsymbol{z}$,パラメータを$\boldsymbol{\theta}$, 観測データ集合を$D$として変分下限$F[q(\boldsymbol{z}, \boldsymbol{\theta})]$を最大化するような${q}(\boldsymbol{\theta})$を求めていきたいと思います.最大化する変分下限$F[q(\boldsymbol{z},\boldsymbol{\theta)}]$を具体的に書くと以下の式になります.

```math
 F[q(\boldsymbol{z},\boldsymbol{\theta)}]=\int q(\boldsymbol{z},\boldsymbol{\theta})\log \frac{p(D,\boldsymbol{\theta})}{q(\boldsymbol{z},\boldsymbol{\theta})} d \boldsymbol{\theta}
```

### 平均場近似
変分ベイス法で重要な概念として、平均場近似があります.推定したい$q(\boldsymbol{z},\boldsymbol{\theta)}$は潜在変数$\boldsymbol{z}$とパラメータ$\boldsymbol{\theta}$の同時分布になっているのでパラメータ$\boldsymbol{\theta}$と潜在変数の$\boldsymbol{z}$の全ての組み合わせを考慮しています.変分ベイズ法では各確率変数は独立であると仮定します.各確率変数に独立性の仮定を置いて推定したい$q(\boldsymbol{z},\boldsymbol{\theta)}$を複数の確率分布の積で表現することを平均場近似と言います.どこまで独立であると仮定するかは問題によりけりですが、今回は潜在変数$\boldsymbol{z}$とパラメータ$\boldsymbol{\theta}$が独立であると仮定して解いていきます.つまり、本記事では$q(\boldsymbol{z},\boldsymbol{\theta})$は以下の式のようにと分解できると仮定します.

```math
q(\boldsymbol{z},\boldsymbol{\theta})=q_1(\boldsymbol{z})q_2(\boldsymbol{\theta})
```




### 変分法
平均場近似によって変分下限$F[q(\boldsymbol{z},\boldsymbol{\theta)}]$は以下のようになります.

```math
F[q_1(\boldsymbol{z}),q_2(\boldsymbol{\theta})]=\int{q_1(\boldsymbol{z})q_2(\boldsymbol{\theta})\log \frac{p(D,\theta)}{q_1(\boldsymbol{z})q_2(\boldsymbol{\theta})} d \boldsymbol{\theta}}d \boldsymbol{z}
```

$F[q_1(\boldsymbol{z}),q_2(\boldsymbol{\theta})]$を最大化するような$q_1(\boldsymbol{z})$,$q_2(\boldsymbol{\theta})$を推定するために$q_1(\boldsymbol{z})$と$q_2(\boldsymbol{\theta})$で微分します.


一般に汎関数$F[q,z]=\int f(q,z)d\boldsymbol{z}$の極値を与える$q(z)$は以下のオイラーラグランジュ方程式をとけばいいことが知られています.

```math
 \frac{\partial f(q,z)}{\partial q(\boldsymbol{z})}=0
```

また、$F[q,z]=\int f(q,q^{'},z)d\boldsymbol{z}$などの場合は解き方が複雑になりますが、変分ベイズ法ではこのような複雑なケースは基本的に出てこないので割愛します.

変分下限$F[q_1(\boldsymbol{z}),q_2(\boldsymbol{\theta})]$を最大化するような$q_1(\boldsymbol{z})$,$q_2(\boldsymbol{\theta})$を変分法を使って解いていきます.本記事では$q_1(\boldsymbol{z})$を変分で求めていきます.$q_2(\boldsymbol{\theta})$の場合も同様なので割愛します.$q_1(\boldsymbol{z})$を求める場合は$f(q_1,\boldsymbol{z})$を以下のように定義します.

```math
 f(q_1,z)=-\int q_1(\boldsymbol{z})q_2(\boldsymbol{\theta}) \log \frac{q_1(\boldsymbol{z})q_2(\boldsymbol{\theta})}{p(D,\boldsymbol{z},\theta)} d\boldsymbol{\theta}
```

変分下限$F[q_1(\boldsymbol{z}),q_2(\boldsymbol{\theta})]$は$f(q_1,z)$を用いて以下の式のようになります.

```math
 F[q_1(\boldsymbol{z}),q_2(\boldsymbol{\theta})]=\int f(q_1,z) d\boldsymbol{z}
```


$F[q_1(\boldsymbol{z}),q_2(\boldsymbol{\theta})]$の$q_1(\boldsymbol{z})$の極値はオイラーラグランジュ方程式により

```math
 \frac{\partial f(q_{1},z)}{\partial q_{1}(\boldsymbol{z})}=0
```

を計算すればいいことになります.

### 頑張って計算
頑張って先ほどの式を計算していきます.

```math
\begin{align}
 \frac{\partial f(q_{1},z)}{\partial q_{1}(\boldsymbol{z})}&=-\frac{\partial}{\partial q_{1}(\boldsymbol{z})}\int q_1(\boldsymbol{z})q_2(\boldsymbol{\theta}) \log \frac{q_1(\boldsymbol{z})q_2(\boldsymbol{\theta})}{p(D,\boldsymbol{z},\boldsymbol{\theta})} d\boldsymbol{\theta}\\
&=-\int q_2(\boldsymbol{\theta}) \log \frac{q_1(\boldsymbol{z})q_2(\boldsymbol{\theta})}{p(D,\boldsymbol{z},\boldsymbol{\theta})} d\boldsymbol{\theta}-\int q_{1}(\boldsymbol{z})q_2(\boldsymbol{\theta})\frac{p(D,\boldsymbol{z},\boldsymbol{\theta})}{q_1(\boldsymbol{z})q_2(\boldsymbol{\theta})}\frac{q_{2}(\boldsymbol{\theta})}{p(D,\boldsymbol{z},\boldsymbol{\theta})} d\boldsymbol{\theta}\\
&=-\int q_2(\boldsymbol{\theta}) \log \frac{q_1(\boldsymbol{z})q_2(\boldsymbol{\theta})}{p(D,\boldsymbol{z},\boldsymbol{\theta})} d\boldsymbol{\theta}+\int q_{2}(\boldsymbol{\theta})d\boldsymbol{\theta}\\
&=-\log q_{1}(\boldsymbol{z}) \int q_{2}(\boldsymbol{\theta})d\boldsymbol{\theta}-\int q_{2}(\theta)\log \frac{q_{2}(\theta)}{p(D,\boldsymbol{z},\boldsymbol{\theta})}d\boldsymbol{\theta}-1
\end{align}
```

$\frac{\partial f(q_{1},\boldsymbol{z})}{\partial q_{1}(\boldsymbol{z})}=0$と置いて式を整理していきます.

```math
\begin{align}
\log q_{1}(\boldsymbol{z})&= \int q_{2}(\boldsymbol{\theta})\log \frac{p(D,z,\boldsymbol{\theta})}{q_{2}(\boldsymbol{\theta})}d\boldsymbol{\theta}-1\\
&=\int q_{2}(\boldsymbol{\theta})\log p(D,z,\boldsymbol{\theta})d\boldsymbol{\theta}+\int q_{2}(\boldsymbol{\theta}) \log q_{2}(\boldsymbol{\theta})d\boldsymbol{\theta}-1
\end{align}
```

対数を外すと$q_1(\boldsymbol{z})$は以下の式のようになり,$q_{1}(\boldsymbol{z})$を求めることができました.

```math
\begin{align}
q_{1}(\boldsymbol{z})
&=\exp \left(\int q_{2}(\theta)\log p(D,z,\theta)d\boldsymbol{\theta}\right) \exp\left(\int q_{2}(\boldsymbol{\theta}) \log q_{2}(\boldsymbol{\theta})d\boldsymbol{\theta}\right)\exp(-1)\\
&\propto \exp \left(\int q_{2}(\theta)\log p(D,z,\theta)d\boldsymbol{\theta}\right)
\end{align}
```

$q_{2}(\theta)$の場合も同様に計算できて、以下のような式になります.

```math
q_{2}(\boldsymbol{\theta}) \propto \exp \left(\int q_{1}(\boldsymbol{z})\log p(D,\boldsymbol{z},\theta)d\boldsymbol{z}\right)
```

$q_1(\boldsymbol{z})$と$q_{2}(\theta)$の結果をまとめると以下のような式になります.

```math
\begin{align}
q_{1}(\boldsymbol{z})&\propto \exp \left(\int q_{2}(\theta)\log p(D,\boldsymbol{z},\theta)d\boldsymbol{\theta}\right)\\
q_{2}(\boldsymbol{\theta})
&\propto  \exp \left(\int q_{1}(\boldsymbol{z})\log p(D,\boldsymbol{z},\boldsymbol{\theta})d\boldsymbol{z}\right)
\end{align}
```

## アルゴリズム
最後に変分ベイズ法のアルゴリズムを書いて終わりにします.

1. $q_1(\boldsymbol{z})$と$q_{2}(\boldsymbol{\theta})$を初期化
2. $q_1(\boldsymbol{z})$を使って$q_{2}(\boldsymbol{\theta})$を更新
3. $q_{2}(\boldsymbol{\theta})$を使って$q_1(\boldsymbol{z})$を更新

2,3を収束するまで繰り返す.非常に簡単ですが、変分ベイズ法のアルゴリズムでした.




## 参考資料
トピックモデルによる統計的潜在意味解析読み会 https://www.slideshare.net/tn1031/33-50249999
作って遊ぶ機械学習変分ベイズ編http://machine-learning.hatenablog.com/entry/2016/01/23/130017
変分ベイズの定式化 http://www.nejicolab.com/akip/?tag=機械学習
Variational Bayesian inference https://kaybrodersen.github.io/talks/Brodersen_2013_03_22.pdf
変分ベイズの復習
https://qiita.com/kento1109/items/99a6bcbf18c1119d127c
ガウス過程と機械学習
トピックモデルによる統計的潜在意味解析