论文信息：Zhang, Y., Xu, Y., & Yang, B. (2021). Lévy Gradient Descent: Augmented Random Search for Geophysical Inverse Problems. Surveys in Geophysics, 42, 1109–1132. https://doi.org/10.1007/s10712-021-09644-6

地球物理反演的难点在于：一方面，位场、电磁、地震波等观测的解析都依赖对地下模型的"正演"，而反演往往是非线性、欠定（ill-posed）的；另一方面，如何评估反演结果的可信度——即给出模型参数的不确定度——始终是实际工作中的难题。本文以一作张壹为代表提出了一种新的随机搜索算法 Lévy 梯度下降（L-GD），它把自然界中常见的**莱维飞行（Lévy flights）**思想引入梯度下降框架，试图同时获得传统两类反演方法的优点：梯度法的效率，以及随机/贝叶斯方法的全局收敛性与不确定度估计。下面围绕这篇 Surveys in Geophysics 论文，梳理其动机、方法、数值实验与结论。

研究背景：确定性反演与概率反演的取舍

把地下空间离散成网格，构造"数据拟合 + 模型约束"的目标函数，反演问题就转化为一个最优化问题。按求解算法，主流方法大体可分两类：

  • 确定性反演：基于目标函数对模型参数的梯度的算法，如经典梯度下降、共轭梯度（CG）、L-BFGS 等拟牛顿类方法。它们对连续凸问题极其高效，广泛用于大规模 2D/3D 反演（重力、磁法、MT、地震走时等）。但面对非线性、病态的非凸问题，这类算法容易收敛到局部极小值（或鞍点），结果强烈依赖初始模型；更重要的是，它们通常不提供模型参数的不确定度。若要估计误差，只能在线性化、局部子空间或修改正则化约束等近似下进行，统计有效性常受到质疑。
  • 概率反演（贝叶斯反演）：把模型参数视为概率分布，用 MCMC 等采样技术探索后验分布。它能天然给出误差分布，但计算量往往比梯度法大一到两个数量级，且随参数个数急剧增长，难以直接用于大规模 3D 反演；通常也不直接给出一个"最优模型"。

于是自然产生一个设想：能否设计一种方法，既像梯度法那样高效、可直接用于大规模反演，又具备随机/贝叶斯方法的全局收敛性与不确定度估计能力？本文提出的 L-GD 正是对这一问题的尝试。

方法：用 Lévy 飞行驱动梯度下降

L-GD 的基本思想非常直观：与经典梯度下降一样，每次迭代都沿目标函数梯度给出的下降方向搜索；唯一的区别在于，它不使用固定（或逐渐减小的）步长，而是从 Lévy 分布中随机抽取步长。

莱维飞行是一类特殊的随机游走：路径由大量小步长与少量大步长组成（图1）。这种运动模式广泛存在于动物觅食、流体动力学、光的输运、人类出行乃至地震时空分布等领域。与著名的布朗运动（步长服从正态分布、覆盖面积小）相比，莱维飞行在步数相同时能覆盖大得多的区域，因而对稀疏目标的搜索效率更高。

关键的性质在于：莱维分布的最大值是无穷大的。这意味着即使偶然陷入局部极小值，L-GD 在有限次搜索内也有概率通过一个很大的步长"跳出"局部极值——这正是其全局收敛性的来源。同时，由于搜索过程本质上是一个马尔可夫过程，随着搜索步增长，它会收敛到与初始模型无关的某种后验分布，因此搜索路径的统计量可以用来评估反演模型的不确定度，而不需要额外的后验采样。

算法对外输出三类结果：平均模型 m_mean、标准差 m_sd（即不确定度）以及历史最优模型 m_best（使目标函数最小的那个）；并设有停滞检测（stagnation test）以避免梯度消失导致算法卡在鞍点。为保证不同量纲、不同尺度的参数能同时反演，算法对参数空间做了归一化（尺度因子 p）。

莱维步长采用Mantegna 算法生成，其中指数因子 𝛽 ∈ (1,2) 控制分布形态（𝛽 越大分布越快趋近高斯过程）。数值统计实验表明：𝛽 越大，最大和平均步长越小、路径越稳定；𝛽 越小，越容易产生超大步长，对参数空间的穿越效率越高（表1中最大步长的均值跨了约三个数量级）。因此 𝛽 就像一把"旋钮"，用来在路径稳定与搜索强度之间权衡。加上一个步长尺度因子 𝛼（默认约 0.01–0.1），用于整体缩放步长、调节收敛速度。

数值实验：非凸问题上的全局收敛与效率

为检验算法，作者先用一个人造的二维非凸势函数 U 作测试：它由若干个负二维高斯叠加而成，含两个局部极小值和一个全局极小值 U_min = U(75, 20)。将 L-GD 与 L-BFGS、模拟退火（SA）、MCMC（Metropolis–Hastings）四种算法从同一初值 m₀=(25,20) 出发对比。

四种算法搜索路径对比
图1　四种算法在非凸势函数上的搜索路径（引自论文 Figure 5）

(a) L-BFGS 陷入局部极小值；(b)(c) SA 与 MCMC 的路径在极值附近大范围散布、效率低；(d) L-GD 路径基本沿梯度方向，并依次找到全部三个极值（含全局最小值）。虚线为搜索路径，全局极小值 U_min=U(75,20)=−1.57826e−3。

