変分推論 〜ガウス混合モデルへの適用〜

はじめに

今回は、ベイズ推論の1手法である変分推論について解説する。具体例としてガウス混合モデルを用いたクラスタリングを行い、実際の計算結果を示す。

変分推論

観測値をX=\{x_1,\dots,x_N\}、潜在変数とパラメータを合わせてZ=\{z_1,\dots,z_N\}とおく。ただし、x_n,z_n\in \mathbb{R}^Dである。このとき、同時確率分布 p(X,Z)を考え、ベイズの定理を適用すると

(1)   \begin{equation*} p(Z|X)=\frac{p(X|Z)p(Z)}{p(X)} \end{equation*}

を得る。ベイズ推論の目的は事後確率p(Z|X)を求めることである。一般に、p(Z|X)を厳密に求めることはできないので、近似関数q(Z)を考える。次式\mathcal{L}（ELBO:Evidence Lower Bound）が最大になるようにq(Z)を最適化する手法を変分推論と呼ぶ。

(2)   \begin{equation*} \mathcal{L}=\int dZ\;q(Z)\ln{\left(\frac{p(X,Z)}{q(Z)}\right)} \end{equation*}

実は、ELBOの最大化とKLダイバージェンス

(3)   \begin{equation*} {\rm KL}(q||p)=-\int dZ\;q(Z)\ln{\left(\frac{p(Z|X)}{q(Z)}\right)} \end{equation*}

を最小化することは等価である。従って、KLダイバージェンスの最小化を用いて変分推論を定式化することもできるが、今回はELBOの最大化を利用して定式化する。

平均場近似

上で個々のデータの集まりとしてZ(=\{z_1,\cdots,z_N\})を定義したが、ここでは複数のパラメータ（と潜在変数）の集まりとしてZ(=\{Z_1,\cdots,Z_M\})を考える。q(Z)=q(Z_1,\cdots,Z_M)を近似するとき最も良く使われる手法が次の平均場近似である。

(4)   \begin{equation*} q(Z)=\prod_{m=1}^{M}q(Z_m) \end{equation*}

この近似は、Z_m間の相関を無視し、それぞれが独立に生成されると仮定したものである。最適な各関数 q(Z_i)を求めるには、\mathcal{L}をq(Z_i)で汎関数微分すれば良い。このことを念頭において式(4)を式(2)に代入すると

(5)   \begin{equation*} \mathcal{L}=\int dZ_i\;q_i\left\{\int \left(\prod_{m\neq i}^{M}dZ_m\;q_m\right)\ln{p(X,Z)}\right\}-\int dZ_i\;q_i\ln{q_i}+C \end{equation*}

を得る。煩雑さを避けるためq_i=q(Z_i)とした。Cはq_iに依存しない全ての項をまとめたものである。いま次の期待値

(6)   \begin{equation*} \mathbb{E}_{m\neq i}\left[\ln{p(X,Z)}\right]=\int \left(\prod_{m\neq i}^{M}dZ_m\;q_m\right)\ln{p(X,Z)} \end{equation*}

を定義すると

(7)   \begin{equation*} \mathcal{L}=\int dZ_i\;q_i\bigl[ \mathbb{E}_{m\neq i}\left[\ln{p(X,Z)}\right]-\ln{q_i} \bigr]+C \end{equation*}

と書くことができる。上式をq_iで汎関数微分したものを0と置くことにより次式を得る。

(8)   \begin{equation*} \ln{q_i}=\mathbb{E}_{m\neq i}\left[\ln{p(X,Z)}\right] \end{equation*}

これが近似関数q_iを求める式である。i=1,\cdots,Mであるから次の反復アルゴリズムで計算することになる。

  1. \{q_2,q_3,\cdots,q_M\}を初期化する。
  2. 式(8)からq_1を求める。
  3. 今求めたq_1を含む\{q_1,q_3,\cdots,q_M\}使い、式(8)からq_2を求める。
  4. この手順を繰り返してq_Mまで計算する。
  5. 手順2へ。
  6. 各q_iの変動がなくなるまで繰り返す。

ガウス混合モデル

混合モデルとは複数の確率分布（クラスタ）を用いてひとつの確率分布を表現するモデルである。具体的に説明するため、下図を考える。

