Principal Component Analysis (PCA)

Cover image is a PCA (Principal Component Analysis) visualization generated by ChatGPT. The prompt used was: “A detailed 2D scatter plot illustrating Principal Component Analysis (PCA). The data points are grouped into multiple clusters with different colors, representing various categories. The X and Y axes are labeled ‘PC1’ and ‘PC2’, with grid lines for reference. A semi-transparent convex hull encases each cluster to emphasize grouping. The background is clean white with subtle shadows. The title ‘PCA Visualization’ is prominently displayed at the top.” 。

Introduction

Principal Component Analysis (PCA) is a commonly used dimension reduction method. It transforms an originally linearly dependent dataset into a set of linearly independent new variables through an orthogonal transformation. This process minimizes the loss of information when projecting the dataset onto these new variables. These new variables are ranked based on their variance, with the one having the highest variance called the first principal component, the second highest called the second principal component, and so on.

PCA is widely used to map high-dimensional data into lower-dimensional space while preserving as much of the original data’s key features and information as possible.

PCA Principles: Maximum Separability and Minimum Reconstruction Error

PCA optimizes a dataset using two key concepts: maximum separability and minimum reconstruction error. Maximum separability ensures that when the data is projected onto a lower-dimensional space, the projected points remain as spread out as possible, preserving the variance in the data. Minimum reconstruction error minimizes the sum of squared distances between each original data point and its corresponding projected point, similar to the least squares method in regression analysis.

Example 1

Suppose we have the following data points:

$$ \mathbf{X} = (X_1, X_2, \cdots, X_p)^\top, $$

where $X_j = (x_{1j}, x_{2j}, \cdots, x_{nj})$, with $i = 1, \cdots, n$, and $j = 1, \cdots, p$, and $\mathbf{X} \in \mathbb{R}^{n \times p}$.

Maximum Separability

From the data points in Example 1, we can draw the M and N line segments as follows.

From the above data points, we can observe that the projections onto the M line segment are more spread out compared to those on the N line segment, meaning that more information is retained.

From a statistical perspective, this can be expressed as follows: the total variance of the data points along variable M is greater than that along variable N. Therefore, we tend to select the variable with the highest total variance as the first principal component (PC1), such as M in Example 1, followed by the second principal component (PC2), and so on.

Minimum Reconstruction Error

From the M and N line segments drawn in Maximum Separability, we can also observe that the sum of squared distances between the data points and their projections is smaller for the M line segment than for the N line segment.

In general, we prioritize the variable with the smallest sum of squared projection distances as the first principal component. In other words, the variable with the minimum projection error (sum of squared Euclidean distances between original points and projected points) is chosen as the first principal component, followed by subsequent components.

By the Pythagorean theorem, we know that when the data centroid is fixed and all data points are equidistant from the centroid, the longer the projection onto a given line segment, the shorter the distance from the projected point to the centroid. Conversely, the shorter the projection, the longer the distance from the projected point to the centroid.

This means that PCA selects directions that simultaneously maximize variance and minimize projection error, demonstrating that Maximum Separability and Minimum Reconstruction Error are fundamentally equivalent.

Centering

Centering transforms the data so that the coordinate origin aligns with the mean of the dataset without altering its distribution. This step is crucial for simplifying matrix operations and facilitating orthogonal transformations in PCA.

In simple terms, given the same data points as in Example 1, we can compute the mean of all data points to determine the original data centroid.

If $\overline{x}_j$ represents the mean of $X_j$, then:

$$ \overline{x}_j = \frac{1}{n} \sum_{i = 1}^n x_{ij}, \quad j = 1, \cdots, p, $$

By subtracting the mean from the original matrix, we obtain the centered data matrix, effectively shifting the data centroid to the origin.

$$ \mathbf{X}_c = \mathbf{X} - \overline{\mathbf{X}} = \begin{pmatrix} x_{11} & x_{12} & \cdots & x_{1p} \\ x_{21} & x_{22} & \cdots & x_{2p} \\ \vdots & \vdots & \ddots & \vdots \\ x_{n1} & x_{n2} & \cdots & x_{np} \end{pmatrix} - \begin{pmatrix} \overline{x}_1 & \overline{x}_2 & \cdots & \overline{x}_p \\ \overline{x}_1 & \overline{x}_2 & \cdots & \overline{x}_p \\ \vdots & \vdots & \ddots & \vdots \\ \overline{x}_1 & \overline{x}_2 & \cdots & \overline{x}_p \end{pmatrix} $$

