高斯混合模型（GMM）

最近在实际工作中用到了高斯混合模型（Gaussian Mixture Model），遂写一篇笔记来整理记录相关知识点，以便复查巩固。

1. 高斯模型(Gaussian model，GM）

简单回顾一下本科概率论讲过的高斯模型。
　　高斯模型是一种常用的变量分布模型，又称正态分布，在数理统计领域有着广泛的应用。
当样本数据 X 是一维数据（Univariate）时，高斯分布遵从下方概率密度函数（Probability Density Function）(下文简称pdf)如下： $p(x|\theta)= \frac{1}{\sqrt{2π}\alpha}exp(-\frac{(x-u)^2}{2\alpha^2})$ 其中 $u$ 为数据均值（期望）， $\alpha$ 为数据标准差（Standard deviation）。
当样本数据 X 是多维数据（Multivariate）时，高斯分布pdf为：
$p(x|\theta)= \frac{1}{(2π)^{\frac{D}{2}}|\sum|^{\frac{1}{2}}}exp(-\frac{(x-u)^T\sum^{-1}(x-u)}{2})$ 其中， $u$ 为数据均值（期望）， $\sum$ 为协方差（Covariance），描述各维变量之间的相关度，D 为数据维度。

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

高斯混合模型可以看作是由 K 个单高斯模型组合而成的模型，这 K 个子模型是混合模型的隐变量（Hidden variable）。一般来说，一个混合模型可以使用任何概率分布，这里使用高斯混合模型是因为高斯分布具备很好的数学性质以及良好的计算性能。

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

先来看一组数据。

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

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

它通过求解两个高斯模型，并通过一定的权重将两个高斯模型融合成一个模型，即最终的混合高斯模型。这个混合高斯模型可以产生这样的样本。
　　更一般化的描述为：假设混合高斯模型由K个高斯模型组成（即数据包含K个类），则GMM的概率密度函数如下：其中

$p(x|\theta)=phi(x|\theta_k)$ 是第k个高斯模型的概率密度函数，可以看成选定第k个模型后，该模型产生x的概率；
$\alpha_k$ 是第k个高斯模型的权重，称作选择第k个模型的先验概率，且满足 $\sum_{k=1}^{K}\alpha_k=1$ 。

所以，混合高斯模型并不是什么新奇的东西，它的本质就是融合几个单高斯模型，来使得模型更加复杂，从而产生更复杂的样本。理论上，如果某个混合高斯模型融合的高斯模型个数足够多，它们之间的权重设定得足够合理，这个混合模型可以拟合任意分布的样本。

3. 模型参数学习

对于单高斯模型，我们可以用最大似然法（Maximum likelihood）估算参数 $\theta$ 的值 $\theta = argmax_{\theta}L(\theta)$ 这里我们假设了每个数据点都是独立的（Independent），似然函数由概率密度函数（PDF）给出。
$L(\theta)=\prod_{j=1}^{N}p(x_j|\theta)$ 由于每个点发生的概率都很小，乘积会变得极其小，不利于计算和观察，因此通常我们用 Maximum Log-Likelihood 来计算（因为 Log 函数具备单调性，不会改变极值的位置，同时在 0-1 之间输入值很小的变化可以引起输出值相对较大的变动）： $logL(\theta)=\sum_{j=1}^{N}logp(x_j|\theta)$ 对其进行求导并令导数为0，所求出的参数就是最佳的高斯分布对应的参数。
　　所以最大化似然函数的意义就是：通过使得样本集的联合概率最大来对参数进行估计，从而选择最佳的分布模型。
对于高斯混合模型，Log-Likelihood 函数是： $logL(\theta)=\sum_{j=1}^{N}logp(x_j|\theta)=\sum_{j=1}^{N}log\sum_{k=1}^{K}\alpha_k\phi(x|\theta_k)$ 如何计算高斯混合模型的参数呢？这里我们无法像单高斯模型那样使用最大似然法来求导求得使 likelihood 最大的参数，因为对于每个观测数据点来说，事先并不知道它是属于哪个子分布的（hidden variable），因此 log 里面还有求和，对于每个子模型都有未知的 $\alpha_k, u_k, \delta_k$ ，直接求导无法计算。需要通过迭代的方法求解。

4. EM算法

EM 算法是一种迭代算法，1977 年由 Dempster 等人总结提出，用于含有隐变量（Hidden variable）的概率模型参数的最大似然估计。

每次迭代包含两个步骤：

E-step：求期望 $E(\gamma_{jk}|X, \theta)$ for all $j = 1,2,3,...,N$
M-step：求极大，计算新一轮迭代的模型参数

这里不具体介绍一般性的 EM 算法，（通过 Jensen 不等式得出似然函数的下界 Lower bound，通过极大化下界做到极大化似然函数，有log(E(x))>=E(log(x))），只介绍怎么在高斯混合模型里应用从来推算出模型参数。
通过 EM 迭代更新高斯混合模型参数的方法（我们有样本数据 $x_1,x_2,...,X_N$ 和一个有 $K$ 个子模型的高斯混合模型，想要推算出这个高斯混合模型的最佳参数）：

首先初始化参数
E-step：依据当前参数，计算每个数据 $j$ 来自子模型 $k$ 的可能性 $\gamma_{jk} = \frac{\alpha_k\phi(x_j|\theta_k)}{\sum_{k=1}^K\alpha_k\phi(x_j|\theta_k)}$
M-step：计算新一轮迭代的模型参数 $u_k=\frac{\sum_j^N(\gamma_{jk}x_j)}{\sum_j^N(\gamma_{jk})}, \sum_k = \frac{\sum_j^N\gamma_{jk}(x_j-u_k)(x_j-u_k)^T}{\sum_j^N(\gamma_{jk})}, \alpha_k = \frac{\sum_j^N(\gamma_{jk})}{N}$
重复计算 E-step 和 M-step 直至收敛 $||\theta_{i+1} - \theta_i|| < \epsilon$ .

至此，我们就找到了高斯混合模型的参数。需要注意的是，EM 算法具备收敛性，但并不保证找到全局最大值，有可能找到局部最大值。解决方法是初始化几次不同的参数进行迭代，取结果最好的那次。