第 03 章：潜在空间优化算法：高斯过程贝叶斯优化与梯度上升

“如果我们已经修好了一条将离散分子映射到连续流形的高速公路，
那么下一个核心命题便是：我们应当如何在高速公路上导航，以最少的实验代价，精准命中性能最高的未知物质？”

1. 引言：连续流形上的导航艺术

在第 02 章中，我们建立了一个连续、低维且可微的潜在流形 $\mathcal{Z} \subset \mathbb{R}^d$ 。每一个潜变量向量 $z \in \mathcal{Z}$ 都可以通过生成解码器 $p_\theta(x \mid z)$ 还原为一个离散的分子拓扑或周期性无机晶格。

现在，逆向设计问题正式转化为一个在连续空间上的函数极值寻优问题：

z^* = \arg\max_{z \in \mathcal{Z}} f(p_\theta(x \mid z))

其中 $f(x)$ 是我们渴望优化的目标物理性能（如太阳能电池吸光系数、受体结合常数、固态离子电导率）。

然而，在材料科学与计算化学中，评估一个候选结构的真实性能 $f(x)$ 是极其昂贵且耗时的：

如果调用实验湿法合成，每个样本的验证周期长达数天甚至数周，金钱成本数千美元；
即使调用第一性原理从头算量子力学（如杂化泛函 DFT 或多体微扰 $GW$ ），单个样本的计算耗时也在数小时到数十小时之间。

因此，优化算法的“样本效率（Sample Efficiency）”与“不确定性量化能力（Uncertainty Quantification）”成为了决定逆向设计成败的生命线。

本章将深入剖析两大核心优化算法——可微代理表面梯度上升法与高斯过程贝叶斯优化（GP-BO），揭示其数理本质、致命陷阱与前沿应对策略。

2. 连续代理表面的梯度上升法及其对抗陷阱

2.1 梯度上升的基础形式

最直观的想法是在隐空间 $\mathcal{Z}$ 上训练一个端到端的可微性质代理回归器 $\hat{y} = g_\psi(z) \approx f(p_\theta(x \mid z))$ 。由于神经网络 $g_\psi$ 是处处可导的，我们可以直接计算关于潜在坐标 $z$ 的雅可比梯度向量，并通过梯度上升法进行多步迭代：

z_{t+1} = z_t + \eta \nabla_z g_\psi(z_t)

               代理性质曲面 g(z) 极大值山峰
                      ▲
                     / \
                    /   \  ◄── 沿着梯度上升 ∇_z g(z)
                   /     \
    起始先导点 z_0 ●───────

2.2 致命陷阱：对抗假象（Adversarial Hallucination / OOD Exploitation）

在理想情况下，梯度上升似乎能毫不费力地把我们引向性能最佳的晶体构型。但在实际计算中，无约束的梯度上升往往会迅速掉入极度危险的“对抗假象陷阱”：

模型黑客行为（Surrogate Hacking）：神经网络代理模型 $g_\psi(z)$ 仅仅在训练集数据覆盖的高维球壳（半径 $R \approx \sqrt{d}$ ）附近具有较好的插值泛化能力；
脱离真实物理支撑流形：在远离训练分布的数据稀疏盲区（Out-of-Distribution, OOD），高维多项式或非线性激活函数很容易产生剧烈的虚假极值尖峰（Pathological Extrapolation Spikes）。梯度上升算法会像寻觅到漏洞的黑客一样，以极快的速度将潜变量 $z$ 强行推向半径巨大（ $\|z\| \gg \sqrt{d}$ ）的极端空洞区域。
解码产物崩溃：当把这个被代理模型打了极高虚假分数的 $z^*$ 喂给解码器 $p_\theta(x \mid z^*)$ 时，解码器由于从未在如此偏远的极端区域训练过，会输出严重的化学键断裂、原子核严重重叠（Nuclear Fusion Collapse）或电荷极度不守恒的“怪胎结构”。

2.3 破局方案：流形支撑惩罚项（Manifold Prior Regularization）

