7.4 Singular Value Decomposition

Print preview

Section 7.4 Singular Value Decomposition

Handout 7.4 Singular Value Decomposition

Objectives: Topics

The Singular Value Decomposition (SVD)

🔗
Applications of the SVD

🔗

🔗

Objectives: Goals

Compute the SVD of a rectangular matrix.

🔗
Use the SVD to:
- construct bases for the four fundamental subspaces,
  
  🔗
- construct spectral decompositions.
  
  🔗
🔗
🔗

🔗

Theorem 7.37. Singular Value Decomposition.

Let \(A\) be an \(m\times n\) matrix with rank \(r\text{.}\) Then \(A\) can be written as

🔗

\begin{equation*} A = U\Sigma V^T = \begin{bmatrix} \vec u_1 \amp \cdots \amp \vec u_m \end{bmatrix} \begin{bmatrix} \sigma_1 \amp 0 \amp 0 \amp \cdots \amp 0 \\ 0 \amp \sigma_2 \amp 0 \amp \cdots \amp 0 \\ \vdots \amp \vdots \amp \ddots \amp \vdots \amp \vdots \\ 0 \amp 0 \amp \cdots \amp \sigma_r \amp 0 \\ 0 \amp 0 \amp \cdots \amp 0 \amp 0 \end{bmatrix}\begin{bmatrix} \vec v_1^T \\ \vdots \\ \vec v_n^T \end{bmatrix} \end{equation*}

where

🔗

\(U\) is an \(m\times m\) orthogonal matrix,

🔗
\(V\) is an \(n\times n\) orthogonal matrix,

🔗
\(\Sigma\) is an \(m\times n\) diagonal matrix with diagonal entries \(\sigma_1,\dots,\sigma_r\text{,}\) where \(\sigma_1 \geq \sigma_2 \geq \cdots \geq \sigma_r \gt 0\text{.}\)

🔗

🔗

Note that Singular Value Decomposition applies to any matrix, not just symmetric or square matrices.

🔗

Example 7.38.

Consider the the linear transformation \(T(\vec x) = A\vec x\text{,}\) where

\begin{equation*} A = \begin{bmatrix} 3 \amp -1 \\ -1 \amp 3 \end{bmatrix}\text{.} \end{equation*}

Note that \(A\) has an orthogonal diagonalization given by

\begin{equation*} A = PDP^T = \frac{1}{\sqrt{2}} \begin{bmatrix}1 \amp 1\\ -1 \amp 1\end{bmatrix} \begin{bmatrix}4 \amp 0\\ 0 \amp 2\end{bmatrix}\dfrac{1}{\sqrt{2}}\begin{bmatrix} 1 \amp -1 \\ 1 \amp 1 \end{bmatrix}\text{.} \end{equation*}

Also, note that the orthogonal diagonalization of \(A\) satisfies the properties of a Singular Value Decomposition. The parts below investigate how multiplication by \(A\) acts on \(\mathbb{R}^2\) geometrically.

🔗

(a)

Can you describe geometrically what multiplication by \(A\) does to an arbitrary vector \(\vec x\text{?}\)

🔗

(b)

Instead of determining what \(A\) does in one fell swoop, let us investigate what each factor of \(A =PDP^T\) does. First, determine what multiplication by \(P^T\) does to an arbitrary vector \(\vec x\text{?}\)

🔗

(c)

Next, determine what multiplication by \(D\) does to an arbitrary vector \(\vec x\text{?}\)

🔗

(d)

Now, determine what multiplication by \(P\) does to an arbitrary vector \(\vec x\text{?}\)

🔗

(e)

Combining parts (b) - (d) above, give a geometric description of what multiplication by \(A\) does to an arbitrary vector \(\vec x\text{.}\)

🔗

Recall from Theorem 6.34 that multiplication by an orthogonal matrix does not change the angles between vectors nor the vector’s length. Therefore, if \(U\) is an orthogonal matrix, we have that the linear transformation \(T(\vec x) = U\vec x\) must simply be a rotation and/or reflection.

🔗

