斯坦福大學機器學習，EM算法求解高斯混合模型

本文轉載自查看原文 2017-06-06 01:33 3967 高斯混合模型/ 斯坦福大學/ 機器學習/ 聚類算法/ 人工智能/ EM算法/ 高斯分量/ 建模

斯坦福大學機器學習，EM算法求解高斯混合模型。一種高斯混合模型算法的改進方法---將聚類算法與傳統高斯混合模型結合起來的建模方法, 並同時提出的運用距離加權的矢量量化方法獲取初始值，並采用衡量相似度的方法來融合高斯分量。從對比結果可以看出，基於聚類的高斯混合模型的說話人識別相對於傳統的高斯混合模型在識別率上有所提高。

------------------------------

高斯模型有單高斯模型（SGM）和混合高斯模型（GMM）兩種。

（1）單高斯模型：

為簡單起見，閾值t的選取一般靠經驗值來設定。通常意義下，我們一般取t=0.7-0.75之間。

二維情況如下所示：

（2）混合高斯模型：

對於(b)圖所示的情況，很明顯，單高斯模型是無法解決的。為了解決這個問題，人們提出了高斯混合模型（GMM），顧名思義，就是數據可以看作是從數個高斯分布中生成出來的。雖然我們可以用不同的分布來隨意地構造 XX Mixture Model ，但是 GMM是最為流行。另外，Mixture Model 本身其實也是可以變得任意復雜的，通過增加 Model 的個數，我們可以任意地逼近任何連續的概率密分布。

每個 GMM 由 K 個 Gaussian 分布組成，每個 Gaussian 稱為一個“Component”，這些 Component 線性加成在一起就組成了 GMM 的概率密度函數：

（1）

其中，πk表示選中這個component部分的概率，我們也稱其為加權系數。

根據上面的式子，如果我們要從 GMM 的分布中隨機地取一個點的話，實際上可以分為兩步：

（1）首先隨機地在這 K 個 Component 之中選一個，每個 Component 被選中的概率實際上就是它的系數 πk，選中了 Component 之后，再單獨地考慮從這個 Component 的分布中選取一個點就可以了──這里已經回到了普通的 Gaussian 分布，轉化為了已知的問題。假設現在有 N 個數據點，我們認為這些數據點由某個GMM模型產生，現在我們要需要確定 πk,μk,σk 這些參數。很自然的，我們想到利用最大似然估計來確定這些參數，GMM的似然函數如下：

（2）

在最大似然估計里面，由於我們的目的是把乘積的形式分解為求和的形式，即在等式的左右兩邊加上一個log函數，但是由上文博客里的（2）式可以看出，轉化為log后，還有log(a+b)的形式，因此，要進一步求解。

我們采用EM算法，分布迭代求解最大值：

EM算法的步驟這里不作詳細的介紹，可以參見博客：

http://blog.pluskid.org/?p=39

函數返回的 Px 是一個的矩陣，對於每一個，我們只要取該矩陣第行中最大的那個概率值所對應的那個 Component 為

所屬的 cluster 就可以實現一個完整的聚類方法了。

------------------------------

EM算法（Expection-Maximizationalgorithm，EM）是一種迭代算法，通過E步和M步兩大迭代步驟，每次迭代都使極大似然函數增加。但是，由於初始值的不同，可能會使似然函數陷入局部最優。辜麗川老師和其夫人發表的論文：基於分裂EM算法的GMM參數估計改進了這一缺陷。下面來談談EM算法以及其在求解高斯混合模型中的作用。

一、高斯混合模型（Gaussian MixtureModel， GMM）

高斯判別分析模型，利用參數估計的方法用於解決二分類問題。下面介紹GMM，它是對高斯判別模型的一個推廣，也能借此引入EM算法。

