E M \mathrm{EM} EM算法
E
M
\mathrm{EM}
EM算法是一种迭代算法，用于含有隐变量
(hidden
variable)
\text{(hidden\;variable)}
(hiddenvariable)的概率模型参数的极大似然估计,其每次迭代分为两步：
E
E
E步，求期望；
M
M
M步，求极大；所以这一算法也称为期望极大算法
(expectation
maximization
algorithm)
\text{(expectation\;maximization\;algorithm)}
(expectationmaximizationalgorithm)。
概率模型有时既含有观测变量，优含有隐变量或潜在变量。如果概率模型的变量都是观测变量，那么给定数据，可以用极大似然估计法或贝叶斯法估计模型参数。当含有隐变量时，就不能简单地使用这些估计方法。
E
M
\mathrm{EM}
EM算法就是含有隐变量的概率模型参数的极大似然估计法，或极大后验概率估计法。
例
为了便于理解首先引入一个例子
（三硬币模型） 假设有
3
3
3枚硬币，分别记作
A
,
B
,
C
A,B,C
A,B,C，这些硬币正面出现的概率分别是
π
,
p
,
q
\pi,p,q
π,p,q.进行如下掷硬币试验：先掷硬币
A
A
A，根据其结果选出硬币
B
B
B或
C
C
C，正面选硬币
B
B
B，反面选硬币
C
C
C，掷硬币的结果，出现正面记作
1
1
1，出现反面记作
0
0
0；独立地重复次试验(这里，
n
=
10
n=10
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假设只能观测到掷硬币的结果，不能观测掷硬币的过程，问如何估计三硬币正面出现的概率，即三硬币模型的参数。
约定：条件概率用
(
∣
)
(|)
(∣)表示，参数应该用
(
;
)
(;)
(;)表示，联合概率用
(
,
)
(,)
(,)表示
解：三硬币模型可以写作
P
(
y
;
θ
)
=
∑
z
P
(
y
,
z
;
θ
)
=
∑
z
P
(
z
;
θ
)
P
(
y
∣
z
;
θ
)
=
π
p
y
(
1
−
p
)
1
−
y
+
(
1
−
π
)
q
y
(
1
−
q
)
1
−
y
\begin{aligned} P(y;\theta) &=\sum_{z} P(y, z; \theta)=\sum_{z} P(z ; \theta) P(y | z; \theta) \\ &=\pi p^{y}(1-p)^{1-y}+(1-\pi) q^{y}(1-q)^{1-y} \end{aligned}
P(y;θ)=z∑P(y,z;θ)=z∑P(z;θ)P(y∣z;θ)=πpy(1−p)1−y+(1−π)qy(1−q)1−y这里，随机变量
y
y
y是观测变量，表示一次试验观测的结果是
1
1
1或
0
0
0；随机变量
z
z
z是隐变量，表示为观测到的掷硬币
A
A
A的结果；
θ
=
(
π
,
p
,
q
)
\theta=(\pi, p, q)
θ=(π,p,q)是模型参数。随机变量
y
y
y的数据可以观测，随机变量
z
z
z的数据不可观测。
将观测数据表示为
Y
=
(
Y
1
,
Y
2
,
⋯
,
Y
n
)
T
Y=\left(Y_{1}, Y_{2}, \cdots, Y_{n}\right)^{\mathrm{T}}
Y=(Y1,Y2,⋯,Yn)T，未观测数据表示为
Z
=
(
Z
1
,
Z
2
,
⋯
,
Z
n
)
T
Z=\left(Z_{1}, Z_{2}, \cdots, Z_{n}\right)^{\mathrm{T}}
Z=(Z1,Z2,⋯,Zn)T，则观测数据的似然函数为
P
(
Y
;
θ
)
=
∑
Z
P
(
Z
;
θ
)
P
(
Y
∣
Z
;
θ
)
P(Y;\theta)=\sum_{Z} P(Z;\theta)P(Y|Z;\theta)
P(Y;θ)=Z∑P(Z;θ)P(Y∣Z;θ)即
P
(
Y
;
θ
)
=
∏
j
=
1
n
[
π
p
y
j
(
1
−
p
)
1
−
y
j
+
(
1
−
π
)
q
y
j
(
1
−
q
)
1
−
y
j
]
P(Y;\theta)=\prod_{j=1}^{n}\left[\pi p^{y_{j}}(1-p)^{1-y_{j}}+(1-\pi) q^{y_{j}}(1-q)^{1-y_{j}}\right]
P(Y;θ)=j=1∏n[πpyj(1−p)1−yj+(1−π)qyj(1−q)1−yj]求模型参数
θ
=
(
π
,
p
,
q
)
\theta=(\pi, p, q)
θ=(π,p,q)的极大似然估计，即
θ
^
=
arg
max
θ
log
P
(
Y
;
θ
)
=
arg
max
θ
∑
j
=
1
n
log
(
π
p
y
j
(
1
−
p
)
1
−
y
j
+
(
1
−
π
)
q
y
j
(
1
−
q
)
1
−
y
j
)
\hat{{\theta}}=\arg \max _{{\theta}} \log P(Y; {\theta})=\arg \max _{{\theta}}\sum_{j=1}^{n}\log(\pi p^{y_{j}}(1-p)^{1-y_{j}}+(1-\pi) q^{y_{j}}(1-q)^{1-y_{j}})
θ^=argθmaxlogP(Y;θ)=argθmaxj=1∑nlog(πpyj(1−p)1−yj+(1−π)qyj(1−q)1−yj)这个问题没有解析解，只有通过迭代的方法求解，
E
M
\mathrm{EM}
EM算法就是可以用于求解这个问题的一种迭代算法。
E
M
\mathrm{EM}
EM算法首先选取参数的初值，记作
θ
(
0
)
=
(
π
(
0
)
,
p
(
0
)
,
q
(
0
)
)
\theta^{(0)}=\left(\pi^{(0)}, p^{(0)}, q^{(0)}\right)
θ(0)=(π(0),p(0),q(0))，然后通过下面步骤迭代计算参数的估计值，直至收敛为止，第
i
i
i次迭代参数的估计值为
θ
(
i
)
=
(
π
(
i
)
,
p
(
i
)
,
q
(
i
)
)
\theta^{(i)}=\left(\pi^{(i)}, p^{(i)}, q^{(i)}\right)
θ(i)=(π(i),p(i),q(i))，
E
M
\mathrm{EM}
EM算法的第
i
+
1
i+1
i+1次迭代如下。
E
E
E步：计算在模型参数
π
(
i
)
,
p
(
i
)
,
q
(
i
)
\pi^{(i)}, p^{(i)}, q^{(i)}
π(i),p(i),q(i)下观测数据
y
i
y_i
yi来自掷硬币
B
B
B的概率
μ
(
i
+
1
)
=
π
(
i
)
(
p
(
i
)
)
y
j
(
1
−
p
(
i
)
)
1
−
y
j
π
(
i
)
(
p
(
i
)
)
y
j
(
1
−
p
(
i
)
)
1
−
y
j
+
(
1
−
π
(
i
)
)
(
q
(
i
)
)
y
j
(
1
−
q
(
i
)
)
1
−
y
j
\mu^{(i+1)}=\frac{\pi^{(i)}\left(p^{(i)}\right)^{y_{j}}\left(1-p^{(i)}\right)^{1-y_{j}}}{\pi^{(i)}\left(p^{(i)}\right)^{y_{j}}\left(1-p^{(i)}\right)^{1-y_{j}}+\left(1-\pi^{(i)}\right)\left(q^{(i)}\right)^{y_{j}}\left(1-q^{(i)}\right)^{1-y_{j}}}
μ(i+1)=π(i)(p(i))yj(1−p(i))1−yj+(1−π(i))(q(i))yj(1−q(i))1−yjπ(i)(p(i))yj(1−p(i))1−yj
M
M
M步：计算模型参数的新估计值
π
(
i
+
1
)
=
1
n
∑
j
=
1
n
μ
j
(
i
+
1
)
\pi^{(i+1)}=\frac{1}{n} \sum_{j=1}^{n} \mu_{j}^{(i+1)}
π(i+1)=n1j=1∑nμj(i+1)
p
(
i
+
1
)
=
∑
j
=
1
n
μ
j
(
i
+
1
)
y
j
∑
j
=
1
n
μ
j
(
i
+
1
)
p^{(i+1)}=\frac{\sum_{j=1}^{n} \mu_{j}^{(i+1)} y_{j}}{\sum_{j=1}^{n} \mu_{j}^{(i+1)}}
p(i+1)=∑j=1nμj(i+1)∑j=1nμj(i+1)yj
q
(
i
+
1
)
=
∑
j
=
1
n
(
1
−
μ
j
(
i
+
1
)
)
y
j
∑
j
=
1
n
(
1
−
μ
j
(
i
+
1
)
)
q^{(i+1)}=\frac{\sum_{j=1}^{n}\left(1-\mu_{j}^{(i+1)}\right) y_{j}}{\sum_{j=1}^{n}\left(1-\mu_{j}^{(i+1)}\right)}
q(i+1)=∑j=1n(1−μj(i+1))∑j=1n(1−μj(i+1))yj数字计算，假设模型参数的初值取为
π
(
0
)
=
0.5
,
p
(
0
)
=
0.5
,
q
(
0
)
=
0.5
\pi^{(0)}=0.5, \quad p^{(0)}=0.5, \quad q^{(0)}=0.5
π(0)=0.5,p(0)=0.5,q(0)=0.5由上面
E
E
E步公式，对
y
j
=
1
y_{j}=1
yj=1与
y
j
=
0
y_{j}=0
yj=0均有
μ
j
(
1
)
=
0.5
\mu_{j}^{(1)}=0.5
μj(1)=0.5
由
M
M
M步公式有
π
(
1
)
=
0.5
,
p
(
1
)
=
0.6
,
q
(
1
)
=
0.6
\pi^{(1)}=0.5, \quad p^{(1)}=0.6, \quad q^{(1)}=0.6
π(1)=0.5,p(1)=0.6,q(1)=0.6
再由
E
E
E步公式得到
μ
j
(
2
)
=
0.5
,
j
=
1
,
2
,
⋯
,
10
\mu_{j}^{(2)}=0.5, \quad j=1,2, \cdots, 10
μj(2)=0.5,j=1,2,⋯,10
继续迭代得，
π
(
2
)
=
0.5
,
p
(
2
)
=
0.6
,
q
(
2
)
=
0.6
\pi^{(2)}=0.5, \quad p^{(2)}=0.6, \quad q^{(2)}=0.6
π(2)=0.5,p(2)=0.6,q(2)=0.6
于是得到模型参数的极大似然估计
π
^
=
0.5
,
p
^
=
0.6
,
q
^
=
0.6
\hat{\pi}=0.5, \quad \hat{p}=0.6, \quad \hat{q}=0.6
π^=0.5,p^=0.6,q^=0.6
如果初值
π
(
0
)
=
0.4
,
p
(
0
)
=
0.6
,
q
(
0
)
=
0.7
\pi^{(0)}=0.4, \quad p^{(0)}=0.6, \quad q^{(0)}=0.7
π(0)=0.4,p(0)=0.6,q(0)=0.7，那么得到的模型参数的极大似然估计
π
^
=
0.4064
,
p
^
=
0.5368
,
q
^
=
0.6432
\hat{\pi}=0.4064, \quad \hat{p}=0.5368, \quad \hat{q}=0.6432
π^=0.4064,p^=0.5368,q^=0.6432，
E
M
\mathrm{EM}
EM算法与初值的选择有关，选择不同的初值可能得到不同的参数估计值。
一般地，用
Y
Y
Y表示观测随机变量的数据，
Z
Z
Z表示隐随机变量的数据。
Y
Y
Y和
Z
Z
Z连在一起称为完全数据(complete-data),观测数据
Y
Y
Y又称为不完全数据(incomplete-data).假设给定观测数据
Y
Y
Y，其概率分布是
P
(
Y
;
θ
)
P(Y;\theta)
P(Y;θ),其中
θ
\theta
θ是需要估计的模型参数，那么不完全数据
Y
Y
Y的似然函数是
P
(
Y
;
θ
)
P(Y;\theta)
P(Y;θ)，对数似然函数是
L
(
θ
)
=
log
P
(
Y
;
θ
)
L(\theta)=\log P(Y;\theta)
L(θ)=logP(Y;θ);假设
Y
Y
Y和
Z
Z
Z的联合概率分布是
P
(
Y
,
Z
;
θ
)
P(Y, Z;\theta)
P(Y,Z;θ),则完全数据的对数似然函数是
log
P
(
Y
,
Z
;
θ
)
\log P(Y, Z;\theta)
logP(Y,Z;θ)
E
M
\mathrm{EM}
EM算法通过迭代求
L
(
θ
)
=
log
P
(
Y
;
θ
)
L(\theta)=\log P(Y;\theta)
L(θ)=logP(Y;θ)的极大似然估计。每次迭代包含两步：
E
E
E步，求期望；
M
M
M步，求极大化。
E
M
\mathrm{EM}
EM算法
输入：观测变量数据
Y
Y
Y，隐变量数据
Z
Z
Z，联合分布
P
(
Y
,
Z
;
θ
)
P(Y, Z;\theta)
P(Y,Z;θ)，条件分布
P
(
Z
∣
Y
;
θ
)
P(Z | Y;\theta)
P(Z∣Y;θ)
输出：模型参数
θ
\theta
θ
(
1
)
(1)
(1)选择参数的初值，开始迭代
(
2
)
(2)
(2)
E
E
E步:记
θ
(
i
)
\theta^{(i)}
θ(i)为第
i
i
i次迭代参数
θ
\theta
θ的估计值，在第
i
+
1
i+1
i+1次迭代的
E
E
E步，计算
Q
(
θ
,
θ
(
i
)
)
=
E
Z
[
log
P
(
Y
,
Z
;
θ
)
∣
Y
;
θ
(
i
)
]
=
∑
Z
P
(
Z
∣
Y
;
θ
(
i
)
)
log
P
(
Y
,
Z
;
θ
)
\begin{aligned} Q\left(\theta, \theta^{(i)}\right) &=E_{Z}\left[\log P(Y, Z;\theta) | Y;\theta^{(i)}\right] \\ &=\sum_{Z} P\left(Z | Y;\theta^{(i)}\right) \log P(Y, Z;\theta) \end{aligned}
Q(θ,θ(i))=EZ[logP(Y,Z;θ)∣Y;θ(i)]=Z∑P(Z∣Y;θ(i))logP(Y,Z;θ)这里，
P
(
Z
∣
Y
;
θ
(
i
)
)
P\left(Z | Y;\theta^{(i)}\right)
P(Z∣Y;θ(i))是在给定观测数据
Y
Y
Y和当前的参数估计
θ
(
i
)
\theta^{(i)}
θ(i)下隐变量数据
Z
Z
Z的条件概率分布。
(
3
)
(3)
(3)
M
M
M步：求使
Q
(
θ
,
θ
(
i
)
)
Q\left(\theta, \theta^{(i)}\right)
Q(θ,θ(i))极大化的
θ
\theta
θ，确定第
i
+
1
i+1
i+1次迭代的参数的估计值
θ
(
i
+
1
)
\theta^{(i+1)}
θ(i+1)
θ
(
i
+
1
)
=
arg
max
θ
Q
(
θ
,
θ
(
i
)
)
\theta^{(i+1)}=\arg \max _{\theta} Q\left(\theta, \theta^{(i)}\right)
θ(i+1)=argθmaxQ(θ,θ(i))
(
4
)
(4)
(4)重复步骤
(
2
)
(2)
(2)和
(
3
)
(3)
(3)直到收敛
上式
Q
(
θ
,
θ
(
i
)
)
Q\left(\theta, \theta^{(i)}\right)
Q(θ,θ(i))是
E
M
\mathrm{EM}
EM算法的核心，称为
Q
Q
Q函数
定义(
Q
Q
Q函数) 完全数据的对数似然函数
log
P
(
Y
,
Z
;
θ
)
\log P(Y, Z;\theta)
logP(Y,Z;θ)关于在给定观测数据
Y
Y
Y和当前的参数估计
θ
(
i
)
\theta^{(i)}
θ(i)下对未观测数据的条件概率分布
P
(
Z
∣
Y
;
θ
(
i
)
)
P\left(Z | Y;\theta^{(i)}\right)
P(Z∣Y;θ(i))的期望称为
Q
Q
Q函数，即
Q
(
θ
,
θ
(
i
)
)
=
E
z
[
log
P
(
Y
,
Z
;
θ
)
∣
Y
;
θ
(
i
)
]
(
1
)
Q\left(\theta, \theta^{(i)}\right)=E_{z}\left[\log P(Y, Z;\theta) | Y;\theta^{(i)}\right] \qquad(1)
Q(θ,θ(i))=Ez[logP(Y,Z;θ)∣Y;θ(i)](1)
E
M
\mathrm{EM}
EM算法的几点说明：
步骤
(
1
)
(1)
(1) 参数的初值可以任意选择，但需注意
E
M
\mathrm{EM}
EM算法对初值是敏感的步骤
步骤
(
2
)
(2)
(2)
E
E
E步求
Q
(
θ
,
θ
(
i
)
)
Q\left(\theta, \theta^{(i)}\right)
Q(θ,θ(i))，
Q
Q
Q函数中
Z
Z
Z是未观测数据，
Y
Y
Y是观测数据，
Q
(
θ
,
θ
(
i
)
)
Q\left(\theta, \theta^{(i)}\right)
Q(θ,θ(i))的第
1
1
1个变元表示要极大化的参数，第
2
2
2个变元表示参数的当前估计值，每次迭代实际在求
Q
Q
Q函数及其极大
步骤
(
3
)
(3)
(3)
M
M
M步求
Q
(
θ
,
θ
(
i
)
)
Q\left(\theta, \theta^{(i)}\right)
Q(θ,θ(i))的极大,得到
θ
(
i
+
1
)
\theta^{(i+1)}
θ(i+1)，完成一次迭代
θ
(
i
)
→
θ
(
i
+
1
)
\theta^{(i)} \rightarrow \theta^{(i+1)}
θ(i)→θ(i+1)，后面将证明每次迭代使似然函数增大或达到局部极值
步骤
(
4
)
(4)
(4) 给出停止迭代的条件，一般是对较小的正数
ε
1
,
ε
2
\varepsilon_{1}, \varepsilon_{2}
ε1,ε2,若满足
∥
θ
(
i
+
1
)
−
θ
(
i
)
∥
<
ε
1
\left\|\theta^{(i+1)}-\theta^{(i)}\right\|<\varepsilon_{1}
∥∥θ(i+1)−θ(i)∥∥<ε1或
∥
Q
(
θ
(
i
+
1
)
,
θ
(
i
)
)
−
Q
(
θ
(
i
)
,
θ
(
i
)
)
∥
<
ε
2
\left\|Q\left(\theta^{(i+1)}, \theta^{(i)}\right)-Q\left(\theta^{(i)}, \theta^{(i)}\right)\right\|<\varepsilon_{2}
∥∥Q(θ(i+1),θ(i))−Q(θ(i),θ(i))∥∥<ε2
则停止迭代
E
M
\mathrm{EM}
EM算法的导出
面对一个含有隐变量的概率模型，目标是极大化观测数据(不完全数据)
Y
Y
Y关于参数
θ
\theta
θ的对数似然函数，即极大化
L
(
θ
)
=
log
P
(
Y
;
θ
)
=
log
∑
Z
P
(
Y
,
Z
;
θ
)
=
log
(
∑
Z
P
(
Y
∣
Z
;
θ
)
P
(
Z
;
θ
)
)
(
2
)
\begin{aligned} L(\theta) &=\log P(Y ; \theta)=\log \sum_{Z} P(Y, Z ; \theta) \\ &=\log \left(\sum_{Z} P(Y | Z;\theta) P(Z ; \theta)\right) \end{aligned}\qquad(2)
L(θ)=logP(Y;θ)=logZ∑P(Y,Z;θ)=log(Z∑P(Y∣Z;θ)P(Z;θ))(2)极大化的主要困难是上式中有未观测数据并有包含和(或积分)的对数。
E
M
\mathrm{EM}
EM算法是通过迭代逐步近似极大化
L
(
θ
)
L(\theta)
L(θ)的,假设在第
i
i
i次迭代后
θ
\theta
θ的估计值是
θ
(
i
)
\theta^{(i)}
θ(i)，我们希望新估计值
θ
\theta
θ能使
L
(
θ
)
L(\theta)
L(θ)增加，即
L
(
θ
)
>
L
(
θ
(
i
)
)
L(\theta)>L\left(\theta^{(i)}\right)
L(θ)>L(θ(i))，并逐步达到极大值，为此，考虑两者的差：
L
(
θ
)
−
L
(
θ
(
i
)
)
=
log
(
∑
Z
P
(
Y
∣
Z
;
θ
)
P
(
Z
;
θ
)
)
−
log
P
(
Y
;
θ
(
i
)
)
L(\theta)-L\left(\theta^{(i)}\right)=\log \left(\sum_{\mathbf{Z}} P(Y | Z;\theta) P(Z ; \theta)\right)-\log P\left(Y ; \theta^{(i)}\right)
L(θ)−L(θ(i))=log(Z∑P(Y∣Z;θ)P(Z;θ))−logP(Y;θ(i))利用
J
e
n
s
e
n
Jensen
Jensen不等式得到其下界
L
(
θ
)
−
L
(
θ
(
i
)
)
=
log
(
∑
z
P
(
Y
∣
Z
;
θ
(
i
)
)
P
(
Y
∣
Z
;
θ
)
P
(
Z
;
θ
)
P
(
Y
∣
Z
;
θ
(
i
)
)
)
−
log
P
(
Y
;
θ
(
i
)
)
⩾
∑
Z
P
(
Z
∣
Y
;
θ
(
i
)
)
log
P
(
Y
∣
Z
;
θ
)
P
(
Z
;
θ
)
P
(
Z
∣
Y
;
θ
(
i
)
)
−
1
∗
log
P
(
Y
;
θ
(
i
)
)
=
∑
Z
P
(
Z
∣
Y
;
θ
(
i
)
)
log
P
(
Y
∣
Z
;
θ
)
P
(
Z
;
θ
)
P
(
Z
∣
Y
;
θ
(
i
)
)
−
∑
Z
P
(
Z
∣
Y
;
θ
(
i
)
)
log
P
(
Y
;
θ
(
i
)
)
=
∑
Z
P
(
Z
∣
Y
;
θ
(
i
)
)
log
P
(
Y
∣
Z
;
θ
)
P
(
Z
;
θ
)
P
(
Z
∣
Y
;
θ
(
i
)
)
P
(
Y
;
θ
(
i
)
)
\begin{aligned} L(\theta)-L\left(\theta^{(i)}\right) &=\log \left(\sum_{z} P\left(Y | Z;\theta^{(i)}\right) \frac{P(Y | Z;\theta) P(Z ; \theta)}{P\left(Y | Z;\theta^{(i)}\right)}\right)-\log P\left(Y ; \theta^{(i)}\right) \\ & \geqslant \sum_{Z} P\left(Z | Y;\theta^{(i)}\right) \log \frac{P(Y | Z; \theta) P(Z ; \theta)}{P\left(Z | Y;\theta^{(i)}\right)}-1*\log P\left(Y ; \theta^{(i)}\right) \\& = \sum_{Z} P\left(Z | Y; \theta^{(i)}\right) \log \frac{P(Y | Z; \theta) P(Z ; \theta)}{P\left(Z | Y; \theta^{(i)}\right)}-\sum_{Z} P\left(Z | Y; \theta^{(i)}\right)\log P\left(Y ; \theta^{(i)}\right) \\&=\sum_{Z} P\left(Z | Y; \theta^{(i)}\right) \log \frac{P(Y | Z; \theta) P(Z ; \theta)}{P\left(Z | Y; \theta^{(i)}\right)P\left(Y ; \theta^{(i)}\right)} \end{aligned}
L(θ)−L(θ(i))=log(z∑P(Y∣Z;θ(i))P(Y∣Z;θ(i))P(Y∣Z;θ)P(Z;θ))−logP(Y;θ(i))⩾Z∑P(Z∣Y;θ(i))logP(Z∣Y;θ(i))P(Y∣Z;θ)P(Z;θ)−1∗logP(Y;θ(i))=Z∑P(Z∣Y;θ(i))logP(Z∣Y;θ(i))P(Y∣Z;θ)P(Z;θ)−Z∑P(Z∣Y;θ(i))logP(Y;θ(i))=Z∑P(Z∣Y;θ(i))logP(Z∣Y;θ(i))P(Y;θ(i))P(Y∣Z;θ)P(Z;θ)这里
∑
Z
P
(
Z
∣
Y
;
θ
(
i
)
)
=
1
\sum_{Z} P\left(Z | Y ; \theta^{(i)}\right)=1
∑ZP(Z∣Y;θ(i))=1令
B
(
θ
,
θ
(
i
)
)
≜
L
(
θ
(
i
)
)
+
∑
Z
P
(
Z
∣
Y
;
θ
(
i
)
)
log
P
(
Y
∣
Z
;
θ
)
P
(
Z
;
θ
)
P
(
Z
∣
Y
;
θ
(
i
)
)
P
(
Y
;
θ
(
i
)
)
(
3
)
B\left(\theta, \theta^{(i)}\right) \triangleq L\left(\theta^{(i)}\right)+\sum_{Z} P\left(Z | Y; \theta^{(i)}\right) \log \frac{P(Y | Z; \theta) P(Z ; \theta)}{P\left(Z | Y;\theta^{(i)}\right) P\left(Y;\theta^{(i)}\right)} \qquad(3)
B(θ,θ(i))≜L(θ(i))+Z∑P(Z∣Y;θ(i))logP(Z∣Y;θ(i))P(Y;θ(i))P(Y∣Z;θ)P(Z;θ)(3)则
L
(
θ
)
⩾
B
(
θ
,
θ
(
i
)
)
(
4
)
L(\theta) \geqslant B\left(\theta, \theta^{(i)}\right) \qquad(4)
L(θ)⩾B(θ,θ(i))(4)即函数
B
(
θ
,
θ
(
i
)
)
B\left(\theta, \theta^{(i)}\right)
B(θ,θ(i))是
L
(
θ
)
L(\theta)
L(θ)的一个下界，从而由上式
(
3
)
(3)
(3)得
L
(
θ
(
i
)
)
=
B
(
θ
(
i
)
,
θ
(
i
)
)
(
5
)
L\left(\theta^{(i)}\right)=B\left(\theta^{(i)}, \theta^{(i)}\right) \qquad(5)
L(θ(i))=B(θ(i),θ(i))(5)因此，任何可以使
B
(
θ
,
θ
(
i
)
)
B\left(\theta, \theta^{(i)}\right)
B(θ,θ(i))增大的
θ
\theta
θ,也可以使
L
(
θ
)
L(\theta)
L(θ)增大，为了使
L
(
θ
)
L(\theta)
L(θ)有尽可能地增长，选择
θ
(
i
+
1
)
\theta^{(i+1)}
θ(i+1)使
B
(
θ
,
θ
(
i
)
)
B\left(\theta, \theta^{(i)}\right)
B(θ,θ(i))达到极大，即
θ
(
i
+
1
)
=
arg
max
θ
B
(
θ
,
θ
(
i
)
)
(
6
)
\theta^{(i+1)}=\arg \max _{\theta} B\left(\theta, \theta^{(i)}\right)\qquad(6)
θ(i+1)=argθmaxB(θ,θ(i))(6)
现求
θ
(
i
+
1
)
\theta^{(i+1)}
θ(i+1)的表达式。省去对
θ
\theta
θ的极大化而言是常数的项，由式
(
6
)
(6)
(6)、式
(
3
)
(3)
(3)及
θ
(
i
+
1
)
=
arg
max
θ
Q
(
θ
,
θ
(
i
)
)
\theta^{(i+1)}=\arg \max _{\theta} Q\left(\theta, \theta^{(i)}\right)
θ(i+1)=argmaxθQ(θ,θ(i))，有
θ
(
i
+
1
)
=
arg
max
θ
(
L
(
θ
(
i
)
)
+
∑
z
P
(
Z
∣
Y
;
θ
(
i
)
)
log
P
(
Y
∣
Z
;
θ
)
P
(
Z
;
θ
)
P
(
Z
∣
Y
;
θ
(
i
)
)
P
(
Y
;
θ
(
i
)
)
)
=
arg
max
θ
(
∑
z
P
(
Z
∣
Y
;
θ
(
i
)
)
log
(
P
(
Y
∣
Z
;
θ
)
P
(
Z
;
θ
)
)
)
=
arg
max
θ
(
∑
Z
P
(
Z
∣
Y
;
θ
(
i
)
)
log
P
(
Y
,
Z
;
θ
)
)
=
arg
max
θ
Q
(
θ
,
θ
(
i
)
)
\begin{aligned} \theta^{(i+1)} &=\arg \max _{\theta}\left(L\left(\theta^{(i)}\right)+\sum_{z} P\left(Z | Y;\theta^{(i)}\right) \log \frac{P(Y | Z; \theta) P(Z ; \theta)}{P\left(Z | Y;\theta^{(i)}\right) P\left(Y ; \theta^{(i)}\right)}\right) \\ &=\arg \max _{\theta}\left(\sum_{z} P\left(Z | Y; \theta^{(i)}\right) \log (P(Y | Z; \theta) P(Z ; \theta))\right) \\ &=\arg \max _{\theta}\left(\sum_{Z} P\left(Z | Y; \theta^{(i)}\right) \log P(Y, Z ; \theta)\right) \\ &=\arg \max _{\theta} Q\left(\theta, \theta^{(i)}\right) \end{aligned}
θ(i+1)=argθmax(L(θ(i))+z∑P(Z∣Y;θ(i))logP(Z∣Y;θ(i))P(Y;θ(i))P(Y∣Z;θ)P(Z;θ))=argθmax(z∑P(Z∣Y;θ(i))log(P(Y∣Z;θ)P(Z;θ)))=argθmax(Z∑P(Z∣Y;θ(i))logP(Y,Z;θ))=argθmaxQ(θ,θ(i))注：
θ
(
i
+
1
)
=
arg
max
θ
(
L
(
θ
(
i
)
)
+
∑
z
P
(
Z
∣
Y
;
θ
(
i
)
)
log
P
(
Y
∣
Z
;
θ
)
P
(
Z
;
θ
)
P
(
Z
∣
Y
;
θ
(
i
)
)
P
(
Y
;
θ
(
i
)
)
)
=
arg
max
θ
(
L
(
θ
(
i
)
)
+
∑
z
P
(
Z
∣
Y
;
θ
(
i
)
)
log
P
(
Y
∣
Z
;
θ
)
P
(
Z
;
θ
)
−
∑
z
P
(
Z
∣
Y
;
θ
(
i
)
)
log
(
P
(
Z
∣
Y
;
θ
(
i
)
)
P
(
Y
;
θ
(
i
)
)
)
)
\theta^{(i+1)}=\arg \max _{\theta}\left(L\left(\theta^{(i)}\right)+\sum_{z} P\left(Z | Y ; \theta^{(i)}\right) \log \frac{P(Y | Z ; \theta) P(Z ; \theta)}{P\left(Z | Y ; \theta^{(i)}\right) P\left(Y ; \theta^{(i)}\right)}\right) \\=\arg \max _{\theta}\left(L(\theta^{(i)})+ \sum_{z}P\left(Z | Y ; \theta^{(i)}\right) \log P(Y | Z ; \theta) P(Z ; \theta)-\sum_{z} P\left(Z | Y ; \theta^{(i)}\right)\log\left(P\left(Z | Y ; \theta^{(i)}\right)P\left(Y ; \theta^{(i)}\right) \right)\right)
θ(i+1)=argθmax(L(θ(i))+z∑P(Z∣Y;θ(i))logP(Z∣Y;θ(i))P(Y;θ(i))P(Y∣Z;θ)P(Z;θ))=argθmax(L(θ(i))+z∑P(Z∣Y;θ(i))logP(Y∣Z;θ)P(Z;θ)−z∑P(Z∣Y;θ(i))log(P(Z∣Y;θ(i))P(Y;θ(i))))式
(
7
)
(7)
(7)等价于
E
M
\mathrm{EM}
EM算法的一次迭代，即求
Q
Q
Q函数及其极大化，
E
M
\mathrm{EM}
EM算法是通过不断求解下界的极大化逼近求解对数似然函数极大化的算法
E M \mathrm{EM} EM算法的直观解释：图中上方曲线为 L ( θ ) L(\theta) L(θ),下方曲线为 B ( θ , θ ( i ) ) B\left(\theta, \theta^{(i)}\right) B(θ,θ(i))，由式 ( 4 ) (4) (4), B ( θ , θ ( i ) ) B\left(\theta, \theta^{(i)}\right) B(θ,θ(i))为对数似然函数 L ( θ ) L(\theta) L(θ)的下界，由式 ( 5 ) (5) (5)，两个函数在点 θ = θ ( i ) \theta=\theta^{(i)} θ=θ(i)处相等，由式 ( 6 ) (6) (6), ( 7 ) (7) (7)，算法 E M \mathrm{EM} EM找到下一个点 θ ( i + 1 ) \theta^{(i+1)} θ(i+1)使函数 B ( θ , θ ( i ) ) B\left(\theta, \theta^{(i)}\right) B(θ,θ(i))极大化，也使函数 Q ( θ , θ ( i ) ) Q\left(\theta, \theta^{(i)}\right) Q(θ,θ(i))极大化，这时由于 L ( θ ) ⩾ B ( θ , θ ( i ) ) L(\theta) \geqslant B\left(\theta, \theta^{(i)}\right) L(θ)⩾B(θ,θ(i))，函数 B ( θ , θ ( i ) ) B\left(\theta, \theta^{(i)}\right) B(θ,θ(i))的增加，保证对数似然函数 L ( θ ) L(\theta) L(θ)在每次迭代中也是增加的， E M \mathrm{EM} EM算法在点 θ ( i + 1 ) \theta^{(i+1)} θ(i+1)重新计算 Q Q Q函数值，进行下一次迭代，在这个过程中，对数似然函数 L ( θ ) L(\theta) L(θ)不断增大，从图中可以推断出 E M \mathrm{EM} EM算法不能保证找到全局最优值
E
M
\mathrm{EM}
EM算法的收敛性
定理1 设
P
(
Y
;
θ
)
P(Y;\theta)
P(Y;θ)是观测数据的似然函数，
θ
(
i
)
(
i
=
1
,
2
,
⋯
)
\theta^{(i)}(i=1,2, \cdots)
θ(i)(i=1,2,⋯)为
E
M
\mathrm{EM}
EM算法得到的参数估计序列，
P
(
Y
;
θ
(
i
)
)
(
i
=
1
,
2
,
⋯
)
P\left(Y;\theta^{(i)}\right)(i=1,2, \cdots)
P(Y;θ(i))(i=1,2,⋯)为对应的似然函数序列，
P
(
Y
;
θ
(
i
)
)
P\left(Y;\theta^{(i)}\right)
P(Y;θ(i))是单调递增的，即
P
(
Y
;
θ
(
i
+
1
)
)
⩾
P
(
Y
;
θ
(
i
)
)
(
8
)
P\left(Y;\theta^{(i+1)}\right) \geqslant P\left(Y;\theta^{(i)}\right)\qquad(8)
P(Y;θ(i+1))⩾P(Y;θ(i))(8)证明：由于
P
(
Y
;
θ
)
=
P
(
Y
,
Z
;
θ
)
P
(
Z
∣
Y
;
θ
)
P(Y;\theta)=\frac{P(Y, Z;\theta)}{P(Z | Y;\theta)}
P(Y;θ)=P(Z∣Y;θ)P(Y,Z;θ)取对数有
log
P
(
Y
;
θ
)
=
log
P
(
Y
,
Z
;
θ
)
−
log
P
(
Z
∣
Y
;
θ
)
\log P(Y;\theta)=\log P(Y, Z;\theta)-\log P(Z | Y;\theta)
logP(Y;θ)=logP(Y,Z;θ)−logP(Z∣Y;θ)
由式
(
1
)
(1)
(1)
Q
(
θ
,
θ
(
i
)
)
=
∑
Z
log
P
(
Y
,
Z
;
θ
)
P
(
Z
∣
Y
;
θ
(
i
)
)
Q\left(\theta, \theta^{(i)}\right)=\sum_{Z} \log P(Y, Z;\theta) P\left(Z | Y;\theta^{(i)}\right)
Q(θ,θ(i))=Z∑logP(Y,Z;θ)P(Z∣Y;θ(i))令
H
(
θ
,
θ
(
i
)
)
=
∑
Z
log
P
(
Z
∣
Y
;
θ
)
P
(
Z
∣
Y
;
θ
(
i
)
)
(
9
)
H\left(\theta, \theta^{(i)}\right)=\sum_{Z} \log P(Z | Y;\theta) P\left(Z | Y;\theta^{(i)}\right)\qquad(9)
H(θ,θ(i))=Z∑logP(Z∣Y;θ)P(Z∣Y;θ(i))(9)于是对数似然函数可以写成
log
P
(
Y
;
θ
)
=
Q
(
θ
,
θ
(
i
)
)
−
H
(
θ
,
θ
(
i
)
)
(
10
)
\log P(Y;\theta)=Q\left(\theta, \theta^{(i)}\right)-H\left(\theta, \theta^{(i)}\right)\qquad(10)
logP(Y;θ)=Q(θ,θ(i))−H(θ,θ(i))(10)在式
(
10
)
(10)
(10)中分别
θ
\theta
θ取为
θ
(
i
)
\theta^{(i)}
θ(i)和
θ
(
i
+
1
)
\theta^{(i+1)}
θ(i+1)并相减，有
log
P
(
Y
;
θ
(
i
+
1
)
)
−
log
P
(
Y
;
θ
(
i
)
)
\log P\left(Y ;\theta^{(i+1)}\right)-\log P\left(Y ;\theta^{(i)}\right)
logP(Y;θ(i+1))−logP(Y;θ(i))
=
[
Q
(
θ
(
i
+
1
)
,
θ
(
i
)
)
−
Q
(
θ
(
i
)
,
θ
(
i
)
)
]
−
[
H
(
θ
(
i
+
1
)
,
θ
(
i
)
)
−
H
(
θ
(
i
)
,
θ
(
i
)
)
]
(
11
)
=\left[Q\left(\theta^{(i+1)}, \theta^{(i)}\right)-Q\left(\theta^{(i)}, \theta^{(i)}\right)\right]-\left[H\left(\theta^{(i+1)}, \theta^{(i)}\right)-H\left(\theta^{(i)}, \theta^{(i)}\right)\right] \qquad(11)
=[Q(θ(i+1),θ(i))−Q(θ(i),θ(i))]−[H(θ(i+1),θ(i))−H(θ(i),θ(i))](11)为证式
(
8
)
(8)
(8)，只需证明
(
11
)
(11)
(11)右端是非负的，式
(
11
)
(11)
(11)右端的第1项，由于
θ
(
i
+
1
)
\theta^{(i+1)}
θ(i+1)使
Q
(
θ
,
θ
(
i
)
)
Q\left(\theta, \theta^{(i)}\right)
Q(θ,θ(i))达到极大，所以有
Q
(
θ
(
i
+
1
)
,
θ
(
i
)
)
−
Q
(
θ
(
i
)
,
θ
(
i
)
)
⩾
0
(
12
)
Q\left(\theta^{(i+1)}, \theta^{(i)}\right)-Q\left(\theta^{(i)}, \theta^{(i)}\right) \geqslant 0 \qquad(12)
Q(θ(i+1),θ(i))−Q(θ(i),θ(i))⩾0(12)其第2项，由式
(
9
)
(9)
(9)可得：
H
(
θ
(
i
+
1
)
,
θ
(
i
)
)
−
H
(
θ
(
i
)
,
θ
(
i
)
)
=
∑
Z
(
log
P
(
Z
∣
Y
,
;
t
h
e
t
a
(
i
+
1
)
)
P
(
Z
∣
Y
;
θ
(
i
)
)
)
P
(
Z
∣
Y
;
θ
(
i
)
)
⩽
log
(
∑
z
P
(
Z
∣
Y
;
θ
(
i
+
1
)
)
P
(
Z
∣
Y
;
θ
(
i
)
)
P
(
Z
∣
Y
;
θ
(
i
)
)
)
=
log
[
∑
z
(
P
(
Z
∣
Y
;
θ
(
i
+
1
)
)
]
=
0
(
13
)
H\left(\theta^{(i+1)}, \theta^{(i)}\right)-H\left(\theta^{(i)}, \theta^{(i)}\right)\\ =\sum_{Z}\left(\log \frac{P\left(Z | Y,;theta^{(i+1)}\right)}{P\left(Z | Y; \theta^{(i)}\right)}\right) P\left(Z | Y;\theta^{(i)}\right)\\ \leqslant \log \left(\sum_{z} \frac{P\left(Z | Y; \theta^{(i+1)}\right)}{P\left(Z | Y; \theta^{(i)}\right)} P\left(Z | Y;\theta^{(i)}\right)\right) \\=\log\left[\sum_z( P\left(Z | Y;\theta^{(i+1)}\right )\right]=0 \qquad(13)
H(θ(i+1),θ(i))−H(θ(i),θ(i))=Z∑(logP(Z∣Y;θ(i))P(Z∣Y,;theta(i+1)))P(Z∣Y;θ(i))⩽log(z∑P(Z∣Y;θ(i))P(Z∣Y;θ(i+1))P(Z∣Y;θ(i)))=log[z∑(P(Z∣Y;θ(i+1))]=0(13)这里的不等号由
J
e
n
s
e
n
Jensen
Jensen不等式得到
由式
(
12
)
(12)
(12)和
(
13
)
(13)
(13)知
(
11
)
(11)
(11)右端是非负的。
定理2 设
L
(
θ
)
=
log
P
(
Y
;
θ
)
L(\theta)=\log P(Y ; \theta)
L(θ)=logP(Y;θ)为观测数据的对数似然函数，
θ
(
i
)
(
i
=
1
,
2
,
⋯
)
\theta^{(i)}(i=1,2, \cdots)
θ(i)(i=1,2,⋯)为
EM
\text{EM}
EM算法得到的参数估计，
L
(
θ
(
i
)
)
(
i
=
1
,
2
,
⋯
)
L\left(\theta^{(i)}\right)(i=1,2, \cdots)
L(θ(i))(i=1,2,⋯)为对应的对数似然函数序列。
(
1
)
(1)
(1)如果
P
(
Y
;
θ
)
P(Y ;\theta)
P(Y;θ)有上界，则
L
(
θ
(
i
)
)
=
log
P
(
Y
;
θ
(
i
)
)
L\left(\theta^{(i)}\right)=\log P\left(Y; \theta^{(i)}\right)
L(θ(i))=logP(Y;θ(i))收敛到某一值
L
∗
L^*
L∗
(
2
)
(2)
(2)在函数
Q
(
θ
,
θ
′
)
Q(\theta, \theta^{'})
Q(θ,θ′)与
L
(
θ
)
L(\theta)
L(θ)满足一定条件下，由
EM
\text{EM}
EM算法得到的参数估计序列
θ
(
i
)
\theta^{(i)}
θ(i)的收敛值
θ
∗
\theta^*
θ∗是
L
(
θ
)
L(\theta)
L(θ)的稳定点。
E
M
\mathrm{EM}
EM算法在高斯混合模型学习中的应用
定义(高斯混合模型) 高斯混合模型是指具有如下形式的概率分布模型:
P
(
y
;
θ
)
=
∑
k
=
1
K
α
k
ϕ
(
y
∣
θ
k
)
(
14
)
P(y ; \theta)=\sum_{k=1}^{K} \alpha_{k} \phi\left(y | \theta_{k}\right)\qquad(14)
P(y;θ)=k=1∑Kαkϕ(y∣θk)(14)其中，
α
k
\alpha_{k}
αk是系数，
α
k
⩾
0
,
∑
k
=
1
K
α
k
=
1
\alpha_{k} \geqslant 0, \quad \sum_{k=1}^{K} \alpha_{k}=1
αk⩾0,∑k=1Kαk=1；
ϕ
(
y
∣
θ
k
)
\phi\left(y | \theta_{k}\right)
ϕ(y∣θk)是高斯分布密度，
θ
k
=
(
μ
k
,
σ
k
2
)
\theta_{k}=\left(\mu_{k}, \sigma_{k}^{2}\right)
θk=(μk,σk2)
ϕ
(
y
;
θ
k
)
=
1
2
π
σ
k
exp
(
−
(
y
−
μ
k
)
2
2
σ
k
2
)
(
15
)
\phi\left(y ; \theta_{k}\right)=\frac{1}{\sqrt{2 \pi} \sigma_{k}} \exp \left(-\frac{\left(y-\mu_{k}\right)^{2}}{2 \sigma_{k}^{2}}\right)\qquad(15)
ϕ(y;θk)=2πσk1exp(−2σk2(y−μk)2)(15)称为第
k
k
k个分模型。
一个混合模型可以由任意概率分布密度代替
(
15
)
(15)
(15)中的高斯分布.
高斯混合模型参数估计的
E
M
\mathrm{EM}
EM算法
假设观测数据
y
1
,
y
2
,
⋯
,
y
N
y_{1}, y_{2}, \cdots, y_{N}
y1,y2,⋯,yN由高斯混合模型生成，
P
(
y
∣
θ
)
=
∑
k
=
1
K
α
k
ϕ
(
y
∣
θ
k
)
P(y | \theta)=\sum_{k=1}^{K} \alpha_{k} \phi\left(y | \theta_{k}\right)
P(y∣θ)=k=1∑Kαkϕ(y∣θk)其中，
θ
=
(
α
1
,
α
2
,
⋯
,
α
K
;
θ
1
,
θ
2
,
⋯
,
θ
K
)
\theta=\left(\alpha_{1}, \alpha_{2}, \cdots, \alpha_{K} ; \theta_{1}, \theta_{2}, \cdots, \theta_{K}\right)
θ=(α1,α2,⋯,αK;θ1,θ2,⋯,θK)。我们用
E
M
\mathrm{EM}
EM算法估计高斯混合模型的参数
θ
\theta
θ
1 明确隐变量，写出完全数据的对数似然函数
可以设想观测数据
y
j
,
j
=
1
,
2
,
⋯
,
N
y_{j}, \quad j=1,2, \cdots, N
yj,j=1,2,⋯,N，是这样产生的：首先依概率
α
k
\alpha_k
αk选择第
k
k
k个高斯分布模型
ϕ
(
y
∣
θ
k
)
\phi\left(y | \theta_{k}\right)
ϕ(y∣θk)；然后依第
k
k
k个分模型的概率分布
ϕ
(
y
∣
θ
k
)
\phi\left(y | \theta_{k}\right)
ϕ(y∣θk)生成观测数据
y
j
y_j
yj，这时观测数据
y
j
,
j
=
1
,
2
,
⋯
,
N
y_{j}, \quad j=1,2, \cdots, N
yj,j=1,2,⋯,N是已知的；反映观测数据
y
j
y_j
yj来自第
k
k
k个分模型的数据是未知的，
k
=
1
,
2
,
⋯
,
K
k=1,2, \cdots, K
k=1,2,⋯,K,以隐变量
γ
j
k
\gamma_{j k}
γjk表示，其定义如下：当第
j
j
j个观测来自第
k
k
k个分模型
γ
j
k
=
1
\gamma_{j k}=1
γjk=1 ，否则
γ
j
k
=
0
,
j
=
1
,
2
,
⋯
,
N
;
k
=
1
,
2
,
⋯
,
K
(
17
)
\gamma_{j k}=0, \quad j=1,2, \cdots, N ; \quad k=1,2, \cdots, K\qquad(17)
γjk=0,j=1,2,⋯,N;k=1,2,⋯,K(17)
γ
j
k
\gamma_{j k}
γjk是
0
−
1
0-1
0−1随机变量。
有了观测数据
y
j
y_j
yj及未观测数据
γ
j
k
\gamma_{j k}
γjk，那么完全数据是
(
y
j
,
γ
1
,
γ
j
2
,
⋯
,
γ
j
K
)
,
j
=
1
,
2
,
⋯
,
N
\left(y_{j}, \gamma_{1}, \gamma_{j 2}, \cdots, \gamma_{j K}\right), \quad j=1,2, \cdots, N
(yj,γ1,γj2,⋯,γjK),j=1,2,⋯,N于是，可以写出完全数据的似然函数：
P
(
y
,
γ
∣
θ
)
=
∏
j
=
1
N
P
(
y
j
,
γ
j
,
γ
j
2
,
⋯
,
γ
j
K
;
θ
)
=
∏
k
=
1
K
∏
j
=
1
N
[
α
k
ϕ
(
y
j
∣
θ
k
)
]
γ
A
=
∏
k
=
1
K
α
k
n
∏
j
=
1
N
[
ϕ
(
y
j
∣
θ
k
)
]
γ
A
=
∏
k
=
1
K
α
k
n
∏
j
=
1
N
[
1
2
π
σ
k
exp
(
−
(
y
j
−
μ
k
)
2
2
σ
k
2
)
]
γ
n
\begin{aligned} P(y, \gamma | \theta) &=\prod_{j=1}^{N} P\left(y_{j}, \gamma_{j}, \gamma_{j 2}, \cdots, \gamma_{j K} ; \theta\right) \\ &=\prod_{k=1}^{K} \prod_{j=1}^{N}\left[\alpha_{k} \phi\left(y_{j} | \theta_{k}\right)\right]^{\gamma_{A}} \\ &=\prod_{k=1}^{K} \alpha_{k}^{n} \prod_{j=1}^{N}\left[\phi\left(y_{j} | \theta_{k}\right)\right]^{\gamma_{A}} \\ &=\prod_{k=1}^{K} \alpha_{k}^{n} \prod_{j=1}^{N}\left[\frac{1}{\sqrt{2 \pi} \sigma_{k}} \exp \left(-\frac{\left(y_{j}-\mu_{k}\right)^{2}}{2 \sigma_{k}^{2}}\right)\right]^{\gamma_{n}} \end{aligned}
P(y,γ∣θ)=j=1∏NP(yj,γj,γj2,⋯,γjK;θ)=k=1∏Kj=1∏N[αkϕ(yj∣θk)]γA=k=1∏Kαknj=1∏N[ϕ(yj∣θk)]γA=k=1∏Kαknj=1∏N[2πσk1exp(−2σk2(yj−μk)2)]γn
式中
n
k
=
∑
j
=
1
N
γ
j
k
,
∑
k
=
1
K
n
k
=
N
n_{k}=\sum_{j=1}^{N} \gamma_{j k}, \quad \sum_{k=1}^{K} n_{k}=N
nk=∑j=1Nγjk,∑k=1Knk=N那么，完全数据的对数似然函数为
log
P
(
y
,
γ
;
θ
)
=
∑
k
=
1
K
(
n
k
log
α
k
+
∑
j
=
1
N
γ
j
k
[
log
(
1
2
π
)
−
log
σ
k
−
1
2
σ
k
2
(
y
j
−
μ
k
)
2
]
)
\log P(y, \gamma ; \theta)=\sum_{k=1}^{K} \left(n_{k} \log \alpha_{k}+\sum_{j=1}^{N} \gamma_{j k}\left[\log \left(\frac{1}{\sqrt{2 \pi}}\right)-\log \sigma_{k}-\frac{1}{2 \sigma_{k}^{2}}\left(y_{j}-\mu_{k}\right)^{2}\right]\right)
logP(y,γ;θ)=k=1∑K(nklogαk+j=1∑Nγjk[log(2π1)−logσk−2σk21(yj−μk)2])
2
EM
\text{EM}
EM算法的
E
E
E步：确定
Q
Q
Q函数
Q
(
θ
,
θ
(
i
)
)
=
E
[
log
P
(
y
,
γ
∣
θ
)
∣
y
,
θ
(
i
)
]
=
E
{
∑
k
=
1
K
(
n
k
log
α
k
+
∑
j
=
1
N
γ
j
k
[
log
(
1
2
π
)
−
log
σ
k
−
1
2
σ
k
2
(
y
j
−
μ
k
)
2
]
)
}
=
∑
k
=
1
K
{
∑
j
=
1
N
(
E
γ
j
k
)
log
α
k
+
∑
j
=
1
N
(
E
γ
j
k
)
[
log
(
1
2
π
)
−
log
σ
k
−
1
2
σ
k
2
(
y
j
−
μ
k
)
2
]
}
(
18
)
\begin{aligned} Q\left(\theta, \theta^{(i)}\right) &=E\left[\log P(y, \gamma | \theta) | y, \theta^{(i)}\right] \\ &=E\left\{\sum_{k=1}^{K} \left(n_{k} \log \alpha_{k}+\sum_{j=1}^{N} \gamma_{j k}\left[\log \left(\frac{1}{\sqrt{2 \pi}}\right)-\log \sigma_{k}-\frac{1}{2 \sigma_{k}^{2}}\left(y_{j}-\mu_{k}\right)^{2}\right]\right)\right\} \\ &=\sum_{k=1}^{K}\left\{\sum_{j=1}^{N}\left(E \gamma_{j k}\right) \log \alpha_{k}+\sum_{j=1}^{N}\left(E \gamma_{j k}\right)\left[\log \left(\frac{1}{\sqrt{2 \pi}}\right)-\log \sigma_{k}-\frac{1}{2 \sigma_{k}^{2}}\left(y_{j}-\mu_{k}\right)^{2}\right]\right\} \end{aligned} \qquad(18)
Q(θ,θ(i))=E[logP(y,γ∣θ)∣y,θ(i)]=E{k=1∑K(nklogαk+j=1∑Nγjk[log(2π1)−logσk−2σk21(yj−μk)2])}=k=1∑K{j=1∑N(Eγjk)logαk+j=1∑N(Eγjk)[log(2π1)−logσk−2σk21(yj−μk)2]}(18)
这里需要计算
E
(
γ
j
k
∣
y
;
θ
)
E\left(\gamma_{j k} | y;\theta\right)
E(γjk∣y;θ),记为
γ
^
j
k
\hat{\gamma}_{j k}
γ^jk
γ
^
j
k
=
E
(
γ
k
∣
y
;
θ
)
=
P
(
γ
j
k
=
1
∣
y
;
θ
)
=
P
(
γ
j
k
=
1
,
y
j
;
θ
)
∑
k
=
1
K
P
(
γ
j
k
=
1
,
y
j
;
θ
)
=
P
(
y
j
∣
γ
j
k
=
1
;
θ
)
P
(
γ
j
k
=
1
;
θ
)
∑
k
=
1
K
P
(
y
j
∣
γ
j
k
=
1
;
θ
)
P
(
γ
j
k
=
1
;
θ
)
=
α
k
ϕ
(
y
j
;
θ
k
)
∑
k
=
1
K
α
k
ϕ
(
y
j
;
θ
k
)
,
j
=
1
,
2
,
⋯
,
N
;
k
=
1
,
2
,
⋯
,
K
\begin{aligned} \hat{\gamma}_{j k} &=E\left(\gamma_{k} | y;\theta\right)=P\left(\gamma_{j k}=1 | y; \theta\right) \\ &=\frac{P\left(\gamma_{j k}=1, y_{j} ; \theta\right)}{\sum_{k=1}^{K} P\left(\gamma_{j k}=1, y_{j} ; \theta\right)} \\ &=\frac{P\left(y_{j} | \gamma_{j k}=1;\theta\right) P\left(\gamma_{j k}=1 ; \theta\right)}{\sum_{k=1}^{K} P\left(y_{j} | \gamma_{j k}=1; \theta\right) P\left(\gamma_{j k}=1 ; \theta\right)} \\ &=\frac{\alpha_{k} \phi\left(y_{j} ; \theta_{k}\right)}{\sum_{k=1}^{K} \alpha_{k} \phi\left(y_{j} ; \theta_{k}\right)}, \quad j=1,2, \cdots, N ; \quad k=1,2, \cdots, K \end{aligned}
γ^jk=E(γk∣y;θ)=P(γjk=1∣y;θ)=∑k=1KP(γjk=1,yj;θ)P(γjk=1,yj;θ)=∑k=1KP(yj∣γjk=1;θ)P(γjk=1;θ)P(yj∣γjk=1;θ)P(γjk=1;θ)=∑k=1Kαkϕ(yj;θk)αkϕ(yj;θk),j=1,2,⋯,N;k=1,2,⋯,K
γ
^
j
k
\hat{\gamma}_{j k}
γ^jk是在当前模型参数下第
j
j
j个观测数据来自第
k
k
k个分模型的概率，称为分模型
k
k
k对观测数据
y
j
y_j
yj的响应度
将
γ
^
j
k
=
E
γ
j
k
\hat{\gamma}_{j k}=E \gamma_{j k}
γ^jk=Eγjk及
n
k
=
∑
j
=
1
N
E
γ
j
k
n_{k}=\sum_{j=1}^{N} E \gamma_{j k}
nk=∑j=1NEγjk代入式
(
18
)
(18)
(18)即得
Q
(
θ
,
θ
(
i
)
)
=
∑
k
=
1
K
(
n
k
log
α
k
+
∑
k
=
1
N
γ
^
j
k
[
log
(
1
2
π
)
−
log
σ
k
−
1
2
σ
k
2
(
y
j
−
μ
k
)
2
]
)
(
19
)
Q\left(\theta, \theta^{(i)}\right)=\sum_{k=1}^{K} \left(n_{k} \log \alpha_{k}+\sum_{k=1}^{N} \hat{\gamma}_{j k}\left[\log \left(\frac{1}{\sqrt{2 \pi}}\right)-\log \sigma_{k}-\frac{1}{2 \sigma_{k}^{2}}\left(y_{j}-\mu_{k}\right)^{2}\right]\right)\qquad(19)
Q(θ,θ(i))=k=1∑K(nklogαk+k=1∑Nγ^jk[log(2π1)−logσk−2σk21(yj−μk)2])(19)
3 确定
EM
\text{EM}
EM算法的
M
M
M步
迭代的
M
M
M步是求函数
Q
(
θ
,
θ
(
i
)
)
Q\left(\theta, \theta^{(i)}\right)
Q(θ,θ(i))对
θ
\theta
θ的极大值，即求新一轮迭代的模型参数：
θ
(
i
+
1
)
=
arg
max
θ
Q
(
θ
,
θ
(
i
)
)
\theta^{(i+1)}=\arg \max _{\theta} Q\left(\theta, \theta^{(i)}\right)
θ(i+1)=argθmaxQ(θ,θ(i))用
μ
^
k
,
σ
^
k
2
\hat{\mu}_{k}, \hat{\sigma}_{k}^{2}
μ^k,σ^k2及
α
^
k
,
k
=
1
,
2
,
⋯
,
K
\hat{\alpha}_{k}, k=1,2, \cdots, K
α^k,k=1,2,⋯,K，表示
θ
(
i
+
1
)
\theta^{(i+1)}
θ(i+1)的各参数，求
μ
^
k
,
σ
^
k
2
\hat{\mu}_{k}, \hat{\sigma}_{k}^{2}
μ^k,σ^k2只需将式
(
19
)
(19)
(19)分别对
μ
k
,
σ
k
2
\mu_{k}, \sigma_{k}^{2}
μk,σk2求偏导数并令其为
0
0
0，即可得到；求
α
^
k
\hat{\alpha}_{k}
α^k是在
∑
k
=
1
K
α
k
=
1
\sum_{k=1}^{K} \alpha_{k}=1
∑k=1Kαk=1条件下求偏导数并令其为
0
0
0得到的，结果如下：
μ
^
k
=
∑
j
=
1
N
γ
^
j
k
y
j
∑
j
=
1
N
γ
^
j
k
,
k
=
1
,
2
,
⋯
,
K
(
20
)
\hat{\mu}_{k}=\frac{\sum_{j=1}^{N} \hat{\gamma}_{j k} y_{j}}{\sum_{j=1}^{N} \hat{\gamma}_{j k}}, \quad k=1,2, \cdots, K\qquad(20)
μ^k=∑j=1Nγ^jk∑j=1Nγ^jkyj,k=1,2,⋯,K(20)
σ
k
2
^
=
∑
j
=
1
N
γ
^
k
(
y
j
−
μ
k
)
2
∑
j
=
1
N
γ
^
k
,
k
=
1
,
2
,
⋯
,
K
(
21
)
\hat{\sigma_{k}^{2}}=\frac{\sum_{j=1}^{N} \hat{\gamma}_{k}\left(y_{j}-\mu_{k}\right)^{2}}{\sum_{j=1}^{N} \hat{\gamma}_{k}}, \quad k=1,2, \cdots, K \qquad(21)
σk2^=∑j=1Nγ^k∑j=1Nγ^k(yj−μk)2,k=1,2,⋯,K(21)
α
^
k
=
n
k
N
=
∑
j
=
1
N
γ
^
j
k
N
,
k
=
1
,
2
,
⋯
,
K
(
22
)
\hat{\alpha}_{k}=\frac{n_{k}}{N}=\frac{\sum_{j=1}^{N} \hat{\gamma}_{j k}}{N}, \quad k=1,2, \cdots, K\qquad(22)
α^k=Nnk=N∑j=1Nγ^jk,k=1,2,⋯,K(22)重复以上计算，直到对数似然函数值不再有明显的变化为止
《统计学习方法》李航
扫一扫
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
