GMM与EM共舞

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)}),它的 K 个聚类中心为 (C_k),现在要求概率值 (p( x in C_k | x))

求解的方法很简单,根据贝叶斯公式:(p(a|b)=frac{p(b|a)p(a)}{p(b)}),我们可以得出:

[p( x in C_k | x)=frac{p(C_k)p( x| xin C_k)}{p( x)} ]

因此,对于每个聚类中心 (C_k),由于分母 (p( x)) 都是相同的,我们只需要计算 (p(C_k)p( x| xin 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_{ xin S(k)}{ x} \ Sigma_k=frac{1}{N_k}sum_{ xin 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)}) 的似然函数:

[egin{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_iin C_k|x_i)} ag{3} ]

这个公式是把所有样本属于 (C_k) 的概率求平均后,作为 (C_k) 这个聚类中心(或者说这个高斯模型)的出现概率。

[mu_k=frac{sum_{i=1}^m{ xP(x_iin C_k|x_i)}}{sum_{i=1}^m{P(x_iin C_k|x_i)}} ag{4} ]

这个求均值的公式,跟单个高斯模型不同的地方在于,我们用的是加权平均。因为每个样本点都有一定的概率属于聚类中心 (C_k),所以,每个样本对 (C_k) 对应的高斯模型的均值也会产生一定的作用,只是由于 (P(x_iin C_k|x_i)) 的值不同,因此这种作用也会有显著差别。

[Sigma_k=frac{sum_{i=1}^m{P(x_iin C_k|x_i)(x_i-mu_k)(x_i-mu_k)^T}}{sum_{i=1}^m{P(x_iin C_k|x_i)}} ag{5} ]

类似地,协方差也是用加权平均求出来的。

(公式 (3) (4) (5) 其实是从极大似然函数推出来的,在周志华老师的西瓜书PRML书中都有详细推导,不过这里我只想给出感性的认识)

不过,以上公式都是基于 (P( x in C_k| x)=frac{p(C_k)p( x| xin 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 算法用一种迭代的策略来优化初始参数,并逐步收敛到最优解。关于这个算法的具体内容,我想进一步深入理解后,再用一篇文章好好写一下。

参考

原文地址:https://www.cnblogs.com/jermmyhsu/p/8251019.html