Standardization

Standardization is used to eliminate the impact of differences in scale and other factors across variables, allowing the data to be evaluated on the same scale (mean = 0, standard deviation = 1). The process of standardization is similar to centering, but it also includes dividing by the standard deviation. Like centering, the standardized dataset retains the same ranking and distribution as the original data. Generally, data preprocessing involves centering first, followed by deciding whether standardization is necessary.

We can also standardize the data points in Example 1. Let $\sigma^2_j$ represent the variance of $X_j$. From the definition of variance, we can express it as:

$$ \sigma^2_j = \frac{1}{n} \sum_{i = 1}^n (x_{ij} - \overline{x}_j)^2, \quad j = 1, \cdots, p. $$

Therefore, the result of standardizing the dataset is

$$ \mathbf{Z} = \frac{\mathbf{X}_c}{\mathbf{\sigma}} = \frac{\mathbf{X}-\overline{\mathbf{X}}}{\mathbf{\sigma}}, $$

where

$$ \mathbf{\sigma} = \begin{pmatrix} \sigma_1 & \sigma_1 & \cdots & \sigma_1 \\ \sigma_2 & \sigma_2 & \cdots & \sigma_2 \\ \vdots & \vdots & \ddots & \vdots \\ \sigma_p & \sigma_p & \cdots & \sigma_p \end{pmatrix}. $$

Data Point Projection

After applying data preprocessing methods such as centering and standardization, the next step is to project the data points onto new variables.

Suppose we have a data vector $x_i$, and we want to project $x_i$ onto $v$, resulting in $w$. From the projection angle $\theta$, we know that the length of $w$ is given by $\|x_i\| \|v\| \cos(\theta)$, where $\|x_i\|$ represents the length of $x_i$. If $v$ is a unit vector, then the length of $w$ can be rewritten as $\|x_i\| \cos(\theta)$.

The above is the projection formula for a single data point.

Now, let’s consider the projection of multiple data points in the dataset. Suppose we use the data points from Example 1, represented as $\mathbf{X} = (X_1, X_2, \cdots, X_p)^\top$, where $X_j = (x_{1j}, x_{2j}, \cdots, x_{nj})$.

Let $\mathbf{\delta} = (\delta_1, \delta_2, \cdots, \delta_p)^\top$ be the projection direction (principal component vector) of $\mathbf{X}$. Then, the projection of all data points can be written as

$$ \mathbf{\mathcal{P}} = \mathbf{\delta}^\top \mathbf{X} = (\delta_1, \delta_2, \cdots, \delta_p) \begin{pmatrix} X_1 \\ X_2 \\ \vdots \\ X_p \end{pmatrix} = \sum_{j=1}^p \delta_j X_j, $$

$\mathbf{\mathcal{P}} = (\mathcal{P}_1, \mathcal{P}_2, \cdots, \mathcal{P}_p)^\top$

where $\mathbf{\delta}$ must satisfy $\sum_{j=1}^p \delta_j^2 = 1$ to ensure that the projected vector is a unit vector.

Variance Calculation

How can we find a suitable projection direction $\mathcal{P}_i$ such that we achieve maximum separability along that projection vector? Given the following conditions

$$ \max_{\{\delta:\|\delta\| = 1\}} Var(\mathbf{\mathcal{P}}) = \max_{\{\delta:\|\delta\| = 1\}} Var(\mathbf{\delta}^\top \mathbf{X}) = \max_{\{\delta:\|\delta\| = 1\}} \mathbf{\delta}^\top Var(\mathbf{X}) \mathbf{\delta}, $$

To find the maximum variance, we must first determine the appropriate $\mathbf{\delta}$. Before doing so, we need to understand the following theorem:

