GMM，即高斯混合模型（Gaussian Mixture Model），簡單地講，就是將多個高斯模型混合起來，作為一個新的模型，這樣就可以綜合運用多模型的表達能力。EM，指的是均值最大化算法（expectation-maximization），它是一種估計模型參數的策略，在 GMM 這類算法中應用廣泛，因此，有時候人們又喜歡把 GMM 這類可以用 EM 算法求解的模型稱為 EM 算法家族。

這篇文章會簡單提一下 GMM 模型的內容，最主要的，還是講一下 EM 算法如何應用到 GMM 模型的參數估計上。

高斯混合模型

什麽是 GMM

GMM 可以認為是 K-means 算法的升級版。在 K-means 中，我們會先計算出幾個聚類中心，然後根據數據點與聚類中心的距離，直接將數據點歸類到最近的聚類中心。這種做法其實很“硬”，因為有很多邊緣點屬於兩個聚類中心的概率可能相差不大，如果一股腦就直接將它歸到某一個中心，實在是太粗暴了。而 GMM 不同於 K-means 的地方就在於，它除了給出聚類中心外，還能告訴你每個點歸屬於某個聚類中心的概率，因此，GMM 又被稱作 soft assignment。

首先，還是給出 GMM 模型的公式：
\[ p( x)=\sum_{k=1}^K{\pi_k N( x|\mu_k, \Sigma_k)} \]
其中，我們規定，\(\sum_{k=1}^K{\pi_k}=1\)。可以看出，GMM 就是將幾個高斯模型線性組合起來，人們習慣上把這裏面的各個高斯模型稱為 Component。其中，\(\pi_k\) 表示每個模型的占比，或者說數據屬於模型 k 的概率，這個值越大，說明聚集在這個模型內的數據越多。

為什麽要用這種模型組合的方式呢？我們知道，高斯模型一般成橢圓狀（二維）或橢球狀（三維），可以把這個橢圓或橢球認為是一種聚類的形狀，而圓心或球心則是聚類中心（對應高斯函數的參數 \(\mu\)

其實，這種組合模型的思路可以應用到很多模型上，比如：泊松模型。而由於高斯模型本身一些良好的性質，因此 GMM 這種模型被用的比較多。

前面說到，GMM 本質上是一種聚類算法，那麽，如果已知一個 GMM 模型，現在給定一個點，我們要怎麽知道這個點屬於哪個聚類中心呢？更具體一點說，怎麽知道這個點屬於每個聚類中心的概率是多少？

用數學的語言表達就是，已知一個 GMM 模型： \(p( x)=\sum_{k=1}^K{\pi_k N( x|\mu_k, \Sigma_k)}\)

求解的方法很簡單，根據貝葉斯公式：\(p(a|b)=\frac{p(b|a)p(a)}{p(b)}\)，我們可以得出：
\[ p( x \in C_k | x)=\frac{p(C_k)p( x| x\in C_k)}{p( x)} \]
因此，對於每個聚類中心 \(C_k\)，由於分母 \(p( x)\) 都是相同的，我們只需要計算 \(p(C_k)p( x| x\in C_k)=\pi_k N( x|\mu_k, \Sigma_k)\) 即可。得到的值就是數據點 $ x$ 屬於 \(C_k\) 的概率，至於具體要將 $ x$ 歸類到哪個中心，可以根據具體情況決定，比如將概率最大的作為歸屬的聚類中心。這一點也是 GMM 優於 K-means 的地方，前者是通過概率的方式來決定歸屬，因此提供了更加豐富的信息。

參數估計

不過，GMM 模型最難的地方在於，如何根據一堆數據點估計出模型的參數？

GMM 需要確定的參數有三類：

1. 高斯模型的個數 \(K\)，這個參數跟 K-means 的 \(K\) 一樣，需要人工事先設定，\(K\) 越大，聚類的粒度也越細；
2. \(\pi_k\)， 每個 Component 的概率分量，或者說在總樣本中的占比；
3. \(\mu_k\)、\(\Sigma_k\)，各個 Component 的參數。

如果樣本所屬分類已知（即已知 \(x\) 屬於哪個 \(C_k\)），那 GMM 的參數就很容易確定了。首先，參數 \(K\) 就一目了然得到了。然後，假設樣本容量為 \(N\)，歸屬於聚類中心 \(C_k\) 的樣本數量為 \(N_k\)，歸屬每個 \(C_k\) 的樣本集合為 \(S(k)\)，可以用以下公式求出其他參數：
\[ \pi_k=\frac{N_k}{N} \\mu_k=\frac{1}{N_k}\sum_{ x\in S(k)}{ x} \\Sigma_k=\frac{1}{N_k}\sum_{ x\in S(k)}{( x-\mu_k)( x-\mu_k)^T} \]
其實，這跟一個高斯模型的情況是一樣的，只不過要依葫蘆畫瓢求出 \(K\) 個。

但如果樣本的分類事先不知道，又該怎麽辦呢？首先，由於 \(K\) 這個值是需要人工確定的，所以這裏暫時假設 \(K\) 已經知道了。現在，我們要預測 \(K\) 個高斯模型的概率分量 \(\pi_k\) 以及每個模型各自的參數 \(\mu_k\) 和 \(\Sigma_k\)。

