変分ベイズ法

量子力学における近似法については、「変分法 (量子力学)」を参照してください。

変分ベイズ法は、ベイズ推論や機械学習において生じる扱いにくい積分を近似するための手法群です。これらの手法は、観測変数（通常は「データ」と呼ばれる）と未知のパラメータ、そして潜在変数から構成される複雑な統計モデルにおいて用いられます。これらのモデルでは、3種類の確率変数間の様々な関係がグラフィカルモデルによって記述される可能性があります。ベイズ推論では典型的に、パラメータと潜在変数は「観測されない変数」としてグループ化されます。変分ベイズ法は主に以下の2つの目的で用いられます。

1. 観測されない変数の事後確率に対する解析的近似値を提供し、これらの変数に対して統計的推論を行うため。
2. 観測データの周辺尤度（証拠と呼ばれることもある）の下限値を導出する（つまり、観測されていない変数に対して周辺化を行ったモデルを与えられたデータの周辺確率）。これは通常、モデル選択を行う際に用いられ、一般的な考え方としては、あるモデルの周辺尤度が高いほど、そのモデルがデータに適合していることを示し、したがって、そのモデルがデータを生成した確率が高いことを示している。（ベイズ因子の記事も参照のこと。）

前者の目的（事後確率の近似）において、変分ベイズはモンテカルロサンプリング法、特にギブスサンプリングなどのマルコフ連鎖モンテカルロ法の代替手段であり、直接評価したりサンプリングしたりすることが難しい複雑な分布に対して、完全なベイズ的アプローチで統計的推論を行います。特に、モンテカルロ法はサンプルセットを用いて正確な事後確率の数値近似値を提供するのに対し、変分ベイズは事後確率の近似値に対して局所最適かつ正確な解析解を提供します。

変分ベイズは、期待値最大化（EM）アルゴリズムの拡張版と見ることができます。EMアルゴリズムは、各パラメータの最も可能性の高い単一の値の最大尤度（ML）または最大事後確率（MAP）推定から、パラメータと潜在変数の事後分布全体（の近似値）を計算する完全ベイズ推定へと進化しました。変分ベイズはEMと同様に、最適なパラメータ値の集合を求め、解析的に解くことのできない、相互に依存した一連の方程式に基づいて、EMと同様の交代構造を持ちます。

多くの応用において、変分ベイズはギブスサンプリング法と同等の精度を持つ解をより高速に生成します。しかし、パラメータを反復的に更新するための方程式を導出するには、同等のギブスサンプリング方程式を導出する場合と比較して、多くの場合、膨大な作業量が必要になります。これは、概念的には非常に単純なモデルであっても当てはまります。以下では、パラメータが2つだけで潜在変数を持たない基本的な非階層モデルの例で示します。

数学的導出

[編集]

問題

[編集]

変分推論では、あるデータが与えられた場合の観測されない変数の集合の事後分布は、いわゆる変分分布によって近似される。${\displaystyle \mathbf {Z} =\{Z_{1}\dots Z_{n}\}}$ ${\displaystyle \mathbf {X} }$ ${\displaystyle Q(\mathbf {Z} ):}$

${\displaystyle P(\mathbf {Z} \mid \mathbf {X} )\approx Q(\mathbf {Z} ).}$

分布は、真の事後分布に類似するように意図されて選択された、より単純な形式の分布の族（例えば、ガウス分布の族）に属するように制限されます。 ${\displaystyle Q(\mathbf {Z} )}$ ${\displaystyle P(\mathbf {Z} \mid \mathbf {X} )}$ ${\displaystyle Q(\mathbf {Z} )}$ ${\displaystyle P(\mathbf {Z} \mid \mathbf {X} )}$

類似度（または非類似度）は非類似度関数で測定されるため、推論はを最小化する分布を選択することによって実行されます。 ${\displaystyle d(Q;P)}$ ${\displaystyle Q(\mathbf {Z} )}$ ${\displaystyle d(Q;P)}$

KLダイバージェンス

[編集]

最も一般的な変分ベイズ法では、類似度関数としてQとPのカルバック・ライブラー距離（KL距離）を用いる。この選択により、この最小化が扱いやすくなる。KL距離は以下のように定義される。

${\displaystyle D_{\mathrm {KL} }(Q\parallel P)\triangleq \sum _{\mathbf {Z} }Q(\mathbf {Z} )\log {\frac {Q(\mathbf {Z} )}{P(\mathbf {Z} \mid \mathbf {X} )}}.}$

QとPが予想とは逆になっていることに注意してください。この逆KLダイバージェンスの使用法は、概念的には期待値最大化アルゴリズムに似ています。（KLダイバージェンスを逆向きに使用すると、期待値伝播アルゴリズムが生成されます。）

扱いにくさ

[編集]

変分法は通常、以下の近似値を求めるために使用されます。

