mathematics//optimization//nonlinear least squares

Nonlinear least squares is the problem of choosing parameters that minimize a sum of squared residuals when the residuals depend nonlinearly on those parameters, and it is the engine of graph SLAM, camera calibration, bundle adjustment in photogrammetry and the fitting of physical models to test data. Each residual \(r_i(\theta)\) is the gap between what the model predicts for a measurement and what was measured.


Nonlinear least squares is the problem of choosing parameters that minimize a sum of squared residuals when the residuals depend nonlinearly on those parameters, and it is the engine of graph SLAM, camera calibration, bundle adjustment in photogrammetry and the fitting of physical models to test data. Each residual ri(θ)r_i(\theta)ri​(θ) is the gap between what the model predicts for a measurement and what was measured.

min⁡θ  ∑i∥ri(θ)∥2\min_{\theta}\; \sum_{i} \big\lVert r_i(\theta) \big\rVert^2θmin​i∑​​ri​(θ)​2

Fitting a motor's thrust curve T=aω2+bωT=a\omega^2+b\omegaT=aω2+bω is linear in aaa and bbb, so ordinary least squares solves it in one step. Fitting the cooling of a motor housing, T(t)=T∞+(T0−T∞)e−t/τT(t)=T_\infty+(T_0-T_\infty)e^{-t/\tau}T(t)=T∞​+(T0​−T∞​)e−t/τ, is nonlinear in the time constant τ\tauτ, and no single step will do. The Gauss-Newton method linearizes every residual around the current guess, r(θ+δ)≈r(θ)+Jδr(\theta+\delta)\approx r(\theta)+J\deltar(θ+δ)≈r(θ)+Jδ with JJJ the Jacobian, solves that linear least-squares problem for the correction δ\deltaδ, and repeats. It is Newton's method with the Hessian approximated by J⊤JJ^{\top}JJ⊤J, so no second derivatives are needed. The Levenberg-Marquardt algorithm adds a damping term that shifts the step toward plain gradient descent whenever the linearization proves poor, which makes it the robust default.

With independent Gaussian errors, minimizing the squared residuals is maximum likelihood estimation, and weighting each residual by the inverse of its covariance keeps it so when sensors differ. That is why the same solver serves calibration, mapping and model fitting: each is a likelihood written as a sum of squares.

At scale the structure is sparse. A SLAM graph has thousands of poses but each measurement touches only two of them, so the Jacobian is mostly zeros, and solvers built for it (g2o, GTSAM, Ceres) handle problems with hundreds of thousands of variables (SLAM).

The problem is non-convex, so the starting guess matters: a fit started far from the truth can converge to a wrong local minimum, and a single gross outlier, squared, can drag the whole solution. Robust losses that grow slower than the square, or a first pass that discards outliers (RANSAC), protect it.

The data must excite the parameters. A thrust test run only between 6000 and 7000 rpm makes ω2\omega^2ω2 and ω\omegaω nearly proportional, J⊤JJ^{\top}JJ⊤J becomes nearly singular (condition number), and the fit returns huge coefficients of opposite sign; the fix is in the test (sweep the whole range), never in the solver.