概率图模型

概率图模型

机器学习最重要的任务是根据已观察到的证据（例如训练样本）对感兴趣的未知变量（例如类别标记）进行估计和推测。

概率模型（probabilistic model）提供了一种描述框架，将学习任务归结为计算变量的概率分布，在概率模型中，利用已知变量推测未知变量的分布称为推断（inference），其核心是如何基于可观测变量推测出未知变量的条件分布：

• 生成式：计算联合分布 $P(Y, R, O)$
• 判别式：计算条件分布 $P(Y, R \mid O)$

其中

• $Y$ 为关心的变量集合（Yield）；
• $O$ 为可观测变量集合（Observation）；
• $R$ 为其他变量集合（Rest）。

给定一组观测变量值，「推断」就是要由 $P(Y, R, O)$ 或 $P(Y, R \mid O)$ 得到条件概率分布 $P(Y \mid O)$。

直接利用概率求和规则消去 $R$ 不可行，时间和空间复杂度为指数级别 $O\left(2^{|Y| + |R|}\right)$。需要一套能简洁紧凑地表达变量间关系的工具。

概率图模型（probabilistic graphical model）是一种用图来表达变量相关关系的概率模型。

图模型提供了一种描述框架：

• 结点：随机变量（集合）
• 边：变量间的概率相关关系

分类：

• 有向图：使用有向无环图表示变量之间的依赖关系
  • 称为「有向图模型」或「贝叶斯网」（Bayesian Network）
• 无向图：使用无向图表示变量间的相关关系
  • 称为「无向图模型」或「马尔可夫网」（Markov Network）

图模型的两种表示

隐马尔可夫模型

隐马尔可夫模型（Hidden Markov Model, HMM）是一种有向图模型，用于描述序列数据的生成过程。

组成：

• 状态变量 $\left\lbrace y_1, \dots, y_n \right\rbrace$，通常假定是隐藏的，不可被观测的
  • 取值范围为 $\mathcal{Y} = \left\lbrace s_1, s_2, \dots, s_N \right\rbrace$，通常是有 $N$ 个可能取值的离散空间。
  • 亦称为「隐变量」（hidden variable）
• 观测变量 $\left\lbrace x_1, \dots, x_n \right\rbrace$，表示第 $i$ 时刻的观测值集合
  • 观测变量可以为离散或连续型，本节只讨论离散型，设为 $\mathcal{X} = \left\lbrace o_1, \dots, o_M \right\rbrace$。

隐马尔科夫模型 $t$ 时刻的状态 $y_t$ 仅依赖于 $t-1$ 时刻的状态 $y_{t-1}$，与其余 $t - 2$ 个状态无关。

于是

$$P(x_1, y_1, \dots, x_n, y_n) = P(y_1) P(x_1 \mid y_1) \prod_{i=2}^{n} P(y_i \mid y_{i-1}) P(x_i \mid y_i)$$

确定一个 HMM 需要三组参数 $\lambda = [\mathbf{A}, \mathbf{B}, \bm{\pi}]$：

• 状态转移概率 $\mathbf{A} = [a_{ij}]_{N \times N}$：模型在各个状态间的转换的概率
  • $a_{ij} = P(y_{y+1} = s_{j} \mid y_t = s_i)$ 表示在任意时刻 $t$，若状态为 $s_i$，则下一时刻状态为 $s_j$ 的概率
• 输出观测概率 $\mathbf{B} = [b_{ij}]_{N \times M}$：模型根据当前状态获得各个观测值的概率
  • $b_{ij} = P(x_t = o_j \mid y_t = s_i)$ 表示在任意时刻 $t$，若状态为 $s_i$，则观测值 $o_{j}$ 被获取的概率
• 初始状态概率 $\bm{\pi} = (\pi_1, \dots, \pi_N)$：模型在初始时刻各个状态出现的概率
  • $\pi_i = P(y_1 = s_i)$ 表示在初始时刻，状态为 $s_i$ 的概率

通过指定状态空间 $\mathcal{Y}$、观测空间 $\mathcal{X}$ 以及上述三组参数，就能确定一个 HMM，通常用其参数 $\lambda$ 来指代。

给定 HMM $\lambda$，它按如下过程产生观测序列 $\left\lbrace x_1, \dots, x_n \right\rbrace$：

1. 设置 $t = 1$，根据初始状态概率 $\bm{\pi}$ 选择初始状态 $y_1$；
2. 根据状态 $y_t$ 和输出观测概率 $\mathbf{B}$ 选择观测变量取值 $x_t$；
3. 根据状态 $y_t$ 和状态转移概率 $\mathbf{A}$ 转移模型状态，即确定下一时刻的状态 $y_{t+1}$；
4. 若 $t < n$，设置 $t = t + 1$，并返回步骤 2；否则结束。

实际应用中，常关注 HMM 的三个基本问题：

1. 给定模型 $\lambda$，如何有效计算其产生观测序列 $\mathbf{x} = \left\lbrace x_1, \dots, x_n \right\rbrace$ 的概率 $P(\mathbf{x} \mid \lambda)$？
   • 如何评估模型与观测序列之间的匹配程度？
   • 应用：根据以往的观测序列预测当前时刻最有可能的观测值 $x_n$
2. 给定模型 $\lambda$ 和观测序列 $\mathbf{x}$，如何找到与此观测序列最匹配的状态序列 $\mathbf{y} = \left\lbrace y_1, \dots, y_n \right\rbrace$？
   • 如何根据观测序列推断出隐藏的模型状态？
   • 应用：语音识别
3. 给定观测序列 $\mathbf{x}$，如何调整模型参数 $\lambda$ 使得该序列出现的概率 $P(\mathbf{x} \mid \lambda)$ 最大？
   • 如何训练模型使其能最好地描述观测数据？
   • 应用：根据数据学习参数（模型训练）

