46、贝叶斯核方法中的推理与高斯过程

贝叶斯核方法中的推理与高斯过程

1. 边缘化与超参数

1.1 边缘化

在概率计算中，若已知 $p(y \mid Y)$ 和 $p(Y)$ 的关系，即 $p(y \mid Y) = \frac{p(y, Y)}{p(Y)}$，其中 $p(Y)$ 可通过对 $p(y, Y)$ 关于 $y$ 积分得到，即 $p(Y) = \int p(y, Y)dy$，那么就可以得到 $p(y \mid Y) = \frac{p(y, Y)}{\int p(y, Y)dy}$，这个过程被称为边缘化。

1.2 超参数

在实际应用中，我们可能会对一些参数的权衡产生疑问，例如 $|f|^2$ 和 $|\nabla_x f|^2$ 的权衡是否正确。此时，我们可以引入超参数 $\lambda$ 来调整，即 $-\ln p(f) = c + \lambda |f|^2 + \gamma |\nabla_x f|^2$，其中 $\gamma > 0$，这实际上改变了我们对函数平坦性假设的尺度。

同时，我们可能并不总是知道似然函数的确切形式，但对从 $f(x)$ 得到 $y$ 的过程中涉及的加性噪声量有一个大致的概念。因此，我们需要一种方法来编码我们对似然和先验概率指定中使用的超参数的不确定性。

我们用 $\theta$ 表示特定情况下所需的所有超参数向量，得到 $p(f \mid \theta) = p(f \mid \theta)p(\theta)$，进而 $p(f) = \int p(f \mid \theta)d\theta = \int p(f \mid \theta)p(\theta)d\theta$。我们称 $p(\theta)$ 为超先验，因为它是对先验 $p(f)$ 的先验假设。

理论上，我们可以对假设 $f$ 进行积分，以得到超参数 $\theta$ 的后验分布 $p(\theta \mid Z) \propto p(Z \mid \theta)p(\theta) = p(\theta) \int p(Z \mid f)p(f \mid \theta)df$，并使用它来得到 $p(f \mid Z) = \int p(f \mid \theta, Z)p(\theta \mid Z)d\theta$。但通常情况下，这些积分的解析解是不可行的，我们必须采用近似方法。

2. 推理方法

2.1 最大后验近似（MAP）

在大多数情况下，对 $p(f \mid X, Y)$ 的积分（如 (16.15) 中 $f$ 的期望）在计算上是难以处理的。因此，我们通常使用近似技术来进行预测。一种流行的近似方法是用后验分布模式下的被积函数值替换对 $p(f \mid X, Y)$ 的积分，即 $p(y \mid X, Y, x) \approx p(y \mid f_{MAP}, x)$，其中 $f_{MAP} = \arg\max_f p(f \mid X, Y)$。

在实践中，我们通过最小化负对数后验来获得 $f_{MAP}$，即 $f_{MAP} = \arg\min_f {-\ln p(f \mid Z)} = \arg\min_f {-\ln p(Z \mid f) - \ln p(f)}$。这种方法的一个额外优点是我们完全避免了归一化问题，因为该公式不依赖于 $p(Z)$。如果所有涉及的分布恰好都是高斯分布，由于均值和模式重合，这种近似是合理的。

对于超参数 $\theta$ 的积分，由于超先验 $p(\theta)$ 的存在，我们也可能需要进行近似。这种情况比计算 MAP 估计更常见，因为通过对 $\theta$ 进行积分得到的复杂先验分布 $p(f)$ 可能会使后续的积分步骤变得难以处理。因此，我们选择 $\theta_{MAP} = \arg\max_{\theta} p(\theta \mid Z) = \arg\min_{\theta} {-\ln p(Z \mid \theta) - \ln p(\theta)}$，这一过程有时被称为 MAP2 估计。

2.2 后验分布的参数近似

