Tiny Volt24. Dimension reduction
The SVD chapter showed that any matrix rotates, scales, and rotates again, with the directions given by the eigenvectors of and . The PSD chapter showed that any matrix of the form is symmetric positive semi-definite, so its eigenvalues are non-negative and its eigenvectors are perpendicular to each other.
One of the more common applications of SVD is dimension reduction. It is used in an algorithm which allows storing -dimensional points of a dataset as -dimensional points where . The key lies in figuring out what matrix to take the SVD of. But before we get into those details, let's start with the problem first as usual.
24.1 Understanding the problem
Let's say you are given a dataset of 20 items. Each item is a 2-dimensional point. Plot . Here's how these points were generated:
import numpy as np
xs = np.linspace(0, 5, 20) - 2.5
ys = 3 + 0.25 * xs + np.random.normal(0, 0.1, size=xs.shape)
You may have noticed that these points roughly lie in a straight line. Now picture this scenario: what if (shown in teal)? How could we use this information?
- In this scenario, the first observation we make is that even though the points are 2-dimensional, they lie in a 1-dimensional subspace. A vector corresponding to this subspace has the same direction as the direction in which these points are arranged.
- In this 1-dimensional space, each point can be described by only one number. Thus we could save these 20 items by using only 20 numbers instead of 40 numbers (1 for the x-value and 1 for the y-value of each point).
Now, even though the points in our dataset are not arranged perfectly in a straight line, we can still approximate the dataset by projecting the points on a straight line (aka a 1D subspace). By making this assumption, we can reduce the dimensionality of the points from 2 to 1. And we don't lose much information if this approximation is reasonable. This is the idea behind dimension reduction. In general:
We can approximate a set of -dimensional points by projecting it on an -dimensional subspace where . Then each point can be represented by an -dimensional vector instead of an -dimensional vector.
In general, we can represent the dataset by an -by- matrix where each column is an -dimensional point and there are such points in the dataset:
24.2 An idea for a solution
Now the question is: how do we find this subspace? In other words, how do we find the vectors that describe this subspace?
Clearly, we want to project each of these points on a subspace. In our case, this subspace is a 1D subspace. But which subspace is suitable for this? Let's visualize this for the possible 1D subspaces in the 2D plane. Draw alongside the dataset points. Then : the teal points are the projections, and each faint line joins a point to its projection. Now and watch how the projections spread out or bunch up. Finally : this is where the projections stay most spread out and the joining lines are shortest, so the least information is lost.
As you may have noticed, the best 1D subspace for projecting the 2D points on is the one which aligns with the direction in which the points vary the most. This direction is roughly given by the vector . In general, here is the strategy for finding the optimal subspace:
- Find the direction in which the data varies the most. This becomes the first 1D subspace you want to project each -dimensional point on.
- Find the direction perpendicular to in which the data varies the most. This becomes the second 1D subspace on which you want to project each of your points.
- Find the direction perpendicular to both in which the data varies the most. This becomes the third 1D subspace to project each point on.
- ...
- Find the direction perpendicular to the directions in which the data varies the most. This becomes the -th 1D subspace to project each point on.
Note that each of these vectors is -dimensional. Once you have the matrix of all such directions , getting the projections on all of these subspaces is as easy as one matrix-matrix multiplication :
The -th column of this matrix holds the projections of the point on the top subspaces in which the dataset varies the most. Also note that this new dataset matrix is now -by- (as compared to -by- before).
24.3 Covariance matrix
Now the question is: how can we find these directions ? The answer lies in the covariance matrix. It is an -by- matrix where the entry at the -th row and -th column, , is the covariance of the -th and -th dimensions of the features.
24.3.1 What is covariance?
Covariance is a function that takes two random variables as inputs. Its output is the expected value of . Roughly, it is a measure of:
- on average, how often and are larger or smaller than their respective means at the same time
- on average, by what amount and vary from their respective means at the same time.
Note that for any random variable , and .
In our example, the points in the dataset are 2-dimensional. The two features are the x-coordinate and the y-coordinate of a point. The covariance matrix in this case looks like:
We can calculate these values from the dataset (with points):
Here and .
However there is a much easier way to calculate the entire covariance matrix at once. First we need to center the dataset:
D_mean = [column_vector - mean_vector for column_vector in D]
To be more precise, if we assign , we have:
The actual code to do so using numpy is:
D = np.stack([xs, ys], axis=-1).T
D_mean = D - D.mean(axis=-1, keepdims=True)
Let's see what centering does. Plot once more, then in teal: the same cloud, shifted so its mean sits at the origin. Once we have , the covariance matrix is , which for our dataset works out to . You can make sure that this is indeed equal to the matrix of and entries above. To see that the covariance matrix has already captured the direction in which the features vary the most, draw : both point along the tilt of the point cloud.
24.4 Covariance matrix is symmetric positive semi-definite
Let's take a detour from the problem and talk about the covariance matrix a bit. This matrix is obviously symmetric. But it is also PSD: any matrix of the form is symmetric PSD, and if is a symmetric PSD matrix then for any non-negative scalar , the matrix is also symmetric PSD. The covariance matrix is exactly of this form.
Since the covariance matrix is PSD, its eigenvalues are all non-negative and its eigenvectors are all perpendicular to each other.
24.5 Complete solution
Now let's come back to the problem. Here comes the most useful insight:
- The eigenvector of corresponding to the largest eigenvalue is the direction in which the features vary the most: .
- The eigenvector corresponding to the second largest eigenvalue is the direction perpendicular to in which the features vary the most: .
- ...
- The eigenvector corresponding to the -th largest eigenvalue is the direction perpendicular to in which the features vary the most: .
In this example dataset, the covariance matrix is and its unit eigenvectors are the columns of with the corresponding eigenvalues and .
Once we have found these vectors which correspond to the 1D subspaces with maximum variance, we can simply take the first of these vectors to create the matrix above and calculate the projections of the -dimensional points on these subspaces. In our example, each point can be approximated by taking its dot product with the first eigenvector . These dot products look like:
This is the same dataset as we saw earlier but with its dimension reduced from 2 to 1.
In the visual below, we project the points on a 1D subspace, but this time we also show the subspaces corresponding to the eigenvectors of the covariance matrix. Draw in rose. Then as before. Now , then : the best approximation is obtained when the input points are projected on the eigenvector with the largest eigenvalue.
The code for the entire process looks like:
num_features = 10
num_points = 100
# this would be the actual dataset you want to reduce the dimensions of
D = np.random.normal(0, 1, size=(num_features, num_points))
D_mean = D - D.mean(axis=-1, keepdims=True)
C = (D_mean @ D_mean.T) / (D.shape[1] - 1)
U, _, _ = np.linalg.svd(C)
num_features_reduced = 5
D_reduced_dim = U[:, :num_features_reduced].T @ D_mean
Note the last line: the SVD of the (symmetric PSD) covariance matrix hands us its eigenvectors in , sorted by decreasing eigenvalue, and one matrix multiplication projects the whole dataset on the top of them.