1. 项目概述:当符号回归遇见高维物理流形
在理论物理的前沿,尤其是在全息对偶(AdS/CFT)的研究中,一个核心挑战是理解共形场论(CFT)中连续参数族的存在性,即所谓的“共形流形”。从引力角度看,这对应于超引力标量势中平坦方向(flat directions)的寻找。然而,超引力势通常涉及数十个标量场,其表达式极其复杂,传统的解析方法往往在符号计算的汪洋大海中搁浅。这就好比试图在一片由高维多项式构成的茂密丛林中,仅凭纸笔绘制出所有可能的平坦小径,其计算复杂度是指数级增长的。
近年来,机器学习,特别是符号回归(Symbolic Regression),为这类问题带来了曙光。符号回归的目标不是拟合一个黑箱模型,而是从数据中直接“发现”潜在的、可解释的数学表达式。然而,当目标流形由高维空间中的多项式约束定义时,传统的符号回归方法(如遗传编程或AIFeynman)常常力不从心。它们要么搜索效率低下,要么难以同时处理多个相互关联的约束方程。
这正是我们工作的切入点。我们开发并应用了一种基于退火序贯蒙特卡洛(Annealed Sequential Monte Carlo, ASMC)采样器的符号回归算法,专门用于从高维数值数据中发现定义几何轨迹的多项式约束。本文将详细拆解这一混合数值-符号方法,从超引力背景设定、数据采样策略,到ASMC算法的核心原理与实现细节,最后展示我们如何应用该方法,在一个五维标量子空间中,成功发现了描述一个三维共形流形的多项式方程组。无论你是对理论物理中的数值方法感兴趣的研究者,还是希望将符号回归应用于复杂工程问题的工程师,这篇文章都将为你提供一个从理论到代码的完整实操指南。
2. 理论背景与问题建模:超引力中的平坦方向
2.1 全息对偶与共形流形
AdS/CFT对偶是理论物理中一个深刻而强大的框架,它将反德西特(AdS)时空中的引力理论与边界上的共形场论(CFT)联系起来。在这个框架下,超引力理论作为弦论低能有效描述,其AdS解对应于边界CFT的固定点。一个关键问题是:这些CFT是孤立的,还是属于一个连续的族?如果是后者,连接不同CFT的连续参数空间就称为共形流形。
从引力侧看,边界理论的共形流形对偶于AdS解的一个连续族,其中AdS因子保持不变。如果存在对这些解的一致截断(consistent truncation),那么共形流形也可以被识别为截断后标量势中的平坦方向。沿着这些方向,标量场构型连续变化,而宇宙学常数(即AdS半径的平方的倒数)保持不变。寻找这些平坦方向,等价于求解标量势梯度为零(∇V = 0)的方程组。
2.2 具体超引力模型与计算挑战
我们工作的具体舞台是三维N=8(半极大)规范超引力理论的一个特定截断。该理论描述了六维N=(1,1)超引力在S^3上的紧化,其标量场参数化了陪集空间SO(8,4)/(SO(8)×SO(4))。经过适当的参数化(详见原始论文附录A),我们得到了一个包含13个标量场(记作向量X)的势能函数V(X)。这个势能V的表达式(论文公式2.5)极其冗长,包含了矩阵的迹、行列式以及复杂的张量收缩。
注意:直接对13个变量的非线性方程组∇V=0进行全局解析求解,即使用上最强大的符号计算软件(如Mathematica),也几乎是不可能的。表达式会迅速膨胀到无法简化或求解的地步。但这并不意味着简单的解不存在,只是传统的“暴力”符号方法找不到它们。
我们的策略是转向数值与符号混合的方法。初步的数值探索表明,注意力可以聚焦在五个标量子集上:x = (x₁, x₂, x₄, x̃₈, x₁₀)。为了确保在这个子空间中找到的解同样是完整理论的解,我们必须采用一个关键技巧:先求导,后置零。即,我们先计算完整13维势能V对所有13个变量的梯度∇V,然后再将其他8个变量y = (x₃, x₅, x₆, x̃₇, x̃₉, x₁₁, x₁₂, x̃₁₃)设为零。这样得到关于x的5个方程,其解自动满足完整理论的运动方程。
核心问题因此被重新表述为:从一个无法直接解析求解的复杂系统∇V(X)=0出发,通过数值方法采样其解空间(即x空间中满足∇V|_y=0 ≈ 0的点),然后利用符号回归技术,从这些数值点云中“反推”出定义该解空间的、相对简单的多项式约束方程。
3. 符号回归算法核心:退火序贯蒙特卡洛(ASMC)采样器
面对从高维点云中发现多项式约束的任务,我们需要一个能在巨大、离散且崎岖的符号空间中进行高效搜索的算法。传统的MCMC方法(如Metropolis-Hastings)在这里容易陷入局部最优。我们的解决方案是ASMC,它融合了序贯蒙特卡洛(SMC)的粒子滤波思想和退火重要性采样(AIS)的渐进优化策略。
3.1 算法总览:从粒子演化到多项式发现
ASMC的核心思想是将“寻找湮灭多项式”的问题,转化为从一个特定概率分布中采样的问题。我们定义多项式空间E_pol上的一个概率密度π(z),它正比于exp(-L(z)),其中L(z)是一个损失函数,当多项式z在数据点上的取值接近零时,L(z)很小。这样,从π(z)中采样得到的高概率多项式,就是对我们想要的“湮灭多项式”(即z(x⁽ⁱ⁾)≈0)的良好近似。
然而,直接从π(z)采样是困难的。ASMC通过构造一系列中间分布{π_n},从一个人为设定的、易于采样的初始分布π_0(例如均匀分布)开始,逐渐“变形”到目标分布π。这个过程由一系列逆温度参数{β_n}控制,其中β_0=0,β_nepochs=1。中间分布定义为: π_n(z) ∝ π_0(z) exp(-β_n L(z))
我们可以把β_n想象成一个“专注度”参数。开始时(β=0),分布π_0完全忽略数据,对所有多项式一视同仁。随着β逐渐增加到1,分布π_n越来越“专注”于那些能使损失L(z)变小的多项式,最终在β=1时,我们得到的目标分布π就高度集中于那些能很好拟合数据的多项式上。
ASMC通过维护一组带权重的“粒子”{z⁽ᵏ⁾_n, w⁽ᵏ⁾_n}来近似每个中间分布π_n。每个粒子代表一个候选多项式。算法迭代地进行以下步骤:
- 重加权:根据从π_{n-1}到π_n的变化,更新每个粒子的权重。
- 重采样:当粒子权重退化(即少数粒子占据绝大部分权重)时,根据权重重新采样粒子,淘汰弱粒子,复制强粒子,以保持探索活力。
- 变异:对每个粒子应用马尔可夫链蒙特卡洛(MCMC)移动,在其邻域内探索新的多项式。
最终,当算法结束时,我们得到的一组粒子{z⁽ᵏ⁾}就是从目标分布π的近似采样,其中权重高的粒子对应的多项式,就是我们要找的、能较好描述数据约束的候选表达式。
3.2 算法细节拆解:粒子、移动与平衡
3.2.1 多项式空间与损失函数定义
首先,我们需要界定搜索空间。我们限制搜索最大次数为maxdegree、最多包含maxmon个单项式的多项式。给定变量数n_var,可能单项式的总数是组合数C(maxdegree+n_var, maxdegree)。这个空间虽然是离散的,但依然非常庞大。
我们采用的损失函数L(z)是拟合精度与模型简洁性的权衡: L(z) = Σ_{i=1}^{n_points} [z(x⁽ⁱ⁾)]² + λ Σ_m |c_m|
- 第一项(数据拟合项):多项式在所有数据点上取值的平方和。理想情况下,��于完美的湮灭多项式,此项为零。
- 第二项(L1正则化项):多项式所有系数c_m的绝对值之和。λ是正则化强度超参数(通常~O(10³))。这一项至关重要,它防止算法找到“平凡解”——即所有系数为零的多项式(显然在所有点上都为零,但毫无意义)。L1正则化也倾向于产生稀疏解(即系数为零的项多),这有助于发现更简洁的表达式。
3.2.2 MCMC移动策略:如何“扰动”一个多项式
变异步骤的核心是定义一个前向传播核q(z_n | z_{n-1}),它描述了如何从当前多项式z_{n-1}产生一个新多项式z_n。我们设计了三种基本移动操作,每次迭代随机选择一种执行:
- 系数扰动:随机选择多项式中的一个系数c_m,为其加上一个高斯噪声:c_m → c_m + δ,其中δ ~ N(0, σ²)。例如,
2x₁x₂ + 3x₂² → 2.1x₁x₂ + 3x₂²。这是最细微的调整,用于对已有良好多项式进行微调。 - 变量乘法:随机选择多项式中的一个单项式,将其乘以一个随机选择的变量。例如,
2x₁x₂ + 3x₂² → 2x₁x₂² + 3x₂²。这增加了多项式的复杂度(次数)。 - 变量除法:随机选择多项式中的一个单项式,如果可能,将其除以一个构成它的变量。例如,
2x₁x₂ + 3x₂² → 2x₁ + 3x₂²。这降低了多项式的复杂度。
每种操作被选中的概率由超参数p_shift,p_multiply,p_divide控制。在实际中,系数扰动应占主导(如80%),以保证局部精细搜索;乘法和除法操作则用于跳出局部最优,探索不同复杂度的模型。
3.2.3 退火调度与重采样策略
退火调度(β_n序列)的选择至关重要。如果β增长太快,系统可能“淬火”过快,陷入局部最优;如果增长太慢,则计算成本高昂。一个常见策略是自适应调度:根据当前粒子集的加权有效样本大小(ESS)来调整β的增长幅度,确保粒子多样性不会过快丧失。
重采样是防止粒子退化的关键。我们监控有效样本大小: ESS_n = ( Σ_{k=1}^{n_particles} (w̃⁽ᵏ⁾_n)² )^{-1} 当ESS_n低于阈值(通常设为n_particles/2)时,执行重采样。重采样后,所有粒子权重重置为均匀的1/n_particles,因为分布信息已经体现在粒子的新分布中了。
实操心得:在实现中,重采样的频率需要仔细权衡。过于频繁的重采样会减少粒子多样性,可能导致早熟收敛;而重采样不足则会让计算资源浪费在权重极低的粒子上。一个实用的技巧是,仅在ESS连续几次低于阈值时才触发重采样,并记录不同β阶段找到的“精英粒子”,在最后进行汇总分析。
4. 数值实验全流程:从数据生成到解析发现
有了理论框架和算法,接下来我们进入实战环节,看看如何一步步从超引力势出发,最终得到描述流形的多项式方程。
4.1 第一步:梯度下降采样解流形
我们的起点是超引力势V(X)。目标是找到满足∇V|_y=0 ≈ 0的点x。
- 初始化:在五维超立方体[-2, 2]^5内,均匀随机生成n_points = 10^5个初始点。选择范围超出[-1,1]至关重要,这能防止后续符号回归算法倾向于产生在[-1,1]区间内偶然拟合良好、但外推性差的高次多项式。
- 损失函数定义:我们使用自动微分(如TensorFlow或Jax)计算完整梯度∇V(X),然后将y变量置零,得到关于x的梯度。损失函数定义为所有数据点梯度范数的平方和:L = Σ_i ||∇V(X⁽ⁱ⁾)|_{y⁽ⁱ⁾=0}||²。
- 优化过程:采用Adam优化器。我们发现,周期性重启优化器能显著加速收敛。具体调度如下:
- 初始学习率α=10^{-2},运行250轮。
- 重启优化器,保持α=10^{-2},再运行250轮(总500轮)。
- 再次重启,α=10^{-2},运行至总750轮。
- 重启,将α降至10^{-3},运行至总1000轮。
- 最后重启,α=10^{-4},运行最后500轮至总1500轮。
这种“热身重启”策略类似于学习率计划,但通过重置优化器的内部动量状态,能更有效地逃离平坦区域或次优盆地。如图2所示,每次重启后损失函数都迎来新一轮的快速下降。
- 结果验证:梯度下降完成后,我们检查结果。图3(a)显示,超过93%的数据点梯度范数低于10^{-4},超过99%低于10^{-3},表明我们成功采样到了势能的平坦方向。同时,图3(b)显示所有点的势能值V都收敛到原点值V(0)=-4附近(误差≲10^{-4}),确认这些点确实位于与原始AdS₃×S³解连续相连的流形上。
注意事项:梯度下降采样存在一个潜在风险:优化过程可能使点聚集在某些“吸引力强”的盆地中,而漏掉其他平坦方向。为了缓解此问题,一是需要生成足够多的初始点(我们用了10^5个),二是初始范围要足够大。此外,可以尝试结合随机噪声注入或不同优化器来增加探索的随机性。
4.2 第二步:流形维度与局部结构分析
在将数据扔给符号回归算法之前,我们需要对点云的几何结构有个初步认识。这有助于我们设定符号回归的预期(例如,我们期望找到几个约束方程?)。
- 局部主成分分析(Local PCA):我们不是对整个数据集进行全局PCA,而是在每个数据点的小邻域内进行PCA。计算该邻域内点的协方差矩阵,并分析其特征值谱。如果流形是d维的,那么我们应该观察到d个显著大于零的特征值(对应流形切空间的方向),以及(5-d)个接近零的特征值(对应法空间方向)。
- 聚类与维度统计:对所有数据点邻域进行Local PCA后,统计每个点估计出的局部维度。通过直方图可以发现,绝大多数点的局部维度估计值聚集在3附近。这表明,我们采样的点云很可能嵌入在一个5维空间中的3维流形上。
- 可视化与相关性检查:绘制所有变量两两之间的散点图(如图4所示)。我们可以直观地看到一些变量对(如x₄/x̃₈, x₄/x₁₀)之间存在明显的非线性相关结构,这暗示了约束方程的存在。而另一些(如x₁/x₂)显示的结构,后续分析表明可能是梯度下降采样引入的人为假象。
结论:数值分析强烈提示,存在一个3维的连续族解。因此,我们的目标是从5个坐标(x₁, x₂, x₄, x̃₈, x₁₀)中,找到2个独立的多项式约束方程,将它们限制到一个3维曲面上。
4.3 第三步:ASMC符号回归实战
现在,我们将准备好的数据(10^5个五维点)输入ASMC算法,寻找湮灭多项式。
4.3.1 算法初始化与参数设置
- 粒子数
n_particles: 通常设置在1000-5000之间。更多的粒子能探索更广的空间,但计算成本更高。 - 多项式空间:设定
maxdegree=4(根据物理直觉,势能V是标量场的函数,其平坦方向方程很可能来自低次多项式),maxmon=10(限制多项式复杂度)。 - 损失函数:λ设为2000,以有效惩罚非稀疏解。
- 退火调度:采用几何增长β_n = β_0 * (β_final/β_0)^{n/n_epochs},并辅以ESS自适应调整。
- MCMC移动概率:
p_shift=0.8,p_multiply=0.1,p_divide=0.1。 - 变异步长σ:初��可设为0.1,并随β增加而逐渐减小,以实现从粗搜索到精调。
4.3.2 典型运行分析与后处理
一次ASMC运行结束后,我们得到一批带权重的粒子(多项式)。权重最高的多项式是我们的一批候选解。然而,由于MCMC移动的随机性,这些多项式的系数可能尚未精确优化到使损失函数绝对最小。
因此,需要一个后处理(精调)步骤:对每个候选多项式z,固定其单项式结构(即哪些变量的多少次幂组合出现),仅对其系数c_m执行一个快速的梯度下降优化,以最小化Σ_i [z(x⁽ⁱ⁾)]²。这一步计算量很小,但能显著提升多项式的精度。
经过精调后,我们筛选出那些在数据集上均方根误差(RMSE)低于某个阈值(例如10^{-6})的多项式。然后对这些多项式进行符号简化(使用如SymPy的simplify函数),并检查它们之间的代数关系。
4.3.3 结果统计与解析发现
在我们的实验中,ASMC算法成功发现了8个不同的多项式约束(记为P1到P8),它们都在数据点上近乎为零。并非所有约束都是独立的。通过计算这些多项式理想(ideal)的Gröbner基或进行符号消元,我们最终将它们约化到2个独立的二次多项式方程。
例如,我们可能发现形式如下的约束(此为示意): P₁: x₄² + x̃₈² - x₁₀² ≈ 0 P₂: x₁x̃₈ + x₂x₄ - C ≈ 0 其中C是某个常数。
这些方程共同定义了一个5维空间中的3维代数簇。求解这个方程组(可能还需要考虑物理上的正定性等条件),就能得到流形的显式参数化,或者至少是隐式定义。这比原始的、包含13个变量的庞杂方程组∇V=0要简洁明了得多。
与AIFeynman的对比:我们在附录B中与先进的符号回归工具AIFeynman进行了对比。AIFeynman在寻找单一、复杂的全局表达式时表现出色,但它本质上是一个序列化的、递归的算法,难以同时发现多个相互关联的、定义高维几何对象的约束方程。而我们的ASMC方法通过粒子群同时探索多个候选解,并且其损失函数和正则化项天然适合于发现一组稀疏的多项式,因此在处理此类“多约束、高维几何”问题时更具优势。
5. 发现的新解及其物理意义
通过上述数值-符号混合流程,我们不仅找到了描述流形的多项式约束,更重要的是,这些约束直接对应着超引力理论的新经典解。
每一个满足约束方程组的具体坐标点x,代入原始的标量场参数化公式,就给出了三维超引力理论的一个AdS₃真空解。通过一致截断的 uplift 公式,这些解可以提升到六维,对应于变形的AdS₃ × M₃几何,其中M₃是三维球面S³的一种变形。这些解是连续的,共享相同的宇宙学常数,因此它们构成了边界二维共形场论(CFT₂)中一个共形流形的引力对偶。
我们的方法系统地揭示了这一连续族的存在,并提供了其解析描述(尽管是隐式的)。这超越了之前工作中通过猜测或对称性找到的孤立解或特定子族。它证明,在看似复杂的超引力势中,确实隐藏着连续的平坦方向,这为全息对偶中寻找共形流形提供了一个可重复、可扩展的计算框架。
6. 实操心得、挑战与扩展方向
6.1 踩过的坑与关键技巧
- 数据范围的艺术:符号回归极度依赖数据点的分布。如果所有数据都集中在[-1,1]区间,算法极易找到一些在该区间内巧合成立的高次多项式,这些多项式缺乏泛化能力,也不是我们想要的本质约束。务必让采样点覆盖足够大的范围(例如[-2,2]),以迫使算法发现更鲁棒、更本质的代数关系。
- 正则化系数λ的权衡:λ太小,算法容易收敛到全零系数的平凡解;λ太大,则会过度惩罚非零系数,可能丢失重要的高次项。需要根据数据规模和预期多项式的稀疏程度进行调优。一个策略是从较大的λ开始,如果发现所有系数都被压制,再逐步减小。
- 退火速度与粒子多样性:退火太快(β增长快)会导致早熟收敛,所有粒子迅速聚集到同一个局部最优解附近。退火太慢则效率低下。监控有效样本大小(ESS)是关键。如果ESS下降过快,应放缓退火速度或增加重采样阈值。
- MCMC移动的接受率:接受率是衡量探索效率的指标。理想接受率通常在20%-50%之间。如果接受率过低,说明提议的移动(变异)太大,粒子难以“爬升”到概率更高的区域;如果接受率过高,则说明移动太小,探索不够。可以通过动态调整变异步长σ来管理接受率。
- 后处理不可或缺:ASMC输出的多项式系数通常不是最优的。一个快速的、基于梯度的系数精调步骤,能以极小的计算成本将多项式的拟合误差降低几个数量级。
6.2 算法扩展与应用前景
- 处理更高维与更复杂流形:当前方法针对的是由多项式方程定义的代数流形。对于更一般的解析流形,可以考虑将基函数从单项式扩展到更一般的函数(如三角函数、指数函数),但这会极大增加搜索空间的复杂度。一种折中方案是分阶段进行:先用多项式发现主要结构,再在残差上使用更复杂的基函数。
- 与神经网络结合:可以使用神经网络(如物理信息神经网络PINNs)来更高效、更全局地采样高维势能的临界点集(即∇V=0),然后将这些点作为ASMC的输入。神经网络的全局逼近能力可以弥补梯度下降可能陷入局部盆地的不足。
- 应用于其他物理问题:这套方法论具有普适性。任何可以表述为“在复杂高维函数中寻找满足特定条件(如梯度为零、值为常数)的子流形”的问题,都可以尝试此方法。例如:
- 弦论景观:寻找膜势能中的平坦方向或慢滚区域。
- 凝聚态物理:从量子多体系统的数值数据中寻找序参量的约束关系。
- 流体力学:从模拟数据中发现守恒律或本构关系的近似解析形式。
6.3 代码实现要点
对于想要复现或应用此方法的研究者,这里给出一些核心代码结构的提示:
import jax.numpy as jnp import jax from functools import partial class ASMCSymbolicRegressor: def __init__(self, n_vars, max_degree, max_monomials, n_particles): self.n_vars = n_vars self.max_degree = max_degree self.monomial_lib = self._build_monomial_library() # 构建单项式库 self.n_particles = n_particles self.particles = self._initialize_particles(max_monomials) # 随机初始化粒子 def _loss(self, coeffs, data): """计算一个多项式(由其系数coeffs定义)在数据data上的损失。""" # coeffs: 形状 (n_monomials,),对应monomial_lib中的单项式 # data: 形状 (n_points, n_vars) poly_vals = jnp.dot(self.monomial_lib(data), coeffs) # 计算多项式值 data_fit_term = jnp.sum(poly_vals**2) regularization = self.lambda_reg * jnp.sum(jnp.abs(coeffs)) return data_fit_term + regularization def _mcmc_move(self, particle, beta, key): """对单个粒子执行一次MCMC移动。""" move_type = random.choice(['perturb', 'multiply', 'divide']) if move_type == 'perturb': # 随机扰动一个系数 idx = random.randint(0, len(particle.coeffs)-1) new_coeffs = particle.coeffs.at[idx].add(random.normal(key) * self.sigma) new_monomials = particle.monomials elif move_type == 'multiply': # 为随机一个单项式乘以一个随机变量 # ... 实现逻辑,注意次数不能超过max_degree pass # ... 其他移动 # 计算接受概率 old_loss = self._loss(particle.coeffs, self.data) new_loss = self._loss(new_coeffs, self.data) log_accept_ratio = -beta * (new_loss - old_loss) # 简化的Metropolis-Hastings # 加上提议分布比的对数(对于对称提议,此项常为零) accept = jnp.log(random.uniform(key)) < log_accept_ratio return jax.lax.cond(accept, lambda: new_particle, lambda: particle) def anneal_step(self, particles, weights, beta_old, beta_new, key): """执行一轮退火步骤:重加权、重采样、变异。""" # 1. 重加权 incremental_weights = jnp.exp(-(beta_new - beta_old) * losses) new_weights = weights * incremental_weights new_weights /= jnp.sum(new_weights) # 2. 检查并执行重采样 ess = 1.0 / jnp.sum(new_weights**2) if ess < self.n_particles / 2: indices = resample_multinomial(new_weights, key) particles = particles[indices] new_weights = jnp.ones(self.n_particles) / self.n_particles # 3. 变异所有粒子 keys = random.split(key, self.n_particles) new_particles = jax.vmap(self._mcmc_move)(particles, beta_new, keys) return new_particles, new_weights def fit(self, data, n_epochs): """主训练循环。""" beta_schedule = jnp.linspace(0, 1, n_epochs+1)[1:] # 退火计划 for epoch, beta_new in enumerate(beta_schedule): self.particles, self.weights = self.anneal_step( self.particles, self.weights, self.beta, beta_new, random.PRNGKey(epoch) ) self.beta = beta_new # 返回精调后的最佳多项式 return self._refine_and_select_best()这个框架提供了ASMC符号回归的核心逻辑。在实际应用中,需要仔细实现单项式库的构建、各种MCMC移动操作及其对应的雅可比行列式(用于计算接受率中的提议分布比),以及高效的重采样算法(如系统重采样)。
最后,我想强调的是,这项工作的价值不仅在于为一个特定的超引力模型找到了新的共形流形,更在于展示了一条结合数值采样与符号推理来解决复杂理论物理问题的新路径。它降低了从高维数值数据中提取解析洞察的门槛,让研究者能够更直接地“询问”数据背后的数学结构。随着算法和计算工具的进一步发展,这类混合方法有望在理论物理、应用数学乃至工程领域的更多“数据丰富但解析困难”的场景中发挥重要作用。