GMM及EM算法

標簽（空格分隔）： 機器學習

---

前言：

• EM(Exception Maximizition) -- 期望最大化算法，用於含有隱變量的概率模型參數的極大似然估計；
• GMM(Gaussian Mixture Model) -- 高斯混合模型，是一種多個高斯分布混合在一起的模型，主要應用EM算法估計其參數；
• 本篇博客首先從簡單的k-means算法給出EM算法的迭代形式，然后用GMM的求解過程給出EM算法的宏觀認識；最后給出EM的標准形式，並分析EM算法為什么收斂。

K-Means Clustering

目標函數(損失函數)

EM過程

1. 隨機初始化樣本中心點，即均值\(\mu_{k}\).
2. E step: 固定均值\(\mu_{k}\)，用\(r_{n}k\)最小化目標函數J。
3. M step: 固定\(r_{n}k\)，用\(\mu_{k}\)最小化目標函數J。
4. 不斷迭代，直至收斂。

Gaussian Mixture Model

將高斯混合模型看成是高斯分量的簡單線性疊加，目標是提供一類比單獨高斯分布更強大的概率模型。

• 假設有K個高斯分布構成了數據X，則其概率分布為：
  \(p(x) = \sum_{k=1}^{K}\pi_{k}N(x|\mu_{k},\Sigma_{k})\).
  其中，數據點x屬於第k個高斯分布的概率為\(\pi_{k}\)，則有$0\leq \pi_{k} \leq 1, \sum_{k=1}^{K}\pi_{k} = 1 $
• x的邊緣概率分布為\(p(x) = \sum_{k=1}^{K}\pi_{k}N(x|\mu_{k},\Sigma_{k})\)

對數似然函數為：

\(\ln P(X|\pi,\mu,\Sigma) = \sum_{n=1}^{N}\ln\{\sum_{k=1}^{K}\pi_{k}N(x_{n}|\mu_{k},\Sigma_{k})\}\)
其中，\(N(x_{n}|\mu_{k},\Sigma_{k}) =\frac{1}{(2\pi)^{\frac{D}{2}}|\Sigma|^{\frac{1}{2}}}exp\{-\frac{1}{2}(x_{n}-\mu_{k})^{T}\Sigma_{k}^{-1}(x_{n}-\mu_{k})\}\)

關於\(\mu_{k}\)使對數似然函數最大：

\(\frac{\partial P(X|\pi_{k},\mu_{k},\Sigma_{k})}{\partial\mu_{k}} = \sum_{n=1}^{N}\frac{\pi_{k}N(x_{n}|\mu_{k},\Sigma_{k})(-\frac{1}{2}*2*\Sigma_{k}^{-1}(x_{n}-\mu_{k}))}{\sum_{j}\pi_{j}N(x_{n}|\mu_{j},\Sigma_{j})} \\ = -\sum_{n=1}^{N}\frac{\pi_{k}N(x_{n}|\mu_{k},\Sigma_{k})}{\sum_{j}\pi_{j}N(x_{n}|\mu_{j},\Sigma_{j})}\Sigma_{k}^{-1}(x_{n}-\mu_{k}) = 0\)
兩邊同時乘以\(\Sigma_{k}\),
記\(\gamma(z_{n}k) = \frac{\pi_{k}N(x_{n}|\mu_{k},\Sigma_{k})}{\sum_{j}\pi_{j}N(x_{n}|\mu_{j},\Sigma_{j})}\)得：

\(\mu_{k} = \frac{1}{\sum_{n=1}^{N}\gamma(z_{n}k)}\sum_{n=1}^{N}\gamma(z_{n}k)·x_{n}\)

關於\(\Sigma_{k}\)使對數似然最大：

\(\Sigma_{k} = \frac{1}{\sum_{n=1}^{N}\gamma(z_{n}k)}\sum_{n=1}^{N}\gamma(z_{n}k)(x_{n}-\mu_{k})(x_{n}-\mu_{k})^{T}\)

關於混合系數\(\pi_{k}\)使對數似然函數最大：

