next up previous
Next: Bibliography

混合因子分析モデルのベイズ解析
宇津木明男(P)    熊谷徹
生命工学工業技術研究所

Bayesian Analysis of Mixtures of Factor Analyzers
Akio Utsugi and Toru Kumagai
National Institute of Bioscience and Human-Technology
email: utsugi@nibh.go.jp, http://www.aist.go.jp/NIBH/~b0616/

Abstract We introduce conjugate priors on the parameters of a mixture of factor analyzers and then construct a Gibbs sampling algorithm. We also present a deterministic algorithm for the estimation of the parameters using modes instead of samples from conditional posteriors in the Gibbs sampling. These algorithms are compared on simulation.


1. はじめに

教師無し学習には、潜在クラスの獲得を目的とするものと、潜在部分空間の 獲得を目的とするものの2種類がある。 確率モデルとしては、前者には混合分布モデルが、後者には 因子分析モデルが使われる事が多い。 混合因子分析モデル(MFA: Mixture of Factor Anlayzers)は、 両者の率直な融合モデルで、パタン認識や画像圧縮等に利用されている [1][2]。

EMアルゴリズムによるパラメータの最尤推定が行なわれる事が多いが、 尤度が有界ではないので、最尤推定値が得られる保証は無い。 また、尤度の発散が防げたとしても、データの次元にほぼ比例して パラメータ数が増加するので、安定した推定値を得るのが難しい場合がある。 このような場合によく使われる対策は、パラメータの正則化である。 確率モデルでは、パラメータに事前確率を導入する事が、正則化に当たる。 以下では、共役分布による事前確率を導入し、 パラメータの事後分布に従うサンプル系列を生成するギブスサンプリングの アルゴリズムを導出する。また、決定論的な近似アルゴリズムも導出し、 ギブスサンプリングの結果と比較する。


2. MFAのデータ生成モデル

n個のデータ点 $\mbox{\boldmath$x$ }_i, (i=1, \ldots, n)$から成るデータを行列 $\mbox{\boldmath$X$ }\;(p \times n) = [\mbox{\boldmath$x$ }_1, \ldots, \mbox{\boldmath$x$ }_n]$で表す。 MFAでは、各データ点がm個の因子分析(FA)ユニットの内の一つから生成されると仮 定する。 データ点とFAユニットの対応を表す2値変数 $\mbox{\boldmath$Z$ }\; (m \times n)=[z_{ki}]$と、 各FAユニットの因子スコア $Y=\{\mbox{\boldmath$Y$ }_k\; (q_k \times n)=[\mbox{\boldmath$y$ }_{k1}, \ldots,
\mbox{\boldmath$y$ }_{kn}]; k=1, \ldots, m\}$を使って、データの生成に関する以下のような確率モデ ルを考える。

   \begin{displaymath}\mbox{\boldmath$X$ }\vert Y, \mbox{\boldmath$Z$ }, \Lambda, \...
...box{\boldmath$\mu$ }_k, \mbox{\boldmath$\Psi$ }_k )\}^{z_{ki}}
\end{displaymath}


   \begin{displaymath}Y\vert\mbox{\boldmath$Z$ }
\sim
\prod_{i=1}^n \prod_{k=1}^m
\...
...math$y$ }_{ki}\vert{\bf {0}}, \mbox{\boldmath$I$ })\}^{z_{ki}}
\end{displaymath}


   \begin{displaymath}\mbox{\boldmath$Z$ }\vert\tau
\sim
\prod_{i=1}^n \prod_{k=1}^m \tau_k^{z_{ki}}
\end{displaymath}

ここで、 $N(\cdot\vert\mbox{\boldmath$\mu$ }, \mbox{\boldmath$\Sigma$ })$は平均 $\mbox{\boldmath$\mu$ }$, 共 分散 $\mbox{\boldmath$\Sigma$ }$の正規分布の密度関数を表す。 Y $\mbox{\boldmath$Z$ }$は観測されない潜在変数と見なされる。 MFAモデルのパラメータ$\Theta$は、 FAユニットの選択確率 $\tau=\{\tau_k; k=1, \ldots, m\}$, 中心位置 $\mu=\{\mbox{\boldmath$\mu$ }_k \; (p \times 1); k=1, \ldots, m\}$, 因子負荷 $\Lambda=\{\mbox{\boldmath$\Lambda$ }_k \; (p \times q_k)=[\mbox{\boldmath$\lambda$ }_{k1}, \ldots, \mbox{\boldmath$\lambda$ }_{kq_k}]; k=1, \ldots,
m\}$ および独自性 $\Psi=\{\mbox{\boldmath$\Psi$ }_k \; (p \times p)=\mathop{\rm diag}(\psi_{k1}, \ldots,
\psi_{kp}); k=1, \ldots, m\}$の4種類から成る。 パラメータ$\Theta$の最尤推定のためのEMアルゴリズムは、 Ghahramaniら[2]によって求められている。


