Minimizing Cost Functions using Iterative Least Squares

Question

I am currently trying to use iterative least squares to solve a system, $y = Hx + v$ where $y$ is a vector of observations, $H$ is the design matrix, and $v$ is the observation error.

From my understanding, iterative least squares tries to minimize a cost $J$ where $J = \min (Hx-y)^2$.

My pseudocode for this is

while ||x^k-x^{k-1}||
   H = h(x^{k-1})
   r = f(x^{k-1}) - y
   dx = pinv(H)*r
   x^{k} = x^{k-1} + dx
end

where $H$ is the Jacobian and the output of $f(x)$ is the expected observations.

Can $r$ be called the cost of the system because it is the residual/difference between $f(x)$ at estimate $x^{k-1}$ and $y$? Should it be squared in an implementation of iterative least squares?

Answer 1

The linear least squares problem $\min \|Hx-y\|^2$ is equivalent to solving $H^THx = H^Ty$, known as the normal equation. $H^T$ denotes the transpose, or adjoint matrix. There are many ways to solve this, depending on your problem, but the conjugate gradient method is a good general method.

Note that the pseudo-inverse (pinv) function you're using will calculate $(H^TH)^{-1}H^T$, so

x=pinv(H)*y

This means your iterative method would have converged in the first iteration. If performance is not an issue and $H$ is reasonably small, that is the easiest way to solve your problem, but for larger $H$ iterative methods will often be significantly faster. This is particularly the case if $H$ is a sparse matrix, or the matrix-vector product $Hr$ can be calculated via a fast method such as the Fourier transform.

As hardmath mentioned, $r$ is the residual vector and $J=\|r\|^2$.