EM算法（期望最大化）——应用：GMM

GMM模型简介

GMM（Gaussian Mixture Model）也叫高斯混合模型。我们（1）可以把它看做是高斯分量的简单线性叠加，其目标是提供一种比单独的高斯分布（GSM，Gaussian Single Model）更为强大的概率模型；（2）也可以利用离散隐变量来描述GMM，并从EM算法层面给出GMM模型的一种优雅解法。
首先，给出高斯混合模型的概率公式（考虑单样本）：
$p(x|\pi,\mu,\Sigma) = \sum_k{\pi_k\mathcal{N}(x|\mu_k,\Sigma_k)}$
在第（1）种理解中，我们认为$\pi_k$是各单独高斯分布$\mathcal{N}(x|\mu_k,\Sigma_k)$所占权重。这里，我们将其作为一个隐变量$z$的先验分布的概率值。对于该离散隐变量$z$的先验分布，有如下描述：

离散变量$z$是一个K维的二值随机变量，即$z=(z_1,z_2,...,z_K)$，其中$z_k=\{0,1\}$且$\sum_kz_k=1, k=1,...,K$，$z$的概率分布可表示为$p(z_k=1|\pi) = \pi_k$其中，$\pi=(\pi_1,\pi_2,...,\pi_K) ,\sum_k\pi_k = 1$。

而原先被认为是单独高斯分布的$\mathcal{N}(x|\mu_k,\Sigma_k)$则被认为是基于先验分布下的条件高斯分布$p(x|z_k=1,\mu_k,\Sigma_k) = \mathcal{N}(x|\mu_k,\Sigma_k)$。

GMM与EM算法的联系

借鉴贝努利分布和狄里克雷分布中表示方法（1-of-K表示法），我们换一种方式对上述分布进行重新表示：

（1）先验分布：$p(z|\pi) =\prod_k{p(z_k=1)}^{z_k}=\prod_k{\pi_k}^{z_k}$
（2）条件高斯分布：$p(x|z,\mu,\Sigma)=\prod_k[p(x|z_k=1)]^{z_k}=\prod_k{[\mathcal{N}(x|\mu_k,\Sigma_k)]^{z_k}}$

所以，高斯混合模型GMM（单样本）可表示为：
$p(x|\pi,\mu,\Sigma) = \sum_k{\pi_k\mathcal{N}(x|\mu_k,\Sigma_k)} = \sum_k{(\prod_k{\pi_k}^{z_k})(\prod_k\mathcal{N}(x|\mu_k,\Sigma_k)}^{z_k}) = \sum_zp(z|\pi)p(x|z,\mu,\Sigma)= \sum_z{p(x,z|\pi,\mu,\Sigma)}$

假设，我们有样本集$\{x_n,n=1,..,N\}$，则边缘概率$p(X|\pi,\mu,\Sigma)=\prod_np(x_n|\pi,\mu,\Sigma) = \prod_n\sum_{z_n}p(z_n|\pi)p(x_n|z_n,\mu,\Sigma)$。概率模型中存在隐变量$\{z_n\},n=1,...,N$，其中$z_n=(z_{n_1},z_{n_2},...,z_{n_K})$。也就是说，对于每一个样本$x_n$，都存在一组隐变量$z_n$。在GMM中，我们的目标同样是$argmax_{(\pi,\mu,\Sigma)}p(X|\pi,\mu,\Sigma)$，但因为$lnp(X|\pi,\mu,\Sigma)=\sum_nln[\sum_{z_n}p(z_n|\pi)p(x_n|z_n,\mu,\Sigma)]$，求和符号在内部，对其求导比较困难。

GMM利用EM算法求解

这与EM算法需要解决的问题完全契合，所以，我们自然而然地想到利用EM算法求解GMM问题。
（1）E步骤：固定参数$\theta^{old} = \{\pi,\mu,\Sigma\}$，计算条件概率分布$q(Z) = p(Z|X,\theta^{old})$，从而获取联合概率分布的对数期望$Q (\theta,\theta^{old})= E[p(X,Z|\theta^{old})]$；值得注意的是：在GMM中，对于每一个样本点$x_n$都有相应的隐变量$z_n=(z_{n1},z_{n2},...,z_{nK})$，所以我们没法给出统一的$p(Z|X,\theta^{old})$表达式，但我们可以用概率集合$\{p(z_{nk}=1|x_n,\pi,\mu,\Sigma)\}$来间接表达。

（2）M步骤：固定$q(Z) = p(Z|X,\theta^{old})$，${argmax}_{\theta}Q(\theta,\theta^{old})$，其中$\theta = \{\pi,\mu,\Sigma\}$。

