[Python Introduction: Principal Component Analysis Part #3] How many principal components should be kept? Deciding with explained variance ratio and scree plots

Should we use all 10 principal components?

Up until the last article, we understood how PCA works. We learned that by performing eigenvalue decomposition on the covariance matrix, we can find the directions of the principal components and the variance in those directions.

However, this raises a new question.

If the original data is 10-dimensional, we get 10 principal components. So, do we use all 10 of them?

...That wouldn't be dimensionality reduction.

On the other hand, using only one might result in too much information loss.

"In the end, how many should I keep?"

In this article, we will learn the tools to answer this question. The keywords are "explained variance ratio" and "scree plot".

🎯 What you will learn in this article

• The meaning and calculation method of explained variance ratio and cumulative explained variance ratio

• How to draw and interpret a scree plot

• Three criteria for determining the number of principal components (Elbow method, cumulative explained variance ratio, and Kaiser criterion)

• Actually compressing Iris data into 2 dimensions and visualizing it

What is the explained variance ratio? Let's think about it using sales analysis

To understand the explained variance ratio, let's consider a business example.

Suppose there are 10 variables that affect a company's sales, such as "advertising expenses," "number of sales staff," "number of products," "location," and "price."

Applying PCA gives us 10 principal components.

If the first principal component explains 40% of the sales fluctuation, and the second principal component explains 25% of it, then these two alone can explain 65% of the total. The remaining eight only explain 35% combined.

The indicator that represents "what percentage of the total variance each principal component explains" is the explained variance ratio.

Written as a mathematical formula, it looks like this.

$$
\text{Contribution ratio of the k-th principal component} = \frac{\lambda_k}{\sum_{i=1}^{p} \lambda_i}
$$

• Numerator: Eigenvalue of the k-th principal component (variance in that direction)

• Denominator: Sum of all eigenvalues (total variance)

In other words, the larger the eigenvalue of a principal component, the more information it contains.

Preparation: Using the Iris dataset

This time, we will use the classic machine learning dataset, "Iris".

For three types of irises (Setosa, Versicolor, and Virginica), four features have been measured.

The code is as follows.

import numpy as np
import pandas as pd
from sklearn.datasets import load_iris
from sklearn.preprocessing import StandardScaler
import matplotlib.pyplot as plt
import seaborn as sns
import japanize_matplotlib # 日本語表示

np.random.seed(42)

# データ読み込み
iris = load_iris()
X = iris.data
y = iris.target
feature_names = iris.feature_names
target_names = iris.target_names

# DataFrameで確認
df_iris = pd.DataFrame(X, columns=feature_names)
df_iris['species'] = [target_names[i] for i in y]

print(f"サンプル数: {len(df_iris)}")
print(f"特徴量: {feature_names}")
print(f"品種: {list(target_names)}")
print()
print(df_iris.head())

The execution results are as follows.

This is a 150-sample, 4-dimensional dataset. Our goal this time is to compress this into 2 dimensions and visualize how the three species are distributed.

Data Standardization

Before PCA, we standardize the data (convert to mean 0 and standard deviation 1). This is essential preprocessing when variables have different scales.

The code is as follows.

scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)
​
print("標準化後の平均:", X_scaled.mean(axis=0).round(2))
print("標準化後の標準偏差:", X_scaled.std(axis=0).round(2))

The execution results are as follows.

Eigenvalue decomposition of the covariance matrix

We perform eigenvalue decomposition on the standardized data.

The code is as follows.

# 共分散行列を計算
n = len(X_scaled)
cov_matrix = (X_scaled.T @ X_scaled) / (n - 1)

# 固有値分解
eigenvalues, eigenvectors = np.linalg.eig(cov_matrix)
​
# 固有値が大きい順にソート
idx = np.argsort(eigenvalues)[::-1]
eigenvalues = eigenvalues[idx].real
eigenvectors = eigenvectors[:, idx].real
​
print("固有値（各主成分の分散）:")
for i, ev in enumerate(eigenvalues):
    print(f"  第{i+1}主成分: {ev:.4f}")

The execution results are as follows.

Four eigenvalues are displayed. The 1st principal component is the largest, and they decrease in order.

📊 Calculating contribution ratios and cumulative contribution ratios

Now, it is time to calculate the contribution ratios.

The code is below.

# 寄与率を計算
total_variance = np.sum(eigenvalues)
explained_ratio = eigenvalues / total_variance

# 累積寄与率を計算
cumulative_ratio = np.cumsum(explained_ratio)

print(f"{'主成分':<8} {'固有値':>10} {'寄与率':>10} {'累積寄与率':>10}")
print("-" * 45)
for i in range(len(eigenvalues)):
    print(
        f"PC{i+1:<6} {eigenvalues[i]:>10.4f} "
        f"{explained_ratio[i]*100:>9.2f}% "
        f"{cumulative_ratio[i]*100:>9.2f}%")

