主成分分析 PCA

本章问题: 你有 50 个特征描述一个用户 (年龄、收入、点击次数、收藏数、登录频率...), 画图只能画 2 维。怎么办? PCA 是答案——找到"数据散得最开"的方向, 把高维数据投影过去, 信息损失最小。

1. 为什么需要降维?

实际工作中高维数据到处都是:

场景	维度	问题
用户画像	50-500 维	画不出图, 模型训练慢
图像	1024×1024 像素	太多冗余 (相邻像素相似)
基因表达	2 万维	维度远大于样本数
文本 TF-IDF	词表大小 (万级)	稀疏 + 噪声多

降维的3 个核心目的:

可视化 — 降到 2D/3D 画散点图
去噪 — 小的特征值方向往往是噪声
加速后续模型 — 维度降 10 倍, 训练快 10 倍

💡 关键直觉: 数据虽然有 50 维, 但很多维度是相关的 (比如"收入"和"消费金额"高度相关)。真正"独立的信息"可能只有 5-10 个。PCA 就是帮你找到这 5-10 个。

2. PCA 的几何直觉

2.1 找"数据散得最开"的方向

想象 2D 平面上有一堆点, 长成一个椭圆:

把所有点投影到 PC1 这一条线上, 就完成了 2D → 1D 降维。投影后的方差 = PC1 对应的特征值 (越大越好)。

🎯 PCA 的核心目标: 找一个方向, 让数据投影过去之后方差最大。方差大 = 信息多。

2.2 为什么是方差最大 = 信息最多?

方差 = 数据"分散程度" = 信息量。

数据全聚在一起 (方差 0) → 没信息, 全是同一个点
数据散得很开 (方差大) → 区分度高, 信息多

但方差大≠区分度好 (分类时还需要看类别, PCA 是无监督, 不看 y)。

2.3 多个主成分

n 维数据 → 最多 n 个主成分, 按方差从大到小排:

PC1 (方差最大)
PC2 (方差次大, 与 PC1 正交)
PC3 (方差第三, 与 PC1、PC2 正交)
...

通常前 2-3 个主成分就解释了 80%+ 方差, 后面的就是噪声。

import numpy as np
import matplotlib.pyplot as plt

# 生成椭圆分布的 2D 数据
np.random.seed(42)
mean = [0, 0]
cov = [[3, 1.5], [1.5, 1]]  # 协方差矩阵
X = np.random.multivariate_normal(mean, cov, 300)

# 计算主成分方向 (= 协方差矩阵的特征向量)
eig_vals, eig_vecs = np.linalg.eigh(np.cov(X.T))
# eig_vecs[:, -1] 是最大特征值对应的特征向量 (PC1)
pc1 = eig_vecs[:, -1]
pc2 = eig_vecs[:, -2]

fig, ax = plt.subplots(figsize=(8, 6))
ax.scatter(X[:, 0], X[:, 1], alpha=0.5, s=20)

# 画主成分方向 (从原点出发)
for i, (vec, val) in enumerate(zip([pc1, pc2], [eig_vals[-1], eig_vals[-2]])):
    ax.arrow(0, 0, vec[0]*3, vec[1]*3, head_width=0.15,
             color='red' if i == 0 else 'blue', linewidth=2,
             label=f'PC{i+1} (λ={val:.2f})')
ax.axhline(0, color='gray', lw=0.5)
ax.axvline(0, color='gray', lw=0.5)
ax.set_title('2D 数据的两个主成分方向')
ax.legend(); ax.grid(alpha=0.3)
plt.show()

3. 数学: 从几何直觉到公式

3.1 协方差矩阵回顾

对于 d 维数据 X (n 行样本 × d 列特征), 先中心化 (每列减均值):

协方差矩阵是 d × d 矩阵:

Σ[i][j] = 特征 i 和特征 j 的协方差
Σ[i][i] = 特征 i 的方差
Σ 对称、半正定

3.2 特征值分解

PCA 的数学核心: 对 Σ 做特征值分解:

v_i = 第 i 个主成分方向 (特征向量, 单位长度)
λ_i = 第 i 个主成分的方差 (特征值)

把所有 v_i 按 λ 从大到小排成矩阵 W (d × k, 取前 k 列), 降维就是:

📌 PCA = 找协方差矩阵的特征向量, 投影到方差最大的 k 个方向上。

3.3 为什么特征向量就是方差最大的方向?

数学证明 (拉格朗日乘子):

目标: 找 w, 最大化 w^T Σ w, 约束 ‖w‖ = 1

构造 L = w^T Σ w - λ(w^T w - 1) 对 w 求导: 2Σw - 2λw = 0 → Σw = λw ✓

所以每个特征向量都是一个驻点, λ 越大, 该方向的方差越大。λ 是 w 方向上的方差。

4. PCA 的完整步骤

4.1 算法流程

