EM算法
一般地,用 Y Y Y表示观测随机变量的数据, Z Z Z表示隐随机变的数据。 Y Y Y和 Z Z Z连在一起称为完全数据。假设给定观测数据 Y Y Y,其概率分布是 P ( Y ∣ θ ) P(Y|\theta) P(Y∣θ),其中 θ \theta θ是需要估计地模型参数。
EM
算法通过迭代求
L
(
θ
)
=
l
o
g
P
(
Y
∣
θ
)
L(\theta)=logP(Y|\theta)
L(θ)=logP(Y∣θ)的极大似然估计。每次迭代包含两步:
E
E
E步,求期望;
M
M
M步,求极大化。
1. 算法流程
输 入 : 观 测 变 量 数 据 Y , 隐 变 量 数 据 Z , 联 合 分 布 P ( Y , Z ∣ θ ) 输入:观测变量数据Y,隐变量数据Z,联合分布P(Y,Z|\theta) 输入:观测变量数据Y,隐变量数据Z,联合分布P(Y,Z∣θ),条件分布 P ( Z ∣ Y , θ ) P(Z|Y,\theta) P(Z∣Y,θ)
输 出 : 模 型 参 数 θ 输出:模型参数\theta 输出:模型参数θ
- 选择参数的初值 θ ( 0 ) \theta^{(0)} θ(0),开始迭代
-
E
E
E步:记
θ
(
i
)
\theta^{(i)}
θ(i)为第
i
i
i次迭代参数的估计值,在第
i
+
1
i+1
i+1次迭代的
E
E
E步,计算
Q ( θ , θ ( i ) ) = E z [ l o g P ( Y , Z ∣ θ ) ∣ Y , θ ( i ) ] = ∑ Z l o g P ( Y , Z ∣ θ ) P ( Z ∣ Y , θ ( i ) ) Q(\theta,\theta^{(i)})=E_z[logP(Y,Z|\theta)|Y,\theta^{(i)}]=\sum_ZlogP(Y,Z|\theta)P(Z|Y,\theta^{(i)}) Q(θ,θ(i))=Ez[logP(Y,Z∣θ)∣Y,θ(i)]=Z∑logP(Y,Z∣θ)P(Z∣Y,θ(i))
其中 P ( Z ∣ Y , θ ( i ) ) P(Z|Y,\theta^{(i)}) P(Z∣Y,θ(i))是在给定观测数据 Y Y Y和当前参数估计 θ ( i ) \theta^{(i)} θ(i)下隐变量数据 Z Z Z的条件概率分布 -
M
M
M步:求使
Q
(
θ
,
θ
(
i
)
)
Q(\theta,\theta^{(i)})
Q(θ,θ(i))的极大化的
θ
\theta
θ,确定第
i
+
1
i+1
i+1次迭代的参数估计值
θ
(
i
+
1
)
\theta^{(i+1)}
θ(i+1)
θ ( i + 1 ) = a r g max θ Q ( θ , θ ( i ) ) \theta^{(i+1)}=arg\max_{\theta}Q(\theta,\theta^{(i)}) θ(i+1)=argθmaxQ(θ,θ(i)) - 重复第2步和第3步,直到收敛
函数 Q ( θ , θ ( i ) ) Q(\theta,\theta^{(i)}) Q(θ,θ(i))是EM算法的核心,称为Q函数
完全数据的对数似然函数 l o g P ( Y , Z ∣ θ ) logP(Y,Z|\theta) logP(Y,Z∣θ)关于在给定观测数据 Y Y Y和当前参数 θ ( i ) \theta^{(i)} θ(i)下对未观测数据 Z Z Z的条件概率分布 P ( Z ∣ Y , θ ( i ) ) P(Z|Y,\theta^{(i)}) P(Z∣Y,θ(i))的期望称为Q函数
2. EM算法的导出
对于含有隐变量的概率模型,目标是极大化观测数据
Y
Y
Y关于参数
θ
\theta
θ的对数似然函数,即极大化
L
(
θ
)
=
l
o
g
P
(
Y
∣
θ
)
=
l
o
g
∑
Z
P
(
Y
,
Z
∣
θ
)
=
l
o
g
∑
Z
P
(
Y
∣
Z
,
θ
)
P
(
Z
∣
θ
)
L(\theta)=logP(Y|\theta)=log\sum_ZP(Y,Z|\theta)=log\sum_ZP(Y|Z,\theta)P(Z|\theta)
L(θ)=logP(Y∣θ)=logZ∑P(Y,Z∣θ)=logZ∑P(Y∣Z,θ)P(Z∣θ)
E
M
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(\theta^{(i)})
L(θ)>L(θ(i)),并逐步到达极大值。为此,考虑两者的差:
L ( θ ) − L ( θ ( i ) ) = l o g ∑ Z P ( Y ∣ Z , θ ) P ( Z ∣ θ ) − l o g P ( Y ∣ θ ( i ) ) L(\theta)-L(\theta^{(i)})=log\sum_ZP(Y|Z,\theta)P(Z|\theta)-logP(Y|\theta^{(i)}) L(θ)−L(θ(i))=logZ∑P(Y∣Z,θ)P(Z∣θ)−logP(Y∣θ(i))
利用Jensen
不等式:
l o g ∑ j λ j y i ≥ ∑ j λ j l o g y i , λ j ≥ 0 , ∑ j λ j = 1 log\sum_j\lambda_jy_i \geq \sum_j\lambda_jlogy_i, \quad \lambda_j \geq0,\sum_j\lambda_j=1 logj∑λjyi≥j∑λjlogyi,λj≥0,j∑λj=1
得到下界:
L
(
θ
)
−
L
(
θ
(
i
)
)
=
l
o
g
(
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
P
(
Y
∣
Z
,
θ
)
P
(
Z
∣
θ
)
P
(
Z
∣
Y
,
θ
(
i
)
)
)
−
l
o
g
P
(
Y
∣
θ
(
i
)
)
≥
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
l
o
g
P
(
Y
∣
Z
,
θ
)
P
(
Z
∣
θ
)
P
(
Z
∣
Y
,
θ
(
i
)
)
−
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
l
o
g
P
(
Y
∣
θ
(
i
)
)
=
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
l
o
g
P
(
Y
∣
Z
,
θ
)
P
(
Z
∣
θ
)
P
(
Z
∣
Y
,
θ
(
i
)
)
P
(
Y
∣
θ
(
i
)
)
L(\theta)-L(\theta^{(i)})=log(\sum_ZP(Z|Y,\theta^{(i)})\frac{P(Y|Z,\theta)P(Z|\theta)} {P(Z|Y,\theta^{(i)})})-logP(Y|\theta^{(i)}) \\ \geq \sum_ZP(Z|Y,\theta^{(i)})log\frac{P(Y|Z,\theta)P(Z|\theta)} {P(Z|Y,\theta^{(i)})} - \sum_ZP(Z|Y,\theta^{(i)})logP(Y|\theta^{(i)}) \\ =\sum_ZP(Z|Y,\theta^{(i)})log\frac{P(Y|Z,\theta)P(Z|\theta)} {P(Z|Y,\theta^{(i)})P(Y|\theta^{(i)})} \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad
L(θ)−L(θ(i))=log(Z∑P(Z∣Y,θ(i))P(Z∣Y,θ(i))P(Y∣Z,θ)P(Z∣θ))−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∣θ)
令
B
(
θ
,
θ
(
i
)
)
=
L
(
θ
(
i
)
)
+
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
l
o
g
P
(
Y
∣
Z
,
θ
)
P
(
Z
∣
θ
)
P
(
Z
∣
Y
,
θ
(
i
)
)
P
(
Y
∣
θ
(
i
)
)
B(\theta,\theta^{(i)})=L(\theta^{(i)})+\sum_ZP(Z|Y,\theta^{(i)})log\frac{P(Y|Z,\theta)P(Z|\theta)} {P(Z|Y,\theta^{(i)})P(Y|\theta^{(i)})}
B(θ,θ(i))=L(θ(i))+Z∑P(Z∣Y,θ(i))logP(Z∣Y,θ(i))P(Y∣θ(i))P(Y∣Z,θ)P(Z∣θ)
则
L ( θ ) ≥ B ( θ , θ ( i ) ) L(\theta) \geq B(\theta, \theta^{(i)}) L(θ)≥B(θ,θ(i))
因此,任何可以使
B
(
θ
,
θ
(
i
)
)
B(\theta,\theta^{(i)})
B(θ,θ(i))增大的
θ
\theta
θ,也可以使
L
(
θ
)
L(\theta)
L(θ)增大。为了使
L
(
θ
)
L(\theta)
L(θ)有尽可能大的增大,选择
θ
(
i
+
1
)
\theta^{(i+1)}
θ(i+1)使
B
(
θ
,
θ
(
i
)
)
B(\theta,\theta^{(i)})
B(θ,θ(i))达到极大,即
θ
(
i
+
1
)
=
a
r
g
max
θ
B
(
θ
,
θ
(
i
)
)
\theta^{(i+1)}=arg\max_\theta B(\theta,\theta^{(i)})
θ(i+1)=argθmaxB(θ,θ(i))
θ
(
i
+
1
)
=
a
r
g
max
θ
(
L
(
θ
(
i
)
)
+
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
l
o
g
P
(
Y
∣
Z
,
θ
)
P
(
Z
∣
θ
)
P
(
Z
∣
Y
,
θ
(
i
)
)
P
(
Y
∣
θ
(
i
)
)
)
=
a
r
g
max
θ
(
∑
Z
P
(
Z
∣
Y
,
θ
(
i
)
)
l
o
g
(
P
(
Y
∣
Z
,
θ
)
P
(
Z
∣
θ
)
)
)
=
a
r
g
max
θ
(
∑
Z
P
(
Z
∣
Y
,
θ
i
)
l
o
g
P
(
Y
,
Z
∣
θ
)
)
=
a
r
g
max
θ
Q
(
θ
,
θ
(
i
)
)
\theta^{(i+1)}=arg\max_\theta (L(\theta^{(i)})+\sum_ZP(Z|Y,\theta^{(i)})log\frac{P(Y|Z,\theta)P(Z|\theta)} {P(Z|Y,\theta^{(i)})P(Y|\theta^{(i)})}) \\ = arg\max_\theta (\sum_ZP(Z|Y,\theta^{(i)})log(P(Y|Z,\theta)P(Z|\theta))) \quad \quad \quad \\ = arg\max_\theta (\sum_ZP(Z|Y,\theta^{i})logP(Y,Z|\theta))=arg\max_\theta Q(\theta, \theta^{(i)})
θ(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))
3. EM算法在高斯混合模型学习中的应用
3.1 高斯混合模型
高斯混合模型是指具有如下形式的概率分布模型:
P
(
y
∣
θ
)
=
∑
k
=
1
K
α
k
ϕ
(
y
∣
θ
k
)
P(y|\theta)=\sum_{k=1}^K\alpha_k\phi(y|\theta_k)
P(y∣θ)=k=1∑Kαkϕ(y∣θk)
其中,
α
k
\alpha_k
αk是系数,
α
k
≥
0
,
∑
k
=
1
K
α
k
=
1
;
ϕ
(
y
∣
θ
k
)
\alpha_k \geq0, \sum_{k=1}^K\alpha_k=1;\phi(y|\theta_k)
αk≥0,∑k=1Kαk=1;ϕ(y∣θk)是高斯密度函数,
θ
k
=
(
μ
k
,
σ
k
)
\theta_k=(\mu_k,\sigma_k)
θk=(μk,σk)
ϕ
(
y
∣
θ
k
)
=
1
2
π
σ
k
e
x
p
(
−
(
y
−
μ
k
)
2
2
σ
k
2
)
\phi(y|\theta_k)=\frac {1} {\sqrt {2 \pi}\sigma_k}exp(-\frac {(y-\mu_k)^2} {2\sigma_k^2})
ϕ(y∣θk)=2πσk1exp(−2σk2(y−μk)2)
称为第
k
k
k个模型
3.2 高斯混合模型参数估计的EM算法
3.2.1 问题描述
假设观测数据
y
1
,
y
2
,
.
.
.
,
y
N
y_1,y_2,...,y_N
y1,y2,...,yN由高斯混合模型生成
P
(
y
∣
θ
)
=
∑
k
=
1
K
α
k
ϕ
(
y
∣
θ
k
)
P(y|\theta)=\sum_{k=1}^K\alpha_k\phi(y|\theta_k)
P(y∣θ)=k=1∑Kαkϕ(y∣θk)
其中,
θ
=
(
α
1
,
α
2
,
.
.
.
,
α
K
;
θ
1
,
θ
2
,
.
.
.
,
θ
K
)
\theta=(\alpha_1,\alpha_2,...,\alpha_K;\theta_1,\theta_2,...,\theta_K)
θ=(α1,α2,...,αK;θ1,θ2,...,θK)是我们需要利用
E
M
EM
EM算法估计的参数
3.2.2 明确隐变量,写出完全数据的对数似然函数
观测数据 y j , j = 1 , 2 , . . . , N y_j,j=1,2,...,N yj,j=1,2,...,N是先以概率 α k \alpha_k αk选择第 k k k个高斯分布模型 ϕ ( y ∣ θ k ) \phi(y|\theta_k) ϕ(y∣θk),然后通过这个高斯分布模型生成观测数据 y j y_j yj
观测数据
y
j
y_j
yj是已知的,反映观测数据
y
j
y_j
yj来自第
k
k
k个高斯分布模型是未知的,以隐变量
γ
j
k
\gamma_{jk}
γjk表示,定义如下:
γ
j
k
=
{
1
,
第j个观测数据来自第k个高斯分布模型
0
,
other
j
=
1
,
2
,
.
.
.
,
N
;
k
=
1
,
2
,
.
.
.
,
K
\gamma_{jk}= \begin{cases} 1, & \text{第j个观测数据来自第k个高斯分布模型} \\ 0, & \text{other} \end{cases} \\ j=1,2,...,N; \quad k=1,2,...,K \quad \quad \quad \quad \quad \quad
γjk={1,0,第j个观测数据来自第k个高斯分布模型otherj=1,2,...,N;k=1,2,...,K
有了观测数据
y
j
y_j
yj和未观测数据
γ
j
k
\gamma_{jk}
γjk,那么完全数据是
(
y
j
,
γ
j
1
,
γ
j
2
,
.
.
.
,
γ
j
K
)
(y_j,\gamma_{j1},\gamma_{j2},...,\gamma_{jK})
(yj,γj1,γj2,...,γjK)
完全数据的似然函数:
P
(
y
,
γ
∣
θ
)
=
∏
j
=
1
N
P
(
y
j
,
γ
j
1
,
γ
j
2
,
.
.
.
,
γ
j
K
∣
θ
)
=
∏
j
=
1
N
∏
k
=
1
K
[
α
k
ϕ
(
y
j
∣
θ
)
]
γ
j
k
=
∏
k
=
1
K
α
k
n
k
∏
j
=
1
N
[
ϕ
(
y
j
∣
k
)
]
γ
j
k
=
∏
k
=
1
K
α
k
n
k
∏
j
=
1
N
[
1
2
π
σ
k
e
x
p
(
−
(
y
−
μ
k
)
2
2
σ
k
2
)
]
γ
j
k
P(y,\gamma|\theta)=\prod_{j=1}^NP(y_j,\gamma_{j1},\gamma_{j2},...,\gamma_{jK}|\theta)\\=\prod_{j=1}^N\prod_{k=1}^K[\alpha_k\phi(y_j|\theta)]^{\gamma_{jk}}\\=\prod_{k=1}^K\alpha_k^{n_k}\prod_{j=1}^N[\phi(y_j|k)]^{\gamma_{jk}}\\=\prod_{k=1}^K\alpha_k^{n_k}\prod_{j=1}^N[\frac {1} {\sqrt {2 \pi}\sigma_k}exp(-\frac {(y-\mu_k)^2} {2\sigma_k^2})]^{\gamma_{jk}}
P(y,γ∣θ)=j=1∏NP(yj,γj1,γj2,...,γjK∣θ)=j=1∏Nk=1∏K[αkϕ(yj∣θ)]γjk=k=1∏Kαknkj=1∏N[ϕ(yj∣k)]γjk=k=1∏Kαknkj=1∏N[2πσk1exp(−2σk2(y−μk)2)]γjk
其中,
n
k
=
∑
j
=
1
N
γ
j
k
,
∑
k
=
1
K
n
k
=
N
n_k=\sum_{j=1}^N\gamma_{jk}, \sum_{k=1}^Kn_k=N
nk=∑j=1Nγjk,∑k=1Knk=N
那么完全数据的对数似然函数为:
l
o
g
P
(
y
,
γ
∣
θ
)
=
∑
k
=
1
K
{
n
k
l
o
g
α
k
+
∑
j
=
1
N
γ
j
k
[
l
o
g
(
1
2
π
)
−
l
o
g
σ
k
−
1
2
σ
k
2
(
y
j
−
μ
k
)
2
]
}
logP(y,\gamma|\theta)=\sum_{k=1}^K\{n_klog\alpha_k+\sum_{j=1}^N\gamma_{jk}[log(\frac {1} {\sqrt {2\pi}})-log\sigma_k-\frac{1}{2\sigma_k^2}(y_j-\mu_k)^2]\}
logP(y,γ∣θ)=k=1∑K{nklogαk+j=1∑Nγjk[log(2π1)−logσk−2σk21(yj−μk)2]}
3.2.3 EM算法的E步:确定Q函数
Q ( θ , θ ( i ) ) = E γ [ l o g P ( y , γ ∣ θ ) ∣ y , θ ( i ) ] = E γ { ∑ k = 1 K { n k l o g α k + ∑ j = 1 N γ j k [ l o g ( 1 2 π ) − l o g σ k − 1 2 σ k 2 ( y j − μ k ) 2 ] } } = ∑ k = 1 K { n k l o g α k + ∑ j = 1 N ( E γ j k ) [ l o g ( 1 2 π ) − l o g σ k − 1 2 σ k 2 ( y j − μ k ) 2 ] } Q(\theta,\theta^{(i)})=E_\gamma[logP(y,\gamma|\theta)|y,\theta^{(i)}]\\=E_\gamma\{\sum_{k=1}^K\{n_klog\alpha_k+\sum_{j=1}^N\gamma_{jk}[log(\frac {1} {\sqrt {2\pi}})-log\sigma_k-\frac{1}{2\sigma_k^2}(y_j-\mu_k)^2]\}\}\\=\sum_{k=1}^K\{n_klog\alpha_k+\sum_{j=1}^N(E\gamma_{jk})[log(\frac {1} {\sqrt {2\pi}})-log\sigma_k-\frac{1}{2\sigma_k^2}(y_j-\mu_k)^2]\} 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{nklogαk+j=1∑N(Eγjk)[log(2π1)−logσk−2σk21(yj−μk)2]}
其中
E
(
γ
r
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
,
.
.
.
,
N
;
k
=
1
,
.
.
,
K
E(\gamma_{rk}|y,\theta)=P(\gamma_{jk}=1|y,\theta)=\frac {P(\gamma_{jk}=1,y_j|\theta)} {\sum_{k=1}^KP(\gamma_{jk}=1,y_j|\theta)}\\=\frac {P(y_j|\gamma_{jk}=1,\theta)P(\gamma_{jk}=1|\theta)}{\sum_{k=1}^KP(y_j|\gamma_{jk}=1,\theta)P(\gamma_{jk}=1|\theta)}\\=\frac {\alpha_k\phi(y_j|\theta_k)} {\sum_{k=1}^K\alpha_k\phi(y_j|\theta_k)},\quad j=1,...,N;k=1,..,K
E(γrk∣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,...,N;k=1,..,K
E ( γ j k ) E(\gamma_{jk}) E(γjk)表示在当前模型参数下第 j j j个观测数据来自第 k k k个高斯分布模型的概率,称为第 k k k个高斯分布模型对观测数据 y j y_j yj的响应度。
3.2.4 确定EM算法的M步
θ
(
i
+
1
)
=
a
r
g
max
θ
Q
(
θ
,
θ
(
i
)
)
\theta^{(i+1)}=arg\max_\theta Q(\theta,\theta^{(i)})
θ(i+1)=argθmaxQ(θ,θ(i))
分别对
μ
k
,
σ
k
2
\mu_k,\sigma_k^2
μk,σk2求偏导并令其为0
∂
Q
(
θ
,
θ
(
i
)
)
∂
μ
k
=
∑
j
=
1
N
E
(
γ
j
k
)
y
i
−
μ
k
σ
k
2
=
0
∑
j
=
1
N
E
(
γ
j
k
)
(
y
j
−
μ
k
)
=
0
μ
^
k
=
∑
j
=
1
N
E
(
γ
j
k
)
y
j
∑
j
=
1
N
E
(
γ
j
k
)
\frac {\partial Q(\theta,\theta^{(i)})}{\partial \mu_k}=\sum_{j=1}^N E(\gamma_{jk})\frac{y_i-\mu_k}{\sigma_k^2}=0\\ \sum_{j=1}^NE(\gamma_{jk})(y_j-\mu_k)=0 \\ \hat \mu_k=\frac {\sum_{j=1}^NE(\gamma_{jk})y_j} {\sum_{j=1}^NE(\gamma_{jk})}
∂μk∂Q(θ,θ(i))=j=1∑NE(γjk)σk2yi−μk=0j=1∑NE(γjk)(yj−μk)=0μ^k=∑j=1NE(γjk)∑j=1NE(γjk)yj
∂ Q ( θ , θ ( i ) ) ∂ σ k 2 = ∑ j = 1 N E ( γ j k ) [ − 1 2 σ k 2 + ( y j − μ k ) 2 2 ( σ k 2 ) 2 ] = 0 ∑ j = 1 N E ( γ j k ) [ ( y j − μ k ) 2 − σ k 2 ] = 0 σ ^ k 2 = ∑ j = 1 N E ( γ j k ) ( y j − μ k ) 2 ∑ j = 1 N E ( γ j k ) \frac {\partial Q(\theta,\theta^{(i)})}{\partial \sigma_k^2}=\sum_{j=1}^N E(\gamma_{jk})[-\frac {1} {2\sigma_k^2}+\frac {(y_j-\mu_k)^2} {2(\sigma_k^2)^2}]=0\\\sum_{j=1}^N E(\gamma_{jk})[(y_j-\mu_k)^2-\sigma_k^2]=0 \\ \hat \sigma_k^2=\frac {\sum_{j=1}^N E(\gamma_{jk})(y_j-\mu_k)^2} {\sum_{j=1}^N E(\gamma_{jk})} ∂σk2∂Q(θ,θ(i))=j=1∑NE(γjk)[−2σk21+2(σk2)2(yj−μk)2]=0j=1∑NE(γjk)[(yj−μk)2−σk2]=0σ^k2=∑j=1NE(γjk)∑j=1NE(γjk)(yj−μk)2
在
∑
k
=
1
K
α
k
=
1
\sum_{k=1}^K\alpha_k=1
∑k=1Kαk=1的条件下对
α
k
\alpha_k
αk求偏导并令其为0。采用拉格朗日乘子法,有:
L
(
θ
)
=
Q
(
θ
,
θ
(
i
)
)
+
λ
(
∑
k
=
1
K
α
k
−
1
)
L(\theta)=Q(\theta,\theta^{(i)})+\lambda(\sum_{k=1}^K\alpha_k-1)
L(θ)=Q(θ,θ(i))+λ(k=1∑Kαk−1)
∂ L ( θ ) ∂ α k = n k α k + λ = 0 ⟺ α ^ k = − n k λ ∑ k = 1 K ( n k + α k λ ) = 0 N + λ = 0 ⟺ λ = − N α ^ k = n k N = ∑ j = 1 N E ( γ j k ) N \frac {\partial L(\theta)}{\partial \alpha_k}=\frac{n_k}{\alpha_k}+\lambda=0 \iff \hat \alpha_k=-\frac{n_k}{\lambda}\\ \sum_{k=1}^K(n_k+\alpha_k\lambda)=0\\N+\lambda=0 \iff \lambda=-N \\ \hat \alpha_k=\frac {n_k} {N} = \frac {\sum_{j=1}^N E(\gamma_{jk})} {N} ∂αk∂L(θ)=αknk+λ=0⟺α^k=−λnkk=1∑K(nk+αkλ)=0N+λ=0⟺λ=−Nα^k=Nnk=N∑j=1NE(γjk)
3.3 高斯混合模型参数估计的EM算法
输 入 : 观 测 数 据 y 1 , y 2 , . . . , y N , 高 斯 混 合 模 型 输入:观测数据y_1,y_2,...,y_N,高斯混合模型 输入:观测数据y1,y2,...,yN,高斯混合模型
输 出 : 高 斯 混 合 模 型 参 数 输出:高斯混合模型参数 输出:高斯混合模型参数
- 取参数的初始值开始迭代
-
E
E
E步:依据当前模型参数,计算高斯分布模型
k
k
k对观测数据
y
j
y_j
yj的响应度
E ( γ r k ∣ y , θ ) = α k ϕ ( y j ∣ θ k ) ∑ k = 1 K α k ϕ ( y j ∣ θ k ) E(\gamma_{rk}|y,\theta)=\frac {\alpha_k\phi(y_j|\theta_k)} {\sum_{k=1}^K\alpha_k\phi(y_j|\theta_k)} E(γrk∣y,θ)=∑k=1Kαkϕ(yj∣θk)αkϕ(yj∣θk) -
M
M
M步:计算新一轮迭代的模型参数
μ ^ k = ∑ j = 1 N E ( γ j k ) y j ∑ j = 1 N E ( γ j k ) σ ^ k 2 = ∑ j = 1 N E ( γ j k ) ( y j − μ k ) 2 ∑ j = 1 N E ( γ j k ) α ^ k = ∑ j = 1 N E ( γ j k ) N \hat \mu_k=\frac {\sum_{j=1}^NE(\gamma_{jk})y_j} {\sum_{j=1}^NE(\gamma_{jk})} \\ \hat \sigma_k^2=\frac {\sum_{j=1}^N E(\gamma_{jk})(y_j-\mu_k)^2} {\sum_{j=1}^N E(\gamma_{jk})} \\ \hat \alpha_k= \frac {\sum_{j=1}^N E(\gamma_{jk})} {N} μ^k=∑j=1NE(γjk)∑j=1NE(γjk)yjσ^k2=∑j=1NE(γjk)∑j=1NE(γjk)(yj−μk)2α^k=N∑j=1NE(γjk) - 重复第2步和第3步,直到收敛