zehua

菜单

Administrator
发布于 2026-09-09 / 0 阅读
0
0

EM 算法

#课程笔记

期望最大化方法（Expectation-Maximization）

EM算法最经典的两个应用是：

  • 高斯混合模型 GMM 的参数估计: 聚类或拟合数据分布
  • 隐式马尔可夫模型 HMM 的参数估计: 时间序列数据，比如语音识别、天气预测等

EM算法还能应用在一些更细分的应用:

  • 学习固定模型的最佳组合方式
  • 复合 Dirichlet 分布的参数估计： Dirichlet分布通常用于表示多项分布的参数（比如投票、市场份额等），EM可以用来估计它的参数。
  • 分离叠加信号： 用EM算法把混合信号中不同来源的信号拆分出来。

在实际使用EM算法时，可能会遇到收敛速度慢、容易陷入局部最优等问题，因此产生了一些变体，后面会介绍

1. The Expectation-Maximization Method

EM算法是用来处理隐藏变量或不完全数据的一种优化迭代算法，通常用于概率模型中寻找某个参数\theta的最大似然估计，找到使 数据最有可能发生 的模型参数。

最简单例子:

把一天的温度记录下来，用 x \in \mathbb{R}^{24} 表示，这是我们得到的数据，我们要估计的隐藏变量是 季节 \theta（夏季、秋季、冬季、春季），因此我们要根据温度分布 p(x|\theta) 来推断当前的季节\theta 。 通俗来讲，我们已经有了当前季节下的温度数据了，我们默认现实中已经发生的事就是概率最大的事，因此我们反过来找哪一种 \theta 能让这个概率最大。

这个情况还好说，但是如一旦问题复杂化，我们并没有完整的 24 小时温度数据 x，我们只能观察到一天的平均温度 y ，所以我们怎么估计季节 \theta ？这就是一个典型的 仅部分数据可见的问题，就要靠 EM 算法了

工作原理:

EM算法的核心是两步循环：

1. E步（Expectation, 期望步）： 在当前的参数\theta下，推测隐藏变量（或者缺失数据）是什么，比如我们推测一天的完整温度 x 应该是什么样子的

2. M步（Maximization, 最大化步）： 基于刚才的推测结果x，重新调整模型参数\theta，让观测数据y在模型中更“可能”发生。 得到刚刚推测后的 x，可以更好的预测季节 \theta。

这两步不断交替进行，直到模型收敛，也就是\theta变化不大为止。

简单来说，EM算法的目的是找到一个最优的参数 \hat{\theta}，让现有的观测数据（例如平均温度y）最可能发生。

1.1 The EM Algorithm

使用EM算法的条件：

  • 要有已知的部分观测数据 y （一天的平均温度）
  • 希望有完整数据x （24小时中每小时的温度）（缺失）
  • 我们已知:
    • p(x|\theta)：在参数\theta （当前季节）下，完整数据x的概率分布。
    • p(y|\theta)：在参数\theta下，观测数据y的概率分布。

目标： 通过最大化p(y|\theta)（观测数据的似然函数）来找到最优参数\theta，记作\hat{\theta}_{MLE}。

因此我们有如下最大似然估计 MLE 公式

\hat{\theta}_{MLE} = \arg\max_{\theta \in \Omega} p(y|\theta)\tag{1.1}

我们将 p(y|\theta) 最大化来找到最优的 \theta，但是直接最大化p(y|\theta)很复杂，因此变成了最大化对数似然函数：

\hat{\theta}_{MLE} = \arg\max_{\theta \in \Omega} \log p(y|\theta) \tag{1.2}
  • 对数函数是单调递增的，两种公式的解相同
  • 对数可以把原来复杂的连乘式化简为加法式，计算更方便。

但是即使这样，直接求解 (1.2) 的最大值也很难，特别是当x是隐藏变量时。于是，EM算法通过两步迭代解决问题：

  • 猜测隐藏数据x的分布；
  • 用这个分布去优化参数\theta，使得似然值尽量变大。