除了用后验分布的模式替换 $p(f \mid Z)$ 外，我们还可以采用更复杂的近似方法。一种改进是使用正态分布 $N(\mu, \Sigma^{-1})$，其中均值 $\mu$ 与 $p(f \mid Z)$ 的模式重合，方差 $\Sigma$ 由 $-\ln p(f_{MAP} \mid Z)$ 的二阶导数确定，即 $f \mid Z \sim N(E[f \mid Z], \Sigma^{-1})$，其中 $\Sigma = \left[-\nabla^2_f \ln p(f \mid Z)\right]_{E[f \mid Z]}$。这种方法的优点是积分仍然是可处理的，而且正态分布是所有具有有界方差的分布中信息最少（熵最大）的分布。

当单个高斯分布不足以捕捉 $p(y \mid X, Y, x)$ 的重要属性时，我们可以使用更复杂的参数模型 $q(f)$（如高斯密度的混合）来改进对 (16.15) 的近似。一种常见的策略是采用变分方法。

2.3 与正则化风险泛函的联系

观察 (16.23) 可以发现，负对数后验 $-\ln p(f \mid Z)$ 与正则化风险泛函 $R_{reg}[f]$ 有相似之处。两者都由两项之和组成：一项依赖于 $f$ 和 $Z$（经验风险 $R_{emp}[f]$ 和 $-\ln p(Z \mid f)$），另一项与 $Z$ 无关（正则化项 $\Omega[f]$ 和 $-\ln p(f)$）。

具体来说，如果我们形式上设置 $mR_{emp}[f] = \sum_{i=1}^m c(x_i, y_i, f(x_i)) \sim -\ln p(Z \mid f)$ 和 $\Omega[f] \sim -\ln p(f)$，我们可以得到正则化风险最小化的贝叶斯解释，即作为 MAP 估计，反之亦然。

正则化风险泛函中的损失函数相当于概率设置中的负对数似然。例如，平方损失对应于假设数据中添加了正态噪声。而正则化项的选择与贝叶斯估计中函数类的先验选择具有相同的影响。

需要注意的是，正则化风险最小化和贝叶斯方法之间的对应关系仅适用于通过最大化对数后验来获得 MAP 解的算法。

2.4 符号翻译

为了清晰起见，我们列出了贝叶斯估计和风险泛函方法中相应量的对应关系，如下表所示：
| 风险泛函方法 | 贝叶斯方法 |
| — | — |
| 经验风险 $\sum_{i=1}^m c(x_i, y_i, f(x_i))$ | 负对数似然 $-\sum_{i=1}^m \ln p(y_i \mid f(x_i))$ |
| 正则化 $\Omega[f]$ | 负对数先验 $-\ln p(f)$ |
| 正则化风险 $R_{emp}[f] + \lambda\Omega[f]$ | 负对数后验 $-\ln p(Z \mid f) - \ln p(f)$ |
| 风险最小化器 | MAP 估计 |

3. 高斯过程

3.1 相关观测

当我们观察函数值 $f(x_i)$ 在位置 $x_i$ 处的值时，自然会假设这些值是相关的，这取决于它们的位置 $x_i$。实际上，我们假设 $f(x_i)$ 形成一个具有均值 $\mu$ 和协方差矩阵 $K$ 的正态分布。

以两个观测值 $f(x_1)$ 和 $f(x_2)$ 为例，假设均值 $\mu = (0, 0)$，协方差矩阵 $K = \begin{bmatrix}1 & \frac{3}{4} \ \frac{3}{4} & \frac{3}{4}\end{bmatrix}$。当我们观察到 $f(x_1)$ 时，我们可以得到 $f(x_2)$ 的条件密度 $p(f(x_2) \mid f(x_1)) = \frac{p(f(x_1), f(x_2))}{p(f(x_1))}$。知道条件密度后，$f(x_2)$ 的均值不再为 0，方差也会减小。

3.2 定义和基本概念

假设我们有一个在点 $x_1, \cdots, x_m$ 处的观测值 $t_i$ 的分布，我们假设这些观测值是由高斯过程生成的。