The execution results are below.

Looking at the results, the cumulative contribution ratio reaches approximately 95% with just the 1st and 2nd principal components.

This means that "even if you reduce 4-dimensional data to 2 dimensions, 95% of the information can be retained." That is an amazing compression ratio!

📈 Drawing a scree plot

A scree plot is a graph of eigenvalues or contribution ratios in the order of the principal components.

"Scree" refers to the debris accumulated on a mountain slope. It was named because the shape of the graph resembles a slope of debris.

The code is below.

fig, axes = plt.subplots(2, 1, figsize=(10, 10))

# 上: 固有値のスクリープロット
ax = axes[0]
x = np.arange(1, len(eigenvalues) + 1)
ax.bar(
    x, eigenvalues, 
    alpha=0.7, color='steelblue'
)
ax.plot(
    x, eigenvalues, 
    'ro-', markersize=8, linewidth=2
)
ax.set_xlabel('主成分')
ax.set_ylabel('固有値')
ax.set_title('スクリープロット（固有値）')
ax.set_xticks(x)
ax.set_xticklabels([f'PC{i}' for i in x])
ax.grid(True, axis='y', alpha=0.3)

# 下: 寄与率と累積寄与率
ax = axes[1]
ax.bar(
    x, explained_ratio * 100, 
    alpha=0.7, color='steelblue', 
    label='寄与率'
)
ax.plot(
    x, cumulative_ratio * 100, 
    'ro-', markersize=8, linewidth=2, 
    label='累積寄与率'
)
ax.axhline(
    y=80, 
    color='green', linestyle='--', linewidth=2, 
    label='80%ライン'
)
ax.set_xlabel('主成分')
ax.set_ylabel('寄与率 (%)')
ax.set_title('寄与率と累積寄与率')
ax.set_xticks(x)
ax.set_xticklabels([f'PC{i}' for i in x])
ax.set_ylim(0, 105)
ax.legend(loc='center right')
ax.grid(True, axis='y', alpha=0.3)

plt.tight_layout()
plt.show()

The execution results are below.

In the graph above, the eigenvalues of the 1st and 2nd principal components are large, and they suddenly become smaller from the 3rd onwards.

In the graph below, you can see at a glance the point where the cumulative contribution ratio exceeds 80%.

🔍 3 criteria for deciding the number of principal components

1. Elbow method

This is a method of looking at the shape of the scree plot and cutting it off at the position of the "elbow" where the curve suddenly flattens out.

In the graph from earlier, the slope becomes suddenly gentler between the 2nd and 3rd principal components. We cut it off at this "elbow" position and adopt the 1st and 2nd principal components.

Pros: Visually easy to understand Cons: The elbow may not always be clearly visible

2. Judgment based on cumulative contribution ratio

This is the most commonly used method. You decide "how much information you want to retain" and adopt the minimum number of principal components that exceeds that percentage.

Commonly used benchmarks are as follows:

• Visualization/Simple analysis: Around 80%

• Preprocessing for analysis: Around 90%

• Minimizing information loss: Around 95%

In this example, since two principal components reach approximately 95%, if the criterion is to "retain 95%," the number of principal components can be determined as 2.

3. Kaiser criterion (eigenvalue > 1)

For standardized data, there is a rule of thumb that "only principal components with eigenvalues greater than 1 should be adopted."

The idea is that principal components with eigenvalues less than 1 contain less information than a single original variable.

In this case, the eigenvalue of the first principal component is well above 1, while the second is less than 1. By this criterion, it is a fairly strict judgment to use only the first principal component.

Which one is correct?

There is no single correct answer.

Decide based on the purpose of the analysis, the nature of the data, and the balance with subsequent processing.

In this case, both the elbow method and the cumulative contribution ratio (95%) support that "two principal components are sufficient," so it is reasonable to adopt two principal components.

🎨 Principal component scores: Viewing data in a new coordinate system

Once the directions of the principal components are known, each data point is projected onto those directions. The coordinates after this projection are called principal component scores.

Written as a mathematical formula, it looks like this:

$$
Z = X \cdot V
$$

• X: Centered data

• V: Matrix of eigenvectors

• Z: Principal component scores (data after dimensionality reduction)

The code is as follows.

# 上位2つの固有ベクトルを取得
W = eigenvectors[:, :2]  # 4×2の行列
​
# 主成分得点を計算
Z = X_scaled @ W  # 150×2の行列
​
print("主成分得点の形状:", Z.shape)

The execution results are as follows.

Compressed from 4 dimensions to 2 dimensions!

Visualizing the 3 varieties with a scatter plot

Finally, we will plot the compressed data on a scatter plot.

The code is as follows.

plt.figure(figsize=(10, 8))

