Alas, life is not so simple as our model is very much non-linear and this will mean that we cannot simply solve for the \(p_k\) values that provide the minimum for \(\chi ^2\), but suppose that we have some approximate \(p_k\) values that are close to the minimum for \(\chi ^2\) then we can define \(p_k + \delta _k\) as the values that are at the minimum and we can use a Taylor expansion to write a linear equation for our model: \[M(\mathbf {p}+\delta ,x) = M(\mathbf {p}+\delta ,x) + \sum _j \delta _j \frac {\partial M}{\partial p_j} + \sum _j \frac {\delta _j^2}{2} \frac {\partial ^2 M}{\partial p_j^2} + ... \] we can then we can follow the previous approach. However in the Gauss-Newton method, we will take the derivatives of \(\chi ^2\) with respect to \(p_k\) rather than \(p_k + \delta _k\) and we will approximate the Taylor expansion by only considering the first derivatives terms.
The math becomes slightly more complex than before and now we must know the approximate \(p_k\) values, but we can still solve for the \(\delta _k\) values:
\[\chi ^2 = \sum _j w_j[y_{obs,j}-M(\mathbf {p},x_j) - \sum _i \delta _i \frac {\partial M(\textbf {p},x_j)}{\partial p_i}]^2\] which gives us the coupled equations:
which can be cast, as before, into a form where we can use linear algebra to solve for the \(\delta _k\) values:
We can cast this in two forms, \(\textbf {H} \textbf {p} = \textbf {b}\), where \(H_{ik} = \sum _j w_j \frac {\partial M(\textbf {p},x_j)}{\partial p_i}\frac {\partial M(\mathbf {p},x_j)}{\partial p_k}\) and \(b_k = \sum _j w_j[y_{obs,j}-M(\mathbf {p},x_j)]\frac {\partial M(\mathbf {p},x_j)}{\partial p_k}\)
but there is also a more computationally efficient (though not memory efficient for large problems) where we have an equation \(\textbf {D}^T \textbf {c} = \textbf {D}^T\textbf {D}\delta \) where now \((\textbf {D}^T\textbf {D})^{-1} \textbf {D}^T \textbf {c} = \delta \) and \(D_{jk} = \sqrt {w_j} \frac {\partial M(\textbf {p},x_j)}{\partial p_k}\) and \(c_j = \sqrt {w_j}[y_{obs,j}-M(\mathbf {p},x_j)]\)
There are now several methodology changes from the linear case. We must have approximate values for the parameters. If these values are too far from the correct ones (where correct means the ones that minimize \(\chi ^2\)) then the approximations will fail. We are now solving for shifts, \(\delta _k\), to apply to the parameter values where \(p_k + \delta _k\) will be closer to the correct value than \(p_k\), but still will not be exactly at the minimum. This is why Gauss-Newton minimization is iterative. It does not give us the right value, but if we repeat it, setting our new \(p_k\) value to the previous \(p_k + \delta _k\) result, each time we repeat the computation we get closer to the minimum for \(\chi ^2\). Finally, note that now we use the partial derivatives in this equation.