详解EM算法与混合高斯模型(Gaussian mixture model, GMM)

 　　最近在看晓川老（shi）师(shu)的博士论文，接触了混合高斯模型（Gaussian mixture model, GMM）和EM(Expectation Maximization)算法，不禁被论文中庞大的数学公式所吓退。本文通过查阅相关资料，在复杂巧妙的推理公式中融入了自己的理解，详细梳理了混合高斯模型和EM算法。

1 单高斯模型(Gaussian single model, GSM）

 　　简单回顾一下概率论讲过的高斯模型。
　　高斯模型是一种常用的变量分布模型，在数理统计领域有着广泛的应用（……好吧读了这么多年书没白费，教科书般的话语已植入骨髓）。一维高斯分布的概率密度函数如下：
$$f(x)=\frac{1}{\sqrt{2\pi }\sigma }\exp (-\frac{{(x-\mu )}^2}{2{\sigma }^2})\text{\hspace{1em}(1)}$$　
　　 $\mu$和${\sigma }^2$分别是高斯分布的均值和方差。
　　譬如将男生身高视为变量X, 假设男生的身高服从高斯分布，则$X\sim N(\mu ,{\sigma }^2)$，女生亦如此。只是男女生身高分布可能具有不同的均值和方差。图1是从谷歌图库中搜索到的男女生身高分布图，来源不清，个人觉得男生的均值身高虚高……四个记号分别表示0~3$\sigma$准则。

图1 男女生身高分布差异

 　　多维变量$X=(x_1,x_2,...x_n)$的联合概率密度函数为：
$$f(X)=\frac{1}{{(2\pi )}^{d/2}{\left|\Sigma \right|}^{1/2}}\exp [-\frac{1}{2}(X-u{)}^T{\Sigma }^{-1}(X-u)],X=(x_1,x_2...x_n)\text{\hspace{1em}(2)}$$
　　其中：
　　d：变量维度。对于二维高斯分布，有d=2;
　　 $u=\left(\begin{array}{l}{u}_1\\ u_2\\ ...\\ u_n\end{array}\right)$：各维变量的均值；
　　 $\Sigma$：协方差矩阵，描述各维变量之间的相关度。对于二维高斯分布，有：
$$\Sigma =\left[\begin{array}{cc}{\delta }_{11}& {\delta }_{12}\\ {\delta }_{21}& {\delta }_{22}\end{array}\right]\text{\hspace{1em}(3)}$$

图2 二维高斯数据分布

 　　图2是二维高斯分布产生的数据示例，参数设定为：$u=\left(\begin{array}{l}0\\ 0\end{array}\right),\Sigma =\left[\begin{array}{cc}1& 0.8\\ 0.8& 5\end{array}\right]$。关于二维高斯分布的参数设定对为高斯曲面的影响，这篇文章二维高斯分布（Two-dimensional Gaussian distribution）的参数分析有提及，主要是为了下文理解混合高斯分布做铺垫。服从二维高斯分布的数据主要集中在一个椭圆内部，服从三维的数据集中在一个椭球内部。

2 混合高斯模型（Gaussian mixture model, GMM）

2.1 为什么要有混合高斯模型

先来看一组数据。　
　　　　　
　　

图3 混合高斯分布所产生数据

 　　如果我们假设这组数据是由某个高斯分布产生的，利用极大似然估计（后文还会提及）对这个高斯分布做参数估计，得到一个最佳的高斯分布模型如下。

图4 用单高斯模型对样本作分析不合理示意

 　　有什么问题吗？一般来讲越靠近椭圆的中心样本出现的概率越大，这是由概率密度函数决定的，但是这个高斯分布的椭圆中心的样本量却极少。显然样本服从单高斯分布的假设并不合理。单高斯模型无法产生这样的样本。
　　实际上，这是用两个不同的高斯分布模型产生的数据。

图5 混合高斯模型对样本作分析示意

 　　正当单高斯模型抓耳挠腮的时候，混合高斯模型就大摇大摆地进场了。它通过求解两个高斯模型，并通过一定的权重将两个高斯模型融合成一个模型，即最终的混合高斯模型。这个混合高斯模型可以产生这样的样本。
　　更一般化的描述为：假设混合高斯模型由K个高斯模型组成（即数据包含K个类），则GMM的概率密度函数如下：
$$p(x)=\sum _{k=1}^{K}p(k)p(x|k)=\sum _{k=1}^K{\pi }_{k}N(x|u_k,{\Sigma }_k)\text{\hspace{1em}(3)}$$
　　其中，$p(x|k)=N(x|u_k,{\Sigma }_k)$是第k个高斯模型的概率密度函数，可以看成选定第k个模型后，该模型产生x的概率；$p(k)={\pi }_k$是第k个高斯模型的权重，称作选择第k个模型的先验概率，且满足$\sum _{k=1}^K{\pi }_k=1$。
　　所以，混合高斯模型并不是什么新奇的东西，它的本质就是融合几个单高斯模型，来使得模型更加复杂，从而产生更复杂的样本。理论上，如果某个混合高斯模型融合的高斯模型个数足够多，它们之间的权重设定得足够合理，这个混合模型可以拟合任意分布的样本。

2.2 直观上理解混合高斯模型

 　　下面通过几张图片来帮助理解混合高斯模型。
　　首先从简单的一维混合高斯模型说起。
　

图6 一维混合高斯模型

 　　在图6中，y1,y2和y3分别表示三个一维高斯模型，他们的参数设定如图所示。y4表示将三个模型的概率密度函数直接相加，注意的是这并不是一个混合高斯模型，因为不满足$\sum _{k=1}^K{\pi }_k=1$的条件。而y5和y6分别是由三个相同的高斯模型融合生成的不同混合模型。由此可见，调整权重将极大影响混合模型的概率密度函数曲线。另一方面也可以直观地理解混合高斯模型可以更好地拟合样本的原因：它有更复杂更多变的概率密度函数曲线。理论上，混合高斯模型的概率密度函数曲线可以是任意形状的非线性函数。
　　下面再给出一个二维空间3个高斯模型混合的例子。
　

(a) 3个类别高斯分布截面轮廓线

(b) 混合高斯分布截面轮廓线

© 二维混合高斯分布概率密度函数图

图7 二维混合高斯模型

 　　(a) 图表示的是3个高斯模型的截面轮廓图，3个模型的权重系数已在图中注明，由截面轮廓图可知3个模型之间存在很多重叠区域。其实这也正是混合高斯模型所希望的。因为如果它们之间的重叠区域较少，那么生成的混合高斯模型一般较为简单，难以生成较为复杂的样本。
　　设定好了3个高斯模型和它们之间的权重系数之后，就可以确定二维混合高斯分布的概率密度函数曲面，如图©所示。图(b)是对于图©概率密度曲面的截面轮廓线。从图7也可以看出，整个混合高斯分布曲面相对比于单高斯分布曲面已经异常复杂。实际上，通过调整混合高斯分布的系数$(\pi ,\mu ,\Sigma )$，可以使得图©的概率密度曲面去拟合任意的三维曲面，从而采样生成所需要的数据样本。

3 极大似然估计(Maximum Likehood Estimate, MLE)（最大化对数似然函数）

● 最大化对数似然函数（log-likelihood function）的意义

 　　首先直观化地解释一下最大化对数似然函数要解决的是什么问题。
　　假设我们采样得到一组样本$y_t$,而且我们知道变量Y服从高斯分布（本文只提及高斯分布，其他变量分布模型类似），数学形式表示为$Y\sim N(\mu ,\Sigma )$。采样的样本如图8所示，我们的目的就是找到一个合适的高斯分布（也就是确定高斯分布的参数$\mu ,\Sigma$），使得这个高斯分布能产生这组样本的可能性尽可能大。
　

图8 最大化似然函数的意义

 　　那怎么找到这个合适的高斯分布呢（在图8中的表示就是1~4哪个分布较为合适）？这时候似然函数就闪亮登场了。
　　似然函数数学化：设有样本集$Y=y_1,y_2...y_N$。$p(y_n|\mu ,\Sigma )$是高斯分布的概率分布函数，表示变量$Y=y_n$的概率。假设样本的抽样是独立的，那么我们同时抽到这N个样本的概率是抽到每个样本概率的乘积，也就是样本集Y的联合概率。此联合概率即为似然函数：
$$L(\mu ,\Sigma )=L(y_1,y_2...y_N;\mu ,\Sigma )=\prod _{n=1}^{N}p(y_n;\mu ,\Sigma )\text{\hspace{1em}(4)}$$
　　对式子(4)进行求导并令导数为0（即最大化似然函数，一般还会先转化为对数似然函数再最大化），所求出的参数就是最佳的高斯分布对应的参数。
　　所以最大化似然函数的意义就是：通过使得样本集的联合概率最大来对参数进行估计，从而选择最佳的分布模型。
　　对于图8产生的样本用最大化似然函数的方法，最终可以得到序号1对应的高斯分布模型是最佳的模型。

4 EM算法（最大化Q函数）

4.1 为什么要有EM算法（EM算法与极大似然估计分别适用于什么问题）

● 尝试用极大似然估计的方法来解GMM模型

 　　解GMM模型，实际上就是确定GMM模型的参数$(\mu ,\Sigma ,\pi )$，使得由这组参数确定的GMM模型最有可能产生采样的样本。
　　先试试看用极大似然估计的方法来解GMM模型会出现什么样的问题。
　　如第3小节所述，要利用极大似然估计求解模型最重要的一步就是求出似然函数，即样本集出现的联合概率。而对于混合高斯模型，如何求解某个样本$y_t$的概率？显然我们得先知道这个样本来源于哪一类高斯模型，然后求这个高斯模型生成这个样本的概率$p(y_t)$。
　　但是问题来了：我们只有样本。不知道样本到底来源于哪一类的高斯模型。那么如何求解样本的生成概率$p(y_t)$？
　　先引入一个隐变量$\gamma$。它是一个K维二值随机变量，在它的K维取值中只有某个特定的元素${\gamma }_k$的取值为1，其它元素的取值为0。实际上，隐变量描述的就是：每一次采样，选择第k个高斯模型的概率，故有：
$$p({\gamma }_k=1)={\pi }_k\text{\hspace{1em}(5)}$$
　　当给定了$\gamma$的一个特定的值之后（也就是知道了这个样本从哪一个高斯模型进行采样），可以得到样本y的条件分布是一个高斯分布，满足：
$$p(y|{\gamma }_k=1)=N(y|{\mu }_k,{\Sigma }_k)\text{\hspace{1em}(6)}$$
　　而实际上，每个样本到底是从这K个高斯模型中哪个模型进行采样的，是都有可能的。故样本y的概率为：
(7) p ( y ) = ∑ γ p ( γ ) p ( y ∣ γ ) = ∑ k = 1 K π k N ( y ∣ μ k , Σ k ) p(y) = \sum\nolimits_\gamma {p(\gamma )} p(y|\gamma ){\rm{ = }}\sum\limits_{ {\rm{k}} = 1}^K { {\pi _k}N(y|{\mu _k},{\Sigma _k})} \tag{7} p(y)=∑γ​p(γ)p(y∣γ)=k=1∑K​πk​N(y∣μ