${\displaystyle P(\mathbf {Z} \mid \mathbf {X} )={\frac {P(\mathbf {X} \mid \mathbf {Z} )P(\mathbf {Z} )}{P(\mathbf {X} )}}={\frac {P(\mathbf {X} \mid \mathbf {Z} )P(\mathbf {Z} )}{\int _{\mathbf {Z} }P(\mathbf {X} ,\mathbf {Z} ')\,d\mathbf {Z} '}}}$

分母における上の周辺化を計算することは、例えば の探索空間が組み合わせ的に大きいため、典型的には困難です。そこで、 を用いて近似値を求めます。 ${\displaystyle \mathbf {Z} }$ ${\displaystyle P(\mathbf {X} )}$ ${\displaystyle \mathbf {Z} }$ ${\displaystyle Q(\mathbf {Z} )\approx P(\mathbf {Z} \mid \mathbf {X} )}$

証拠の下限

[編集]

主要記事:証拠の下限

とすると、上記のKLダイバージェンスは次のようにも書ける。 ${\displaystyle P(\mathbf {Z} \mid \mathbf {X} )={\frac {P(\mathbf {X} ,\mathbf {Z} )}{P(\mathbf {X} )}}}$

${\displaystyle {\begin{array}{rl}D_{\mathrm {KL} }(Q\parallel P)&=\sum _{\mathbf {Z} }Q(\mathbf {Z} )\left[\log {\frac {Q(\mathbf {Z} )}{P(\mathbf {Z} ,\mathbf {X} )}}+\log P(\mathbf {X} )\right]\\&=\sum _{\mathbf {Z} }Q(\mathbf {Z} )\left[\log Q(\mathbf {Z} )-\log P(\mathbf {Z} ,\mathbf {X} )\right]+\sum _{\mathbf {Z} }Q(\mathbf {Z} )\left[\log P(\mathbf {X} )\right]\end{array}}}$

は定数であり、は分布な ので、${\displaystyle P(\mathbf {X} )}$ ${\displaystyle \mathbf {Z} }$ ${\displaystyle \sum _{\mathbf {Z} }Q(\mathbf {Z} )=1}$ ${\displaystyle Q(\mathbf {Z} )}$

${\displaystyle D_{\mathrm {KL} }(Q\parallel P)=\sum _{\mathbf {Z} }Q(\mathbf {Z} )\left[\log Q(\mathbf {Z} )-\log P(\mathbf {Z} ,\mathbf {X} )\right]+\log P(\mathbf {X} )}$

これは、離散確率変数の期待値の定義によれば、次のように書ける。

${\displaystyle D_{\mathrm {KL} }(Q\parallel P)=\mathbb {E} _{\mathbf {Q} }\left[\log Q(\mathbf {Z} )-\log P(\mathbf {Z} ,\mathbf {X} )\right]+\log P(\mathbf {X} )}$

これを並べ替えると

${\displaystyle {\begin{array}{rl}\log P(\mathbf {X} )&=D_{\mathrm {KL} }(Q\parallel P)-\mathbb {E} _{\mathbf {Q} }\left[\log Q(\mathbf {Z} )-\log P(\mathbf {Z} ,\mathbf {X} )\right]\\&=D_{\mathrm {KL} }(Q\parallel P)+{\mathcal {L}}(Q)\end{array}}}$

対数エビデンスは に関して固定されているため、最終項を最大化するとから への KL ダイバージェンスが最小化されます。 を適切に選択することで、を計算し、 を最大化することが扱いやすくなります。したがって、事後分布 の解析的近似と、対数エビデンスの下限値（KL ダイバージェンスが非負であるため）の両方が得られます。 ${\displaystyle \log P(\mathbf {X} )}$ ${\displaystyle Q}$ ${\displaystyle {\mathcal {L}}(Q)}$ ${\displaystyle Q}$ ${\displaystyle P}$ ${\displaystyle Q}$ ${\displaystyle {\mathcal {L}}(Q)}$ ${\displaystyle Q}$ ${\displaystyle P(\mathbf {Z} \mid \mathbf {X} )}$ ${\displaystyle {\mathcal {L}}(Q)}$ ${\displaystyle \log P(\mathbf {X} )}$

下限は、熱力学的自由エネルギーと同様に、負のエネルギーにエントロピーを加えた値として表すこともできるため、（負の）変分自由エネルギーと呼ばれます。この用語は、データの対数証拠の下限（最悪のケース）であることを強調するために、証拠下限（Evidence Lower Bound ） （ELBOと略される）とも呼ばれます。 ${\displaystyle {\mathcal {L}}(Q)}$ ${\displaystyle \operatorname {E} _{Q}[\log P(\mathbf {Z} ,\mathbf {X} )]}$ ${\displaystyle Q}$ ${\displaystyle {\mathcal {L}}(Q)}$

証明

[編集]

KLダイバージェンスが特別なケースであるブレグマンダイバージェンスの一般化されたピタゴラスの定理により、次のことが示される: ^{[ 1 ] [ 2 ]}

ブレグマンダイバージェンスに対する一般化ピタゴラスの定理^{[ 2 ]}

${\displaystyle D_{\mathrm {KL} }(Q\parallel P)\geq D_{\mathrm {KL} }(Q\parallel Q^{*})+D_{\mathrm {KL} }(Q^{*}\parallel P),\forall Q^{*}\in {\mathcal {C}}}$

ここで は凸集合であり、次の条件を満たす場合等式が成り立ちます。 ${\displaystyle {\mathcal {C}}}$

${\displaystyle Q=Q^{*}\triangleq \arg \min _{Q\in {\mathcal {C}}}D_{\mathrm {KL} }(Q\parallel P).}$

この場合、グローバル最小値は次のように求められる。^{[ 1 ]}${\displaystyle Q^{*}(\mathbf {Z} )=q^{*}(\mathbf {Z} _{1}\mid \mathbf {Z} _{2})q^{*}(\mathbf {Z} _{2})=q^{*}(\mathbf {Z} _{2}\mid \mathbf {Z} _{1})q^{*}(\mathbf {Z} _{1}),}$ ${\displaystyle \mathbf {Z} =\{\mathbf {Z_{1}} ,\mathbf {Z_{2}} \},}$

${\displaystyle {\begin{array}{rl}q^{*}(\mathbf {Z} _{2})&={\frac {P(\mathbf {X} )}{\zeta (\mathbf {X} )}}{\frac {P(\mathbf {Z} _{2}\mid \mathbf {X} )}{\exp(D_{\mathrm {KL} }(q^{*}(\mathbf {Z} _{1}\mid \mathbf {Z} _{2})\parallel P(\mathbf {Z} _{1}\mid \mathbf {Z} _{2},\mathbf {X} )))}}\\&={\frac {1}{\zeta (\mathbf {X} )}}\exp \mathbb {E} _{q^{*}(\mathbf {Z} _{1}\mid \mathbf {Z} _{2})}\left(\log {\frac {P(\mathbf {Z} ,\mathbf {X} )}{q^{*}(\mathbf {Z} _{1}\mid \mathbf {Z} _{2})}}\right),\end{array}}}$

ここで正規化定数は次の通りです。

${\displaystyle {\begin{array}{rl}\zeta (\mathbf {X} )&=P(\mathbf {X} )\int _{\mathbf {Z} _{2}}{\frac {P(\mathbf {Z} _{2}\mid \mathbf {X} )}{\exp(D_{\mathrm {KL} }(q^{*}(\mathbf {Z} _{1}\mid \mathbf {Z} _{2})\parallel P(\mathbf {Z} _{1}\mid \mathbf {Z} _{2},\mathbf {X} )))}}\\&=\int _{\mathbf {Z} _{2}}\exp \mathbb {E} _{q^{*}(\mathbf {Z} _{1}\mid \mathbf {Z} _{2})}\left(\log {\frac {P(\mathbf {Z} ,\mathbf {X} )}{q^{*}(\mathbf {Z} _{1}\mid \mathbf {Z} _{2})}}\right).\end{array}}}$

この用語は実際には証拠下限値（ELBO ）^と呼ばれることが多い。これは[ ^{1 ]}上記の通りである。 ${\displaystyle \zeta (\mathbf {X} )}$ ${\displaystyle P(\mathbf {X} )\geq \zeta (\mathbf {X} )=\exp({\mathcal {L}}(Q^{*}))}$

と の役割を入れ替えることで、真のモデルの周辺分布の近似値とをそれぞれ反復的に計算することができます。この反復法は単調収束することが保証されていますが、^{[ 1 ]}収束した はの局所的最小値にすぎません。 ${\displaystyle \mathbf {Z} _{1}}$ ${\displaystyle \mathbf {Z} _{2},}$ ${\displaystyle q^{*}(\mathbf {Z} _{1})}$ ${\displaystyle q^{*}(\mathbf {Z} _{2})}$ ${\displaystyle P(\mathbf {Z} _{1}\mid \mathbf {X} )}$ ${\displaystyle P(\mathbf {Z} _{2}\mid \mathbf {X} ),}$ ${\displaystyle Q^{*}}$ ${\displaystyle D_{\mathrm {KL} }(Q\parallel P)}$

制約空間が独立空間内に限定されている場合、つまり上記の反復スキームは以下に示すようにいわゆる平均場近似になります。 ${\displaystyle {\mathcal {C}}}$ ${\displaystyle q^{*}(\mathbf {Z} _{1}\mid \mathbf {Z} _{2})=q^{*}(\mathbf {Z_{1}} ),}$ ${\displaystyle Q^{*}(\mathbf {Z} )=q^{*}(\mathbf {Z} _{1})q^{*}(\mathbf {Z} _{2}),}$

平均場近似

[編集]

変分分布は通常、潜在変数のある分割、すなわち潜在変数のある分割に対して因数分解すると仮定される。 ${\displaystyle Q(\mathbf {Z} )}$ ${\displaystyle \mathbf {Z} }$ ${\displaystyle \mathbf {Z} _{1}\dots \mathbf {Z} _{M}}$

${\displaystyle Q(\mathbf {Z} )=\prod _{i=1}^{M}q_{i}(\mathbf {Z} _{i}\mid \mathbf {X} )}$

変分法（「変分ベイズ」と呼ばれる）を用いると、各因子に対する「最良の」分布（前述のKLダイバージェンスを最小化する分布の観点から）は次式を満たすことが示される。 ^{[ 3 ]}${\displaystyle q_{j}^{*}}$ ${\displaystyle q_{j}}$

${\displaystyle q_{j}^{*}(\mathbf {Z} _{j}\mid \mathbf {X} )={\frac {e^{\operatorname {E} _{q_{-j}^{*}}[\ln p(\mathbf {Z} ,\mathbf {X} )]}}{\int e^{\operatorname {E} _{q_{-j}^{*}}[\ln p(\mathbf {Z} ,\mathbf {X} )]}\,d\mathbf {Z} _{j}}}}$

ここで、はデータと潜在変数の結合確率の対数の期待値であり、分割に含まれないすべての変数についてとられたものである。分布の導出については、^{[ 4 ]}の補題4.1を参照のこと。 ${\displaystyle \operatorname {E} _{q_{-j}^{*}}[\ln p(\mathbf {Z} ,\mathbf {X} )]}$ ${\displaystyle q^{*}}$ ${\displaystyle q_{j}^{*}(\mathbf {Z} _{j}\mid \mathbf {X} )}$

実際には、通常は対数で作業します。つまり、

${\displaystyle \ln q_{j}^{*}(\mathbf {Z} _{j}\mid \mathbf {X} )=\operatorname {E} _{q_{-j}^{*}}[\ln p(\mathbf {Z} ,\mathbf {X} )]+{\text{constant}}}$

上記の式の定数は、正規化定数（上記の式の の分母）に関連しており、式の残りの部分は通常、既知の分布の種類（ガウス、ガンマなど）として認識できるため、通常は検査によって復元されます。 ${\displaystyle q_{j}^{*}}$

期待値の特性を用いると、式は通常、潜在変数上の事前分布の固定された超パラメータと、現在のパーティションにない潜在変数（すなわち、 に含まれない潜在変数）の期待値（および場合によっては分散 などの高次モーメント）の関数に簡略化できる。これにより、1つのパーティションの変数上の分布のパラメータと他のパーティションの変数の期待値との間に循環依存関係が生じる。これは当然、EM（期待値最大化アルゴリズム）によく似た反復アルゴリズムを示唆しており、このアルゴリズムでは、潜在変数の期待値（および場合によっては高次モーメント）が何らかの方法（おそらくランダム）で初期化され、次に各分布のパラメータが期待値の現在の値を用いて順に計算され、その後、新たに計算された分布の期待値が計算されたパラメータに従って適切に設定される。この種のアルゴリズムは に収束することが保証されている。^{[ 5 ]}${\displaystyle \operatorname {E} _{q_{-j}^{*}}[\ln p(\mathbf {Z} ,\mathbf {X} )]}$ ${\displaystyle \mathbf {Z} _{j}}$

言い換えれば、変数の各パーティションについて、パーティションの変数に対する分布の式を簡略化し、問題の変数に対する分布の機能的依存性を調べることによって、分布の族を通常は決定できます（これにより定数の値が決定されます）。分布のパラメータの式は、事前分布のハイパーパラメータ（既知の定数）で表されますが、他のパーティションの変数の関数の期待値でも表されます。通常、これらの期待値は、変数自体の期待値（つまり、平均 ）の関数に簡略化できます。場合によっては、変数の二乗の期待値（変数の分散に関連します）や、より高次のべき乗の期待値（つまり、高次のモーメント）も表示されます。ほとんどの場合、他の変数の分布は既知の族からのものであり、関連する期待値の式を調べることができます。ただし、これらの式は分布のパラメータに依存し、パラメータは他の変数に関する期待値に依存します。その結果、各変数の分布のパラメータの式は、変数間に相互非線形依存関係を持つ一連の方程式として表すことができます。通常、この方程式系を直接解くことはできません。しかし、前述のように、これらの依存関係は単純な反復アルゴリズムを示唆しており、ほとんどの場合収束することが保証されています。例を挙げてこのプロセスをより明確にします。

変分推論の双対性公式

[編集]

双対公式による座標上昇変分推論アルゴリズムの図解^{[ 4 ]}

次の定理は変分推論の双対性公式と呼ばれています。^{[ 4 ]}これは変分ベイズ法で使用される変分分布のいくつかの重要な性質を説明しています。

定理2つの確率空間 とを考える。そしてとなる共通の支配確率測度が存在すると仮定する。 を満たす任意の実数値確率変数を とおく。このとき、次の等式が成立する。 ${\displaystyle (\Theta ,{\mathcal {F}},P)}$ ${\displaystyle (\Theta ,{\mathcal {F}},Q)}$ ${\displaystyle Q\ll P}$ ${\displaystyle \lambda }$ ${\displaystyle P\ll \lambda }$ ${\displaystyle Q\ll \lambda }$ ${\displaystyle h}$ ${\displaystyle (\Theta ,{\mathcal {F}},P)}$ ${\displaystyle h\in L_{1}(P)}$

${\displaystyle \log E_{P}[\exp h]={\text{sup}}_{Q\ll P}\{E_{Q}[h]-D_{\text{KL}}(Q\parallel P)\}.}$

さらに、右辺の上限は、次の式が成り立つ 場合にのみ達成される。

${\displaystyle {\frac {q(\theta )}{p(\theta )}}={\frac {\exp h(\theta )}{E_{P}[\exp h]}},}$

は確率測度 に関してほぼ確実に成り立ちます。ここで、 と はそれぞれ確率測度 と のラドン・ニコディム微分を表します。 ${\displaystyle Q}$ ${\displaystyle p(\theta )=dP/d\lambda }$ ${\displaystyle q(\theta )=dQ/d\lambda }$ ${\displaystyle P}$ ${\displaystyle Q}$ ${\displaystyle \lambda }$

基本的な例

[編集]

ガウス分布からのiid観測値の集合から成り、平均と分散が不明な単純な非階層的ベイズモデルを考えてみましょう。^{[ 6 ]}以下では、変分ベイズ法の仕組みを説明するために、このモデルを詳細に検討します。

数学的な便宜上、以下の例では、分散そのものではなく、精度、つまり分散の逆数（または多変量ガウス分布の場合は共分散行列の逆数）を基準とします。（理論的な観点からは、精度と分散は1対1で対応しているため、両者は 同等です。）

数学モデル

[編集]

未知の平均値と精度に共役事前分布を適用します。つまり、平均値はガウス分布に従い、精度はガンマ分布に従います。言い換えると、 ${\displaystyle \mu }$ ${\displaystyle \tau }$

${\displaystyle {\begin{aligned}\tau &\sim \operatorname {Gamma} (a_{0},b_{0})\\\mu |\tau &\sim {\mathcal {N}}(\mu _{0},(\lambda _{0}\tau )^{-1})\\\{x_{1},\dots ,x_{N}\}&\sim {\mathcal {N}}(\mu ,\tau ^{-1})\\N&={\text{number of data points}}\end{aligned}}}$

事前分布におけるハイパーパラメータ と は固定値で与えられます。これらを小さな正の値に設定することで、およびの事前分布が不明であることを示す広い事前分布を与えることができます。 ${\displaystyle \mu _{0},\lambda _{0},a_{0}}$ ${\displaystyle b_{0}}$ ${\displaystyle \mu }$ ${\displaystyle \tau }$

与えられたデータポイントからパラメータの事後分布を推定することが目的です。${\displaystyle N}$ ${\displaystyle \mathbf {X} =\{x_{1},\ldots ,x_{N}\}}$ ${\displaystyle q(\mu ,\tau )=p(\mu ,\tau \mid x_{1},\ldots ,x_{N})}$ ${\displaystyle \mu }$ ${\displaystyle \tau .}$

結合確率

[編集]

すべての変数の結合確率は次のように書き直すことができる。

${\displaystyle p(\mathbf {X} ,\mu ,\tau )=p(\mathbf {X} \mid \mu ,\tau )p(\mu \mid \tau )p(\tau )}$

個々の要因は

${\displaystyle {\begin{aligned}p(\mathbf {X} \mid \mu ,\tau )&=\prod _{n=1}^{N}{\mathcal {N}}(x_{n}\mid \mu ,\tau ^{-1})\\p(\mu \mid \tau )&={\mathcal {N}}\left(\mu \mid \mu _{0},(\lambda _{0}\tau )^{-1}\right)\\p(\tau )&=\operatorname {Gamma} (\tau \mid a_{0},b_{0})\end{aligned}}}$

どこ

${\displaystyle {\begin{aligned}{\mathcal {N}}(x\mid \mu ,\sigma ^{2})&={\frac {1}{\sqrt {2\pi \sigma ^{2}}}}e^{\frac {-(x-\mu )^{2}}{2\sigma ^{2}}}\\\operatorname {Gamma} (\tau \mid a,b)&={\frac {1}{\Gamma (a)}}b^{a}\tau ^{a-1}e^{-b\tau }\end{aligned}}}$

因数分解近似

[編集]

と仮定します。つまり、事後分布はとについて独立な因子に分解されます。この種の仮定は変分ベイズ法の基礎となっています。真の事後分布は実際にはこのようには分解されません（実際、この単純なケースでは、ガウスガンマ分布であることが知られています）。したがって、得られる結果は近似値となります。 ${\displaystyle q(\mu ,\tau )=q(\mu )q(\tau )}$ ${\displaystyle \mu }$ ${\displaystyle \tau }$

q ( μ )の導出

[編集]

それから

${\displaystyle {\begin{aligned}\ln q_{\mu }^{*}(\mu )&=\operatorname {E} _{\tau }\left[\ln p(\mathbf {X} \mid \mu ,\tau )+\ln p(\mu \mid \tau )+\ln p(\tau )\right]+C\\&=\operatorname {E} _{\tau }\left[\ln p(\mathbf {X} \mid \mu ,\tau )\right]+\operatorname {E} _{\tau }\left[\ln p(\mu \mid \tau )\right]+\operatorname {E} _{\tau }\left[\ln p(\tau )\right]+C\\&=\operatorname {E} _{\tau }\left[\ln \prod _{n=1}^{N}{\mathcal {N}}\left(x_{n}\mid \mu ,\tau ^{-1}\right)\right]+\operatorname {E} _{\tau }\left[\ln {\mathcal {N}}\left(\mu \mid \mu _{0},(\lambda _{0}\tau )^{-1}\right)\right]+C_{2}\\&=\operatorname {E} _{\tau }\left[\ln \prod _{n=1}^{N}{\sqrt {\frac {\tau }{2\pi }}}e^{-{\frac {(x_{n}-\mu )^{2}\tau }{2}}}\right]+\operatorname {E} _{\tau }\left[\ln {\sqrt {\frac {\lambda _{0}\tau }{2\pi }}}e^{-{\frac {(\mu -\mu _{0})^{2}\lambda _{0}\tau }{2}}}\right]+C_{2}\\&=\operatorname {E} _{\tau }\left[\sum _{n=1}^{N}\left({\frac {1}{2}}(\ln \tau -\ln 2\pi )-{\frac {(x_{n}-\mu )^{2}\tau }{2}}\right)\right]+\operatorname {E} _{\tau }\left[{\frac {1}{2}}(\ln \lambda _{0}+\ln \tau -\ln 2\pi )-{\frac {(\mu -\mu _{0})^{2}\lambda _{0}\tau }{2}}\right]+C_{2}\\&=\operatorname {E} _{\tau }\left[\sum _{n=1}^{N}-{\frac {(x_{n}-\mu )^{2}\tau }{2}}\right]+\operatorname {E} _{\tau }\left[-{\frac {(\mu -\mu _{0})^{2}\lambda _{0}\tau }{2}}\right]+\operatorname {E} _{\tau }\left[\sum _{n=1}^{N}{\frac {1}{2}}(\ln \tau -\ln 2\pi )\right]+\operatorname {E} _{\tau }\left[{\frac {1}{2}}(\ln \lambda _{0}+\ln \tau -\ln 2\pi )\right]+C_{2}\\&=\operatorname {E} _{\tau }\left[\sum _{n=1}^{N}-{\frac {(x_{n}-\mu )^{2}\tau }{2}}\right]+\operatorname {E} _{\tau }\left[-{\frac {(\mu -\mu _{0})^{2}\lambda _{0}\tau }{2}}\right]+C_{3}\\&=-{\frac {\operatorname {E} _{\tau }[\tau ]}{2}}\left\{\sum _{n=1}^{N}(x_{n}-\mu )^{2}+\lambda _{0}(\mu -\mu _{0})^{2}\right\}+C_{3}\end{aligned}}}$

上記の導出において、、、はに関して定数である値を指します。項 はの関数ではなく、 の値に関わらず同じ値を持つことに注意してください。したがって、3行目では、これを最後の定数項に吸収することができます。7行目でも同じことを行います。 ${\displaystyle C}$ ${\displaystyle C_{2}}$ ${\displaystyle C_{3}}$ ${\displaystyle \mu }$ ${\displaystyle \operatorname {E} _{\tau }[\ln p(\tau )]}$ ${\displaystyle \mu }$ ${\displaystyle \mu }$

最後の行は の単純な2次多項式です。これは の対数なので、自体はガウス分布であることがわかります。 ${\displaystyle \mu }$ ${\displaystyle q_{\mu }^{*}(\mu )}$ ${\displaystyle q_{\mu }^{*}(\mu )}$

ある程度の面倒な計算（中括弧内の四角形を展開し、 と を含む項を分離してグループ化し、上で平方完成させる）を行うと、ガウス分布のパラメータを導くことができます。 ${\displaystyle \mu }$ ${\displaystyle \mu ^{2}}$ ${\displaystyle \mu }$

${\displaystyle {\begin{aligned}\ln q_{\mu }^{*}(\mu )&=-{\frac {\operatorname {E} _{\tau }[\tau ]}{2}}\left\{\sum _{n=1}^{N}(x_{n}-\mu )^{2}+\lambda _{0}(\mu -\mu _{0})^{2}\right\}+C_{3}\\&=-{\frac {\operatorname {E} _{\tau }[\tau ]}{2}}\left\{\sum _{n=1}^{N}(x_{n}^{2}-2x_{n}\mu +\mu ^{2})+\lambda _{0}(\mu ^{2}-2\mu _{0}\mu +\mu _{0}^{2})\right\}+C_{3}\\&=-{\frac {\operatorname {E} _{\tau }[\tau ]}{2}}\left\{\left(\sum _{n=1}^{N}x_{n}^{2}\right)-2\left(\sum _{n=1}^{N}x_{n}\right)\mu +\left(\sum _{n=1}^{N}\mu ^{2}\right)+\lambda _{0}\mu ^{2}-2\lambda _{0}\mu _{0}\mu +\lambda _{0}\mu _{0}^{2}\right\}+C_{3}\\&=-{\frac {\operatorname {E} _{\tau }[\tau ]}{2}}\left\{(\lambda _{0}+N)\mu ^{2}-2\left(\lambda _{0}\mu _{0}+\sum _{n=1}^{N}x_{n}\right)\mu +\left(\sum _{n=1}^{N}x_{n}^{2}\right)+\lambda _{0}\mu _{0}^{2}\right\}+C_{3}\\&=-{\frac {\operatorname {E} _{\tau }[\tau ]}{2}}\left\{(\lambda _{0}+N)\mu ^{2}-2\left(\lambda _{0}\mu _{0}+\sum _{n=1}^{N}x_{n}\right)\mu \right\}+C_{4}\\&=-{\frac {\operatorname {E} _{\tau }[\tau ]}{2}}\left\{(\lambda _{0}+N)\mu ^{2}-2\left({\frac {\lambda _{0}\mu _{0}+\sum _{n=1}^{N}x_{n}}{\lambda _{0}+N}}\right)(\lambda _{0}+N)\mu \right\}+C_{4}\\&=-{\frac {\operatorname {E} _{\tau }[\tau ]}{2}}\left\{(\lambda _{0}+N)\left(\mu ^{2}-2\left({\frac {\lambda _{0}\mu _{0}+\sum _{n=1}^{N}x_{n}}{\lambda _{0}+N}}\right)\mu \right)\right\}+C_{4}\\&=-{\frac {\operatorname {E} _{\tau }[\tau ]}{2}}\left\{(\lambda _{0}+N)\left(\mu ^{2}-2\left({\frac {\lambda _{0}\mu _{0}+\sum _{n=1}^{N}x_{n}}{\lambda _{0}+N}}\right)\mu +\left({\frac {\lambda _{0}\mu _{0}+\sum _{n=1}^{N}x_{n}}{\lambda _{0}+N}}\right)^{2}-\left({\frac {\lambda _{0}\mu _{0}+\sum _{n=1}^{N}x_{n}}{\lambda _{0}+N}}\right)^{2}\right)\right\}+C_{4}\\&=-{\frac {\operatorname {E} _{\tau }[\tau ]}{2}}\left\{(\lambda _{0}+N)\left(\mu ^{2}-2\left({\frac {\lambda _{0}\mu _{0}+\sum _{n=1}^{N}x_{n}}{\lambda _{0}+N}}\right)\mu +\left({\frac {\lambda _{0}\mu _{0}+\sum _{n=1}^{N}x_{n}}{\lambda _{0}+N}}\right)^{2}\right)\right\}+C_{5}\\&=-{\frac {\operatorname {E} _{\tau }[\tau ]}{2}}\left\{(\lambda _{0}+N)\left(\mu -{\frac {\lambda _{0}\mu _{0}+\sum _{n=1}^{N}x_{n}}{\lambda _{0}+N}}\right)^{2}\right\}+C_{5}\\&=-{\frac {1}{2}}(\lambda _{0}+N)\operatorname {E} _{\tau }[\tau ]\left(\mu -{\frac {\lambda _{0}\mu _{0}+\sum _{n=1}^{N}x_{n}}{\lambda _{0}+N}}\right)^{2}+C_{5}\end{aligned}}}$

上記の手順はすべて、 2 つの二次方程式の和の公式を使用することで短縮できることに注意してください。

言い換えると：

${\displaystyle {\begin{aligned}q_{\mu }^{*}(\mu )&\sim {\mathcal {N}}(\mu \mid \mu _{N},\lambda _{N}^{-1})\\\mu _{N}&={\frac {\lambda _{0}\mu _{0}+N{\bar {x}}}{\lambda _{0}+N}}\\\lambda _{N}&=(\lambda _{0}+N)\operatorname {E} _{\tau }[\tau ]\\{\bar {x}}&={\frac {1}{N}}\sum _{n=1}^{N}x_{n}\end{aligned}}}$

q(τ)の導出

[編集]

の導出は上記と同様ですが、簡潔にするために詳細の一部は省略します。 ${\displaystyle q_{\tau }^{*}(\tau )}$

${\displaystyle {\begin{aligned}\ln q_{\tau }^{*}(\tau )&=\operatorname {E} _{\mu }[\ln p(\mathbf {X} \mid \mu ,\tau )+\ln p(\mu \mid \tau )]+\ln p(\tau )+{\text{constant}}\\&=(a_{0}-1)\ln \tau -b_{0}\tau +{\frac {1}{2}}\ln \tau +{\frac {N}{2}}\ln \tau -{\frac {\tau }{2}}\operatorname {E} _{\mu }\left[\sum _{n=1}^{N}(x_{n}-\mu )^{2}+\lambda _{0}(\mu -\mu _{0})^{2}\right]+{\text{constant}}\end{aligned}}}$

両辺をべき乗すると、ガンマ分布であることがわかります。具体的には、 ${\displaystyle q_{\tau }^{*}(\tau )}$

${\displaystyle {\begin{aligned}q_{\tau }^{*}(\tau )&\sim \operatorname {Gamma} (\tau \mid a_{N},b_{N})\\a_{N}&=a_{0}+{\frac {N+1}{2}}\\b_{N}&=b_{0}+{\frac {1}{2}}\operatorname {E} _{\mu }\left[\sum _{n=1}^{N}(x_{n}-\mu )^{2}+\lambda _{0}(\mu -\mu _{0})^{2}\right]\end{aligned}}}$

パラメータを計算するアルゴリズム

[編集]

前のセクションの結論を要約してみましょう。

${\displaystyle {\begin{aligned}q_{\mu }^{*}(\mu )&\sim {\mathcal {N}}(\mu \mid \mu _{N},\lambda _{N}^{-1})\\\mu _{N}&={\frac {\lambda _{0}\mu _{0}+N{\bar {x}}}{\lambda _{0}+N}}\\\lambda _{N}&=(\lambda _{0}+N)\operatorname {E} _{\tau }[\tau ]\\{\bar {x}}&={\frac {1}{N}}\sum _{n=1}^{N}x_{n}\end{aligned}}}$

そして

${\displaystyle {\begin{aligned}q_{\tau }^{*}(\tau )&\sim \operatorname {Gamma} (\tau \mid a_{N},b_{N})\\a_{N}&=a_{0}+{\frac {N+1}{2}}\\b_{N}&=b_{0}+{\frac {1}{2}}\operatorname {E} _{\mu }\left[\sum _{n=1}^{N}(x_{n}-\mu )^{2}+\lambda _{0}(\mu -\mu _{0})^{2}\right]\end{aligned}}}$

いずれの場合も、一方の変数に関する分布のパラメータは、もう一方の変数に関する期待値に依存します。ガウス分布とガンマ分布のモーメントの期待値を求める標準的な公式を用いて、期待値を展開することができます。

${\displaystyle {\begin{aligned}\operatorname {E} [\tau \mid a_{N},b_{N}]&={\frac {a_{N}}{b_{N}}}\\\operatorname {E} \left[\mu \mid \mu _{N},\lambda _{N}^{-1}\right]&=\mu _{N}\\\operatorname {E} \left[X^{2}\right]&=\operatorname {Var} (X)+(\operatorname {E} [X])^{2}\\\operatorname {E} \left[\mu ^{2}\mid \mu _{N},\lambda _{N}^{-1}\right]&=\lambda _{N}^{-1}+\mu _{N}^{2}\end{aligned}}}$

これらの式を上記の式に適用するのはほとんどの場合簡単ですが、 の式の場合はさらに手間がかかります。 ${\displaystyle b_{N}}$

${\displaystyle {\begin{aligned}b_{N}&=b_{0}+{\frac {1}{2}}\operatorname {E} _{\mu }\left[\sum _{n=1}^{N}(x_{n}-\mu )^{2}+\lambda _{0}(\mu -\mu _{0})^{2}\right]\\&=b_{0}+{\frac {1}{2}}\operatorname {E} _{\mu }\left[(\lambda _{0}+N)\mu ^{2}-2\left(\lambda _{0}\mu _{0}+\sum _{n=1}^{N}x_{n}\right)\mu +\left(\sum _{n=1}^{N}x_{n}^{2}\right)+\lambda _{0}\mu _{0}^{2}\right]\\&=b_{0}+{\frac {1}{2}}\left[(\lambda _{0}+N)\operatorname {E} _{\mu }[\mu ^{2}]-2\left(\lambda _{0}\mu _{0}+\sum _{n=1}^{N}x_{n}\right)\operatorname {E} _{\mu }[\mu ]+\left(\sum _{n=1}^{N}x_{n}^{2}\right)+\lambda _{0}\mu _{0}^{2}\right]\\&=b_{0}+{\frac {1}{2}}\left[(\lambda _{0}+N)\left(\lambda _{N}^{-1}+\mu _{N}^{2}\right)-2\left(\lambda _{0}\mu _{0}+\sum _{n=1}^{N}x_{n}\right)\mu _{N}+\left(\sum _{n=1}^{N}x_{n}^{2}\right)+\lambda _{0}\mu _{0}^{2}\right]\\\end{aligned}}}$

次に、期待値なしでパラメータ方程式を次のように記述します。

${\displaystyle {\begin{aligned}\mu _{N}&={\frac {\lambda _{0}\mu _{0}+N{\bar {x}}}{\lambda _{0}+N}}\\\lambda _{N}&=(\lambda _{0}+N){\frac {a_{N}}{b_{N}}}\\{\bar {x}}&={\frac {1}{N}}\sum _{n=1}^{N}x_{n}\\a_{N}&=a_{0}+{\frac {N+1}{2}}\\b_{N}&=b_{0}+{\frac {1}{2}}\left[(\lambda _{0}+N)\left(\lambda _{N}^{-1}+\mu _{N}^{2}\right)-2\left(\lambda _{0}\mu _{0}+\sum _{n=1}^{N}x_{n}\right)\mu _{N}+\left(\sum _{n=1}^{N}x_{n}^{2}\right)+\lambda _{0}\mu _{0}^{2}\right]\end{aligned}}}$

と の式の間には循環依存関係があることに注意してください。これは自然にEMのようなアルゴリズムを示唆しています。 ${\displaystyle \lambda _{N}}$ ${\displaystyle b_{N}}$

1. 計算してこれらの値を使用して計算し、${\displaystyle \sum _{n=1}^{N}x_{n}}$ ${\displaystyle \sum _{n=1}^{N}x_{n}^{2}.}$ ${\displaystyle \mu _{N}}$ ${\displaystyle a_{N}.}$
2. 任意の値に初期化します。${\displaystyle \lambda _{N}}$
3. の現在の値と他のパラメータの既知の値を使用して、 を計算します。${\displaystyle \lambda _{N},}$ ${\displaystyle b_{N}}$
4. の現在の値と他のパラメータの既知の値を使用して、 を計算します。${\displaystyle b_{N},}$ ${\displaystyle \lambda _{N}}$
5. 収束するまで（つまり、どちらの値も少ししか変化しなくなるまで）最後の 2 つの手順を繰り返します。

次に、事後パラメータの近似分布のハイパーパラメータの値を取得し、これを使用して事後分布の必要なプロパティ（平均と分散、95% 最高密度領域（全確率の 95% を含む最小の間隔）など）を計算できます。

このアルゴリズムは局所的最大値に収束することが保証されていることがわかります。

また、事後分布は対応する事前分布と同じ形になることにも注意してください。これは仮定していません。私たちが行った唯一の仮定は、分布が因数分解され、分布の形が自然に従うというものでした。事後分布が事前分布と同じ形になるという事実は偶然ではなく、事前分布が指数族に属する場合、つまり標準的な分布のほとんどが指数族に属する場合の一般的な結果であることがわかります（以下を参照）。

さらなる議論

[編集]

ステップバイステップのレシピ

[編集]

上記の例は、与えられたベイジアンネットワークにおける事後確率密度に対する変分ベイジアン近似を導出する方法を示しています。

1. ネットワークをグラフィカルモデルで記述し、観測変数（データ）と非観測変数（パラメータと潜在変数）およびそれらの条件付き確率分布を特定します。変分ベイズ法は事後確率の近似値を構築します。この近似値は、因子分解分布、すなわち非観測変数の互いに素な部分集合上の2つ以上の独立した分布の積であるという基本的な特性を持ちます。${\displaystyle \mathbf {X} }$ ${\displaystyle {\boldsymbol {\Theta }}}$ ${\displaystyle \mathbf {Z} }$ ${\displaystyle p(\mathbf {Z} ,{\boldsymbol {\Theta }}\mid \mathbf {X} )}$
2. 観測されない変数を2つ以上のサブセットに分割し、それらから独立因子を導出します。これを行うための普遍的な手順はありません。サブセットが多すぎると近似値が悪くなり、少なすぎると変分ベイズ法全体が扱いにくくなります。通常、最初の分割はパラメータと潜在変数を分離することです。多くの場合、これだけで扱いやすい結果が得られます。分割を と仮定します。${\displaystyle \mathbf {Z} _{1},\ldots ,\mathbf {Z} _{M}}$
3. 与えられたパーティションに対して、基本方程式を使用して最も近似する分布の式を書きます。${\displaystyle \mathbf {Z} _{j}}$ ${\displaystyle q_{j}^{*}(\mathbf {Z} _{j}\mid \mathbf {X} )}$ ${\displaystyle \ln q_{j}^{*}(\mathbf {Z} _{j}\mid \mathbf {X} )=\operatorname {E} _{i\neq j}[\ln p(\mathbf {Z} ,\mathbf {X} )]+{\text{constant}}}$
4. グラフィカルモデルを用いて結合確率分布の式を記入してください。どの変数も含まない条件付き分布は無視できます。それらは定数項に折り込まれます。${\displaystyle \mathbf {Z} _{j}}$
5. 上記の例に従って、式を簡略化し、期待値演算子を適用します。理想的には、これは に含まれない変数の基本関数の期待値（例えば、1次または2次モーメントの生のモーメント、対数の期待値など）に簡略化されるはずです。変分ベイズ法が適切に機能するためには、これらの期待値は一般に、これらの変数の分布のパラメータおよび／またはハイパーパラメータの関数として解析的に表現できる必要があります。いずれの場合も、これらの期待値項は、現在のパーティション内の変数に関して定数です。${\displaystyle \mathbf {Z} _{j}}$
6. 現在のパーティション内の変数に関する式の関数形は、分布の種類を示します。特に、式をべき乗すると、分布の確率密度関数(PDF) が生成されます（または少なくとも、正規化定数が不明な場合に、それに比例する何かが生成されます）。この方法全体を扱いやすくするためには、関数形が既知の分布に属するものと認識できる必要があります。式を既知の分布のPDFと一致する形式に変換するには、かなりの数学的操作が必要になる場合があります。これが可能であれば、正規化定数は定義により復元でき、式の適切な部分を抽出することで、既知の分布のパラメータの式を導出できます。
7. すべての期待値を現在のパーティションに含まれない変数の関数で解析的に置き換えることができ、PDF を既知の分布との識別を可能にする形式にすると、結果は、他のパーティション内の変数のパラメータの関数として最適パラメータの値を表す一連の方程式になります。
8. この手順をすべてのパーティションに適用できる場合、結果はすべてのパラメータの最適値を指定する相互にリンクされた方程式のセットになります。
9. 次に、期待値最大化（EM）型の手順を適用し、各パラメータの初期値を選択し、一連のステップを反復します。各ステップでは、方程式を循環的に処理し、各パラメータを順に更新します。この手順は収束することが保証されています。

最も重要なポイント

[編集]

数学的な操作が多岐にわたるため、全体像を見失いがちです。重要なのは以下の点です。

1. 変分ベイズの考え方は、データが与えられた場合、観測されない変数（パラメータと潜在変数）のセットの事後確率の解析的近似を構築することです。つまり、解の形式はギブス サンプリングなどの他のベイズ推論方法と似ています。つまり、変数についてわかっていることすべてを記述しようとする分布です。他のベイズ法と同様に（ただし、期待最大化（EM）法や他の最大尤度法などとは異なり）、観測されない変数の両方のタイプ（パラメータと潜在変数）は同じ、つまりランダム変数として扱われます。その後、変数の推定値は、分布の平均を計算して単一の点推定値を取得したり、信頼区間や最高密度領域などを導出するなどの標準的なベイズの方法で導出できます。
2. 「解析的近似」とは、事後分布について式を記述できることを意味します。この式は一般に、よく知られた確率分布の積で構成され、各分布は観測されない変数の集合を因数分解します（つまり、観測データが与えられた場合、他の変数から条件付きで独立です）。この式は真の事後分布ではなく、その近似値です。特に、平均や分散といった観測されない変数の最小モーメントにおいては、一般的にかなり近い値を示します。
3. すべての数学的操作の結果は、(1) 因子を構成する確率分布の恒等式、および (2) これらの分布のパラメータに対する相互に依存する式です。これらのパラメータの実際の値は、EM法によく似た交代反復手順を通じて数値的に計算されます。

期待最大化（EM）と比較して

[編集]

変分ベイズ法（VB）はしばしば期待最大化法（EM）と比較されます。実際の数値計算手順は非常に似ており、どちらも最適なパラメータ値へと収束していく反復手順を交互に繰り返していくものです。それぞれの手順を導出するための初期ステップも漠然と似ており、どちらも確率密度の式から始まり、かなりの量の数学的操作を伴います。

しかし、いくつかの違いがあります。最も重要なのは、何が計算されるかということです。

• EMは、「パラメータ」として分類できる確率変数の事後分布の点推定値を計算しますが、これは潜在変数の実際の事後分布の推定値のみとなります（少なくとも「ソフトEM」においては、また多くの場合、潜在変数が離散的な場合のみ）。計算される点推定値はこれらのパラメータの最頻値であり、それ以外の情報は得られません。
• 一方、VB は、すべての変数（パラメータと潜在変数の両方）の実際の事後分布の推定値を計算します。点推定値を導出する必要がある場合、ベイズ推論では通常のように、モードではなく平均が一般的に使用されます。これに付随して、VB で計算されるパラメータは、 EM のパラメータと同じ重要性を持ちません。EM は、ベイズ ネットワーク自体のパラメータの最適値を計算します。VB は、ベイズ ネットワークのパラメータと潜在変数を近似するために使用される分布のパラメータの最適値を計算します。たとえば、一般的なガウス混合モデルには、混合コンポーネントのそれぞれの平均と分散のパラメータがあります。EM はこれらのパラメータの最適値を直接推定します。ただし、VB は、最初にこれらのパラメータに分布を適合させます（通常は事前分布の形式、たとえば正規尺度の逆ガンマ分布）。次に、この事前分布のパラメータの値、つまり基本的にハイパーパラメータを計算します。この場合、VB は、コンポーネントの平均と分散の結合分布を記述する正規スケール逆ガンマ分布の 4 つのパラメータの最適推定値を計算します。

より複雑な例

[編集]

プレート記法を用いたベイズ混合ガウスモデル。小さな四角は固定パラメータ、大きな円はランダム変数を示します。塗りつぶされた図形は既知の値を示します。[K] はサイズKのベクトル、[ D , D ] はサイズD × Dの行列を意味します。Kのみの場合は、K 個の結果を持つカテゴリ変数を意味します。zから始まり、最後に横棒で終わる波線はスイッチを示します。この変数の値は、他の入力変数に対して、サイズKの可能な値の配列からどの値を使用するかを選択します。

以下のように記述されるベイズガウス混合モデルを想像してください。^{[ 3 ]}

${\displaystyle {\begin{aligned}\mathbf {\pi } &\sim \operatorname {SymDir} (K,\alpha _{0})\\\mathbf {\Lambda } _{i=1\dots K}&\sim {\mathcal {W}}(\mathbf {W} _{0},\nu _{0})\\\mathbf {\mu } _{i=1\dots K}&\sim {\mathcal {N}}(\mathbf {\mu } _{0},(\beta _{0}\mathbf {\Lambda } _{i})^{-1})\\\mathbf {z} [i=1\dots N]&\sim \operatorname {Mult} (1,\mathbf {\pi } )\\\mathbf {x} _{i=1\dots N}&\sim {\mathcal {N}}(\mathbf {\mu } _{z_{i}},{\mathbf {\Lambda } _{z_{i}}}^{-1})\\K&={\text{number of mixing components}}\\N&={\text{number of data points}}\end{aligned}}}$

注記：

• SymDir() は、各成分のハイパーパラメータが に設定された次元 の対称ディリクレ分布です。ディリクレ分布は、カテゴリ分布または多項分布の共役事前分布です。${\displaystyle K}$ ${\displaystyle \alpha _{0}}$
• ${\displaystyle {\mathcal {W}}()}$ はウィシャート分布であり、多変量ガウス分布の精度行列（逆共分散行列）の共役事前分布です。
• Mult() は、単一の観測値に対する多項分布です（カテゴリ分布に相当します）。状態空間は「K 個の中の 1 個」の表現、つまり、要素の 1 つが 1（観測値の ID を表す）で、他のすべての要素が 0 である -次元ベクトルです。${\displaystyle K}$
• ${\displaystyle {\mathcal {N}}()}$ はガウス分布であり、この場合は具体的には多変量ガウス分布です。

上記の変数の解釈は次のとおりです。

• ${\displaystyle \mathbf {X} =\{\mathbf {x} _{1},\dots ,\mathbf {x} _{N}\}}$ は、多変量ガウス分布に従って分布する次元ベクトルであるデータポイントの集合です。${\displaystyle N}$ ${\displaystyle D}$
• ${\displaystyle \mathbf {Z} =\{\mathbf {z} _{1},\dots ,\mathbf {z} _{N}\}}$ は、データポイントごとに1つずつの潜在変数の集合であり、対応するデータポイントがどの混合成分に属するかを指定します。これには、前述のように、 の成分を持つ「K 個の中の 1 つ」のベクトル表現が使用されます。${\displaystyle z_{nk}}$ ${\displaystyle k=1\dots K}$
• ${\displaystyle \mathbf {\pi } }$ 混合成分の混合比率です。${\displaystyle K}$
• ${\displaystyle \mathbf {\mu } _{i=1\dots K}}$ 各混合成分に関連付けられたパラメータ（平均と精度）を指定します。${\displaystyle \mathbf {\Lambda } _{i=1\dots K}}$

すべての変数の結合確率は次のように書き直すことができる。

${\displaystyle p(\mathbf {X} ,\mathbf {Z} ,\mathbf {\pi } ,\mathbf {\mu } ,\mathbf {\Lambda } )=p(\mathbf {X} \mid \mathbf {Z} ,\mathbf {\mu } ,\mathbf {\Lambda } )p(\mathbf {Z} \mid \mathbf {\pi } )p(\mathbf {\pi } )p(\mathbf {\mu } \mid \mathbf {\Lambda } )p(\mathbf {\Lambda } )}$

個々の要因は

${\displaystyle {\begin{aligned}p(\mathbf {X} \mid \mathbf {Z} ,\mathbf {\mu } ,\mathbf {\Lambda } )&=\prod _{n=1}^{N}\prod _{k=1}^{K}{\mathcal {N}}(\mathbf {x} _{n}\mid \mathbf {\mu } _{k},\mathbf {\Lambda } _{k}^{-1})^{z_{nk}}\\p(\mathbf {Z} \mid \mathbf {\pi } )&=\prod _{n=1}^{N}\prod _{k=1}^{K}\pi _{k}^{z_{nk}}\\p(\mathbf {\pi } )&={\frac {\Gamma (K\alpha _{0})}{\Gamma (\alpha _{0})^{K}}}\prod _{k=1}^{K}\pi _{k}^{\alpha _{0}-1}\\p(\mathbf {\mu } \mid \mathbf {\Lambda } )&=\prod _{k=1}^{K}{\mathcal {N}}(\mathbf {\mu } _{k}\mid \mathbf {\mu } _{0},(\beta _{0}\mathbf {\Lambda } _{k})^{-1})\\p(\mathbf {\Lambda } )&=\prod _{k=1}^{K}{\mathcal {W}}(\mathbf {\Lambda } _{k}\mid \mathbf {W} _{0},\nu _{0})\end{aligned}}}$

どこ

${\displaystyle {\begin{aligned}{\mathcal {N}}(\mathbf {x} \mid \mathbf {\mu } ,\mathbf {\Sigma } )&={\frac {1}{(2\pi )^{D/2}}}{\frac {1}{|\mathbf {\Sigma } |^{1/2}}}\exp \left\{-{\frac {1}{2}}(\mathbf {x} -\mathbf {\mu } )^{\rm {T}}\mathbf {\Sigma } ^{-1}(\mathbf {x} -\mathbf {\mu } )\right\}\\{\mathcal {W}}(\mathbf {\Lambda } \mid \mathbf {W} ,\nu )&=B(\mathbf {W} ,\nu )|\mathbf {\Lambda } |^{(\nu -D-1)/2}\exp \left(-{\frac {1}{2}}\operatorname {Tr} (\mathbf {W} ^{-1}\mathbf {\Lambda } )\right)\\B(\mathbf {W} ,\nu )&=|\mathbf {W} |^{-\nu /2}\left\{2^{\nu D/2}\pi ^{D(D-1)/4}\prod _{i=1}^{D}\Gamma \left({\frac {\nu +1-i}{2}}\right)\right\}^{-1}\\D&={\text{dimensionality of each data point}}\end{aligned}}}$

と仮定します。 ${\displaystyle q(\mathbf {Z} ,\mathbf {\pi } ,\mathbf {\mu } ,\mathbf {\Lambda } )=q(\mathbf {Z} )q(\mathbf {\pi } ,\mathbf {\mu } ,\mathbf {\Lambda } )}$

その後^{[ 3 ]}

${\displaystyle {\begin{aligned}\ln q^{*}(\mathbf {Z} )&=\operatorname {E} _{\mathbf {\pi } ,\mathbf {\mu } ,\mathbf {\Lambda } }[\ln p(\mathbf {X} ,\mathbf {Z} ,\mathbf {\pi } ,\mathbf {\mu } ,\mathbf {\Lambda } )]+{\text{constant}}\\&=\operatorname {E} _{\mathbf {\pi } }[\ln p(\mathbf {Z} \mid \mathbf {\pi } )]+\operatorname {E} _{\mathbf {\mu } ,\mathbf {\Lambda } }[\ln p(\mathbf {X} \mid \mathbf {Z} ,\mathbf {\mu } ,\mathbf {\Lambda } )]+{\text{constant}}\\&=\sum _{n=1}^{N}\sum _{k=1}^{K}z_{nk}\ln \rho _{nk}+{\text{constant}}\end{aligned}}}$

ここで定義したのは

${\displaystyle \ln \rho _{nk}=\operatorname {E} [\ln \pi _{k}]+{\frac {1}{2}}\operatorname {E} [\ln |\mathbf {\Lambda } _{k}|]-{\frac {D}{2}}\ln(2\pi )-{\frac {1}{2}}\operatorname {E} _{\mathbf {\mu } _{k},\mathbf {\Lambda } _{k}}[(\mathbf {x} _{n}-\mathbf {\mu } _{k})^{\rm {T}}\mathbf {\Lambda } _{k}(\mathbf {x} _{n}-\mathbf {\mu } _{k})]}$

利回り の式の両辺を累乗する${\displaystyle \ln q^{*}(\mathbf {Z} )}$

${\displaystyle q^{*}(\mathbf {Z} )\propto \prod _{n=1}^{N}\prod _{k=1}^{K}\rho _{nk}^{z_{nk}}}$

これを正規化することを要求すると、のすべての値にわたって合計が1になることが要求され、 ${\displaystyle \rho _{nk}}$ ${\displaystyle k}$

${\displaystyle q^{*}(\mathbf {Z} )=\prod _{n=1}^{N}\prod _{k=1}^{K}r_{nk}^{z_{nk}}}$

どこ

${\displaystyle r_{nk}={\frac {\rho _{nk}}{\sum _{j=1}^{K}\rho _{nj}}}}$

言い換えると、は単一観測多項分布、および各個体 上の因子の積であり、のパラメータを持つ単一観測多項分布として分布します。 ${\displaystyle q^{*}(\mathbf {Z} )}$ ${\displaystyle \mathbf {z} _{n}}$ ${\displaystyle r_{nk}}$ ${\displaystyle k=1\dots K}$

さらに、我々は

${\displaystyle \operatorname {E} [z_{nk}]=r_{nk}\,}$

これはカテゴリ分布の標準的な結果です。

ここで、因子 について考えると、上で指定したガウス混合モデルを定義するグラフィカル モデルの構造により、 それが に自動的に因数分解されることに注意してください。${\displaystyle q(\mathbf {\pi } ,\mathbf {\mu } ,\mathbf {\Lambda } )}$ ${\displaystyle q(\mathbf {\pi } )\prod _{k=1}^{K}q(\mathbf {\mu } _{k},\mathbf {\Lambda } _{k})}$

それから、

${\displaystyle {\begin{aligned}\ln q^{*}(\mathbf {\pi } )&=\ln p(\mathbf {\pi } )+\operatorname {E} _{\mathbf {Z} }[\ln p(\mathbf {Z} \mid \mathbf {\pi } )]+{\text{constant}}\\&=(\alpha _{0}-1)\sum _{k=1}^{K}\ln \pi _{k}+\sum _{n=1}^{N}\sum _{k=1}^{K}r_{nk}\ln \pi _{k}+{\text{constant}}\end{aligned}}}$

両辺の指数をとると、ディリクレ分布として認識される。${\displaystyle q^{*}(\mathbf {\pi } )}$

${\displaystyle q^{*}(\mathbf {\pi } )\sim \operatorname {Dir} (\mathbf {\alpha } )\,}$

どこ

${\displaystyle \alpha _{k}=\alpha _{0}+N_{k}\,}$

どこ

${\displaystyle N_{k}=\sum _{n=1}^{N}r_{nk}\,}$

ついに

${\displaystyle \ln q^{*}(\mathbf {\mu } _{k},\mathbf {\Lambda } _{k})=\ln p(\mathbf {\mu } _{k},\mathbf {\Lambda } _{k})+\sum _{n=1}^{N}\operatorname {E} [z_{nk}]\ln {\mathcal {N}}(\mathbf {x} _{n}\mid \mathbf {\mu } _{k},\mathbf {\Lambda } _{k}^{-1})+{\text{constant}}}$

とを含む項をグループ化して読み取ると、結果は次式で表される ガウス・ウィシャート分布となる。${\displaystyle \mathbf {\mu } _{k}}$ ${\displaystyle \mathbf {\Lambda } _{k}}$

${\displaystyle q^{*}(\mathbf {\mu } _{k},\mathbf {\Lambda } _{k})={\mathcal {N}}(\mathbf {\mu } _{k}\mid \mathbf {m} _{k},(\beta _{k}\mathbf {\Lambda } _{k})^{-1}){\mathcal {W}}(\mathbf {\Lambda } _{k}\mid \mathbf {W} _{k},\nu _{k})}$

定義を踏まえると

${\displaystyle {\begin{aligned}\beta _{k}&=\beta _{0}+N_{k}\\\mathbf {m} _{k}&={\frac {1}{\beta _{k}}}(\beta _{0}\mathbf {\mu } _{0}+N_{k}{\bar {\mathbf {x} }}_{k})\\\mathbf {W} _{k}^{-1}&=\mathbf {W} _{0}^{-1}+N_{k}\mathbf {S} _{k}+{\frac {\beta _{0}N_{k}}{\beta _{0}+N_{k}}}({\bar {\mathbf {x} }}_{k}-\mathbf {\mu } _{0})({\bar {\mathbf {x} }}_{k}-\mathbf {\mu } _{0})^{\rm {T}}\\\nu _{k}&=\nu _{0}+N_{k}\\N_{k}&=\sum _{n=1}^{N}r_{nk}\\{\bar {\mathbf {x} }}_{k}&={\frac {1}{N_{k}}}\sum _{n=1}^{N}r_{nk}\mathbf {x} _{n}\\\mathbf {S} _{k}&={\frac {1}{N_{k}}}\sum _{n=1}^{N}r_{nk}(\mathbf {x} _{n}-{\bar {\mathbf {x} }}_{k})(\mathbf {x} _{n}-{\bar {\mathbf {x} }}_{k})^{\rm {T}}\end{aligned}}}$

最後に、これらの関数は の値を必要とし、 は を利用し、 は、、に基づいて定義されることに注意してください。これらの期待値が取られる分布が決定されたので、それらの公式を導くことができます。 ${\displaystyle r_{nk}}$ ${\displaystyle \rho _{nk}}$ ${\displaystyle \operatorname {E} [\ln \pi _{k}]}$ ${\displaystyle \operatorname {E} [\ln |\mathbf {\Lambda } _{k}|]}$ ${\displaystyle \operatorname {E} _{\mathbf {\mu } _{k},\mathbf {\Lambda } _{k}}[(\mathbf {x} _{n}-\mathbf {\mu } _{k})^{\rm {T}}\mathbf {\Lambda } _{k}(\mathbf {x} _{n}-\mathbf {\mu } _{k})]}$

${\displaystyle {\begin{aligned}\operatorname {E} _{\mathbf {\mu } _{k},\mathbf {\Lambda } _{k}}[(\mathbf {x} _{n}-\mathbf {\mu } _{k})^{\rm {T}}\mathbf {\Lambda } _{k}(\mathbf {x} _{n}-\mathbf {\mu } _{k})]&=D\beta _{k}^{-1}+\nu _{k}(\mathbf {x} _{n}-\mathbf {m} _{k})^{\rm {T}}\mathbf {W} _{k}(\mathbf {x} _{n}-\mathbf {m} _{k})\\\ln {\widetilde {\Lambda }}_{k}&\equiv \operatorname {E} [\ln |\mathbf {\Lambda } _{k}|]=\sum _{i=1}^{D}\psi \left({\frac {\nu _{k}+1-i}{2}}\right)+D\ln 2+\ln |\mathbf {W} _{k}|\\\ln {\widetilde {\pi }}_{k}&\equiv \operatorname {E} \left[\ln |\pi _{k}|\right]=\psi (\alpha _{k})-\psi \left(\sum _{i=1}^{K}\alpha _{i}\right)\end{aligned}}}$

これらの結果は、

${\displaystyle r_{nk}\propto {\widetilde {\pi }}_{k}{\widetilde {\Lambda }}_{k}^{1/2}\exp \left\{-{\frac {D}{2\beta _{k}}}-{\frac {\nu _{k}}{2}}(\mathbf {x} _{n}-\mathbf {m} _{k})^{\rm {T}}\mathbf {W} _{k}(\mathbf {x} _{n}-\mathbf {m} _{k})\right\}}$

対応する値の合計が 1 になるように 正規化することで、これらを比例値から絶対値に変換できます。${\displaystyle k}$

ご了承ください：

1. 変数およびのパラメータ、、およびの更新方程式は、統計量、、およびに依存し、これらの統計量は、に依存します。${\displaystyle \beta _{k}}$ ${\displaystyle \mathbf {m} _{k}}$ ${\displaystyle \mathbf {W} _{k}}$ ${\displaystyle \nu _{k}}$ ${\displaystyle \mathbf {\mu } _{k}}$ ${\displaystyle \mathbf {\Lambda } _{k}}$ ${\displaystyle N_{k}}$ ${\displaystyle {\bar {\mathbf {x} }}_{k}}$ ${\displaystyle \mathbf {S} _{k}}$ ${\displaystyle r_{nk}}$
2. 変数のパラメータの更新方程式は統計量 に依存し、統計量は​​ に依存します。${\displaystyle \alpha _{1\dots K}}$ ${\displaystyle \mathbf {\pi } }$ ${\displaystyle N_{k}}$ ${\displaystyle r_{nk}}$
3. の更新方程式は、、、に直接循環依存関係があり、また、およびを介して、、、およびに間接的に循環依存関係があります。${\displaystyle r_{nk}}$ ${\displaystyle \beta _{k}}$ ${\displaystyle \mathbf {m} _{k}}$ ${\displaystyle \mathbf {W} _{k}}$ ${\displaystyle \nu _{k}}$ ${\displaystyle \mathbf {W} _{k}}$ ${\displaystyle \nu _{k}}$ ${\displaystyle \alpha _{1\dots K}}$ ${\displaystyle {\widetilde {\pi }}_{k}}$ ${\displaystyle {\widetilde {\Lambda }}_{k}}$

これは、次の 2 つのステップを交互に実行する反復的な手順を示唆しています。

1. 他のすべてのパラメータの現在の値を使用して値を計算する E ステップ。${\displaystyle r_{nk}}$
2. の新しい値を使用して、他のすべてのパラメータの新しい値を計算するM ステップ。${\displaystyle r_{nk}}$

これらのステップは、ガウス混合モデルのパラメータに対する最大尤度解または最大事後確率（MAP）解を導出する標準的なEMアルゴリズムと密接に対応していることに注意してください。Eステップの役割は、データが与えられた場合の潜在変数の事後確率、すなわち統計量、、の計算と密接に対応しています。統計量、、の計算は、データ全体にわたる対応する「ソフトカウント」統計量の計算と密接に対応しています。また、これらの統計量を用いてパラメータの新しい値を計算することは、ガウス混合モデルにおける通常のEMでソフトカウントを用いて新しいパラメータ値を計算することに密接に対応しています。 ${\displaystyle r_{nk}}$ ${\displaystyle p(\mathbf {Z} \mid \mathbf {X} )}$ ${\displaystyle N_{k}}$ ${\displaystyle {\bar {\mathbf {x} }}_{k}}$ ${\displaystyle \mathbf {S} _{k}}$

指数分布族

[編集]

前の例では、観測されない変数上の分布が「パラメータ」上の分布と「潜在データ」上の分布に分解されると仮定すると、各変数について導出された「最良」分布は、その変数の対応する事前分布と同じ族に属することに注目してください。これは、指数族から導出されるすべての事前分布に当てはまる一般的な結果です。

参照

[編集]

• 変分メッセージパッシング: 変分ベイズ推論のためのモジュラーアルゴリズム。
• 変分オートエンコーダ: 確率的グラフィカルモデルと変分ベイズ法のファミリーに属する人工ニューラル ネットワーク。
• 期待最大化アルゴリズム: 変分ベイズ推論の特殊なケースに対応する関連アプローチ。
• 一般化フィルタリング: 非線形状態空間モデル用の変分フィルタリング方式。
• 変分法: 関数の最大化または最小化を扱う数学的分析の分野。
• 最大エントロピー識別：これは、追加の大きなマージン制約を導入し、考慮することを可能にする変分推論の枠組みである^{[ 7 ]}

参考文献

[編集]

1. Tran, Viet Hung (2018). 「情報幾何学によるコピュラ変分ベイズ推論」. arXiv : 1803.10998 [ cs.IT ].
2. Adamčík, Martin (2014). 「ブレグマンダイバージェンスの情報幾何学とマルチエキスパート推論への応用」 .エントロピー. 16 (12): 6338– 6381. Bibcode : 2014Entrp..16.6338A . doi : 10.3390/e16126338 .
3. Nguyen, Duy (2023年8月15日). 「変分ベイズ法の詳細な入門ノート」 . doi : 10.2139/ssrn.4541076 . SSRN 4541076. 2023年8月15日閲覧。 
4. Lee, Se Yoon (2021). 「ギブスサンプラーと座標上昇変分推論：集合論的レビュー」. Communications in Statistics - Theory and Methods . 51 (6): 1– 21. arXiv : 2008.01006 . doi : 10.1080/03610926.2021.1921214 . S2CID 220935477 . 
5. Boyd, Stephen P.; Vandenberghe, Lieven (2004).凸最適化(PDF) . Cambridge University Press. ISBN  978-0-521-83378-3. 2011年10月15日閲覧。
6. ビショップ、クリストファー・M. (2006). 「第10章」.パターン認識と機械学習. シュプリンガー. ISBN 978-0-387-31073-2。
7. Sotirios P. Chatzis, “ Infinite Markov-Switching Maximum Entropy Discrimination Machines ,” Proc. 30th International Conference on Machine Learning (ICML). Journal of Machine Learning Research: Workshop and Conference Proceedings, vol. 28, no. 3, pp. 729–737, 2013年6月.

外部リンク

[編集]

• オンライン教科書「Information Theory, Inference, and Learning Algorithms Archived 2017-05-12 at the Wayback Machine」、David JC MacKay著では、変分法の入門書が紹介されています (p. 422)。
• 変分ベイズに関するチュートリアル。Fox, C. および Roberts, S. 2012. 人工知能レビュー、doi : 10.1007/s10462-011-9236-8。
• 変分ベイズ リポジトリ2003 年までの近似ベイズ学習のための変分法の使用に関連する研究論文、ソフトウェア、およびリンクのリポジトリです。
• MJ Beal 著の「近似ベイズ推論の変分アルゴリズム」には、EM と変分ベイズ EM の比較と、変分ベイズ HMM を含むいくつかのモデルの導出が含まれています。
• より数学的に詳細な説明を読む前に、Jason Eisner 著の「High-Level Explaination of Variational Inference」を読む価値があるかもしれません。
• 情報幾何学によるコピュラ変分ベイズ推論（PDF）（Tran, VH 2018）。この論文は主に学生向けに書かれています。ブレグマン情報量を用いて、変分ベイズ推論は真のモデルを任意に相関する（コピュラ）分布空間に一般化したピタゴラス射影に過ぎず、独立空間はその特別なケースに過ぎないことを示しています。
• 変分ベイズ法の詳細な入門書 ノート. Nguyen, D. 2023