具体做法为：
（1）E步骤：求解$p(Z|X,\theta^{old}) = \{p(z_{nk}=1|x_n,\pi^{old},\mu^{old},\Sigma^{old})\}$：
给定样本$\{x_n,n=1,...,N\}$且样本之间相互独立，$p(z_{nk}=1|x_n,\pi,\mu,\Sigma) = \frac{p(z_{nk}=1,x_n|\pi,\mu,\Sigma) }{p(x_n|\pi,\mu,\Sigma) } = \frac{p(z_{nk}=1|\pi)p(x_n|z_{nk}=1,\mu_k,\Sigma_k)}{p(x_n|\pi,\mu,\Sigma) }= \frac{p(z_{nk}=1|\pi)p(x_n|z_{nk}=1,\mu_k,\Sigma_k)}{\sum_j{p(z_{nj}=1|\pi)p(x_n|z_{nj}=1,\mu_j,\Sigma_j)}}$
令$\gamma_{nk} = p(z_{nk}=1|x_n,\pi,\mu,\Sigma) = \frac{\pi_k\mathcal{N}(x_n|\mu_k,\Sigma_k)}{\sum_j{\pi_j\mathcal{N}(x_n|\mu_j,\Sigma_j)}}$表示后验概率集合中的一个值，且对于单个样本$x_n$有$\sum_k\gamma_{nk}=1$。

（2）M步骤：${argmax}_{\theta}Q(\theta,\theta^{old})$，其中$\theta = \{\pi,\mu,\Sigma\}$：
$Q(\theta,\theta^{old}) = \sum_Zp(Z|X,\theta^{old})ln[p(X,Z|\theta)] = \sum_{z_{nk}}p(z_{nk}=1|x_n,\pi^{old},\mu^{old},\Sigma^{old})ln[p(x_n,z_{nk}=1|\pi,\mu,\Sigma)] = \\\sum_n\sum_k\gamma_{nk}ln[\pi_k\mathcal{N}(x_n|\mu_k,\Sigma_k)] = \sum_n\sum_k\gamma_{nk}[ln\pi_k+ln\mathcal{N}(x_n|\mu_k,\Sigma_k)]$

利用$\mu_k,\Sigma_k$分别对$Q(\theta,\theta^{old})$求导：
$\frac{\partial{Q}}{\partial{\mu_k}} = \frac{\partial}{\partial{\mu_k}}({\sum_n\sum_k\gamma_{nk}[ln\pi_k+ln\mathcal{N}(x_n|\mu_k,\Sigma_k)]}) = \frac{\partial}{\partial{\mu_k}}({\sum_n\sum_k\gamma_{nk}(x_n-\mu_k)^T\Sigma_k^{-1}(x_n-\mu_k)}) = 2\sum_n\gamma_{nk}\Sigma_k^{-1}(x_n-\mu_k) = 0$
$\Rightarrow 2\sum_n\gamma_{nk}\Sigma_k\Sigma_k^{-1}(x_n-\mu_k) = \Sigma_k·0$
$\Rightarrow \mu_k = \frac{\sum_n\gamma_{nk}·x_n}{\sum_n\gamma_{nk}}$

$\frac{\partial{Q}}{\partial{\Sigma_k}} = \frac{\partial}{\partial{\Sigma_k}}({\sum_n\sum_k\gamma_{nk}[ln\pi_k+ln\mathcal{N}(x_n|\mu_k,\Sigma_k)]}) = \frac{\partial}{\partial{\Sigma_k}}({\sum_n\sum_k\gamma_{nk}(x_n-\mu_k)^T\Sigma_k^{-1}(x_n-\mu_k)}) = 0$
$\Rightarrow \Sigma_k =\frac {\sum_n\gamma_{nk}(x_n-\mu_k)(x_n-\mu_k)^T}{\sum_n\gamma_{nk}}$
对于 $\pi_k$ 则需要考虑有限制条件$\sum_k\pi_k=1$，引入拉格朗日乘数法，最大化量为$Q(\theta,\theta^{old})+\lambda(\sum_k\pi_k-1)$

$\Rightarrow \frac{\partial}{\partial{\pi_k}}[Q(\theta,\theta^{old})+\lambda(\sum_k\pi_k-1)] = \frac{\partial}{\partial{\pi_k}}[{\sum_n\sum_k\gamma_{nk}(x_n-\mu_k)^T\Sigma_k^{-1}(x_n-\mu_k)}+\lambda(\sum_k\pi_k-1)] = \sum_n\gamma_{nk}\frac1{\pi_k}+\lambda = 0$
两边同时乘以$\pi_k$，并对等式求和$\sum_k$：

$\sum_n\sum_k\gamma_{nk}+\lambda\sum_k\pi_k = 0$
$\because \sum_k\gamma_{nk}=1,\sum_k\pi_k=1$
$\therefore \lambda = -N$

$\Rightarrow \pi_k = \frac{\sum_n\gamma_{nk}}{N}$

令$N_k =\sum_n\gamma_{nk}$，则有

$\mu_k = \frac1{N_k}{\sum_n\gamma_{nk}·x_n}$
$\Sigma_k =\frac1{N_k} {\sum_n\gamma_{nk}(x_n-\mu_k)(x_n-\mu_k)^T}$
$\pi_k = \frac{N_k}{N}$

相关文章：
EM算法（期望最大化）——理论部分http://blog.csdn.net/tingyue_/article/details/70474042
EM算法（期望最大化）——从EM算法角度理解K-Means与GMM的区别http://blog.csdn.net/tingyue_/article/details/70739671