According to the theorem above, we can determine that the projection direction $\mathbf{\delta}$ is given by the eigenvector $\gamma_1$ corresponding to the largest eigenvalue $\lambda_1$ of the covariance matrix $\Sigma = Var(\mathbf{X}) \in \mathbb{R}^{p \times p}$.

Principal Component Analysis

First Principal Component

Returning to the initial example, let the dataset be $\mathbf{X}$. Based on the theorem mentioned in Variance Calculation, we can determine that the first principal component (PC1) is

$$ Y_1 = \gamma_1^\top \mathbf{X}, $$

where $Y_1 \in \mathbb{R}^p$ has the largest variance after projection

$$ \argmax_{\{\gamma_1:\|\gamma_1\| = 1\}} \mathbf{\gamma_1}^\top Var(\mathbf{X}) \mathbf{\gamma_1}, $$

which means that $\gamma_1$ satisfies

$$ \Sigma \gamma_1 = \lambda_1 \gamma_1. $$

From the dataset $\mathbf{X}$, we know that $E(\mathbf{X}) = \overline{\mathbf{X}}$ and $Var(\mathbf{X}) = \Sigma$. By diagonalizing the covariance matrix $\Sigma$, we get

$$ \Sigma = \Gamma \Lambda \Gamma^\top, $$

Thus, PC1 can also be written as

$$ Y_1 = \Gamma^\top_1 (\mathbf{X}-\overline{\mathbf{X}}) = \Gamma^\top_1 \mathbf{X}_{c} \in \mathbb{R}^n. $$

Second Principal Component

Now that we have determined PC1 $Y_1$ and its projection vector $\mathcal{P}_1$, we can similarly find the second principal component (PC2). It is important to note that the projection vector of PC2 must be orthogonal to PC1, meaning

$$ \mathcal{P}_2^\top \mathcal{P}_1 = 0, $$

The variance of PC2 is given by

$$ \argmax_{\left\{\gamma_2:\|\gamma_2\| = 1, \gamma_2^\top \gamma_1 = 0 \right\}} \mathbf{\gamma_2}^\top Var(\mathbf{X}) \mathbf{\gamma_2}. $$

If $\gamma_2$ satisfies these conditions, it defines PC2 as

$$ Y_2 = \Gamma^\top_2 (\mathbf{X}-\overline{\mathbf{X}}) = \Gamma^\top_2 \mathbf{X}_{c} \in \mathbb{R}^n. $$

Other Principal Components

To summarize the above, for the $q$-th principal component, we have

$$ Y_q = \Gamma^\top_q (\mathbf{X}-\overline{\mathbf{X}}) = \Gamma^\top_q \mathbf{X}_{c} \in \mathbb{R}^n, $$

and its projection vector $\mathcal{P}_q$ must be orthogonal to $\mathcal{P}_{q-1}$ (i.e., $\mathcal{P}_q^\top \mathcal{P}_{q-1} = 0$). The total number of principal components does not exceed the number of features in the dataset, $p$.

For the principal component $\mathbf{Y} = \Gamma^\top \mathbf{X}_{c}$, the following property holds:

Python Example

Dataset Introduction

The dataset used in this example is the Boston housing data. This dataset contains housing price information from 1978 for the city of Boston, Massachusetts. It includes 506 observations for each census district in the Boston metropolitan area.

Variable Description

According to the book Applied Multivariate Statistical Analysis, this dataset contains 506 observations and 14 variables. The description for each variable is as follows:

Dataset

Here’s a Python code to load and view the dataset for PCA analysis.

1
2
3
4
5
6
7
8
9

import pandas as pd

# Read the file
file_path = "bostonh.dat"
df = pd.read_csv(file_path, sep = r"\s+", header = None)  # r"\s+" means split by one or more spaces
df.columns = [f"X{i+1}" for i in range(df.shape[1])]

# View the dataset
print(df.head())

The partial information of the dataset is as follows.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13

        X1    X2    X3  X4     X5     X6    X7      X8  X9    X10   X11  \