定义 ：设 $t(x)$ 是由 $x \in \mathcal{X}$（$\mathcal{X}$ 是任意索引集）参数化的随机过程。如果对于任何 $m \geq 1$ 和 $x_1, \cdots, x_m \in \mathcal{X}$，随机变量 $(t(x_1), \cdots, t(x_m))$ 是正态分布的，则称 $t(x)$ 是高斯过程。

我们用 $k(x, x’)$ 表示生成协方差矩阵 $K = \text{cov}[t(x_1), \cdots, t(x_m)]$ 的函数，用 $\mu$ 表示分布的均值，且 $K_{ij} = k(x_i, x_j)$，则 $(t(x_1), \cdots, t(x_m)) \sim N(\mu, K)$。

函数 $k(x, x’)$ 是明确定义的、对称的，并且矩阵 $K$ 是正定的。因此，$k$ 是一个可接受的核。

在标准设置中，$\mu = 0$，这意味着我们对估计的特定值没有先验知识，但假设较小的值更受青睐。对于给定的一组 $(t(x_1), \cdots, t(x_m)) =: t$，先验密度函数为 $p(t) = (2\pi)^{-\frac{m}{2}} (\det K)^{-\frac{1}{2}} \exp\left(-\frac{1}{2}t^T K^{-1}t\right)$。

通过简单的代换 $t = K\alpha$，我们有 $\alpha \sim N(0, K^{-1})$，并且 $p(\alpha) = (2\pi)^{-\frac{m}{2}} (\det K)^{-\frac{1}{2}} \exp\left(-\frac{1}{2}\alpha^T K\alpha\right)$。取对数后，我们发现这与正则化框架中的 $\Omega[f]$ 相同，这表明高斯过程先验与使用再生核希尔伯特空间框架的估计器之间存在联系。

3.3 简单假设

高斯过程先验认为低复杂度的假设对应于向量 $y$，使得 $y^T K^{-1}y$ 较小。特别是对于具有大特征值 $\lambda_i$ 的 $K$ 的（归一化）特征向量 $v_i$，因为 $Kv_i = \lambda_i v_i$ 意味着 $v_i^T K^{-1}v_i = \lambda_i^{-1}$。这意味着 $K$ 的谱和特征系统代表了一种实际手段，可以观察特定先验对估计平滑度的影响。

例如，对于高斯协方差核 $k(x, x’) = \exp\left(-\frac{|x - x’|^2}{2\sigma^2}\right)$（其中 $\sigma = 1$），在 $[-5, 5]$ 上的均匀分布假设下，我们可以得到一些简单的基础假设函数，这些函数与傅里叶分解有相似之处，表明该核强烈偏好缓慢振荡的函数。

3.4 回归

在回归估计中，我们通常假设在生成 $t(x_i)$ 的过程中添加了噪声，即 $y_i = t(x_i) + \epsilon_i$，其中 $\epsilon_i$ 是独立的零均值随机变量。为了简化符号，我们假设所有 $\epsilon_i$ 都来自相同的分布 $\epsilon \sim p(\epsilon)$。

我们可以得到似然 $p(y \mid t(x))$ 为 $p(y_i \mid t(x_i)) = p(y_i - t(x_i))$，后验分布为 $p(t \mid y) \propto p(y \mid t)p(t) = \left(\prod_{i} p(y_i - t(x_i))\right) \frac{1}{\sqrt{(2\pi)^m \det K}} \exp\left(-\frac{1}{2}t^T K^{-1}t\right)$。

在高斯过程回归中，一个流行的选择是假设加性正态噪声 $\epsilon_i \sim N(0, \sigma^2)$。这样做有几个优点：所有涉及的分布在推理过程中都保持正态，我们可以计算精确的解；我们可以通过简单的矩阵求逆找到分布的模式。