假設樣本集為 $\{x^{(1)},x^{(2)},...,x^{(m)}\}$ 並且樣本和標簽滿足聯合分布 $p(x^{(i)},z^{(i)})=p(x^{(i)}|z^{(i)})p(z^{(i)})$ 。這里： $z^{(i)}$ 服從多項式分布，即 $z^{(i)}\sim Multinomial(\phi)$ （ $\phi_{j}=p(z^{(i)}=j)$ ， $\phi_{j}\ge 0$ ， $\underset{j=1}{\overset{k}{\sum}}\phi_{j}=1$ ），且 $z\in\{1,2,...,k\}$ ；在 $z^{(i)}$ 給定的情況下， $x^{(i)}$ 服從正態分布，即 $x^{(i)}|z^{(i)}=j\sim N(\mu_{j},\Sigma_{j})$ 。這樣的模型稱為高斯混合模型。

該模型的似然函數為：

$\begin{aligned}l(\phi,\mu,\Sigma)&=\underset{i=1}{\overset{m}{\sum}}log\; p(x^{(i)};\phi,\mu,\Sigma)\\ &= \underset{i=1}{\overset{m}{\sum}}log \underset{z^{(i)}=1}{\overset{k}{\sum}}p(x^{(i)}|z^{(i)};\mu,\Sigma)p(z^{(i)};\phi) \end{aligned}$

如果直接令 $l(\phi,\mu,\Sigma)$ 的各變量偏導為0，試圖分別求出各參數，我們會發現根本無法求解。但如果變量 $z^{(i)}$ 是已知的，求解便容易許多，上面的似然函數可以表示為：

$l(\phi,\mu,\Sigma)=\underset{i=1}{\overset{m}{\sum}}log\;p(x^{(i)}|z^{(i)};\mu,\Sigma)+log\;p(z^{(i)};\phi)$

利用偏導求解上述式，可分別得到參數 $\phi,\mu,\Sigma$ 的值：

$\phi_{j}=\frac{1}{m}\underset{i=1}{\overset{m}{\sum}}\#\{z^{(i)}=j\}$

