3
$\begingroup$

Suppose I have a system:

$$ Ax = b $$

where $A$ is a $m$ by $n$ matrix which is less than full rank (neither full column nor row rank). In my particular case $m<n$.

I'd like a combination of a minimum norm solution (for elements of $x$ which are not determined and a least squares solution for those that are determined or over determined).

I know I can approximate this with Tikhonov regularization:

$$ \mathop{\text{minimize}}\limits_x \|Ax-b\|^2 + \alpha\|x\|^2. $$

But this is not ideal as the regularizer tugs on all elements of $x$ not just those left undetermined. If I make $\alpha$ too large then my $x$ at determined or over determined elements will not be optimal (w.r.t. $\|Ax-b\|^2$). If $\alpha$ is too small, I imagine I run into conditioning issues.

Is there no way to formulate this correctly using some sort null space/QR decomposition?

Update: See accepted answer, but for completeness here's a full recipe: \begin{align} Ax &= b\\ AE &= QR &\text{rank-revealing QR decomposition of A}\\ AE &= Q\left(\begin{array}{c}I\\0\end{array}\right)\left(\begin{array}{cc}R_1&R_2\end{array}\right) & \text{only keep $r$ (rank) many rows of $R$}\\ \left(\begin{array}{c}R_1^T\\R_2^T\end{array}\right)F &= VT_\text{full} &\text{rank-revealing QR decomposition of $\left(\begin{array}{c}R_1^T\\R_2^T\end{array}\right)$} \\ \left(\begin{array}{c}R_1^T\\R_2^T\end{array}\right)F &= V\left(\begin{array}{c}I\\0\end{array}\right)T & \text{only keep $r$ many rows of $T$}\\ AE &= Q\left(\begin{array}{c}I\\0\end{array}\right)FT^T\left(\begin{array}{cc}I&0\end{array}\right)V^T \\ x &= A^+b\\ x &= EV\left(\begin{array}{c}I\\0\end{array}\right)T^{-T}F^T\left(\begin{array}{cc}I&0\end{array}\right)Q^Tb \end{align}

Improve this question
$\endgroup$
5
  • $\begingroup$ Preemptive strike: Here's an example of an under and over determined system: $$\left(\begin{array}{ccc}1&0&0\\1&0&0\\0&1&0\end{array}\right)\left(\begin{array}x_1\\x_2\\x_3\end{array}\right) = \left(\begin{array}{c}4\\5\\6\end{array}\right)$$ $\endgroup$
    – Alec Jacobson
    Commented May 6, 2014 at 22:12
  • $\begingroup$ Does SVD not give you what you want? $\endgroup$
    – Igor Khavkine
    Commented May 7, 2014 at 8:54
  • $\begingroup$ My matrices are very large and sparse, so complete SVD is not really an option (I assume). $\endgroup$
    – Alec Jacobson
    Commented May 7, 2014 at 14:32
  • $\begingroup$ The SVD will give the answer, but it's not always the minimal amount of work to get a desired datum. They are implementations for sparse cases. stats.stackexchange.com/questions/41259/… I expect that using one of them is probably better than trying to implement one's own minimalist approach (with QR rather than full SVD). $\endgroup$
    – Conder
    Commented May 8, 2014 at 4:06
  • $\begingroup$ I'm happy to implement my own approach with QR if it avoids a dense matrix of singular vectors... $\endgroup$
    – Alec Jacobson
    Commented May 8, 2014 at 12:00

1 Answer 1

Reset to default
4
$\begingroup$

The Moore-Penrose pseudoinverse is probably what you're looking for. The pseudoinverse solution $A^+b$ is the smallest norm $x$ such that $\|Ax-b\|_2$ is minimized. It can be computed using QR decomposition although you have to use rank-revealing QR when $A$ does not have full column rank. A more complete explanation is given here.

Improve this answer
$\endgroup$
4
  • $\begingroup$ So if I have $x = A^+b$ and $A^+=(A^TA)^{-1}A^T$ and the reduced QR decomposition of $A$ is $A=Q_1R_1$, then I'll want to solve using $R_1^TR_1x = A^Tb$ apply back-,forward- substitution with $R_1$? But if $A$ doesn't have full column rank, won't $R_1$ be rank deficient? Or what do you mean by "reduced QR"? $\endgroup$
    – Alec Jacobson
    Commented May 6, 2014 at 22:55
  • $\begingroup$ @mangledorf Sorry, I think the right term is rank-revealing QR decomposition. This presentation explains the method and its application to this problem. $\endgroup$
    – Alex Becker
    Commented May 6, 2014 at 23:23
  • $\begingroup$ I think you do want the pseudoinverse. An alternative way to compute, take random (orthogonal) right multiplicators $T$ until from QR you decompose $AT=QR$ where $R$ is of the special form $\Bigl({U|B\atop Z}\Bigr)$ where $Z$ is numericaly zero (possibly empty) and $U$ is invertible and upper triangular. When $B$ is empty (full column rank), you are done (take $U^{-1}$, suitably zero-extended). Else you can recurse (only once) to compute the pseudo inverse of $(U|B)$ by taking the transpose before and after, and if the same parameters are used for zero-detection, this has full column rank. $\endgroup$
    – Conder
    Commented May 7, 2014 at 4:25
  • $\begingroup$ @Conder can you elaborate? Is this the same rank-revealing QR that AlexBecker writes of? How can $B$ be empty and full column rank? What does suitably zero-extended mean? Does $B$ participate in the final solve? $\endgroup$
    – Alec Jacobson
    Commented May 8, 2014 at 2:03

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .