An In-Depth Guide to PCA
Principal Component Analysis (PCA) is an unsupervised learning algorithm. Before diving into PCA itself, it's helpful to understand the broader idea of unsupervised learning.
The General Idea Behind Unsupervised Learning
In unsupervised learning, the goal is to build a data model from feature vectors. A loss function in this setting measures how implausible a data point is under the model. A small loss means the point fits the model well; a large loss suggests the point is unusual or potentially an outlier.
You might have seen the K-means clustering algorithm. It's a classic example of an unsupervised learning algorithm. In K-means clustering, the goal is to learn a set of model parameters, often called cluster centers, so that every data point lies close to at least one of them.
So essentially, the cluster centers are representatives of the data. The loss function in K-means clustering is defined as:
$$l_\Theta(x) = \underset{i=1,\ldots,k}{\min} \| x - \theta_i \|^2$$If this loss is large, the point $x$ isn't close to any of the cluster centers and is thus considered an outlier.
Essentially, a loss function in unsupervised learning characterizes what the data looks like
Linear Algebra, Subspaces, and Distance to a Subspace
Suppose $\theta_1, \ldots, \theta_r$ are vectors in $\mathrm{R}^d$, then the set of all linear combinations of $\theta_1, \ldots, \theta_r$ is called a subspace $\mathcal{S}$ of $\mathrm{R}^d$. Every point in $\mathcal{S}$ has the form $\Theta a$, where
$$\Theta = \begin{bmatrix} \mid & \ldots & \mid \\ \theta_1 & \ldots & \theta_r \\ \mid & \ldots & \mid \\ \end{bmatrix}$$is a $d \times r$ matrix whose columns are the vectors $\theta_1, \ldots, \theta_r$. To see why, understand that the matrix-vector product $\Theta a$ is, by definition, a linear combination of the columns of $\Theta$:
$$\Theta a = \begin{bmatrix} \mid & \ldots & \mid \\ \theta_1 & \ldots & \theta_r \\ \mid & \ldots & \mid \end{bmatrix} \begin{bmatrix} a_1 \\ \vdots \\ a_r \end{bmatrix} = a_1 \begin{bmatrix} \mid \\ \theta_1 \\ \mid \end{bmatrix} + \cdots + a_r \begin{bmatrix} \mid \\ \theta_r \\ \mid \end{bmatrix}$$The right-hand term is clearly a linear combination of vectors $\theta_1, \ldots, \theta_r$ with coefficients coming from the vector $a$. Thus, by definition, $\Theta a \in \mathcal{S}$ for every $a \in \mathrm{R}^r$. So we can represent $\mathcal{S}$ as:
$$\mathcal{S} = \left\{ \Theta a: a \in \mathrm{R}^r \right\}$$The distance of a point $x$ to the subspace $\mathcal{S}$ is the minimum distance from $x$ to any point within that subspace.
To build some geometric intuition, a subspace spanned by two vectors in $\mathrm{R}^3$ (which forms a plane) and the distance of a point $x$ to this subspace is illustrated below:
Let's see how to compute this distance. We want to find a point $\hat{x}$ in $\mathcal{S}$ that is closest to $x$. Mathematically, this can be represented as:
$$\hat{x} = \underset{x^\prime \in \mathcal{S}}{\arg\min} \| x - x^\prime \|_2^2$$We have seen that every point $x^\prime \in \mathcal{S}$ can be represented as $x^\prime = \Theta a$ for some vector $a \in \mathrm{R}^r$. So now our objective is to find an optimal vector $a^*$ that minimize this quantity:
$$a^* = \underset{a \in \mathrm{R}^r}{\arg\min} \| x - \Theta a \|_2^2$$Does this expression look familiar? It's the standard least-squares problem. To stir up your memory, the vector $y$ that minimizes
$$\| Ay - b \|_2^2$$is given by
$$y = A^\dagger b,$$where $A^\dagger$ is the pseudo-inverse of $A$. If $A$ has linearly independent columns, then the formula for $A^\dagger$ is given by:
$$A^\dagger = (A^TA)^{-1}A^T$$In our case, $A$ is $\Theta$, $b$ is $x$, and the vector of coefficients $y$ is $a$. So, the optimal vector $a^*$ we're looking for is given by:
$$a^* = \Theta^\dagger x$$To find the actual closest point $\hat{x}$, we just plug $a^*$ back in:
$$\hat{x} = \Theta a^* = \Theta \Theta^\dagger x$$Therefore, the distance from $x$ to the subspace $\mathcal{S}$ is simply the distance between $x$ and this closest point $\hat{x}$:
$$\begin{align*} \operatorname{dist}(x, \mathcal{S}) &= \| x - \hat{x} \|_2 \\[4pt] &= \| x - \Theta \Theta^\dagger x \|_2^2 \end{align*}$$Introduction to PCA
PCA is an unsupervised learning algorithm. Its underlying data model assumes that every data point $x$ lies close to a subspace $\mathcal{S}$ spanned by a list of $r$ vectors, $\theta_1, \ldots, \theta_r$. These vectors are called the principal components
Therefore, the PCA loss function for a single point $x$ is its squared distance to this subspace, which we derived in the previous section:
$$\begin{align*} l_\Theta(x) &= \operatorname{dist}(x, \mathcal{S})^2 \\[4pt] &= \| x - \Theta \Theta^\dagger x \|_2^2 \end{align*}$$Given a data matrix $X \in \mathrm{R}^{n \times d}$, where each row $x_i^T$ is a data point in $\mathrm{R}^d$, the goal of PCA is to find $r$ vectors $\theta_1, \ldots, \theta_r$ such that every data point lies close to the subspace $\mathcal{S}$ spanned by these vectors. We typically choose $r \ll d$.
A geometric intuition is illustrated below:
In the next section, we'll discuss how to find these $r$ vectors
How to Compute PCA
Our goal is to find the matrix $\Theta$ (whose columns define the subspace) that minimizes the total loss over the entire dataset $X$. This objective is the sum of squared distances for all $n$ data points:
$$\begin{align*} L(\Theta) &= \sum_{i=1}^n l_\Theta(x_i) \\[4pt] &= \sum_{i=1}^n \| x_i - \Theta \Theta^\dagger x_i \|_2^2 \end{align*}$$The subspace $\mathcal{S}$ can be defined by an infinite number of different bases. To find a unique solution, we add a constraint: we require the basis vectors $\theta_1, \ldots, \theta_r$ to be orthonormal.
The key benefit of this constraint is that it allows us to replace the pseudo-inverse $\Theta^\dagger$ with the simple transpose $\Theta^T$. This is because for a matrix $\Theta$ with orthonormal columns, the product $\Theta^T \Theta = I$. This simplifies the pseudo-inverse formula:
$$\begin{align*} \Theta^\dagger &= (\Theta^T \Theta)^{-1} \Theta^T \\[4pt] &= (I)^{-1} \Theta^T \\[4pt] &= \Theta^T \end{align*}$$Substituting this into our objective, it becomes:
$$\begin{align*} L(\Theta) = \sum_{i=1}^n \| x_i - \Theta \Theta^T x_i \|_2^2 \end{align*}$$Now, we want to express this sum in a compact matrix form. The data matrix $X$ is composed of row vectors $x_i^T$:
$$X = \begin{bmatrix} - & x_1^T & - \\ & \vdots & \\ - & x_n^T & - \end{bmatrix}$$Using the fact that the squared norm of a vector is equal to the squared norm of its transpose, we can rewrite our objective:
$$\begin{align*} L(\Theta) &= \sum_{i=1}^n \| x_i - \Theta \Theta^T x_i \|_2^2 \\[2pt] &= \sum_{i=1}^n \| (x_i - \Theta \Theta^T x_i)^T \|_2^2 \\[2pt] &= \sum_{i=1}^n \| x_i^T - (\Theta \Theta^T x_i)^T \|_2^2 \\[2pt] &= \sum_{i{=1}}^n \| x_i^T - x_i^T \Theta \Theta^T \|_2^2 \end{align*}$$This final expression is the sum of the squared Euclidean norms of the rows of the matrix $X - X \Theta \Theta^T$.
By definition, the squared Frobenius norm of a matrix is the sum of the squared Euclidean norms of its individual rows. Therefore, our objective function can be rewritten as:
$$L(\Theta) = \| X - X \Theta \Theta^T \|_F^2$$Let's define a new matrix $A = X \Theta$. The objective can be rewritten as:
$$\begin{align*} L(\Theta) &= \| X - (X \Theta) \Theta^T \|_F^2 \\[4pt] &= \| X - A \Theta^T \|_F^2 \end{align*}$$Consider the matrix product $A \Theta^T$. We have $A = X \Theta$, where $X$ is a tall matrix, and we have assumed $r \ll d$. Thus, the product $A = X \Theta$ is also a tall matrix, as its number of rows is inherited from $X$ (which is large) and its number of columns, $r$, is small. The matrix $\Theta^T$ is a wide matrix because its number of rows is $r$ (small) and its number of columns is $d$ (large). So we can say that $A \Theta^T$ is a tall-wide product.
Note that both $\operatorname{rank}(A) \leq r$ and $\operatorname{rank}(\Theta) \leq r$. This implies $\operatorname{rank}(A \Theta^T) \leq r$. So now our job is to find two low-rank matrices $A$ and $\Theta^T$ that minimize the objective
$$L(A, \Theta) = \| X - A \Theta^T \|_F^2$$Most articles on PCA directly define PCA as minimizing the above objective. However, in this article we have derived it from first principles
Next, I'll discuss two methods to minimize this objective.
SVD Approach
Going forward, we'll rename $A$ as $Z$ and $\Theta^T$ as $W$. So our goal for PCA is to find two low rank matrices, $Z$ and $W$, that minimize the following objective:
$$L(Z, W) = \| X - ZW \|_F^2$$The Eckart–Young theorem states that the closest rank-$r$ matrix approximation to any matrix $X$ is given by its truncated Singular Value Decomposition (SVD).
First, let the full SVD of $X$ be:
$$\begin{align*} X &=U \Sigma V^T \\ &= \sum_{j=1}^d \sigma_j u_j v_j^T \end{align*}$$where $d$ is the rank of $X$.
The Eckart–Young theorem says that the best rank-$r$ approximation $\hat{X}$ (where $r \ll d$) is created by keeping only the top $r$ singular values and their corresponding singular vectors, and chopping off the remaining ones:
$$\hat{X} = \sum_{j=1}^r \sigma_j u_j v_j^T$$Mathematically, this $\hat{X}$ is the solution to the minimization problem:
$$\hat{X} = \underset{\operatorname{rank}(X^\prime) \leq r}{\operatorname{argmin}} \| X - X^\prime \|_F^2$$In matrix form, this truncated SVD is written as $\hat{X} = U_r \Sigma_r V_r^T$.
Our PCA objective is given by
$$L(Z, W) = \| X - ZW \|_F^2$$The theorem tells us this is minimized when we set the product $ZW$ equal to the best rank-$r$ approximation of $X$:
$$ZW = U_r \Sigma_r V_r^T$$Thus, we can obtain a solution for PCA by simply setting $Z = U_r$ and $W = \Sigma_r V_r^T$
Alternating Minimization Technique
The alternating minimization technique finds the matrices $Z$ and $W$ by alternating between the following two steps until convergence. We initialize $W_0$. Then, for each $t = 1, 2, 3, \ldots$ we alternate between:
$$Z^t = \underset{Z}{\text{argmin}} \| X - ZW^{t-1} \|_F^2$$and
$$W^t = \underset{W}{\text{argmin}} \| X - Z^tW \|_F^2$$This technique may look different at first, but you have likely already been exposed to it in K-Means clustering. There, we also perform alternating minimization: we alternate between finding an assignment vector that assigns each data point to its nearest cluster centroid and then finding new cluster centroids based on those assignments, repeating until convergence.
Alright, let's look at each step and think about how we can solve it.
Solving for $W$ (Step 2)
First, let's take a look at the second step. At step $t$, we want to find $W^t$ given a fixed $Z^t$ that minimizes:
$$\| X - Z^tW \|_F^2$$We can rewrite the objective as:
$$\| Z^tW - X \|_F^2$$We recognize this as a standard matrix least squares problem, and the least squares approximate solution for $W$ is given by:
$$W^t = \left( Z^t \right)^\dagger X$$Solving for $Z$ (Step 1)
Now, let's look at the first step. At step $t$, we want to find $Z^t$ given a fixed $W^{t-1}$ that minimizes:
$$\| X - ZW^{t-1} \|_F^2$$Again, using the fact that the Frobenius norm of a matrix is equal to the Frobenius norm of its transpose, this objective can be rewritten. This allows us to re-frame the problem into the standard least-squares form.
$$\begin{align*} \| X - ZW^{t-1} \|_F^2 &= \| ZW^{t-1} - X \|_F^2 \\[4pt] &= \| (ZW^{t-1} - X)^T \|_F^2 \\[4pt] &= \| (W^{t-1})^TZ^T - X^T \|_F^2 \\[4pt] \end{align*}$$This is now a matrix least squares problem for the unknown matrix $Z^T$. The approximate solution is given by:
$$(Z^t)^T = \big( (W^{t-1})^T \big)^\dagger X^T$$Taking the transpose on both sides, and using the fact that the pseudo-inverse and transpose operations can be interchanged, the least squares approximate solution for $Z$ is given by:
$$Z^t = X(W^{t-1})^\dagger$$Finally, we'll cover PCA's most widespread application: dimensionality reduction
PCA for dimensionality reduction
Suppose we have solved for PCA using the methods discussed above and found the matrix $\Theta$. What do we do next?
A simple and widely used idea is to find the embedding of a data point. Consider a $d$-dimensional data point $x$. We want to find its compressed representation.
From our earlier discussion, the projection of $x$ onto the subspace $\mathcal{S}$ is $\hat{x} = \Theta \Theta^T x$. The vector of coefficients that defines this projection is $\Theta^T x$. This $r$-dimensional vector is the embedding.
So, while $x$ is a $d$-dimensional vector, $\Theta^T x$ is an $r$-dimensional vector. Since we have assumed $r \ll d$, the vector $\Theta^T x$ is a compressed representation of $x$. This is the embedding.
Furthermore, this transformation has an approximate isometry property. This means that if two vectors $x_1$ and $x_2$ are close, their corresponding embeddings $\Theta^T x_1$ and $\Theta^T x_2$ are also close. This is a very important property that helps in visualization.
Conclusion
I hope you enjoyed this article. The original EASE paper does not fully establish the context nor does it detail the mathematical derivations. In this article, however, we have discussed the foundational concepts and derived the solution entirely from first principles.