主成分分析代码实战与降维可视化

小象实战讲义 · 多元统计分析

面对高维数据，变量间的相关性常导致信息冗余与分析困难。本节通过Python实战，系统学习主成分分析（PCA）的核心原理与实现流程，掌握如何将高维相关数据压缩为低维互不相关的综合变量，并实现数据的可视化与解释。

💡 核心导读

• 降维思想：理解PCA通过正交线性变换，提取数据主要变异方向，实现信息浓缩与维度压缩。
• 数学原理：掌握协方差矩阵特征分解与主成分提取的对应关系，理解方差最大化准则。
• 实战流程：学习Python中数据标准化、PCA模型拟合、主成分选择与结果解释的完整步骤。
• 可视化与评估：通过碎石图与重建误差，直观评估降维效果。

一、主成分分析核心原理

问题定义与目标

设 $X = (X_1, X_2, \ldots, X_p)^\top$ 为 $p$ 维随机向量，其均值向量为 $\mu = E(X)$，协方差矩阵为 $\Sigma = \text{Cov}(X) \in \mathbb{R}^{p \times p}$。PCA旨在寻找一组新的综合变量 $Y_1, Y_2, \ldots, Y_m$（$m \leq p$），满足：

1. 每个主成分是原始变量的线性组合：$Y_k = T_k^\top X$，其中 $T_k$ 为单位向量（$T_k^\top T_k = 1$）。
2. 主成分间互不相关：$\text{Cov}(Y_i, Y_j) = 0 \ (i \neq j)$。
3. 方差依次最大化：$Y_1$ 具有最大方差，$Y_2$ 在与 $Y_1$ 不相关的条件下方差最大，依此类推。

核心结论与推导

上述优化问题等价于对协方差矩阵 $\Sigma$ 进行特征值分解。设 $\Sigma$ 的特征值从大到小排列为 $\lambda_1 \geq \lambda_2 \geq \cdots \geq \lambda_p \geq 0$，对应的单位正交特征向量为 $T_1, T_2, \ldots, T_p$，则有：

• 第 $k$ 个主成分：$Y_k = T_k^\top X$
• 主成分方差：$\text{Var}(Y_k) = \lambda_k$
• 主成分协方差：$\text{Cov}(Y_i, Y_j) = 0 \ (i \neq j)$

关键指标与选择准则

1. 方差贡献率：衡量单个主成分的重要性。 $$\phi_k = \frac{\lambda_k}{\sum_{i=1}^p \lambda_i}$$
2. 累计方差贡献率：衡量前 $m$ 个主成分对总体信息的保留程度。 $$\psi_m = \frac{\sum_{i=1}^m \lambda_i}{\sum_{i=1}^p \lambda_i}$$
3. 主成分选择常用准则：
   • Kaiser准则：保留特征值大于1的主成分（适用于标准化数据）。
   • 累计贡献率准则：通常保留累计贡献率达到80%-85%的前 $m$ 个主成分。
   • 碎石图拐点法：观察特征值下降曲线的拐点，拐点后的主成分贡献平缓，可舍弃。

二、Python实战：PCA完整流程与可视化

以下代码演示了使用Python进行PCA的完整流程，包括数据模拟、标准化、模型拟合、主成分选择与结果可视化。

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

# 设置随机种子，确保结果可复现
np.random.seed(42)

# 1. 模拟学生成绩数据：50名学生，6门课程（3理+3文）
n_students, n_courses = 50, 6
# 生成两个潜在因子：逻辑思维（理科相关）和语言能力（文科相关）
logic = np.random.randn(n_students, 1) * 15 + 75
language = np.random.randn(n_students, 1) * 12 + 70

# 构建相关成绩：理科成绩与逻辑因子强相关，文科成绩与语言因子强相关，且文理成绩呈一定负相关
science_scores = logic + np.random.randn(n_students, 3) * 8
arts_scores = language - 0.6*logic + np.random.randn(n_students, 3) * 9

# 合并并限制分数在0-100之间
all_scores = np.clip(np.hstack([science_scores, arts_scores]), 0, 100)
columns = ['数学', '物理', '化学', '语文', '历史', '英语']
df = pd.DataFrame(all_scores, columns=columns)
df.index.name = '学生编号'

print("前5名学生成绩示例：")
print(df.head().round(1))
print("\n" + "="*50)

# 2. 数据标准化（消除量纲影响，等价于基于相关系数矩阵进行PCA）
scaler = StandardScaler()
X_scaled = scaler.fit_transform(df)

# 3. 拟合PCA模型
pca = PCA()
pca.fit(X_scaled)

# 4. 提取核心结果
eigenvalues = pca.explained_variance_          # 特征值，即主成分方差
var_ratio = pca.explained_variance_ratio_      # 方差贡献率
cum_var = np.cumsum(var_ratio)                 # 累计贡献率
loadings = pca.components_.T * np.sqrt(pca.explained_variance_)  # 因子载荷矩阵（载荷 = 特征向量 × √特征值）

# 5. 打印结果摘要
print("PCA分析结果摘要")
print("="*50)
print(f"{'主成分':<6} {'特征值':<8} {'贡献率':<8} {'累计贡献率':<10}")
for i in range(n_courses):
    print(f"PC{i+1:<5} {eigenvalues[i]:<8.3f} {var_ratio[i]:<8.3f} {cum_var[i]:<10.3f}")

