Part I · 第 6 章 优化 · 6.6Part I · Ch. 6 Optimization · 6.6

Newton、Gauss-Newton 与 LMNewton, Gauss-Newton & LM

梯度下降只用"坡度",每一步都像在雾里摸索;如果还知道曲率,就能直接猜出谷底在哪里。Newton 法用 Hessian 拟合一个抛物面,一步跳到它的最低点。而机器人里的优化问题大多是最小二乘:误差的平方和。利用它的特殊结构,Gauss-Newton 只需要一阶导数 \(J\),就能得到近似的 Newton 步;加一点阻尼,就成了 Levenberg–Marquardt,几乎所有标定、SLAM 和数值逆运动学求解器的核心。 Gradient descent uses only the "slope", feeling its way as if in fog; if you also know the curvature you can guess where the valley floor is. Newton's method fits a paraboloid using the Hessian and jumps to its lowest point in one step. Most robotics optimization problems are least squares, sums of squared errors. Exploiting that structure, Gauss-Newton gets an approximate Newton step from first derivatives \(J\) alone; add a little damping and you have Levenberg–Marquardt, the core of nearly every calibration, SLAM and numerical inverse-kinematics solver.

约 30 分钟~30 min 前置:6.3 梯度下降 · 6.5 线搜索与信赖域 · 3.3 雅可比 · 1.13 最小二乘Prereqs: 6.3 Gradient descent · 6.5 Line search & trust region · 3.3 Jacobian · 1.13 Least squares 5 个交互5 interactives

    1为什么机器人需要它Why robots need this

    让你的手臂末端到达一个位置需要解逆运动学;让相机和激光雷达对齐需要外参标定;让几万个路标和几千个位姿一致需要光束法平差 (bundle adjustment) 与位姿图优化。这些问题的代价函数都长一个样:很多项"预测减测量"的平方和,即非线性最小二乘。

    用 6.3 节的梯度下降解它们,往往要几千步,而且步长要小心地调。用本节的 Gauss-Newton / LM,常常 5 到 20 步就到了小数点后十位。这就是为什么 Ceres、g2o、GTSAM 这些机器人领域的标准优化库,默认的求解器都是 LM:用一阶导数的代价,换来近似二阶的收敛速度。

    Bringing the arm's gripper to a position needs inverse kinematics; aligning a camera and a lidar needs extrinsic calibration; making tens of thousands of landmarks and thousands of poses consistent needs bundle adjustment and pose-graph optimization. All of these have the same kind of cost: many "prediction minus measurement" terms, squared and summed, that is, nonlinear least squares.

    Solving them with the gradient descent of Lesson 6.3 often takes thousands of steps and a carefully tuned step size. With the Gauss-Newton / LM of this lesson, 5 to 20 steps often reach ten decimal places. That is why the standard robotics libraries Ceres, g2o and GTSAM default to LM: first-derivative cost, nearly second-order convergence speed.

    2直觉:三种"看得更远"的办法Intuition: three ways to see further

    Newton:拟合抛物面。在当前点 \(\x_k\),用函数值、梯度和 Hessian 拟合一个二次函数(泰勒展开到二阶),然后直接跳到这个二次函数的最低点。如果函数本来就是二次的,一步到位。

    Gauss-Newton:拉直残差。最小二乘的代价是 \(\tfrac12\|\mathbf r(\x)\|^2\),里面的 \(\mathbf r\) 是弯的。不去拟合整个代价,而是把残差在 \(\x_k\) 处线性化 \(\mathbf r(\x_k+\mathbf p)\approx\mathbf r+J\mathbf p\),这时代价变成一个关于 \(\mathbf p\) 的普通线性最小二乘(1.13 节),一步解出。它省掉了二阶导数。

    Levenberg–Marquardt:给 Gauss-Newton 拴一根绳。线性化只在附近可信(6.5 节),所以对步长加一个惩罚 \(\lambda\|\mathbf p\|^2\):\(\lambda\) 小,就是 Gauss-Newton;\(\lambda\) 大,就退化成一小步梯度下降。算法根据"实际下降/预测下降"自动调节 \(\lambda\),就是 6.5 节信赖域的另一种写法。

    Newton: fit a paraboloid. At the current point \(\x_k\), use the value, gradient and Hessian to fit a quadratic (the Taylor expansion to second order) and jump straight to that quadratic's lowest point. If the function is itself quadratic, one step does it.

    Gauss-Newton: straighten the residuals. A least-squares cost is \(\tfrac12\|\mathbf r(\x)\|^2\) and the \(\mathbf r\) inside is curved. Instead of fitting the whole cost, linearize the residual at \(\x_k\): \(\mathbf r(\x_k+\mathbf p)\approx\mathbf r+J\mathbf p\). The cost then becomes an ordinary linear least-squares problem in \(\mathbf p\) (Lesson 1.13), solved in one go. It spares you the second derivatives.

    Levenberg–Marquardt: put Gauss-Newton on a leash. The linearization is trustworthy only nearby (Lesson 6.5), so penalize the step with \(\lambda\|\mathbf p\|^2\): small \(\lambda\) is Gauss-Newton; large \(\lambda\) degrades into a short gradient-descent step. The algorithm tunes \(\lambda\) from the ratio "actual over predicted decrease", which is just another way of writing the trust region of Lesson 6.5.

    3动手玩:一维的 Newton 法Play: Newton's method in one dimension

    拖动黄点(当前位置 \(x_k\))。紫色抛物线是 \(f\) 在 \(x_k\) 处的二阶泰勒模型 \(m(x)=f+f'(x-x_k)+\tfrac12f''(x-x_k)^2\),青色竖线落在它的顶点上:这就是 Newton 步 \(x_{k+1}=x_k-f'/f''\)。小箭头是同一点上学习率 0.3 的梯度下降一步,用来对比。按"Newton 一步"让它走,下方的表格实时显示误差。

    Drag the amber dot (the current position \(x_k\)). The purple parabola is the second-order Taylor model \(m(x)=f+f'(x-x_k)+\tfrac12f''(x-x_k)^2\) of \(f\) at \(x_k\), and the cyan vertical line lands on its vertex: that is the Newton step \(x_{k+1}=x_k-f'/f''\). The small arrow is one gradient-descent step with rate 0.3 from the same point, for comparison. Press "Newton step" to walk; the table below shows the error live.

    FIG 1用抛物线拟合再跳到谷底Fit a parabola, jump to its vertex
    函数:\(0.5x^2+0.1x^4\)(好脾气的凸函数);\(\log\cosh x\)(凸,但 Newton 会在远处发散);双谷 \(\tfrac14(x^2-1)^2+0.15x\)(中间 \(f''<0\));\(\tfrac14x^4\)(最小点处 \(f''=0\))。表中 \(e_k=|x_k-x^\star|\),\(x^\star\) 由同一方法从当前起点迭代到收敛得到。 Functions: \(0.5x^2+0.1x^4\) (a well-behaved convex function); \(\log\cosh x\) (convex, yet Newton diverges far from the minimum); the double well \(\tfrac14(x^2-1)^2+0.15x\) (with \(f''<0\) in the middle); \(\tfrac14x^4\) (\(f''=0\) at the minimum). In the table \(e_k=|x_k-x^\star|\), with \(x^\star\) obtained by iterating the same method from the current start to convergence.
    试一试Try this
    • 在"碗 + 四次项"上,从 \(x_0=2\) 出发连按"Newton 一步":误差的位数每步大约翻一倍(二次收敛),不到 7 步就到机器精度。On "bowl + quartic" press "Newton step" repeatedly from \(x_0=2\): the number of correct digits roughly doubles each step (quadratic convergence), reaching machine precision in under 7 steps.
    • 在 \(\log\cosh x\) 上把起点放到 \(|x|>1.2\),按 Newton 一步:它跳得比原来更远,反复迭代就发散。(这个函数是凸的!)On \(\log\cosh x\) start with \(|x|>1.2\) and press Newton: it jumps farther away than before and diverges on repetition. (The function is convex!)
    • 在双谷函数上把 \(x_k\) 放到 \(f''<0\) 的区域(红色背景),按 Newton 一步:它走向抛物线的最高点,\(f\) 反而升高。On the double well put \(x_k\) in a region with \(f''<0\) (red background) and press Newton: it walks to the parabola's highest point and \(f\) goes up.
    • 在 \(\tfrac14x^4\) 上从 \(x_0=2\) 开始迭代:误差每步只乘以 \(2/3\),只有线性收敛,因为最小点处曲率为零。On \(\tfrac14x^4\) iterate from \(x_0=2\): the error is multiplied by \(2/3\) each step, only linear convergence, because the curvature at the minimum is zero.
    核心洞见Key insight

    Newton 法的步长自动正好:曲率大的地方(碗窄)步子短,曲率小的地方(碗宽)步子长,不需要调学习率。代价是它完全相信二阶模型:模型准确时快得惊人,不准时(离得远,或 \(f''\le0\))会走得比梯度下降还糟。所以 6.5 节的线搜索或信赖域对 Newton 法是标配。

    Newton's step length is automatically just right: short where curvature is large (narrow bowl), long where it is small (wide bowl); there is no learning rate to tune. The price is that it fully believes the second-order model: astonishingly fast when the model is accurate, worse than gradient descent when it is not (far away, or \(f''\le0\)). That is why the line search or trust region of Lesson 6.5 is standard equipment for Newton's method.

    4二维:Hessian 椭圆与 Newton 步Two dimensions: the Hessian ellipse and the Newton step

    在二维里,二阶模型 \(m(\mathbf p)=f+\mathbf g^\top\mathbf p+\tfrac12\mathbf p^\top H\mathbf p\) 的等高线是以 Newton 点 \(\x_k-H^{-1}\mathbf g\) 为中心的椭圆(\(H\) 正定时)。Newton 步指向椭圆中心;梯度下降的方向 \(-\mathbf g\) 垂直于等高线,对扁椭圆严重偏离中心(6.3 节的"之"字)。\(H^{-1}\) 恰好把椭圆"压圆",这就是 Newton 法对条件数不敏感的原因。

    拖动起点,看紫色模型等高线与真实等高线(青色);"二次函数"上两者完全重合;"Rosenbrock"里重合只在一小块区域;"双谷"里 \(H\) 会变成不定,模型的"中心"是鞍点,Newton 步反而爬坡。

    In two dimensions the contours of the second-order model \(m(\mathbf p)=f+\mathbf g^\top\mathbf p+\tfrac12\mathbf p^\top H\mathbf p\) are ellipses centered at the Newton point \(\x_k-H^{-1}\mathbf g\) (when \(H\) is positive definite). The Newton step points at the ellipse center, while the gradient direction \(-\mathbf g\) is perpendicular to the contour and badly off-center for flat ellipses (the zig-zag of Lesson 6.3). \(H^{-1}\) "squashes the ellipse round", which is why Newton's method does not care about the condition number.

    Drag the start point and compare the purple model contours with the true ones (cyan): on the "quadratic" they coincide exactly; on "Rosenbrock" only in a small patch; on the "two minima" terrain \(H\) becomes indefinite, the model's "center" is a saddle, and the Newton step climbs.

    FIG 2Newton 模型与 Newton 步The Newton model and the Newton step
    黄色空心点是 \(\x_k\)(可拖动);紫色等高线是二阶模型;紫色实心点是 Newton 点(\(H\succ0\) 时是模型最低点;\(H\) 不定时是鞍点,画为空心方块);青色小箭头是梯度下降方向(长度已归一化)。"Newton 一步"把 \(\x_k\) 移到 Newton 点,留下轨迹。 The hollow amber dot is \(\x_k\) (draggable); the purple contours are the second-order model; the filled purple dot is the Newton point (the model's lowest point when \(H\succ0\); a saddle when \(H\) is indefinite, drawn as a hollow square); the small cyan arrow is the gradient-descent direction (normalized). "Newton step" moves \(\x_k\) to the Newton point, leaving a trail.
    • 在"二次函数"上把条件数拉到 50 以上,按"Newton 一步":不管多扁,一步就落在最小点。On the "quadratic" raise the condition number above 50 and press "Newton step": however flat, it lands on the minimum in one step.
    • 在"双谷"里找到一个起点,使 Newton 步让 \(f\) 升高(读数里会标红)。On "two minima" find a start where the Newton step raises \(f\) (the readout turns red).

    5最小二乘:Gauss-Newton 省掉了什么Least squares: what Gauss-Newton leaves out

    最小二乘 \(f=\tfrac12\|\mathbf r(\x)\|^2\) 的梯度和 Hessian 有特殊结构(推导见"再看数学"):

    The gradient and Hessian of a least-squares cost \(f=\tfrac12\|\mathbf r(\x)\|^2\) have a special structure (derived in "The math"):

    $$ \nabla f=J^\top\mathbf r,\qquad \nabla^2f=\hl{gn2}{\cw{J^\top J}}+\textstyle\sum_i r_i\,\nabla^2r_i $$

    第一项 \(J^\top J\) 只需要一阶导数,而且永远半正定;第二项需要每个残差的二阶导数,笨重,且可能不定。Gauss-Newton 就是扔掉第二项。什么时候这样做是合理的?当残差 \(r_i\) 本身很小(模型拟合得好,解附近 \(\mathbf r\approx0\)),或者残差关于参数几乎是线性的(\(\nabla^2r_i\approx0\))。机器人里的好数据通常满足前者。

    The first term \(J^\top J\) needs only first derivatives and is always positive semi-definite; the second needs the second derivatives of every residual, is bulky and may be indefinite. Gauss-Newton simply drops the second term. When is that justified? When the residuals \(r_i\) are themselves small (the model fits well, \(\mathbf r\approx0\) near the solution), or when the residuals are nearly linear in the parameters (\(\nabla^2r_i\approx0\)). Good robotics data usually satisfies the former.

    下面拟合一条指数曲线 \(y=a\,e^{bt}\)(8 个测量点,真值 \(a=2,\ b=-1.2\),加噪声 \(\sigma\))。四个求解器从同一点出发:GD + 回溯、Newton + 回溯(Hessian 不正定时加 \(\mu I\) 使之正定)、Gauss-Newton + 回溯、LM。拖动 \(\sigma\) 改变残差大小:残差越大,Gauss-Newton 丢掉的第二项越重要,它的优势越小。每条曲线都是真实运行的结果。

    Below we fit an exponential \(y=a\,e^{bt}\) (8 data points, truth \(a=2,\ b=-1.2\), noise \(\sigma\)). Four solvers start from the same point: GD + backtracking, Newton + backtracking (add \(\mu I\) to make the Hessian positive definite when it is not), Gauss-Newton + backtracking and LM. Drag \(\sigma\) to change the residual size: the larger the residuals, the more the dropped second term matters and the smaller Gauss-Newton's advantage. Every curve is a real run.

    FIG 3曲线拟合:四个求解器Curve fitting: four solvers
    左:数据(带描边的中性色圆点)、真实曲线(虚线)和各求解器当前迭代的拟合曲线。右:参数 \((a,b)\) 平面的代价地形与四条迭代路径。下:到最优解的距离 \(\|\x_k-\x^\star\|\)(对数)。读数里"第二项占比"是在最优解处 \(\|\sum r_i\nabla^2r_i\|_F/\|J^\top J\|_F\),即 Gauss-Newton 近似的相对误差。 Left: the data (outlined neutral-coloured dots),the true curve (dashed) and each solver's current fit. Right: the cost landscape in the \((a,b)\) plane with the four paths. Bottom: distance to the optimum \(\|\x_k-\x^\star\|\) (log). The "second-term share" in the readout is \(\|\sum r_i\nabla^2r_i\|_F/\|J^\top J\|_F\) at the optimum, the relative error of the Gauss-Newton approximation.
    • 令 \(\sigma=0\)(零残差问题):Gauss-Newton 与 Newton 一样快(都是 6 步左右),GD 要上百步。Set \(\sigma=0\) (a zero-residual problem): Gauss-Newton is as fast as Newton (about 6 steps each) while GD needs over a hundred.
    • 把 \(\sigma\) 拉到 0.25 左右:Gauss-Newton(和 LM)变慢,Newton 仍然只要 5、6 步。第二项开始起作用了。Raise \(\sigma\) to about 0.25: Gauss-Newton (and LM) slow down while Newton still needs only 5 or 6 steps. The second term has started to matter.

    6LM 的阻尼旋钮 \(\lambda\)The LM damping knob \(\lambda\)

    Levenberg–Marquardt 的一步是解 \((J^\top J+\lambda I)\mathbf p=-J^\top\mathbf r\)。\(\lambda\to0\) 时就是 Gauss-Newton;\(\lambda\to\infty\) 时 \(\mathbf p\approx-J^\top\mathbf r/\lambda\),是沿负梯度的一小步。扫过所有 \(\lambda\),\(\mathbf p(\lambda)\) 画出一条从 Gauss-Newton 点到"原地不动"的曲线,路上每一点都是某个信赖域半径下的最优步(6.5 节)。下面在 Rosenbrock(作为最小二乘,\(\mathbf r=[1-x,\ 10(y-x^2)]\))上拖动起点和 \(\lambda\)。

    One Levenberg–Marquardt step solves \((J^\top J+\lambda I)\mathbf p=-J^\top\mathbf r\). As \(\lambda\to0\) this is Gauss-Newton; as \(\lambda\to\infty\), \(\mathbf p\approx-J^\top\mathbf r/\lambda\), a short step along the negative gradient. Sweeping all \(\lambda\), \(\mathbf p(\lambda)\) traces a curve from the Gauss-Newton point to "stand still", and every point on it is the optimal step for some trust-region radius (Lesson 6.5). Below, on Rosenbrock (as least squares, \(\mathbf r=[1-x,\ 10(y-x^2)]\)), drag the start point and \(\lambda\).

    FIG 4阻尼曲线 \(\mathbf p(\lambda)\)The damping curve \(\mathbf p(\lambda)\)
    f(x+p) 随 λ 的变化f(x+p) vs. λ
    紫色曲线是 \(\mathbf p(\lambda)\)(\(\lambda\) 从 \(10^{-3}\) 到 \(10^{5}\)),黄点是当前 \(\lambda\) 对应的一步;空心方块是 Gauss-Newton 点(\(\lambda=0\));青色短线是负梯度方向。右下图:\(f(\x+\mathbf p(\lambda))\) 的变化,星号是使 \(f\) 最低的 \(\lambda\),虚线是当前的 \(f(\x)\):在虚线以上说明这一步让 \(f\) 变差。 The purple curve is \(\mathbf p(\lambda)\) (\(\lambda\) from \(10^{-3}\) to \(10^{5}\)); the amber dot is the step for the current \(\lambda\); the hollow square is the Gauss-Newton point (\(\lambda=0\)); the short cyan line is the negative-gradient direction. Lower right: \(f(\x+\mathbf p(\lambda))\); the star marks the \(\lambda\) giving the lowest \(f\), and the dashed line is the current \(f(\x)\): above it means the step makes \(f\) worse.
    • 把 \(\lambda\) 拉到最小:黄点与空心方块(Gauss-Newton 点)重合。在 Rosenbrock 起点附近,这一步是好是坏?Pull \(\lambda\) to its minimum: the amber dot coincides with the hollow square (the Gauss-Newton point). Near the Rosenbrock start, is that step good or bad?
    • 把 \(\lambda\) 拉到最大:步子变得极短,方向与青色短线(负梯度)一致。Pull \(\lambda\) to its maximum: the step becomes tiny and lines up with the cyan stub (the negative gradient).
    • 把 \(\lambda\) 调到星号附近,让 \(f(\x+\mathbf p)\) 接近这条曲线上的最小值。再拖动起点:最佳 \(\lambda\) 是不是总在变?(所以 LM 要自适应地调它。)Tune \(\lambda\) near the star so \(f(\x+\mathbf p)\) approaches the minimum along the curve. Drag the start: does the best \(\lambda\) keep changing? (That is why LM adapts it on the fly.)

    7再看数学The math

    提示:悬停在带底色的符号上,FIG 3 里对应的求解器路径会被加粗。Tip: hover the highlighted symbols; the matching solver paths in FIG 3 are thickened.

    Newton 法Newton's method

    把 \(f\) 在 \(\x_k\) 处展开到二阶,\(f(\x_k+\mathbf p)\approx f_k+\mathbf g_k^\top\mathbf p+\tfrac12\mathbf p^\top H_k\mathbf p\)。若 \(H_k\succ0\),右边的最小点由梯度为零给出:

    Expand \(f\) to second order at \(\x_k\): \(f(\x_k+\mathbf p)\approx f_k+\mathbf g_k^\top\mathbf p+\tfrac12\mathbf p^\top H_k\mathbf p\). If \(H_k\succ0\), the minimizer of the right side follows from setting its gradient to zero:

    $$ H_k\,\mathbf p_k=-\mathbf g_k\qquad\Longrightarrow\qquad \hl{newton}{\x_{k+1}=\x_k-H_k^{-1}\mathbf g_k} $$

    性质。(1) 在二次函数上一步到位。(2) 二次收敛:若 \(H(\x^\star)\succ0\),\(\nabla^2f\) 是 Lipschitz 的,且 \(\x_0\) 足够近,则 \(\|\x_{k+1}-\x^\star\|\le C\|\x_k-\x^\star\|^2\),正确的位数每步大约翻倍(FIG 1)。(3) 仿射不变:把变量做可逆线性变换,Newton 迭代不变,所以条件数无关。(4) 代价:每步要算(\(n\times n\) 的)Hessian 并解线性方程组,\(O(n^3)\);\(H\) 不正定时方向不一定是下降方向;离解远时可能发散。常见补救:线搜索或信赖域(6.5 节)、给 \(H\) 加 \(\mu I\) 使其正定、拟牛顿法(BFGS)用梯度差近似 \(H\)。

    Properties. (1) One step on a quadratic. (2) Quadratic convergence: if \(H(\x^\star)\succ0\), \(\nabla^2f\) is Lipschitz and \(\x_0\) is close enough, then \(\|\x_{k+1}-\x^\star\|\le C\|\x_k-\x^\star\|^2\); the number of correct digits roughly doubles each step (FIG 1). (3) Affine invariance: an invertible linear change of variables leaves the Newton iteration unchanged, hence indifference to the condition number. (4) Costs: each step forms the \(n\times n\) Hessian and solves a linear system, \(O(n^3)\); when \(H\) is not positive definite the direction need not be a descent direction; far from the solution it can diverge. Standard remedies: line search or trust region (Lesson 6.5), adding \(\mu I\) to make \(H\) positive definite, and quasi-Newton methods (BFGS) that approximate \(H\) from gradient differences.

    最小二乘与 Gauss-NewtonLeast squares and Gauss-Newton

    设残差 \(\mathbf r:\R^n\to\R^m\),雅可比 \(J=\partial\mathbf r/\partial\x\in\R^{m\times n}\)(3.3 节),\(f=\tfrac12\mathbf r^\top\mathbf r=\tfrac12\sum_ir_i^2\)。逐项求导:\(\partial_jf=\sum_ir_i\,\partial_jr_i\),即 \(\nabla f=J^\top\mathbf r\);再求一次:\(\partial_{jk}f=\sum_i\big(\partial_jr_i\,\partial_kr_i+r_i\,\partial_{jk}r_i\big)\),即

    Let the residual be \(\mathbf r:\R^n\to\R^m\) with Jacobian \(J=\partial\mathbf r/\partial\x\in\R^{m\times n}\) (Lesson 3.3) and \(f=\tfrac12\mathbf r^\top\mathbf r=\tfrac12\sum_ir_i^2\). Differentiating term by term: \(\partial_jf=\sum_ir_i\,\partial_jr_i\), that is \(\nabla f=J^\top\mathbf r\); once more: \(\partial_{jk}f=\sum_i\big(\partial_jr_i\,\partial_kr_i+r_i\,\partial_{jk}r_i\big)\), that is

    $$ \nabla^2f=J^\top J+\sum_{i=1}^mr_i\,\nabla^2r_i $$

    丢掉第二项,Newton 方程 \(H\mathbf p=-\mathbf g\) 变成正规方程:

    Dropping the second term, the Newton equation \(H\mathbf p=-\mathbf g\) becomes the normal equations:

    $$ \hl{gn2}{\cw{(J^\top J)\,\mathbf p=-J^\top\mathbf r}}\qquad\Longleftrightarrow\qquad \mathbf p=\arg\min_{\mathbf p}\ \tfrac12\big\|\mathbf r+J\mathbf p\big\|^2 $$

    右边的等价说明:Gauss-Newton 一步就是对线性化的残差做线性最小二乘(1.13 节),数值上可以用 QR 或 SVD 解,避免把条件数平方。它的性质:(1) 只要 \(J\) 列满秩,\(J^\top J\succ0\),\(\mathbf p\) 一定是下降方向(\(\mathbf g^\top\mathbf p=-\|J\mathbf p\|^2<0\))。(2) 零残差问题(\(\mathbf r(\x^\star)=0\))上,被丢掉的项在解处恰为零,GN 与 Newton 同样二次收敛;残差小时接近二次,残差大或曲率大时,只有线性收敛,甚至不收敛(FIG 3 滑块)。(3) \(J\) 接近秩亏(例如机械臂的奇异位形,FIG 5)时,\(J^\top J\) 近乎奇异,\(\mathbf p\) 巨大而不可靠。

    The equivalence on the right says a Gauss-Newton step is a linear least-squares solve of the linearized residual (Lesson 1.13), which can be done with QR or SVD to avoid squaring the condition number. Properties: (1) as long as \(J\) has full column rank, \(J^\top J\succ0\) and \(\mathbf p\) is a descent direction (\(\mathbf g^\top\mathbf p=-\|J\mathbf p\|^2<0\)). (2) On zero-residual problems (\(\mathbf r(\x^\star)=0\)) the dropped term vanishes at the solution, and GN converges quadratically just like Newton; with small residuals it is nearly quadratic, with large residuals or strong curvature only linear, or not at all (the slider in FIG 3). (3) When \(J\) is close to rank-deficient (for instance an arm at a kinematic singularity, FIG 5), \(J^\top J\) is nearly singular and \(\mathbf p\) is huge and unreliable.

    Levenberg–MarquardtLevenberg–Marquardt

    $$ \hl{lm2}{\ck{(J^\top J+\lambda I)\,\mathbf p=-J^\top\mathbf r}}\ ,\qquad \rho=\frac{f(\x)-f(\x+\mathbf p)}{\tfrac12\mathbf p^\top(\lambda\mathbf p-\mathbf g)}\qquad\Longrightarrow\qquad \begin{cases}\rho>0:\ \text{accept},\ \lambda\leftarrow\lambda\cdot\max\!\big(\tfrac13,\,1-(2\rho-1)^3\big)\\ \rho\le0:\ \text{reject},\ \lambda\leftarrow2\lambda\ (\text{and double the factor on repeated failure})\end{cases} $$

    分母是模型预测的下降量(由 \((J^\top J+\lambda I)\mathbf p=-\mathbf g\) 可得 \(m(0)-m(\mathbf p)=\tfrac12\mathbf p^\top(\lambda\mathbf p-\mathbf g)\)),所以 \(\rho\) 正是 6.5 节的"实际/预测"。上面的更新规则来自 Madsen–Nielsen:\(\rho\) 接近 1 时 \(\lambda\) 缩小到原来的 \(1/3\),\(\rho=\tfrac12\) 时不变,\(\rho\) 接近 0(勉强接受)时几乎翻倍;失败时乘以 \(\nu\)(初值 2),连续失败时 \(\nu\) 也翻倍。这个 \(\lambda\) 与信赖域半径一一对应:\(\lambda\) 越大,\(\|\mathbf p\|\) 越小。Marquardt 的原版把 \(\lambda I\) 换成 \(\lambda\,\mathrm{diag}(J^\top J)\),使阻尼随各参数的尺度缩放,参数单位悬殊时更稳健(6.3 节误区 2)。

    三者关系。Newton:\(H=J^\top J+\sum r_i\nabla^2r_i\),无阻尼;GN:丢掉第二项;LM:GN + \(\lambda I\)。再往后,6.9 节会把同样的 LM 用在由几千个残差块组成的因子图上,只是 \(J\) 变成稀疏矩阵。

    The denominator is the predicted decrease of the model (from \((J^\top J+\lambda I)\mathbf p=-\mathbf g\) one gets \(m(0)-m(\mathbf p)=\tfrac12\mathbf p^\top(\lambda\mathbf p-\mathbf g)\)), so \(\rho\) is exactly the "actual over predicted" of Lesson 6.5. The update rule is from Madsen–Nielsen: for \(\rho\) near 1, \(\lambda\) shrinks to \(1/3\); at \(\rho=\tfrac12\) it is unchanged; for \(\rho\) near 0 (barely accepted) it almost doubles; on failure it is multiplied by \(\nu\) (initially 2), and \(\nu\) itself doubles on repeated failure. This \(\lambda\) corresponds one-to-one to a trust-region radius: the larger \(\lambda\), the smaller \(\|\mathbf p\|\). Marquardt's original replaces \(\lambda I\) with \(\lambda\,\mathrm{diag}(J^\top J)\), so the damping scales with each parameter, which is sturdier when the parameters have wildly different units (Pitfall 2 of Lesson 6.3).

    How the three relate. Newton: \(H=J^\top J+\sum r_i\nabla^2r_i\), undamped; GN: the second term dropped; LM: GN plus \(\lambda I\). Lesson 6.9 later runs the same LM on factor graphs made of thousands of residual blocks, where \(J\) simply becomes a sparse matrix.

    8机器人中的它:奇异位形附近的逆运动学In robotics: inverse kinematics near a singularity

    两连杆臂(连杆长 1.0、0.8)的 IK 就是最小二乘:残差 \(\mathbf r(\mathbf q)=\mathrm{FK}(\mathbf q)-\mathbf t\),雅可比 \(J\) 就是机器人雅可比(3.3 节)。Gauss-Newton 的一步 \(\mathbf p=-J^{-1}\mathbf r\)(这里 \(J\) 是方阵)正是经典的"雅可比逆"IK。但当手臂几乎伸直(\(q_2\approx0\))时,\(J\) 近乎奇异:手端没法沿手臂方向瞬时移动,要求它移动,\(J^{-1}\) 就给出巨大的关节转角。

    下面让手臂从几乎伸直的姿态出发去够一个更近的目标。GN(无阻尼)的第一步就是十几弧度的疯狂转动,之后还要扑腾很多步;LM 因为 \(\lambda\) 限制了步长,平稳地收敛;GD + 回溯也安全,但更慢。拖动目标和关节空间里的起点自己试试。

    IK of the two-link arm (link lengths 1.0 and 0.8) is least squares: the residual \(\mathbf r(\mathbf q)=\mathrm{FK}(\mathbf q)-\mathbf t\), and its Jacobian \(J\) is the robot Jacobian (Lesson 3.3). The Gauss-Newton step \(\mathbf p=-J^{-1}\mathbf r\) (here \(J\) is square) is the classic "Jacobian-inverse" IK. But when the arm is nearly straight (\(q_2\approx0\)), \(J\) is nearly singular: the gripper cannot instantaneously move along the arm, and asking it to makes \(J^{-1}\) return enormous joint rotations.

    Below the arm starts nearly straight and reaches for a nearer target. GN (undamped) makes a wild turn of more than ten radians on its first step and flounders for many more; LM, whose \(\lambda\) limits the step, converges smoothly; GD + backtracking is also safe but slower. Drag the target and the start in joint space to try it yourself.

    FIG 5奇异位形附近的 IK:GN 对 LMIK near a singularity: GN versus LM
    左上:工作空间(拖动目标);手臂画的是 LM 当前的姿态,细线是 GN 的末端轨迹。右上:关节空间中三个求解器的路径(拖动圆圈改变起点)。左下:末端误差 \(\|\mathbf r\|\)(对数)。右下:每次迭代的关节位移 \(\|\mathbf p\|\)(rad,对数),GN 的首步高达十几弧度。LM 的 \(\lambda_0=0.3\),按 \(\rho\) 自适应;GN 无任何保护。 Top left: the workspace (drag the target); the arm drawn is LM's current posture and the thin line is GN's end-effector trail. Top right: the three solvers' paths in joint space (drag the circle to change the start). Bottom left: the end-effector error \(\|\mathbf r\|\) (log). Bottom right: joint displacement \(\|\mathbf p\|\) per iteration (rad, log); GN's first step is over ten radians. LM starts at \(\lambda_0=0.3\) and adapts with \(\rho\); GN has no protection at all.
    • 读出 GN 的第一步有多大(右下图和读数),比较 LM 的第一步。Read how large GN's first step is (bottom right and the readout) and compare LM's.
    • 把关节空间里的起点拖到 \(q_2\approx1.5\)(弯曲的、远离奇异的姿态):GN 的表现是否好多了?Drag the joint-space start to \(q_2\approx1.5\) (a bent posture, far from the singularity): does GN do much better?
    真实系统里In real systems

    这就是为什么工业级 IK(TRAC-IK、KDL 的 ChainIkSolverVel_wdls、MoveIt 的数值求解器)几乎总用阻尼最小二乘:\(\mathbf p=-(J^\top J+\lambda I)^{-1}J^\top\mathbf r\) 或等价的 \(J^\top(JJ^\top+\lambda I)^{-1}\mathbf r\)。在奇异点附近,它牺牲一点精度,换来有界的关节速度,不会"炸"。同样的 LM 也是视觉 SLAM 里 PnP 位姿优化、光束法平差的默认引擎;相机标定(OpenCV 的 calibrateCamera)内部跑的就是 LM。

    That is why industrial IK (TRAC-IK, KDL's ChainIkSolverVel_wdls, MoveIt's numerical solvers) nearly always uses damped least squares: \(\mathbf p=-(J^\top J+\lambda I)^{-1}J^\top\mathbf r\) or the equivalent \(J^\top(JJ^\top+\lambda I)^{-1}\mathbf r\). Near a singularity it trades a little accuracy for bounded joint velocities, and does not "explode". The same LM is also the default engine of PnP pose optimization and bundle adjustment in visual SLAM; camera calibration (OpenCV's calibrateCamera) runs LM internally.

    9常见误区Common pitfalls

    误区 1:Newton 法"更高级",所以总比梯度下降好Pitfall 1: Newton is "more advanced", so always better than gradient descent

    Newton 法每步要算和分解 Hessian,对 \(n\) 为百万级的神经网络完全不可行;在非凸区域它会走向鞍点或山顶(FIG 1、FIG 2);离解远时还可能发散。局部快,全局不一定稳:因此一定要配线搜索或信赖域。

    Newton needs to form and factor the Hessian at every step, which is hopeless for networks with millions of parameters; in non-convex regions it heads for saddles or hilltops (FIG 1, FIG 2) and far from the solution it can diverge. Fast locally is not stable globally: always pair it with a line search or trust region.

    误区 2:Gauss-Newton 总是收敛到最小点Pitfall 2: Gauss-Newton always converges to the minimum

    它丢掉了 \(\sum r_i\nabla^2r_i\)。残差大(模型错了,数据里有离群点)时,GN 慢得多,甚至绕圈。若还没有线搜索,纯 GN 的步长可能大到把你甩出线性化有效的区域。实用的 GN 总是带线搜索或阻尼,即 LM。离群点的处理是另一个故事:用鲁棒核函数 (Huber、Cauchy),在 6.9 节展开。

    It drops \(\sum r_i\nabla^2r_i\). With large residuals (a wrong model, outliers in the data) GN is much slower and may even circle. Without a line search, a pure GN step can be long enough to throw you out of the region where the linearization holds. Practical GN always carries a line search or damping, that is, LM. Handling outliers is another story: robust kernels (Huber, Cauchy), taken up in Lesson 6.9.

    误区 3:直接解正规方程 \((J^\top J)\mathbf p=-J^\top\mathbf r\)Pitfall 3: solving the normal equations \((J^\top J)\mathbf p=-J^\top\mathbf r\) directly

    \(J^\top J\) 的条件数是 \(J\) 的平方(1.12 节:奇异值平方)。\(J\) 的条件数 \(10^4\) 时 \(J^\top J\) 就到了 \(10^8\),双精度也勉强。更稳的做法是对 \(J\) 做 QR 分解,或对增广矩阵 \(\begin{bmatrix}J\\\sqrt\lambda I\end{bmatrix}\) 做最小二乘。LM 的 \(\lambda\) 也同时起到"保证可逆"的作用。

    The condition number of \(J^\top J\) is the square of \(J\)'s (Lesson 1.12: squared singular values). A \(J\) with condition number \(10^4\) gives \(10^8\) for \(J^\top J\), already marginal even in double precision. More stable is a QR factorization of \(J\), or a least-squares solve of the augmented matrix \(\begin{bmatrix}J\\\sqrt\lambda I\end{bmatrix}\). The LM \(\lambda\) also guarantees invertibility.

    误区 4:\(\lambda\) 设成一个固定常数Pitfall 4: setting \(\lambda\) to a fixed constant

    FIG 4 已经表明,最佳 \(\lambda\) 随位置变化几个数量级。固定 \(\lambda\) 的"阻尼最小二乘"只适合对精度要求不高的实时场景(例如速度级 IK 里的 \(\lambda\approx0.01\)–\(0.1\))。需要收敛到高精度的离线问题,要用 \(\rho\) 自适应 \(\lambda\)。

    FIG 4 showed that the best \(\lambda\) varies by orders of magnitude with position. A fixed-\(\lambda\) "damped least squares" suits only real-time settings where accuracy is not critical (for instance \(\lambda\approx0.01\)–\(0.1\) in velocity-level IK). Offline problems needing high accuracy must adapt \(\lambda\) with \(\rho\).

    10小练习Check yourself

    1. 对二次函数 \(f(\x)=\tfrac12\x^\top A\x-\mathbf b^\top\x\)(\(A\succ0\)),从任意点出发,纯 Newton 法需要几步到达最小点?1. For a quadratic \(f(\x)=\tfrac12\x^\top A\x-\mathbf b^\top\x\) (\(A\succ0\)), how many pure Newton steps from an arbitrary start reach the minimum?
    二阶泰勒展开对二次函数是精确的:\(\mathbf g=A\x-\mathbf b\),\(H=A\),\(\x-A^{-1}(A\x-\mathbf b)=A^{-1}\mathbf b\),正是最小点,与起点和 \(\kappa\) 都无关(FIG 2 的"二次函数")。\(n\) 步的说法属于共轭梯度法。The second-order Taylor expansion of a quadratic is exact: \(\mathbf g=A\x-\mathbf b\), \(H=A\), \(\x-A^{-1}(A\x-\mathbf b)=A^{-1}\mathbf b\), exactly the minimizer, independent of the start and of \(\kappa\) (the "quadratic" of FIG 2). The "\(n\) steps" claim belongs to conjugate gradients.
    2. 要用 Gauss-Newton 求 \(x^2=2\) 的根,即最小化 \(f=\tfrac12r^2\),\(r(x)=x^2-2\)。从 \(x_0=2\) 出发,\(x_1\) 是多少?2. Use Gauss-Newton to solve \(x^2=2\), i.e. minimize \(f=\tfrac12r^2\) with \(r(x)=x^2-2\). From \(x_0=2\), what is \(x_1\)?
    \(r(2)=2\),\(J=r'=2x=4\)。GN:\(J^2p=-Jr\Rightarrow p=-r/J=-2/4=-0.5\),\(x_1=1.5\)(这正是求平方根的古老的 Babylonian 方法)。对比完整的 Newton:\(f''=J^2+r\,r''=16+2\cdot2=20\),\(f'=Jr=8\),\(p=-8/20=-0.4\),\(x_1=1.6\)。差别来自被 GN 丢掉的第二项 \(r\,r''=4\)。\(r(2)=2\), \(J=r'=2x=4\). GN: \(J^2p=-Jr\Rightarrow p=-r/J=-2/4=-0.5\), so \(x_1=1.5\) (the ancient Babylonian square-root method). Compare full Newton: \(f''=J^2+r\,r''=16+2\cdot2=20\), \(f'=Jr=8\), \(p=-8/20=-0.4\), \(x_1=1.6\). The difference is the second term \(r\,r''=4\) that GN dropped.
    3. LM 中 \(\lambda\) 取得非常大(例如 \(10^6\))时,步 \(\mathbf p\) 近似是什么?3. In LM with a very large \(\lambda\) (say \(10^6\)), what does the step \(\mathbf p\) approximate?
    \(\lambda\) 远大于 \(J^\top J\) 的特征值时,\((J^\top J+\lambda I)^{-1}\approx I/\lambda\),所以 \(\mathbf p\approx-\mathbf g/\lambda\):方向是负梯度,长度被 \(\lambda\) 压得很短(不是严格的零,所以不是选项 C)。\(\lambda\to0\) 才是 Gauss-Newton。这就是 LM"在 GN 与 GD 之间插值"的含义(FIG 4)。When \(\lambda\) greatly exceeds the eigenvalues of \(J^\top J\), \((J^\top J+\lambda I)^{-1}\approx I/\lambda\) and \(\mathbf p\approx-\mathbf g/\lambda\): direction along the negative gradient, length squashed by \(\lambda\) (not exactly zero, so not option C). \(\lambda\to0\) gives Gauss-Newton. That is what "LM interpolates between GN and GD" means (FIG 4).
    4. 一个机械臂的 IK 求解器在手臂几乎伸直时突然让关节转了好几圈。最可能的原因和对策是?4. An IK solver suddenly spins the joints several turns when the arm is nearly straight. What is the most likely cause, and the remedy?
    伸直时 \(\sin q_2\approx0\),\(J\) 的一个奇异值接近 0,\(J^{-1}\) 对应方向被放大约 \(1/\sigma_{\min}\) 倍。阻尼把 \(1/\sigma\) 换成 \(\sigma/(\sigma^2+\lambda)\),在 \(\sigma\) 很小时不再爆炸(FIG 5)。此例中目标在工作空间内,所以 C 错;也不是学习率的问题(GN 没有学习率)。When straight, \(\sin q_2\approx0\) and one singular value of \(J\) is near 0; \(J^{-1}\) amplifies that direction by about \(1/\sigma_{\min}\). Damping replaces \(1/\sigma\) by \(\sigma/(\sigma^2+\lambda)\), which no longer blows up for tiny \(\sigma\) (FIG 5). Here the target is within the workspace, so C is wrong, and it is not a learning-rate issue (GN has no learning rate).

    11延伸阅读Further reading

    • J. Nocedal & S. Wright, Numerical Optimization, 2nd ed., Ch. 3.3 (Newton), Ch. 10 (least squares: Gauss-Newton, Levenberg–Marquardt) — 本节内容的标准参考,收敛性证明完整The standard reference for this lesson, with full convergence proofs
    • K. Madsen, H. B. Nielsen & O. Tingleff, Methods for Non-Linear Least Squares Problems, DTU lecture notes, 2004 — LM 阻尼更新规则(本节用的 \(1-(2\rho-1)^3\))的出处,写得极清楚Source of the LM damping update (the \(1-(2\rho-1)^3\) used here), very clearly written
    • S. Boyd & L. Vandenberghe, Convex Optimization, Ch. 9.5 — Newton 法、仿射不变性与二次收敛的分析Newton's method, affine invariance and the analysis of quadratic convergence
    • T. D. Barfoot, State Estimation for Robotics, Ch. 4 & 6 — Gauss-Newton 在机器人状态估计和位姿估计中的应用Gauss-Newton for robot state and pose estimation
    • S. R. Buss, Introduction to Inverse Kinematics with Jacobian Transpose, Pseudoinverse and Damped Least Squares methods, 2004 — 阻尼最小二乘 IK 的经典短文The classic short note on damped-least-squares IK