Note that multiplication by a \(m \times n\) diagonal matrix \(\Sigma\) changes the vector from being in \(\mathbb{R}^n\) to being in \(\mathbb{R}^m\text{,}\) and then scales it along the appropriate coordinate axes.

🔗

Therefore, in general, we can use the SVD of a matrix \(A\) to describe how the linear transformation \(T(\vec x) = A\vec x\) acts on the vectors in \(\mathbb{R}^n\text{.}\)

🔗

Visualization of the SVD expressing a matrix as rotations and scalings.

🔗

Note 7.39. Algorithm for Finding the SVD of \(A\).

How do we find the SVD of a matrix \(A\text{?}\)

Start by forming the symmetric matrix \(A^TA\) and finding its eigenvalues \(\lambda_i\text{,}\) with the eigenvalues assigned in nonincreasing order. Set \(\sigma_i = \sqrt{\lambda_i}\text{.}\) These are the singular values of \(A\) and they are the entries on the diagonal of \(\Sigma\text{.}\) Fill any missing diagonal elements with 0.
- Since \(A^TA\) is symmetric, we know from Theorem 7.18 that it has \(n\) real eigenvalues. A bonus fact is that for symmetric matrices of the form \(A^TA\text{,}\) the \(n\) eigenvalues must be nonnegative, so
  
  \begin{equation*} \lambda_1 \ge \lambda_2 \ge \cdots \ge \lambda_n \ge 0. \end{equation*}
  
  🔗
  
  🔗
🔗
🔗
Next, find \(n\) normalized and orthogonal eigenvectors of \(A^TA\) and call them \(\vec v_i\text{,}\) with \(\vec v_i\) corresponding to eigenvalue \(\lambda_i\text{.}\) These vectors form the rows of \(V^T\text{.}\)
- Make sure to find enough orthonormal eigenvectors for each repeated eigenvalue, which we know we can do by Theorem 7.18.
  🔗
  
  🔗
- Alternatively, we can think of the \(\vec v_i\) as the columns of \(V\text{,}\) but since SVD uses \(V^T\text{,}\) it is usually easier to think with respect to \(V^T\text{.}\)
  🔗
  
  🔗
🔗
🔗
Use the formula

\begin{equation*} \vec u_i = \dfrac{1}{\sigma_i}A\vec v_i \end{equation*}

to find the columns of \(U\text{.}\) If the rank of \(A\) is less than \(m\text{,}\) then this formula will not produce enough columns to fill \(U\text{.}\) We must find the remaining columns by extending the set \(\left\{\vec u_1, \dots, \vec u_r\right\}\) to an orthonormal basis of \(\mathbb{R}^m\text{.}\)
- The calculated vectors \(\vec u_i\) will be normalized and orthogonal to one another as long as the original \(\vec v_i\) were.
  🔗
  
  🔗
- If the rank of \(A\) is less than \(n\text{,}\) then some of the eigenvalues found in (1) will be 0, which results in singular values of 0. These singular values can not be used in the formula, since that would result in dividing by 0.
  🔗
  
  🔗
🔗
🔗

🔗

Example 7.40.

Construct the SVD of

🔗

\begin{equation*} A = \begin{bmatrix}3 \amp 8 \\ 0 \amp 3 \end{bmatrix}. \end{equation*}

🔗

Example 7.41.

Construct the SVD of

🔗

\begin{equation*} A=\begin{bmatrix}1 \amp -1\\ -2 \amp 2\\ 2 \amp -2\end{bmatrix}. \end{equation*}

🔗

The SVD has been applied to many modern applications in CS, engineering, and mathematics. Some applications are:

Constructing bases for the four fundamental subspaces

🔗
Forming approximations of a matrix

🔗
Estimating the rank and the condition number of a matrix

🔗
Computing the pseudoinverse of a matrix

🔗
Linear least-squares problems

🔗
Non-linear least-squares problems

🔗
Image compression

🔗
Machine learning and data mining

🔗
Facial Recognition

🔗
Principal component analysis

🔗

🔗

Students are expected to be familiar with the first two items. We will explore these two items below.

🔗

Let \(A = U\Sigma V^T\) be an \(m \times n\) matrix with rank \(r\text{.}\)

🔗

Consider the columns of \(U\text{:}\)

Recall that the formula for calculating the first few columns of \(U\) is

\begin{equation*} \vec u_i = \dfrac{1}{\sigma_i}A\vec v_i\text{.} \end{equation*}

🔗

🔗
It follows that \(\left\{\vec u_1, \dots, \vec u_r\right\}\) is a set of \(\underline{\hspace{2in}}\) vectors that all live in \(\underline{\hspace{1.5in}}\text{.}\)
🔗

🔗
Therefore, the first \(r\) columns of \(U\) form an \(\underline{\hspace{3in}}\) for \(\underline{\hspace{1.5in}}\)
🔗

🔗
Since the remaining columns of \(U\) are created so that all columns of \(U\) form an orthonormal basis for \(\mathbb{R}^m\text{,}\) we see that \(\left\{\vec u_{r+1}, \dots, \vec u_m\right\}\) is a set of \(\underline{\hspace{2in}}\) vectors which live in \(\underline{\hspace{3in}}\text{.}\)
🔗

🔗
Therefore, the last \(m - r\) columns of \(U\) form an \(\underline{\hspace{3in}}\) for \(\underline{\hspace{1.5in}}\text{.}\)
🔗

🔗

Now, consider the rows of \(V^T\text{:}\)

The last \(n-r\) rows of \(V^T\) are eigenvectors of \(A^TA\) which correspond to an eigenvalue of \(\underline{\hspace{1in}}\text{.}\)
🔗

🔗
Therefore, \(\left\{\vec v_{r+1}, \dots, \vec v_n\right\}\) form an \(\underline{\hspace{3in}}\) for \(\underline{\hspace{2in}}\text{.}\)
🔗

🔗
We have that \(\Nul(A^TA) = \Nul(A)\text{.}\) Thus, the last \(n - r\) rows of \(V^T\) form an \(\underline{\hspace{3in}}\) for \(\underline{\hspace{1.5in}}\text{.}\)
🔗

🔗
Since \(V^T\) is an orthogonal matrix, we know that it’s \(n\) rows form an orthonormal basis for \(\mathbb{R}^n\text{.}\) Therefore, \(\left\{\vec v_1, \dots, \vec v_r\right\}\) is a set of \(\underline{\hspace{2in}}\) vectors which live in \(\underline{\hspace{3in}}\text{.}\)
🔗

🔗
It follows that the first \(r\) rows of \(V^T\) form an \(\underline{\hspace{3in}}\) for \(\underline{\hspace{1.5in}}\text{.}\)
🔗

🔗

🔗

Theorem 7.42.

Let \(A\) be an \(m \times n\) matrix with rank \(r\) and with SVD \(A = U\Sigma V^T\text{,}\) where

\begin{equation*} U = \begin{bmatrix} \vec u_1 \amp \cdots \amp \vec u_m \end{bmatrix}, \quad \Sigma = \begin{bmatrix} \sigma_1 \amp 0 \amp 0 \amp \cdots \amp 0 \\ 0 \amp \sigma_2 \amp 0 \amp \cdots \amp 0 \\ \vdots \amp \vdots \amp \ddots \amp \vdots \amp \vdots \\ 0 \amp 0 \amp \cdots \amp \sigma_r \amp 0 \\ 0 \amp 0 \amp \cdots \amp 0 \amp 0 \end{bmatrix}, \quad V^T = \begin{bmatrix} \vec v_1^T \\ \vdots \\ \vec v_n^T \end{bmatrix} \end{equation*}

Then,

\(A\vec v_i = \sigma_i \vec u_i\) for \(i = 1, \cdots, r\text{.}\)
🔗

🔗
\(\left\{\vec v_1, \dots, \vec v_r\right\}\) is an orthonormal basis for \(\Row(A)\text{.}\)
🔗

🔗
\(\left\{\vec u_1, \dots, \vec u_r\right\}\) is an orthonormal basis for \(\Col(A)\text{.}\)
🔗

🔗
\(\left\{\vec v_{r+1}, \dots, \vec v_n\right\}\) is an orthonormal basis for \(\Nul(A)\text{.}\)
🔗

🔗
\(\left\{\vec u_{r+1}, \dots, \vec u_m\right\}\) is an orthonormal basis for \(\Nul(A^T)\text{.}\)
🔗

🔗

🔗

A graphic depicting the relationship between the rows of V transpose and the columns of U.

🔗

Example 7.43.

Let \(A = U\Sigma V^T\) for

\begin{equation*} U = \begin{bmatrix} 0 \amp 0 \amp 1 \amp 0 \\ 0 \amp 1 \amp 0 \amp 0 \\ 0 \amp 0 \amp 0 \amp -1 \\ 1 \amp 0 \amp 0 \amp 0 \end{bmatrix}, \quad \Sigma = \begin{bmatrix} 4 \amp 0 \amp 0 \amp 0 \amp 0 \\ 0 \amp 3 \amp 0 \amp 0 \amp 0 \\ 0 \amp 0 \amp \sqrt{5} \amp 0 \amp 0 \\ 0 \amp 0 \amp 0 \amp 0 \amp 0 \end{bmatrix}, \quad V^T = \begin{bmatrix} 0 \amp 1 \amp 0 \amp 0 \amp 0 \\ 0 \amp 0 \amp 1 \amp 0 \amp 0 \\ \sqrt{0.2} \amp 0 \amp 0 \amp 0 \amp \sqrt{0.8} \\ 0 \amp 0 \amp 0 \amp 1 \amp 0 \\ -\sqrt{0.8} \amp 0 \amp 0 \amp 0 \amp \sqrt{0.2} \end{bmatrix}\text{.} \end{equation*}

Find an orthonormal basis for each of the four subspaces.

🔗

We can create the spectral decomposition for a general matrix \(A\) in a similar way to Theorem 7.19.

🔗

Theorem 7.44.

Let \(A = U\Sigma V^T\) be the SVD of a rank \(r\) matrix, with the columns of \(U\) given by \(\vec u_i\text{,}\) the entries on the diagonal of \(\Sigma\) given by \(\sigma_i\text{,}\) and the rows of \(V^T\) given by \(\vec v_i\text{.}\) Then, the spectral decomposition of \(A\) is given by

\begin{equation*} A=\sum_{i=1}^r \sigma_i\,\vec u_i\vec v_i^T. \end{equation*}

🔗

Example 7.45.

Find the spectral decomposition of \(A = \begin{bmatrix}3 \amp 8 \\ 0 \amp 3 \end{bmatrix}\text{.}\) Use the decomposition to identify the best rank 1 approximation of \(A\text{.}\)

🔗

We end this section, and this course, by investigating the connection between Spectral Decomposition and Image Compression. Consider an image in grayscale consisting of pixels in a \(880 \times 1100\) grid. Since each pixel is in grayscale, we can assign it an integer value between 0 and 255 depending on how gray it is, where 0 is black and 255 is white. We can form a \(880 \times 1100\) matrix where each entry is the grayscale value of the associated pixel. In general, the rank of this matrix could be quite high.

🔗

For example, suppose we have selected an image and created the associated matrix and it has rank \(776\text{.}\) One way to store this image in a slightly compressed manner is to store the \(776\) basis vectors for \(\Col(A)\text{,}\) and then store a list of linear combinations describing how to use these \(776\) vectors to recreate all \(880\) columns (the dimensions of an image are reversed, first horizontal, then vertical).

🔗

Question 7.46.

Could we compress the file size further? Suppose we were okay with a little error in our image. Perhaps just a hair of fuzziness that could be corrected open reopening the compressed file. How could we go about doing this?

🔗

The idea is to find a good enough approximation of the image. We can use Spectral Decomposition to find the best rank \(j\) approximation for \(j = 1, \dots, r\text{.}\) We can save tons of storage by letting \(j\) be small, but then our image will be super blurry. We can have a very accurate image by letting \(j\) be large, but then the compression’s filesize is not much smaller than the original’s filesize. We need to find a \(j\) which balances these, a good approximation with small filesize.

🔗

Here is the best rank 1 approximation for the image. Can you tell what the image is?