EM算法原理分析

EM算法主要用于含有隐藏变量的参数估计问题。
在将EM算法之前,先讲一下Jensen不等式。
定理:假设f是一个凸函数,X是随机变量,即:

E[f(X)]f(EX) E [ f ( X ) ] ≥ f ( E X )

此外,如果f是严格凸的,当且仅当 X=E[X]=E[f(X)]=f(EX) X = E [ X ] = 常 数 ( 不 再 是 随 机 变 量 ) 时 E [ f ( X ) ] = f ( E X ) .
不理解的可以看下面的图:
这里写图片描述
是不是一目了然?简单解释一下:假设X是一个随机变量,有0.5的概率落在a点,有0.5的概率落在b点,因此X的期望 E[X] E [ X ] 便落在a,b 的中点处。根据f是凸函数,我们可以在图上画出 f(a),f(b),f(E[X]) f ( a ) , f ( b ) , f ( E [ X ] ) 的位置,而 E[f(X)] E [ f ( X ) ] 则落在 f(a),f(b) f ( a ) , f ( b ) 的中点处。
由上图可知,因为f是凸函数,所以有 E[f(X)]f(EX) E [ f ( X ) ] ≥ f ( E X ) 。同理,如果f是凹函数,则有 E[f(X)]f(EX) E [ f ( X ) ] ≤ f ( E X )
EM算法
假设我们有m个独立样本(独立性假设) {x(1),...,x(m)} { x ( 1 ) , . . . , x ( m ) } ,给定以下似然函数:
l(θ)=i=1mlog p(x;θ)=i=1mlogzp(x,z;θ) l ( θ ) = ∑ i = 1 m l o g   p ( x ; θ ) = ∑ i = 1 m l o g ∑ z p ( x , z ; θ )

我们希望求出模型 p(x,z) p ( x , z ) 的参数 θ θ . 然而,由于存在隐藏变量 z z ,θ的求解是很困难的,如果能够提前得到 z z ,那么最大似然估计将变得简单起来。(请记住这一点,因为后面的EM算法的E步其实就相当于给z做了一个先验假设,然后再做优化)
对于每个样本i,假设Qi是关于z的分布( zQi(z)=1,Qi(z)0 ∑ z Q i ( z ) = 1 , Q i ( z ) ≥ 0 ),因此可得到下列不等式:
l(θ)=i=1mlog p(x(i);θ)=i=1mlogz(i)p(x(i),z(i);θ)=i=1mlogz(i)Qi(z(i))p(x(i),z(i);θ)Qi(z(i))i=1mz(i)Qi(z(i))logp(x(i),z(i);θ)Qi(z(i))(1)(2)(3)(4) (1) l ( θ ) = ∑ i = 1 m l o g   p ( x ( i ) ; θ ) (2) = ∑ i = 1 m l o g ∑ z ( i ) p ( x ( i ) , z ( i ) ; θ ) (3) = ∑ i = 1 m l o g ∑ z ( i ) Q i ( z ( i ) ) p ( x ( i ) , z ( i ) ; θ ) Q i ( z ( i ) ) (4) ≥ ∑ i = 1 m ∑ z ( i ) Q i ( z ( i ) ) l o g p ( x ( i ) , z ( i ) ; θ ) Q i ( z ( i ) )

