聚類之高斯混合模型（Gaussian Mixture Model）【轉】

時間 2019-12-08

標籤高斯混合模型 gaussian mixture model 简体版

原文原文鏈接

k-means應該是原來級別的聚類方法了，這整理下一個使用後驗機率準確評測其精度的方法—高斯混合模型。算法

咱們談到了用 k-means 進行聚類的方法，此次咱們來講一下另外一個很流行的算法：Gaussian Mixture Model (GMM)。事實上，GMM 和 k-means 很像，不過 GMM 是學習出一些機率密度函數來（因此 GMM 除了用在 clustering 上以外，還常常被用於 density estimation ），簡單地說，k-means 的結果是每一個數據點被 assign 到其中某一個 cluster 了，而 GMM 則給出這些數據點被 assign 到每一個 cluster 的機率，又稱做 soft assignment 。dom

得出一個機率有不少好處，由於它的信息量比簡單的一個結果要多，好比，我能夠把這個機率轉換爲一個 score ，表示算法對本身得出的這個結果的把握。也許我能夠對同一個任務，用多個方法獲得結果，最後選取「把握」最大的那個結果；另外一個很常見的方法是在諸如疾病診斷之類的場所，機器對於那些很容易分辨的狀況（患病或者不患病的機率很高）能夠自動區分，而對於那種很難分辨的狀況，好比，49% 的機率患病，51% 的機率正常，若是僅僅簡單地使用 50% 的閾值將患者診斷爲「正常」的話，風險是很是大的，所以，在機器對本身的結果把握很小的狀況下，會「拒絕發表評論」，而把這個任務留給有經驗的醫生去解決。機器學習

廢話說了一堆，不過，在回到 GMM 以前，咱們再稍微扯幾句。咱們知道，無論是機器仍是人，學習的過程均可以看做是一種「概括」的過程，在概括的時候你須要有一些假設的前提條件，例如，當你被告知水裏遊的那個傢伙是魚以後，你使用「在一樣的地方生活的是同一種東西」這相似的假設，概括出「在水裏遊的都是魚」這樣一個結論。固然這個過程是徹底「本能」的，若是不仔細去想，你也不會了解本身是怎樣「認識魚」的。另外一個值得注意的地方是這樣的假設並不老是徹底正確的，甚至能夠說老是會有這樣那樣的缺陷的，所以你有可能會把蝦、龜、甚至是潛水員當作魚。也許你以爲能夠經過修改前提假設來解決這個問題，例如，基於「生活在一樣的地方而且穿着一樣衣服的是同一種東西」這個假設，你得出結論：在水裏有而且身上長有鱗片的是魚。但是這樣仍是有問題，由於有些沒有長鱗片的魚如今又被你排除在外了。函數

在這個問題上，機器學習面臨着和人同樣的問題，在機器學習中，一個學習算法也會有一個前提假設，這裏被稱做「概括偏執 (bias)」（bias 這個英文詞在機器學習和統計裏還有其餘許多的意思）。例如線性迴歸，目的是要找一個函數儘量好地擬合給定的數據點，它的概括偏執就是「知足要求的函數必須是線性函數」。一個沒有概括偏執的學習算法從某種意義上來講毫無用處，就像一個徹底沒有概括能力的人同樣，在第一次看到魚的時候有人告訴他那是魚，下次看到另外一條魚了，他並不知道那也是魚，由於兩條魚總有一些地方不同的，或者就算是同一條魚，在河裏不一樣的地方看到，或者只是看到的時間不同，也會被他認爲是不一樣的，由於他沒法概括，沒法提取主要矛盾、忽略次要因素，只好要求全部的條件都徹底同樣──然而哲學家已經告訴過咱們了：世界上不會有任何樣東西是徹底同樣的，因此這我的即便是有無比強悍的記憶力，也絕學不到任何一點知識。學習

這個問題在機器學習中稱做「過擬合 (Overfitting)」，例如前面的迴歸的問題，若是去掉「線性函數」這個概括偏執，由於對於 N 個點，咱們老是能夠構造一個 N-1 次多項式函數，讓它完美地穿過全部的這 N 個點，或者若是我用任何大於 N-1 次的多項式函數的話，我甚至能夠構造出無窮多個知足條件的函數出來。若是假定特定領域裏的問題所給定的數據個數老是有個上限的話，我能夠取一個足夠大的 N ，從而獲得一個（或者無窮多個）「超級函數」，可以 fit 這個領域內全部的問題。然而這個（或者這無窮多個）「超級函數」有用嗎？只要咱們注意到學習的目的（一般）不是解釋現有的事物，而是從中概括出知識，並能應用到新的事物上，結果就顯而易見了。ui