由於所有混合系數的和為1，即滿足\(\sum_{k=1}^{K}\pi_{k} = 1\),則原來的最大似然問題變為帶約束條件的似然函數最大化，使用Lagrange乘子法：
\(max \ln P(X|\pi,\mu,\Sigma) + \lambda(\sum_{n=1}^{N}\pi_{k} - 1)\)
對\(\pi_{k}\)求偏導，令導數為0：
$0 = \sum_{n=1}^{N}\frac{N(x_{n}|\mu_{k},\Sigma_{k})}{\sum_{j}\pi_{j}N(x_{n}|\mu_{j},\Sigma_{j})} + \lambda $

兩邊同乘\(\pi_{k}\),得
\(0 = \sum_{n=1}^{N}\frac{\pi_{k}N(x_{n}|\mu_{k},\Sigma_{k})}{\sum_{j}\pi_{j}N(x_{n}|\mu_{j},\Sigma_{j})} + \lambda·\pi_{k}\ \ \ \ (*)\)

對k = 1,..,K求和，得：
\(\sum_{n=1}^{N}\sum_{k=1}^{K}\gamma(z_{n}k) + \lambda = 0\),

所以\(\lambda = -N\)，代入(*)式，得：
\(\pi_{k} = \frac{\sum_{n=1}^{N}\gamma(z_{n}k)}{N}\)

Exception Maximizition

期望最大化EM算法的標准形式：

• E-step：來自第j個組份的概率
• \(x_{j}^{i} = Q_{i}(z_{i} = j) = p(z^{i} = j | x^{i};\Phi, \mu, \Sigma)\)
• S-step：估計每個組份的參數。
• \(\sum_{i=1}^{m}\sum_{z^{(i)}}Q_{i}(z^{(i)})log\frac{P(x^{i},z^{i};\Phi,\mu\Sigma)}{Q_{i}(z^{(i)})}\)

EM process:

Initilization:

First, choose some initial values for the means \(\mu\{\mu_{1}, \mu_{2},...,\mu_{K}\}\), covariances \(\Sigma\{\Sigma_{1},\Sigma_{2},...,\Sigma_{K}\}\), and mixing coefficients \(\pi\{\pi_{1},\pi_{2},...,\pi_{K}\}\).

E-Step:

Use the current values for the parameters \(\mu, \Sigma, \pi\) to evaluate the posterior probabilities, or responsibilities \(\gamma(z_{n}k) = \frac{\pi_{k}N(x_{n}|\mu_{k},\Sigma_{k})}{\sum_{j}\pi_{j}N(x_{n}|\mu_{j},\Sigma_{j})}\);

M-step:

Re-estimate the means \(\mu\), covariances \(\Sigma\), and mixing coefficients \(\pi\) using the results:
\(\mu_{k} = \frac{1}{\sum_{n=1}^{N}\gamma(z_{n}k)}\sum_{n=1}^{N}\gamma(z_{n}k)·x_{n}\)

\(\Sigma_{k} = \frac{1}{\sum_{n=1}^{N}\gamma(z_{n}k)}\sum_{n=1}^{N}\gamma(z_{n}k)(x_{n}-\mu_{k})(x_{n}-\mu_{k})^{T}\)

\(\pi_{k} = \frac{\sum_{n=1}^{N}\gamma(z_{n}k)}{N}\)

Evaluate the log likelihood

\(\ln P(X|\pi,\mu,\Sigma) = \sum_{n=1}^{N}\ln\{\sum_{k=1}^{K}\pi_{k}N(x_{n}|\mu_{k},\Sigma_{k})\}\)
and check for convergence of either the parameters or the log likelihood.

隱變量的視角解釋GMM

有時模型內既存在觀測變量(observable variable)，又存在隱變量(latent variable). 當存在隱變量時，直接使用極大似然法或者貝葉斯估計模型參數就比較困難。而EM算法的目標是找到具有潛在變量的模型的最大似然解。
圖模型：

• 記所有的觀測數據的集合為\(X = {x_{1}, x_{2}, ..., x_{n}}\), 所有隱變量的集合為\(Z = {z_{1}, z_{2}, ..., z_{k}}\),模型所有參數的集合為\(\theta\).
  則對數似然函數為：
  \(\ln p(X|\theta) = \ln\{\sum_{z} p(X, Z|\theta) \}\)；
