# Blown-up iterates in Gauss-Newton method

I am working on a non-linear least squares problem with standard form, in which I need to calibrate a parameter vector $$\Theta$$ to a set of inputs $$\mathbf{x}$$ and outputs $$\mathbf{y}$$:

\begin{align} & \min_\Theta ||R(\mathbf{x},\mathbf{y};\Theta)||_2^2 \\[4pt] &\quad \text{st: } R(\mathbf{x},\mathbf{y};\Theta):=\mathbf{y}-f(\mathbf{x};\Theta) \end{align}

More specifically, I want to fit the 2-dimensional vector $$\Theta=(\theta_1,\theta_2)$$ to $$\mathbf{y}=\{y_1,y_2\}$$ based on inputs $$\mathbf{x}=\{x_1,x_2\}$$. To do so, I have implemented a Gauss-Newton algorithm. The function $$f$$ is non-linear and requires Monte Carlo methods to evaluate (I am fixing a given simulation seed to avoid noise during the fitting procedure).

I have started by writing a 1-dimensional Newton algorithm for benchmarking purposes. I have then tried to fit $$\theta\in\{\theta_1,\theta_2\}$$ to one point $$y\in\{y_1,y_2\}$$ $$-$$ so yielding four fitting programs in total. All four optimizations above have worked, meaning satisfactory convergence in few iterations.

Next, I have written a Gauss-Newton routine to fit $$\Theta$$ to $$\mathbf{y}$$, and this is where my problems start. My Gauss-Newton update at the first iteration is way too large compared to what I was obtaining when fitting individual parameters. For example:

• When fitting $$\theta_1$$ to $$y_1$$ using Newton, on the first iteration I was obtaining a Newton update approximately equal to $$0.002$$. However, under the Gauss-Newton algorithm, I get an update 3 orders of magnitude larger for $$\theta_1$$, i.e. approximately $$-4.268$$.
• Similarly, when fitting $$\theta_2$$ to $$y_2$$ using Newton, on the first iteration of the Newton algorithm I get an update of approximately $$2.6$$, whereas under Gauss-Newton I get again 3 orders of magnitude more, $$1616.5$$.

I've done the following checks:

• My Gauss-Newton algorithm yields the same results than my Newton one when applied to individual parameter/data-point fitting.
• The four entries on the Jacobian matrix from the first iteration of Gauss-Newton, coincide with the derivatives I was obtaining on the first iteration of each one of the four single parameter-target problems above.
• Using the solutions from the individual Newton fits as initializers for $$\Theta^{(0)}$$.

Below is the Jacobian I am getting on the first iteration to 7 decimals, as well as the initial state of the squared residuals $$R^2_0$$: $$J_{\text{iter}_1}=\begin{pmatrix} 4.998\times 10^{-4} & 1.320\times 10^{-6} \\ 1.520\times 10^{-3} & 4.021\times 10^{-5} \end{pmatrix} \qquad R_0^2 =\begin{pmatrix} 1.189\times 10^{-6}\\ 1.048\times 10^{-6} \end{pmatrix}$$ such that the Gauss-Newton update is: $$\Theta^{(1)} = \Theta^{(0)}-J_{\text{iter}_1}^{-1}R_0^2$$

I am unsure what the problem is here. Starts to look like there is a fundamental issue with applying the Gauss-Newton method, but for the problem I am tackling this algorithm should be applicable.

Any pointers as to what the problem might be or what to look for to continue debugging are appreciated.

• Why is the question being downvoted? – Daneel Olivaw Apr 1 at 16:31
• The Gauss-Newton update would normally be $\Theta^{1}=\Theta^{0}-({J^{0}}^{T}J^{0})^{-1}{J^{0}}^{T}R^{0}$. – Brian Borchers Apr 1 at 16:54