colors = ['#E74C3C', '#2ECC71', '#3498DB']
for i, (species, color) in enumerate(zip(target_names, colors)):
    mask = y == i
    plt.scatter(
        Z[mask, 0], Z[mask, 1],
        c=color, label=species,
        s=80, alpha=0.7, edgecolors='white'
    )

plt.xlabel(f'第1主成分 (寄与率: {explained_ratio[0]*100:.1f}%)')
plt.ylabel(f'第2主成分 (寄与率: {explained_ratio[1]*100:.1f}%)')
plt.title('Irisデータセットの主成分分析')
plt.legend(title='品種')
plt.grid(True, alpha=0.3)
plt.axhline(y=0, color='gray', linestyle='-', linewidth=0.5)
plt.axvline(x=0, color='gray', linestyle='-', linewidth=0.5)
plt.show()

The execution results are as follows.

The 3 varieties are separated quite clearly!

In particular, Setosa (red) is completely separated from the other two varieties. Versicolor and Virginica can also be generally distinguished, although there is some overlap.

We cannot directly see 4-dimensional data, but by compressing it to 2 dimensions, the structure of the data has become visible.

📊 Principal component loadings: Which variables are effective?

Let's look at how much each feature contributes to the principal components. This is called principal component loading.

$$
\text{loading} = \sqrt{\lambda_k} \cdot v_k
$$

This is the eigenvector (direction) multiplied by the standard deviation (spread).

The code is as follows.

# 主成分負荷量を計算
loadings = pd.DataFrame(
    eigenvectors[:, :2] * np.sqrt(eigenvalues[:2]),
    index=feature_names,
    columns=['PC1', 'PC2']
)
​
print("主成分負荷量:")
print(loadings.round(3))

The execution results are as follows.

• PC1: The absolute values of the loadings for petal length, petal width, and sepal length are large → An axis representing "overall flower size"

• PC2: The absolute value of the loading for sepal width is large → An axis representing "sepal width"

Biplot: Visualizing samples and variables simultaneously

A graph that plots both principal component scores (sample positions) and principal component loadings (variable contributions) is called a biplot.

The code is below.

plt.figure(figsize=(10, 8))
​
# サンプルをプロット
for i, (species, color) in enumerate(zip(target_names, colors)):
    mask = y == i
    plt.scatter(Z[mask, 0], Z[mask, 1],
                c=color, label=species,
                s=80, alpha=0.7, edgecolors='white')
​
# 変数の矢印をプロット
scale = 2  # 矢印の見やすさ調整
loadings_scaled = loadings * scale
for i in range(len(loadings)):
    plt.arrow(0, 0, 
              loadings_scaled.iloc[i, 0], 
              loadings_scaled.iloc[i, 1],
              color='black', alpha=0.6, head_width=0.05)
    plt.text(loadings_scaled.iloc[i, 0] * 1.1,
             loadings_scaled.iloc[i, 1] * 1.1,
             loadings.index[i], fontsize=9)
​
plt.xlabel(f'第1主成分 (寄与率: {explained_ratio[0]*100:.1f}%)')
plt.ylabel(f'第2主成分 (寄与率: {explained_ratio[1]*100:.1f}%)')
plt.title('バイプロット')
plt.legend(title='品種')
plt.grid(True, alpha=0.3)
plt.axhline(y=0, color='gray', linestyle='-', linewidth=0.5)
plt.axvline(x=0, color='gray', linestyle='-', linewidth=0.5)
plt.show()

The execution results are below.

How to read the arrows is as follows.

• Long arrows: That variable is well explained by the principal component

• Direction of the arrows: Positive/negative influence on the principal component

• Angle between arrows: Acute angle means positive correlation, obtuse angle means negative correlation

Summary

In this session, we learned how to decide "how many principal components to keep."

Contribution ratio represents what percentage of the total variance each principal component explains. It is calculated by dividing the eigenvalue by the sum of all eigenvalues.

Cumulative contribution ratio represents what percentage of the total information can be retained by the top k principal components.

Scree plot is a tool for visually grasping the importance of principal components by graphing eigenvalues or contribution ratios.

Criteria for deciding the number of principal components include cumulative contribution ratio (80-95%), Kaiser criterion (eigenvalue > 1), and the elbow method.

Principal component scores are the coordinates obtained by projecting data onto the principal component directions. This becomes the data after dimensionality reduction.

Principal component loadings represent how much each variable contributes to the principal component.

📚 Referenced sources

• Sales Analytics - Introduction to Dimensionality Reduction with Python Part 3

Next Episode Preview

So far, we have understood the mechanics of PCA through manual calculations.

However, in practice, we do not perform these calculations every time. Using a library, it can be done in just a few lines.

Next time, we will learn how to quickly execute PCA using scikit-learn. Because you have the knowledge of manual calculations, you will be able to understand the library's output.

Thank you for reading until the end! Clicking the "Like" button would be encouraging 🙌