输入: X (n × d), 目标维度 k
1. 中心化: X̃ = X - mean(X)        ← 每列减均值
2. (可选) 标准化: X̃ / std(X)       ← 各特征量纲不同时必做
3. 算协方差矩阵: Σ = (1/(n-1)) X̃ᵀ X̃
4. 特征值分解: Σ = V Λ Vᵀ
5. 取 V 的前 k 列: W = V[:, :k]
6. 投影: Z = X̃ W                  ← n × k, 即降维结果
7. (可选) 反投影重建: X̂ = Z Wᵀ + mean(X)

4.2 标准化必须吗?

必须做的情况:

特征量纲不同 (比如"年龄 0-100"和"收入 0-100万")
否则方差大的特征 (收入) 会主导主成分

不用做:

所有特征同量纲 (比如都是像素值 0-255)

from sklearn.preprocessing import StandardScaler
X_scaled = StandardScaler().fit_transform(X)

4.3 怎么选 k (降多少维)?

两种常用方法:

方法 1: 累计方差贡献率 ≥ 95%

方法 2: 画 scree plot (碎石图), 找"肘部"

特征值从大到小画图, 找曲线急剧变平的位置 = 肘部, 肘部之前的主成分保留。

from sklearn.decomposition import PCA
import matplotlib.pyplot as plt

pca = PCA().fit(X_scaled)
explained = pca.explained_variance_ratio_

# 累计方差贡献率
cumsum = np.cumsum(explained)
fig, axes = plt.subplots(1, 2, figsize=(12, 4))

# 左: 每个主成分贡献率
axes[0].bar(range(1, len(explained)+1), explained)
axes[0].set_xlabel('主成分'); axes[0].set_ylabel('方差贡献率')
axes[0].set_title('每个主成分的解释方差')

# 右: 累计贡献率
axes[1].plot(range(1, len(cumsum)+1), cumsum, 'o-')
axes[1].axhline(0.95, color='red', ls='--', label='95% 阈值')
axes[1].set_xlabel('主成分'); axes[1].set_ylabel('累计贡献率')
axes[1].set_title('累计方差贡献率'); axes[1].legend()
axes[1].grid(alpha=0.3)
plt.tight_layout()
plt.show()

5. sklearn 上手: IRIS 4 维 → 2 维可视化

IRIS 数据 4 维 (花萼长/宽、花瓣长/宽), 人眼看不出 4 维结构。用 PCA 降到 2D:

from sklearn.datasets import load_iris
from sklearn.decomposition import PCA
from sklearn.preprocessing import StandardScaler
import matplotlib.pyplot as plt

# 1. 加载数据
iris = load_iris()
X, y = iris.data, iris.target  # X: (150, 4), y: (150,)

# 2. 标准化 (IRIS 各特征量纲差异大, 必须做)
X_scaled = StandardScaler().fit_transform(X)

# 3. PCA 降到 2 维
pca = PCA(n_components=2)
X_pca = pca.fit_transform(X_scaled)  # (150, 2)

print(f"PC1 解释方差: {pca.explained_variance_ratio_[0]:.2%}")
print(f"PC2 解释方差: {pca.explained_variance_ratio_[1]:.2%}")
print(f"累计: {(pca.explained_variance_ratio_.sum()):.2%}")
# PC1 解释方差: 72.96%
# PC2 解释方差: 22.85%
# 累计: 95.81%  ← 2 维就保留了 95.8% 信息!

# 4. 画散点图
fig, ax = plt.subplots(figsize=(8, 6))
for label, name, color in [(0, '山鸢尾', 'red'), (1, '变色鸢尾', 'green'), (2, '维吉尼亚鸢尾', 'blue')]:
    mask = y == label
    ax.scatter(X_pca[mask, 0], X_pca[mask, 1], label=name, s=50, alpha=0.7)
ax.set_xlabel(f'PC1 ({pca.explained_variance_ratio_[0]:.1%})')
ax.set_ylabel(f'PC2 ({pca.explained_variance_ratio_[1]:.1%})')
ax.set_title('IRIS 4 维 → 2 维 (PCA)')
ax.legend(); ax.grid(alpha=0.3)
plt.show()

结果: 三类花在 2D 图上清晰分开! 而且只用了 2 个主成分, 保留了 95.8% 信息。这就是 PCA 的威力。

6. PCA 跟 SVD 的关系

PCA 还有第二种实现方式, 工业界更常用 (数值稳定 + 处理稀疏矩阵):

# 方式 1: sklearn 默认用 SVD 实现 PCA
from sklearn.decomposition import PCA
pca = PCA(n_components=2)
Z = pca.fit_transform(X_scaled)

# 方式 2: 直接用 scipy 的 SVD
import numpy as np
U, S, Vt = np.linalg.svd(X_scaled, full_matrices=False)
Z = X_scaled @ Vt[:2].T  # 等价

数学等价性: 对中心化后的 X 算 SVD, 右奇异向量 V 就是主成分方向。

# 验证 PCA == SVD 方向
from sklearn.decomposition import PCA
import numpy as np

X_centered = X - X.mean(axis=0)
pca = PCA(n_components=2).fit(X_centered)
U, S, Vt = np.linalg.svd(X_centered, full_matrices=False)

# PCA 的 components_ 就是 V 的前 k 列
print("PCA components_:", pca.components_.shape)
print("SVD Vt 前 2 行:", Vt[:2].shape)
# 方向相同 (可能差符号, 但代表同一根轴)