经过代换 $t = K\alpha$ 并取对数后，忽略与 $\alpha$ 无关的项，我们得到 $\ln p(\alpha \mid y) = -\frac{1}{2\sigma^2} |y - K\alpha|^2 - \frac{1}{2}\alpha^T K\alpha + c$。这显然是一个正态分布，因此模式和均值重合，MAP 近似 (16.22) 变得精确。通过最大化 (16.46) 得到 $\alpha = (K + \sigma^2I)^{-1}y$。

知道 $\alpha$ 后，我们可以预测新位置 $x$ 处的 $y$。对于加性高斯噪声，我们可以通过计算得到 $y(x_{m+1})$ 的方差 $Var y_{m+1} = k(x_{m+1}, x_{m+1}) + \sigma^2 - k^T (K + \sigma^2I)^{-1}k$，从而知道 $y(x_{m+1}) \sim N(K(K + \sigma^2I)^{-1}y, k(x_{m+1}, x_{m+1}) + \sigma^2 - k^T (K + \sigma^2I)^{-1}k)$。

3.5 分类

为了简单起见，我们仅限于二分类问题。我们使用逻辑回归来建模条件概率 $P(y = 1 \mid x)$ 和 $P(y = -1 \mid x)$。

高斯过程分类模型的后验密度为 $p(f \mid Z) \propto p(Y \mid X, t) = \prod_{i=1}^m p(y_i \mid t(x_i)) \exp\left(-\frac{1}{2}t^T K^{-1}t\right)$。通过变换 $t = K\alpha$，负对数后验密度变为 $-\ln p(f \mid X, Y) = \sum_{i=1}^m -\ln p(y_i \mid t(x_i)) + \frac{1}{2}\alpha^T K\alpha$。

如果我们采用 16.2.1 节中的 MAP 方法，推理通过搜索密度 $p(f \mid X, Y)$ 的模式来进行，因此估计问题可以简化为一个非线性函数最小化问题。

3.6 调整高斯过程的超参数

通常，我们事先不知道加性噪声的确切量或协方差核的具体形式。为了解决这个问题，我们需要使用 16.1.4 节中的超参数形式。

对于回归问题，假设加性高斯噪声，(16.48) 可以写成 $p(y \mid \theta) = \frac{1}{\sqrt{(2\pi)^m \det(K + \sigma^2I)}} \exp\left(-\frac{1}{2}y^T (K + \sigma^2)^{-1}y\right)$，其中 $K$ 和 $\sigma$ 是 $\theta$ 的函数。

MAP2 估计的基本思想是通过最大化 $p(y \mid \theta)p(\theta)$ 来最大化 $p(\theta \mid y)$。在实践中，这可以通过梯度上升或二阶方法来实现。我们需要计算 (16.55) 关于 $\theta$ 的梯度，具体表达式为 $\nabla_{\theta} [-\ln p(y \mid \theta)p(\theta)] = \frac{1}{2} \nabla_{\theta} (\ln\det Q) - \frac{1}{2} \nabla_{\theta} \left[y^T Q^{-1}y\right] - \nabla_{\theta} \ln p(\theta) = -\frac{1}{2} \text{tr}\left[Q^{-1} \nabla_{\theta} Q\right] + \frac{1}{2}y^T Q^{-1} (\nabla_{\theta} Q) Q^{-1}y - \nabla_{\theta} \ln p(\theta)$，其中 $Q = K + \sigma^2I$。

如果我们假设一个平坦的超先验 $p(\theta) = \text{const}$，则对 $\theta$ 的优化就变成了在 $-\ln p(y \mid \theta)$ 上的梯度下降。但计算 (16.57) 在数值上仍然非常昂贵，因为它涉及到 $Q$（一个 $m \times m$ 矩阵）的求逆。

下面是高斯过程推理和超参数调整的流程图：

graph TD;
    A[开始] --> B[确定高斯过程模型];
    B --> C[获取观测数据];
    C --> D[选择推理方法];
    D --> E{是否需要调整超参数};
    E -- 是 --> F[设置超参数先验];
    F --> G[计算似然和后验];
    G --> H[使用 MAP 或 MAP2 估计];
    H --> I[更新超参数];
    I --> J[进行预测];
    E -- 否 --> K[直接进行预测];
    J --> L[结束];
    K --> L;

