Suppose I want to solve an overdetermined linear least squares problem

$$ x = \operatorname*{argmin}_{x\in\mathbb{R}^n} \| Ax - b\|^2 $$

where $A \in \mathbb{R}^{m\times n}$ has full column rank. I make no assumptions that this system is nearly consistent.

Suppose I know a have access to a matrix $P\in\mathbb{R}^{n\times m}$ which is a good approximation to the Moore–Penrose pseudoinverse $P \approx A^+$. I only have access to $P$ through black-box subroutines to compute $y \mapsto Py$ and $z\mapsto P^*z$. Can I use $P$ to compute $x$ rapidly by some kind of iterative scheme?

My instinct is that this is not possible. All of the work I'm preconditioning for least squares I know of either use a right preconditioner (i.e., assume I can apply maps $R$ and $R^{-1}$ such that $AR^{-1}$ is well-conditioned) or use an (approximate) QR factorization of $A$. It makes sense why left preconditioning may not be possible for this kind of problem, but I haven't seen any authoritative source state that this type of left-preconditioning is not possible/not known to be possible.

  • $\begingroup$ From your minimization problem you can derive: $A^TAx = A^Tb$. Multiply by $P^T$ on the left and substitute $x = Pz$ then: $P^TA^TAPz = P^TA^Tb$. Apply the conjugate gradients solver for the normal equations (CGNR). There should also be an explicit version of CGNR that accepts a preconditioner and saves some operations. $\endgroup$
    – lightxbulb
    Oct 11, 2022 at 9:59
  • $\begingroup$ Thanks for the answer @lightxbulb! Do you have a reference for this approach? My instinct was that this shouldn't work since $AP$ is singular, but testing it out in Matlab (with LSQR instead of CGNR) indeed seemed to work $\endgroup$
    – eepperly16
    Oct 11, 2022 at 15:07
  • $\begingroup$ LSQR is fine, it's a numerically more stable alternative to CGNR, but LSQR requires a bit more computation per iteration. You can see Saad's book for CGNR. Regarding singularity, $A^TA$ is full rank in your case, but even if it were singular if you pick an initial guess of zero with CGNR you will end up in the $A^{+}b$ solution, if you aren't using a preconditioner. If you use a left and right preconditioner then you essentially solve $\min_z \|APz - b\|^2$. If $x^* = A^{+}b$ is in the range of $P$ there should be no issues, otherwise you get a solution on a subspace of the original problem. $\endgroup$
    – lightxbulb
    Oct 11, 2022 at 17:12
  • $\begingroup$ @lightxbulb My point was that $AP$ is singular, not that $A^\top A$ is, so we're in agreement there. Is there a reference for this preconditioning strategy of solving $\min_z \|APz - b\|$? I haven't seen it in the literature $\endgroup$
    – eepperly16
    Oct 11, 2022 at 17:24
  • $\begingroup$ I don't know about literature, but solving $\min_z \|Ax - b\|$ with preconditioners $P, P^T$ is equivalent to solving $\min_z\|APz - b\|$ since the two result in the same normal equations $ P^TA^TAPz = P^TA^Tb$. So the only question is whether the solution $x^*$ to $A^TAx = A^Tb$ matches $Pz^*$ where $z^*$ is the solution of $P^TA^TAPz = P^TA^Tb$. The latter depends on your preconditioner, one requirement is that $x^*$ be in the range of $P$. I am not certain whether an additional condition is required to make the solution unique, but it shouldn't be too hard to derive. $\endgroup$
    – lightxbulb
    Oct 11, 2022 at 17:34


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.