0  0.00632  18.0  2.31   0  0.538  6.575  65.2  4.0900   1  296.0  15.3   
1  0.02731   0.0  7.07   0  0.469  6.421  78.9  4.9671   2  242.0  17.8   
2  0.02729   0.0  7.07   0  0.469  7.185  61.1  4.9671   2  242.0  17.8   
3  0.03237   0.0  2.18   0  0.458  6.998  45.8  6.0622   3  222.0  18.7   
4  0.06905   0.0  2.18   0  0.458  7.147  54.2  6.0622   3  222.0  18.7   

      X12   X13   X14  
0  396.90  4.98  24.0  
1  396.90  9.14  21.6  
2  392.83  4.03  34.7  
3  394.63  2.94  33.4  
4  396.90  5.33  36.2  

At the same time, use the following code to check the data and whether there are any missing values.

1

print(f"NAs:\n{df.isna().sum()}")

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16

NAs:
X1     0
X2     0
X3     0
X4     0
X5     0
X6     0
X7     0
X8     0
X9     0
X10    0
X11    0
X12    0
X13    0
X14    0
dtype: int64

From the above results, it can be seen that there are no missing values in this dataset.

Define Variables and Standardization

Next, we define the feature variables x and the target variable y. By defining the feature and target variables, after performing PCA, the results (PCs) can be used for further analysis with other models, such as Random Forests. If only PCA is needed, all variables can be defined as feature variables.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13

from sklearn.decomposition import PCA
from sklearn.preprocessing import StandardScaler
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd

# Define Feature Variables and Target Variable
x = df.drop("X14", axis = 1)  # Define all variables except the price as feature variables
y = df["X14"]  # Define the house price as the target variable

# Standardization
scaler = StandardScaler()
x_scaled = pd.DataFrame(scaler.fit_transform(x), columns = x.columns, index=x.index)

PCA

Next, we will perform Principal Component Analysis. Here, we will use scikit-learn’s PCA() function, and the n_components parameter of the PCA() function allows us to decide how much variance to retain. By not specifying the n_components parameter, it means we will retain all components (maximum variance).

1
2
3
4
5
6
7

pca = PCA()  # Retain all PCs
# pca = PCA(n_components = 0.95)  # Set to retain 95% of the variance

print(f"Shape of the data before dimensionality reduction: {x_scaled.shape}")
x_reduced = pd.DataFrame(pca.fit_transform(x_scaled), index = x.index)
print(f"Shape of the data after dimensionality reduction: {x_reduced.shape}")
print(f"Number of principal components selected: {pca.n_components_}")

1
2
3

Shape of the data before dimensionality reduction: (506, 13)
Shape of the data after dimensionality reduction: (506, 13)
Number of principal components selected: 13

Here, we can also intuitively observe that the total number of principal components does not exceed the number of features in the dataset.

Next, we will output the PCA component matrix and explained variance, and display them in a DataFrame for easier viewing.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18

# Create a DataFrame for PCA components, including explained variance and most important features
pca_components_df = pd.DataFrame(
    pca.components_,  # PCA components matrix
    columns = x.columns,  # Original feature names
    index = [f"PC{i+1}" for i in range(pca.n_components_)]  # Naming PCA components
)

# Calculate explained variance
pca_components_df["Explained Variance"] = pca.explained_variance_ratio_

# Calculate cumulative explained variance
pca_components_df["Cumulative Explained Variance"] = pca.explained_variance_ratio_.cumsum()

# Identify the most important feature for each PCA component (feature with the highest absolute weight)
pca_components_df["Most Important Feature"] = pca_components_df.iloc[:, :-3].abs().idxmax(axis = 1)

# Display the result
print(pca_components_df)

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44

            X1        X2        X3        X4        X5        X6        X7  \
PC1   0.250951 -0.256315  0.346672  0.005042  0.342852 -0.189243  0.313671   
PC2  -0.315252 -0.323313  0.112493  0.454829  0.219116  0.149332  0.311978   
PC3   0.246566  0.295858 -0.015946  0.289781  0.120964  0.593961 -0.017675   
PC4   0.061771  0.128712  0.017146  0.815941 -0.128226 -0.280592 -0.175206   
PC5  -0.082157 -0.320617  0.007811 -0.086531 -0.136854  0.423447 -0.016691   
PC6   0.219660  0.323388  0.076138 -0.167490  0.152983 -0.059267  0.071709   
PC7   0.777607 -0.274996 -0.339576  0.074136 -0.199635  0.063940  0.116011   
PC8  -0.153350  0.402680 -0.173932  0.024662 -0.080121  0.326752  0.600823   
PC9   0.260390  0.358137  0.644416 -0.013728 -0.018522  0.047898 -0.067562   
PC10 -0.019369 -0.267527  0.363532  0.006182 -0.231056  0.431420 -0.362779   
PC11  0.109644 -0.262756  0.303169 -0.013927 -0.111319 -0.053162  0.459159   
PC12  0.086761 -0.071425 -0.113200 -0.003983  0.804323  0.152873 -0.211936   
PC13  0.045952 -0.080919 -0.251077  0.035922  0.043630  0.045567 -0.038551   

            X8        X9       X10       X11       X12       X13  \
PC1  -0.321544  0.319793  0.338469  0.204942 -0.202973  0.309760   
PC2  -0.349070 -0.271521 -0.239454 -0.305897  0.238559 -0.074322   
PC3  -0.049736  0.287255  0.220744 -0.323446 -0.300146 -0.267000   
PC4   0.215436  0.132350  0.103335  0.282622  0.168498  0.069414   
PC5  -0.098592  0.204132  0.130461  0.584002  0.345607 -0.394561   
PC6  -0.023439  0.143194  0.192934 -0.273153  0.803455  0.053216   
PC7  -0.103900 -0.137943 -0.314887  0.002324  0.070295  0.087011   
PC8   0.121812 -0.080358 -0.082774  0.317884  0.004923  0.424353   
PC9  -0.153291 -0.470891 -0.176563  0.254428 -0.044898 -0.195221   
PC10  0.171213 -0.021909  0.035168 -0.153430  0.096515  0.600711   
PC11  0.695693 -0.036544  0.104836 -0.174505 -0.019275 -0.271382   
PC12  0.390941 -0.107026 -0.215191  0.209599  0.041723  0.055226   
PC13 -0.018299 -0.633490  0.720233  0.023398 -0.004463  0.024432   

      Explained Variance  Cumulative Explained Variance Most Important Feature  
PC1             0.471296                       0.471296                     X3  
PC2             0.110252                       0.581548                     X4  
PC3             0.095586                       0.677134                     X6  
PC4             0.065967                       0.743101                     X4  
PC5             0.064217                       0.807318                    X11  
PC6             0.050570                       0.857888                    X12  
PC7             0.041181                       0.899069                     X1  
PC8             0.030469                       0.929538                     X7  
PC9             0.021303                       0.950841                     X3  
PC10            0.016941                       0.967783                     X6  
PC11            0.014309                       0.982091                     X8  
PC12            0.013023                       0.995115                     X5  
PC13            0.004885                       1.000000                    X10

From the above results, we can observe that the explained variance of PC1 is four times that of PC2, highlighting the dominant role of PC1 in explaining the dataset’s variability. Upon further inspection of the cumulative explained variance, we notice that as more principal components are added, the rate of cumulative variance explanation slows down. This allows us to make informed decisions on how many principal components to retain after dimensionality reduction, balancing effective explanation of the dataset while avoiding excessive dimensions that lead to redundant information and computational burden. Therefore, based on these results, we can choose to keep only the first few principal components to efficiently explain the data, without retaining all of them.

Below is a plot of the remaining cumulative variance, allowing us to more intuitively observe the remaining cumulative variance and decide how many principal components to retain.

Use the following code to find all the points after the principal component transformation (PC transformation):

1
2
3
4
5

# Find all the points after the principal component (PC) transformation
pc_df = pd.DataFrame()
for i in range(pca.components_.shape[0]):
    pc_df[f"PC{i+1}"] = x.values @ pca.components_[i]
print(pc_df.head())

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20

         PC1        PC2        PC3        PC4         PC5         PC6  \