两步交替执行，直到找到最优的\theta。

EM算法的五个步骤

将EM算法传统的两步（E步和M步）进一步拆解为五个步骤：

Step 1: 初始化参数

  • 设m = 0，初始估计参数\theta^{(m)}。
  • 通常通过随机初始化\theta来开始，或者根据经验值猜测一个初始值。

Step 2: 条件分布计算（E步的一部分）

  • 基于当前参数\theta^{(m)}，计算隐藏数据x的条件分布：
p(x | y, \theta^{(m)}).
  • 这一步的本质是：猜测隐藏数据x可能是什么。

Step 3: 构建Q函数（E步的另一部分）

  • 根据上一步的结果，计算Q函数，即对数似然的期望：
\begin{align*} Q(\theta \mid \theta^{(m)}) &= \int_{\mathcal{X}(y)} \log p(x \mid \theta) p(x \mid y, \theta^{(m)}) \, dx \\ &= \mathbb{E}_{X \mid y, \theta^{(m)}}[\log p(X \mid \theta)] \end{align*}
  • 这是EM算法的核心，Q函数用于量化在当前参数\theta下，完整数据x的对数似然。

Step 4: 最大化Q函数（M步）

  • 通过最大化Q(\theta|\theta^{(m)})找到新的参数估计\theta^{(m+1)}：
\theta^{(m+1)} = \arg\max_{\theta \in \Omega} Q(\theta|\theta^{(m)}).
  • 根据上一步对x的推测，更新参数\theta。

Step 5: 检查收敛条件

  • 判断新的\theta^{(m+1)}是否与之前的\theta^{(m)}足够接近，或者对数似然的变化是否小于某个阈值\epsilon：
\|\theta^{(m+1)} - \theta^{(m)}\| < \epsilon.
  • 如果满足收敛条件，算法终止；否则回到Step 2。

E步： 用当前参数\theta^{(m)}计算Q函数，猜测隐藏数据x

M步： 最大化Q函数，更新\theta^{(m+1)}

即使我们无法保证能找到全局最优解，但是 EM 保证每次迭代的似然值p(y|\theta) 都不会变差。如果p(y|\theta)有多个局部最大值，可能需要多次初始化。

1.2 EM算法的简化变体

在标准EM算法中，我们在E步会计算完整数据x的条件分布，即p(x|y, \theta^{(m)})，然后在M步中用这个条件分布去优化\theta。尽管理论很丰满，但是实际计算很骨感，条件分布涉及积分或者期望的计算。因此为了简化省略复杂数学过程，我们引入简单的变种。

我们不再计算完整的条件分布p(x|y, \theta^{(m)})，而是直接找到一个单一的最有可能的隐藏数据x^{(m)}，然后基于这个点继续优化参数\theta。步骤如下：

  • E-like步（E类步）：

    x^{(m)} = \arg\max_{x \in \mathcal{X}(y)} p(x|y, \theta^{(m)}).
    • 给定当前的参数\theta^{(m)}，在所有可能的x中，找到让p(x|y, \theta^{(m)})最大的那个点。
    • x^{(m)} 就是我们对隐藏数据x的最佳猜测。

  • M-like步（M类步）：

    \theta^{(m+1)} = \arg\max_{\theta \in \Omega} p(x^{(m)}|\theta).
    • 这一步是：在x^{(m)}固定的情况下，寻找能够最大化p(x^{(m)}|\theta)的参数\theta。

简化之处:

  • 标准EM中，E步是计算条件分布的整个期望，而这里的E-like步只是找一个“最优点”x^{(m)}

这种简化方法有两个名字：

  • 点估计EM（Point-estimate EM）： 因为它只寻找隐藏数据x的一个“点估计”x^{(m)}，而不是分布的期望。

  • 分类EM（Classification EM）： 因为它在某些情况下会直接将观测数据分到一个具体的类别，比如后面提到的K-means算法。

K-means聚类

K-means 是点估计EM的一个实际应用场景，我们可以用它来帮助理解这一变体的思想。