最后一步怎么得来的呢?其实就是用到了Jensen不等式。特别的, f(x)=log x f ( x ) = l o g   x 是一个凹函数,因为 f′′(x)=1x2<0 f ″ ( x ) = − 1 x 2 < 0 .因此有 E[f(x)]f(E(x)) E [ f ( x ) ] ≤ f ( E ( x ) ) ,其中自变量x为 p(x(i),z(i);θ)Qi(z(i)) p ( x ( i ) , z ( i ) ; θ ) Q i ( z ( i ) ) ,代入得:
f(Ez(i)Qi[p(x(i),z(i);θ)Qi(z(i))])Ez(i)Qi[f(p(x(i),z(i);θ)Qi(z(i)))] f ( E z ( i ) ∼ Q i [ p ( x ( i ) , z ( i ) ; θ ) Q i ( z ( i ) ) ] ) ≥ E z ( i ) ∼ Q i [ f ( p ( x ( i ) , z ( i ) ; θ ) Q i ( z ( i ) ) ) ]
.
综上便可得到上文所述不等式。
那么,不等式什么时候取等号呢?其实上文的定理已经提到了,当自变量为常数时等号成立,对应到我们得不等式中,即:
p(x(i),z(i);θ)Qi(z(i))=c p ( x ( i ) , z ( i ) ; θ ) Q i ( z ( i ) ) = c
.
事实上,我们知道 zQi(z)=1 ∑ z Q i ( z ) = 1 ,因此我们可以得到下面得推导:
Qi(z(i))=p(x(i),z(i);θ)zp(x(i),z;θ)=p(x(i),z(i);θ)p(x(i);θ)=p(z(i)|x(i);θ)(5)(6)(7) (5) Q i ( z ( i ) ) = p ( x ( i ) , z ( i ) ; θ ) ∑ z p ( x ( i ) , z ; θ ) (6) = p ( x ( i ) , z ( i ) ; θ ) p ( x ( i ) ; θ ) (7) = p ( z ( i ) | x ( i ) ; θ )

也就是说,我们可以简单设置 Qi Q i 为在参数 θ θ 下给定 x(i) x ( i ) 时,关于 z(i) z ( i ) 的后验分布。
因此,我们可以得到EM算法的迭代过程如下:
循环以下两步直到收敛{
(E-step)对于每个样本i,
Qi(z(i)):=p(z(i)|x(i);θ). Q i ( z ( i ) ) := p ( z ( i ) | x ( i ) ; θ ) .

(M-step)
θ:=argmaxθi=1mz(i)Qi(z(i))logp(x(i),z(i);θ)Qi(z(i)). θ := a r g m a x θ ∑ i = 1 m ∑ z ( i ) Q i ( z ( i ) ) l o g p ( x ( i ) , z ( i ) ; θ ) Q i ( z ( i ) ) .

}
那么我们怎么知道EM算法是否收敛呢?我们假设 θ(t)θ(t+1) θ ( t ) 和 θ ( t + 1 ) 为迭代过程中的参数,那么我们只要证明 l(θ(t))l(θ(t+1)) l ( θ ( t ) ) ≤ l ( θ ( t + 1 ) ) ,那么就可以得到EM算法是在不断优化,直至收敛。顺着这个思想,我们假设 Qi(z(i)):=p(z(i)|x(i);θ) Q i ( z ( i ) ) := p ( z ( i ) | x ( i ) ; θ ) ,此时
l(θ(t))=i=1mz(i)Qi(z(i))logp(x(i),z(i);θ(t))Qi(z(i)) l ( θ ( t ) ) = ∑ i = 1 m ∑ z ( i ) Q i ( z ( i ) ) l o g p ( x ( i ) , z ( i ) ; θ ( t ) ) Q i ( z ( i ) )

参数 θ(t+1) θ ( t + 1 ) 通过最大化等式右边的式子获得,因此:
θ(t+1)i=1mz(i)Qi(z(i))logp(x(i),z(i);θ(t+1))Qi(z(i))i=1mz(i)Qi(z(i))logp(x(i),z(i);θ(t))Qi(z(i))=l(θ(t))(8)(9)(10) (8) θ ( t + 1 ) ≥ ∑ i = 1 m ∑ z ( i ) Q i ( z ( i ) ) l o g p ( x ( i ) , z ( i ) ; θ ( t + 1 ) ) Q i ( z ( i ) ) (9) ≥ ∑ i = 1 m ∑ z ( i ) Q i ( z ( i ) ) l o g p ( x ( i ) , z ( i ) ; θ ( t ) ) Q i ( z ( i ) ) (10) = l ( θ ( t ) )

θ θ θ(t+1) θ ( t + 1 ) 时, p(x(i),z(i);θ(t+1))Q(t)i(z(i)) p ( x ( i ) , z ( i ) ; θ ( t + 1 ) ) Q i ( t ) ( z ( i ) ) 不一定为常数了,所以等号不一定成立,因此上述第一个式子为大于等于。
至于第二个不等式,由EM算法的M步可知, θ(t+1) θ ( t + 1 ) 是通过最大化上一步的函数值得到的,即:

θ(t+1):=argmaxθi=1mz(i)Qi(z(i))logp(x(i),z(i);θ(t))Qi(z(i)). θ ( t + 1 ) := a r g m a x θ ∑ i = 1 m ∑ z ( i ) Q i ( z ( i ) ) l o g p ( x ( i ) , z ( i ) ; θ ( t ) ) Q i ( z ( i ) ) .

再把 θ(t+1) θ ( t + 1 ) 迭代回去,得到的函数值肯定会大于等于上一步的函数值,因此第二个不等式成立。
综上,通过EM算法,我们总可以得到 l(θ(t+1))l(θ(t)) l ( θ ( t + 1 ) ) ≥ l ( θ ( t ) ) ,从而不断优化,直到收敛,收敛条件是函数值增长小于等于阈值(阈值自己设定)时,停止迭代。

二 高斯混合模型(Gaussian Misture Model, GMM)
EM算法的一个重要应用就是高斯混合模型的参数估计。
高斯混合模型(Gaussian Misture Model, GMM)是指具有如下形式的概率分布模型:

p(y|θ)=j=1kϕjp(y|θj) p ( y | θ ) = ∑ j = 1 k ϕ j p ( y | θ j )

其中, ϕj ϕ j 是系数, ϕj0,kj=1ϕj=1 ϕ j ≥ 0 , ∑ j = 1 k ϕ j = 1 ; p(y|θj) p ( y | θ j ) 是高斯分布密度, θj=(μj,σ2j)=((μj,Σj) θ j = ( μ j , σ j 2 ) = ( ( μ j , Σ j ) ,
p(y|θj)=1(2π)12σjexp((yμj)22σ2j) p ( y | θ j ) = 1 ( 2 π ) 1 2 σ j e x p ( − ( y − μ j ) 2 2 σ j 2 )

称为第j个分模型。
一般混合模型可以由任意概率分布密度代替上式中的高斯分布密度,我们这里只介绍最常用的高斯混合模型。

E-step:计算

w(i)j=Qi(z(i)=j)=P(z(i)=j|x(i);ϕ,μ,Σ). w j ( i ) = Q i ( z ( i ) = j ) = P ( z ( i ) = j | x ( i ) ; ϕ , μ , Σ ) .

w(i)j w j ( i ) 是针对第i个样本,在参数为 ϕ,μ,Σ ϕ , μ , Σ 已知样本特征 x(i) x ( i ) 的情况下,属于第j个分模型的概率。
M-step:最大化以下式子优化参数 ϕ,μ,Σ ϕ , μ , Σ

L=i=1mz(i)Qi(z(i))logp(x(i),z(i);ϕ,μ,Σ)Qi(z(i))=i=1mz(i)Qi(z(i)=j)logp(x(i)|z(i)=j;μ,Σ)p(z(i)=j;ϕ)Qi(z(i))=j=i=1mz(i)w(i)jlog1(2π)12|Σj|12exp(12(x(i)μj)TΣ1j(x(i)μj))ϕjw(i)j(18)(19)(20) (18) L = ∑ i = 1 m ∑ z ( i ) Q i ( z ( i ) ) l o g p ( x ( i ) , z ( i ) ; ϕ , μ , Σ ) Q i ( z ( i ) ) (19) = ∑ i = 1 m ∑ z ( i ) Q i ( z ( i ) = j ) l o g p ( x ( i ) | z ( i ) = j ; μ , Σ ) p ( z ( i ) = j ; ϕ ) Q i ( z ( i ) ) = j (20) = ∑ i = 1 m ∑ z ( i ) w j ( i ) l o g 1 ( 2 π ) 1 2 | Σ j | 1 2 e x p ( − 1 2 ( x ( i ) − μ j ) T Σ j − 1 ( x ( i ) − μ j ) ) ⋅ ϕ j w j ( i )

首先我们关于 μl μ l 最大化以上式子。将L对 μl μ l 求导,得到:

Lμl=μli=1mj=1kw(i)jlog1(2π)12|Σj|12exp(12(x(i)μj)TΣ1j(x(i)μj))ϕjw(i)j=μli=1mj=1kw(i)j12(x(i)μj)TΣ1j(x(i)μj)=12i=1mw(i)lμl2μTlΣ1lx(i)μTlΣ1lμl=i=1mw(i)l(Σ1lx(i)Σ1lμl)(21)(22)(23)(24) (21) ∂ L ∂ μ l = ∇ μ l ∑ i = 1 m ∑ j = 1 k w j ( i ) l o g 1 ( 2 π ) 1 2 | Σ j | 1 2 e x p ( − 1 2 ( x ( i ) − μ j ) T Σ j − 1 ( x ( i ) − μ j ) ) ⋅ ϕ j w j ( i ) (22) = ∇ μ l ∑ i = 1 m ∑ j = 1 k w j ( i ) 1 2 ( x ( i ) − μ j ) T Σ j − 1 ( x ( i ) − μ j ) (23) = 1 2 ∑ i = 1 m w l ( i ) ∇ μ l 2 μ l T Σ l − 1 x ( i ) − μ l T Σ l − 1 μ l (24) = ∑ i = 1 m w l ( i ) ( Σ l − 1 x ( i ) − Σ l − 1 μ l )

令导数等于零,可得到 μl μ l 的更新规则如下:

μl:=mi=1w(i)lx(i)mi=1w(i)l. μ l := ∑ i = 1 m w l ( i ) x ( i ) ∑ i = 1 m w l ( i ) .

至于 Σ Σ 的更新跟 μl μ l 类似,不再赘述。下面讲一下 ϕ ϕ 的更新。
通过观察式子,我们可以把无关变量去掉,得到:
L(ϕ)=i=1mj=1kw(i)jlogϕj. L ( ϕ ) = ∑ i = 1 m ∑ j = 1 k w j ( i ) l o g ϕ j .

另一方面,因为 ϕj=p(z(i)=j;ϕ) ϕ j = p ( z ( i ) = j ; ϕ ) ,所以有约束条件 kj=1ϕj=1 ∑ j = 1 k ϕ j = 1 .因此,我们使用拉格朗日乘子 β β 将有约束问题转换成无约束问题,如下:
L(ϕ)=i=1mj=1kw(i)jlogϕj+β(j=1kϕj1) L ( ϕ ) = ∑ i = 1 m ∑ j = 1 k w j ( i ) l o g ϕ j + β ( ∑ j = 1 k ϕ j − 1 )

值得注意的是,这里并没有把约束条件 ϕj>0 ϕ j > 0 加上,这是为什么呢?别急,下文会提到。
对以上式子求导,得到:
L(ϕ)ϕj=i=1mw(i)jϕj+β ∂ L ( ϕ ) ∂ ϕ j = ∑ i = 1 m w j ( i ) ϕ j + β

令导数等于零,可得到 ϕj ϕ j 的更新规则如下:
ϕj:=mi=1w(i)jβ. ϕ j := ∑ i = 1 m w j ( i ) − β .

使用约束条件 kj=1ϕj=1 ∑ j = 1 k ϕ j = 1 ,我们可以得到 β=mi=1kj=1w(i)j=mi=11=m(使w(i)j=Qi(z(i)=j),kj=1w(i)j=1) − β = ∑ i = 1 m ∑ j = 1 k w j ( i ) = ∑ i = 1 m 1 = m ( 使 用 条 件 w j ( i ) = Q i ( z ( i ) = j ) , 从 而 ∑ j = 1 k w j ( i ) = 1 ) ,因此,我们可以进一步化简得到:
ϕj:=1mi=1mw(i)j. ϕ j := 1 m ∑ i = 1 m w j ( i ) .

我们可以看到, ϕj ϕ j 恒大于零,默认满足约束条件 ϕj>0 ϕ j > 0

再简单说明一下我理解的EM算法与Kmeans算法的联系与区别:
联系:Kmeans算法可以看作EM算法的一个特例,Kmeans中的簇即为EM算法中的隐藏变量;
区别:Kmeans中每一个数据点都只属于一个簇中,属于硬分隔;
而EM算法使用后验概率的方法,相当于一个数据点分到每一个簇都有一个概率,概率和为1.

参考:吴恩达CS229 Lecture notes “The EM algorithm”
《统计学习方法》(李航著)

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值