• 但是由於對數似然內部存在求和符號，所以一般情況下邊緣概率分布\(p(X|\theta)\)一定不是指數族函數，從而使得對數似然的形式變復雜了；
• 將\(\{X, Z\}\)數據集稱為完全(complete)數據集，稱實際觀測數據集\(\{X\}\)為不完全(incomplete)數據集.一般我們所擁有的數據集只有不完全數據集。
• E-step: 使用當前的參數\(\theta^{old}\)計算隱變量的后驗概率\(p(Z|Z, \theta^{old})\)，然后計算完全數據(X,Z)對數似然對於參數值\(\theta\)的期望：
  \(Q(\theta, \theta^{old}) = \sum_{z}p(Z|X,\theta^{old})\ln p(X, Z|\theta)\);
• M-step: 通過最大化Q函數更新模型參數\(\theta\):
  \(\theta^{new} = argmax_{\theta}Q(\theta, \theta^{old})\)

引入隱變量z

用1-of-K的形式表示K維的二元隨機變量z，\(z_{k} \in \{0, 1\} 並且 \sum_{k=1}^{K} = 1\)，則\(p(z_{k} = 1) = \pi_{k}\),

• z的概率分布可以表示為：\(p(z) = \Pi_{k=1}^{K}\pi_{k}^{z_{k}}\);
• 在給定隱變量z的條件下x的條件概率分布為：\(p(x|z) = \Pi_{k=1}^{K}N(x|\mu_{k}, \Sigma_{k})^{z_{k}}\)
• 完全數據(X,Z)的似然函數為：\(p(X,Z|\pi,\mu,\Sigma) = P(X|Z,\pi,\mu,\Sigma)·P(Z|\pi,\mu,\Sigma)\\ = \Pi_{n=1}^{N}\Pi_{k=1}^{K}\pi_{k}^{z_{n}k}N(x_{n}|\mu_{k},\Sigma_{k})^{z_{n}k}\)
• 完全數據(X,Z)取對數，得到對數似然：
  \(\ln P(X,Z|\pi,\mu,\Sigma) = \sum_{n=1}^{N}\sum_{k=1}^{K}z_{nk}(\ln\pi_{k} + \ln N(x_{n}|\mu_{k}, \Sigma_{k})) \ \ \ \ \ (1)\)
• 不完全數據(X)的對數似然為：
  \(\ln P(X|\pi,\mu,\Sigma) = \sum_{n=1}^{N}\ln\{\sum_{k=1}^{K}\pi_{k}N(x_{n}|\mu_{k},\Sigma_{k})\} \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ (2)\)

對比完全數據的對數似然函數(1)式和不完全數據的對數似然(2)式，可以發現，對於(1)式，log運算放到了里面，而高斯分布是指數族函數，所以大大簡化了計算量。

• 計算隱變量z的條件概率：\(p(Z|X,\pi,\mu,\Sigma) = \frac{p(Z,X|\pi,\mu,\Sigma)}{p(X|\pi,\mu,\Sigma)} \\ = \Pi_{n=1}^{N}\Pi_{k=1}^{K}[\pi_{k}N(x_{n}|\mu_{k},\Sigma_{k})]^{z_{n}k}\)
• 對\(z_{n}k\)求期望：
  

Exception Maximizition Algorithm

Input: Observed variables X, latent variables Z, parameters \(\theta\), and joint distribution \(p(X,Z|\theta)\).
Goal: maxizimize the likelihood function \(p(X|\theta)\).
Output: parameters \(\theta\).

1. Choose an initial setting for the parameters \(\theta^{old}\).
2. E step Calculate \(Q(\theta, \theta^{old}) = \sum_{z}p(Z|X,\theta^{old})\ln p(X,Z|\theta)\).
3. M step Evaluate \(\theta^{new}\) given by \(\theta^{new} = argmax_{\theta}Q(\theta, \theta^{old})\).
4. Check for convergence of either the log likelihood or the parameter values.

