机器学习EM算法原理及推导
在机器学习与统计推断中,我们经常会遇到“缺失数据”或“潜在变量”(latent variables)的情形:样本并非完全可观测,而部分信息被隐藏或丢失。这种情况下,直接对观测数据做极大似然估计(Maximum Likelihood Estimation, MLE)往往非常困难,因其对数似然函数通常是非凸、甚至没有封闭解的高维积分。期望—最大化(Expectation–Maximization, EM)算法正是在这种背景下应运而生。EM算法以迭代的方式交替执行“期望步”(E步)和“最大化步”(M步),通过引入辅助函数将难以直接优化的对数似然问题,转化为一系列更易处理的子问题,从而逐步逼近最优解。本文将以通俗且系统的角度,详细剖析EM算法的原理、推导、实现步骤,以及在高斯混合模型等典型案例中的应用,并探讨其优缺点与常见改进方法。
EM算法的核心思想
EM算法的核心在于引入“完整数据”(complete data)的概念。假设观测数据为 X X X,潜在变量为 Z Z Z,模型参数为 θ \theta θ。如果我们可以同时观测 ( X , Z ) (X,Z) (X,Z),则完整数据的对数似然函数
log p ( X , Z ∣ θ ) \log p(X,Z \mid \theta) logp(X,Z∣θ)
通常比只依赖观测数据的对数似然
log p ( X ∣ θ ) = log ∫ p ( X , Z ∣ θ ) d Z \log p(X \mid \theta) = \log \int p(X,Z \mid \theta)\,\mathrm{d}Z logp(X∣θ)=log∫p(X,Z∣θ)dZ
更容易优化。EM算法便是通过以下两步迭代实现参数估计:
-
E步(Expectation): 在当前参数 θ ( t ) \theta^{(t)} θ(t) 下,计算潜在变量 Z Z Z 的后验分布对数似然的期望,即构造辅助函数
Q ( θ ∣ θ ( t ) ) = E Z ∣ X , θ ( t ) [ log p ( X , Z ∣ θ ) ] . Q(\theta \mid \theta^{(t)}) = \mathbb{E}_{Z\mid X,\,\theta^{(t)}}\bigl[\log p(X,Z \mid \theta)\bigr]. Q(θ∣θ(t))=EZ∣X,θ(t)[logp(X,Z∣θ)].
-
M步(Maximization): 固定辅助函数中的期望部分,将 Q ( θ ∣ θ ( t ) ) Q(\theta \mid \theta^{(t)}) Q(θ∣θ(t)) 关于 θ \theta θ 最大化,获得更新参数
θ ( t + 1 ) = arg max θ Q ( θ ∣ θ ( t ) ) . \theta^{(t+1)} = \arg\max_{\theta} Q(\theta \mid \theta^{(t)}). θ(t+1)=argθmaxQ(θ∣θ(t)).
通过反复交替 E 步和 M 步,EM 算法能够保证观测数据的对数似然不会下降,直至收敛到某个局部极值。
EM算法的数学推导
-
对数似然的下界
利用Jensen不等式,我们可以得到:log p ( X ∣ θ ) = log ∫ p ( X , Z ∣ θ ) d Z = log ∫ q ( Z ) p ( X , Z ∣ θ ) q ( Z ) d Z ≥ ∫ q ( Z ) log p ( X , Z ∣ θ ) q ( Z ) d Z , \log p(X \mid \theta) = \log \int p(X,Z \mid \theta)\,\mathrm{d}Z = \log \int q(Z)\,\frac{p(X,Z \mid \theta)}{q(Z)}\,\mathrm{d}Z \ge \int q(Z) \log \frac{p(X,Z \mid \theta)}{q(Z)}\,\mathrm{d}Z, logp(X∣θ)=log∫p(X,Z∣θ)dZ=log∫q(Z)q(Z)p(X,Z∣θ)dZ≥∫q(Z)logq(Z)p(X,Z∣θ)dZ,
对任意分布 q ( Z ) q(Z) q(Z) 均成立。令 q ( Z ) = p ( Z ∣ X , θ ( t ) ) q(Z) = p(Z\mid X,\theta^{(t)}) q(Z)=p(Z∣X,θ(t)),即可得到对数似然的一个“下界”(lower bound),记为
L ( q , θ ) = E q [ log p ( X , Z ∣ θ ) ] − E q [ log q ( Z ) ] . \mathcal{L}(q,\theta) = \mathbb{E}_{q}\bigl[\log p(X,Z\mid\theta)\bigr] - \mathbb{E}_{q}\bigl[\log q(Z)\bigr]. L(q,θ)=Eq[logp(X,Z∣θ)]−Eq[logq(Z)].
-
E步:构造Q函数
在 E 步中,我们固定 θ ( t ) \theta^{(t)} θ(t),选择最佳 q ( Z ) q(Z) q(Z);通过变分推导可证最佳解为真实后验 p ( Z ∣ X , θ ( t ) ) p(Z\mid X,\theta^{(t)}) p(Z∣X,θ(t))。此时下界与对数似然之差为后验的熵项,但对参数更新无影响,于是我们只需最大化第一项:Q ( θ ∣ θ ( t ) ) = E Z ∣ X , θ ( t ) [ log p ( X , Z ∣ θ ) ] . Q(\theta\mid\theta^{(t)}) = \mathbb{E}_{Z\mid X,\theta^{(t)}}[\log p(X,Z\mid\theta)]. Q(θ∣θ(t))=EZ∣X,θ(t)[logp(X,Z∣θ)].
-
M步:更新参数
在 M 步中,固定 Q Q Q 中的期望计算结果,求解θ ( t + 1 ) = arg max θ Q ( θ ∣ θ ( t ) ) . \theta^{(t+1)} = \arg\max_\theta Q(\theta \mid \theta^{(t)}). θ(t+1)=argθmaxQ(θ∣θ(t)).
由于 Q Q Q 一般比原始对数似然更简单(常为凸函数或可分解项),常可得到解析解或快速数值解。
经过上述两步交替,EM算法维护了对数似然的单调递增性:
log p ( X ∣ θ ( t + 1 ) ) ≥ Q ( θ ( t + 1 ) ∣ θ ( t ) ) ≥ Q ( θ ( t ) ∣ θ ( t ) ) = log p ( X ∣ θ ( t ) ) . \log p(X \mid \theta^{(t+1)}) \ge Q(\theta^{(t+1)}\mid\theta^{(t)}) \ge Q(\theta^{(t)}\mid\theta^{(t)}) = \log p(X \mid \theta^{(t)}). logp(X∣θ(t+1))≥Q(θ(t+1)∣θ(t))≥Q(θ(t)∣θ(t))=logp(X∣θ(t)).
因此,算法能保证不退步地爬升,直至局部收敛。
EM算法的伪代码
输入:观测数据 X,初始参数 θ^(0),最大迭代次数 T,收敛阈值 ε
输出:估计参数 θ*for t = 0,1,2,…,T-1 doE步:计算 Q(θ | θ^(t)) = E_{Z|X,θ^(t)}[log p(X,Z | θ)]M步:θ^(t+1) = arg max_θ Q(θ | θ^(t))若 |θ^(t+1) - θ^(t)| < ε,则跳出
end for返回 θ* = θ^(t+1)
应用示例:高斯混合模型
高斯混合模型(Gaussian Mixture Model, GMM)是 EM 算法最常见的应用之一。假设数据集由 K K K 个不同高斯分布生成,模型参数包括各分布的混合系数 π k \pi_k πk、均值 μ k \mu_k μk 和协方差矩阵 Σ k \Sigma_k Σk。
-
E步:计算每个样本属于第 k k k 个高斯成分的后验概率(责任度)
γ i k = π k ( t ) N ( x i ∣ μ k ( t ) , Σ k ( t ) ) ∑ j = 1 K π j ( t ) N ( x i ∣ μ j ( t ) , Σ j ( t ) ) . \gamma_{ik} = \frac{\pi_k^{(t)} \mathcal{N}(x_i \mid \mu_k^{(t)},\Sigma_k^{(t)})}{\sum_{j=1}^K \pi_j^{(t)} \mathcal{N}(x_i \mid \mu_j^{(t)},\Sigma_j^{(t)})}. γik=∑j=1Kπj(t)N(xi∣μj(t),Σj(t))πk(t)N(xi∣μk(t),Σk(t)).
-
M步:基于责任度更新参数
N k = ∑ i = 1 n γ i k , π k ( t + 1 ) = N k n , μ k ( t + 1 ) = 1 N k ∑ i = 1 n γ i k x i , Σ k ( t + 1 ) = 1 N k ∑ i = 1 n γ i k ( x i − μ k ( t + 1 ) ) ( x i − μ k ( t + 1 ) ) ⊤ . N_k = \sum_{i=1}^n \gamma_{ik},\quad \pi_k^{(t+1)} = \frac{N_k}{n},\quad \mu_k^{(t+1)} = \frac{1}{N_k}\sum_{i=1}^n \gamma_{ik} x_i,\quad \Sigma_k^{(t+1)} = \frac{1}{N_k}\sum_{i=1}^n \gamma_{ik}(x_i - \mu_k^{(t+1)})(x_i - \mu_k^{(t+1)})^\top. Nk=i=1∑nγik,πk(t+1)=nNk,μk(t+1)=Nk1i=1∑nγikxi,Σk(t+1)=Nk1i=1∑nγik(xi−μk(t+1))(xi−μk(t+1))⊤.
通过上述迭代,GMM 能够在数据中自动识别聚类中心与分布形态,并给出每个点所属各聚类的“软分配”概率。
收敛性、优缺点与改进
-
收敛性:EM 保证对数似然单调递增,但仅能收敛到局部极大值或鞍点;全局最优解无法保证。
-
优点:
- 将难以直接优化的积分式对数似然,分解为更易处理的子问题;
- 代码实现简单,E步与M步常有解析解;
- 在诸多统计模型(GMM、隐马尔可夫模型、潜在狄利克雷分配等)中表现出色。
-
缺点:
- 对初值敏感,易陷入次优解;
- 收敛速度较慢(线性收敛),在高维或大规模数据时迭代次数多;
- 需手动指定潜在变量个数(如GMM的聚类数)。
-
常见改进:
- 随机EM(Stochastic EM):在 E 步中对子集数据进行采样,减少计算量;
- 变分EM(Variational EM):引入可参数化的变分分布,适用于后验分布无解析解的情形;
- 加速收敛:结合牛顿法或拟牛顿法(如EM+Fisher信息矩阵)提高 M 步的收敛率;
- 初始化策略:多次随机初始化或使用 K-means 聚类结果作为初始参数,以降低陷入局部最优的风险。
EM算法以其清晰的“E步—M步”框架,成功地将复杂的潜变量极大似然估计问题拆解为一系列易处理的子问题,并在众多模型中得到了广泛应用。尽管存在局部极值与收敛速度缓慢等局限,通过变分推断、加速方法和更合理的初始化策略,EM算法依然在现代大规模机器学习任务中保持竞争力。未来,随着深度学习与概率图模型的不断融合,EM及其变种将继续演进,服务于更为复杂的隐含结构挖掘与不确定性建模。