结果一目了然：L-BFGS 陷入局部极小（其对初始模型的依赖在非凸问题上暴露无遗）；SA 与 MCMC 的路径在极值附近大范围散布、靠近极小值的效率低；而 L-GD 的路径基本沿目标函数的梯度前进，不仅高效逼近极值，还依次探到了全部三个极小值。定量上（表3），L-GD 首次抵达全局极值的步数 N_first 仅 866，比 SA（8085）快约 9 倍、比 MCMC（4429）快约 5 倍；在给定总搜索次数下，L-GD 的成功搜索次数 N_k 也最多。

为何如此高效？作者给出了一个很妙的视角：传统随机搜索靠"概率性接受使目标函数变大的步"来跳出局部极值；L-GD 则不同——它的步长大多相近，仅偶有极大步。极大步天然承担了"扩大搜索区域/逃离局部极值"的角色，统计上相当于其他方法中的拒绝采样，但 L-GD 把这类大跳全部接受了下来，不浪费算力。相比无梯度方法，L-GD 的搜索方向由梯度引导，因此在高维问题中搜索方向的组合数远小于纯随机方法，有望规模扩展到大型反演。

数值实验：重力与地震走时反演

除了非凸测试函数，作者还将 L-GD 用于两类真实尺度的地球物理反演，并与 L-BFGS 对比。

重力界面反演（3D）：用球坐标下由三角面片描述的地形界面，反演一密度界面的起伏（半径 r=1e5 m，密度 1.0 g/cm³）。构造了两组正演模型——一组观测多于未知量（超定），一组欠定且噪声更强。超定情形下，L-GD 反演的 mean/best 地形与 L-BFGS 的结果几乎一致，误差多在 ±5 m 内；更重要的是，L-GD 额外给出了平均地形的标准差（图2 中叠加的等值线），且其分布与真实误差基本吻合（误差主要集中于中央高区）。欠定且含强噪声时，L-BFGS 在中央高区出现明显错误（最高超过 580 m，而真实约为 420 m），而 L-GD 反演结果仍接近真实模型，误差多在 ±25 m 内。

重力界面反演结果
图2　重力界面反演结果（引自论文 Figure 8，超定情形）

(a)(b) L-GD 反演的平均地形与最优地形；(c) 平均地形与真实地形之差，叠加等值线为其估计标准差（不确定度），与真实误差分布基本一致；(d) L-BFGS 反演结果。两法结果几乎一致，但只有 L-GD 能同时给出误差估计。

地震首波走时反演（2D）：用两个三角形网格（一个贴合异常体做正演，一个用于反演）对合成速度结构作反演，走时正演用非结构网格上的快速行进法（FMM）。由于问题欠定，目标函数中加入了对模型平滑性的正则化（权重 𝜇=5，使问题变为强凸）。结果 L-GD 反演的 mean/best 慢度结构均与 L-BFGS 非常接近，高慢度异常体的位置与倾角得到较好恢复。估计标准差（SD）的分布大体与真实误差重合，但在高慢度异常体处，估计 SD 集中在异常体本身而非其两端。作者据此强调了一个重要提醒：L-GD 给出的 SD 本质上是搜索过程所得模型参数的统计量，只有当反演模型足够接近真实解时，它才代表真实误差分布，因此在多局部极小值的情形下其不确定性估计精度会下降——它更多是定性的参考，而非严格的定量误差。

结论与讨论

综合来看，L-GD 在两类反演（重力、地震走时，分属线性/非线性）中都能以可接受的耗时得到与梯度法相当的精度，同时具备：

  1. 全局收敛性：非凸问题上明显优于纯梯度法，且对初始模型不敏感；
  2. 较高的搜索效率：由于方向由梯度引导，所需搜索步数远少于无梯度随机方法，且所需步数更多取决于问题的非线性程度而非未知参数个数，因此可扩展到大型反演；
  3. 顺带的不确定度估计：通过对搜索路径的统计给出模型参数的估计误差分布，虽本质上是定性的，但对评估反演结果、指导实际工作很有价值。

当然，方法的局限也需正视：对凸问题，L-GD 的收敛步数通常仍远多于 L-BFGS 等纯梯度法；估计标准差与真实误差的一致性依赖于反演模型的准确性与问题的非线性程度。作者的下一步设想颇为开放：因为 L-GD 实际只需要梯度的方向信息而非精确梯度，因此有望通过近似梯度或次梯度把它推广到更广的、甚至不连续可定义的问题上，并尝试用于大地电磁、地震波形等更多反演场景。

对读者而言，这是一篇把"自然界高效觅食策略"数学化并落地到地球物理反演的方法论文，其"全局搜索 + 降本 + 顺带误差估计"的组合，对我们处理非线性、病态反演时如何平衡效率与可信度，很有启发。

延伸阅读与引用

  • Pavlyukevich, I. (2007). Lévy flights, non-local search and simulated annealing. Journal of Computational Physics, 226(2), 1830–1844.
  • Mantegna, R. N., & Stanley, H. E. (1994). Stochastic process with ultraslow convergence to a Gaussian: the truncated Lévy flight. Physical Review Letters, 73(22), 2946–2949.
  • Humphries, N. E., et al. (2010). Environmental context explains Lévy and Brownian movement patterns of marine predators. Nature, 465(7301), 1066–1069.
  • Zhang, Y., Mooney, W. D., Chen, C., & Du, J. (2019). Interface inversion of gravitational data using spherical triangular tessellation. Geophysical Journal International, 217(1), 703–713.

本文为对该论文的中文解读，图片均引自论文原文。