GMM與K-Means的關系

K-Means對數據點進行了硬分配，即每個數據點只屬於唯一的聚類；
而GMM應用EM算法，進行了軟分配。

• 考慮一個GMM，設所有混合分量的協方差矩陣為\(\epsilon I\),則
  \(p(x|\mu_{k},\Sigma_{k}) = \frac{1}{(2\pi\epsilon)^{\frac{D}{2}}}exp\{-\frac{1}{2\epsilon}||x - \mu_{k}||^{2} \}\)
• 一個數據點屬於一個分量的responsibility為：
  
• 當\(\epsilon \to 0\), 分子分母中只有當某項\(||x_{n} - \mu_{j}||^{2} \to 0\)時，整個式子的極限才存在。因此在這種極限情況下，GMM會將一個數據點分配給距它最近的那個中心點所在的簇類中。
• 當\(\epsilon \to 0\),完全數據的對數似然函數為：
  
• 這就退化為K-means的損失函數（目標函數）。

因此，K-means沒有估計聚類的協方差，只是估計了聚類的均值。當GMM的方差取極限時，就退化成了k-means算法。
K-means是GMM的方差取為0的極限情況。

Exception Maximizition

期望最大化EM算法的標准形式：

• E-step：來自第j個組份的概率
• \(x_{j}^{i} = Q_{i}(z_{i} = j) = p(z^{i} = j | x^{i};\Phi, \mu, \Sigma)\)
• S-step：估計每個組份的參數。
• \(\sum_{i=1}^{m}\sum_{z^{(i)}}Q_{i}(z^{(i)})log\frac{P(x^{i},z^{i};\Phi,\mu\Sigma)}{Q_{i}(z^{(i)})}\)

為什么EM算法是可以收斂的

EM算法,是尋找具有潛在變量的概率模型的最大似然解的一種通用的方法.所以，在含有隱變量z的情況下，我們的目標是：最大化似然函數：
\(p(X|\theta) = \sum_{z}p(X, Z|\theta)\)

• 直接優化\(p(X|\theta)\)是困難的，但是優化完整數據的似然\(p(X, Z|\theta)\)就容易的多。

所以，我們引入隱變量的分布\(q(z)\)，將似然函數取對數，總能化簡成如下形式：

因為KL散度總是大於等於0，所以\(L(q,\theta)\)可以看做是對數似然的下界。因此EM算法的思想就是：找到似然函數的一個下界，用這個簡單的下界來逼近最終的對數似然。

• E step: 保持\(\theta\)固定，\(L(q,\theta^{old})\)關於q(z)的最大化。實際是求對z的期望。如下圖所示，期望為：當\(KL(q||p) = 0\)(即p = q)時,\(L(q,\theta^{old})\)與對數似然相等。
  

• M step: 保持\(q(z)\)固定，更新參數\(\theta\)，使得\(L(q,\theta^{old})\)關於\(\theta\)求最大值。對每個參數求偏導，然后令偏導為0. 由於\(KL(q||p) > 0\), 因此\(\ln p(X|\theta)\)增加量一定大於其下界\(L(q,\theta)\)的增大量。

• 通過E step 和 M step不斷迭代的過程，如果目標不完全數據的對數似然\(\ln p(X|\theta)\)存在最大值,則EM算法一定能夠找到。
• EM算法對參數初值的選擇是敏感的，不同的參數往往會收斂到不同的結果。因此在實際應用中，應該設置幾組參數，然后選擇一組最好的結果。

關於EM算法導出的幾點說明

• 在[李航,2012]中，使用\(\ln(\theta) - \in(\theta^{old})\)推導出似然函數的下界。然后每次迭代用下界逼近似然的極大值；
• 在PRML中，將完全數據的似然函數\(L(Z,X|\theta)\)拆分為 下界 + \(KL(q||p)\)的形式，依然是使用下界逼近目標函數的極大值。

[參考文獻]

1 M. Jordan, J. Kleinberg, ect. Pattern Recognition and Machine Learning. 2006
2 李航,統計學習方法,清華大學出版社,2012.