Announcements¶

• Homework:
  • Optional bonus homework out, due Sunday
    • Will add up to 5% to your final grade
• Final exam: Monday, May 8 from 12-2pm
• Upcoming office hours:
  • Today: Peer tutor Daniel Cho from 12:30-3:30pm in CCDS 16th floor
  • Today: Abhishek Tiwari from 4-5pm on CCDS 13th floor
• Reading
  • Aggarwal Sections 8.1-8.2

Recap¶

Every matrix $A$ has a singular value decomposition. We can take the SVD of a rank-$r$ matrix $A$,

and create a rank-$k$ approximation by modifying SVD so that

1. $U$ is $m\times k$
2. $V$ is $n\times k$
3. The matrix $\Sigma$ is a $k\times k$ diagonal matrix, whose diagonal values are ${\sigma }_1\ge {\sigma }_2\ge \cdots \ge {\sigma }_k>0$.

Low-effective-rank is common¶

There are two helpful interpretations:

1. Common Patterns
2. Latent Factors

1. Finding common patterns in your data

2. Identifying potential latent factors

Lecture 34: Principal Component Analysis¶

[This lecture is based on Prof. Crovella's CS 132 and CS 506 lecture notes, and Lior Pachter's blog.]

34.1 The Pseudoinverse¶

Consider the case where we are working with the reduced SVD of $A$:

In the reduced SVD, $\Sigma$ is invertible (it is a diagonal matrix with all positive entries on the diagonal).

Using this decomposition we can define an important matrix corresponding to $A$.

This matrix $A^{+}$ is called the pseudoinverse of $A$.

(Sometimes called the Moore-Penrose pseudoinverse).

Obviously, $A$ cannot have an inverse, because it is not even square (let alone invertible) in general.

So why is $A^{+}$ called the pseudoinverse?

Let's go back to our favorite equation, $Ax=b$, specifically in the case where there are no solutions.

In that case, we can find least-squares solutions by finding $\hat{x}$ such that $A\hat{x}$ is the projection of $b$ onto $\mathrm{Col}A$.

And, if $A^{T}A$ is invertible, that $\hat{x}$ is given by $$\hat{x}=(A^{T}A{)}^{-1}{A}^{T}b$$

But, what if $A^{T}A$ is not invertible?

There are still least-square solutions, but now there are an infinite number.

What if we just want to find one of them?

Let's use the pseudoinverse:

Then: $$A\hat{x}=A{A}^{+}b$$

Now, $U$ is an orthonormal basis for $\mathrm{Col}A$.

And, $U^{T}b$ are the coefficients of the projection of $b$ onto each column of $U$, since the columns are unit length.

So, $U{U}^{T}b$ is the projection of $b$ onto $\mathrm{Col}A$.

So, $\hat{x}=A^{+}b$ is a least squares solution of $Ax=b$,

even when $A^{T}A$ is not invertible,

ie, this formula works for any $A$.

Remember, any $A$ has an SVD, and so any $A$ has a pseudoinverse!

34.2 Dimensionality Reduction and PCA¶

In today's lecture we will explore a new way to use SVD: as a way to transform our data objects.

As a reminder, here is what the SVD looks like:

Notice that $U$ contains a row for each object.

In a sense we have transformed objects from an $n$ dimensional space to a $k$ dimensional space, where $k$ is (probably much) smaller than $n$.

This is an example of dimensionality reduction.

When we take our data to be the rows of $U\Sigma$ instead of rows of $A$, we are reducing the dimension of our data from $n$ dimensions to $k$ dimensions.

This suggests a question:

Is there an optimal transformation of the data into $k$ dimensions?

What would that mean?

Here is a natural criterion for the "best" $k$-dimensional transformation:

Find the $k$-dimensional hyperplane that is "closest" to the points.

More precisely:

Given n points in $R^n$, find the hyperplane (affine space) of dimension $k$ with the property that the squared distance of the points to their orthogonal projection onto the hyperplane is minimized.

This sounds like an appealing criterion.

But it also turns out to have a strong statistical guarantee.

In fact, this criterion is a transformation that captures the maximum variance in the data.

That is, the resulting $k$-dimensional dataset is the one with maximum variance.

Let's see why this is the case.

Centroid and Variance¶

First, let's recall the idea of a dataset centroid.

Given a $n\times d$ data matrix $X$ with observations on the rows (as always):

then

In words: the centroid -- ie, the mean vector -- is the average over the rows.

It is the "center of mass" of the dataset.

Next, we define the sample variance of a dataset as:

where $x_j^T$ is row $j$ of $X$.

In other words, the sample variance of the set of points is the average squared distance from each point to the centroid.

