EM 算法
EM 算法(Expectation-Maximization)是处理含隐变量概率模型参数估计的通用框架,在 motif discovery、聚类和 HMM 训练中有广泛应用。
- 重点理解 E 步和 M 步的交替迭代过程
- EM 算法不保证全局最优,但保证单调改进似然
EM 算法(Expectation-Maximization,期望最大化算法)是一种用于寻找含隐变量(latent variable)概率模型参数最大似然估计的迭代算法。
它的核心思想是:
- 当直接优化似然函数困难时,引入隐变量
- 交替进行 E 步(计算期望)和 M 步(最大化期望)
- 逐步改进参数估计,直到收敛
在生物信息学中,EM 算法的价值体现在:
- Motif discovery:MEME 等工具使用 EM 算法寻找转录因子结合位点
- 聚类:Gaussian Mixture Model 的参数估计
- HMM 训练:Baum-Welch 算法是 EM 算法在 HMM 上的特例
- 缺失数据处理:处理观测数据不完整的情况
EM 算法提供了一个统一的框架,处理那些"如果能知道隐变量,问题就很简单"的场景。
给定观测数据 和隐变量 ,我们希望找到参数 使得对数似然最大化:
由于求和在对数内部,直接优化很困难。
EM 的洞察
Section titled “EM 的洞察”如果我们知道隐变量 ,那么完整数据的对数似然为:
这个函数通常更容易优化。EM 算法的思想是:
- 当不知道 时,用当前参数估计 的分布
- 基于这个分布,计算完整数据似然的期望
- 最大化这个期望得到新参数
- E 步(Expectation)
- 基于当前参数 $ heta^{(t)}$,计算隐变量 $Z$ 的后验分布 $P(Z|X, heta^{(t)})$。
- M 步(Maximization)
- 基于 E 步得到的分布,计算 Q 函数并最大化,得到新参数 $ heta^{(t+1)}$。
- Q 函数
- 完整数据对数似然在隐变量分布下的期望:$Q( heta, heta^{(t)}) = E_{Z|X, heta^{(t)}}[log P(X, Z; heta)]$。
通用 EM 算法
Section titled “通用 EM 算法”算法 1:通用 EM 算法
输入:观测数据 X,初始参数 θ^(0),收敛阈值 ε输出:估计参数 θ
1. t = 02. repeat: // E步:计算Q函数 a. Q(θ, θ^(t)) = E_\{Z|X,θ^(t)\}[log P(X, Z; θ)]
// M步:最大化Q函数 b. θ^(t+1) = argmax_θ Q(θ, θ^(t))
c. t = t + 13. until |ℓ(θ^(t)) - ℓ(θ^(t-1))| < ε4. return θ^(t)时间复杂度:取决于具体模型,通常每次迭代 ,其中 为数据量, 为参数维度
收敛性:保证单调增加似然,但不保证全局最优
在 Motif Discovery 中的应用
Section titled “在 Motif Discovery 中的应用”给定一组序列,寻找一个共同的 motif 模式:
- 每个序列中 motif 的位置是隐变量
- motif 的 PWM(位置权重矩阵)是需要估计的参数
EM for Motif Discovery
Section titled “EM for Motif Discovery”算法 2:EM 算法用于 motif discovery
输入:序列集合 S = \{s_1, s_2, ..., s_n\},motif长度 L输出:motif的PWM M
1. 初始化:随机选择每个序列中motif的起始位置2. repeat: // E步:计算每个序列每个位置是motif起始位置的概率 a. for each 序列 s_i: for each 可能起始位置 j: P(位置 j 是motif | s_i, M) ∝ P(s_i[j:j+L] | M) · P(j)
// M步:基于概率更新PWM b. for each motif位置 k = 1 to L: for each 字符 c ∈ \{A, C, G, T\}: M[k][c] = (1 + Σ_i Σ_j P(位置 j 是motif | s_i, M) · I(s_i[j+k-1] = c)) / (4 + Σ_i Σ_j P(...))3. until 收敛4. return M时间复杂度:每次迭代 ,其中 为平均序列长度
在 Gaussian Mixture Model 中的应用
Section titled “在 Gaussian Mixture Model 中的应用”给定数据点,假设它们来自多个高斯分布的混合:
- 每个数据点来自哪个高斯分布是隐变量
- 每个高斯分布的均值、方差和混合权重是需要估计的参数
EM for GMM
Section titled “EM for GMM”算法 3:EM 算法用于 GMM
输入:数据 X = \{x_1, x_2, ..., x_n\},混合分量数 K输出:每个高斯分量的参数 \{μ_k, Σ_k, π_k\}
1. 初始化:随机初始化 μ_k, Σ_k, π_k2. repeat: // E步:计算每个数据点来自每个分量的概率 a. for each 数据点 x_i: for each 分量 k: γ_ik = P(z_i = k | x_i, θ) = π_k · N(x_i | μ_k, Σ_k) / Σ_j π_j · N(x_i | μ_j, Σ_j)
// M步:更新参数 b. N_k = Σ_i γ_ik c. π_k = N_k / n d. μ_k = Σ_i γ_ik · x_i / N_k e. Σ_k = Σ_i γ_ik · (x_i - μ_k)(x_i - μ_k)^T / N_k3. until 收敛4. return \{μ_k, Σ_k, π_k\}时间复杂度:每次迭代 ,其中 为数据维度
简单的硬币抛掷问题
Section titled “简单的硬币抛掷问题”假设有两枚硬币 A 和 B,但不知道哪枚是哪枚。我们进行了 5 轮实验,每轮抛掷 10 次:
轮次: 1 2 3 4 5结果: H T T T H (5 heads, 5 tails) T H T H T (5 heads, 5 tails) H H T H T (7 heads, 3 tails) H T H T T (4 heads, 6 tails) T H H H T (6 heads, 4 tails)每轮随机选择硬币 A 或 B,但我们不知道选了哪枚。
初始化:假设 ,
第 1 轮迭代:
E 步:计算每轮使用硬币 A 的概率
- 第 1 轮:
- 第 2 轮:
- ...
M 步:基于概率更新 和
复杂度与收敛性
Section titled “复杂度与收敛性”- 单调性:每次迭代保证不降低对数似然
- 收敛到局部最优:可能收敛到局部最优,不保证全局最优
- 收敛速度:通常是线性收敛,在某些情况下可以是二次收敛
初始化敏感性
Section titled “初始化敏感性”EM 算法对初始化敏感:
- 不同的初始参数可能收敛到不同的局部最优
- 实践中常用多次随机初始化,选择似然最高的结果
- Hard EM:在 E 步中选择最可能的隐变量值(而非分布)
- Generalized EM:在 M 步只部分改进 Q 函数(而非完全最大化)
- Stochastic EM:在 E 步中使用采样而非精确计算
在生物信息学中的位置
Section titled “在生物信息学中的位置”EM 算法在生物信息学中的核心应用:
- Motif discovery:MEME 等工具的核心算法
- 聚类:单细胞数据聚类、基因表达模式识别
- HMM 训练:Baum-Welch 算法是 EM 算法的特例
- 基因预测:某些基因预测工具使用 EM 估计参数
与其他算法的关系
Section titled “与其他算法的关系”| 算法 | 与 EM 的关系 | 应用场景 | |------|-----------|---------| | Baum-Welch | EM 在 HMM 上的特例 | HMM 参数估计 | | K-means | Hard EM 在 GMM 上的特例 | 聚类 | | Gibbs sampling | MCMC 方法,可用于近似 E 步 | 复杂模型的近似推断 | | Variational inference | 确定性近似,可用于加速 E 步 | 大规模模型的近似推断 |
- Dempster, A. P., Laird, N. M., & Rubin, D. B. (1977). "Maximum likelihood from incomplete data via the EM algorithm"
- Bishop, C. M. Pattern Recognition and Machine Learning, Chapter 9
- Bailey, T. L., & Elkan, C. (1994). "Fitting a mixture model by expectation maximization to discover motifs in biopolymers"