为了杜绝梯度上升脱轨越界，必须在目标函数中引入物理先验密度惩罚项，强制限制优化轨迹始终留在先验高斯分布的致密流形之内：

\max_z \mathcal{J}(z) = g_\psi(z) + \lambda_{\text{prior}} \log p(z) = g_\psi(z) - \frac{\lambda_{\text{prior}}}{2} \|z\|^2

相应的梯度更新规则变为带有**高斯收缩阻尼项（Weight Decay Damping）**的修正更新：

z_{t+1} = (1 - \eta \lambda_{\text{prior}}) z_t + \eta \nabla_z g_\psi(z_t)

这一约束有效拉住了试图脱缰而去的梯度粒子，确保生成的结构在物理有效性与性能提升之间维持良性平衡。

3. 高斯过程贝叶斯优化（Gaussian Process Bayesian Optimization, GP-BO）

相比于对未知空间盲目自信的简单神经网络，**高斯过程贝叶斯优化（GP-BO）**以其严密的不确定性推断和无与伦比的超高样本效率，成为了逆向设计领域最为经典且广受信赖的黄金武器。

图 3-1：高斯过程贝叶斯优化在材料潜在空间中的代理建模与主动探索机制示意图

图 3-2：高斯过程后验均值、置信区间（±1.96σ）与期望改进采集函数（EI）的主动学习极值决策验证

3.1 高斯过程（Gaussian Process, GP）的严密数学推导