この図は以下のことを表現している。

  • 混合比率\piが事前分布p(\pi)から生成される。ここで、\pi=(\pi_1,\cdots,\pi_K)である。Kはクラスタの数を表す。\pi_{k}は0以上1以下の実数値であり、\sum_{k=1}^{K}\pi_k=1を満たす。
  • 各クラスタ(k=1,\cdots,K)に対するパラメータ\theta_kが事前分布p(\theta_k)から生成される。\Theta=(\theta_1,\dots,\theta_K)
  • x_nがどのクラスタに割り振られるかを決める変数s_nが\piにより選ばれる。ここで、 s_n=(s_{n1},\cdots,s_{nK})である。s_{nk}は0か1の値をとる整数値であり、\sum_{k=1}^{K}s_{nk}=1を満たす。S=(s_1,\cdots,s_N)
  • s_nにより選択されたk番目の確率分布p(x_n|\theta_k)からx_nが生成される。

    • 上図を式で表すと次式となる。

    (9)   \begin{eqnarray*} p(X,S,\Theta,\pi) &=&p(X|S,\Theta)p(S|\pi)p(\Theta)p(\pi) \nonumber\\ &=&\Bigl\{\prod_{n=1}^{N}p(x_n|s_n,\Theta)p(s_n|\pi)\Bigr\} \Bigl\{\prod_{k=1}^{K}p(\theta_k)\Bigr\}p(\pi) \end{eqnarray*}

    例えば、クラスタ（p(x_n|\theta_k)）をガウス分布で表す場合、複数のガウス分布を適当な比率（\pi_{k}）で混ぜ合わせて目的の分布を近似することになる（下図）。このモデルをガウス混合モデル（Gaussian Mixture Model：GMM）と呼ぶ。

    ガウス混合モデルの変分推論

    ガウス混合モデルに平均場近似を適用する。出発式となるのは式(8)で示した次式である。

    (10)   \begin{equation*} \ln{q(Z_i)}=\mathbb{E}_{m\neq i}\left[\ln{p(X,Z)}\right] \end{equation*}

    Zに相当する変数は(S,\Theta,\pi)である。ガウス分布を考えるので \Theta=(\mu,\Lambda)と置ける。ここで\mu=(\mu_1,\cdots,\mu_K)は平均、\Lambda=(\Lambda_1,\cdots,\Lambda_K)は精度行列である。いま観測値x_nの次元をDとしているので\mu_kの次元もDである。また、行列\Lambda_kはD\times Dのサイズを持つ。

    式(10)より次の式を得る。

    (11)   \begin{eqnarray*} \ln{q(S)}&=&\mathbb{E}_{q(\mu,\Lambda,\pi)}\left[\ln{p(X,S,\mu,\Lambda,\pi)}\right]\\ \ln{q(\mu,\Lambda,\pi)}&=&\mathbb{E}_{q(S)}\left[\ln{p(X,S,\mu,\Lambda,\pi)}\right] \end{eqnarray*}

    この2式を辛抱強く変形していくと、事後確率 q(s_n),q(\pi),q(\mu_k),q(\Lambda_k)についての式を得る。詳細はPRMLや須山本を見てほしい。これら4式はお互いに依存しているので既に述べたように反復して解くことになる。

    計算結果

    観測値の次元Dを2としてPythonによる実装を行った。今回はベイズ推論関係のライブラリを使わず、フルスクラッチで実装した。参照した文献は、須山本のp144から始まる「4.4 ガウス混合モデルにおける推論」である。私のソースコードはここにある。最初に適当なガウス分布を用いて観測値を作成する（下図参照）。

    これを入力値として、K=3として変分推論を行う。結果は以下の通り。

    各クラスタを3色で示した。綺麗に分割できている。また、観測値以外の点についても、どのクラスタに属するかを確率として予測することができる。つまり、1点ごとに確率が求まるので、全平面を塗り分けることができる。今回は1点つき100回確率を求めその平均値を用いて色を決定した。図に見るように境界部分では色が混ざり連続的にもう一方の色へ変化していくことが分かる。このような変化を捉えることができるのがベイス推論の大きな特徴である。

    別の観測値を考える。

    推論結果は以下の通り。

    赤の境界近傍にある観測点の1つが紫色となっている。これは青である確率と赤で確率がほぼ互角となったためであろう。この図の注目すべき点は、右上の領域が青で塗られている点である。青の観測値は斜め方向に分散が大きいため、このように予測されたのである。単純に距離だけで判断するのであれば赤や緑で塗られるはずである。

    以下、数値計算する上での注意点を列記する。

  • 入力データは、平均0、標準偏差1となるように正規化した。この前処理により安定して収束するようになる。
  • 反復計算を行う過程で、式内に現れる半正定値行列の行列式が負になることがある。式変形は正しい（たぶん）ので純粋に数値計算時の誤差が引き起こす現象と思われる。負になった場合は、その時点で再度その行列を初期化した。
  • クラスタ中心の初期位置は、あらかじめK-Meansクラスタリングにより求めておくと安定して収束するようになる。乱数による初期化を行うと収束する場合もあるし、しない場合もあり、安定しない。

    • まとめ

    今回は、ベイス推論の1手法である変分推論を紹介し、クラスタリングに適用した事例を示した。ベイズ推論自体は、回帰や分類問題など幅広く応用できる機械学習の手法である。深層学習に触発されたAIブームの中、ベイズ推論もまた注目されている。

    参考文献

  • PRML
  • 須山本

    • どちらも私のバイブルである。PRMLについては演習問題もおすすめです。

    Kumada Seiya

    Kumada Seiya

    仕事であろうとなかろうと勉強し続ける、その結果”中身”を知ったエンジニアになれる

    最近の記事

    • 関連記事
    • おすすめ記事
    • 特集記事

    アーカイブ

    カテゴリー

    PAGE TOP