综上所述，贝叶斯核方法中的推理和高斯过程为我们提供了强大的工具，用于处理复杂的概率建模和预测问题。通过合理选择近似方法和调整超参数，我们可以在不同的应用场景中获得较好的性能。

4. 推理方法与高斯过程的应用总结

4.1 推理方法的选择

在实际应用中，选择合适的推理方法至关重要。以下是对不同推理方法的总结和适用场景分析：
| 推理方法 | 优点 | 缺点 | 适用场景 |
| — | — | — | — |
| 最大后验近似（MAP） | 避免归一化问题，计算相对简单；在高斯分布情况下近似合理 | 对于多峰分布可能效果不佳 | 当后验分布集中在模式附近，且分布接近高斯分布时 |
| 后验分布的参数近似（如高斯近似） | 积分可处理，正态分布信息最少 | 单个高斯分布可能无法捕捉复杂分布特征 | 后验分布相对简单，可由正态分布较好近似时 |
| 变分近似 | 可使用更复杂的参数模型改进近似 | 技术细节复杂 | 后验分布复杂，单个高斯无法满足需求时 |

4.2 高斯过程的应用步骤

高斯过程在回归和分类问题中都有广泛应用，以下是其一般应用步骤：
1. 确定高斯过程模型 ：
- 选择合适的协方差函数 $k(x, x’)$，如高斯协方差核 $k(x, x’) = \exp\left(-\frac{|x - x’|^2}{2\sigma^2}\right)$。
- 确定均值 $\mu$，通常设为 0。
2. 获取观测数据 ：收集在点 $x_1, \cdots, x_m$ 处的观测值 $y_i$。
3. 选择推理方法 ：根据问题的复杂度和数据特点，选择 MAP 近似、参数近似或其他方法。
4. 处理超参数 ：
- 若需要调整超参数，设置超参数先验 $p(\theta)$。
- 使用 MAP2 估计等方法更新超参数。
5. 进行预测 ：
- 对于回归问题，计算新位置 $x$ 处的预测值和方差。
- 对于分类问题，通过搜索后验密度的模式进行推理。

4.3 超参数调整的注意事项

在调整高斯过程的超参数时，需要注意以下几点：
- 计算复杂度 ：计算梯度和 Hessian 涉及矩阵求逆，对于大规模数据，计算量可能非常大。可以考虑使用稀疏贪婪近似等方法降低计算复杂度。
- 超先验的选择 ：平坦超先验 $p(\theta) = \text{const}$ 可以简化计算，但可能会影响估计的准确性。在实际应用中，需要根据问题的特点选择合适的超先验。
- 收敛性 ：使用梯度上升或二阶方法进行优化时，需要确保算法的收敛性。可以通过调整学习率、迭代次数等参数来提高收敛速度和稳定性。

5. 案例分析

5.1 回归案例

假设我们有一个回归问题，数据是在区间 $[-5, 5]$ 上均匀分布的点，并且添加了高斯噪声。我们使用高斯过程进行回归估计，具体步骤如下：
1. 确定模型 ：选择高斯协方差核 $k(x, x’) = \exp\left(-\frac{|x - x’|^2}{2}\right)$，均值 $\mu = 0$。
2. 获取数据 ：生成 200 个均匀分布的点 $x_i$，并计算对应的 $y_i = t(x_i) + \epsilon_i$，其中 $\epsilon_i \sim N(0, 0.1^2)$。
3. 选择推理方法 ：使用 MAP 近似，假设噪声为加性高斯噪声。
4. 调整超参数 ：使用 MAP2 估计调整噪声方差 $\sigma^2$。
5. 进行预测 ：计算新位置 $x$ 处的预测值和方差。

以下是 Python 代码示例：

import numpy as np
import matplotlib.pyplot as plt

