In nonlinear least squares, why does Gauss-Newton approximate the Hessian by J^T J?
answer
- the objective is a sum of squares
- gradient is J transpose times r
- the dropped term is weighted by residuals
- positive semidefinite for free
- small residuals make it nearly exact
basics
~20 sFor a sum of squared residuals the exact Hessian is J^T J plus a term weighted by the residuals themselves. Gauss-Newton drops that second term, so it needs only first derivatives and always gets a positive semidefinite curvature matrix.
solid answer
~50 sWrite the objective as `f(theta) = 0.5 * sum r_i(theta)^2` with Jacobian `J` of the residual vector. Then the gradient is `J^T r` and the exact Hessian is `J^T J + sum r_i * (second-derivative matrix of r_i)`. Gauss-Newton keeps only `J^T J` and solves `J^T J * delta = -J^T r`. Two things are bought by that: you never need second derivatives of the model, and `J^T J` is positive semidefinite by construction, so the step cannot be an ascent direction when the matrix is invertible. The approximation is good when the residuals at the solution are small or the model is nearly linear in the parameters, because those are exactly the cases where the dropped term is negligible. On large-residual or strongly curved problems the approximation degrades and convergence slows from near-quadratic to linear or worse; adding a damping term `lambda * I` to `J^T J` is the usual fix.
go deeper
Know that least-squares fitting has special structure: the residual Jacobian alone gives both the gradient and a usable curvature matrix, so no second derivatives are required.
Be able to write the objective, its gradient J^T r, and the split of the exact Hessian into J^T J plus a residual-weighted term, and say which term is dropped and why.
Diagnose a fit that stalls or oscillates: check whether residuals stay large at the optimum, whether the Jacobian has lost rank, and whether damping is needed to keep steps trustworthy.
Judge when a curve-fitting problem deserves a specialised least-squares method versus a general optimizer, and what the model-misspecification signal in persistent large residuals is telling you about the problem itself.
## The setting Nonlinear least squares fits a model to data by minimising a sum of squared residuals: ``` f(theta) = 0.5 * sum_i r_i(theta)^2, r_i(theta) = y_i - model(x_i, theta) ``` Here `theta` is the parameter vector of length p and there are n residuals. Stack the residuals into a vector `r` and their first derivatives into the n-by-p Jacobian `J`, where `J_ij = d r_i / d theta_j`. ## Exact derivatives Differentiating once gives the gradient ``` grad f = J^T r ``` Differentiating again produces two pieces: ``` Hessian f = J^T J + sum_i r_i * G_i ``` where `G_i` is the p-by-p matrix of second derivatives of the single residual `r_i`. The first piece depends only on first derivatives of the model. The second piece requires n separate second-derivative matrices and is weighted by the residuals. ## The approximation and what it buys Gauss-Newton simply deletes the second piece and uses `J^T J` as the curvature matrix. The step solves ``` J^T J * delta = -J^T r ``` which is recognisably a linear least-squares problem in `delta`: at each iteration you linearise the model about the current parameters and fit the linearised model to the current residuals. Three practical benefits follow. 1. **No second derivatives.** Only `J` is needed. For a model with many parameters, supplying or computing n second-derivative matrices is often the difference between feasible and not. 2. **Guaranteed curvature sign.** `J^T J` is positive semidefinite for any `J`, because `v^T J^T J v = ||J v||^2 >= 0`. When it is invertible it is positive definite, and then the step has a strictly negative inner product with the gradient — a descent direction, with no eigenvalue repair needed. The exact Hessian carries no such guarantee: a large residual multiplied by a strongly curved residual function can make it indefinite. 3. **Scale of the work.** Forming `J^T J` costs on the order of n p^2 and solving the p-by-p system on the order of p^3, independent of how nasty the model's second derivatives are. ## When the dropped term matters The neglected sum is `sum_i r_i * G_i`. It is small in two regimes: - **Small-residual problems.** If the model fits well, the `r_i` near the solution are close to zero and multiply the curvature matrices down to nothing. This is the classic case — a well-specified curve fit to data with modest noise — and there Gauss-Newton converges nearly as fast as a full Newton method. - **Nearly linear models.** If each residual is close to linear in the parameters, each `G_i` is close to zero regardless of residual size. It is *not* small when the model is badly misspecified (large residuals persist at the optimum) or when the model is strongly curved in the parameters. Then the curvature matrix is systematically wrong, the asymptotic rate drops to linear, and full steps can overshoot enough to increase the objective. ## Standard safeguards Two issues appear in practice. If `J` is rank deficient or nearly so — parameters that are not identifiable from the data, or redundant directions — `J^T J` is singular or badly conditioned and the step explodes. And even with a well-conditioned matrix the full step assumes the linearisation holds over its whole length. The usual response is a damped variant: solve `(J^T J + lambda * I) delta = -J^T r` with a positive `lambda` adjusted between iterations. Large `lambda` makes the step short and biased toward the steepest-descent direction, which is safe far from the solution; small `lambda` recovers the fast Gauss-Newton behaviour near it. This damped form is the Levenberg-Marquardt algorithm and is the default choice for practical curve fitting. A line search that requires the objective to actually decrease serves a similar purpose. ## Relationship to full Newton Gauss-Newton is best understood as Newton's method with a structured, cheap surrogate for the Hessian, exploiting the fact that a sum-of-squares objective hands you most of its curvature for free through the Jacobian. It is not a general-purpose optimizer: applied to an objective that is not a sum of squares there is no `J^T J` to form. Where it applies, it is usually the first thing to try before reaching for a full second-order method.
- When does the Gauss-Newton approximation break down badly?When residuals at the solution stay large and the residual functions are strongly curved in the parameters, because then the dropped term `sum r_i * G_i` is a real part of the curvature. Symptoms are steps that increase the objective and convergence that slows to linear. A misspecified model is the usual cause, so the fix is often to reconsider the model rather than to tune the optimizer.
- What happens if the Jacobian is rank deficient at the current parameters?`J^T J` becomes singular or nearly so, the linear system has no unique solution, and the computed step blows up along the poorly determined directions. It signals parameters that the data cannot separately identify. The practical fix is to add a positive multiple of the identity to `J^T J`, which bounds the step length and makes the system solvable.
- How does adding lambda times the identity to J^T J change the step?It shrinks the step and rotates it toward the negative gradient direction; as lambda grows the move becomes short and cautious, and as lambda goes to zero it returns to the pure Gauss-Newton step. Adaptive schemes raise lambda when a trial step fails to reduce the objective and lower it when steps succeed, giving robustness far away and speed nearby.
saying these in an interview costs you the question
- Says Gauss-Newton uses the exact Hessian
- Claims it works for any objective, not just sums of squares
- Thinks J^T J can be indefinite
- Ignores that large residuals invalidate the approximation
- Confuses the Jacobian of residuals with the Hessian