沒有概括偏執或者概括偏執太寬泛會致使 Overfitting ，然而另外一個極端──限制過大的概括偏執也是有問題的：若是數據自己並非線性的，強行用線性函數去作迴歸一般並不能獲得好結果。難點正在於在這之間尋找一個平衡點。不過人在這裏相對於（如今的）機器來講有一個很大的優點：人一般不會孤立地用某一個獨立的系統和模型去處理問題，一我的天天都會從各個來源獲取大量的信息，而且經過各類手段進行整合處理，概括所得的全部知識最終得以統一地存儲起來，並能有機地組合起來去解決特定的問題。這裏的「有機」這個詞頗有意思，搞理論的人總能提出各類各樣的模型，而且這些模型都有嚴格的理論基礎保證能達到指望的目的，然而絕大多數模型都會有那麼一些「參數」（例如 K-means 中的 k ），一般沒有理論來講明參數取哪一個值更好，而模型實際的效果卻一般和參數是否取到最優值有很大的關係，我以爲，在這裏「有機」不妨看做是全部模型的參數已經自動地取到了最優值。另外，雖然進展不大，可是人們也一直都指望在計算機領域也創建起一個統一的知識系統（例如語意網就是這樣一個嘗試）。atom

廢話終於說完了，回到 GMM 。按照咱們前面的討論，做爲一個流行的算法，GMM 確定有它本身的一個至關體面的概括偏執了。其實它的假設很是簡單，顧名思義，Gaussian Mixture Model ，就是假設數據服從 Mixture Gaussian Distribution ，換句話說，數據能夠看做是從數個 Gaussian Distribution 中生成出來的。實際上，咱們在 K-means 和 K-medoids 兩篇文章中用到的那個例子就是由三個 Gaussian 分佈從隨機選取出來的。實際上，從中心極限定理能夠看出，Gaussian 分佈（也叫作正態 (Normal) 分佈）這個假設實際上是比較合理的，除此以外，Gaussian 分佈在計算上也有一些很好的性質，因此，雖然咱們能夠用不一樣的分佈來隨意地構造 XX Mixture Model ，可是仍是 GMM 最爲流行。另外，Mixture Model 自己其實也是能夠變得任意複雜的，經過增長 Model 的個數，咱們能夠任意地逼近任何連續的機率密分佈。spa

每一個 GMM 由 K 個 Gaussian 分佈組成，每一個 Gaussian 稱爲一個「Component」，這些 Component 線性加成在一塊兒就組成了 GMM 的機率密度函數： .net

根據上面的式子，若是咱們要從 GMM 的分佈中隨機地取一個點的話，實際上能夠分爲兩步：首先隨機地在這 K 個 Component 之中選一個，每一個 Component 被選中的機率實際上就是它的係數 $π_{k}$ scala

那麼如何用 GMM 來作 clustering 呢？其實很簡單，如今咱們有了數據，假定它們是由 GMM 生成出來的，那麼咱們只要根據數據推出 GMM 的機率分佈來就能夠了，而後 GMM 的 K 個 Component 實際上就對應了 K 個 cluster 了。根據數據來推算機率密度一般被稱做 density estimation ，特別地，當咱們在已知（或假定）了機率密度函數的形式，而要估計其中的參數的過程被稱做「參數估計」。

如今假設咱們有 N 個數據點，並假設它們服從某個分佈（記做 p(x) ），如今要肯定裏面的一些參數的值，例如，在 GMM 中，咱們就須要肯定 $π_{k}$

下面讓咱們來看一看 GMM 的 log-likelihood function ：

\sum i = 1 N log {\sum k = 1 K π k N (x i | μ k, Σ k)}

$π_{k}$

估計數據由每一個 Component 生成的機率（並非每一個 Component 被選中的機率）：對於每一個數據 x_i 來講，它由第 k 個 Component 生成的機率爲
$γ (i, k) = π k N ( x i | μ k , Σ k ) \sum K j = 1 π j$
因爲式子裏的 $μ_{k} 和 Σ_{k}$
估計每一個 Component 的參數：如今咱們假設上一步中獲得的 $γ (i, k)$ $μ k Σ k = 1 N k \sum i = 1 N γ ( i , k ) x i = 1$ $γ (i, k)$

重複迭代前面兩步，直到似然函數的值收斂爲止。
固然，上面給出的只是比較「直觀」的解釋，想看嚴格的推到過程的話，能夠參考 Pattern Recognition and Machine Learning 這本書的第九章。有了實際的步驟，再實現起來就很簡單了。Matlab 代碼以下：

