Singular Value Decomposition and its various interpretations and applications (with interactive examples)

Published: 23 Mar 2024

💡 This is a long read, but don't miss the interactive examples I have sprinkled throughout that might help you with learning!

Introducing SVD using linear transformations

Most of the time when I am about to learn a new technical topic, I would visit Wikiepdia for an introduction first. And so... here is what Wikipedia has to say about Singular Value Decomposition.

📖 "In linear algebra, the singular value decomposition (SVD) is a factorization of a real or complex matrix into a rotation, followed by a rescaling, followed by another rotation." - Wikipedia

In other words, we only need 3 matrices to form any complex matrix. If we visualize this complex matrix as a linear transformation, we can see it as:

Let us see this in action - step through the interactive visualization to transform the black square to fit the skewed grey shape with just the rotation > scale > rotation operations.

As illustrated above, we rotated the shape by $R_1$ = 45°, then we scaled it by $S$ = (1.5, 0.1, 1.0), then we rotated it again by $R_2$ = 45°. This is represented by the following expression:

Let us rewrite this as the key expression that is commonly used to describe the SVD of a matrix:

When we compare the expression $U\Sigma V^T$ with the expression $R_{2}S{R}_1$, where $U=R_2$, $V^T=R_1$, and $\Sigma =S$, we see that:

• as $R_1$ and $R_2$ are rotation matrices, they are orthonormal. Similarly, $U$ and $V$ are orthonormal matrices too, i.e. they are made up of linearly independent unit vectors that are perpendicular to each other, represented by $\left\{u_1\mathrm{...}{u}_n\right\}$ and $\left\{v_1\mathrm{...}{v}_n\right\}$. Their determinants are also equal to 1, as rotation matrices don't change the scale of the space they transform.
• $\Sigma$ is a diagonal matrix that represents "scaling". It is the matrix that determines how the matrix A will transform the scale of the space along various dimensions.

Here's another way to look at this - let $v$ be a vector that we want to transform using the linear transformation $A$. This is simply shown as:

$Av=x$

If we split the resulting vector $x$ as a scalar value with a unit vector, we shall get:

$Av=\sigma u$

Imagine that now instead of just transforming 1 vector, we are going to transform a series of vectors using this expression:

which can be re-arranged as

$AV=U\Sigma$

$A=U\Sigma V^T$

But why does SVD matter? Because it exposes some interesting properties that can be used to solve a variety of problems!

Dimensionality Reduction and Data Compression

Consider this rank 1 matrix

This matrix is of rank = 1 because it only has 1 linearly independent column vector. By decomposing this bigger matrix to 2 smaller ones $u$ and $v^T$, we represent the same matrix with less data (16 vs just 8 numbers!). This becomes significant when we have a bigger, more complex matrix.

Hence the idea is - when we factorize a complex matrix, we can compress the amount of data required to represent it! In this case, we have achieved lossless compression.

Now imagine we have a higher ranked matrix A. We cannot simply reduce it into 1 pair of row matrix and column matrix. But what we can do is to approximate A by using a sum of ranked 1 matrices like so:

$A={\sigma }_1{u}_1{v_1}^T+{\sigma }_2{u}_2{v_2}^T+\mathrm{...\; ...}+{\sigma }_n{u}_n{v_n}^T$

$=\left[\begin{array}{cccc}{u}_1& u_2& \mathrm{...}& u_n\end{array}\right]\left[\begin{array}{cccc}{\sigma }_1& & & \\ & {\sigma }_2& & \\ & & ...& \\ & & & {\sigma }_n\end{array}\right]\left[\begin{array}{c}{v_1}^T\\ {v_2}^T\\ {\mathrm{...}}^{}\\ {v_n}^T\end{array}\right]$

$=U\Sigma V^T$

Again, we arrive at the standard SVD expression.

We see that $\Sigma$ is a diagonal matrix consisting of scalar values "σ". These are known as singular values. Each of them represents the "scale of influence" of its respective pair of left and right singular vectors $u$ and $v$ (i.e. a "singular vector pair").

In $\Sigma$, these singular values σ are arranged in descending order of magnitude. This means that the singular vector pairs earlier in the sequence has a higher influence in approximating the complex matrix $A$.

This is a helpful property of SVD, as that means we can remove singular vector pairs that has the least influence in approximating $A$ by turning their corresponding singular values σ to 0. This helps us achieve lossy compression with minimal information loss!

In order to solidify this concept further, let us look at this interactive example that uses the same skewed parallelogram we had. Use the slider to play with different values of ${\sigma }_2$ and see how it changes the shape.

${\sigma }_1$ {{s3Sigma[0][0]}}

${\sigma }_2$ {{s3Sigma[1][1]}}

$A=\left[\begin{array}{cc}{u}_1& u_2\end{array}\right]\left[\begin{array}{cc}\mathrm{<msub><mi>\sigma </mi><mi>1</mi></msub>}& 0\\ 0& \mathrm{<msub><mi>\sigma </mi><mi>2</mi></msub>}\end{array}\right]\left[\begin{array}{c}{v_1}^T\\ {v_2}^T\end{array}\right]$

$=\left[\begin{array}{cc}\mathrm{cos(45°)}& \mathrm{-sin(45°)}\\ \mathrm{sin(45°)}& \mathrm{cos(45°)}\end{array}\right]\left[\begin{array}{cc}\mathrm{\{\{s3Sigma[0][0]\}\}}& 0.00\\ 0.00& \mathrm{\{\{s3Sigma[1][1]\}\}}\end{array}\right]\left[\begin{array}{cc}\mathrm{cos(45°)}& \mathrm{-sin(45°)}\\ \mathrm{sin(45°)}& \mathrm{cos(45°)}\end{array}\right]$

$=\left[\begin{array}{cc}\mathrm{\{\{s3A[0][0]\}\}}& \mathrm{\{\{s3A[0][1]\}\}}\\ \mathrm{\{\{s3A[1][0]\}\}}& \mathrm{\{\{s3A[0][1]\}\}}\end{array}\right]$

We can see that when we adjust the values for ${\sigma }_1$ and ${\sigma }_2$, we change the scaling factor of their respective dimensions. As the parallelogram has the least variance along the dimension represented by ${\sigma }_2$, when ${\sigma }_2=0$, we compress the parallelogram into a line that is relatively close to the original shape (i.e. having the least information loss).

The fact that the 2D parallelogram has been turned into a 1D line is also of significance! This shows the idea of dimensionality reduction - by turning ${\sigma }_2=0$, we have turned $\Sigma$ from being of rank = 2 to being of rank = 1, which in turn changes the rank of the linear transformation A from being of rank = 2 to being of rank = 1 as well (i.e. both matrices are now made up of only 1 linearly independent column vector)! This shows that the rank of $\Sigma$ actually exposes the rank of A. Cool!

Let us now visualize this concept of lossy compression via dimensionality reduction using another example. Adjust the "quality" slider to see how it compresses the image!

Quality {{s4Slider}}

In the image matrix, we can tell that out of 10 column vectors, there are only 5 that are linearly independent. This is because the left half of the image is the same as the right half.

Hence $\Sigma$ and the image matrix are both of rank = 5 (the last 5 singular values σ are 0), i.e. the last 5 singular vector pairs do not tell us any new patterns or information about the image, and so there is no change to the image when we adjust the "quality" slider between 50% - 100%. But the change to the image quality gets increasingly visible as the affected singular value gets bigger. The first 5 singular vector pairs have a higher influence in approximating our image matrix.

Finding hidden correlations and winning the Netflix Prize

NETFLIX
Prize

On 2 October 2006, Netflix held a competition with a grand prize of US $1,000,000 for anyone who is able to beat Netflix's own algorithm for predicting user ratings for films, based on previous user ratings. Many teams that participated realized that SVD is the basis for the winning algorithm.

Not only did dimensionality reduction become a useful feature (given that the teams were given very large data matrices), SVD could also be used to find hidden patterns or correlations between the rows and columns of a data matrix!

💡 This kind of statistical inference from large data sets is one of the methods for machine learning! I believe this is one of the methods used by various online shoppnig and social media platforms for their recommendation engines.

As you might be able to tell based on intuition from the previous example, each singular vector pair describes the relationship of some unique "characteristics" of the image, and the strength of each relationship is shown in its corresponding singular value.

Another way to put this is that in $A=U\Sigma V^T$:

• $V$ describes the unique characteristics among the columns of $A$,
• $U$ describes the unique characteristics among the rows of $A$, and
• $\Sigma$ describes the strength of the relationship between the unique column-wise and row-wise characteristics described by the corresponding singular vector pair in $U$ and $V$.

Let us look at how this is reflected mathematically:

Let $A=U\Sigma V^T$, and hence $A^T=V{\Sigma }^T{U}^T$

$A^{T}A=V{\Sigma }^T{U}^{T}U\Sigma V^T$

$A^{T}A=V{\Sigma }^2{V}^T$

$(A^{T}A)V=V{\Sigma }^2$

Here $A^{T}A$ is what is known as a covariance matrix (it describes the dot product, and hence the directional relationship) of each column in $A$ with every other columns in $A$.

Also, if you haven't already noticed, $(A^{T}A)V=V{\Sigma }^2$ actually follows the equation for finding eigenvectors and eigenvalues. You should see that the eigenvectors are the column vectors of $V$, and the eigenvalues are the squares of the singular values in $\Sigma$!

This means that each column vector $v$ describes the column-wise characteristics of the data matrix, while the square of each singular value ${\sigma }^2$ describes the strength of that characteristic.

We can use the same steps to derive that $(A{A}^T)U=U{\Sigma }^2$, and hence $U$ describes the row-wise characteristics of the data matrix.

Returning to the "Netflix Prize", remember all we had in our data matrix are each user's ratings for each film. We do not know anything else about the users or about the films (e.g. genre). But with SVD, we can sniff out the "eigen-concept" that describes the relationship between the user and the film, and the strength of that relationship.

Imagine the following data matrix of user ratings, in which you have not rated for any shows yet. Please give a rating of 1-10 for the shows "Star Wars" and "Twilight", and we would recommend you your next watch, be it "Star Trek" or "Lovesick", even though we do not know anything else about the shows and our users beyond just the ratings they have given.

Go ahead, give it a try!

Your next recommended watch: {{userRecommendation}}

💡 Note that we can get more accurate results with more data.

The recommendation for your next watch is discovered via SVD. We factorize the data matrix into 3 separate matrices represented by $U\Sigma V^T$. Based on your inputs, we have shown in the tables below the matrices showing the relationships between users, movies and concepts.

$\Sigma$ - "concept strength" matrix

This is the "concept strength" matrix that is represented by $\Sigma$. Right now, we do not know what each "concept" might mean, except that it refers to some "property" of the movie (e.g. perhaps "genre")? We can get a better guess when we compare this with the $U$ and $V$ matrices.

For now, let us apply dimensionality reduction by removing concepts that have the least significance (i.e. lowest singular values). In our case, as shown in the table, we shall only keep concepts 1 and 2.

$V^T$ - concept to movie matrix

This is the "concept to movie" matrix that is represented by $V^T$. We can see that Twilight and Harry Potter weighs heavily on Concept 1, while Star Wars and Star Trek weighs heavily on Concept 2.

Note that we do not know what Concept 1 or Concept 2 represents. It could represent genre or popularity, or even both!

$U$ - user to concept matrix

This is the "user to concept" matrix that is represented by $U$. Interestingly, we can see that FantasyLove and LoveSickBoy share the same taste for movies belonging to Concept 1, whereas GeekGal likes movies belonging to Concept 2.

And you? {{userConcept}}

Solving linear regression with SVD

You have probably learnt about drawing "best-fit lines" (i.e. linear regression) in Math/Science class in Secondary/middle school. This is an important concept as it helps us build a predictive model that describes the relationship between a dependent variable and its independent variables (for e.g., predicting housing prices based on proximity distances and land area, or stock prices from past prices and trading volume, etc.). As you can see, what we are about to cover has applications in Machine Learning.

Finding the best-fit line analytically was easy, as we could simply use human-sight to estimate a line that has the least average distance from every data point. However, it is impossible for a computer to do that (unless you make this into a Computer Vision problem). The computer can only calculate the best-fit line numerically (this is also known as solving the "linear least squares problem"). And as it turns out, we can solve this problem with, you guess it, SVD!

Here is an example of SVD being used to solve a linear least squares problem. Click/tap on the graph to add new data points and watch the program find the best-fit line!

How is this done? First, let's get a hint from looking at the points you've just drawn in the interactive example - have you noticed that we are actually drawing a line along the direction with the biggest variance among the data points (i.e. the direction with the biggest spread)?

Taking a look at what we learn from "Dimensionality Reduction", we know that each singular value calculated from SVD shows the variance of points along each principle direction. So, let's form our data matrix (let's call it $A$) from the data points you drew, and then work towards applying SVD on this data matrix.

$A=\left[\begin{array}{cc}{x}_1& y_1\\ x_2& y_2\\ ...& ...\\ x_n& y_n\end{array}\right]$

Remember that we are able to use SVD to analyze the rotational and scalar components of a data matrix, but it does not account for "translation". So we need to remove the translational component from our data by centering its data points around their mean. This will give us a new data matrix $\mathrm{A\text{'}}$.

$\mathrm{A\text{'}}=\left[\begin{array}{cc}{x}_1-\bar{x}& y_1-\bar{y}\\ x_2-\bar{x}& y_2-\bar{y}\\ ...& ...\\ x_n-\bar{x}& y_n-\bar{y}\end{array}\right]$

where $\bar{x}=\sum _{i=1}^n\frac{x_i}{n}$ and $\bar{y}=\sum _{i=1}^n\frac{y_i}{n}$

For clarity, here is a code snippet:

// Finding the mean data point
let meanX = 0, meanY = 0;
for(let i = 0; i < data.length ; i++) {
    meanX += data[i][0];
    meanY += data[i][1];
}
meanX /= data.length;
meanY /= data.length;

// Building the new data matrix A' where points are centered around the mean data point
const centeredData = [];
for (let i = 0; i < data.length; i++) {
    centeredData.push([data[i][0] - meanX, data[i][1] - meanY]);
}
            

Following that, we decompose this data matrix using SVD to get $U\Sigma V^T$.

$\mathrm{A\text{'}}=U\Sigma V^T$
$=$ you have not inserted enough points in the interactive example $=\left[\begin{array}{cc}\mathrm{\{\{s4U[0][0]\}\}}& \mathrm{\{\{s4U[0][1]\}\}}\\ \mathrm{\{\{s4U[1][0]\}\}}& \mathrm{\{\{s4U[1][1]\}\}}\end{array}\right]\left[\begin{array}{cc}\mathrm{\{\{s4Q[0][0]\}\}}& \\ & \mathrm{\{\{s4Q[0][1]\}\}}\end{array}\right]\left[\begin{array}{cc}\mathrm{\{\{s4V[0][0]\}\}}& \mathrm{\{\{s4V[1][0]\}\}}\\ \mathrm{\{\{s4V[1][0]\}\}}& \mathrm{\{\{s4V[1][1]\}\}}\end{array}\right]$

We want to select the right singular vector in $V$ that has the highest corresponding singular value in $\Sigma$ (highlighted in red), because it is the direction vector / eigenvector along which the data points has the biggest variance!

Now, we can simply define the best-fit line using the equation:

$y-\bar{y}=\left(\frac{v_y}{v_x}\right)(x-\bar{x}>)$

The steps we have described is also known as Principle Component Analysis.

Solving linear equations with SVD

Let's look at solving the following system of linear equations:

$\{\begin{array}{c}a+7b+3c=0\\ 2a+4b+c=0\\ 4a+8b+6c=0\end{array}$

This can be forumulated in the following form:

$Ax=0⟹\left[\begin{array}{ccc}1& 7& 3\\ 2& 4& 1\\ 4& 8& 6\end{array}\right]\left[\begin{array}{c}a\\ b\\ c\end{array}\right]=\left[\begin{array}{c}0\\ 0\\ 0\end{array}\right]$

This is equivalent to saying that we are looking for the nullspace of A.

Since A is full-ranked and square, and hence invertible:

$x=A^{\mathrm{-1}}0$
$x=0$

That is easy! We see that our linear system only has a trivial solution, i.e., the zero vector is the only vector in the nullspace of A.

However, what if A is not full-ranked and hence not invertible? Consider the following linear system:

$\{\begin{array}{c}a+b+\mathrm{0.25}c=0\\ 2a+4b+c=0\\ 4a+8b+2c=0\end{array}$

where $A=\left[\begin{array}{ccc}1& 1& \mathrm{0.25}\\ 2& 4& 1\\ 4& 8& 2\end{array}\right]$

Since A is not full-ranked (the last 2 columns are linearly dependent), A is not invertible. How then do we calculate the nullspace of A?

The answer - we shall use SVD on A, and the nullspace of A shall be the set of right singular vectors that has zero as their singular values!

But how did we arrive at this conclusion? Let's break this down together.

$A=U\Sigma V^T$

$A=\left[\begin{array}{ccc}\{\{ffsNullU[0][0]\}\}& \{\{ffsNullU[0][1]\}\}& \{\{ffsNullU[0][2]\}\}\\ \{\{ffsNullU[1][0]\}\}& \{\{ffsNullU[1][1]\}\}& \{\{ffsNullU[1][2]\}\}\\ \{\{ffsNullU[2][0]\}\}& \{\{ffsNullU[2][1]\}\}& \{\{ffsNullU[2][2]\}\}\end{array}\right]\left[\begin{array}{ccc}\{\{ffsNullQ[0][0]\}\}& & \\ & \{\{ffsNullQ[0][1]\}\}& \\ & & \{\{ffsNullQ[0][2]\}\}\end{array}\right]\left[\begin{array}{ccc}\{\{ffsNullV[0][0]\}\}& \{\{ffsNullV[1][0]\}\}& \{\{ffsNullV[2][0]\}\}\\ \{\{ffsNullV[0][1]\}\}& \{\{ffsNullV[1][1]\}\}& \{\{ffsNullV[2][1]\}\}\\ \{\{ffsNullV[0][2]\}\}& \{\{ffsNullV[1][2]\}\}& \{\{ffsNullV[2][2]\}\}\end{array}\right]$

$A=\left[\begin{array}{ccc}{u}_1& u_2& u_3\end{array}\right]\left[\begin{array}{cccc}{\sigma }_1& & & \\ & {\sigma }_2& & \\ & & & 0\end{array}\right]\left[\begin{array}{c}{v_1}^T\\ {v_2}^T\\ {v_3}^T\end{array}\right]$

$A\left[\begin{array}{ccc}{v}_1& v_2& v_3\end{array}\right]=\left[\begin{array}{ccc}{u}_1& u_2& u_3\end{array}\right]\left[\begin{array}{cccc}{\sigma }_1& & & \\ & {\sigma }_2& & \\ & & & 0\end{array}\right]$

$A\left[\begin{array}{ccc}{v}_1& v_2& v_3\end{array}\right]=\left[\begin{array}{ccc}{\sigma }_1{u}_1& {\sigma }_2{u}_2& 0.u_3\end{array}\right]$

$A{v}_3=0$

We see that the nullspace of A = span { $v_3$ }, because $v_3$ is the only vector in the SVD that can become the zero vector when transformed by A. This shows that the right singular vectors, with corresponding singular values that are 0, form the nullspace of A.

Alright, now let us look at another similar problem.

Solve $Ax=c$, where $c\ne 0$.

We can solve for $x$ using SVD.

$Ax=U\Sigma V^{T}x$

$Ax=\left[\begin{array}{ccc}{u}_1& u_2& u_3\end{array}\right]\left[\begin{array}{cccc}{\sigma }_1& & & \\ & {\sigma }_2& & \\ & & & 0\end{array}\right]\left[\begin{array}{c}{v_1}^T\\ {v_2}^T\\ {v_3}^T\end{array}\right]x$

$Ax=\left[\begin{array}{ccc}{u}_1& u_2& u_3\end{array}\right]\left[\begin{array}{}{\sigma }_1{v}_1.x\\ {\sigma }_2{v}_2.x\\ 0v_3.x\end{array}\right]$

$Ax=a{u}_1+b{u}_2$ , where a, b ∈ ℝ

This shows that $Ax=\mathrm{span}\{u_1,u_2\}$

And since we knew that

$A\left[\begin{array}{ccc}{v}_1& v_2& v_3\end{array}\right]=\left[\begin{array}{ccc}{\sigma }_1{u}_1& {\sigma }_2{u}_2& 0.u_3\end{array}\right]$

...which means...

$A{v}_1={\sigma }_1{u}_1$ and

$A{v}_2={\sigma }_2{u}_2$

From this, we can infer that $x=''span''\{v_1,v_2\}$, (i.e., spanning the right singular vectors where the singular value is not zero) because any linear combination of $v_1$ and $v_2$, when transformed by $A$, will fall on the plane with its basis described by $''span''\{u_1,u_2\}$.

The Moore-Penrose Pseudoinverse

Lastly, there is another way we can use SVD to solve linear systems, by deriving what we call the Moore-Penrose Pseudoinverse - $A^{+}$.

Again, let us consider a linear system of the form:

$Ax=b$

We can solve for $x$ using the form:

$x=A^{+}b=V{\Sigma }^{+}{U}^{T}b$

where $A^{+}=V{\Sigma }^{+}{U}^T$ is the pseudoinverse.

To calculate the pseudoinverse $A^{+}$, we can obtain $U$ and $V$ via SVD, and obtain ${\Sigma }^{+}$ by replacing each non-zero singular value $\sigma$ in $\Sigma$ with it's reciprocal (i.e. $1$/$\sigma$).

If it happens that $A$ is invertible, you will find that $A^{\mathrm{-1}}=A^{+}$

However, if $A$ is overdetermined (i.e., height > width), you'll get a least-squares solution, and if $A$ is underdetermined (i.e., width > height) you'll get a minimum-norm solution, both helpful in estimating a solution for $x$.

Calculating the SVD of a matrix

After all the talk about various interpretation and application of SVD, we have not even considered how we may calculate the SVD of a matrix! Of course, if you are writing for software, there are already math libraries that do this for any language of your choosing.

• Python - numpy.linalg.svd
• Matlab - svd(A)
• Javascript - SVD-JS by Danilo Salvati, as explained in "Singular Value Decomposition and Least Squares Solutions. By G.H. Golub et al."
• C# - MathNet.Numerics

We will not be going into the details of the algorithms for computing the SVD of a given matrix. Perhaps that can be a separate blog post of its own.

Thank you for reading! It takes time to create content such as this. If you'd like to support free and open education, consider dropping a tip!