EM算法与高斯混合模型学习总结

       

          我们现在考虑一类具有隐变量的统计推断问题，用数学的语言就是说:

         1) 我们的总体样本空间$\mathcal{S}\in\mathbb{R}^{M}\times\mathbb{R}^{L}$, 其分布的概率密度函数是$P(x,y\mid \theta)$，其中$x\in\mathbb{R}^{M},y\in\mathbb{R}^{L}$, $\theta$是未知的参数。

         2) 我们在抽样的时候，由于条件限制，对每一个样本$(x,y)$，$x\in\mathbb{R}^{M},y\in\mathbb{R}^{L}$，后L个变量无法观测到，也就是向量$y$无法获悉，只能观测到$x$。我们称前M个变量为可观测变量，x为可观测向量，后L个变量为不可观测变量，y为不可观测向量。

 

     问题1：

         我们有数据集合$D=\lbrace x_{1},...,x_{N} \rbrace$, 其中$x_{i}\in\mathbb{R}^{M}$为第i个样本可观测向量，样本之间相互独立，试推断可能的参数$\theta$。 

 

EM原算法:

     

     为了方便起见，我们记大写字母$X=(x_{1},...,x_{N})$表示所有的样本的可观测向量的总体，而$Y=(y_{1},...,y_{N})$表示所有的样本的不可观测向量的总体, 我们自然地将概率密度表示为紧凑的形式： 

                   $$P(X,Y\mid\theta)\triangleq P((x_{1},y_{1}),...,(x_{N},y_{N})\mid\theta)=\prod_{i=1}^{N}P((x_{1},y_{1}),...,(x_{N},y_{N})\mid\theta)$$

        这个时候自然的，我们的似然函数就是：

                     \begin{equation}P(X\mid\theta)=\int P(X,Y\mid\theta)dY\end{equation}

 当$Y$各个分量连续，或者是：

                    \begin{equation}P(X\mid\theta)=\sum_{Y}P(X,Y\mid\theta)\end{equation}

当$Y$各个分量离散，或者是还有混合的情形，Y部分分量离散部分连续。为了记号方便起见我们假定Y各个分量均离散，其实各种情形下计算过程基本无差异。

 

        现在我们的目标是使得似然函数$P(X\mid \theta)$尽量大。注意到一般来说由于表达式(2)中求和号的存在，直接求似然函数的极值比较困难，现在我们想办法设计一种迭代算法使得极大似然函数尽量大。假设我们第$i$步迭代已经得到参数$\theta^{(i)}$, 现在希望迭代更新使得似然函数尽量大。注意到：

 

\begin{split}&\frac{P(X\mid\theta)}{P(X\mid\theta^{(i)})}\newline =&\sum_{Y}\frac{P(X,Y\mid\theta)}{P(X\mid\theta^{(i)})}\newline =&\sum_{Y}\frac{P(X,Y\mid \theta)P(Y\mid X,\theta^{(i)})}{P(X,Y\mid \theta^{(i)})}\end{split}

        所以两边取$log$并且运用Jensen不等式我们得到：

\begin{split}\log(\frac{P(X\mid\theta)}{P(X\mid\theta^{(i)})})=&\log(\sum_{Y}\frac{P(X,Y\mid \theta)P(Y\mid X,\theta^{(i)})}{P(X,Y\mid \theta^{(i)})})\newline \geq &\sum_{Y}P(Y\mid X,\theta^{(i)})\log\frac{P(X,Y\mid \theta)}{P(X,Y\mid \theta^{(i)})}\newline =&\sum_{Y}[\log P(X,Y\mid \theta)]P(Y\mid X,\theta^{(i)})-\sum_{Y}[\log P(X,Y\mid \theta^{(i)})]P(Y\mid X,\theta^{(i)})\newline =&E_{Y}(\log P(X,Y\mid \theta)\big | \ X,\theta^{(i)})-E_{Y}(\log P(X,Y\mid \theta^{(i)})\big | \ X,\theta^{(i)})\end{split} 

 

我们令辅助函数:

\begin{equation}Q(\theta,\theta^{(i)})\triangleq E_{Y}(\log P(X,Y\mid \theta)\big | \ X,\theta^{(i)})\end{equation}

 

所以有：

\begin{equation}\log P(X\mid \theta)-\log P(X\mid \theta^{(i)})\geq Q(\theta,\theta^{(i)})-Q(\theta^{(i)},\theta^{(i)})\end{equation} 

 

为了使得 $\log p(X\mid \theta)$显著增加，我们只需要选取一个:

\begin{equation}\theta^{(i+1)}\triangleq\mathop{\text{argmax}}\limits_{\theta} Q(\theta,\theta^{(i)})\end{equation}

       从上述推导我们看到，我们首先得显式的算出函数$Q$，它由一个求期望(Expectation)的计算给出，然后下一步是求Q的极大值使得Q极大化(Maximization)，所以我们的算法也会分为两步，一步是求期望(E步)，再求极大值（M步），这种方法于是被称为EM算法。