3. パラメータの事前確率分布の導入

共役分布を利用して、次のようなパラメータの事前確率を導入する。

   \begin{displaymath}\tau\vert\gamma \sim D(\mbox{\boldmath$\tau$ }\vert \gamma)
\end{displaymath}

$D(\cdot\vert \gamma)$はディレクレ分布の密度関数である。

   \begin{displaymath}\mu\vert\Psi^{-1}, \bar{\mbox{\boldmath$\mu$ }}, \alpha_1
\si...
...box{\boldmath$\mu$ }}, \alpha_1^{-1}\mbox{\boldmath$\Psi$ }_k)
\end{displaymath}


   \begin{displaymath}\Lambda\vert\Psi^{-1}, \alpha_2
\sim
\prod_{k=1}^m \prod_{j=...
... }_{kj}\vert{\bf {0}}, \alpha_2^{-1}\mbox{\boldmath$\Psi$ }_k)
\end{displaymath}

ただし、 $\bar{\mbox{\boldmath$\mu$ }}$は、簡単のためデータの平均 $\sum_i^n \mbox{\boldmath$x$ }_i/n$によって推定す る事にする。ここでは、データを中心化することで、 $\bar{\mbox{\boldmath$\mu$ }}=0$とする。

\begin{displaymath}\Psi^{-1}\vert\delta, \beta
\sim
\prod_{k=1}^m \prod_{j=1}^p
G(\psi_{kj}^{-1}\vert\delta, \beta)
\end{displaymath}

$G(\cdot\vert\delta, \beta)$は平均 $\delta/\beta$, 分散 $\delta/\beta^2$となるようにパラメタライズされたガンマ分布の密度関数である。 $H=\{\alpha_1, \alpha_2, \beta, \delta, \gamma\}$は、ハイパーパラメータ と呼ばれる。しかし、以下では、 $Y, \mbox{\boldmath$Z$ }, \Theta, H$をまとめて パラメータと呼ぶ事にする。


4. ギブスサンプリング

ベイズ推論では、パラメータの事後分布 $f(Y, \mbox{\boldmath$Z$ }, \Theta, H\vert X)$に基づいて統計的推論を行なう。 ベイズの定理を使って、形式的に事後分布の関数形を求める事はできるが、 それ以上の計算は困難である。このような場合、事後分布に従うパラメータ のサンプル系列を生成し、これに基づいて近似的な推論計算を行なう事 が良くやられる。しかし、事後分布から直接サンプル系列を生成する 事も難しい。そこで、ギブスサンプリングを使う事にする。

まず、パラメータをブロックに分け、各ブロック内のパラメー タの事後分布を他のパラメータで条件付けて求める。 この条件付き事後分布からのサンプル生成を交互に繰り返す事 によって得られるパラメータ系列の分布は、 ある条件の下でパラメータの事後分布に近付いていくことが示されている。

潜在変数 $Y, \mbox{\boldmath$Z$ }$の条件つき事後分布は、 Ghahramaniら[2]によって求められているので、ここでは、$\Theta$Hの条件つき事後分布を求めよう。データの生成確率とパラメータの 事前確率から、$\Theta$の条件つき事後分布が以下のように 求められる。

\begin{displaymath}\mbox{\boldmath$\tau$ }\vert\mbox{\boldmath$X$ }, \mbox{\bold...
...D(\mbox{\boldmath$\tau$ }\vert n_1+\gamma, \ldots, n_k+\gamma)
\end{displaymath}


$\displaystyle {\vec(\tilde{\mbox{\boldmath$\Lambda$ }}_k) \vert\mbox{\boldmath$X$ }, Y, \mbox{\boldmath$Z$ }, \mbox{\boldmath$\Psi$ }_k^{-1}, \alpha_1, \alpha_2}$
  $\textstyle \sim$ $\displaystyle N(\vec(\tilde{\mbox{\boldmath$\Lambda$ }}_k) \vert\vec[\mbox{\boldmath$C$ }_{XY\;k}(\mbox{\boldmath$C$ }_{YY\;k}+\mbox{\boldmath$A$ }_k)^{-1}],$  
    $\displaystyle (\mbox{\boldmath$C$ }_{YY\;k} + \mbox{\boldmath$A$ }_k)^{-1} \otimes \mbox{\boldmath$\Psi$ }_k )$  


$\displaystyle {\mbox{\boldmath$\Psi$ }_k^{-1}\vert\mbox{\boldmath$X$ }, Y, \mbo...
...$Z$ }, \tilde{\mbox{\boldmath$\Lambda$ }}_k,
\alpha_1, \alpha_2, \beta, \delta}$
  $\textstyle \sim$ $\displaystyle G(\mbox{\boldmath$\Psi$ }_k^{-1} \vert\; \frac{1}{2}(n_k+q_k+2\delta+1), \;
\frac{1}{2}
\mathop{\rm diag}[\mbox{\boldmath$C$ }_{XX\;k}$  
    $\displaystyle -2 \mbox{\boldmath$C$ }_{XY\;k} \tilde{\mbox{\boldmath$\Lambda$ }...
... }_k)\tilde{\mbox{\boldmath$\Lambda$ }}_k^{\prime}+2\beta\mbox{\boldmath$I$ }])$  