马尔可夫随机场

马尔可夫随机场（Markov Random Field, MRF）是一种无向图模型，用于描述多个随机变量之间的相关关系。

图模型表示：

• 结点表示变量（集），边表示依赖关系。
• 有一组势函数（potential function），亦称为因子（factor），是定义在变量子集上的非负函数，主要用于定义概率分布函数。

对于图中结点的一个子集，若其中任意两结点间都有边连接，则称该结点子集为团（clique）。若在一个团中加入任意一个结点都不再是团，则称该团为极大团（maximal clique）。

MRF 使用基于极大团的因子，多个变量之间的联合概率分布可基于团分解为多个因子的乘积，每个因子仅与一个团相关。

具体来说，对于 $n$ 个变量 $\mathbf{x} = \left\lbrace x_1, \dots, x_n \right\rbrace$，所有团构成的集合为 $\mathcal{C}$，与团 $Q \in \mathcal{C}$ 对应的变量集合记为 $\mathbf{x}_Q$，则联合概率 $P(\mathbf{x})$ 定义为

$$P(\mathbf{x}) = \dfrac{1}{Z} \prod_{Q \in \mathcal{C}} \psi_Q(\mathbf{x}_Q)$$

其中

• $\psi_Q$ 为与团 $Q$ 对应的势函数，用于对团 $Q$ 中的变量关系进行建模;
• $Z = \sum_{\mathbf{x}}\prod_{Q \in \mathcal{C}}\psi_{Q}(\mathbf{x}_Q)$ 为规范化因子。
  • 实际应用中 $Z$ 的精确计算通常很困难，但许多任务并不需要获得其精确值。

若变量个数较多，团的数目会很多，上式会有很多乘积项，会给计算带来负担。

注意到若团 $Q$ 不是极大团，则其必被一个极大团 $Q^{*}$ 所包含，即 $\mathbf{x}_{Q} \subseteq \mathbf{x}_{Q^{*}}$，这意味着变量 $\mathbf{x}_{Q}$ 之间的关系不仅体现在 $\psi_Q$ 中，还体现在 $\psi_{Q^{*}}$ 中。于是联合概率可基于极大团来定义：

$$P(\mathbf{x}) = \dfrac{1}{Z^{*}} \prod_{Q \in \mathcal{C}^{*}} \psi_Q(\mathbf{x}_Q)$$

其中

• $\mathcal{C}^{*}$ 为所有极大团的集合；
• $Z^{*} = \sum_{\mathbf{x}}\prod_{Q \in \mathcal{C}^{*}}\psi_{Q}(\mathbf{x}_Q)$ 为规范化因子。

MRF 可以得到「条件独立性」。如下图所示，若从结点集 $A$ 中的结点到结点集 $B$ 中的结点，都必须经过结点集 $C$ 中的结点，则称结点集 $A, B$ 被结点集 $C$ 分离，$C$ 为「分离集」（separating set）。

对 MRF 有：

• 「全局马尔可夫性」（global Markov property）：给定两个变量子集的分离集，则这两个变量子集条件独立。

在上图中，若三个结点集对应的变量集分别为 $\mathbf{x}_A, \mathbf{x}_B, \mathbf{x}_C$，则 $\mathbf{x}_A, \mathbf{x}_B$ 在给定 $\mathbf{x}_C$ 下条件独立，记为 $\mathbf{x}_A \perp \mathbf{x}_B \mid \mathbf{x}_C$。

证明

为了便于讨论，简化为单变量 $x_A, x_B, x_C$，联合概率为

$$P(x_A, x_B, x_C) = \dfrac{1}{Z} \psi_{AC}(x_A, x_C) \psi_{BC}(x_B, x_C)$$

条件概率定义可得

$$\begin{aligned} P(x_A, x_B \mid x_C) &= \dfrac{P(x_A, x_B, x_C)}{P(x_C)} \\ &= \dfrac{P(x_A, x_B, x_C)}{\sum_{x_{A'}}\sum_{x_{B'}}P(x_{A'}, x_{B'}, x_C)} \\ &= \dfrac{\dfrac{1}{Z} \psi_{AC}(x_A, x_C) \psi_{BC}(x_B, x_C)}{\sum_{x_{A'}}\sum_{x_{B'}}\dfrac{1}{Z} \psi_{AC}(x_{A'}, x_C) \psi_{BC}(x_{B'}, x_C)} \\ &= \dfrac{\psi_{AC}(x_A, x_C)}{\sum_{x_{A'}}\psi_{AC}(x_{A'}, x_C)} \dfrac{\psi_{BC}(x_B, x_C)}{\sum_{x_{B'}}\psi_{BC}(x_{B'}, x_C)} \end{aligned}$$

而

$$\begin{aligned} P(x_A \mid x_C) &= \dfrac{P(x_A, x_C)}{P(x_C)} \\ &= \dfrac{\sum_{x_{B'}}P(x_A, x_{B'}, x_C)}{\sum_{x_{A'}}\sum_{x_{B'}}P(x_{A'}, x_{B'}, x_C)} \\ &= \dfrac{\sum_{x_{B'}}\dfrac{1}{Z} \psi_{AC}(x_A, x_C) \psi_{BC}(x_{B'}, x_C)}{\sum_{x_{A'}}\sum_{x_{B'}}\dfrac{1}{Z} \psi_{AC}(x_{A'}, x_C) \psi_{BC}(x_{B'}, x_C)} \\ &= \dfrac{\psi_{AC}(x_A, x_C)}{\sum_{x_{A'}}\psi_{AC}(x_{A'}, x_C)} \end{aligned}$$

同理可得 $P(x_B \mid x_C)$，于是

$$P(x_A, x_B \mid x_C) = P(x_A \mid x_C) P(x_B \mid x_C)$$

即 $x_A \perp x_B \mid x_C$。

