本课程来自深度之眼,部分截图来自课程视频以及李航老师的《统计学习方法》第二版。
公式输入请参考: 在线Latex公式
前言
任务简介:理解EM算法的思想和E步、M步的求解过程。
详细说明:第9章介绍了EM算法,EM算法用于含有隐变量的概率模型的参数估计。EM算法不是一个具体的分类或回归算法,而是广泛用于含有隐变量的模型的求解问题。通过学习第1节,掌握EM算法E步和M步的求解过程;通过学习第2节,需要掌握在高斯混合模型中如何用EM算法估计参数。
学习目标:
0.导读视频。
1.通过例题9.1掌握EM算法E步和M步的求解过程。
2.了解EM算法求解如何用从最大化观测数据似然函数导出。
3.掌握高斯混合模型如何用EM算法估计参数。
例9.1 三硬币模型
假设有3枚硬币,分别记作A,B,C。这些硬币正面出现的概率分别是
π
,
p
和
q
π,p和q
π,p和q。进行如下掷硬币试验:先掷硬币A,根据其结果选出硬币B或硬币C,正面选硬币B,反面选硬币C;然后掷选出的硬币,掷硬币的结果,出现正面记作1,出现反面记作0;独立地重复n次试验(这里,n=10),观测结果如下:
1
,
1
,
0
,
1
,
0
,
0
,
1
,
0
,
1
,
1
1,1,0,1,0,0,1,0,1,1
1,1,0,1,0,0,1,0,1,1
假设只能观测到掷硬币的结果,不能观测掷硬币的过程。问如何估计三硬币正面出现的概率,即三硬币模型的参数,使得出现上面观测结果的概率最大。
虽然我们不知道模型具体的参数,但是我们可以看到一串序列,这串序列叫做:观测数据。那么不可观测数据也叫隐变量。例如第一个结果是1,也就是得到正面,但是这个正面我们不知道是因为A正面导致抛B得到的正面还是A背面导致抛C得到的正面。也就是这个例子中A的结果就是隐变量。如果我们假设A抛出后的结果为
z
z
z,最后的结果为
y
y
y,参数为
θ
\theta
θ,则根据似然估计可以写出最后丢出结果为正面
y
=
1
y=1
y=1的概率:
p
(
y
=
1
∣
θ
)
=
p
(
z
=
1
∣
θ
)
p
(
y
=
1
∣
z
=
1
,
θ
)
+
p
(
z
=
0
∣
θ
)
p
(
y
=
1
∣
z
=
0
,
θ
)
p(y=1|\theta)=p(z=1|\theta)p(y=1|z=1,\theta)+p(z=0|\theta)p(y=1|z=0,\theta)
p(y=1∣θ)=p(z=1∣θ)p(y=1∣z=1,θ)+p(z=0∣θ)p(y=1∣z=0,θ)
在本例中,做一次实验得到结果的概率可以写成:
p
(
y
∣
θ
)
=
∑
z
p
(
z
∣
θ
)
p
(
y
∣
z
,
θ
)
=
∑
z
)
p
(
y
,
z
∣
θ
)
=
π
p
y
(
1
−
p
)
1
−
y
+
(
1
−
π
)
q
y
(
1
−
q
)
1
−
y
p(y|\theta)=\sum_zp(z|\theta)p(y|z,\theta)=\sum_z)p(y,z|\theta)\\ =\pi p^y(1-p)^{1-y}+(1-\pi)q^y(1-q)^{1-y}
p(y∣θ)=z∑p(z∣θ)p(y∣z,θ)=z∑)p(y,z∣θ)=πpy(1−p)1−y+(1−π)qy(1−q)1−y
那么做N次实验的结果可以写成:
∏
i
=
1
N
p
(
y
i
∣
θ
)
=
∏
i
=
1
N
[
∑
z
p
(
y
i
,
z
∣
θ
)
]
=
∏
i
=
1
N
[
π
p
y
i
(
1
−
p
)
1
−
y
i
+
(
1
−
π
)
q
y
i
(
1
−
q
)
1
−
y
i
]
\prod_{i=1}^Np(y_i|\theta)=\prod_{i=1}^N\left[\sum_zp(y_i,z|\theta)\right]\\ =\prod_{i=1}^N\left[\pi p^{y_i}(1-p)^{1-y_i}+(1-\pi)q^{y_i}(1-q)^{1-y_i}\right]
i=1∏Np(yi∣θ)=i=1∏N[z∑p(yi,z∣θ)]=i=1∏N[πpyi(1−p)1−yi+(1−π)qyi(1−q)1−yi]
这里要最大化概率,就可以加log,连乘变连加。然后
p
(
y
i
,
z
∣
θ
)
p(y_i,z|\theta)
p(yi,z∣θ)里面有两个要估计的参数,就要用EM算法来求。
EM算法
输入:观测变量数据Y,隐变量数据Z,联合分布
P
(
Y
,
Z
∣
θ
)
P(Y,Z|\theta)
P(Y,Z∣θ),条件分布
P
(
Z
∣
Y
,
θ
)
P(Z|Y,\theta)
P(Z∣Y,θ)
输出:模型参数
θ
\theta
θ
(1)选择参数的初值
θ
(
0
)
\theta^{(0)}
θ(0),开始迭代;
(2)E步:记
θ
(
i
)
\theta^{(i)}
θ(i)为第
i
i
i次迭代参数
θ
\theta
θ的估计值,在第
i
+
1
i+1
i+1次迭代的E步,计算
Q
(
θ
,
θ
(
i
)
)
=
E
Z
[
log
P
(
Y
,
Z
∣
θ
)
]
∣
Y
,
θ
(
i
)
=
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
log
P
(
Y
,
Z
∣
θ
)
(1)
Q(\theta,\theta^{(i)})=E_Z[\log P(Y,Z|\theta)]|Y,\theta^{(i)}\\ =\sum_ZP(Z|Y,\theta^{(i)})\log P(Y,Z|\theta)\tag1
Q(θ,θ(i))=EZ[logP(Y,Z∣θ)]∣Y,θ(i)=Z∑P(Z∣Y,θ(i))logP(Y,Z∣θ)(1)
(3)M步:求使
Q
(
θ
,
θ
(
i
)
)
Q(\theta,\theta^{(i)})
Q(θ,θ(i))极大化的
θ
\theta
θ,确定第
i
+
1
i+1
i+1次迭代的参数的估计值
θ
(
i
+
1
)
\theta^{(i+1)}
θ(i+1)
θ
(
i
)
=
a
r
g
max
θ
Q
(
θ
,
θ
(
i
)
)
\theta^{(i)}=arg\underset{\theta}{\max}Q(\theta,\theta^{(i)})
θ(i)=argθmaxQ(θ,θ(i))
(4)重复第(2)步和第(3)步,直到收敛。
·式(1)的函数
Q
(
θ
,
θ
(
i
)
)
Q(\theta,\theta^{(i)})
Q(θ,θ(i))是EM算法的核心,称为
Q
Q
Q函数(
Q
Q
Q function)。
这里的核心思想是:EM算法分为E步和M步,E步是求期望,M步是求极大。
EM算法的导出
上节叙述了EM算法。为什么EM算法能近似实现对观测数据的极大似然估计呢?下面通过近似求解观测数据的对数似然函数的极大化问题来导出EM算法,由此可以清楚地看出EM算法的作用。
关于参数
θ
\theta
θ的似然函数可以写为:
L
(
θ
)
=
ln
P
(
Y
∣
θ
)
=
ln
∑
Z
P
(
Y
,
Z
∣
θ
)
=
ln
∑
Z
P
(
Z
∣
θ
)
P
(
Y
∣
Z
,
θ
)
(2)
L(\theta)=\ln P(Y|\theta)=\ln\sum_ZP(Y,Z|\theta)=\ln\sum_ZP(Z|\theta)P(Y|Z,\theta)\tag2
L(θ)=lnP(Y∣θ)=lnZ∑P(Y,Z∣θ)=lnZ∑P(Z∣θ)P(Y∣Z,θ)(2)
直接对这个似然函数求极大是很难的,里面包含有未观察到的数据
Z
Z
Z,因此换一个思路,我们希望找到的参数
θ
\theta
θ能使得似然函数
L
(
θ
)
L(\theta)
L(θ)增加,也就是
L
(
θ
)
>
L
(
θ
(
i
)
)
L(\theta)>L(\theta^{(i)})
L(θ)>L(θ(i)),这样经过若干次迭代后,似然函数就出现了极值,因此考虑
L
(
θ
)
−
L
(
θ
(
i
)
)
L(\theta)-L(\theta^{(i)})
L(θ)−L(θ(i)),看其是否有下界。
把上面似然函数的公式2带入:
L
(
θ
)
−
L
(
θ
(
i
)
)
=
ln
(
∑
Z
P
(
Z
∣
θ
)
P
(
Y
∣
Z
,
θ
)
)
−
ln
P
(
Y
∣
θ
(
i
)
)
L(\theta)-L(\theta^{(i)})=\ln\left(\sum_ZP(Z|\theta)P(Y|Z,\theta)\right)-\ln P(Y|\theta^{(i)})
L(θ)−L(θ(i))=ln(∑ZP(Z∣θ)P(Y∣Z,θ))−lnP(Y∣θ(i))
第一项乘一个除一个项:
=
ln
(
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
P
(
Z
∣
θ
)
P
(
Y
∣
Z
,
θ
)
P
(
Z
∣
Y
,
θ
(
i
)
)
)
−
ln
P
(
Y
∣
θ
(
i
)
)
=\ln\left(\sum_ZP(Z|Y,\theta^{(i)})\cfrac{P(Z|\theta)P(Y|Z,\theta)}{P(Z|Y,\theta^{(i)})}\right)-\ln P(Y|\theta^{(i)})
=ln(Z∑P(Z∣Y,θ(i))P(Z∣Y,θ(i))P(Z∣θ)P(Y∣Z,θ))−lnP(Y∣θ(i))
加的这个项不是乱弄的,整个第一项可以满足Jensen不等式(Jensen inequality):
log
∑
j
λ
j
y
j
≥
∑
j
λ
j
log
y
j
\log \sum_j\lambda_jyj\ge\sum_j\lambda_j\log y_j
logj∑λjyj≥j∑λjlogyj
这里的
λ
\lambda
λ就是
P
(
Z
∣
Y
,
θ
(
i
)
)
P(Z|Y,\theta^{(i)})
P(Z∣Y,θ(i)),而且jensen不等式满足
∑
j
λ
j
=
1
\sum_j\lambda_j=1
∑jλj=1,也就是
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
=
1
\sum_ZP(Z|Y,\theta^{(i)})=1
∑ZP(Z∣Y,θ(i))=1
因此上式大于等于下面的式子:
≥
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
ln
(
P
(
Z
∣
θ
)
P
(
Y
∣
Z
,
θ
)
P
(
Z
∣
Y
,
θ
(
i
)
)
)
−
ln
P
(
Y
∣
θ
(
i
)
)
\ge \sum_ZP(Z|Y,\theta^{(i)})\ln\left(\cfrac{P(Z|\theta)P(Y|Z,\theta)}{P(Z|Y,\theta^{(i)})}\right)-\ln P(Y|\theta^{(i)})
≥Z∑P(Z∣Y,θ(i))ln(P(Z∣Y,θ(i))P(Z∣θ)P(Y∣Z,θ))−lnP(Y∣θ(i))
后面这项乘一个值为1的项:
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
=
1
\sum_ZP(Z|Y,\theta^{(i)})=1
∑ZP(Z∣Y,θ(i))=1,上式等于:
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
ln
(
P
(
Z
∣
θ
)
P
(
Y
∣
Z
,
θ
)
P
(
Z
∣
Y
,
θ
(
i
)
)
)
−
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
ln
P
(
Y
∣
θ
(
i
)
)
=
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
[
ln
(
P
(
Z
∣
θ
)
P
(
Y
∣
Z
,
θ
)
P
(
Z
∣
Y
,
θ
(
i
)
)
)
−
ln
P
(
Y
∣
θ
(
i
)
)
]
=
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
ln
(
P
(
Z
∣
θ
)
P
(
Y
∣
Z
,
θ
)
P
(
Z
∣
Y
,
θ
(
i
)
)
P
(
Y
∣
θ
(
i
)
)
)
\sum_ZP(Z|Y,\theta^{(i)})\ln\left(\cfrac{P(Z|\theta)P(Y|Z,\theta)}{P(Z|Y,\theta^{(i)})}\right)-\sum_ZP(Z|Y,\theta^{(i)})\ln P(Y|\theta^{(i)})\\ =\sum_ZP(Z|Y,\theta^{(i)})\left [\ln\left(\cfrac{P(Z|\theta)P(Y|Z,\theta)}{P(Z|Y,\theta^{(i)})}\right)-\ln P(Y|\theta^{(i)})\right]\\ =\sum_ZP(Z|Y,\theta^{(i)})\ln\left(\cfrac{P(Z|\theta)P(Y|Z,\theta)}{P(Z|Y,\theta^{(i)})P(Y|\theta^{(i)})}\right)
Z∑P(Z∣Y,θ(i))ln(P(Z∣Y,θ(i))P(Z∣θ)P(Y∣Z,θ))−Z∑P(Z∣Y,θ(i))lnP(Y∣θ(i))=Z∑P(Z∣Y,θ(i))[ln(P(Z∣Y,θ(i))P(Z∣θ)P(Y∣Z,θ))−lnP(Y∣θ(i))]=Z∑P(Z∣Y,θ(i))ln(P(Z∣Y,θ(i))P(Y∣θ(i))P(Z∣θ)P(Y∣Z,θ))
这个是下界,我们想让下界最大,可以看到求最大的时候分母
P
(
Z
∣
Y
,
θ
(
i
)
)
P
(
Y
∣
θ
(
i
)
)
P(Z|Y,\theta^{(i)})P(Y|\theta^{(i)})
P(Z∣Y,θ(i))P(Y∣θ(i))是常数,可以忽略,分子按条件概率简写,然后可以写成:
L
(
θ
)
−
L
(
θ
(
i
)
)
≥
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
ln
P
(
Y
,
Z
∣
θ
)
L(\theta)-L(\theta^{(i)})\ge \sum_ZP(Z|Y,\theta^{(i)})\ln P(Y,Z|\theta)
L(θ)−L(θ(i))≥Z∑P(Z∣Y,θ(i))lnP(Y,Z∣θ)
令
B
(
θ
,
θ
(
i
)
)
=
L
(
θ
(
i
)
)
+
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
ln
P
(
Y
,
Z
∣
θ
)
B(\theta,\theta^{(i)})=L(\theta^{(i)})+\sum_ZP(Z|Y,\theta^{(i)})\ln P(Y,Z|\theta)
B(θ,θ(i))=L(θ(i))+Z∑P(Z∣Y,θ(i))lnP(Y,Z∣θ)
求B的最大的时候
L
(
θ
(
i
)
)
L(\theta^{(i)})
L(θ(i))是常数,可以忽略,因此:
θ
(
i
+
1
)
=
a
r
g
max
θ
B
(
θ
,
θ
(
i
)
)
=
a
r
g
max
θ
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
ln
P
(
Y
,
Z
∣
θ
)
=
a
r
g
max
θ
Q
(
θ
,
θ
(
i
)
)
\begin{aligned} \theta^{(i+1)} &=arg\underset{\theta}{\max}B(\theta,\theta^{(i)}) \\ &= arg\underset{\theta}{\max}\sum_ZP(Z|Y,\theta^{(i)})\ln P(Y,Z|\theta)\\ &= arg\underset{\theta}{\max}Q(\theta,\theta^{(i)}) \end{aligned}
θ(i+1)=argθmaxB(θ,θ(i))=argθmaxZ∑P(Z∣Y,θ(i))lnP(Y,Z∣θ)=argθmaxQ(θ,θ(i))
小结
EM算法通过迭代求解观测数据的对数似然函数的极大化,实现极大似然估计。
EM的极大似然估计中包含两个步骤:E步求期望和M步求极大。
为了每次都能够极大化L,需保证每个步骤中
L
(
θ
)
−
L
(
θ
(
i
)
)
>
0
L(\theta)-L(\theta^{(i)})>0
L(θ)−L(θ(i))>0。
通过找到
L
(
θ
)
−
L
(
θ
(
i
)
)
L(\theta)-L(\theta^{(i)})
L(θ)−L(θ(i))的下界,不断提高该下界即可。