Editing Gauss–Newton algorithm (section)

== Description ==
Given <math>m</math> functions <math>\textbf{r} = (r_1, \ldots, r_m)</math> (often called residuals) of <math>n</math> variables <math>\boldsymbol{\beta} = (\beta_1, \ldots \beta_n),</math> with <math>m \geq n,</math> the Gauss–Newton algorithm [[iterative method|iteratively]] finds the value of <math>\beta</math> that minimize the sum of squares<ref name=ab>Björck (1996)</ref>
<math display="block"> S(\boldsymbol \beta) = \sum_{i=1}^m r_i(\boldsymbol \beta)^{2}.</math>

Starting with an initial guess <math>\boldsymbol \beta^{(0)}</math> for the minimum, the method proceeds by the iterations
<math display="block"> \boldsymbol \beta^{(s+1)} = \boldsymbol \beta^{(s)} - \left(\mathbf{J_r}^\operatorname{T} \mathbf{J_r} \right)^{-1} \mathbf{J_r}^\operatorname{T} \mathbf{r}\left(\boldsymbol \beta^{(s)}\right), </math>

where, if '''r''' and '''''β''''' are [[column vectors]], the entries of the [[Jacobian matrix]] are
<math display="block"> \left(\mathbf{J_r}\right)_{ij} = \frac{\partial r_i \left(\boldsymbol \beta^{(s)}\right)}{\partial \beta_j},</math>

and the symbol <math>^\operatorname{T}</math> denotes the [[matrix transpose]].

At each iteration, the update <math>\Delta = \boldsymbol \beta^{(s+1)} - \boldsymbol \beta^{(s)}</math> can be found by rearranging the previous equation in the following two steps:

*<math>\Delta = -\left(\mathbf{J_r}^\operatorname{T} \mathbf{J_r} \right)^{-1} \mathbf{J_r}^\operatorname{T} \mathbf{r}\left(\boldsymbol \beta^{(s)}\right)</math>
*<math>\mathbf{J_r}^\operatorname{T} \mathbf{J_r} \Delta = -\mathbf{J_r}^\operatorname{T} \mathbf{r}\left(\boldsymbol \beta^{(s)}\right) </math> 

With substitutions <math display="inline">A = \mathbf{J_r}^\operatorname{T} \mathbf{J_r} </math>, <math>\mathbf{b} = -\mathbf{J_r}^\operatorname{T} \mathbf{r}\left(\boldsymbol \beta^{(s)}\right) </math>, and <math>\mathbf {x} = \Delta </math>, this turns into the conventional matrix equation of form <math>A\mathbf {x} = \mathbf {b} </math>, which can then be solved in a variety of methods (see [[#Notes|Notes]]).

If {{math|1=''m'' = ''n''}}, the iteration simplifies to

<math display="block"> \boldsymbol \beta^{(s+1)} = \boldsymbol \beta^{(s)} - \left(\mathbf{J_r}\right)^{-1} \mathbf{r}\left(\boldsymbol \beta^{(s)}\right),</math>

which is a direct generalization of [[Newton's method]] in one dimension.

In data fitting, where the goal is to find the parameters <math>\boldsymbol{\beta}</math> such that a given model function <math> \mathbf{f}(\mathbf{x}, \boldsymbol{\beta}) </math> best fits some data points <math> (x_i, y_i) </math>, the functions <math> r_i </math>are the [[residual (statistics)|residuals]]:
<math display="block">r_i(\boldsymbol \beta) = y_i - f\left(x_i, \boldsymbol \beta\right).</math>

Then, the Gauss–Newton method can be expressed in terms of the Jacobian <math> \mathbf{J_f} = -\mathbf{J_r} </math> of the function <math> \mathbf{f} </math> as
<math display="block"> \boldsymbol \beta^{(s+1)} = \boldsymbol \beta^{(s)} + \left(\mathbf{J_f}^\operatorname{T} \mathbf{J_f} \right)^{-1} \mathbf{J_f}^\operatorname{T} \mathbf{r}\left(\boldsymbol \beta^{(s)}\right). </math>

Note that <math>\left(\mathbf{J_f}^\operatorname{T} \mathbf{J_f}\right)^{-1} \mathbf{J_f}^\operatorname{T}</math> is the left [[Moore–Penrose pseudoinverse|pseudoinverse]] of <math>\mathbf{J_f}</math>.