K-means 的基本问题

  1. 输入： n个观测数据点y_1, y_2, \dots, y_n，每个点是d维向量，构成矩阵y = [y_1, y_2, \dots, y_n]^T。
  2. 目标： 将这些点分成 k 个簇（cluster），并估计每个簇的“中心点”，也就是k个簇的中心位置\theta = {\theta_1, \dots, \theta_k}。
  3. 隐藏数据： 每个点属于哪个簇是“隐藏的”，我们需要通过算法来确定。

K-means 的步骤

  1. 初始化：

    • 随机初始化 k 个簇中心\theta_1^{(0)}, \dots, \theta_k^{(0)}。
    • 这相当于 EM 算法的“Step 1”。

  2. E-like步：分配簇（E类步）

    • 对于每个点 y_i，分配到最近的簇中心：

      z_i^{(m)} = \arg\min_{j=1,\dots,k} \|y_i - \theta_j^{(m)}\|^2.

    • z_i^{(m)} 表示第 i 个点属于簇 j。

    • 在 EM 的语言中，这相当于寻找隐藏数据x^{(m)}。

  3. M-like步：更新簇中心（M类步）

    • 计算每个簇的新中心点：

      \theta_j^{(m+1)} = \frac{1}{n_j} \sum_{i: z_i^{(m)} = j} y_i.

    • n_j 是簇 j 的点数。

    • 这一步就是通过最大化似然函数来更新参数\theta。

  4. 重复迭代：

    • 重复上述两个步骤，直到簇中心不再显著变化。

K-means 和 EM 的区别

  1. 硬分配 vs. 概率分配：

    • K-means 是“硬分配”，每个点 y_i 只能属于一个簇 j。
    • EM 则是“软分配”，用概率分布p(x|y,\theta)来表示一个点属于每个簇的可能性。

  2. 隐藏数据的处理：

    • K-means 找的是一个明确的“点估计”，类似于简化版的EM。
    • EM 算法计算的是隐藏数据的整个条件分布p(x|y, \theta)，更精确，但计算量大。

标准 EM 算法可能在某些情况下失败，原因是对数似然函数可能出现奇点问题。举个例子：

  • 假如我们用 EM 算法来学习一个高斯混合模型（GMM），其中有 10 个高斯分布（components）。在某些情况下，可能会发生以下问题：
    • 某一个高斯分布只被一个数据点“分配”。
    • 这个分布的协方差（covariance）会被估计为 0，导致计算无法继续（因为协方差矩阵不可逆）。

这是 EM 算法的一个潜在问题：当数据分布稀疏或模型过于复杂时，算法可能陷入不合理的解。

因此我们将先验信息引入 EM 算法，加入关于参数 \theta 的先验信息，让模型能够避免陷入那些极端情况。

先验（Prior）： 是对参数 \theta 的一种预先假设或偏好。

举例：

  • 在高斯混合模型中，先验信息可以是：所有的协方差矩阵都应该有一定的最小值，避免出现“零协方差”。
  • 或者可以假设 \theta（如高斯的均值）遵循一个正态分布，而不是完全自由选择。

用一个概率分布 p(\theta) 表示先验信息。如果我们假设 p(\theta) 是均匀分布，那其实就相当于对 \theta 的约束是平等的。

从最大似然到最大后验

在传统 EM 算法中，我们的目标是找到参数 \theta，使得观测数据 y 的似然函数 p(y|\theta) 最大：

\hat{\theta}_{MLE} = \arg\max_{\theta \in \Omega} \log p(y|\theta)

但当引入先验信息 p(\theta) 后，目标就变成了最大后验估计（MAP）：

\hat{\theta}_{MAP} = \arg\max_{\theta \in \Omega} \log p(\theta|y).

通过贝叶斯公式，我们可以写成：

\log p(\theta|y) = \log p(y|\theta) + \log p(\theta) - \log p(y).

其中：

  • \log p(y|\theta) 是似然函数，来自数据本身。
  • \log p(\theta) 是先验信息，来自我们对参数的“偏好”。
  • \log p(y) 是归一化常数，与优化无关，可以忽略。