最簡單也最容易想到的方法是極大似然估計。假設有 m 個樣本，首先，寫出 \(p( x)=\sum_{k=1}^K{\pi_k N( x|\mu_k, \Sigma_k)}\) 的似然函數：
\[ \begin{eqnarray} \ln{[\prod_{i=1}^m p( x_i)]}&=&\ln{[\prod_{i=1}^m{\sum_{k=1}^K{\pi_k N( x|\mu_k, \Sigma_k)}}]} \&=&\sum_{i=1}^m{\ln{[\sum_{k=1}^K{\pi_k N( x|\mu_k, \Sigma_k)]}}} \\end{eqnarray} \]
不過，這個對數函數卻出奇的復雜，直接求導數的方法很難求出 \(\mu_k\) 和 \(\Sigma_k\)，因此，我們只能換用其他方式來求解。而這就是 EM 算法發揮作用的地方。

均值最大化算法 EM

K-means 的啟示

在正式開講 EM 之前，我們先回憶一下，K-means 是怎麽求出聚類中心的。其實，總共分三步進行：

1. 隨機初始化 K 個聚類中心的位置；
2. 將所有樣本點，按照跟各個聚類中心的距離進行歸類；
3. 根據樣本重新歸類的結果，更新聚類中心的位置，重復步驟 2 直到收斂（即聚類中心重新調整的幅度小於某個閾值）。

既然 GMM 本身也屬於一種聚類算法，那麽，我們能不能用 K-means 的思路來求出 GMM 的參數呢？答案當然是可以的。

不過，在這之前，我們需要先知道 GMM 的幾個參數（\(\pi_k\)，\(\mu_k\)，\(\Sigma_k\)）要怎麽計算。假設我們已經知道了後驗概率 \(P( x \in C_k| x)\)，則可以根據以下公式計算參數（其中，m 表示樣本數量）：
\[ \pi_k=\frac{1}{n}\sum_{i=1}^m{P( x_i\in C_k|x_i)} \tag{3} \]
這個公式是把所有樣本屬於 \(C_k\) 的概率求平均後，作為 \(C_k\) 這個聚類中心（或者說這個高斯模型）的出現概率。
\[ \mu_k=\frac{\sum_{i=1}^m{ xP(x_i\in C_k|x_i)}}{\sum_{i=1}^m{P(x_i\in C_k|x_i)}} \tag{4} \]
這個求均值的公式，跟單個高斯模型不同的地方在於，我們用的是加權平均。因為每個樣本點都有一定的概率屬於聚類中心 \(C_k\)，所以，每個樣本對 \(C_k\) 對應的高斯模型的均值也會產生一定的作用，只是由於 \(P(x_i\in C_k|x_i)\) 的值不同，因此這種作用也會有顯著差別。
\[ \Sigma_k=\frac{\sum_{i=1}^m{P(x_i\in C_k|x_i)(x_i-\mu_k)(x_i-\mu_k)^T}}{\sum_{i=1}^m{P(x_i\in C_k|x_i)}} \tag{5} \]
類似地，協方差也是用加權平均求出來的。

（公式 (3) (4) (5) 其實是從極大似然函數推出來的，在周誌華老師的西瓜書和PRML書中都有詳細推導，不過這裏我只想給出感性的認識）

不過，以上公式都是基於 \(P( x \in C_k| x)=\frac{p(C_k)p( x| x\in C_k)}{p( x)}\)計算出來的，而這個公式本身又需要知道 \(P(C_k)\)（即 \(\pi_k\)）等參數，這就陷入一個雞生蛋還是蛋生雞的怪圈。

但是，借助 K-means 的思路，我們可以事先隨機初始化這些參數，然後計算出 \(P( x \in C_k| x)\)，再用它更新 GMM 參數，然後再用更新的模型計算 \(P( x \in C_k| x)\)，如此叠代下去，總有收斂的時候，這樣，我們不就可以像 K-means 一樣計算出參數了嗎？！

下面，我們就仿照 K-means 的方法，給出叠代計算 GMM 參數的步驟：

1. 隨機初始化各個高斯模型的參數；
2. 根據參數，計算 \(P( x \in C_k| x)\)，這一步其實是計算出每一個樣本歸屬於每一個聚類中心的概率；
3. 根據第 2 步計算得到的 \(P( x \in C_k| x)\)，按照公式 (3) (4) (5) 重新計算 GMM 參數，並重復步驟 2 直到收斂。

EM 算法

其實，上面仿照 K-means 算法的計算步驟，就是 EM 算法的核心部分了。

EM 算法主要分為 E 和 M 兩步：

1. E 指的是 Expectation，即計算均值的過程，對應的是上面的步驟 2，這一步主要是計算每個樣本歸屬的聚類中心；
2. M 指的是 Maximum，即對參數的最大似然估計，對應的是上面的步驟 3。我前面也說了，公式 (3) (4) (5) 計算參數的公式是用最大似然函數推出來的，所以，這一步其實是在根據步驟 2 的分類結果，重新用最大似然函數來估計參數。

下面這幅圖是從西瓜書上截下來的，是 EM 算法求解 GMM 參數的完整過程。

圖中有很多公式標記，可能要參考原書才看得懂，不過，它的流程和我之前給出的 3 個步驟是一致。另外，算法的停止條件可以是達到最大叠代次數，或者是似然函數（公式 (2)）的增長小於某個閾值。

好了，本文到這裏就不再深入下去了。EM 算法博大精深，吳軍老師在《數學之美》稱它為上帝的算法，可見這個算法的強大之處。EM 算法可以應用的場合非常多，從本文給出的 GMM 例子來看，它其實很類似梯度下降算法，在給定的目標函數很復雜、難以求解時，EM 算法用一種叠代的策略來優化初始參數，並逐步收斂到最優解。關於這個算法的具體內容，我想進一步深入理解後，再用一篇文章好好寫一下。

參考

• 漫談 Clustering (3): Gaussian Mixture Model
• Mixture of Gaussian clustering
• EM及高斯混合模型
• 機器學習西瓜書
• PRML

GMM與EM共舞