1.Orthogonal matrices

The fundamental theorem of linear algebra partitions the domain and codomain of a linear mapping $𝒇:U\to V$. For real vectors spaces $U={ℝ}^n$, $V={ℝ}^m$ the partition properties are stated in terms of spaces of the associated matrix $𝑨$ as

${\displaystyle \begin{array}{cccc}C\left(𝑨\right)\oplus N\left({𝑨}^T\right)={ℝ}^m& C\left(𝑨\right)\perp N\left({𝑨}^T\right)& C\left({𝑨}^T\right)\oplus N\left(𝑨\right)={ℝ}^n& C\left({𝑨}^T\right)\perp N\left(𝑨\right)\end{array}.}$

The dimension of the column and row spaces $r=\dim C\left(𝑨\right)=\dim C\left({𝑨}^T\right)$ is the rank of the matrix, $n-r$ is the nullity of $𝑨$, and $m-r$ is the nullity of $A^T$. A infinite number of bases could be defined for the domain and codomain. It is of great theoretical and practical interest to define bases with properties that facilitate insight or computation.

The above partitions of the domain and codomain are orthogonal, and suggest searching for orthogonal bases within these subspaces. Introduce a matrix representation for the bases

${\displaystyle 𝑼=\left[\begin{array}{cccc}{𝒖}_1& {𝒖}_2& \dots & {𝒖}_m\end{array}\right]\in {ℝ}^{m\times m},𝑽=\left[\begin{array}{cccc}{𝒗}_1& {𝒗}_2& \dots & {𝒗}_n\end{array}\right]\in {ℝ}^{n\times n},}$

with $C\left(𝑼\right)={ℝ}^m$ and $C\left(𝑽\right)={ℝ}^n$. Orthogonality between columns ${𝒖}_i$, ${𝒖}_j$ for $i\ne j$ is expressed as ${𝒖}_i^T{𝒖}_j=0$. For $i=j$, the inner product is positive ${𝒖}_i^T{𝒖}_i>0$, and since scaling of the columns of $𝑼$ preserves the spanning property $C\left(𝑼\right)={ℝ}^m$, it is convenient to impose ${𝒖}_i^T{𝒖}_i=1$. Such behavior is concisely expressed as a matrix product

${\displaystyle {𝑼}^T𝑼={𝑰}_m,}$

with ${𝑰}_m$ the identity matrix in ${ℝ}^m$. Expanded in terms of the column vectors of $𝑼$ the first equality is

${\displaystyle {\left[\begin{array}{cccc}{𝒖}_1& {𝒖}_2& \dots & {𝒖}_m\end{array}\right]}^T\left[\begin{array}{cccc}{𝒖}_1& {𝒖}_2& \dots & {𝒖}_m\end{array}\right]=\left[\begin{array}{c}{𝒖}_1^T\\ {𝒖}_2^T\\ ⋮\\ {𝒖}_m^T\end{array}\right]\left[\begin{array}{cccc}{𝒖}_1& {𝒖}_2& \dots & {𝒖}_m\end{array}\right]=\left[\begin{array}{cccc}{𝒖}_1^T{𝒖}_1& {𝒖}_1^T{𝒖}_2& \dots & {𝒖}_1^T{𝒖}_m\\ {𝒖}_2^T{𝒖}_1& {𝒖}_2^T{𝒖}_2& \dots & {𝒖}_2^T{𝒖}_m\\ ⋮& ⋮& \ddots & ⋮\\ {𝒖}_m^T{𝒖}_1& {𝒖}_m^T{𝒖}_2& \dots & {𝒖}_m^T{𝒖}_m\end{array}\right]={𝑰}_m.}$

It is useful to determine if a matrix $𝑿$ exists such that $𝑼𝑿={𝑰}_m$, or

${\displaystyle 𝑼𝑿=𝑼\left[\begin{array}{cccc}{𝒙}_1& {𝒙}_2& \dots & {𝒙}_m\end{array}\right]=\left[\begin{array}{cccc}{𝒆}_1& {𝒆}_2& \dots & {𝒆}_m\end{array}\right].}$

The columns of $𝑿$ are the coordinates of the column vectors of ${𝑰}_m$ in the basis $𝑼$, and can readily be determined

${\displaystyle 𝑼{𝒙}_j={𝒆}_j\Rightarrow {𝑼}^T𝑼{𝒙}_j={𝑼}^T{𝒆}_j\Rightarrow {𝑰}_m{𝒙}_j=\left[\begin{array}{c}{𝒖}_1^T\\ {𝒖}_2^T\\ ⋮\\ {𝒖}_m^T\end{array}\right]{𝒆}_j\Rightarrow {𝒙}_j={\left({𝑼}^T\right)}_j,}$