因此，最大后验估计的目标变成了：

\hat{\theta}_{MAP} = \arg\max_{\theta \in \Omega} \big( \log p(y|\theta) + \log p(\theta) \big). \tag{1.3}
  • 标准 EM 优化的是 \log p(y|\theta)，即“数据的可能性”。
  • MAP EM 优化的是 \log p(y|\theta) + \log p(\theta)，即“数据可能性 + 先验信息”。

修改后的算法流程：

在标准 EM 中：

  1. E步： 计算 Q 函数，即期望对数似然：

    Q(\theta|\theta^{(m)}) = \mathbb{E}_{X|Y, \theta^{(m)}}[\log p(X|\theta)].

  2. M步： 最大化 Q(\theta|\theta^{(m)})，找到新的参数\theta^{(m+1)}。

在 MAP EM 中：

  • E步：不变。 仍然计算 Q 函数，和标准 EM 相同：

    Q(\theta|\theta^{(m)}) = \mathbb{E}_{X|Y, \theta^{(m)}}[\log p(X|\theta)].

  • M步：修改为最大化后验概率：
    加入先验项 \log p(\theta)，变为：

    \theta^{(m+1)} = \arg\max_{\theta \in \Omega} \big( Q(\theta|\theta^{(m)}) + \log p(\theta) \big).

加入先验的好处：

  • 避免奇点问题：
    • 比如在 GMM 的例子中，通过给 \theta 一个合理的先验（如限制协方差矩阵的最小值），可以避免奇异解。
  • 引入额外信息：
    • 如果我们对参数 \theta 有特定的知识（例如高斯分布的均值或方差范围），可以通过先验表达这种信息，从而得到更合理的结果。

与标准 EM 的关系：

  • 如果先验 p(\theta) 是均匀分布，MAP 就退化为标准 EM，因为 \log p(\theta) 是常数，不会影响优化目标。

1.4 完整数据

完整数据是 EM 算法中的一个核心概念，因为 EM 的目标是通过最大化完整数据的对数似然函数来优化模型参数 \theta。

完整数据 x 的定义需要满足以下两个实用要求：

  1. 在给定完整数据 x 的情况下，优化 p(x|\theta)（以找到参数 \theta）应该是简单的。这确保了 EM 算法的 M 步可以容易计算。
  2. x 通常由观测数据 y 和一些隐藏或缺失的数据构成。

从理论上，完整数据 X 需要满足一个马尔科夫关系：\theta \to X \to Y。

  • 意思是：参数 \theta 和观测数据 Y 之间的关系完全通过完整数据 X 来描述。
  • 用公式表示：p(y|x, \theta) = p(y|x)，也就是说，给定 x 后，y 与 \theta 是条件独立的。

当观测数据 Y 是完整数据 X 的函数时（例如，Y 是从 X 的某些部分生成的），情况会更简单。

  • T(X) 是一个确定性函数： X \to Y 的映射关系是确定的。
  • 马尔科夫关系自动满足： 因为 Y 是 X 的函数，所以 p(y|x, \theta) = p(y|x) 始终成立。

一个例子是：假设 X 是一天24小时的完整温度记录，而 Y 是这些温度的平均值。这里 Y 是通过一个确定函数（取平均）从 X 生成的。

1.4.1 EM 在缺失数据问题中的应用

GMM（高斯混合模型）：

  • 在 GMM 中，Z 是每个数据点属于某个高斯分布的分类标签（即“隐藏变量”）。
  • E步：计算 p(z|y, \theta^{(m)})，即每个数据点属于每个高斯分布的概率。
  • M步：优化高斯分布的均值、协方差和混合权重。

HMM（隐马尔可夫模型）：

  • 在 HMM 中，Z 是隐藏的状态序列。
  • E步：计算给定观测数据 Y 的状态序列分布 p(z|y, \theta^{(m)})。
  • M步：优化状态转移概率、发射概率等参数。