ここで、 $\tilde{\mbox{\boldmath$\Lambda$ }}_k = [\mbox{\boldmath$\mu$ }_k, \mbox{\boldmath$\Lambda$ }_k]$, $n_k = \sum_{i=1}^n z_{ki}$, $\mbox{\boldmath$C$ }_{XX\;k}= \sum_{i=1}^n
z_{ki}\mbox{\boldmath$x$ }_i\mbox{\b...
...ki} \tilde{\mbox{\boldmath$y$ }}_{ki}\tilde{\mbox{\boldmath$y$ }}_{ki}^{\prime}$, $\tilde{\mbox{\boldmath$y$ }}_{ki}=[1, \mbox{\boldmath$y$ }_{ki}^{\prime}]'$, $\mbox{\boldmath$A$ }_k=\mathop{\rm diag}(\alpha_1, \alpha_2{\bf 1}_{q_k}^{\prime})$, である。

また、ハイパーパラメータ $\alpha_1, \alpha_2, \beta$にプライア

\begin{displaymath}h \sim G(h\vert d_h, s_h), \; (h=\alpha_1, \alpha_2, \beta)
\end{displaymath}

を導入する事で、 $\alpha_1, \alpha_2, \beta$の条件つき事後分布が

   \begin{displaymath}\alpha_1 \vert \mbox{\boldmath$X$ }, \Theta, d_{\alpha 1}, s_...
...ldmath$\Psi$ }_k^{-1} \mbox{\boldmath$\mu$ }_k +s_{\alpha 1} )
\end{displaymath}


   \begin{displaymath}\alpha_2\vert\mbox{\boldmath$X$ }, \Theta, d_{\alpha 2}, s_{\...
...th$\Psi$ }_k^{-1} \mbox{\boldmath$\Lambda$ }_k) +s_{\alpha 2})
\end{displaymath}


   \begin{displaymath}\beta\vert\mbox{\boldmath$X$ }, \Theta, \delta, d_\beta, s_\b...
...k=1}^m\mathop{\rm tr}(\mbox{\boldmath$\Psi$ }_k^{-1})+s_\beta)
\end{displaymath}

として求まる。ここで、 $q = \sum_k^m q_k$である。 残りのハイパーパラメータは、良く使われる値に固定する。 これで、ギブスサンプリングが完成した。


5. 決定論的なベイズ推定アルゴリズム

条件付き事後確率からランダムサンプルを得る代わり に、条件付き事後モードを使う事によって決定論的なベイズ推定アルゴリズムが得られる。 特に、ハイパーパラメータHを固定した場合は、これは$\Theta$の最大事後確率(MAP)推定の EMアルゴリズムと一致する。Hを固定しない場合は、 ハイパーパラメータ探索付きのMAP推定EMアルゴリズムと見なす事ができる。 決定論的アルゴリズムは、定常状態への到達が容易に確認できるため、 ギブスサンプリングに比べて早い段階で学習を打ち切る事ができる。


6. シミュレーション

図1は、MFAモデルに従って人工的に生成したデータに対する、決定論的アルゴリ ズムによるパラメータ推定の結果である。 FAユニットの中心と選択確率が正方形の位置と大きさによって、 因子負荷と独自性が線分によって、それぞれ表されている。 ギブスサンプリングの生成するパラメータのサンプル系列は、定常状態に達した後は、 このパラメータ配置の周辺を揺れ動く。 図2はハイパーパラメータ$\alpha_1$$\alpha_2$の時間変化のグラフである。決定論的アルゴリズムは ギブスサンプリングの中心的傾向を良く模倣している事が分かる。


図1:決定論的アルゴリズムによる推定値

n=500, p=2, q_k=1, m=5


図2:ハイパーパラメータの時間変化



 
next up previous
Next: Bibliography

1999-09-30