1. 从“猜”到“算”为什么我们需要牛顿法在工程优化、机器学习模型训练乃至金融衍生品定价中我们常常会遇到一个核心问题寻找一个函数的最小值或最大值。比如你想调整一个神经网络的参数让它的预测误差最小或者你想找到一个机械结构的最优设计使其在满足强度要求下重量最轻。这本质上都是在求解一个优化问题。最朴素的想法是什么是“猜”或者叫“试错”。比如你想找函数f(x) x²的最小值你可能会在数轴上随机选几个点比如x1x-2x0.5然后比较它们的函数值。你会发现x0.5时f(x)0.25比较小于是你可能会在0.5附近再试试。这个过程效率很低尤其是在参数维度成百上千的高维空间里无异于大海捞针。于是我们引入了“梯度”的概念。梯度指向了函数值增长最快的方向那么它的反方向自然就是函数值下降最快的方向。这就是梯度下降法的核心思想沿着当前点的负梯度方向走一小步不断迭代期望能走到最低点。这个方法直观、简单是深度学习的基石之一。但梯度下降法有个明显的缺点它只利用了函数的一阶信息梯度就像一个蒙着眼睛下山的人只知道脚下哪个方向最陡但不知道这个陡坡会持续多久也不知道前面是不是悬崖。因此它常常要走很多“之”字形的弯路收敛速度慢尤其是在接近最优解谷底的平坦区域步长会变得非常小进展缓慢。有没有更“聪明”的方法有的。如果我们不仅能知道脚下的坡度一阶导数还能知道坡度的变化率二阶导数即曲率我们就能对函数的局部形状有一个更准确的“预测”。我们知道一个光滑函数在一点附近可以用泰勒展开来近似。梯度下降法只用了泰勒展开的一阶项而牛顿法则豪爽地用到了二阶项。它通过构建一个局部二次模型来近似原函数并直接跳到这个二次模型的极小值点。理论上这能让迭代次数大幅减少达到“二次收敛”的速度——这意味着每迭代一步误差的平方会缩小一个常数倍是一种非常快的收敛。然而天下没有免费的午餐。牛顿法虽然快但代价高昂它需要计算目标函数的海森矩阵Hessian Matrix即所有二阶偏导数组成的矩阵。对于有n个参数的函数海森矩阵是一个n×n的矩阵。计算它本身就需要O(n²)的量级更别提后续还需要求解一个n维的线性方程组计算复杂度约O(n³)。当n很大时现代机器学习动辄数百万甚至上亿参数这完全是不可承受之重。此外牛顿法要求海森矩阵必须是正定的才能保证找到的是极小值点但这个条件在实际中并不总能满足。于是数学家们开始思考我们能否找到一种方法既保留牛顿法快速收敛的精髓又避免计算和存储庞大的海森矩阵这个问题的答案就是一系列拟牛顿法的诞生。它们的目标是构造一个海森矩阵或其逆矩阵的近似矩阵这个近似矩阵能随着迭代不断更新以较低的成本模拟出海森矩阵所包含的曲率信息。所以从梯度下降一阶到牛顿法二阶再到拟牛顿法近似二阶这条技术演进路径的核心驱动力始终是在计算效率与收敛速度之间寻找最佳平衡点。接下来我们就深入这个“平衡的艺术”内部看看各种方法是如何具体工作的以及在实际中我们该如何选择和驾驭它们。2. 牛顿法理想很丰满现实很骨感牛顿法为我们描绘了一个快速找到最优解的蓝图但当我们真正动手实现时会发现蓝图和施工图之间隔着不少坑。理解这些细节是决定我们能否用好牛顿法以及何时该转向拟牛顿法的关键。2.1 核心推导从泰勒展开到迭代公式我们从一个无约束优化问题开始min f(x)其中x是一个n维向量f是二阶可微的函数。在当前位置x_k我们对f(x)做二阶泰勒展开f(x) ≈ f(x_k) ∇f(x_k)^T (x - x_k) 1/2 (x - x_k)^T H(x_k) (x - x_k)其中∇f(x_k)是梯度向量一阶导数H(x_k)是海森矩阵二阶导数。牛顿法的思想是直接找到这个二次近似函数的极小值点作为下一个迭代点。对于二次函数其极小值点满足导数为零。我们对上面的近似式关于x求导并令其为零∇f(x_k) H(x_k)(x - x_k) 0如果海森矩阵H(x_k)是可逆的我们就能解出下一步的迭代点x_{k1} x_k - [H(x_k)]^{-1} ∇f(x_k)这就是经典的牛顿法迭代公式。其中- [H(x_k)]^{-1} ∇f(x_k)被称为牛顿方向和牛顿步长的乘积。它与梯度下降方向-∇f(x_k)有着本质区别梯度下降方向只考虑了局部最速下降而牛顿方向还通过海森矩阵的逆进行了“缩放”和“旋转”使得搜索方向能直接指向当前二次模型的最低点。注意这里有一个关键点。在纯牛顿法中步长是隐含在[H]^{-1}∇f这个量里的我们通常认为步长为1。但在实际应用中为了确保全局收敛性即从任意初始点出发都能收敛我们常常会在牛顿方向前加一个步长因子α_k通过线搜索来确定即x_{k1} x_k - α_k [H(x_k)]^{-1} ∇f(x_k)。这被称为阻尼牛顿法或线搜索牛顿法它比纯牛顿法更鲁棒。2.2 海森矩阵计算、存储与正定性之殇牛顿法最大的瓶颈毫无悬念地落在了海森矩阵H上。1. 计算成本高昂 对于由用户显式给出表达式的函数我们可以通过符号微分或自动微分来精确计算二阶导数。但对于复杂函数尤其是那些由大量代码组成的函数如一个深度神经网络的前向传播计算完整的海森矩阵在计算上是灾难性的。其复杂度是O(n²)对于n10000这就是一亿个二阶导数。2. 存储成本高昂 海森矩阵是n×n的需要O(n²)的存储空间。在n很大的情况下这可能会直接耗尽内存。例如一个百万参数n1e6的模型其海森矩阵理论上需要1e12个元素的存储空间这远远超出了普通计算机的容量。3. 正定性要求苛刻 牛顿法要求海森矩阵H(x_k)是正定矩阵这样才能保证牛顿方向是下降方向即- [H]^{-1}∇f与负梯度-∇f的点积大于0。但在非凸函数的优化中特别是在鞍点附近海森矩阵很可能不是正定的甚至是不定的。此时牛顿方向可能指向一个上升方向导致算法发散。4. 求解线性方程组的开销 即使我们有了海森矩阵每一步迭代还需要求解线性方程组H(x_k) d -∇f(x_k)来得到搜索方向d。对于稠密矩阵求解的典型复杂度是O(n³)这比计算海森矩阵本身还要昂贵。正是这些“骨感”的现实限制了标准牛顿法在大规模优化问题中的应用。它更像一个理论上的“黄金标准”揭示了利用二阶信息所能达到的速度上限但在实践中我们往往需要它的各种“平价替代品”。2.3 实战中的变通与技巧尽管有诸多限制牛顿法在小规模、黑塞矩阵易于计算且性质良好的问题上依然是无冕之王。在实际使用中我们会有一些变通使用自动微分框架对于定义在像 PyTorch、JAX 这样的框架中的函数可以利用其高阶自动微分功能来计算海森向量积Hessian-vector product而无需显式构造整个矩阵。这对于某些需要海森矩阵信息的算法如信任域法很有用但依然无法完全规避O(n²)级别的信息需求。处理非正定问题当遇到非正定海森矩阵时常见的修正方法有Levenberg-Marquardt 修正求解(H μI) d -∇f其中μ是一个正数I是单位矩阵。通过增加μ可以迫使矩阵(H μI)正定。这实质上是将牛顿法与梯度下降法做了折中μ很大时方向趋近于负梯度方向μ很小时方向趋近于牛顿方向。修改 Cholesky 分解在尝试对H进行 Cholesky 分解H LL^T要求H正定失败时可以动态地增加对角线元素直到分解成功。这被称为“修正乔列斯基分解”。聚焦“牛顿-高斯”法在非线性最小二乘问题中形式为min Σ [r_i(x)]²目标函数的海森矩阵有一个特殊的近似H ≈ J^T J其中J是残差r_i(x)的雅可比矩阵。这个近似矩阵总是半正定的避免了正定性问题。使用这个近似的牛顿法特称为高斯-牛顿法。如果再加上 Levenberg-Marquardt 修正来处理J^T J可能奇异的情况就成了非常强大的Levenberg-Marquardt 算法在曲线拟合、Bundle Adjustment 等问题中应用极广。尽管如此对于绝大多数大规模机器学习问题标准牛顿法依然遥不可及。这时拟牛顿法的舞台灯光便亮了起来。3. 拟牛顿法用一阶信息“拼凑”二阶智能拟牛顿法的核心哲学非常巧妙既然直接计算海森矩阵H或其逆B H^{-1}太贵那我们就不直接算它。我们利用迭代过程中“免费”获得的一阶信息梯度来逐步构建和更新一个对H或B的近似矩阵。这个近似矩阵需要满足一些由微分中值定理导出的基本性质。3.1 拟牛顿条件连接两点的桥梁设第k步我们有一个对海森矩阵逆B_k ≈ [H(x_k)]^{-1}的近似。我们走到了新点x_{k1} x_k s_k其中s_k是位移并计算了新点的梯度g_{k1} ∇f(x_{k1})。记梯度变化为y_k g_{k1} - g_k。对于光滑函数由中值定理可知存在某个ξ使得y_k H(ξ) s_k。如果我们假设海森矩阵在x_k和x_{k1}之间变化不大那么一个合理的期望是我们的近似逆海森矩阵B_{k1}应该满足B_{k1} y_k ≈ s_k或者如果我们近似的是海森矩阵本身H_k ≈ H(x_k)则应满足H_{k1} s_k ≈ y_k这两个等式被称为拟牛顿条件或割线方程。它意味着用我们更新的矩阵作用在梯度变化上应该能近似地还原出位移。这是拟牛顿法更新公式需要满足的核心约束。3.2 百花齐放的更新公式DFP、BFGS 与 L-BFGS如何在满足拟牛顿条件的前提下从B_k更新到B_{k1}不同的更新规则催生了不同的算法。目标是让更新尽可能“小”即B_{k1}与B_{k1}在某种度量下差异最小同时满足拟牛顿条件。1. DFP 方法Davidon–Fletcher–Powell这是历史上第一个拟牛顿法。它直接对逆海森矩阵B进行更新。其更新公式为B_{k1}^{DFP} (I - ρ_k s_k y_k^T) B_k (I - ρ_k y_k s_k^T) ρ_k s_k s_k^T其中ρ_k 1 / (y_k^T s_k)。 DFP 公式具有很好的数学性质但在实际数值计算中有时会表现出数值不稳定性尤其是在线搜索不精确的情况下。2. BFGS 方法Broyden–Fletcher–Goldfarb–ShannoBFGS 可以说是当今最著名、应用最广泛的拟牛顿法。与 DFP 对偶它选择直接更新海森矩阵H的近似H_k然后通过求逆或求解方程组来得到搜索方向。其更新公式为H_{k1}^{BFGS} H_k - (H_k s_k s_k^T H_k) / (s_k^T H_k s_k) (y_k y_k^T) / (y_k^T s_k)对应的逆矩阵的更新公式为B_{k1}^{BFGS} (I - ρ_k s_k y_k^T) B_k (I - ρ_k y_k s_k^T) ρ_k s_k s_k^T等等这个B的更新公式看起来和 DFP 一模一样仔细观察中间项DFP 公式中是B_k而 BFGS 的逆矩阵更新公式中中间项是(I - ρ_k y_k s_k^T) B_k (I - ρ_k s_k y_k^T)这里需要澄清上面给出的B_{k1}^{BFGS}公式实际上是BFGS 方法应用于逆海森矩阵B的更新公式。而经典的 BFGS 是先更新H再求逆。通过 Sherman–Morrison–Woodbury 公式可以推导出B的更新公式为B_{k1} B_k \frac{(s_k - B_k y_k)(s_k - B_k y_k)^T}{(s_k - B_k y_k)^T y_k}但更常见的是另一种等价的、数值更稳定的形式。为了避免混淆我们记住核心BFGS 通常指更新海森矩阵H的那个版本并且它比 DFP 在实践中更稳定、性能更好因此成为了事实上的标准。3. L-BFGS大规模优化的救星BFGS 虽然避免了计算海森矩阵但它仍需存储和更新一个n×n的稠密矩阵H_k或B_k对于大规模问题n很大存储依然是瓶颈。L-BFGSLimited-memory BFGS解决了这个问题。它的思想非常聪明不存储完整的n×n矩阵而是只保存最近m步通常m在 5 到 20 之间的位移s_k和梯度变化y_k这两个向量序列。在需要计算搜索方向d -B_k g_k时它利用这两组向量通过一个巧妙的“双循环递归”算法直接计算出矩阵向量积B_k g_k而无需显式构造B_k。 这个算法的计算复杂度是O(mn)存储开销是O(mn)。当m远小于n时这相比O(n²)是一个巨大的节约。L-BFGS 在保持 BFGS 优良收敛性质的同时完美适配了大规模问题成为机器学习、深度学习优化库如 SciPy, PyTorch, TensorFlow中的标配优化器之一。为了更直观地对比这几种方法的核心特点我整理了下表特性梯度下降法牛顿法BFGSL-BFGS利用信息一阶梯度一阶梯度 二阶海森矩阵一阶梯度近似二阶一阶梯度近似二阶收敛速度线性收敛局部二次收敛超线性收敛超线性收敛每步计算成本O(n)O(n²)~O(n³)O(n²)O(mn)m很小存储成本O(n)O(n²)O(n²)O(mn)关键瓶颈收敛慢需调步长海森矩阵计算/存储/求逆稠密矩阵存储 (n²)需要选择记忆长度m适用场景大规模、初始阶段、简单模型小规模、海森矩阵易得且正定中等规模问题大规模优化问题主流选择3.3 拟牛顿法的实战心得与坑点在实际项目中应用 L-BFGS我积累了一些在官方文档里不一定写得那么直白的经验初始矩阵H0的选择L-BFGS 需要提供一个初始的逆海森矩阵近似B0。通常的做法是将其设为标量矩阵γI。一个经典且有效的启发式方法是取γ (y_k^T s_k) / (y_k^T y_k)即用最近一步的信息来估计初始尺度。在 SciPy 或 PyTorch 中这通常对应initial_hess_inv或line_search_fn相关的参数设置。记忆长度m不是越大越好增大m保留了更多的曲率信息理论上能让近似更准收敛更快。但更大的m也意味着更多的存储和每步O(mn)计算中的常数因子增大。对于许多问题m取 10-20 已经能获得很好的效果继续增加带来的收益很小。有时对于非常病态的问题适当增加m比如到 50可能有帮助。线搜索的重要性拟牛顿法包括 L-BFGS通常假设使用的步长α_k1即纯牛顿步是合适的。但在实际中为了确保全局收敛性和稳定性往往需要搭配一个线搜索过程。一个强 Wolfe 条件的线搜索对于 L-BFGS 的稳定运行至关重要。它要求步长不仅能使函数值充分下降Armijo 条件还要使梯度的变化满足一定关系这有助于保证生成的s_k和y_k满足y_k^T s_k 0这是维持近似矩阵正定性的关键。如果你发现 L-BFGS 不稳定或发散首先检查你的线搜索是否足够稳健。梯度的精度拟牛顿法通过梯度差分y_k来估计曲率。如果梯度计算本身就不精确例如使用数值差分近似或者在随机优化中梯度噪声很大那么构建的曲率信息就是垃圾算法性能会严重下降甚至失败。确保你的梯度计算是精确的使用自动微分或噪声是可控的在随机优化中采用大批次。不适合所有问题对于高度非凸、存在大量平坦区域或悬崖梯度突变的函数拟牛顿法可能不如一些自适应学习率的梯度下降法如 Adam表现好。Adam 等算法为每个参数维护不同的学习率对噪声和病态条件更鲁棒在深度学习训练的早期阶段和复杂损失地形中常常是更安全的选择。4. 超越经典现代优化中的牛顿类方法变体牛顿法和拟牛顿法的思想影响深远催生了许多适应特定场景或解决特定痛点的变体算法。了解它们能帮助我们在面对具体问题时有更丰富的工具箱。4.1 共轭梯度法另一种利用历史信息的方式共轭梯度法最初是为求解大型线性方程组Axb其中A对称正定而设计的。在优化语境下对于二次函数f(x) 1/2 x^T A x - b^T x其海森矩阵是常数矩阵A。共轭梯度法通过构造一系列相互共轭的搜索方向可以在最多n步内找到极小值。将其推广到非线性优化就得到了非线性共轭梯度法如 Fletcher-Reeves, Polak-Ribière 等方法。它不需要存储矩阵只需要当前和上一步的梯度信息存储开销仅为O(n)。它的核心思想是新的搜索方向d_{k1}由当前负梯度-g_{k1}和上一步方向d_k的线性组合构成即d_{k1} -g_{k1} β_k d_k。不同的β_k计算公式对应不同的变体。共轭梯度法可以看作是一种特殊的、记忆长度仅为 1 的拟牛顿法。它比最速下降法快但通常比 BFGS 慢。它的优势在于极低的内存占用适用于那些n极大以至于连存储几组向量如 L-BFGS 需要的O(mn)都困难的问题例如在某些偏微分方程数值求解产生的大规模优化中。4.2 随机拟牛顿法拥抱数据噪声在机器学习中目标函数通常是大量样本损失函数的和f(x) (1/N) Σ_{i1}^N f_i(x)。计算完整梯度全批量的成本是O(Nn)当N很大时不可行。随机梯度下降SGD通过使用单个或小批次样本的梯度来估计真实梯度将每次迭代成本降为O(Bn)其中B是批次大小。那么能否将拟牛顿法的思想与随机梯度结合起来这就是随机拟牛顿法的研究领域。其挑战在于梯度噪声小批次估计的梯度g_k噪声很大导致根据y_k g_{k1} - g_k估计的曲率信息y_k被噪声严重污染。方差控制需要设计机制来降低海森矩阵近似的方差。常见的思路包括重叠采样计算s_k和y_k时使用更大的批次甚至全批量以获取更准确的曲率信息。但这增加了计算成本。正则化更新在更新公式中引入正则化项防止近似矩阵因噪声而变得病态。方差缩减技术结合 SVRG、SAGA 等方差缩减方法先得到一个低方差的梯度估计再用它来更新拟牛顿矩阵。随机拟牛顿法如 oLBFGS在逻辑回归、线性支持向量机等中等规模问题上表现出色能够显著加快收敛。但对于非常深的神经网络由于其损失曲面高度非凸且噪声复杂自适应学习率方法Adam 等的鲁棒性和简便性使其更受欢迎随机拟牛顿法在这里的应用相对较少。4.3 海森向量积不求矩阵只求作用在许多二阶优化算法中我们并不需要完整的海森矩阵只需要知道海森矩阵乘以一个任意向量v的结果即海森向量积Hv。例如在牛顿共轭梯度法中我们需要反复计算Hp其中p是共轭方向来迭代求解牛顿方程Hd -g。神奇的是对于由自动微分框架定义的函数计算Hv的成本与计算梯度g的成本是同一个量级O(n)远低于计算完整的HO(n²)。这可以通过“前向模式自动微分”或更常用的“反向模式自动微分结合梯度点乘”即先计算梯度g ∇f(x)再对g·v求导来实现。在 PyTorch 中你可以使用torch.autograd.functional.vhp(vector-Hessian product) 来计算在 JAX 中则有jax.jvp或jax.vjp等组合。利用海森向量积我们可以实现牛顿共轭梯度法。它内部使用共轭梯度法来迭代求解牛顿方程Hd -g。由于共轭梯度法本身只需要矩阵向量积Hp而不需要显式的H这就完美地结合了起来。这种方法适用于那些海森矩阵稀疏或可以通过快速矩阵向量积来求解的问题在某些中等规模的机器学习模型训练中它能比 L-BFGS 收敛得更快但实现起来也更复杂。5. 如何选择从理论到实践的决策指南面对这么多方法在实际项目中到底该怎么选没有放之四海而皆准的答案但可以遵循一个清晰的决策流程。第一步评估问题规模与特性参数数量nn很小1000几乎可以随意选择。标准牛顿法如果海森矩阵易得、BFGS 都是优秀选择。可以尝试牛顿法看能否快速收敛。n中等1000 ~ 10^5L-BFGS 是主力。其O(mn)的存储和计算成本可以接受。n很大10^5 尤其是 10^6内存成为首要约束。首先考虑 L-BFGSm设小点如5。如果内存依然紧张考虑非线性共轭梯度法存储O(n)。对于深度学习首选通常是自适应学习率的梯度方法Adam AdamW因为它们对初始学习率不敏感在早期阶段进展快。目标函数性质凸/强凸拟牛顿法L-BFGS和共轭梯度法表现会很好通常能稳定收敛到全局最优。高度非凸如神经网络损失曲面复杂存在大量鞍点、平坦区。此时自适应学习率方法Adam的鲁棒性往往更好。它们能自动调整每个参数的学习率更容易逃离鞍点和平坦区。牛顿类方法在非凸区域可能因为海森矩阵不定而失败。实践中Adam 及其变体如 AdamW 带权重衰减修正是训练深度网络的事实标准起步器。最小二乘形式优先考虑Levenberg-Marquardt或高斯-牛顿法它们是为这类问题量身定做的。随机性/大数据目标函数是大量样本的求和。对于中等规模凸问题如逻辑回归可以尝试随机拟牛顿法oSVRG-LBFGS。对于大规模非凸问题深度学习Adam或带动量的SGD是更主流、更稳定的选择。第二步考量计算与存储资源梯度计算成本如果计算一次完整梯度的代价极高例如需要运行一次昂贵的仿真那么采用收敛速度更快的二阶或拟二阶方法如L-BFGS是值得的因为它能用更少的迭代次数即更少的梯度计算达到所需精度。内存限制这是最硬的约束。务必检查O(n²)的矩阵是否能存下。存不下就立刻排除标准牛顿法和BFGS。在O(mn)的 L-BFGS 和O(n)的共轭梯度法或 SGD 之间做选择。第三步实现复杂性与调参成本易用性Adam胜出。它几乎不需要调学习率默认的 0.001 在很多问题上都能工作对参数初始化不敏感是“开箱即用”的首选。L-BFGS需要实现或使用一个稳健的线搜索强 Wolfe 条件并且对初始步长、记忆长度m、初始逆海森估计等参数有一定敏感度。虽然 SciPy 等库提供了很好的实现但在随机优化或非凸问题中仍需小心调试。牛顿法实现最复杂需要计算海森矩阵或海森向量积并处理可能出现的非正定问题。除非问题有特殊结构如最小二乘否则通常不作为首选。一个实用的推荐流程默认起点对于大多数未知的深度学习任务先用Adam或AdamW跑起来。它快速、鲁棒能给你一个不错的基线。中期调优如果模型收敛后你想进一步微调以期达到可能更优的解或者训练一个凸/较浅的模型可以尝试切换到L-BFGS。使用 Adam 训练得到的参数作为 L-BFGS 的初始点并启用强 Wolfe 线搜索。这常常能带来最后一点精度的提升。特殊问题对于已知是凸且中等规模的问题如逻辑回归、线性 SVM、某些物理仿真直接使用L-BFGS它通常比一阶方法快得多。资源极限当参数多到 L-BFGS 的内存也成问题时回到带动量的 SGD或Adam它们是处理超大规模参数的基石。在我个人的很多项目中尤其是在训练传统的机器学习模型如 GBDT 不适合的线性模型、矩阵分解时L-BFGS 是我的首选它的收敛速度令人满意。而在探索新的神经网络结构时Adam 则提供了无与伦比的稳定性和开发速度。理解牛顿法与拟牛顿法不仅让你多了一些工具更重要的是让你明白了优化算法设计中的权衡艺术——如何在信息、计算、存储之间取得平衡从而在面对具体问题时能做出更明智、更高效的选择。