在许多应用中（如 GMM 和 HMM 模型），完整数据 X 是由两部分构成的：

  1. 观测数据 Y： 这是你能够看到的部分。
  2. 隐藏数据或缺失数据 Z： 这是你无法直接观测到的部分。

完整数据可以写成 X = (Y, Z)。这意味着 EM 的任务可以分解为：

  • 推测 Z（隐藏数据）。
  • 基于 Z 和 Y 来优化参数 \theta。

Q 函数的改写

\begin{align} Q(\theta \mid \theta^{(m)}) &= \int_{\mathcal{X}} \log p(x \mid \theta) p(x \mid y, \theta^{(m)}) \, dx \notag \\ &= \int_{\mathcal{X}} \log p(y, z \mid \theta) p(y, z \mid y, \theta^{(m)}) \, dx \notag \\ &= \int_{\mathcal{Z}} \log p(y, z \mid \theta) p(z \mid y, \theta^{(m)}) \, dz \notag \\ &= \mathbb{E}_{Z \mid y, \theta^{(m)}}[\log p(y, Z \mid \theta)]. \tag{1.5} \end{align}

推导过程:

完整数据的对数似然：

Q(\theta|\theta^{(m)}) = \int_{\mathcal{X}} \log p(x|\theta) p(x|y, \theta^{(m)}) dx.
  • 这里 p(x|\theta) 是完整数据的概率分布。
  • p(x|y, \theta^{(m)}) 是当前参数 \theta^{(m)} 下完整数据 x 的条件分布。

将完整数据分解为 (Y, Z)：
完整数据 x 可以写成 (y, z)，其中 y 是观测数据，z 是缺失数据。所以：

Q(\theta|\theta^{(m)}) = \int_{\mathcal{X}} \log p(y, z|\theta) p(y, z|y, \theta^{(m)}) dx.

化简积分域：
在观测数据 y 给定的情况下，完整数据的随机部分只剩下缺失数据 z，所以积分可以只对 z 进行：

Q(\theta|\theta^{(m)}) = \int_{\mathcal{Z}} \log p(y, z|\theta) p(z|y, \theta^{(m)}) dz.

分解联合概率：
p(y, z|\theta) 可以分解为 p(y|z, \theta) \cdot p(z|\theta)：

Q(\theta|\theta^{(m)}) = \int_{\mathcal{Z}} \log \big[ p(y|z, \theta) \cdot p(z|\theta) \big] p(z|y, \theta^{(m)}) dz.

对数展开：
将对数展开为两部分后，得到：

Q(\theta|\theta^{(m)}) = \int_{\mathcal{Z}} \big[ \log p(y|z, \theta) + \log p(z|\theta) \big] p(z|y, \theta^{(m)}) dz.

期望形式：
使用期望的符号表示，最终的 Q 函数为：

Q(\theta|\theta^{(m)}) = \mathbb{E}_{Z|Y, \theta^{(m)}}[\log p(Y, Z|\theta)].

1.4.2 独立同分布样本的EM算法

对于许多常见的应用，例如学习高斯混合模型（GMM）或隐马尔可夫模型（HMM），完整数据X是一组n个独立同分布（i.i.d.）的随机向量，X = [X_1 \ X_2 \ \dots \ X_n]^\top，且第i个观测样本y_i仅是X_i的函数。那么，以下命题对于将Q函数分解为一个求和形式是非常有用的：

命题 1.1

假设对于所有x \in \mathcal{X}^n和\theta \in \Omega，有p(x|\theta) = \prod_{i=1}^n p(x_i|\theta)，并且对于所有i = 1, \dots, n，满足马尔可夫关系\theta \to X_i \to Y_i，即：

p(y_i \mid x, y_1, \dots, y_{i-1}, y_{i+1}, \dots, y_n, \theta) = p(y_i \mid x_i),

如公式(1.6)所示。

因此，Q(\theta|\theta^{(m)})可以被分解为：

Q(\theta|\theta^{(m)}) = \sum_{i=1}^n Q_i(\theta|\theta^{(m)}),