# 6. 根据准则确定主成分个数
# Kaiser准则
n_kaiser = int(np.sum(eigenvalues > 1.0))
# 累计贡献率≥80%准则
n_cum = np.argmax(cum_var >= 0.80) + 1
n_comp = max(n_kaiser, n_cum)  # 取两者中较大的，确保信息保留充分

print(f"\n主成分选择建议：")
print(f"  Kaiser准则(λ>1)：保留 {n_kaiser} 个主成分")
print(f"  累计贡献率≥80%：保留 {n_cum} 个主成分")
print(f"  最终建议保留主成分数：{n_comp} 个")
print(f"  前{n_comp}个主成分累计解释方差：{cum_var[n_comp-1]:.3f}")

# 7. 碎石图可视化
plt.figure(figsize=(8, 4))
plt.plot(range(1, n_courses+1), eigenvalues, 'bo-', linewidth=2, markersize=8, label='特征值')
plt.axhline(y=1, color='r', linestyle='--', label='Kaiser准则 (λ=1)')
plt.axvline(x=n_comp, color='g', linestyle=':', label=f'建议保留数 ({n_comp})')
plt.xlabel('主成分序号')
plt.ylabel('特征值（方差）')
plt.title('PCA碎石图 - 特征值下降曲线')
plt.legend()
plt.grid(alpha=0.3)
plt.tight_layout()
plt.show()

# 8. 计算降维后的重建误差（评估信息损失）
pca_selected = PCA(n_components=n_comp)
X_transformed = pca_selected.fit_transform(X_scaled)
X_reconstructed = pca_selected.inverse_transform(X_transformed)
reconstruction_mse = np.mean((X_scaled - X_reconstructed) ** 2)
print(f"\n使用前{n_comp}个主成分进行降维重建的均方误差(MSE)：{reconstruction_mse:.4f}")

代码输出与结果解读

运行上述代码，你将得到类似如下的输出和可视化图形：

前5名学生成绩示例：
        数学    物理    化学    语文    历史    英语
学生编号
0     71.1  79.1  79.7  13.1  32.7  43.5
1     66.5  71.6  76.2  30.9   7.9  17.3
2     99.8  86.1  86.8  22.4   4.7  15.0
3     97.2  82.5  97.6  25.6  10.3  18.1
4     72.0  91.2  69.9  10.3  30.3  37.2

==================================================
PCA分析结果摘要
==================================================
主成分  特征值    贡献率   累计贡献率
PC1    3.748    0.612    0.612
PC2    1.188    0.194    0.806
PC3    0.382    0.062    0.869
PC4    0.318    0.052    0.921
PC5    0.303    0.050    0.970
PC6    0.183    0.030    1.000

主成分选择建议：
  Kaiser准则(λ>1)：保留 2 个主成分
  累计贡献率≥80%：保留 2 个主成分
  最终建议保留主成分数：2 个
  前2个主成分累计解释方差：0.806

使用前2个主成分进行降维重建的均方误差(MSE)：0.1938

结果解读：

1. 特征值与贡献率：第一主成分（PC1）特征值最大（3.748），单独解释了约61.2%的总方差。第二主成分（PC2）解释了19.4%的方差。
2. 主成分个数：
   • Kaiser准则建议保留2个主成分（特征值>1）。
   • 累计贡献率≥80%准则建议保留2个主成分（累计80.6%）。
   • 综合建议保留2个主成分，以平衡简约性与信息保留度。
3. 碎石图：图形清晰显示，特征值在第二个主成分后下降曲线变得平缓，印证了保留2个主成分的合理性。
4. 重建误差：使用前2个主成分重建数据，均方误差仅为0.1938，表明降维导致的信息损失很小。

📝 动手练一练

题目：对某数据集进行PCA，得到特征值列表为：[3.5, 1.8, 0.8, 0.4, 0.3, 0.2]。

1. 根据Kaiser准则，应保留几个主成分？
2. 若要求累计贡献率不低于85%，应保留几个主成分？

参考答案：

1. Kaiser准则：特征值大于1的有3.5和1.8，共2个，因此保留 2个 主成分。
2. 累计贡献率：总特征值和为 3.5+1.8+0.8+0.4+0.3+0.2 = 7.0。
   • 前2个累计贡献率：(3.5+1.8)/7.0 ≈ 0.757 (75.7%)，未达到85%。
   • 前3个累计贡献率：(3.5+1.8+0.8)/7.0 ≈ 0.871 (87.1%)，达到要求。
   • 因此，为满足累计贡献率≥85%，应保留 3个 主成分。

本章小结

本节系统讲解了主成分分析的理论核心与Python实战应用。核心要点回顾：

1. 本质：PCA是一种基于数据协方差矩阵/相关系数矩阵特征分解的线性、正交降维方法。
2. 关键：主成分方向是特征向量，其方差（重要性）由对应的特征值衡量。
3. 流程：数据标准化 → 计算协方差/相关阵 → 特征值分解 → 根据贡献率或碎石图确定主成分数 → 解释主成分意义。
4. 选择：主成分个数需权衡信息保留（累计贡献率高）与模型简约（主成分数少），常用Kaiser准则与累计贡献率准则结合判断。

行动清单：

1. 使用sklearn.decomposition.PCA对你熟悉的一个多维数据集（如鸢尾花数据集sklearn.datasets.load_iris）进行降维。
2. 绘制碎石图，并分别应用Kaiser准则和累计贡献率≥80%准则确定主成分个数。
3. 尝试只保留前两个主成分，将数据投影到二维平面进行可视化，观察类别分离情况。

— 小象教研组