什么是EM算法
EM算法属于似然思想下,对于模型参数更新的方法。具体而言,我们以神经网络为例,我们记这个神经网络的所有参数为 θ \theta θ,可以观测到的因变量为 y y y,不可观测的因素(隐变量)记为 x x x。通过EM算法,我们可以将参数 θ \theta θ更新到使得可观测变量 y y y的log likelihood变得最大,或者说最贴合我们的数据分布。
问题描述与化简
我们要优化的函数是:
L ( θ ) = log p ( y ∣ θ ) \mathcal{L}(\theta) = \log p(y|\theta) L(θ)=logp(y∣θ)
但是这个函数我们不太好操作,我们对它做一个变换
L ( θ ) = log ∫ p ( x , y ∣ θ ) d x \begin{aligned} \mathcal{L}(\theta) & =\log \int p(x,y|\theta)dx \\ \end{aligned} L(θ)=log∫p(x,y∣θ)dx
对于任意一个 x x x的分布 q ( x ) q(x) q(x),我们将上式转化为
L ( θ ) = log ∫ q ( x ) p ( x , y ∣ θ ) q ( x ) d x \begin{aligned} \mathcal{L}(\theta)=\log \int q(x) \frac{p(x,y|\theta)}{q(x)}dx \end{aligned} L(θ)=log∫q(x)q(x)p(x,y∣θ)dx
通过Jensen不等式,我们有关系:
L ( θ ) = log ∫ q ( x ) p ( x , y ∣ θ ) q ( x ) d x ≤ ∫ q ( x ) log p ( x , y ∣ θ ) q ( x ) d x ≜ F ( q ( x ) , θ ) \begin{aligned} \mathcal{L}(\theta)=&\log \int q(x) \frac{p(x,y|\theta)}{q(x)}dx \\ \leq& \int q(x)\log \frac{p(x,y|\theta)}{q(x)}dx \triangleq \mathcal{F}(q(x),\theta) \end{aligned} L(θ)=≤log∫q(x)q(x)p(x,y∣θ)dx∫q(x)logq(x)p(x,y∣θ)dx≜F(q(x),θ)
这里我们得到的 F ( q ( x ) , θ ) \mathcal{F}(q(x), \theta) F(q(x),θ)就是我们EM算法需要去优化的函数,后面我们会说明,为什么这个优化是可以接受的。
我们不妨对 F ( q ( x ) , θ ) \mathcal{F}(q(x), \theta) F(q(x),θ)进行一点拆解,
F ( q ( x ) , θ ) = ∫ q ( x ) [ log p ( x , y ∣ θ ) − log q ( x ) ] d x = ∫ q ( x ) log p ( x , y ∣ θ ) − ∫ q ( x ) log q ( x ) d x = ∫ q ( x ) log p ( x , y ∣ θ ) − H ( q ( x ) ) \begin{aligned} \mathcal{F}(q(x), \theta) =& \int q(x)[\log p(x,y|\theta) - \log q(x)]dx \\ =&\int q(x)\log p(x,y|\theta) - \int q(x)\log q(x)dx \\ =&\int q(x)\log p(x,y|\theta) - \mathcal{H}(q(x)) \end{aligned} F(q(x),θ)===∫q(x)[logp(x,y∣θ)−logq(x)]dx∫q(x)logp(x,y∣θ)−∫q(x)logq(x)dx∫q(x)logp(x,y∣θ)−H(q(x))
这里的这个 H ( x ) \mathcal{H}(x) H(x)我们可以看出,其实是对于分布 q ( x ) q(x) q(x)计算的熵(entropy),并且不含参数 θ \theta θ项。
EM算法机制
EM算法分为两步,上面我们得到的目标优化函数 F \mathcal{F} F有两个参数进行优化,一个是隐变量 x x x的分布 q ( x ) q(x) q(x),一个是我们模型的参数 θ \theta θ,EM算法的两步分别就这两个参数进行优化。
E步
这里我们固定已经有的参数 θ \theta θ,可以是随机初始化的,只对分布 q ( x ) q(x) q(x)进行改变。
q ( k ) ( x ) = arg max q ( x ) F ( q ( x ) , θ ( k − 1 ) ) q^{(k)}(x) = \argmax_{q(x)} \mathcal{F}(q(x), \theta^{(k-1)}) q(k)(x)=q(x)argmaxF(q(x),θ(k−1))
这里我们得到的 q ( k ) ( x ) q^{(k)}(x) q(k)(x)就是这一步我们更新到的东西,把它带入 F \mathcal{F} F进行M步。
M步
首先我们固定分布 q ( x ) q(x) q(x),在这一步我们只对 θ \theta θ进行优化:
θ ( k ) = arg max θ F ( q ( k ) ( x ) , θ ) \begin{aligned} \theta^{(k)} = \argmax_{\theta} \mathcal{F}(q^{(k)}(x),\theta) \end{aligned} θ(k)=θargmaxF(q(k)(x),θ)
这里我们得到的 θ ( k ) \theta^{(k)} θ(k)就是这一步我们更新的参数,带入 F \mathcal{F} F后进入E步。
总而言之这俩步骤就形成了一个内循环,我们可以设置迭代的上限或者拟合值的阈值作为条件来跳出优化循环。
为什么EM算法能用
我们由上面的推导可以发现,我们实际操作的函数并不是一开始的log似然函数 L ( θ ) \mathcal{L}(\theta) L(θ),而是它变换后的一个函数,我们将在这里证明为什么对于这个 F \mathcal{F} F进行的优化能够最终对原函数进行改善。
首先我们将进行说明, F \mathcal{F} F实质上是 L \mathcal{L} L的一个下界。
为了证明这个东西,最直观的想法就是把两个函数减一减,看看最后结果的正负号。依据这个思路,我们进行下面的操作,每一步的依据在下面阐述:
L ( θ ) − F ( q ( x ) , θ ) = log p ( y ∣ θ ) − ∫ q ( x ) log p ( x , y ∣ θ ) q ( x ) d x = log p ( y ∣ θ ) − ∫ q ( x ) log p ( x ∣ y , θ ) p ( y ∣ θ ) q ( x ) d x = log p ( y ∣ θ ) − log p ( y ∣ θ ) ∫ q ( x ) d x − ∫ q ( x ) log p ( x ∣ y , θ ) q ( x ) d x = − ∫ q ( x ) log p ( x ∣ y , θ ) q ( x ) d x ≜ K L ( q ( x ) , p ( x ∣ y , θ ) ) \begin{aligned} \mathcal{L}(\theta) - \mathcal{F}(q(x),\theta) =& \log p(y|\theta) - \int q(x) \log \frac{p(x,y|\theta)}{q(x)}dx \\ =&\log p(y|\theta) - \int q(x) \log \frac{p(x|y,\theta)p(y|\theta)}{q(x)}dx \\ =&\log p(y|\theta) - \log p(y|\theta) \int q(x) dx - \int q(x) \log \frac{p(x|y,\theta)}{q(x)}dx \\ =& - \int q(x) \log \frac{p(x|y,\theta)}{q(x)}dx \triangleq \bm{KL}(q(x), p(x|y,\theta)) \end{aligned} L(θ)−F(q(x),θ)====logp(y∣θ)−∫q(x)logq(x)p(x,y∣θ)dxlogp(y∣θ)−∫q(x)logq(x)p(x∣y,θ)p(y∣θ)dxlogp(y∣θ)−logp(y∣θ)∫q(x)dx−∫q(x)logq(x)p(x∣y,θ)dx−∫q(x)logq(x)p(x∣y,θ)dx≜KL(q(x),p(x∣y,θ))
第一行只是把两个东西进行代入,第二行利用了重期望公式,第三行利用log函数把内部乘在一起的分开成加在一起的,第四行是因为第三行第二项的积分是单纯密度函数的积分,最后积出来就是1。
最后我们得到的东西称为 Kullback-Liebler divergence,我们后面讨论这个东西的正负性。
这里我们插入一个lemma,叫做 Gibb’s inequality,是一个不等式,它说明的是:
对于变量x而言,他的任意两个分布 p ( x ) , q ( x ) p(x), q(x) p(x),q(x),均满足:
− ∑ x ∈ X p ( x ) log p ( x ) ≤ − ∑ x ∈ X p ( x ) log ( q ( x ) ) -\sum_{x \in \mathcal{X}} p(x) \log p(x) \leq -\sum_{x \in \mathcal{X}}p(x)\log(q(x)) −x∈X∑p(x)logp(x)≤−x∈X∑p(x)log(q(x))
对于连续型随机变量而言,是这个形式:
− ∫ p ( x ) log p ( x ) d x ≤ − ∫ p ( x ) log q ( x ) d x -\int p(x) \log p(x) dx \leq -\int p(x) \log q(x) dx −∫p(x)logp(x)dx≤−∫p(x)logq(x)dx
并且等号成立当且仅当 p ( x ) = q ( x ) p(x) = q(x) p(x)=q(x)
换言之,在这种求和方式下,一个分布的熵是最小的,不等式左边正是一个分布的熵。
那么利用这个引理,我们可以就可以轻松说明 KL divergence 的正负性的。
K L ( q ( x ) , p ( x ∣ y , θ ) ) = − ∫ q ( x ) log p ( x ∣ y , θ ) q ( x ) d x = − ∫ q ( x ) log p ( x ∣ y , θ ) d x − ( − ∫ q ( x ) log q ( x ) d x ) \begin{aligned} KL(q(x), p(x|y,\theta)) =& - \int q(x) \log \frac{p(x|y,\theta)}{q(x)}dx \\ =& -\int q(x) \log p(x|y,\theta)dx - (-\int q(x)\log q(x)dx) \end{aligned} KL(q(x),p(x∣y,θ))==−∫q(x)logq(x)p(x∣y,θ)dx−∫q(x)logp(x∣y,θ)dx−(−∫q(x)logq(x)dx)
利用 Gibb 不等式,显然得到,KL divergence 是一个非负的函数,并且等号成立当且仅当 q ( x ) = p ( x ∣ y , θ ) q(x) = p(x|y,\theta) q(x)=p(x∣y,θ),即和那个两个变量之间是独立的。
这里,我们完成了整体的第一部分证明:即 F \mathcal{F} F是 L \mathcal{L} L的一个下界。
下面我们利用这个不等关系说明EM算法的合理性,理由同样在下面可以找到:
L ( θ ( k − 1 ) ) = F ( q ( k ) , θ ( k − 1 ) ) ≤ F ( q ( k ) ( x ) , θ ( k ) ) ≤ L ( θ ( k ) ) \mathcal{L}(\theta^{(k-1)}) = \mathcal{F}(q^{(k)},\theta^{(k-1)}) \leq \mathcal{F}(q^{(k)}(x), \theta^{(k)}) \leq \mathcal{L}(\theta^{(k)}) L(θ(k−1))=F(q(k),θ(k−1))≤F(q(k)(x),θ(k))≤L(θ(k))
第一个等号是E步保证的,通过argmax这个分布q,得到函数 F \mathcal{F} F的上界,即令 q ( k ) ( x ) = p ( k − 1 ) ( x ∣ y , θ ) q^{(k)}(x) = p^{(k-1)}(x|y,\theta) q(k)(x)=p(k−1)(x∣y,θ)。
第二个不等号是M步保证的,通过更新参数 θ \theta θ来argmax当前的分布 q ( x ) q(x) q(x)。
第三个不等号是我们最前面提到的Jensen不等式。
由这个不等式,我们得到了重要结论:EM算法永远不会变差。
参考
https://www.cs.cmu.edu/~tom/10-702/Zoubin-702.pdf