0  38.890181  32.935324 -51.873961  92.167019  178.427648  382.499065   
1  33.023432  54.798669 -71.207997  83.321389  176.841195  367.268772   
2  26.538735  48.768409 -67.853639  85.185186  178.071394  362.405072   
3  12.756987  47.783518 -72.338822  86.619758  177.275189  358.347339   
4  15.652406  50.778719 -73.709208  85.656877  177.036592  360.899950   

         PC7        PC8        PC9       PC10       PC11       PC12  \
0 -63.253078  32.971896 -64.423458  25.122328  47.667764 -56.089807   
1 -41.167341  40.135316 -59.975715  26.991812  53.496420 -45.713799   
2 -43.914223  27.501830 -57.556203  30.316431  46.747978 -42.276562   
3 -37.954676  20.639702 -56.427659  32.854127  36.984715 -33.688927   
4 -36.574611  26.755055 -57.546991  31.525144  40.145396 -35.216528   

         PC13  
0  206.962371  
1  167.304566  
2  167.918899  
3  154.655918  
4  154.388828  

PC1 Violin Plot

We can use a violin plot to visualize the distribution of PC1 data.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21

import plotly.express as px

fig = px.violin(pc_df,
                 y = "PC1", 
                 title = "PC1 violin plot",
                 labels = {"PC1": "PC1"},
                 points = "all",
                 box = True  # Display quartiles
                 )

fig.update_traces(
    box_width = 0.1  # Adjust the width of the quartile box to make it smaller
)

fig.update_layout(
    width = 1200,  # Set width  
    height = 600,  # Set height
)

fig.write_html("PC1 violin plot.html")
fig.show()

If you are unable to view the interactive PC1 violin plot or need a full-screen view, please click here to access it.

PC1 and PC2 Scatter Plot

We can use a scatter plot to visualize the distribution of PC1 and PC2 data.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15

import plotly.express as px
fig = px.scatter(pc_df,
                 x = "PC1", 
                 y = "PC2", 
                 title = "PC1 and PC2 scatter plot",
                 labels = {"PC1": "PC1", "PC2": "PC2"}
                 )

fig.update_layout(
    width = 600,   # Set width  
    height = 600,  # Set height
)

fig.write_html("PC1 and PC2 scatter plot.html")
fig.show()

If you are unable to view the interactive PC1 and PC2 scatter plot or need a full-screen view, please click here to access it.

Conclusion

PCA is a useful dimensionality reduction method that maps a dataset from complex $p$ dimensions to a $q$-dimensional space while retaining as much information as possible. Using PCA helps reduce the dimensionality of data, thereby lowering computation and storage costs, and improving model efficiency and accuracy. By removing unnecessary noise and correlations between features, PCA can improve the generalization ability of machine learning models and make the data more interpretable.

Environment

• Operating System: Windows 11 24H2

• Programming Language: Python 3.12.9

Further Learning

• The ipynb file used in this article.

See Also

• Dataset Exploration and Analysis

References

• Wolfgang Karl Härdle & Léopold Simar. (2015). Applied Multivariate Statistical Analysis (4th Edition). Springer. https://link.springer.com/book/10.1007/978-3-662-45171-7

• 动画讲编程. (2022, August 2). 机器学习，什么是PCA降维？一个动画告诉你答案 [Video]. YouTube. Retrieved March 25, 2025. https://youtu.be/FI0H_HStNJU

• Tommy Huang. (2018, April 20). 機器/統計學習:主成分分析(Principal Component Analysis, PCA). Retrieved March 25, 2025. https://chih-sheng-huang821.medium.com/機器-統計學習-主成分分析-principle-component-analysis-pca-58229cd26e71

• Jin-Zhou Li. (2024, January 27). 主成分分析 (Principal Component Analysis, PCA). Retrieved March 25, 2025. https://jinzhou.netlify.app/en/posts/principal-component-analysis/

• Quantnet. (2025, January 27). MVA. GitHub. https://github.com/QuantLet/MVA

• Violin Plots in Python. (No Date). Plotly | Graphing Libraries. Retrieved April 1, 2025. https://plotly.com/python/violin/

• 李航. (2022). 機器學習聖經：最完整的統計學習方法. (First Edition). 深智數位股份有限公司.