现在总结EM算法如下：

                输入：观测变量$X$, 隐含变量$Y$, 联合分布$p(X,Y\mid\theta)$, $\epsilon$

                输出：模型参数$\theta$

                Step1. 选择模型的初始参数$\theta^{(0)}$

                Step2. E步：记$\theta^{(i)}$是第$i$步迭代得到的参数，则在第$i+1$步迭代中，计算函数：

                            \begin{equation}Q(\theta,\theta^{(i)})\triangleq E_{Y}(\log p(X,Y\mid \theta)\big | \ X,\theta^{(i)})\end{equation}

                Step3  M步：更新 $\theta^{(i+1)}\triangleq\mathop{\text{argmax}}\limits_{\theta} Q(\theta,\theta^{(i)})$, 如果 $Q(\theta^{(i+1)},\theta^{(i)})-Q(\theta^{(i)},\theta^{(i)})\geq \epsilon$ 则回到Step2继续，否则停止，输出$\theta^{(i)}$.

 

 

高斯混合模型

     

     所谓的高斯混合模型是对如下的隐变量模型应用EM算法统计推断：

 

         问题2：

          我们研究一类问题1的特例，即样本总体$\mathcal{S}=\mathbb{R}^{M}\times\lbrace 1,2,...,K\rbrace$, 分布密度函数$p(x,k\mid\theta)=\pi_{k}P(x\mid\Sigma_{k},\mu_{k})$, 其中参数$\theta=(\Sigma_{1},...,\Sigma_{K},\mu_{1},...,\mu_{K},\pi_{1},...,\pi_{i})$ 满足对任意的k=1,...,K：

          1)$\Sigma_{k}$为对称正定$M\times M$矩阵；

          2)$\mu_{k}\in\mathbb{R}^{M}$；

          3)$\pi_{k}\geq 0$,$\sum_{k=1}^{K}\pi_{i}=1$,

且我们用$P(x\mid\Sigma,\mu)$表示$M$维正态分布$\mathcal{N}(\mu_{k},\Sigma_{k})$的概率密度函数。现在如果我们有数据集合$D=\lbrace x_{1},...,x_{N} \rbrace$, 其中$x_{i}\in\mathbb{R}^{M}$为第i个样本可观测向量，样本之间相互独立，试推断可能的参数$\theta$。

                                                 

       E步：

         我们对上面的问题用EM算法求解，先开始E步，求函数Q。为了方便起见，对于任意$i=1,...,N$,我们定义函数： $l_{i}: SGL(M)\times \mathbb{R}^{M}\longrightarrow \mathbb{R}:$

$$l_{i}(\Sigma,\mu)=\log[P(x_{i}\mid \Sigma, \mu)],$$ 

对任意的$\Sigma\in SGL(M), \mu\in\mathbb{R}^{M}$, 这里$SGL(M)$表示全体对称可逆的$M\times M$矩阵的集合。

         这时：

\begin{split}Q(\theta,\theta^{(s)})=& E_{Y}(\log P(x_{1},...,x_{N},y_{1},...,y_{N}\mid \theta)\big | x_{1},...,x_{N},\theta^{(s)})\newline =&\sum_{i=1}^{N}E_{Y}(\log P(x_{i},y_{i}\mid \theta)\big | X, \theta^{(s)})\newline =&\sum_{i=1}^{N}E_{y_{i}}(\log P(x_{i},y_{i}\mid \theta)\big | x_{i}, \theta^{(s)})\newline =&\sum_{i=1}^{N}\sum_{k=1}^{K}[\log P(x_{i},k\mid \theta)]P(k\mid x_{i}, \theta^{(s)})\newline =&\sum_{i=1}^{N}\sum_{k=1}^{K}P_{ki}^{(s)}\log [P(x_{i}\mid \Sigma_{k},\mu_{k})\cdot\pi_{k}]\newline =&\sum_{i=1}^{N}\sum_{k=1}^{K}P_{ki}^{(s)}l_{i}(\Sigma_{k},\mu_{k})+\sum_{k=1}^{K}(\sum_{i=1}^{N}P_{ki}^{(s)})\log(\pi_{k})\end{split}

        

       M步：

         现在开始求Q的极大值，先得计算一下偏导数。下面我们先重点计算一下函数$l_{i}$的微分. 这时候 $P(x_{i}\mid \Sigma, \mu)=(2\pi)^{-\frac{M}{2}}\det(\Sigma)^{-\frac{1}{2}}\exp\lbrace -\frac{1}{2}(x_{i}-\mu)^{T}\cdot\Sigma^{-1}\cdot(x_{i}-\mu)\rbrace$, $$l_{i}(\Sigma,\mu)=-\frac{1}{2}\log(\det(\Sigma))-\frac{1}{2}(x_{i}-\mu)^{T}\Sigma^{-1}(x_{i}-\mu)$$.

我们注意到行列式求导公式：$d(\log(\det\Sigma))=\text{tr}(\Sigma^{-1}d\Sigma)$所以容易得到：

 