Consider if we move the points (translate each point by some constant amount). Clearly, the sample variance does not change.

So, let's move the points to be centered on the origin.

The sample variance of the new points $\tilde{X}$ is the same as the old points $X$, but the centroid of the new point set is the origin.

Now that the mean of the points is the zero vector, we can reason geometrically.

Here is a picture to show why the distance-minimizing subspace is variance-maximizing.

In this figure,

• the red point is one example point from $\tilde{X}$,
• the green point is the origin / centroid, and
• the blue point is the $k$-dimensional projection of the red point.

The length of the black line is fixed -- it is the distance of the original point from the origin.

So the squared length of the black line is this point's contribution to the sample variance.

Now, regardless of how we shift the green line, a right triangle is formed because the projection is orthogonal.

So -- by virtue of the Pythagorean Theorem, the blue line squared plus the red line squared equals the black line squared.

Which means that when we shift the subspace (green line) so as to minimize the squared distances to all the example points $\tilde{X}$ (blue lines) we automatically maximize the squared distance of all the resulting blue points to the origin (red lines).

And the squared distance of the blue point from the origin (red dashed line) is its contribution to the new $k$-dimensional sample variance.

In other words, the distance-minimizing projection is the variance-maximizing projection!

Example. Here is a dataset $X$ in $R^2$.

Next, we compute

This translates each point so that the new mean is the origin.

Now, the last step is to find the $k$-dimensional subspace that minimizes the Euclidean distance between the data (red points) and their projection on the subspace.

In other words, what rank $k$ matrix $X^{\left(k\right)}\in R^{n\times k}$ is closest to $\tilde{X}$? We seek:

We know how to find this matrix -- via the SVD!

So for this case, let's construct the best 1-D approximation of the mean-centered data:

This method is called Principal Component Analysis.

In summary, PCA consists of:

1. Mean center the data, and
2. Reduce the dimension of the mean-centered data via SVD.

This is equivalent to projecting the data onto the hyperplane that captures the maximum variance in the data.

It winds up constructing the best low dimensional linear approximation of the data.

What are "principal components"?

These are nothing more than the columns of $U$ (or the rows of $V^T$). Because they capture the direction of maximum variation, they are called "principal" components.

There are many uses of PCA (and SVD).

• Data compression (discussed last time)
• Visualization
• Denoising
• Anomaly Detection

We won't have time to cover anomaly detection today. The basic idea is that rather than trying to find anomalous points in a high-dimensional dataset, we can use PCA to project the points onto a lower $k$-dimensional dataset and look for anomalies there.

34.3 Using PCA for Visualization and Denoising¶

We will study both visualization and denoising in the context of text processing.

As we have seen, a common way to work with documents is to build a matrix whose rows contain all documents and whose columns contain all terms.

You could imagine encoding a dataset as such:

Notice that this representation does not take into account word order or sentence structure. It's an example of a bag of words approach.

The bag of words encoding for a document is to treat the document as a multiset of words. That is, we simply count how many times each word occurs. It is a "bag" because all the order of the words in the document is lost.

Surprisingly, we can still tell a lot about the document even without knowing its word ordering.

Counting the number of times each word occurs in a document yields a vector of term frequencies.

However, simply using the "bag of words" directly has a number of drawbacks. We will adjust our vectors to account for two of them.

1. If a single document contains many words, then we expect all term frequencies to be high, even for relatively unimportant words. To adjust for this, it often makes sense to normalize the frequency vectors to have unit length.

2. If a single word occurs in many documents, then it's probably not that useful, so again let's try to reduce this word's weight. The idea here is that rare words are more informative than common words.

This is described directly in the scikit-learn documentation:

In a large text corpus, some words will be very [frequent] (e.g. “the”, “a”, “is” in English) hence carrying very little meaningful information about the actual contents of the document.

If we were to feed the direct count data directly to a classifier those very frequent terms would shadow the frequencies of rarer yet more interesting terms.

In order to re-weight the count features into floating point values suitable for usage by a classifier it is very common to use the tf–idf transform.

Tf means term-frequency while tf–idf means term-frequency times inverse document-frequency.

This is a originally a term weighting scheme developed for information retrieval (as a ranking function for search engines results), that has also found good use in document classification and clustering.

Hence, the definition of tf-idf is as follows.

First:

Next, if $N$ is the total number of documents in the corpus $D$ then:

where the denominator is the number of documents in which the term $t$ appears.

And finally:

Matrices that represent text documents will commonly contain TF-IDF scores.

Often, terms are correlated -- they appear together in combinations that suggest a certain "concept".

That is, term-document matrices often show low effective rank -- many columns can be approximated as combinations of other columns.

