EM算法是一种迭代算法,用于含有隐变量(hidden variable)的概率模型参数的极大似然估计,或极大后验概率估计。
EM算法的每次迭代由两步组成:
- E步,求期望(expectation);
- M步,求极大(maximization)。
所以这一算法称为期望极大算法(expectation maximization algorithm),简称EM算法。
9.1 EM算法的引入
概率模型有时既含有观测变量(observable variable),又含有隐变量或潜在变量(latent variable)。如果概率模型的变量都是观测变量,那么给定数据,可以直接用极大似然估计发,或贝叶斯估计法估计模型参数。但是,当模型含有隐变量时,就不能简单地使用这些估计方法。EM算法就是含有隐变量的概率模型参数的极大似然估计,或极大后验概率估计。
9.1.1 EM算法
EM算法与初值的选择有关,选择不同的初值可能得到不同的参数估计值。
一般地,用$Y$表示观测随机变量的数据,$Z$表示隐随机变量的数据。$Y$和$Z$连在一起称为完全数据(complete-data),观测数据$Y$又称为不完全数据(incomplete-data)。假设给定观测数据$Y$,其概率分布是$P(Y|\theta)$,其中$\theta$是需要估计的模型参数,那么不完全数据$Y$的似然函数是$P(Y|\theta)$,对数似然函数$L(\theta) = \log P(Y|\theta)$;假设$Y$和$Z$的联合概率分布是$P(Y,Z|\theta)$,那么完全数据的对数似然函数是$\log P(Y,Z|\theta)$。
EM算法通过迭代求$L(\theta) = \log P(Y|\theta)$的极大似然估计。每次迭代都包含两步:E步,求期望;M步,求极大化。
算法 9.1(EM算法)
输入:观测变量数据$Y$,隐变量数据$Z$,联合分布$P(Y,Z|\theta)$,条件分布$P(Z|Y,\theta)$;
输出:模型参数$\theta$。
(1)选择参数的初值$\theta^{(0)}$,开始迭代;
(2)E步:记$\theta^{(i)}$为第$i$次迭代参数$\theta$的估计值,在第$i+1$次迭代的E步,计算
这里,$P(Z|Y,\theta)$是在给定观测数据$Y$和当前的参数估计$\theta^{(i)}$下隐变量数据$Z$的条件概率分布;
(3)M步:求使$Q(\theta,\theta^{(i)})$极大化的$\theta$,确定第$i+1$次迭代的参数的估计值$\theta^{(i+1)}$
(4)重复第(2)步和第(3)步,直到收敛。
式$(9.9)$的函数$Q(\theta,\theta^{(i)})$是EM算法的核心,称为$Q$函数(Q function)。
定义 9.1($Q$函数) 完全数据的对数似然函数$\log P(Y,Z|\theta)$关于在给定观测数据$Y$和当前参数$\theta^{(i)}$下对未观测数据$Z$的条件概率分布$P(Z|Y,\theta^{(i)})$的期望称为$Q$函数,即
关于EM算法的几点说明:
步骤(1) 参数的初值可以任意选择,但需要注意EM算法对初值是敏感的。
步骤(2) E步求$Q(\theta,\theta^{(i)})$。$Q$函数式中$Z$是未观测数据,$Y$是观测数据。注意,$Q(\theta,\theta^{(i)})$的第一个变元表示要极大化的参数,第二个变元表示参数的当前估计值。每次迭代实际在求$Q$函数及其极大。
步骤(3) M步求$Q(\theta,\theta^{(i)})$的极大化,得到$\theta^{(i+1)}$,完成一次迭代$\theta^{(i)} \longrightarrow \theta^{(i+1)}$。
步骤(4)给出停止迭代的条件,一般是对较小的正数$\epsilon_1,\epsilon_2$,若满足
则停止迭代。
9.1.2 EM算法的导出
对一个含有隐变量的概率模型,目标是极大化观测数据(不完全数据)$Y$关于参数$\theta$的对数似然函数,即极大化
这一极大化的主要困难是式$(9.12)$中有未观测数据并有包含和(或积分)的对数。
事实上,EM算法是通过迭代逐步近似极大化$L(\theta)$的。假设在第$i$次迭代后$\theta$的估计值是$\theta^{(i)}$。我们希望新估计值$\theta$能使$L(\theta)$增加,即$L(\theta) > L(\theta^{(i)})$,并逐步达到极大值。为此,考虑两者的差:
利用$Jensen$不等式(Jensen inequality)得到其下界:
令
则
即函数$B(\theta,\theta^{(i)})$是$L(\theta)$的一个下界,而且由式$(9.13)$可知,
因此,任何可以使$B(\theta,\theta^{(i)})$增大的$\theta$,也可以使$L(\theta)$增大。为了使$L(\theta)$有尽可能大的增长,选择$\theta^{(i+1)}$使$B(\theta,\theta^{(i)})$达到极大,即
现在求$\theta^{(i+1)}$的表达式。省去对$\theta$的极大化而言是常数的项,由式$(9.16)$、式$(9.13)$即式$(9.10)$,有
式$(9.17)$等价于EM算法的一次迭代,即求$Q$函数及其极大化。EM算法是通过不断求解下界的极大化逼近求解对数似然函数极大化的算法。
EM算法不能保证找到全局最优值。
9.1.3 EM算法在无监督学习中的应用
EM算法可以用于生成模型的无监督学习。生成模型由联合概率分布$P(X,Y)$表示,可以认为无监督学习训练数据是联合概率分布产生的数据。$X$为观测数据,$Y$为未观测数据。
9.2 EM算法的收敛性
定理 9.1 设$P(Y|\theta)$为观察数据的似然函数,$\theta^{(i)}(i=1,2,…)$为EM算法得到的参数估计序列,$P(Y|\theta^{(i)})(i=1,2,…)$为对应的似然函数序列,则$P(Y|\theta^{(i)})$是单调递增的,即
定理 9.2 设$L(\theta) = \log P(Y|\theta)$为观测数据的对数似然函数,$\theta^{(i)}(i=1,2,…)$为EM算法得到的参数估计序列,$L(\theta^{(i)})(i=1,2,…)$为对应的对数似然函数序列。
(1)如果$P(Y|\theta)$有上界,则$L(\theta^{(i)}) = \log P(Y|\theta^{(i)})$收敛到某一值$L^*$;
(2)在函数$Q(\theta,\theta^{‘})$与$L(\theta)$满足一定条件下,由EM算法得到的参数估计序列$\theta^{(i)}$的收敛值$\theta^{*}$是$L(\theta)$的稳定点。
EM算法的收敛性包含关于对数似然函数序列$L(\theta^{(i)})$的收敛性和关于参数估计序列$\theta^{(i)}$的收敛性两层意思,前者并不蕴含后者。此外,定理只能保证参数估计序列收敛到对数似然函数序列的稳定点,不能保证收敛到极大值点。所以在应用中,初值的选择变得非常重要,常用的办法是选取几个不同的初值进行迭代,然后对得到的各个估计值加以比较,从中选择最好的。
9.3 EM算法在高斯混合模型学习中的应用
EM算法是学习高斯混合模型(Gaussian mixture model)的有效方法。
9.3.1 高斯混合模型
定义 9.2(高斯混合模型) 高斯混合模型是指具有如下形式的概率分布模型:
其中,$\alpha_k$是系数,$\alpha_k \geq 0$,$\sum\limits_{k=1}^{K}\alpha_{k} = 1$;$\phi(y|\theta_k)$是高斯分布密度,$\theta_k = (\mu_k,\sigma_k^2)$,
称为第$k$个分模型。
一般混合模型可以由任意概率分布密度代替式$(9.25)$中的高斯分布密度。
9.3.2 高斯混合模型参数估计的EM算法
算法 9.2(高斯混合模型参数估计的EM算法)
输入:观测数据$y_1,y_2,…,y_N$,高斯混合模型;
输出:高斯混合模型参数。
(1)取参数的初始值开始迭代;
(2)E步:依据当前模型参数,计算分模型$k$对观测数据$y_j$的响应度
(3)M步:计算新一轮迭代的模型参数
(4)重复第(2)步和第(3)步,直到收敛。
9.4 EM算法的推广
EM算法还可以解释为$F$函数(F function)的极大-极大算法(maximization-maximization algorithm),基于这个解释有若干变形与推广,如广义期望极大(generalized expectation maximization,GEM)算法。
9.4.1 $F$函数的极大-极大算法
定义 9.3($F$函数) 假设隐变量数据$Z$的概率分布为$\tilde{P}(Z)$,定义分布$\tilde{P}$与参数$\theta$的函数$F(\tilde{P},\theta)$如下:
称为$F$函数。式中$H(\tilde{P}) = - E_{\tilde{P}}\log \tilde{P}(Z)$是分布$\tilde{P}(Z)$的熵。
引理 9.1 对于固定的$\theta$,存在唯一的分布$\tilde{P}_\theta$极大化$F(\tilde{P},\theta)$,这时$\tilde{P}_\theta$由下式给出:
并且$\tilde{P}_\theta$随$\theta$连续变化。
引理 9.2 若$\tilde{P}_\theta(Z) = P(Z|Y,\theta)$,则
定理 9.3 设$L(\theta) = \log P(Y|\theta)$为观测数据的对数似然函数,$\theta^{(i)},i=1,2,…$,为EM算法得到的参数估计序列,函数$F(\tilde{P},\theta)$由式$(9.33)$定义。如果$F(\tilde{P},\theta)$在$\tilde{P}^$和$\theta^$有局部极大值,那么$L(\theta)$也在$\theta^$有局部极大值。类似地,如果$F(\tilde{P},\theta)$在$\tilde{P}^$和$\theta^$达到全局最大值,那么$L(\theta)$也在$\theta^$达到全局最大值。
定理 9.4 EM算法的一次迭代可由$F$函数的极大-极大算法实现。
设$\theta^{(i)}$为第$i$次迭代参数$\theta$的估计,$\tilde{P}^{(i)}$为第$i$次迭代函数$\tilde{P}$的估计。在第$i+1$次迭代的两步为:
(1)对固定的$\theta^{(i)}$,求$\tilde{P}^{(i+1)}$使$F(\tilde{P},\theta^{(i)})$极大化;
(2)对固定的$\tilde{P}^{(i+1)}$,求$\theta^{(i+1)}$使$F(\tilde{P}^{(i+1)},\theta)$极大化。
9.4.2 GEM算法
算法 9.3(GEM算法1)
输入:观测数据,$F$函数;
输出:模型参数。
(1)初始化参数$\theta^{(0)}$,开始迭代。
(2)第$i+1$次迭代,第一步:记$\theta^{(i)}$为参数$\theta$的估计值,$\tilde{P}^{(i)}$为函数$\tilde{P}$的估计,求$\tilde{P}^{(i+1)}$使$\tilde{P}$极大化$F(\tilde{P},\theta^{(i)})$;
(3)第二步:求$\theta^{(i+1)}$使$F(\tilde{P}^{(i+1)},\theta)$极大化;
(4)重复(2)和(3),直到收敛。
在GEM算法1中,有时求$Q(\theta,\theta^{(i)})$的极大化是很困难的。下面介绍GEM算法2和GEM算法3.
算法 9.4(GEM算法2)
输入:观测数据,$Q$函数;
输出:模型参数。
(1)初始化参数$\theta^{(0)}$,开始迭代。
(2)第$i+1$次迭代,第一步:记$\theta^{(i)}$为参数$\theta$的估计值,计算
(3)第二步:求$\theta^{(i+1)}$使
(4)重复(2)和(3),直到收敛。
当参数的维数为时,可采用一种特殊的GEM算法,它将EM算法的M步分解为次条件极大化,每次只改变参数向量的一个分量,其余分量不改变。
算法 9.5(GEM算法3)
输入:观测数据,$Q$函数;
输出:模型参数。
(1)初始化参数$\theta^{(0)} = \{\theta^{(0)}_1,\theta^{(0)}_2,…,\theta^{(0)}_d\}$,开始迭代。
(2)第$i+1$次迭代,第一步:记$\theta^{(i)} = \{\theta^{(i)}_1,\theta^{(i)}_2,…,\theta^{(i)}_d\}$为参数$\theta = \{\theta_1,\theta_2,…,\theta_d\}$的估计值,计算
(3)第二步:进行$d$次条件极大化:
首先,在$\theta^{(i)}_2,…,\theta^{(i)}_d$保持不变的条件下求使$Q(\theta,\theta^{(i)}) $达到极大的$\theta_1^{(i+1)}$;然后,在$\theta_1 = \theta_1^{(i+1)},\theta_j = \theta_j^{(i)},j=3,4,..,d$的条件下求使$Q(\theta,\theta^{(i)})$达到极大的$\theta_2^{(i+1)}$;如此继续,经过$d$次条件极大化,得到$\theta^{(i+1)} = \{\theta^{(i+1)}_1,\theta^{(i+1)}_2,…,\theta^{(i+1)}_d\}$使得
(4)重复(2)和(3),直到收敛。