由全局马尔可夫性可得两个推论：

• 「局部马尔可夫性」（local Markov property）：给定某变量的邻接变量，则该变量条件独立于其他变量。
  • 令 $V$ 为图的结点集，$n(v)$ 为结点 $v$ 在图上的邻接结点，$n^{*}(v) = n(v) \cup \left\lbrace v \right\rbrace$，有 $\mathbf{x}_v \perp \mathbf{x}_{V \setminus n^{*}(v)} \mid \mathbf{x}_{n(v)}$。
• 「成对马尔可夫性」（pairwise Markov property）：给定所有其他变量，两个非邻接变量条件独立。
  • 令图的结点集和边集分别为 $V, E$，对图中的两个结点 $u, v$，若 $\left\lbrace u, v \right\rbrace \notin E$，则有 $\mathbf{x}_u \perp \mathbf{x}_v \mid \mathbf{x}_{V \setminus \left\lbrace u, v \right\rbrace}$。

势函数 $\psi_Q(\mathbf{x}_Q)$ 的作用是定量刻画变量集 $\mathbf{x}_Q$ 中变量的相关关系，应为非负函数，且在所偏好的变量取值上有较大的函数值。例如：

为了满足非负性，指数函数常被用于定义势函数，即

$$\psi_Q(\mathbf{x}_Q) = \exp \left( -H_Q(\mathbf{x}_Q) \right)$$

其中 $H_Q(\mathbf{x}_Q)$ 是一个定义在变量 $\mathbf{x}_Q$ 上的实值函数，常见形式为

$$H_Q(\mathbf{x}_Q) = \sum_{u, v \in Q,\, u \ne v}\alpha_{uv}x_u x_v + \sum_{v \in Q}\beta_v x_v$$

其中 $\alpha_{uv}, \beta_v$ 是参数。上式第二项仅考虑单结点，第一项则考虑每一对结点的关系。

条件随机场

条件随机场（Conditional Random Field, CRF）是一种判别式无向图模型，可看作是给定观测值的 MRF。CRF 对多个变量给定相应观测值后的条件概率进行建模。

令 $\mathbf{x} = \left\lbrace x_1, \dots, x_n \right\rbrace$ 为观测序列，$\mathbf{y} = \left\lbrace y_1, \dots, y_n \right\rbrace$ 为与之相应的标记序列，则 CRF 的目标是构建条件概率模型 $P(\mathbf{y} \mid \mathbf{x})$。

其中 $\mathbf{y}$ 可以是结构型变量，即其分量之间具有某种相关性。

• 自然语言处理的词性标注任务中，观测数据为语句（单词序列），标记为相应的词性序列，具有线性序列结构
• 在语法分析任务中，输出标记是语法树，具有树形结构

令 $G = \left\langle V, E \right\rangle$ 表示结点与标记变量 $\mathbf{y}$ 中的元素一一对应的无向图，$y_v$ 表示与结点 $v$ 对应的标记变量，$n(v)$ 表示结点 $v$ 的邻接结点，若图 $G$ 的每个变量 $y_v$ 都满足马尔可夫性，即 $P(y_v \mid \mathbf{x}, \mathbf{y}_{V \setminus \left\lbrace v \right\rbrace}) = P(y_v \mid \mathbf{x}, \mathbf{y}_{n(v)})$，则 $(\mathbf{y}, \mathbf{x})$ 构成一个 CRF。

理论上 $G$ 可具有任意结构，现实应用中最常用的是如下图所示的链式结构，即「链式条件随机场」（chain-structured CRF）。

类似 MRF，CRF 也使用势函数和图结构的团来定义条件概率 $P(\mathbf{y} \mid \mathbf{x})$ 。

给定观测序列 $\mathbf{x}$，上图所示的链式 CRF 主要包含两种关于标记变量的团，即单个标记变量 $\left\lbrace y_i \right\rbrace$ 以及相邻的标记变量 $\left\lbrace y_{i-1}, y_i \right\rbrace$。

选择合适的势函数，可得到形如 $P(\mathbf{x}) = \dfrac{1}{Z} \displaystyle \prod_{Q \in \mathcal{C}}\psi_Q(\mathbf{x}_Q)$ 的条件概率定义。在条件概率场中，通过选用指数势函数并引入特征函数（feature function），条件概率被定义为

$$P(\mathbf{y} \mid \mathbf{x}) = \dfrac{1}{Z} \exp\left( \sum_{j} \sum_{i=1}^{n-1}\lambda_{j} t_{j}(y_{i+1}, y_i, \mathbf{x}, i) + \sum_{k}\sum_{i=1}^{n}\mu_{k}s_{k}(y_i, \mathbf{x}, i) \right)$$

其中

• $t_{j}(y_{i+1}, y_i, \mathbf{x}, i)$ 是定义在观测序列的两个相邻标记位置上的转移特征函数（transition feature function），用于刻画相邻标记变量之间的相关关系以及观测序列对它们的影响；
• $s_{k}(y_i, \mathbf{x}, i)$ 是定义在观测序列的标记位置 $i$ 上的状态特征函数（state feature function），用于刻画观测序列对标记变量的影响；
• $\lambda_{j}, \mu_{k}$ 为参数；
• $Z$ 为规范化因子。

要使用 CRF，还需定义合适的特征函数，特征函数通常是实值函数，以刻画数据的一些很可能成立或期望成立的经验特性。例如：

条件随机场和马尔可夫随机场均使用团上的势函数定义概率，两者在形式上没有显著区别；但条件随机场处理的是条件概率，而马尔可夫随机场处理的是联合概率。

图模型推断

基于概率图模型定义的分布，能对目标变量的边际分布（marginal distribution）或某些可观测变量为条件的条件分布进行推断。