\begin{equation}dl_{i}=-\frac{1}{2}tr\lbrace[\Sigma^{-1}-\Sigma^{-1}(x_{i}-\mu)\cdot(x_{i}-\mu)^{T}\Sigma^{-1}]d\Sigma\rbrace+(x_{i}-\mu)^{T}\Sigma^{-1}d\mu\end{equation}

注意到我们要求极值问题：

                                       \begin{split}\text{max }&\ Q(\theta,\theta^{(s)}) \newline \text{s.t. }&\ \pi_{k}\in [0,1],&\ \sum_{k=1}^{K}\pi_{k}=1\end{split} 

只需要求某个模型参数$\theta$以及常数$\lambda\in\mathbb{R}$ 满足：

                                       \begin{equation} dQ=\lambda \sum_{k=1}^{K}d\pi_{k}\end{equation}

 结合(5)我们很容易得到：

$$dQ=-\frac{1}{2}\text{tr}(A_{k}\cdot d\Sigma_{k})+[\sum_{i=1}^{N}P^{(s)}_{ki}x_{i}-(\sum_{i=1}^{N}P_{ki})\mu_{k}]d \mu_{k}+\sum_{k=1}^{K}(\sum_{i=1}^{N}P_{ki}^{(s)})\frac{d \pi_{k}}{\pi_{k}}$$其中矩阵$A_{k}=\Sigma_{k}^{-1}(\sum_{i=1}^{N}P_{ki}^{(s)})-\Sigma_{k}^{-1}[\sum_{i=1}^{N}P_{ki}^{(s)}(x_{i}-\mu_{k})\cdot (x_{i}-\mu_{k})^{T}]\Sigma_{k}^{-1}$ 

再结合（6）我们有：$A_{k}=0$, $\sum_{i=1}^{N}P^{(s)}_{ki}x_{i}-(\sum_{i=1}^{N}P_{ki})\mu_{k}=0$, $\sum_{i=1}^{N}P_{ki}^{(s)}=\lambda\pi_{k}$, 

所以立即得到对任意$k=1,...,K$, 当$\theta$各个参数取如下值的时候$Q(\theta,\theta_{k})$将取极大值：

\begin{equation}\Sigma_{k}=\frac{\sum_{i=1}^{N}P_{ki}^{(s)}(x_{i}-\mu_{k})\cdot(x_{i}-\mu_{k})^{T}}{\sum_{i=1}^{N}P_{ki}^{(s)}},\end{equation}

\begin{equation}\mu_{k}=\frac{\sum_{i=1}^{N}P_{ki}^{(s)}x_{i}}{\sum_{i=1}^{N}P_{ki}^{(s)}},\end{equation}

\begin{equation}\pi_{k}=\frac{\sum_{i=1}^{N}P^{(s)}_{ki}}{N},\end{equation}

 

此时$M$步计算完毕。

 

此时我们可以总结一下Gauss混合模型的算法如下：        

         

---

 

                    输入：数据集$D= \lbrace x_{1},...,x_{N}\rbrace $, 混合类别个数$K$，某停止条件

                    输出：高斯混合模型的参数$\theta=(\Sigma_{1},...,\Sigma_{K},\mu_{1},...,\mu_{K},\pi_{1},...,\pi_{K})$

                    Step1. 初始化某参数$\theta^{(0)}$

                    Step2. 对于$s=1,2,...$执行：

                                     1. E-Step: 计算： \begin{equation}P^{(s)}_{ki}=\frac{P(x_{i}\mid \Sigma_{k}^{(s)},\mu_{k}^{(s)})\pi_{k}^{(s)}}{\sum_{j=1}^{K} P(x_{i}\mid \Sigma_{j}^{(s)},\mu_{j}^{(s)})\pi_{j}^{(s)}}\end{equation}

                                     2. M-Step:  对任意的$k=1,...,K$计算：\begin{equation}\Sigma_{k}^{(s+1)}=\frac{\sum_{i=1}^{N}P_{ki}^{(s)}(x_{i}-\mu_{k})\cdot(x_{i}-\mu_{k})^{T}}{\sum_{i=1}^{N}P_{ki}^{(s)}},\end{equation}

\begin{equation}\mu_{k}^{(s+1)}=\frac{\sum_{i=1}^{N}P_{ki}^{(s)}x_{i}}{\sum_{i=1}^{N}P_{ki}^{(s)}},\end{equation}

                        \begin{equation}\pi_{k}^{(s+1)}=\frac{\sum_{i=1}^{N}P^{(s)}_{ki}}{N},\end{equation}

                                        令$\theta^{(s+1)}=(\Sigma_{1}^{(s+1)},...,\Sigma_{K}^{(s+1)},\mu_{1}^{(s+1)},...,\mu_{K}^{(s+1)},\pi_{1}^{(s+1)},...,\pi_{K}^{(s+1)})$, 

                                        如果满足停止条件，则输出$\theta^{(s)}$

 

参考文献：

【1】李航：《统计学习方法》，北京，清华大学出版社，2012

【2】周志华：《机器学习》，北京，清华大学出版社，2016