📌 实际工程: sklearn 的 PCA 默认走 SVD (scipy 路径), 不显式算协方差矩阵。好处是数值稳定 (不会因矩阵接近奇异而出错), 还能处理 n \ll d 的情况 (如基因数据)。

7. 应用场景

7.1 数据可视化 (最常见)

高维 → 2D/3D 画散点图, 直观看出聚类 / 异常 / 类别分布。

# MNIST 784 维 (28×28 像素) → 2 维
from sklearn.datasets import load_digits
digits = load_digits()
X = digits.data  # (1797, 64) 实际是 8x8 像素

pca = PCA(n_components=2)
X_2d = pca.fit_transform(X)
plt.scatter(X_2d[:, 0], X_2d[:, 1], c=digits.target, cmap='tab10', s=10)
plt.colorbar(); plt.title('MNIST 64 维 → 2 维')
plt.show()

7.2 去噪

PCA 只保留大方差方向, 小方差方向往往是噪声。把数据投影到 top-k 主成分再重建, 等于滤掉了噪声。

# 图像去噪: 把图像 X 看作样本, 算 PCA, 重建
pca = PCA(n_components=50).fit(X_noisy)
X_denoised = pca.inverse_transform(pca.transform(X_noisy))

7.3 特征压缩 + 加速模型

训练前先 PCA 降维, 模型快几倍到几十倍:

from sklearn.pipeline import Pipeline
from sklearn.linear_model import LogisticRegression

pipe = Pipeline([
    ('scaler', StandardScaler()),
    ('pca', PCA(n_components=20)),
    ('clf', LogisticRegression(max_iter=1000))
])
pipe.fit(X_train, y_train)

7.4 人脸识别 (Eigenfaces)

经典应用: 把人脸图像展平, 用 PCA 提取"主要人脸特征", 脸库压缩 100 倍。

8. PCA 的局限

8.1 只能抓线性关系

PCA 是线性降维。如果数据分布是弯曲的 (瑞士卷数据集), PCA 会失效。

from sklearn.datasets import make_swiss_roll
X, color = make_swiss_roll(n_samples=1000, noise=0.1, random_state=42)
# PCA 投影会"压扁"瑞士卷, 丢失局部结构
# 这时候用 t-SNE / UMAP / Isomap (非线性降维)

8.2 不看标签, 可能损失分类信息

PCA 是无监督——它只看 X 的方差, 不看 y。最大方差方向不一定是分类最有用的方向。分类任务更推荐 LDA (Linear Discriminant Analysis)。

8.3 主成分难解释

PC1 可能是"收入 × 0.7 + 年龄 × 0.5 + 消费 × 0.3"的线性组合, 没有明确业务含义。要可解释性, 用特征选择 (而非降维)。

9. 跟其他降维方法的对比

方法	类型	保留什么	适合场景
PCA	线性	全局方差最大	通用降维、去噪、压缩
LDA	线性	类间可分性最大	分类前降维 (有标签)
t-SNE	非线性	局部邻域结构	可视化高维 (2D/3D)
UMAP	非线性	局部 + 全局	可视化, 比 t-SNE 快
Isomap	非线性	流形测地距离	瑞士卷这类弯曲数据

💡 降维方法选择: 看全局 + 通用 → PCA; 分类 → LDA; 仅可视化 → t-SNE/UMAP; 流形 → Isomap。

10. 实战检查清单

做 PCA 之前先回答这 5 个问题:

量纲统一了吗? → 没统一先 StandardScaler
想降多少维? → 累计方差 ≥ 95% 或肘部法
数据是线性结构吗? → 否, 考虑 t-SNE/UMAP
有标签吗? → 有标签且做分类, 考虑 LDA
要可解释吗? → 要的话用特征选择 (SelectKBest) 而非降维

小结

概念	一句话
PCA 目标	找方差最大的投影方向
核心数学	协方差矩阵的特征值分解
必做步骤	中心化 + (量纲不同时) 标准化
选 k 方法	累计方差 ≥ 95% 或碎石图肘部
跟 SVD 关系	等价 (sklearn 默认走 SVD)
典型应用	可视化、去噪、特征压缩
最大局限	只能抓线性关系

🎯 一句话总结: PCA 是一种无监督、线性、降维方法, 通过找数据方差最大的方向把高维数据压到低维, 尽可能保留信息。它是机器学习工程师的瑞士军刀——遇到高维问题, 先试 PCA。

练习

基础: 在 IRIS 数据上, 分别用 PCA 降到 1D/2D/3D, 画图比较三类的可分性。
思考: MNIST (784 维) 降到多少维能保留 95% 方差? 试 n_components=10, 50, 100, 200 跑 LogisticRegression, 看准确率变化。
进阶: 给一张带噪声的人脸图像 (用 fetch_olivetti_faces 加载), 用 PCA 重建, 比较 n_components=10, 50, 200 的去噪效果。
对比: 用 t-SNE 也对 MNIST 降维 2D, 比较跟 PCA 的可视化效果 (哪个类分得更开?)。