对概率图模型，还需确定具体分布的参数，称为「参数估计」或「参数学习」问题，通常使用极大似然估计或最大后验概率估计求解。但若将参数视为待推测的变量，则参数估计过程和推断十分相似，可以「吸收」到推断问题中

假设图模型对应的变量集 $\mathbf{x} = \left\lbrace x_1, \dots, x_n \right\rbrace$ 能分成两个不相交的变量集 $\mathbf{x}_F, \mathbf{x}_E$，推断问题的目标就是计算「边际概率」 $P(\mathbf{x}_F)$ 或「条件概率」 $P(\mathbf{x}_F \mid \mathbf{x}_E)$。

条件概率有

$$P(\mathbf{x}_F \mid \mathbf{x}_E) = \dfrac{P(\mathbf{x}_F, \mathbf{x}_E)}{P(\mathbf{x}_E)} = \dfrac{P(\mathbf{x}_F, \mathbf{x}_E)}{\sum_F P(\mathbf{x}_F, \mathbf{x}_E)}$$

分子为联合分布，分母为边际分布。

联合概率 $P(\mathbf{x}_F, \mathbf{x}_E)$ 可基于图模型获得，推断问题的关键在于高效计算边际分布，即

$$P(\mathbf{x}_E) = \sum_F P(\mathbf{x}_F, \mathbf{x}_E)$$

分类：

• 精确推断方法
  • 计算出目标变量的边际分布或条件分布的精确值
  • 一般情况下，该类方法的计算复杂度随极大团规模增长呈指数增长，适用范围有限
• 近似推断方法
  • 在较低的时间复杂度下获得原问题的近似解
  • 在实际问题中更常用

精确推断

精确推断实质是一种动态规划算法, 利用图模型所描述的条件独立性来削减计算目标概率值所需的计算量。

如下图，计算边缘概率 $P(x_5)$：

按 $\left\lbrace x_1, x_2, x_4, x_3 \right\rbrace$ 的顺序计算加法，有

$$\begin{aligned} P(x_5) &= \sum_{x_4} \sum_{x_3} \sum_{x_2} \sum_{x_1} P(x_1, x_2, x_3, x_4, x_5) \\ &= \sum_{x_4} \sum_{x_3} \sum_{x_2} \sum_{x_1} P(x_1) P(x_2 \mid x_1) P(x_3 \mid x_2) P(x_4 \mid x_3) P(x_5 \mid x_3) \\ &= \sum_{x_3} P(x_5 \mid x_3) \sum_{x_4} P(x_4 \mid x_3) \sum_{x_2} P(x_3 \mid x_2) \sum_{x_1} P(x_2 \mid x_1) P(x_1) \\ &= \sum_{x_3} P(x_5 \mid x_3) \sum_{x_4} P(x_4 \mid x_3) \sum_{x_2} P(x_3 \mid x_2) m_{12}(x_2) \end{aligned}$$

其中 $m_{ij}(x_{j})$ 为求加过程的中间结果，下标 $i$ 表示此项是对 $x_i$ 求加的结果，下标 $j$ 表示此项中的剩下的变量。$m_{ij}(x_{j})$ 是关于 $x_{j}$ 的函数。

不断执行可得

$$\begin{aligned} P(x_5) &= \sum_{x_3} P(x_5 \mid x_3) \sum_{x_4} P(x_4 \mid x_3) m_{23}(x_3)\\ &= \sum_{x_3} P(x_5 \mid x_3) m_{23}(x_3) \sum_{x_4} P(x_4 \mid x_3) \\ &= \sum_{x_3} P(x_5 \mid x_3) m_{23}(x_3) m_{43}(x_3)\\ &= m_{35}(x_5) \end{aligned}$$

上述方法对无向图模型同样适用，对上面相应的无向图，有

$$P(x_1, \dots, x_5) = \dfrac{1}{Z} \psi_{12}(x_1, x_2) \psi_{23}(x_2, x_3) \psi_{34}(x_3, x_4) \psi_{35}(x_3, x_5)$$

边际分布可这样计算：

$$\begin{aligned} P(x_5) &= \dfrac{1}{Z} \sum_{x_3} \psi_{35}(x_3, x_5) \sum_{x_4} \psi_{34}(x_3, x_4) \sum_{x_2} \psi_{23}(x_2, x_3) \sum_{x_1} \psi_{12}(x_1, x_2) \\ &= \dfrac{1}{Z} \sum_{x_3} \psi_{35}(x_3, x_5) \sum_{x_4} \psi_{34}(x_3, x_4) \sum_{x_2} \psi_{23}(x_2, x_3) \sum_{x_1} m_{12}(x_2)\\ &= \dots\\ &= \dfrac{1}{Z} m_{35}(x_5) \end{aligned}$$

显然，通过利用乘法对加法的分配律，变量消去法把多个变量的积的求和问题，转化为对部分变量交替进行求积与求和的问题。这种转化使得每次的求和与求积运算限制在局部，仅与部分变量有关，从而简化了计算。

变量消去法有一个明显的缺点：若需计算多个边际分布，重复使用变量消去法将会造成大量的冗余计算。

信念传播（Belief Propagation）算法将变量消去法中的求和操作看作一个信息传递过程，较好地解决了求解多个边际分布时的重复计算问题。

具体来说，变量消去法通过求和操作

$$m_{ij}(x_{j}) = \sum_{x_i} \psi(x_i, x_{j}) = \prod_{k \in n(i) \setminus j} m_{ki}(x_i)$$

消去变量 $x_i$。而在信念传播算法中，这个操作被看作从 $x_i$ 向 $x_{j}$ 传递了一个消息 $m_{ij}(x_{j})$。

这样上面所描述的变量消去过程就能描述为下面右图所示的消息传递过程：

每次传递消息仅与变量 $x_i$ 及其邻接结点直接相关。换言之，消息传递相关的计算被限制在图的局部进行。