一个高斯过程是由其均值函数 $m(z)$ 与协方差核函数 $k(z, z')$ 唯一定义的无限维连续随机过程：

f(z) \sim \mathcal{GP}\left(m(z), k(z, z')\right)

在材料科学中，最广泛采用的是能够刻画不同平滑度的 Matérn 5/2 核函数：

k_{\text{Matérn5/2}}(z, z') = \sigma_f^2 \left( 1 + \sqrt{5}r + \frac{5}{3}r^2 \right) \exp(-\sqrt{5}r), \quad r = \sqrt{\sum_{i=1}^d \frac{(z_i - z'_i)^2}{l_i^2}}

假设我们已经通过实验或高精度 DFT 真实评测了 $N$ 个潜空间候选点，构成数据集 $\mathcal{D}_{1:N} = \left\{ (z_n, y_n) \right\}_{n=1}^N$ 。对于任意未知的测试点 $z_*$ ，其先验联合高斯分布为：

\begin{bmatrix} \mathbf{y}_{1:N} \\ f(z_*) \end{bmatrix} \sim \mathcal{N}\left( \mathbf{0}, \begin{bmatrix} \mathbf{K} + \sigma_n^2 I & \mathbf{k}_* \\ \mathbf{k}_*^T & k(z_*, z_*) \end{bmatrix} \right)

其中 $\mathbf{K}_{ij} = k(z_i, z_j)$ 为历史样本间的协方差矩阵， $\mathbf{k}_* = [k(z_1, z_*), \dots, k(z_N, z_*)]^T$ 。

利用多元高斯分布的条件化定理（Conditioning Theorem），我们可以解析导出未知点 $z_*$ 的后验分布依然是单变量高斯分布：

f(z_*) \mid \mathcal{D}_{1:N} \sim \mathcal{N}\left(\mu(z_*), \sigma^2(z_*)\right)

后验均值（代表对性质的最佳期望预测）：
$\mu(z_*) = \mathbf{k}_*^T (\mathbf{K} + \sigma_n^2 I)^{-1} \mathbf{y}_{1:N}$
后验方差（严格量化由于数据缺乏带来的认知不确定度）：
$\sigma^2(z_*) = k(z_*, z_*) - \mathbf{k}_*^T (\mathbf{K} + \sigma_n^2 I)^{-1} \mathbf{k}_*$

3.2 采集函数（Acquisition Functions）推导：平衡利用与探索

我们不仅需要知道哪里期望值高，更需要知道哪里值得我们“下注冒险”。采集函数 $\alpha(z)$ 就是将后验均值 $\mu(z)$ 与后验方差 $\sigma(z)$ 综合折算为一个评测价值标量的决策准则：

z_{N+1} = \arg\max_{z \in \mathcal{Z}} \alpha(z)

1. 期望提升准则 (Expected Improvement, EI)

设当前历史已发现的最优性质值为 $y^+ = \max_{1 \le n \le N} y_n$ 。对于新候选点 $z$ ，其相对于历史最佳的提升量为：

I(z) = \max\left(0, f(z) - y^+\right)

期望提升 EI 是改善量在后验高斯分布下的数学期望：

\alpha_{\text{EI}}(z) = \mathbb{E}[I(z)] = \int_{y^+}^\infty (y - y^+) \frac{1}{\sqrt{2\pi}\sigma(z)} \exp\left( -\frac{(y - \mu(z))^2}{2\sigma^2(z)} \right) \mathrm{d}y

利用变量代换与高斯积分，我们可以解析推导出无与伦比的闭式解析解（Closed-form Formula）：

\alpha_{\text{EI}}(z) = \begin{cases} (\mu(z) - y^+ - \xi) \, \Phi(Z) + \sigma(z) \, \phi(Z), & \text{若 } \sigma(z) > 0 \\ 0, & \text{若 } \sigma(z) = 0 \end{cases}

其中：

Z = \frac{\mu(z) - y^+ - \xi}{\sigma(z)}

$\Phi(\cdot)$ 与 $\phi(\cdot)$ 分别为标准正态分布的累积分布函数（CDF）与概率密度函数（PDF）；
$\xi \ge 0$ 为探索调节参数。
物理项深度解构：
- 第一项 $(\mu(z) - y^+) \Phi(Z)$ 促使算法走向均值高的区域——利用（Exploitation）；
- 第二项 $\sigma(z) \phi(Z)$ 促使算法走向不确定性大、从未探索过的盲区——探索（Exploration）！

2. 上置信界准则 (Upper Confidence Bound, UCB)

\alpha_{\text{UCB}}(z) = \mu(z) + \kappa \cdot \sigma(z)

通过调节超参数 $\kappa > 0$ ，科学家可以直接控制算法是以保守稳健为主，还是以开拓发现冷门新骨架为主。

4. 工业级代码实战：端到端潜空间高斯过程贝叶斯优化

以下代码使用纯 PyTorch 构建高斯过程后验更新器与闭式期望提升（EI）采集函数优化器，完成从历史候选池采样、贝叶斯代理表面更新到捕获未知全局极值点的完整可闭环流程：

python

"""
文件名: latent_bayesian_optimization.py
功能: 实现基于高斯过程与期望提升 (Expected Improvement) 的潜在空间贝叶斯优化器。
"""

import math
import torch
import torch.nn as nn


class ExactGPInference:
    """
    高斯过程 (Gaussian Process) 闭式推断引擎
    采用 Matérn 5/2 核函数与噪声方差修正
    """
    def __init__(self, length_scale: float = 1.0, signal_var: float = 1.0, noise_var: float = 1e-4):
        self.length_scale = length_scale
        self.signal_var = signal_var
        self.noise_var = noise_var
        self.Z_train = None
        self.y_train = None
        self.K_inv = None

    def matern52_kernel(self, X1: torch.Tensor, X2: torch.Tensor) -> torch.Tensor:
        """Matérn 5/2 协方差核函数计算"""
        # 计算两两欧氏距离矩阵
        dist = torch.cdist(X1 / self.length_scale, X2 / self.length_scale, p=2)
        sqrt5_r = math.sqrt(5.0) * dist
        k = self.signal_var * (1.0 + sqrt5_r + (5.0 / 3.0) * dist.pow(2)) * torch.exp(-sqrt5_r)
        return k

    def fit(self, Z: torch.Tensor, y: torch.Tensor):
        """记录历史评测点并预计算逆协方差矩阵 (Cholesky 分解)"""
        self.Z_train = Z
        self.y_train = y.view(-1, 1)

        N = Z.shape[0]
        K = self.matern52_kernel(Z, Z) + self.noise_var * torch.eye(N, device=Z.device)
        # 数值稳定的逆矩阵求解 (通过 Cholesky 分解)
        L = torch.linalg.cholesky(K)
        self.K_inv = torch.cholesky_inverse(L)

    def predict(self, Z_test: torch.Tensor):
        """解析推断测试点潜变量的高斯后验均值与方差"""
        K_star = self.matern52_kernel(self.Z_train, Z_test)  # (N_train, N_test)
        K_test = self.matern52_kernel(Z_test, Z_test)        # (N_test, N_test)

        # mu = K_*^T * K^{-1} * y
        mu = torch.matmul(K_star.T, torch.matmul(self.K_inv, self.y_train)).squeeze(-1)

        # var = diag(K_** - K_*^T * K^{-1} * K_*)
        var = torch.diag(K_test - torch.matmul(K_star.T, torch.matmul(self.K_inv, K_star)))
        sigma = torch.sqrt(torch.clamp(var, min=1e-8))

        return mu, sigma


class BayesianOptimizationLoop:
    """
    潜空间贝叶斯优化主控制器
    包含期望提升 (Expected Improvement, EI) 计算与梯度上升极大化采集函数
    """
    def __init__(self, gp_model: ExactGPInference, xi: float = 0.01):
        self.gp = gp_model
        self.xi = xi

    def compute_expected_improvement(self, Z_candidates: torch.Tensor) -> torch.Tensor:
        """计算候选点的 EI 采集分值 (Closed-form Analytical Formulation)"""
        mu, sigma = self.gp.predict(Z_candidates)
        current_best = torch.max(self.gp.y_train)

        diff = mu - current_best - self.xi
        Z_norm = diff / sigma

        # 正态累积分布函数与概率密度函数
        normal = torch.distributions.Normal(0.0, 1.0)
        cdf = normal.cdf(Z_norm)
        pdf = torch.exp(normal.log_prob(Z_norm))

        ei = diff * cdf + sigma * pdf
        # 边界修正: 当不确定度极低且均值低于历史最佳时，EI 归零
        ei = torch.where(sigma < 1e-6, torch.zeros_like(ei), ei)
        return ei

    def suggest_next_candidate(self, latent_dim: int, n_restarts: int = 10,
                               steps: int = 50, lr: float = 0.05) -> torch.Tensor:
        """
        通过多启动梯度上升在潜空间中求解 argmax_z alpha_EI(z)
        """
        best_candidate = None
        highest_ei = -float("inf")

        # 采样多个初始点进行多峰寻优，避免陷入局部鞍点
        init_guesses = torch.randn(n_restarts, latent_dim, requires_grad=True)
        optimizer = torch.optim.Adam([init_guesses], lr=lr)

        for _ in range(steps):
            optimizer.zero_grad()
            # 极大化 EI 等价于极小化 -EI
            ei_values = self.compute_expected_improvement(init_guesses)
            loss = -torch.sum(ei_values)
            loss.backward()
            optimizer.step()

        with torch.no_grad():
            final_eis = self.compute_expected_improvement(init_guesses)
            best_idx = torch.argmax(final_eis)
            best_candidate = init_guesses[best_idx].clone()
            highest_ei = final_eis[best_idx].item()

        return best_candidate, highest_ei


# =====================================================================
# 单元验证模块: 模拟未知物理黑盒函数，验证贝叶斯优化闭环寻优能力
# =====================================================================
if __name__ == "__main__":
    torch.manual_seed(42)

    latent_dim = 8
    # 模拟一个真实的复杂物理响应黑盒目标函数 (例如带隙匹配度)
    def true_physical_blackbox(z: torch.Tensor) -> torch.Tensor:
        # 真实的最优解藏在潜空间 z^* 处，其余区域存在非凸起伏
        z_star = torch.ones(latent_dim) * 0.75
        target_score = 10.0 - torch.norm(z - z_star, dim=-1).pow(2)
        # 叠加微小物理振荡
        target_score += torch.sin(z.sum(dim=-1) * 2.0)
        return target_score

    print(">>> [STEP 1] 初始化低容量历史实验数据集 (样本数 N=8)...")
    z_init = torch.randn(8, latent_dim)
    y_init = true_physical_blackbox(z_init)
    print(f"    初始历史样本中最优性质: {torch.max(y_init):.4f}")

    gp = ExactGPInference(length_scale=1.5, signal_var=2.0, noise_var=1e-3)
    gp.fit(z_init, y_init)

    bo_loop = BayesianOptimizationLoop(gp, xi=0.05)

    print(">>> [STEP 2] 执行 5 轮自适应贝叶斯主动优化闭环迭代...")
    current_Z = z_init
    current_y = y_init

    for iteration in range(1, 6):
        gp.fit(current_Z, current_y)
        next_z, expected_gain = bo_loop.suggest_next_candidate(latent_dim=latent_dim, n_restarts=8)

        # 真实物理黑盒评测新候选
        real_y = true_physical_blackbox(next_z.unsqueeze(0))

        # 更新历史数据集
        current_Z = torch.cat([current_Z, next_z.unsqueeze(0)], dim=0)
        current_y = torch.cat([current_y, real_y], dim=0)

        print(f"    第 {iteration} 轮推荐: 期望增益 EI={expected_gain:.4f}, 真实评测值={real_y.item():.4f}, 迄今最佳={torch.max(current_y):.4f}")

    improvement = torch.max(current_y) - torch.max(y_init)
    assert improvement > 0, "贝叶斯优化未实现性能提升！"
    print(f">>> [SUCCESS] 贝叶斯优化闭环验证全部通过！累计性质突破提升: +{improvement:.4f}")

5. 本章小结

本章系统厘清了在连续潜在流形上实现全局按需逆向设计的数学算法机制：

洞悉梯度上升的致命盲区：揭示了无约束梯度上升极易受代理模型虚假高估绑架，导致潜编码跌入流形空洞产生“对抗假象结构”，明确了流形先验正则化的必要性；
高斯过程的数学优雅：严密推导了基于 Matérn 核的后验均值与认知不确定性解析计算式，展现了 GP 如何在极小样本量下提供严格的“自信度度量”；
期望提升（EI）的自动权衡：深入解构了闭式 EI 采集函数如何同时驱动“在已知好区域深挖（利用）”与“在未知盲区开拓（探索）”。

然而，当我们要处理的结构本身具有极其离散的图构建步骤，或者我们希望以概率直接正比于目标奖励采样出成千上万个截然不同但活性相当的多样化骨架时，纯粹连续流形贝叶斯优化便会遇到算力瓶颈。

下一章我们将迈向更为激进与前沿的离散策略智能体——第 04 章：强化学习与条件引导生成（从 PPO 到 GFlowNets）！

第 03 章：潜在空间优化算法：高斯过程贝叶斯优化与梯度上升 ​

1. 引言：连续流形上的导航艺术 ​

2. 连续代理表面的梯度上升法及其对抗陷阱 ​

2.1 梯度上升的基础形式 ​

2.2 致命陷阱：对抗假象（Adversarial Hallucination / OOD Exploitation） ​

2.3 破局方案：流形支撑惩罚项（Manifold Prior Regularization） ​

3. 高斯过程贝叶斯优化（Gaussian Process Bayesian Optimization, GP-BO） ​

3.1 高斯过程（Gaussian Process, GP）的严密数学推导 ​

3.2 采集函数（Acquisition Functions）推导：平衡利用与探索 ​

1. 期望提升准则 (Expected Improvement, EI) ​

2. 上置信界准则 (Upper Confidence Bound, UCB) ​

4. 工业级代码实战：端到端潜空间高斯过程贝叶斯优化 ​

5. 本章小结 ​