$\mu_{j}=\frac{\underset{i=1}{\overset{m}{\sum}}\#\{z^{(i)}=j\}x^{(i)}}{\underset{i=1}{\overset{m}{\sum}}\#\{z^{(i)}=j\}}$

$\Sigma_{j}=\frac{\underset{i=1}{\overset{m}{\sum}}\#\{z^{(i)}=j\}(x^{(i)}-\mu_{j})(x^{(i)}-\mu_{j})^T}{\underset{i=1}{\overset{m}{\sum}}\#\{z^{(i)}=j\}}$

其中,#{ }為指示函數，表示滿足括號內條件的數目。

那么，變量 $z^{(i)}$ 無法通過觀察直接得到， $z^{(i)}$ 就稱為隱變量，就需要通過EM算法，求解GMM了。下面從Jensen不等式開始，介紹下EM算法：

二、 Jensen不等式（Jensen’s inequality）

引理：如果函數f的定義域為整個實數集，並且對於任意x或 $\overset{\rightarrow}{x}$ 存在 $f''(x)\ge 0$ 或函數的Hessian矩陣 $H\ge0$ ，那么函數f稱為凹函數。或函數的Hessian矩陣H>0，那么函數f為嚴格凹函數。

（存在 $f''(x)\le 0$ 或函數的Hessian矩陣 $H\le0$ ，那么函數f稱為凸函數；如果或函數的Hessian矩陣 H<0，那么函數f為嚴格凸函數。）

定理：如果函數f是凹函數，X為隨機變量，那么：

$E[f(x)]\ge f(EX)$

不幸的是很多人都會講Jensen不等式記混，我們可以通過圖形的方式幫助記憶。下圖中，橫縱坐標軸分別為X和f(X)，f(x)為一個凹函數，a、b分別為變量X的定義域，E[X]為定義域X的期望。圖中清楚的看到各個量的位置和他們間的大小關系。反之，如果函數f是凸函數，X為隨機變量，那么：

$E[f(x)]\le f(EX)$

Jensen不等式等號成立的條件為：，即X為一常數。

三、 EM算法

假設訓練集 $\{x^{(1)},x^{(2)},...,x^{(m)}\}$ 是由m個獨立的樣本構成。我們的目的是要對概率密度函數進行參數估計。它的似然函數為：

$\begin {aligned}l(\theta)&=\underset{i=1}{\overset{m}{\sum}}log\;p(x,\theta)\\ &=\underset{i=1}{\overset{m}{\sum}}log\underset{z}{\sum}p(x,z;\theta) \end{aligned}$

然而僅僅憑借似然函數，無法對參數進行求解。因為這里的隨機變量 $z^{(i)}$ 是未知的。

EM算法提供了一種巧妙的方式，可以通過逐步迭代逼近最大似然值。下面就來介紹下EM算法：

假設對於所有i， $Q_{i}$ 皆為隨機變量 $z^{(i)}$ 的分布函數。即： $\sum_{Z}Q_{i}(z)=1,Q_{i}(z)\ge0$ 。那么：

$\begin{aligned}l(\theta)&=\underset{i}{\sum}log \; p(x^{(i)};\theta)\\ &=\underset{i}{\sum}log\underset{z^{(i)}}{\sum}p(x^{(i)},z^{(i)};\theta)\quad\quad\quad\quad &(1)\\&=\underset{i}{\sum}log\underset{z^{(i)}}{\sum}Q_{i}(z^{(i)})\frac{p(x^{(i)},z^{(i)};\theta)}{Q_{i}(z^{(i)})}&(2)\\ &\ge \underset{i}{\sum}\underset{z^{(i)}}{\sum}Q_{i}(z^{(i)})log\frac{p(x^{(i)},z^{(i)};\theta)}{Q_{i}(z^{(i)})}&(3) \end{aligned}$

其中第（2）步至第（3）步的推導就使用了Jensen不等式。其中：f(x)=log x， $f''(x)=-\frac{1}{x^{2}}<0$ ，因此為凸函數； $\underset{z^{(i)}}{\sum}Q_{i}(z^{(i)})\frac{p(x^{(i)},z^{(i)};\theta)}{Q_{i}(z^{(i)})}$ 表示隨機變量為 $Q_{i}(z^{(i)})$ 概率分布函數為 $\frac{p(x^{(i)},z^{(i)};\theta)}{Q_{i}(z^{(i)})}$ 的期望。因此有：

$f(\underset{z^{(i)}\sim Q_{i}}{E}[\frac{p(x^{(i)},z^{(i)};\theta)}{Q_{i}(z^{(i)})}])\ge \underset{z^{(i)}\sim Q_{i}}{E}[f(\frac{p(x^{(i)},z^{(i)};\theta)}{Q_{i}(z^{(i)})})]$

這樣，對於任意分布 $Q_{i}$ ，（3）都給出了 $l(\theta)$ 的一個下界。如果我們現在通過猜測初始化了一個 $\theta$ 的值，我們希望得到在這個特定的 $\theta$ 下，更緊密的下界，也就是使等號成立。根據Jensen不等式等號成立的條件，當 $\frac{p(x^{(i)},z^{(i)};\theta)}{Q_{i}(z^{(i)})}$ 為一常數時，等號成立。即：

$\frac{p(x^{(i)},z^{(i)};\theta)}{Q_{i}(z^{(i)})} =c$

由上式可得 $p(x^{(i)},z^{(i)};\theta)=cQ_{i}(z^{(i)})$ ，又 $\sum_{Z}Q_{i}(z)=1$ ，因此 $p(x^{(i)},z^{(i)};\theta)=c$ 。再由上式可得：

$\begin {aligned}Q_{i}(z^{(i)})&=\frac{p(x^{(i)},z^{(i)};\theta)}{\sum_{z}{p(x^{(i)},z;\theta)}}\\&=\frac{p(x^{(i)},z^{(i)};\theta)}{{p(x^{(i)};\theta)}}\\&=p(z^{(i)}|x^{(i)};\theta) \end{aligned}$

上述等式最后一步使用了貝葉斯公示。

EM算法有兩個步驟：

（1）通過設置初始化 $\theta$ 值，求出使似然方程最大的 $Q_{i}$ 值，此步驟稱為E-步（E-step）

（2）利用求出的 $Q_{i}$ 值，更新 $\theta$ 。此步驟稱為M-步（M-step）。過程如下：

repeat until convergence{

(E-step) for each i, set

$Q_{i}(z^{(i)}):=p(z^{(i)}|x^{(i)};\theta)$ $Q_{i}(z^{(i)})=p(z^{(i)}|x^{(i)};\theta)$

(M-step) set

$\theta:=arg\,\underset{\theta}{max}\underset{i}{\sum}\underset{z^{(i)}}{\sum}Q_{i}(z^{(i)})log\frac{p(x^{(i)},z^{(i)};\theta)}{Q_{i}(z^{(i)})}$

｝

那么，如何保證EM算法是收斂的呢？下面給予證明：

假設 $\theta^{(t)}$ 和 $\theta^{(t+1)}$ 是EM算法第t次和第t+1次迭代所得到的參數 $\theta$ 的值，如果有 $l(\theta^{(t)})\le l(\theta^{(t+1)})$ ，即每次迭代后似然方程的值都會增大，通過逐步迭代，最終達到最大值。以下是證明：

$\begin{aligned} l(\theta^{(t+1)}) &\ge \underset{i}{\sum}\underset{z^{(i)}}{\sum}Q_{i}^{(t)}(z^{(i)})log\frac{p(x^{(i)},z^{(i)};\theta^{(t+1)})}{Q_{i}^{(t)}(z^{(i)})}\quad\quad&(4)\\ &\ge \underset{i}{\sum}\underset{z^{(i)}}{\sum}Q_{i}^{(t)}(z^{(i)})log\frac{p(x^{(i)},z^{(i)};\theta^{(t)})}{Q_{i}^{(t)}(z^{(i)})}&(5)\\ &=l(\theta^{(t)})&(6) \end {aligned}$

不等式（4）是由不等式（3）得到，對於任意 $Q_{i}$ 和 $\theta$ 值都成立；得到不等式（5）是因為我們需要選擇特定的 $\theta^{(t+1)}$ 使得方程 $arg\,\underset{\theta}{max}\underset{i}{\sum}\underset{z^{(i)}}{\sum}Q_{i}(z^{(i)})log\frac{p(x^{(i)},z^{(i)};\theta)}{Q_{i}(z^{(i)})}$ 在 $\theta^{(t+1)}$ 處的值大於在 $\theta^{(t)}$ 處的值；等式（6）是找到特定的 $\theta^{(t)}$ 的值，使得等號成立。

最后我們通過圖形的方式再更加深入細致的理解EM算法的特點：

$J(Q,\theta)=\underset{i}{\sum}\underset{z^{(i)}}{\sum}Q_{i}(z^{(i)})log\frac{p(x^{(i)},z^{(i)};\theta)}{Q_{i}(z^{(i)})}$

由上文我們知道有這樣的關系： $l(\theta)\ge J(Q,\theta)$ ，EM算法就是不斷最大化這個下界，逐步得到似然函數的最大值。如下圖所示：

首先，初始化 $\theta^{(0)}$ ，調整使得 $J(Q,\theta^{(0)})$ 與 $l(\theta^{(0)})$ 相等，然后求出 $J(Q,\theta^{(0)})$ 使得到最大值的 $\theta^{(1)}$ ；固定 $\theta^{(1)}$ ，調整，使得 $J(Q,\theta^{(1)})$ 與 $l(\theta^{(1)})$ 相等，然后求出使 $J(Q,\theta^{(1)})$ 得到最大值的 $\theta^{(2)}$ ；……；如此循環，使得 $l(\theta)$ 的值不斷上升，直到k次循環后，求出了 $l(\theta)$ 的最大值 $l(\theta^{(k)})$ 。

四、 EM算法應用於混合高斯模型（GMM）

再回到GMM：上文提到由於隱變量 $z^{(i)}$ 的存在，無法直接求解參數，但可以通過EM算法進行求解：

E-Step:

$\omega_{j}^{(i)}=Q_{i}(z^{(i)}=j) =P(z^{(i)}=j|x^{(i)};\phi,\mu,\Sigma)$

M-Step:

$\begin{aligned} &\underset{i=1}{\overset{m}{\sum}}\underset{z^{(i)}}{\sum}Q_{i}(z^{(i)})log \frac{p(x^{(i)},z^{(i)};\phi,\mu,\Sigma)}{Q_{i}(z^{(i)})}\\ = &\underset{i=1}{\overset{m}{\sum}}\underset{j=1}{\overset{k}{\sum}}Q_{i}(z^{(i)})log \frac{p(x^{(i)}|z^{(i)}=j;\mu,\Sigma)p(z^{(i)}=j;\phi)}{Q_{i}(z^{(i)})}\\ = &\underset{i=1}{\overset{m}{\sum}}\underset{j=1}{\overset{k}{\sum}}\omega_{j}^{(i)}log\frac{\frac{1}{(2\pi)^{\frac{n}{2}}\left|\Sigma_{j}\right|^\frac{1}{2}}exp(-\frac{1}{2}(x^{(i)}-\mu_{j})^{T}\Sigma_{j}^{-1}(x^{(i)}-\mu_{j}))\cdot\phi_{j}}{\omega_{j}^{(i)}} \end {aligned}$

（1）參數 $\mu$
對期望 $\mu$ 的每個分量 $\mu_{l}$ 求偏導：

$\begin{aligned}&\triangledown_{\mu_{l}}\underset{i=1}{\overset{m}{\sum}}\underset{j=1}{\overset{k}{\sum}}\omega_{j}^{(i)}log\frac{\frac{1}{(2\pi)^{\frac{n}{2}}\left|\Sigma_{j}\right|^\frac{1}{2}}exp(-\frac{1}{2}(x^{(i)}-\mu_{j})^{T}\Sigma_{j}^{-1}(x^{(i)}-\mu_{j}))\cdot\phi_{j}}{\omega_{j}^{(i)}}\\ =&-\triangledown_{\mu_{l}}\underset{i=1}{\overset{m}{\sum}}\underset{j=1}{\overset{k}{\sum}}\omega_{j}^{(i)}\frac{1}{2}(x^{(i)}-\mu_{j})^{T}\Sigma_{j}^{-1}(x^{(i)}-\mu_{j})\\= &\frac{1}{2}\underset{i=1}{\overset{m}{\sum}}\omega_{l}^{(i)}\Delta_{\mu_{l}}(2\mu_{l}^{T}\Sigma_{l}^{-1}x^{(i)}-\mu_{l}^{T}\Sigma_{l}^{-1}\mu_{l})\\=&\underset{i=1}{\overset{m}{\sum}}\omega_{l}^{(i)}\Sigma_{l}^{-1}(x^{(i)}-\mu_{l}) \end{aligned}$

令上式為0，得：

$\mu_{l}:=\frac{\sum_{i=1}^{m}\omega_{l}^{(i)}x^{(i)}}{\sum_{i=1}^{m}{\omega_{l}^{(i)}}}$

（2）參數

觀察M-Step，可以看到，跟 $\phi_{i}$ 相關的變量僅僅有 $\omega_{j}^{(i)}$ 。因此，我們僅僅需要最大化下面的目標函數：

$max\,\underset{i=1}{\overset{m}{\sum}}\underset{j=1}{\overset{k}{\sum}}\omega_{j}^{(i)}log\phi_{j}$

又由於 $\sum_{j=1}^{k}\phi_{j}=1$ ，為約束條件。因此，可以構造拉格朗日算子求目標函數：

$\mathcal L(\phi)=\underset{i=1}{\overset{m}{\sum}}\underset{j=1}{\overset{k}{\sum}}\omega_{j}^{(i)}log\phi_{j}+\beta(\underset{j=1}{\overset{k}{\sum}}\phi_{j}-1)$

求偏導：

$\frac{\partial}{\partial \phi_{j}} \mathcal L(\phi)=\underset{i=1}{\overset{m}{\sum}}\frac{\omega_{j}{(i)}}{\phi_{j}}+1$

令 $\frac{\partial}{\partial \phi_{j}} \mathcal L(\phi)=0$ 得：

$\phi_{j}:=\frac{\underset{i=1}{\overset{m}{\sum}}\omega_{j}^{(i)}}{-\beta}$

將 $\sum_{j=1}^{k}\phi_{j}=1$ 帶入上式得：

$\begin{aligned} -\beta=\sum_{i=1}^{m}\sum_{j=1}^{k}\omega_{j}^{(i)} &= \sum_{i=1}^{m}1 &=m\end{aligned}$

最后將 $\beta=-m$ 帶入得：

$\phi_{j}:=\frac{\underset{i=1}{\overset{m}{\sum}}\omega_{j}^{(i)}}{m}$

（3）參數 $\Sigma$

$\begin{aligned}&\triangledown _{\Sigma_{j}} \underset{i=1}{\overset{m}{\sum}}\underset{j=1}{\overset{k}{\sum}}\omega_{j}^{(i)}log\frac{\frac{1}{(2\pi)^{\frac{n}{2}}\left|\Sigma_{j}\right|^\frac{1}{2}}exp(-\frac{1}{2}(x^{(i)}-\mu_{j})^{T}\Sigma_{j}^{-1}(x^{(i)}-\mu_{j}))\cdot\phi_{j}}{\omega_{j}^{(i)}}\\ =&\triangledown_{\Sigma_{j}}\underset{i=1}{\overset{m}{\sum}}\underset{j=1}{\overset{k}{\sum}} \omega_{j}^{(i)}[-\frac{n}{2}log2\pi-\frac{1}{2}log\left|\Sigma_{j}\right|-\frac{1}{2}(x^{(i)}-\mu_{j})^{T}\Sigma_{j}^{-1}(x^{(i)}-\mu_{j})+log\phi_{j}-log\omega_{j}^{(i)}]\\=&\underset{i=1}{\overset{m}{\sum}}\underset{j=1}{\overset{m}{\sum}}\omega_{j}^{(i)}\frac{1}{2}\Sigma_{j}^{-T}+\frac{1}{2}(x^{(i)}-\mu_{j})^{T}\Sigma^{-2}(x^{(i)}-\mu_{j}) \end{aligned}$

令上式為零，解得：

$\Sigma_{j}:=\frac{\sum_{i=1}^{m}\omega_{j}^{(i)} (x^{(i)}-\mu_{j})(x^{(i)}-\mu_{j})^{T} }{ \sum_{i=1}^{m}\omega_{j}^{(i)} }$

五、總結

EM算法利用不完全的數據，進行極大似然估計。通過兩步迭代，逐漸逼近最大似然值。而GMM可以利用EM算法進行參數估計。

最后提下辜老師論文的思路：EM模型容易收斂到局部最大值，並且嚴重依賴初試值。傳統的方法即上文中使用的方法是每次迭代過程中，同時更新高斯分布中所有參數，而辜老師的方法是把K個高斯分布中的一個分量，利用奇異值分解的方法將其分裂為兩個高斯分布，並保持其他分量不變的情況下，對共這K+1個高斯分布的權值進行更新，直到符合一定的收斂條件。這樣一來，雖然算法復雜度沒有降低，但每輪只需要更新兩個參數，大大降低了每輪迭代的計算量。

------------------------------

本人微信公眾帳號：心禪道（xinchandao）

本人微信公眾帳號：雙色球預測合買（ssqyuce）

免責聲明！

本站轉載的文章為個人學習借鑒使用，本站對版權不負任何法律責任。如果侵犯了您的隱私權益，請聯系本站郵箱yoyou2525@163.com刪除。

猜您在找 斯坦福大學機器學習筆記及代碼（一） LR 算法總結--斯坦福大學機器學習公開課學習筆記斯坦福大學Andrew Ng教授主講的《機器學習》公開課觀后感 CS229 斯坦福大學機器學習復習材料(數學基礎) - 線性代數機器學習算法總結(六)——EM算法與高斯混合模型機器學習之高斯混合模型及EM算法 [斯坦福大學2014機器學習教程筆記]第四章-正規方程在矩陣不可逆情況下的解決方法（選學）【筆記】邁向人工智能 | Coursera-斯坦福大學-機器學習-吳恩達（已完結）斯坦福大學sql練習（基礎篇）機器學習基礎知識筆記（一）-- 極大似然估計、高斯混合模型與EM算法