Gauss-Newton

Gauss-Newton is the fundamental iterative algorithm for nonlinear least squares — the problem class that virtually every SLAM back-end computation (bundle adjustment, pose graph optimization, PnP refinement, direct image alignment) reduces to. It exploits the special sum-of-squares structure of the cost to get near-second-order convergence while only ever computing first derivatives.

Derivation

We minimize a cost built from a residual vector e(x)Rm\mathbf{e}(\mathbf{x}) \in \mathbb{R}^m over a state xRn\mathbf{x} \in \mathbb{R}^n:

F(x)=12e(x)2F(\mathbf{x}) = \frac{1}{2} \|\mathbf{e}(\mathbf{x})\|^2

Linearize the residual (not the cost) around the current estimate xk\mathbf{x}_k using the Jacobian Jk=e/xJ_k = \partial \mathbf{e} / \partial \mathbf{x} evaluated at xk\mathbf{x}_k:

e(xk+Δx)e(xk)+JkΔx\mathbf{e}(\mathbf{x}_k + \Delta\mathbf{x}) \approx \mathbf{e}(\mathbf{x}_k) + J_k \Delta\mathbf{x}

Substituting gives a quadratic model of the cost in Δx\Delta\mathbf{x}:

F(xk+Δx)12ek2+ΔxTJkTek+12ΔxTJkTJkΔxF(\mathbf{x}_k + \Delta\mathbf{x}) \approx \frac{1}{2}\|\mathbf{e}_k\|^2 + \Delta\mathbf{x}^T J_k^T \mathbf{e}_k + \frac{1}{2} \Delta\mathbf{x}^T J_k^T J_k \Delta\mathbf{x}

Setting the derivative with respect to Δx\Delta\mathbf{x} to zero yields the normal equations:

(JkTJk)Δx=JkTek(J_k^T J_k)\, \Delta\mathbf{x} = -J_k^T \mathbf{e}_k

Solve for Δx\Delta\mathbf{x} (in practice by Cholesky factorization of the sparse JTJJ^T J, never by explicit inversion), update xk+1=xk+Δx\mathbf{x}_{k+1} = \mathbf{x}_k + \Delta\mathbf{x}, re-linearize, repeat until the update or the cost change is negligible.

With measurement covariances, residuals are weighted by information matrices Ω=Σ1\Omega = \Sigma^{-1} and the normal equations become JTΩJΔx=JTΩeJ^T \Omega J \,\Delta\mathbf{x} = -J^T \Omega\, \mathbf{e} — the same algebra, and precisely MAP estimation under Gaussian noise.

Relation to Newton’s method

Newton’s method uses the true Hessian of FF:

2F=JTJ+iei2ei\nabla^2 F = J^T J + \sum_i e_i \, \nabla^2 e_i

Gauss-Newton drops the second term, approximating HJTJH \approx J^T J. This is cheap (no second derivatives) and accurate when residuals are small at the optimum (a good model fitting well) or nearly linear — exactly the regime of a converging SLAM problem. Near the solution the convergence is close to quadratic. The approximation also guarantees H0H \succeq 0, so the computed step is a descent direction whenever HH is nonsingular.

Failure modes

On manifolds

Poses live on SE(3)\mathrm{SE}(3), not Rn\mathbb{R}^n, so the update is applied through the exponential map: parameterize the increment as ξR6\boldsymbol{\xi} \in \mathbb{R}^6, solve the normal equations for ξ\boldsymbol{\xi}, and update TTexp(ξ^)T \leftarrow T \cdot \exp(\hat{\boldsymbol{\xi}}). The Jacobians are taken with respect to this local perturbation. All SLAM solvers (g2o, Ceres, GTSAM) implement Gauss-Newton/LM in this “optimize a local vector, retract onto the manifold” form.

Why it matters for SLAM

Gauss-Newton is the inner loop of SLAM. Bundle adjustment is Gauss-Newton/LM on reprojection errors, where the sparse block structure of JTJJ^T J (poses coupled to points only through observations) is exploited via the Schur complement; pose graph optimization is Gauss-Newton on relative-pose residuals; direct methods (LSD-SLAM, DSO) run it on photometric errors; even ICP’s alignment step is a Gauss-Newton iteration. Reading any SLAM back-end paper requires fluency in the vocabulary defined here — residual, Jacobian, Hessian approximation, normal equations, damping — and most practical debugging (why did my optimization diverge? why is my Hessian singular?) traces back to the assumptions listed above.