在信念传播算法中，一个结点仅在接收到来自其他所有结点的消息后才能向另一个结点发送消息，且结点的边际分布正比于它所接收的消息的乘积，即

$$P(x_i) \propto \prod_{k \in n(i)} m_{ki}(x_i)$$

例如在上图中，结点 $x_3$ 要向 $x_5$ 发消息，必须事先收到来自结点 $x_2, x_4$ 的消息，且传递到 $x_5$ 的消息 $m_{35}(x_5)$ 恰为概率 $P(x_5)$。

若图结构中没有环，则信念传播算法经过两个步骤即可完成所有消息传递，进而能计算所有变量上的边际分布：

1. 指定一个根结点，从所有叶结点开始向根结点传递消息，直到根结点收到所有邻接结点的消息；
2. 从根结点开始向叶结点传递消息，直到所有叶结点均收到消息。

例如在上图中，令 $x_1$ 为根节点，则 $x_4, x_5$ 为叶节点。消息传递过程如下图所示：

此时图的每条边上都有方向不同的两条消息，基于这些消息和 $P(x_i) \propto \prod_{k \in n(i)} m_{ki}(x_i)$ 即可获得所有变量的边际概率。

近似推断

精确推断方法需要很大的计算开销，因此在现实应用中近似推断方法更为常用。

近似推断方法大致可以分为两类：

• 采样法（sampling）：通过使用随机化方法完成近似，如 MCMC 采样
• 变分推断（variational inference)：使用确定性近似完成近似推断

MCMC 采样法

很多任务中，我们关心的并非概率分布本身，而是基于概率分布的期望，并且还能基于期望进一步作出决策。若直接计算或逼近这个期望比推断概率分布更容易，则直接操作无疑将使推断问题的求解更为高效。

采样法基于这个思路，假定目标是计算函数 $f(x)$ 在概率密度函数 $p(x)$ 下的期望

$$\mathbb{E}_p[f] = \int f(x) p (x) \d x$$

则可根据 $p(x)$ 抽取一组样本 $\left\lbrace x_1, \dots, x_N \right\rbrace$，然后计算 $f(x)$ 在这些样本上的均值

$$\hat{f} = \dfrac{1}{N} \sum_{i=1}^{N} f(x_i)$$

以近似目标期望 $\mathbb{E}[f]$。

若样本 $\left\lbrace x_1, \dots, x_N \right\rbrace$ 独立，基于大数定理，这种通过大量采样的方法就能获得较高的近似精度。

因此问题的关键在于「如何采样」。对概率图模型来说，就是如何高效地基于图模型所描述的概率分布来获取样本。

概率图模型中最常用的采样技术是马尔可夫链蒙特卡罗（Markov Chain Monte Carlo, MCMC）方法：

给定连续变量 $x \in X$ 的概率密度函数 $p(x)$，$x$ 在 $A$ 中的概率可计算为

$$P(A) = \int_A p(x) \d x$$

若有函数 $f\colon X \to \R$，则可计算 $f(x)$ 的期望

$$p(f) = \mathbb{E}_p [f(X)] = \int_x f(x) p(x) \d x$$

若 $x$ 不是单变量，而是一个个高维多元变量 $\mathbf{x}$，且服从一个非常复杂的分布，对上式求积分通常很困难。

为此 MCMC 先构造出服从 $p$ 分布的独立同分布随机变量 $\mathbf{x}_1, \dots, \mathbf{x}_N$，然后再得到上式的无偏估计

$$\tilde{p}(f) = \dfrac{1}{N} \sum_{i=1}^{N} f(\mathbf{x}_i)$$

然而当概率密度函数 $p(\mathbf{x})$ 很复杂时，构造服从 $p$ 分布的独立同分布样本也很困难。MCMC 方法的关键就在于构造「平稳分布为 $p$ 的马尔可夫链」来产生样本：

若马尔可夫链运行时间足够长（即收敛到平稳状态），则此时产出的样本 $\mathbf{x}$ 近似服从分布 $p$。

如何判断马尔可夫链达到平稳状态呢？假定平稳马尔可夫链 $T$ 的状态转移概率（即从状态 $\mathbf{x}$ 转移到状态 $\mathbf{x}'$ 的概率）为 $T(\mathbf{x}' \mid \mathbf{x})$，$t$ 时刻状态的分布为 $p(\mathbf{x}^t)$，则若在某个时刻马尔可夫链满足「平稳条件」

$$p(\mathbf{x}^t) T(\mathbf{x}^{t-1} \mid \mathbf{x}^t) = p(\mathbf{x}^{t-1}) T(\mathbf{x}^t \mid \mathbf{x}^{t-1})$$

则 $p(\mathbf{x})$ 是该马尔可夫链的平稳分布，且马尔可夫链在满足该条件时已收敛到平稳状态。

MCMC 方法先设法构造一条马尔可夫链，使其收敛到平稳分布恰为待估计参数的后验分布，然后通过这条马尔可夫链来产生符合后验分布的样本，并基于这些样本来进行估计。

因此，这里的马尔可夫链转移概率的构造至关重要，不同的构造方法将产生不同的 MCMC 算法。

Metropolis-Hastings(MH) 算法是 MCMC 的重要代表。它基于「拒绝采样」（reject sampling）来逼近平稳分布 $p$：

算法每次根据上一轮采样结果 $\mathbf{x}^{t-1}$ 来采样获得候选状态样本 $\mathbf{x}^{*}$，但这个候选样本会以一定概率被「拒绝」。

假定从状态 $\mathbf{x}^{t-1}$ 到状态 $\mathbf{x}^{*}$ 的转移概率为 $Q(\mathbf{x}^{*} \mid \mathbf{x}^{t-1}) A (\mathbf{x}^{*} \mid \mathbf{x}^{t-1})$，其中 $Q(\mathbf{x}^{*} \mid \mathbf{x}^{t-1})$ 是用户给定的先验概率，$A(\mathbf{x}^{*} \mid \mathbf{x}^{t-1})$ 是 $\mathbf{x}^{*}$ 被接受的概率。