其中，

Q_i(\theta|\theta^{(m)}) = \mathbb{E}_{X_i|y_i, \theta^{(m)}}[\log p(X_i|\theta)], \quad i = 1, \dots, n.

这表明Q函数可以按样本分量的条件期望形式逐项分解。

下面证明:

首先，我们证明，对于给定的\theta，集合\{(X_i, Y_i)\}, i = 1, \dots, n中的元素是相互独立的，即：

p(x, y|\theta) = \prod_{i=1}^n p(x_i, y_i|\theta). \tag{1.7}

这种相互独立性成立是因为：

\begin{aligned} p(x, y|\theta) &= p(y_1|y_2, \dots, y_n, x, \theta) \cdots p(y_n|x, \theta)p(x|\theta) \quad \text{（由链式法则）} \\ &= p(y_1|x_1, \theta) \cdots p(y_n|x_n, \theta)p(x|\theta) \quad \text{（根据公式(1.6)，但保留} \theta \text{在条件中）} \\ &= p(y_1|x_1, \theta) \cdots p(y_n|x_n, \theta) \prod_{i=1}^n p(x_i|\theta) \quad \text{（由} X \text{的独立性假设）} \\ &= \prod_{i=1}^n p(y_i|x_i, \theta)p(x_i|\theta) \\ &= \prod_{i=1}^n p(x_i, y_i|\theta) \end{aligned}

接下来，我们证明对于所有 i = 1, \dots, n，有：

p(x_i|y, \theta) = p(x_i|y_i, \theta).

这是因为：

\begin{aligned} p(x_i|y, \theta) &= \frac{p(x_i, y|\theta)}{p(y|\theta)} \quad \text{（由贝叶斯定理）}, \\ &= \frac{\int_{\mathcal{X}^{n-1}} p(x, y|\theta) dx_1 \dots dx_{i-1} dx_{i+1} \dots dx_n}{\int_{\mathcal{X}^n} p(x, y|\theta) dx}\\ &= \frac{\int_{\mathcal{X}^{n-1}} \prod_{j=1}^n p(x_j, y_j|\theta) dx_1 \dots dx_{i-1} dx_{i+1} \dots dx_n}{\int_{\mathcal{X}^n} \prod_{j=1}^n p(x_j, y_j|\theta) dx_1 \dots dx_n} \quad \text{（由公式(1.7)）} \\ &= \frac{p(x_i, y_i|\theta) \prod_{j=1, j \neq i}^n \int_{\mathcal{X}} p(x_j, y_j|\theta) dx_j}{\prod_{j=1}^n \int_{\mathcal{X}} p(x_j, y_j|\theta) dx_j} \\ &= \frac{p(x_i, y_i|\theta) \prod_{j=1, j \neq i}^n p(y_j|\theta)}{\prod_{j=1}^n p(y_j|\theta)} \\ &= \frac{p(x_i, y_i|\theta)}{p(y_i|\theta)} \\ &= p(x_i|y_i, \theta) \end{aligned}

接下来，对于 Q(\theta|\theta^{(m)}) 的推导，我们有：

\begin{aligned} Q(\theta|\theta^{(m)}) &= \mathbb{E}_{X|y, \theta^{(m)}}[\log p(X|\theta)]\\ &= \mathbb{E}_{X|y, \theta^{(m)}} \left[ \log \prod_{i=1}^n p(X_i|\theta) \right] \quad \text{（由于 $X$ 的独立性假设）}\\ &= \mathbb{E}_{X|y, \theta^{(m)}} \left[ \sum_{i=1}^n \log p(X_i|\theta) \right]\\ &= \sum_{i=1}^n \mathbb{E}_{X_i|y_i, \theta^{(m)}}[\log p(X_i|\theta)] \\ &= \sum_{i=1}^n Q_i(\theta|\theta^{(m)}) \end{aligned}

最后一行的成立是基于公式 (1.8)。

决策树
K-means 聚类实验：三种实现的比较

评论