where ${\left({𝑼}^T\right)}_j$ is the $j^{\mathrm{th}}$ column of ${𝑼}^T$, hence $𝑿={𝑼}^T$, leading to

${\displaystyle {𝑼}^T𝑼=𝑰=𝑼{𝑼}^T.}$

Note that the second equality

${\displaystyle \left[\begin{array}{cccc}{𝒖}_1& {𝒖}_2& \dots & {𝒖}_m\end{array}\right]{\left[\begin{array}{cccc}{𝒖}_1& {𝒖}_2& \dots & {𝒖}_m\end{array}\right]}^T=\left[\begin{array}{cccc}{𝒖}_1& {𝒖}_2& \dots & {𝒖}_m\end{array}\right]\left[\begin{array}{c}{𝒖}_1^T\\ {𝒖}_2^T\\ ⋮\\ {𝒖}_m^T\end{array}\right]={𝒖}_1{𝒖}_1^T+{𝒖}_2{𝒖}_2^T+\cdots +{𝒖}_m{𝒖}_m^T=𝑰}$

acts as normalization condition on the matrices ${𝑼}_j={𝒖}_j{𝒖}_j^T$.

Definition. A square matrix $𝑼$ is said to be orthogonal if ${𝑼}^T𝑼=𝑼{𝑼}^T=𝑰$.

2.Intrinsic basis of a linear mapping

Given a linear mapping $𝒇:U\to V$, expressed as $𝒚=𝒇\left(𝒙\right)=𝑨𝒙$, the simplest description of the action of $𝑨$ would be a simple scaling, as exemplified by $𝒈\left(𝒙\right)=a𝒙$ that has as its associated matrix $a𝑰$. Recall that specification of a vector is typically done in terms of the identity matrix $𝒃=𝑰𝒃$, but may be more insightfully given in some other basis $𝑨𝒙=𝑰𝒃$. This suggests that especially useful bases for the domain and codomain would reduce the action of a linear mapping to scaling along orthogonal directions, and evaluate $𝒚=𝑨𝒙$ by first re-expressing $𝒚$ in another basis $𝑼$, $𝑼𝒔=𝑰𝒚$ and re-expressing $𝒙$ in another basis $𝑽$, $𝑽𝒓=𝑰𝒙$. The condition that the linear operator reduces to simple scaling in these new bases is expressed as $s_i={\sigma }_i{r}_i$ for $i=1,\dots ,min(m,n)$, with ${\sigma }_i$ the scaling coefficients along each direction which can be expressed as a matrix vector product $𝒔=𝚺𝒓$, where $𝚺\in {ℝ}^{m\times n}$ is of the same dimensions as $𝑨$ and given by

${\displaystyle 𝚺=\left[\begin{array}{ccccccc}{\sigma }_1& 0& \dots & 0& 0& \dots & 0\\ 0& {\sigma }_2& \dots & 0& 0& \dots & 0\\ ⋮& ⋮& \ddots & 0& ⋮& \ddots & ⋮\\ 0& 0& \dots & {\sigma }_r& 0& \dots & 0\\ 0& 0& \dots & 0& 0& \dots & 0\\ ⋮& ⋮& \ddots & ⋮& ⋮& \ddots & ⋮\\ 0& 0& \dots & 0& 0& \dots & 0\end{array}\right].}$

Imposing the condition that $𝑼,𝑽$ are orthogonal leads to

${\displaystyle 𝑼𝒔=𝒚\Rightarrow 𝒔={𝑼}^T𝒚,𝑽𝒓=𝒙\Rightarrow 𝒓={𝑽}^T𝒙,}$

which can be replaced into $𝒔=𝚺𝒓$ to obtain

${\displaystyle {𝑼}^T𝒚=𝚺{𝑽}^T𝒙\Rightarrow 𝒚=𝑼𝚺{𝑽}^T𝒙.}$

From the above the orthogonal bases $𝑼,𝑽$ and scaling coefficients $𝚺$ that are sought must satisfy $𝑨=𝑼𝚺{𝑽}^T$. The SVD theorem states that the matrix factors $𝑼,𝚺,𝑽$ do indeed exist.

Theorem. Every matrix $𝑨\in {ℝ}^{m\times n}$ has a singular value decomposition (SVD)

${\displaystyle 𝑨=𝑼𝚺{𝑽}^T,}$

with properties:

The scaling coefficients ${\sigma }_j$ are called the singular values of $𝑨$. The columns of $𝑼$ are called the left singular vectors, and those of $𝑽$ are called the right singular vectors. Carrying out computation of the matrix products