When PCA is used for dimensionality reduction of documents, it tends to to extract these "concept" vectors.

The application of PCA to term-document matrices is called Latent Semantic Analysis (LSA).

Among other benefits, LSA can improve the performance of clustering of documents.

This happens because the important concepts are captured in the most significant principal components.

Our Dataset: 20 Newsgroups¶

Once again, let's use as an example the text corpus of 18,000 forum posts in 20 newsgroups. We will take a subset of just 3 newsgroups for this example.

Let's read a few of the forum posts within this dataset.

Basic Clustering¶

As a start, let's use $k$-means clustering using tf-idf scores.

We can use the tokenizer within sklearn to compute $n$-grams for $n=1$ and $2$. (An $n$-gram is a set of $n$ consecutive terms.)

We'll compute a document-term matrix dtm.

Let's first cluster the documents using the raw tf-idf scores. This is without any use of PCA, and so includes lots of noisy or meaningless terms.

Let's evaluate the clusters. We'll assume that the newsgroup the article came from is the 'ground truth.' Remember that we studied two metrics to evaluate how good our clustering is:

• Rand Index: If we are testing an algorithm on data for which we know ground truth, the Rand Index allows us to assess the algorithm’s accuracy.
• Silhouette Score: Even without knowing the ground truth, the Silhouette Score compares the distances between points within a cluster, to the distances between points in different clusters. A higher Silhouette Score relates to a model with “better defined” clusters.

Improvement: Stemming¶

One source of noise that we can eliminate (before we use LSA) comes from word endings.

For example: a Google search on 'run' will return web pages on 'running.'

This is useful, because the difference between 'run' and 'running' in practice is not enough to matter.

The usual solution taken is to simply 'chop off' the part of the word that indicates a variation from the base word. This process is called 'stemming.'

A very good stemmer is the "Snowball" stemmer.

You can read more at http://www.nltk.org and http://www.nltk.org/howto/stem.html.

Installation Note: From a cell you need to call nltk.download() and select the appropriate packages from the interface that appears. In particular you need to download: stopwords from corpora and punkt and snowball_data from models.

Let's stem the data using the Snowball stemmer:

And now let's see how well we can cluster on the stemmed data.

So the Rand Index went from 0.816 to 0.846 as a result of stemming.

Demonstrating PCA¶

OK. Now, let's apply PCA.

Our data matrix is in sparse form.

First, we mean center the data. Note that dtm is a sparse matrix, but when we mean center the data we can reduce the sparsity.

Then we use PCA to reduce the dimension of the mean-centered data.

Note that if you have sparse data, you may want to use scipy.sparse.linalg.svds() and for large data it may be advantageous to use sklearn.decomposition.TruncatedSVD().

The principal components (rows of $V^T$) encode the extracted concepts with terms.

Each LSA concept is a linear combination of words.

The rows of $U$ correspond to documents, which are linear combinations of concepts.

Denoising¶

In order to improve our clustering accuracy, we will exclude the less significant concepts from the documents' feature vectors.

That is, we will choose the leftmost $k$ columns of $U$ and the topmost $k$ rows of $V^T$.

The reduced set of columns of $U$ are our new document encodings, and it is those that we will cluster.

It looks like the first 2 principal components have larger singular values than the rest.

Now let's measure how clustering would work with the first few principal components. Remember the ground truth: the data really come from one of 3 different forums.

Note that we can get good accuracy and coherent clusters with just two principal components.

Visualization¶

That's a good thing, because it means that we can also visualize the data well with the help of PCA.

The challenge of visualization is that the data live in a high dimensional space.

We humans can only look at 2 (or maybe 3) dimensions at a time... but it's often not clear which dimensions to look at.

The idea behind using PCA for visualization is that since low-numbered principal components capture most of the variance in the data, these are the "directions" from which it is most useful to inspect the data.

We saw that the first two principal components were particularly large -- let's start by using them for visualization.

Points in this plot have been labelled with their "true" (aka "ground truth") cluster labels.

Notice how clearly the clusters separate and how coherently they present themselves. This is obvious an excellent visualization that is provided by PCA.

Since this visualization is so clear, we can use it to examine the results of our various clustering methods and get some insight into how they differ.

What happens if we misjudge the number of clusters? Let's try forming 6 clusters instead.

What about the other principal components? Are they useful for visualization?

A common approach is to look at all pairs of (low-numbered) principal components, to look for additional structure in the data.

Let's form 2-dimensional plots using any two of the first five principal components. The colors in the figures below correspond to the ground truth of which newsgroup each post belongs to.

We see that using the first two principal components does indeed perform the best job at clustering.