（Update 2012.07.03：若是你直接把下面的代碼拿去運行了，碰到 covariance 矩陣 singular 的狀況，能夠參見這篇文章。）

function varargout = gmm(X, K_or_centroids) % ============================================================ % Expectation-Maximization iteration implementation of % Gaussian Mixture Model. % % PX = GMM(X, K_OR_CENTROIDS) % [PX MODEL] = GMM(X, K_OR_CENTROIDS) % % - X: N-by-D data matrix. % - K_OR_CENTROIDS: either K indicating the number of % components or a K-by-D matrix indicating the % choosing of the initial K centroids. % % - PX: N-by-K matrix indicating the probability of each % component generating each point. % - MODEL: a structure containing the parameters for a GMM: % MODEL.Miu: a K-by-D matrix. % MODEL.Sigma: a D-by-D-by-K matrix. % MODEL.Pi: a 1-by-K vector. % ============================================================ threshold = 1e-15; [N, D] = size(X); if isscalar(K_or_centroids) K = K_or_centroids; % randomly pick centroids rndp = randperm(N); centroids = X(rndp(1:K), :); else K = size(K_or_centroids, 1); centroids = K_or_centroids; end % initial values [pMiu pPi pSigma] = init_params(); Lprev = -inf; while true Px = calc_prob(); % new value for pGamma pGamma = Px .* repmat(pPi, N, 1); pGamma = pGamma ./ repmat(sum(pGamma, 2), 1, K); % new value for parameters of each Component Nk = sum(pGamma, 1); pMiu = diag(1./Nk) * pGamma' * X; pPi = Nk/N; for kk = 1:K Xshift = X-repmat(pMiu(kk, :), N, 1); pSigma(:, :, kk) = (Xshift' * ... (diag(pGamma(:, kk)) * Xshift)) / Nk(kk); end % check for convergence L = sum(log(Px*pPi')); if L-Lprev < threshold break; end Lprev = L; end if nargout == 1 varargout = {Px}; else model = []; model.Miu = pMiu; model.Sigma = pSigma; model.Pi = pPi; varargout = {Px, model}; end function [pMiu pPi pSigma] = init_params() pMiu = centroids; pPi = zeros(1, K); pSigma = zeros(D, D, K); % hard assign x to each centroids distmat = repmat(sum(X.*X, 2), 1, K) + ... repmat(sum(pMiu.*pMiu, 2)', N, 1) - ... 2*X*pMiu'; [dummy labels] = min(distmat, [], 2); for k=1:K Xk = X(labels == k, :); pPi(k) = size(Xk, 1)/N; pSigma(:, :, k) = cov(Xk); end end function Px = calc_prob() Px = zeros(N, K); for k = 1:K Xshift = X-repmat(pMiu(k, :), N, 1); inv_pSigma = inv(pSigma(:, :, k)); tmp = sum((Xshift*inv_pSigma) .* Xshift, 2); coef = (2*pi)^(-D/2) * sqrt(det(inv_pSigma)); Px(:, k) = coef * exp(-0.5*tmp); end end end

函數返回的 Px 是一個 $N \times K$

相對於以前 K-means 給出的結果，這裏的結果更好一些，左下角的比較稀疏的那個 cluster 有一些點跑得比較遠了。固然，由於這個問題本來就是徹底有 Mixture Gaussian Distribution 生成的數據，GMM （若是能求得全局最優解的話）顯然是能夠對這個問題作到的最好的建模。

另外，從上面的分析中咱們能夠看到 GMM 和 K-means 的迭代求解法其實很是類似（均可以追溯到 EM 算法，下一次會詳細介紹），所以也有和 K-means 一樣的問題──並不能保證老是能取到全局最優，若是運氣比較差，取到很差的初始值，就有可能獲得不好的結果。對於 K-means 的狀況，咱們一般是重複必定次數而後取最好的結果，不過 GMM 每一次迭代的計算量比 K-means 要大許多，一個更流行的作法是先用 K-means （已經重複並取最優值了）獲得一個粗略的結果，而後將其做爲初值（只要將 K-means 所得的 centroids 傳入 gmm 函數便可），再用 GMM 進行細緻迭代。

如咱們最開始所討論的，GMM 所得的結果（Px）不只僅是數據點的 label ，而包含了數據點標記爲每一個 label 的機率，不少時候這其實是很是有用的信息。最後，須要指出的是，GMM 自己只是一個模型，咱們這裏給出的迭代的辦法並非惟一的求解方法。感興趣的同窗能夠自行查找相關資料。

來源：https://blog.csdn.net/u014665013/article/details/78970184