Least Squares: Gauss-Newton

The normal equations work beautifully when the model is linear in its parameters. But many scientific models are not like that.

Think of exponential decay, enzyme kinetics, adsorption isotherms, or chemical equilibrium models. Their parameters may appear inside exponentials, denominators, or other nonlinear expressions. In those cases, the model cannot be written as:

$$\mathbf{f}(\mathbf{X},\mathbf{\beta })=\mathbf{X}\mathbf{\beta }$$

So we need another idea.

The Gauss-Newton method solves nonlinear least squares by repeatedly pretending, locally, that the nonlinear model is linear [1],[2].

The nonlinear least squares problem

The objective is still:

$$\underset{\mathbf{\beta }}{\min }{\Vert \mathbf{y}-\mathbf{f}(\mathbf{X},\mathbf{\beta })\Vert }^2$$

The difference is that $\mathbf{f}$ is now nonlinear in $\mathbf{\beta }$.

Define the residual vector:

$$\mathbf{r}(\mathbf{\beta })=\mathbf{y}-\mathbf{f}(\mathbf{X},\mathbf{\beta })$$

We want to make that residual vector as small as possible.

The local linear approximation

At the current parameter estimate ${\mathbf{\beta }}^{(k)}$, Gauss-Newton approximates the model with a first-order Taylor expansion [1],[2]:

$$\mathbf{f}(\mathbf{X},\mathbf{\beta })\approx \mathbf{f}(\mathbf{X},{\mathbf{\beta }}^{(k)})+\mathbf{J}({\mathbf{\beta }}^{(k)})(\mathbf{\beta }-{\mathbf{\beta }}^{(k)})$$

The matrix $\mathbf{J}$ is the Jacobian. It contains the partial derivatives of each model prediction with respect to each parameter.

If the model has many outputs and several parameters, the Jacobian tells us how every prediction changes when every parameter changes slightly.

Interactive example

The example below fits a nonlinear exponential curve. Each step builds a local linear approximation around the current parameters, solves a small least squares correction, and updates the curve. Watch how the fit improves quickly when the local approximation is good.

The Gauss-Newton step

Instead of solving the nonlinear problem directly, we solve a linear least squares problem for a parameter correction $\mathbf{\delta }$ [1],[2]:

$${\mathbf{J}}^T\mathbf{J}\mathbf{\delta }={\mathbf{J}}^T\mathbf{r}$$

Then we update:

$${\mathbf{\beta }}^{(k+1)}={\mathbf{\beta }}^{(k)}+\mathbf{\delta }$$

So each iteration says:

1. Look at the current model.
2. Build a local linear approximation.
3. Solve a least squares problem for the correction.
4. Update the parameters.

Then repeat.

When it works well

Gauss-Newton works well when the model is not too nonlinear near the solution and the starting guess is already reasonably close. In that situation, the local linear approximation is accurate, and the method can converge quickly [1],[2].

It can struggle when the starting point is poor, the residuals are large, or the model is strongly nonlinear. Then the local linear approximation may point in a direction that is mathematically reasonable but practically too bold.

That is where Levenberg-Marquardt becomes useful.

What to remember

Gauss-Newton is not a completely different philosophy from linear least squares. It is more like a clever loop around it. At each iteration, it turns a nonlinear problem into a temporary linear one.

That is why understanding the normal equations first pays off: they quietly reappear inside nonlinear optimization.

References