${\displaystyle }$ ${\displaystyle 𝑨=\left[\begin{array}{ccccccc}{𝒖}_1& {𝒖}_2& \dots & {𝒖}_r& {𝒖}_{r+1}& \dots & {𝒖}_m\end{array}\right]\left[\begin{array}{ccccccc}{\sigma }_1& 0& \dots & 0& 0& \dots & 0\\ 0& {\sigma }_2& \dots & 0& 0& \dots & 0\\ ⋮& ⋮& \ddots & 0& ⋮& \ddots & ⋮\\ 0& 0& \dots & {\sigma }_r& 0& \dots & 0\\ 0& 0& \dots & 0& 0& \dots & 0\\ ⋮& ⋮& \ddots & ⋮& ⋮& \ddots & ⋮\\ 0& 0& \dots & 0& 0& \dots & 0\end{array}\right]\left[\begin{array}{c}{𝒗}_1^T\\ {𝒗}_2^T\\ ⋮\\ {𝒗}_r^T\\ {𝒗}_{r+1}^T\\ ⋮\\ {𝒗}_n^T\end{array}\right]\Rightarrow }$ ${\displaystyle 𝑨=\left[\begin{array}{ccccccc}{𝒖}_1& {𝒖}_2& \dots & {𝒖}_r& {𝒖}_{r+1}& \dots & {𝒖}_m\end{array}\right]\left[\begin{array}{c}{\sigma }_1{𝒗}_1^T\\ {\sigma }_2{𝒗}_2^T\\ ⋮\\ {\sigma }_r{𝒗}_r^T\\ ⋮\\ 0\end{array}\right]}$

leads to a representation of $𝑨$ as a sum

${\displaystyle 𝑨=\sum _{i=1}^r{\sigma }_i{𝒖}_i{𝒗}_i^T,}$

with $r⩽min(m,n)$. Written out in full, the above sum is

${\displaystyle 𝑨={\sigma }_1{𝒖}_1{𝒗}_1^T+{\sigma }_2{𝒖}_2{𝒗}_2^T+\cdots +{\sigma }_r{𝒖}_r{𝒗}_r^T.}$

Each product ${𝒖}_i{𝒗}_i^T$ is a matrix of rank one, and is called a rank-one update. Truncation of the above sum to $p$ terms leads to an approximation of $𝑨$

${\displaystyle 𝑨\cong {𝑨}_p=\sum _{i=1}^p{\sigma }_i{𝒖}_i{𝒗}_i^T.}$

In very many cases the singular values exhibit rapid decay, ${\sigma }_1\gg {\sigma }_2\gg \cdots$, such that the approximation above is an accurate representation of the matrix $𝑨$ for $p\ll r$.

3.SVD solution of linear algebra problems

The SVD can be used to solve common problems within linear algebra.

Linear systems.

1. Compute the SVD, $𝑼𝚺{𝑽}^T=𝑨$;

2. Find the coordinates of $𝒃$ in the orthogonal basis $𝑼$, $𝒄={𝑼}^T𝒃$;

3. Scale the coordinates of $𝒄$ by the inverse of the singular values $y_i=c_i/{\sigma }_i$, $i=1,\dots ,m$, such that $\Sigma 𝒚=𝒄$ is satisfied;

4. Find the coordinates of $𝒚$ in basis ${𝑽}^T$, $𝒙=𝑽𝒚$.

Least squares.

The pseudo-inverse.

${\displaystyle (𝑼𝚺{𝑽}^T)𝒙=𝒃\Rightarrow 𝑼(𝚺{𝑽}^T𝒙)=𝒃\Rightarrow (𝚺{𝑽}^T𝒙)={𝑼}^T𝒃\Rightarrow {𝑽}^T𝒙={𝚺}^{+}{𝑼}^T𝒃\Rightarrow 𝒙=𝑽{𝚺}^{+}{𝑼}^T𝒃,}$

where the matrix ${𝚺}^{+}\in {ℝ}^{n\times m}$ (notice the inversion of dimensions) is defined as a matrix with elements ${\sigma }_i^{-1}$ on the diagonal, and is called the pseudo-inverse of $𝚺$. Similarly the matrix

${\displaystyle {𝑨}^{+}=𝑽{𝚺}^{+}{𝑼}^T}$

that allows stating the solution of $𝑨𝒙=𝒃$ simply as $𝒙={𝑨}^{+}𝒃$ is called the pseudo-inverse of $𝑨$. Note that in practice ${𝑨}^{+}$ is not explicitly formed. Rather the notation ${𝑨}^{+}$ is simply a concise reference to carrying out steps 1-4 above.