# 生成数据
np.random.seed(0)
x = np.linspace(-5, 5, 200)
t = np.sin(x)
y = t + np.random.normal(0, 0.1, size=x.shape)

# 定义协方差函数
def kernel(x1, x2, sigma=1):
    dist = np.sum(x1**2, 1).reshape(-1, 1) + np.sum(x2**2, 1) - 2 * np.dot(x1, x2.T)
    return np.exp(-dist / (2 * sigma**2))

# 计算后验分布
def posterior(x, y, x_new, sigma_n=0.1, sigma=1):
    K = kernel(x[:, np.newaxis], x[:, np.newaxis], sigma)
    K_s = kernel(x[:, np.newaxis], x_new[:, np.newaxis], sigma)
    K_ss = kernel(x_new[:, np.newaxis], x_new[:, np.newaxis], sigma)
    L = np.linalg.cholesky(K + sigma_n**2 * np.eye(len(x)))
    alpha = np.linalg.solve(L.T, np.linalg.solve(L, y))
    mu = np.dot(K_s.T, alpha)
    v = np.linalg.solve(L, K_s)
    cov = K_ss - np.dot(v.T, v)
    return mu, cov

# 进行预测
x_new = np.linspace(-5, 5, 100)
mu, cov = posterior(x, y, x_new)
std = np.sqrt(np.diag(cov))

# 绘制结果
plt.figure(figsize=(10, 6))
plt.scatter(x, y, c='r', label='Data')
plt.plot(x_new, mu, c='b', label='Mean')
plt.fill_between(x_new, mu - 2 * std, mu + 2 * std, color='gray', alpha=0.2, label='95% CI')
plt.title('Gaussian Process Regression')
plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.show()

5.2 分类案例

在二分类问题中，我们可以使用高斯过程进行分类。假设我们有一个简单的二维数据集，使用逻辑回归和高斯过程进行分类，步骤如下：
1. 确定模型 ：选择合适的协方差函数和均值。
2. 获取数据 ：生成二维数据集，并标记类别。
3. 选择推理方法 ：使用 MAP 方法搜索后验密度的模式。
4. 进行预测 ：计算新样本的分类概率。

以下是一个简化的流程图，展示了高斯过程分类的步骤：

graph TD;
    A[确定高斯过程模型] --> B[获取二维数据集];
    B --> C[选择 MAP 推理方法];
    C --> D[计算后验密度];
    D --> E[搜索后验密度模式];
    E --> F[进行分类预测];

6. 总结与展望

6.1 总结

贝叶斯核方法中的推理和高斯过程为解决复杂的概率建模和预测问题提供了强大的工具。通过边缘化和超参数的引入，我们可以更好地处理不确定性。不同的推理方法（如 MAP 近似、参数近似）适用于不同的场景，在实际应用中需要根据具体情况进行选择。

高斯过程通过假设相关观测和使用协方差函数，能够有效地处理回归和分类问题。在回归中，通过假设加性噪声和使用高斯分布，可以得到精确的解；在分类中，使用逻辑回归和搜索后验密度的模式可以进行有效的分类。

6.2 展望

虽然贝叶斯核方法和高斯过程已经取得了很多成果，但仍有一些挑战和发展方向：
- 计算复杂度 ：在大规模数据和复杂模型下，计算复杂度仍然是一个问题。未来可以研究更高效的算法和近似方法，如稀疏高斯过程、随机变分推理等。
- 模型选择 ：如何选择合适的协方差函数和超参数仍然是一个难题。可以探索自动化的模型选择方法，如贝叶斯模型平均、信息准则等。
- 多模态数据处理 ：随着数据的多样化，如何处理多模态数据（如图像、文本、音频等）是一个重要的研究方向。可以将高斯过程与深度学习等方法相结合，处理更复杂的数据。

总之，贝叶斯核方法和高斯过程在机器学习和统计学领域具有广阔的应用前景，未来的研究将不断推动其发展和应用。