若 $\mathbf{x}^{*}$ 最终收敛到平稳状态，根据上面的平稳条件，有

$$p(\mathbf{x}^{t-1}) Q(\mathbf{x}^{*} \mid \mathbf{x}^{t-1}) A (\mathbf{x}^{*} \mid \mathbf{x}^{t-1}) = p(\mathbf{x}^{*}) Q(\mathbf{x}^{t-1} \mid \mathbf{x}^{*}) A (\mathbf{x}^{t-1} \mid \mathbf{x}^{*})$$

于是为了达到平稳状态，只需将接受率设置为

$$A (\mathbf{x}^{*} \mid \mathbf{x}^{t-1}) = \min \left\lbrace 1, \dfrac{p(\mathbf{x}^{*}) Q(\mathbf{x}^{t-1} \mid \mathbf{x}^{*})}{p(\mathbf{x}^{t-1}) Q(\mathbf{x}^{*} \mid \mathbf{x}^{t-1})} \right\rbrace$$

吉布斯采样（Gibbs sampling）有时被视为 MH 算法的特例，它也使用马尔可夫链获取样本，该马尔可夫链的平稳分布也是采样的目标分布 $p(\mathbf{x})$。

具体来说，假定 $\mathbf{x} = \left\lbrace x_1, \dots, x_N \right\rbrace$，目标分布为 $p(\mathbf{x})$，在初始化 $\mathbf{x}$ 的取值后，通过循环执行以下步骤完成采样：

1. 随机或以某个次序选取某变量 $x_i$；
2. 根据 $\mathbf{x}$ 中除 $x_i$ 外的现有取值，计算条件概率 $p(x_i \mid \mathbf{x}_{\bar{\imath}})$，其中 $\mathbf{x}_{\bar{\imath}} = \left\lbrace x_1, \dots, x_{i-1}, x_{i+1}, \dots, x_N \right\rbrace$；
3. 根据 $p(x_i \mid \mathbf{x}_{\bar{\imath}})$ 对样本 $x_i$ 采样，用采样值代替原值。

变分推断

变分推断通过使用已知简单分布来逼近需推断的复杂分布，并通过限制近似分布的类型，从而得到一种局部最优、但具有确定解的近似后验分布。

下面是概率图模型的一种简洁的表示方法——盘式记法（plate notation）：

• 相互独立的、由相同机制生成的多个变量（独立同分布）被放在一个方框（盘）内，并在方框中标出类似变量重复出现的个数 $N$
• 方框可以嵌套
• 通常用阴影标注出已知的、能观察到的变量

上图表示 $N$ 个变量 $\left\lbrace x_1, \dots, x_N \right\rbrace$ 均依赖于其他变量 $\mathbf{z}$。

所有能观察到的变量 $x$ 的联合分布的概率密度函数为

$$p(\mathbf{x} \mid \Theta) = \prod_{i=1}^{N} \sum_{\mathbf{z}} p(x_i, \mathbf{z} \mid \Theta)$$

所对应的对数似然函数为

$$\ln p(\mathbf{x} \mid \Theta) = \sum_{i=1}^{N} \ln \left[ \sum_{\mathbf{z}} p(x_i, \mathbf{z} \mid \Theta) \right]$$

其中 $\mathbf{x} = \left\lbrace x_1, \dots, x_N \right\rbrace$，而 $\Theta$ 是 $\mathbf{x}, \mathbf{z}$ 服从的分布参数。

推断和学习任务主要是由观察变量 $\mathbf{x}$ 来估计隐变量 $\mathbf{z}$ 和分布参数 $\Theta$，即求解 $p(\mathbf{z} \mid \mathbf{x}, \Theta)$ 和 $\Theta$。

概率模型的参数估计通常以最大化对数似然函数为手。对对数似然函数使用 EM 算法：

• E 步：根据 $t$ 时刻的参数 $\Theta^t$ 对 $p(\mathbf{z} \mid \mathbf{x}, \Theta^t)$ 进行推断，并计算联合似然函数 $p(\mathbf{x}, \mathbf{z} \mid \Theta)$；
• M 步：基于 E 步的结果进行最大化寻优，即对关于变量 $\Theta$ 的函数 $\mathcal{Q}(\Theta; \Theta^t)$ 进行最大化，从而求取

$$\begin{aligned} \Theta^{t+1} &= \argmax_{\Theta} \mathcal{Q}(\Theta; \Theta^t)\\ &= \argmax_{\Theta} \sum_{\mathbf{z}} p(\mathbf{z} \mid \mathbf{x}, \Theta^t) \ln p(\mathbf{x}, \mathbf{z} \mid \Theta) \end{aligned}$$

上式中的 $\mathcal{Q}(\Theta; \Theta^t)$ 实际上是对数联合似然函数 $\ln p(\mathbf{x}, \mathbf{z} \mid \Theta)$ 在分布 $p(\mathbf{z} \mid \mathbf{x}, \Theta^t)$ 下的期望，当这个分布与 $\mathbf{z}$ 的真实后验分布相等时，$\mathcal{Q}(\Theta; \Theta^t)$ 近似于对数似然函数。于是 EM 算法最终可获得稳定的参数 $\Theta$，而隐变量 $\mathbf{z}$ 的分布也能通过该参数获得。

需要注意的是，$p(\mathbf{z} \mid \mathbf{x}, \Theta^t)$ 未必是隐变量 $\mathbf{z}$ 服从的真实分布，而只是一个近似分布。若将这个近似分布用 $q(\mathbf{z})$ 表示，则有

$$\ln p(\mathbf{x}) = \mathcal{L}(q) + \operatorname{KL}(q \parallel p)$$

其中

$$\begin{aligned} \mathcal{L}(q) &= \int q(\mathbf{z}) \ln \dfrac{p(\mathbf{x}, \mathbf{z})}{q(\mathbf{z})}\d \mathbf{z}\\ \operatorname{KL}(q \parallel p) &= - \int q(\mathbf{z}) \ln \dfrac{p(\mathbf{z} \mid \mathbf{x})}{q(\mathbf{z})} \d \mathbf{z} \end{aligned}$$

过程

$$\begin{aligned} \ln p(\mathbf{x}) &= \ln p(\mathbf{x}, \mathbf{z}) - \ln p(\mathbf{z} \mid \mathbf{x})\\ &= \ln \dfrac{p(\mathbf{x}, \mathbf{z})}{q(\mathbf{z})} - \ln \dfrac{p(\mathbf{z} \mid \mathbf{x})}{q(\mathbf{z})} \end{aligned}$$

对上式两边取 $q(\mathbf{z})$ 的期望，有

$$\begin{aligned} \int q(\mathbf{z}) \ln p(\mathbf{x}) \d \mathbf{z} &= \int q(\mathbf{z}) \ln \dfrac{p(\mathbf{x}, \mathbf{z})}{q(\mathbf{z})} \d \mathbf{z} - \int q(\mathbf{z}) \ln \dfrac{p(\mathbf{z} \mid \mathbf{x})}{q(\mathbf{z})} \d \mathbf{z}\\ &= \mathcal{L}(q) - \operatorname{KL}(q \parallel p) \end{aligned}$$

在这个式子中，$\mathcal{L}(q)$ 构成了下界，KL 散度 $\operatorname{KL}(q \parallel q)$ 衡量了近似分布 $q(\mathbf{z})$ 与真实后验分布 $p(\mathbf{z} \mid \mathbf{x})$ 之间的差异。

在现实任务中，E 步对 $p(\mathbf{z} \mid \mathbf{x}, \Theta^t)$ 的推断很可能因 $\mathbf{z}$ 模型复杂而难以进行，此时可借助变分推断：

通常假设 $\mathbf{z}$ 服从分布

$$q(\mathbf{z}) = \prod_{i=1}^{M} q_i (\mathbf{z}_i)$$

即假设复杂的多变量 $\mathbf{z}$ 可拆解为一系列相互独立的多变量 $\mathbf{z}_i$，且可令 $q_i$ 分布相对简单，或有很好的结构，例如假设 $q_i$ 为指数族（exponential family）分布，此时有（过程略）

$$\begin{aligned} \mathcal{L}(q) &= \int q_{j} \ln \tilde{q}(\mathbf{x}, \mathbf{z}_{j}) \d \mathbf{z}_{j} - \int q_{j} \ln q_{j} \d \mathbf{z}_{j} + \text{const} \end{aligned}$$

其中

$$\begin{aligned} \ln \tilde{q}(\mathbf{x}, \mathbf{z}_{j}) &= \mathbb{E}_{i \ne j} \left[ \ln p(\mathbf{x}, \mathbf{z}) \right] + \text{const}\\ \mathbb{E}_{i \ne j}\left[ \ln p(\mathbf{x}, \mathbf{z}) \right] &= \int \ln p(\mathbf{x}, \mathbf{z}) \prod_{i \ne j} q_i(\mathbf{z}_i) \d \mathbf{z}_i \end{aligned}$$

要得到 $q_{j}$，可固定 $q_{i\ne j}$，再对 $\mathcal{L}(q)$ 最大化，可发现 $\mathcal{L}(q) = - \operatorname{KL}(q_{j} \parallel \tilde{p}(\mathbf{x}, \mathbf{z}_{j}))$，即 $q_{j} = \tilde{p}(\mathbf{x}, \mathbf{z}_{j})$ 时 $\mathcal{L}(q)$ 最大。

于是可知变量子集 $\mathbf{z}_{j}$ 所服从的最优分布 $q_{j}^{*}$ 应满足

$$\ln q_{j}^{*}(\mathbf{z}_{j}) = \mathbb{E}_{i \ne j} \left[ \ln p(\mathbf{x}, \mathbf{z}) \right] + \text{const}$$

即

$$q_{j}^{*}(\mathbf{z}_{j}) = \dfrac{\exp\left( \mathbb{E}_{i \ne j}\left[ \ln p(\mathbf{x}, \mathbf{z}) \right] \right) }{\int \exp\left( \mathbb{E}_{i \ne j}\left[ \ln p(\mathbf{x}, \mathbf{z}) \right] \right) \d \mathbf{z}_{j}}$$

因此通过适当分割变量子集 $\mathbf{z}_{j}$ 并选择 $q_i$ 服从的分布，$\mathbb{E}_{i \ne j}\left[ \ln p (\mathbf{x}, \mathbf{z}) \right]$ 往往有闭式解，使得上式能对隐变量高效推断。

由于在对 $\mathbf{z}_{j}$ 所服从的分布 $q_{j}^{*}$ 估计时融合了 $z_{j}$ 以外的 $z_{i \ne j}$ 的信息，这是通过联合似然函数 $\ln p(\mathbf{x}, \mathbf{z})$ 在 $\mathbf{z}_{j}$ 之外的隐变量分布上求期望得到的，因此亦成为平均场（mean field）方法。

在实际应用中，最重要的是考虑如何对隐变量进行拆解，以及假设各变量子集服从何种分布，在此基础之上结合EM算法对概率图模型进行推断和参数估计。

话题模型

话题模型（topic model）是一类生成式有向图模型，主要用来处理离散型的数据集合（如文本集合）。

作为一种非监督产生式模型，话题模型能够有效利用海量数据发现文档集合中隐含的语义。

隐狄里克雷分配模型（Latent Dirichlet Allocation, LDA）是话题模型的典型代表。

LDA 的基本单元：

• 词（word）：词是待处理数据的基本离散单元。
  • 例如在文本处理任务中，一个词就是一个英文单词或有独立意义的中文。
• 文档（document）：文档是待处理的数据对象，它由一组词组成，这些词在文档中是不计顺序的。
  • 例如一篇论文、一个网页等。
  • 这样的表示方式称为词袋（bag-of-words）。
  • 数据对象只要能用词袋描述，就可使用话题模型。
• 话题（topic）：话题表示一个概念，具体表示为一系列相关的词，以及它们在该概念下出现的概率。

形象地说，如下图所示：

• 一个话题就像是一个箱子，里面装着在这个概念下出现概率较高的那些词；
• 假设数据集中一共包含 $K$ 个话题和 $T$ 篇文档；
• 文档中的词来自一个包含 $N$ 个词的词典；
• 用 $T$ 个 $N$ 维向量 $\mathbf{W} = \left\lbrace \bm{w}_1, \dots, \bm{w}_T \right\rbrace$ 表示数据集（即文档集合）；
• $K$ 个 $N$ 维向量 $\bm{\beta}_{k}\, (k = 1, \dots, K)$ 表示话题；
• $\bm{w}_t \in \R^N$ 的第 $n$ 个分量 $w_{t, n}$ 表示文档 $t$ 中词 $n$ 的词频；
• $\bm{\beta}_{k} \in \R^N$ 的第 $n$ 个分量 $\beta_{k, n}$ 表示话题 $k$ 中词 $n$ 的词频。

现实任务中可通过统计文档中出现的词来获得词频向量 $\bm{w}_i$，但通常并不知道这组文档谈论了哪些话题，也不知道每篇文档与哪些话题有关。

LDA 从生成式模型的角度看待文档和话题。具体来说，LDA 认为每篇文档包含多个话题，用向量 $\Theta_t \in \R^{K}$ 表示文档 $t$ 中所包含的每个话题的比例，$\Theta_{t, k}$ 即表示文档 $t$ 中包含话题 $k$ 的比例，进而通过下面的步骤由话题「生成」文档 $t$：

1. 根据参数为 $\bm{\alpha}$ 的狄利克雷分布随机采样一个话题分布 $\Theta_t$；
2. 按如下步骤生成文档中的 $N$ 个词：
   1. 根据 $\Theta_t$ 进行话题指派，得到文档 $t$ 中词 $n$ 的话题 $z_{t, n}$；
   2. 根据指派的话题所对应的词频分布 $\bm{\beta}_{k}$ 随机采样生成词。

上图也演示出根据以上步骤生成文档的过程。这样生成的文档自然地以不同比例包含多个话题（1.），文档中的每个词来自一个话题（2.2.），而这个话题是根据话题比例产生的（2.1.）。

下图描述了 LDA 的变量关系，其中文档中的词频 $w_{t, n}$ 是唯一的已观测变量，它依赖于对这个词进行的话题指派 $z_{t, n}$，以及话题所对应的词频 $\bm{\beta}_{k}$；同时话题指派 $z_{t, n}$ 依赖于文档中的话题分布 $\Theta_t$，$\Theta_t$ 依赖于狄利克雷分布的参数 $\bm{\alpha}$，而话题词频则依赖于参数 $\bm{\eta}$。

于是 LDA 模型对应的概率分布为

$$p(\mathbf{W}, \mathbf{z}, \bm{\beta}, \Theta \mid \bm{\alpha}, \bm{\eta}) = \prod_{t = 1}^T p(\Theta_t \mid \bm{\alpha}) \prod_{k = 1}^K p(\bm{\beta}_k \mid \bm{\eta}) \left(\prod_{n = 1}^N p(w_{t, n} \mid z_{t, n}, \bm{\beta}_{k}) p(z_{t, n} \mid \Theta_t) \right)$$

其中 $p(\Theta_t \mid \bm{\alpha}),\, p(\bm{\beta}_{k} \mid \bm{\eta})$ 通常分别设置为以 $\bm{\alpha}, \bm{\eta}$ 为参数的 $K, N$ 维狄利克雷分布。

狄利克雷分布相关内容略。

模型参数估计部分：

给定训练数据 $\mathbf{W} = \left\lbrace \bm{w}_1, \dots, \bm{w}_T \right\rbrace$，LDA 的模型参数可通过极大似然法估计，即寻找 $\bm{\alpha}, \bm{\eta}$ 以最大化对数似然

$$LL(\bm{\alpha}, \bm{\eta}) = \sum_{t=1}^{T} \ln p(\bm{w}_t \mid \bm{\alpha}, \bm{\eta})$$

由于 $p(\bm{w}_t \mid \bm{\alpha}, \bm{\eta})$ 不易计算，上式难以直接求解，因此实践中常采用「变分法」求取近似解。

模型推断部分：

若模型已知，即 $\bm{\alpha}, \bm{\eta}$ 确定，则根据词频 $w_{t, n}$ 来推断文档集所对应的话题结构（即推断 $\Theta_t,\, \bm{\beta}_{k},\, z_{t, n}$）可通过求解

$$p(\mathbf{z}, \bm{\beta}, \Theta \mid \mathbf{W}, \bm{\alpha}, \bm{\eta}) = \dfrac{p(\mathbf{W}, \mathbf{z}, \bm{\beta}, \Theta \mid \bm{\alpha}, \bm{\eta})}{p(\mathbf{W} \mid \bm{\alpha}, \bm{\eta})}$$

由于分母上的 $p(\mathbf{W} \mid \bm{\alpha}, \bm{\eta})$ 难以获取，上式难以直接求解，事件中常采用「吉布斯采样」或「变分法」进行近似推断。