Four Checks Before Trusting a Pseudoinverse Calculation
Pseudoinverse and truncated-SVD calculations are especially vulnerable to errors that look plausible in notation. A transpose can be put on the wrong factor, a zero singular value can be “reciprocated,” or a least-squares answer can be presented as if it were the only minimizer. Even a correct truncated approximation does not have a single numerical error until the norm is specified. This tutorial develops checks for these situations using the singular value decomposition (SVD).
Let \(\mathbf{A}\in\mathbb{R}^{m\times n}\) have rank \(r\). An SVD writes \(\mathbf{A}=\mathbf{U}\boldsymbol{\Sigma}\mathbf{V}^{\top}\), where \(\mathbf{U}\in\mathbb{R}^{m\times m}\) and \(\mathbf{V}\in\mathbb{R}^{n\times n}\) are orthogonal, and \(\boldsymbol{\Sigma}\in\mathbb{R}^{m\times n}\) has the singular values \(\sigma_1\ge\cdots\ge\sigma_r>0\) on its diagonal and zeros elsewhere. Orthogonal means, for example, \(\mathbf{U}^{\top}\mathbf{U}=\mathbf{I}_m\). The pseudoinverse has the opposite shape to \(\mathbf{A}\): \(\mathbf{A}^{+}\in\mathbb{R}^{n\times m}\).
This definition makes two basic checks immediate: the pseudoinverse has shape \(n\times m\), and it contains no infinite or undefined entries from attempting to divide by zero. The formula also satisfies the four Moore–Penrose equations, as the next result verifies.
Proof. Let \(\mathbf{D}=\boldsymbol{\Sigma}\boldsymbol{\Sigma}^{+}\in\mathbb{R}^{m\times m}\) and \(\mathbf{E}=\boldsymbol{\Sigma}^{+}\boldsymbol{\Sigma}\in\mathbb{R}^{n\times n}\). Each is diagonal: its first \(r\) diagonal entries are \(1\), and all other entries are \(0\). Using \(\mathbf{U}^{\top}\mathbf{U}=\mathbf{I}_m\) and \(\mathbf{V}^{\top}\mathbf{V}=\mathbf{I}_n\), we get \(\mathbf{A}\mathbf{A}^{+}=\mathbf{U}\mathbf{D}\mathbf{U}^{\top}\) and \(\mathbf{A}^{+}\mathbf{A}=\mathbf{V}\mathbf{E}\mathbf{V}^{\top}\). Both are symmetric because \(\mathbf{D}\) and \(\mathbf{E}\) are symmetric. Also \(\mathbf{D}\boldsymbol{\Sigma}=\boldsymbol{\Sigma}\), since \(\mathbf{D}\) retains exactly the nonzero diagonal entries of \(\boldsymbol{\Sigma}\). Thus \(\mathbf{A}\mathbf{A}^{+}\mathbf{A}=\mathbf{U}\mathbf{D}\boldsymbol{\Sigma}\mathbf{V}^{\top}=\mathbf{U}\boldsymbol{\Sigma}\mathbf{V}^{\top}=\mathbf{A}\). Likewise, \(\mathbf{E}\boldsymbol{\Sigma}^{+}=\boldsymbol{\Sigma}^{+}\), so \(\mathbf{A}^{+}\mathbf{A}\mathbf{A}^{+}=\mathbf{V}\mathbf{E}\boldsymbol{\Sigma}^{+}\mathbf{U}^{\top}=\mathbf{V}\boldsymbol{\Sigma}^{+}\mathbf{U}^{\top}=\mathbf{A}^{+}\). \(\square\)
Dimensions: Check the Product Before the Entries
Worked Example: A Rectangular Pseudoinverse
Consider \(\mathbf{A}=\begin{bmatrix}3&0&0\\0&2&0\end{bmatrix}\), a \(2\times3\) matrix. Its singular values are \(3\) and \(2\), and the corresponding SVD can use identity matrices for \(\mathbf{U}\) and \(\mathbf{V}\). The pseudoinverse is therefore
A common shape error is to write \(\mathbf{A}^{+}\) as another \(2\times3\) matrix. Under that incorrect assumption, \(\mathbf{A}\mathbf{A}^{+}\) would attempt to multiply a \(2\times3\) matrix by a \(2\times3\) matrix. This product is undefined because the inner dimensions, \(3\) and \(2\), do not match. With the correct \(3\times2\) pseudoinverse, the products are defined and have different shapes:
Check: \(\mathbf{A}\mathbf{A}^{+}\) is \(2\times2\), while \(\mathbf{A}^{+}\mathbf{A}\) is \(3\times3\); the computed shapes agree with the multiplication rules. Multiplying once more gives \(\mathbf{A}\mathbf{A}^{+}\mathbf{A}=\mathbf{I}_2\mathbf{A}=\mathbf{A}\).
The two products should not be casually identified with each other: they act on different spaces and need not have the same size. In this example, the first is the identity on \(\mathbb{R}^2\); the second is a \(3\times3\) matrix that retains the first two coordinates and sends the third to zero.
Zero Singular Values Stay Zero
Worked Example: A Zero Singular Value
Let \(\mathbf{B}=\begin{bmatrix}4&0&0\\0&0&0\end{bmatrix}\). Its rank is \(1\), with singular values \(4\) and \(0\). The reciprocal of \(4\) is \(1/4\), but the zero singular value contributes a zero, not a reciprocal:
In particular, writing \(1/0\) in the pseudoinverse is not a valid calculation. Check the products directly:
Then \(\mathbf{B}\mathbf{B}^{+}\mathbf{B}=\begin{bmatrix}1&0\\0&0\end{bmatrix}\begin{bmatrix}4&0&0\\0&0&0\end{bmatrix}=\mathbf{B}\), and \(\mathbf{B}^{+}\mathbf{B}\mathbf{B}^{+}=\begin{bmatrix}1&0&0\\0&0&0\\0&0&0\end{bmatrix}\begin{bmatrix}\frac14&0\\0&0\\0&0\end{bmatrix}=\mathbf{B}^{+}\). Both checks agree with the defining equations.
The zero entries encode directions on which the matrix loses information. The pseudoinverse does not attempt to undo that loss by division by zero. Instead, it inverts the action on the singular directions where the singular value is positive and assigns zero to the remaining directions.
Least-Squares Solutions and the Minimum-Norm Choice
A least-squares solution to \(\mathbf{A}\mathbf{x}=\mathbf{b}\) is a vector \(\mathbf{x}\) minimizing the residual norm \(\|\mathbf{A}\mathbf{x}-\mathbf{b}\|\). A minimum-norm least-squares solution is a least-squares solution whose own norm \(\|\mathbf{x}\|\) is as small as possible. These are two different minimization criteria: first minimize the residual, then among all residual-minimizers choose the one with least norm.
Proof. Write \(\beta_i=\mathbf{u}_i^{\top}\mathbf{b}\) for the columns \(\mathbf{u}_i\) of \(\mathbf{U}\), and write \(\mathbf{x}=\sum_{j=1}^{n}\gamma_j\mathbf{v}_j\). Since \(\mathbf{U}\) and \(\mathbf{V}\) are orthogonal, \(\|\mathbf{A}\mathbf{x}-\mathbf{b}\|^2=\sum_{i=1}^{r}(\sigma_i\gamma_i-\beta_i)^2+\sum_{i=r+1}^{m}\beta_i^2\). The second sum is fixed, and every term in the first sum is minimized, with value zero, exactly when \(\gamma_i=\beta_i/\sigma_i\) for \(1\le i\le r\). The remaining coordinates \(\gamma_{r+1},\ldots,\gamma_n\) do not affect the residual. By the definition of \(\mathbf{A}^{+}\), \(\mathbf{A}^{+}\mathbf{b}=\sum_{i=1}^{r}(\beta_i/\sigma_i)\mathbf{v}_i\), so this gives exactly the stated family of minimizers. Because the \(\mathbf{v}_j\) are orthonormal, the squared norm of any such minimizer is \(\sum_{i=1}^{r}(\beta_i/\sigma_i)^2+\sum_{j=r+1}^{n}c_j^2\). This is smallest exactly when every \(c_j=0\), which proves both minimality and uniqueness. If \(r=n\), the second sum is empty and there are no free coordinates, so the same conclusion holds. \(\square\)
Worked Example: Many Least-Squares Solutions, One Minimum-Norm Solution
Take \(\mathbf{C}=\begin{bmatrix}1&1\\0&0\end{bmatrix}\) and \(\mathbf{b}=\begin{bmatrix}1\\1\end{bmatrix}\). The equation \(\mathbf{C}\mathbf{x}=\mathbf{b}\) is inconsistent because its second coordinate would require \(0=1\). A least-squares solution minimizes \(\|\mathbf{C}\mathbf{x}-\mathbf{b}\|\). For \(\mathbf{x}=(x_1,x_2)^{\top}\), the residual is \((x_1+x_2-1,-1)^{\top}\), so its squared norm is \((x_1+x_2-1)^2+1\). The least-squares solutions are exactly those with \(x_1+x_2=1\).
The rank-one SVD has singular value \(\sqrt{2}\), left singular vector \(\mathbf{u}_1=(1,0)^{\top}\), and right singular vector \(\mathbf{v}_1=(1,1)^{\top}/\sqrt{2}\). Hence \(\mathbf{C}^{+}=\frac12\begin{bmatrix}1&0\\1&0\end{bmatrix}\), and
To see why this is the minimum-norm choice, write any least-squares solution as \((\frac12+t,\frac12-t)^{\top}\). Its squared norm is \((\frac12+t)^2+(\frac12-t)^2=\frac12+2t^2\), minimized uniquely at \(t=0\). The residual at the selected vector is \((0,-1)^{\top}\), whose squared norm is \(1\), the minimum found above.
Correction to a common mistake: saying “the least-squares solution is \(\mathbf{C}^{+}\mathbf{b}\)” can incorrectly suggest that there are no others. Here every vector with coordinate sum \(1\) is a least-squares solution; \(\mathbf{C}^{+}\mathbf{b}\) is the unique one with minimum norm.
Truncated SVD Errors Need a Norm
For a matrix of rank \(q\), write its nonzero singular values in descending order and let \(0\le k<q\). Its rank-\(k\) truncated SVD is \(\mathbf{A}_k=\sum_{i=1}^{k}\sigma_i\mathbf{u}_i\mathbf{v}_i^{\top}\), where the sum is zero when \(k=0\). The difference \(\mathbf{A}-\mathbf{A}_k\) is the discarded tail. Two commonly used measures of that difference are the Frobenius norm \(\|\mathbf{M}\|_F=(\sum_{i,j}m_{ij}^2)^{1/2}\) and the operator norm \(\|\mathbf{M}\|_2=\max_{\|\mathbf{x}\|=1}\|\mathbf{M}\mathbf{x}\|\).
Proof. The discarded tail is \(\mathbf{A}-\mathbf{A}_k=\sum_{i=k+1}^{q}\sigma_i\mathbf{u}_i\mathbf{v}_i^{\top}\). The matrices \(\mathbf{u}_i\mathbf{v}_i^{\top}\) are orthonormal under the Frobenius inner product \(\langle\mathbf{X},\mathbf{Y}\rangle_F=\sum_{a,b}x_{ab}y_{ab}\): the inner product of the terms indexed by \(i\) and \(j\) is \((\mathbf{u}_i^{\top}\mathbf{u}_j)(\mathbf{v}_i^{\top}\mathbf{v}_j)\), which is \(1\) when \(i=j\) and \(0\) otherwise. Therefore the squared Frobenius norm of the tail is \(\sum_{i=k+1}^{q}\sigma_i^2\). For the operator norm, write any unit vector as \(\mathbf{x}=\sum_{j=1}^{n}c_j\mathbf{v}_j\), with \(\sum_j c_j^2=1\). Then \(\|(\mathbf{A}-\mathbf{A}_k)\mathbf{x}\|^2=\sum_{i=k+1}^{q}\sigma_i^2c_i^2\le \sigma_{k+1}^2\sum_{i=k+1}^{q}c_i^2\le\sigma_{k+1}^2\). Equality is attained at \(\mathbf{x}=\mathbf{v}_{k+1}\), proving the operator-norm formula. \(\square\)
Worked Example: The Same Truncation, Two Error Values
Let \(\mathbf{D}=\begin{bmatrix}5&0&0\\0&3&0\\0&0&1\end{bmatrix}\) and retain only its first singular component, so \(k=1\). Then \(\mathbf{D}_1=\begin{bmatrix}5&0&0\\0&0&0\\0&0&0\end{bmatrix}\) and \(\mathbf{D}-\mathbf{D}_1=\begin{bmatrix}0&0&0\\0&3&0\\0&0&1\end{bmatrix}\).
In the Frobenius norm, the error is \(\sqrt{0^2+3^2+1^2}=\sqrt{10}\). In the operator norm, it is \(3\), the largest discarded singular value. The two values measure the same residual matrix in different ways; neither is “the error” without a norm being named. Here the operator-norm value is attained on \(\mathbf{e}_2\), while the Frobenius norm combines the squared contributions from both discarded diagonal entries.
Common Pitfalls and a Final Checklist
Wrong shape: If \(\mathbf{A}\) is \(m\times n\), then \(\mathbf{A}^{+}\) is \(n\times m\). A proposed pseudoinverse with the same shape as a rectangular \(\mathbf{A}\) should trigger a dimension check before any multiplication is attempted. Do not call a product “wrong-sized” if it is actually undefined: for two \(2\times3\) matrices, the product is undefined because the inner dimensions do not match.
Wrong treatment of zero: A diagonal entry \(0\) in \(\boldsymbol{\Sigma}\) contributes \(0\) in \(\boldsymbol{\Sigma}^{+}\), not \(1/0\). Only positive singular values are reciprocated.
Wrong least-squares conclusion: The pseudoinverse always gives the minimum-norm least-squares solution, but a rank-deficient matrix can have other least-squares solutions. Those differ from it by vectors in the null space.
Unspecified truncation error: A statement such as “the error after truncation is \(3\)” is incomplete unless it gives the norm. For the preceding example, \(3\) is the operator-norm error, whereas the Frobenius-norm error is \(\sqrt{10}\).
If \(\mathbf{A}\in\mathbb{R}^{m\times n}\), verify \(\mathbf{A}^{+}\in\mathbb{R}^{n\times m}\) and check inner dimensions in every product.
Reciprocate only the positive singular values; keep zero entries zero.
First identify the least-squares minimizers; then distinguish their unique minimum-norm member.
For a truncated SVD, state whether the residual is measured in the Frobenius norm or the operator norm.
Check Your Understanding
For each question, state the dimensions or norm explicitly and verify any proposed formula.
- If \(\mathbf{A}\in\mathbb{R}^{4\times2}\), what are the dimensions of \(\mathbf{A}^{+}\), \(\mathbf{A}\mathbf{A}^{+}\), and \(\mathbf{A}^{+}\mathbf{A}\)?
- A diagonal singular-value matrix has diagonal entries \(6,0,0\). What entries appear on the corresponding diagonal of its pseudoinverse?
- For a rank-deficient matrix, how can all least-squares solutions be written in terms of \(\mathbf{A}^{+}\mathbf{b}\) and \(\operatorname{Nul}\mathbf{A}\)?
- A matrix has singular values \(4,2,1\), and its rank-one truncated SVD is used. What are the Frobenius-norm and operator-norm errors?
- Explain why the phrase “the truncation error is \(2\)” is incomplete without naming a norm.
