主成分分析（PCA）
主成分分析（PCA）：数据降维的艺术大师
引言
想象一下，你有一个装满各种物品的行李箱，但航空公司突然说只能带一半的重量。你会怎么办？聪明的做法是：保留最重要的物品，丢弃不重要的，同时尽可能保持行李箱的"完整功能"。
这就是主成分分析（PCA）要解决的问题！它是数据科学中最重要的降维技术之一。
什么是主成分分析（PCA）？
主成分分析（Principal Component Analysis，PCA）是一种数据降维技术，它能够：
- 将高维数据转换为低维数据
- 保留最重要的信息
- 去除冗余和噪声
- 发现数据中的主要变化方向
生活中的类比
摄影师的视角选择 📸
- 选择最能表现主题的角度
- 一张好照片能展现事物的主要特征
地图的制作 🗺️
- 3D的地球表面 → 2D的平面地图
- 保留最重要的地理信息
- 不同投影方式适合不同用途
简历的撰写 📝
- 丰富的人生经历 → 一页纸的简历
- 突出最重要的技能和经验
- 用最少的信息展现最大的价值
核心思想：寻找数据的"主要方向"
什么是主成分？
主成分就是数据变化最大的方向：
- 第一主成分：数据变化最大的方向
- 第二主成分：与第一主成分垂直，变化第二大的方向
- 第三主成分：与前两个都垂直，变化第三大的方向
- 以此类推…
2D示例：理解主成分
import numpy as np
import matplotlib.pyplot as plt
from sklearn.decomposition import PCA
from sklearn.datasets import make_blobs
# 创建示例数据
np.random.seed(42)
X, _ = make_blobs(n_samples=100, centers=1, n_features=2,
cluster_std=2.0, random_state=42)
# 让数据有明显的方向性
rotation_matrix = np.array([[1, 0.5], [0, 1]])
X = X @ rotation_matrix
# 应用PCA
pca = PCA(n_components=2)
X_pca = pca.fit_transform(X)
# 可视化
fig, axes = plt.subplots(1, 3, figsize=(15, 5))
# 原始数据
axes[0].scatter(X[:, 0], X[:, 1], alpha=0.7)
axes[0].set_title('原始数据')
axes[0].set_xlabel('特征1')
axes[0].set_ylabel('特征2')
axes[0].grid(True)
# 原始数据 + 主成分方向
axes[1].scatter(X[:, 0], X[:, 1], alpha=0.7)
# 绘制主成分方向
mean_point = np.mean(X, axis=0)
for i, (component, variance) in enumerate(zip(pca.components_, pca.explained_variance_)):
axes[1].arrow(mean_point[0], mean_point[1],
component[0] * variance, component[1] * variance,
head_width=0.3, head_length=0.3,
fc=f'C{i+1}', ec=f'C{i+1}',
label=f'主成分{i+1}')
axes[1].set_title('原始数据 + 主成分方向')
axes[1].set_xlabel('特征1')
axes[1].set_ylabel('特征2')
axes[1].legend()
axes[1].grid(True)
# PCA变换后的数据
axes[2].scatter(X_pca[:, 0], X_pca[:, 1], alpha=0.7)
axes[2].set_title('PCA变换后的数据')
axes[2].set_xlabel('第一主成分')
axes[2].set_ylabel('第二主成分')
axes[2].grid(True)
plt.tight_layout()
plt.show()
print("主成分解释的方差比例:", pca.explained_variance_ratio_)
print("累计解释的方差比例:", np.cumsum(pca.explained_variance_ratio_))PCA的数学原理（简化版）
核心目标
PCA要找到一组新的坐标轴（主成分），使得：
- 数据在新轴上的方差最大
- 各个新轴互相垂直（正交）
- 按方差大小排序
步骤分解
- 数据中心化：将数据移动到原点
- 计算协方差矩阵：衡量特征之间的关系
- 特征值分解：找到主要的变化方向
- 选择主成分：保留最重要的几个方向
- 数据变换：将原始数据投影到新的坐标系
def manual_pca(X, n_components):
"""手动实现PCA算法"""
# 1. 数据中心化
X_centered = X - np.mean(X, axis=0)
print("步骤1: 数据中心化完成")
# 2. 计算协方差矩阵
cov_matrix = np.cov(X_centered.T)
print("步骤2: 协方差矩阵计算完成")
print("协方差矩阵:")
print(cov_matrix)
# 3. 特征值分解
eigenvalues, eigenvectors = np.linalg.eigh(cov_matrix)
print("步骤3: 特征值分解完成")
# 4. 按特征值大小排序（降序）
idx = np.argsort(eigenvalues)[::-1]
eigenvalues = eigenvalues[idx]
eigenvectors = eigenvectors[:, idx]
print("特征值（按大小排序）:", eigenvalues)
print("解释的方差比例:", eigenvalues / np.sum(eigenvalues))
# 5. 选择前n_components个主成分
selected_eigenvectors = eigenvectors[:, :n_components]
# 6. 数据变换
X_pca = X_centered @ selected_eigenvectors
return X_pca, selected_eigenvectors, eigenvalues
# 使用手动实现的PCA
X_manual, components, eigenvals = manual_pca(X, n_components=2)
print(f"\n手动PCA结果与sklearn PCA结果的差异:")
print(f"最大差异: {np.max(np.abs(np.abs(X_manual) - np.abs(X_pca))):.10f}")详细代码实现
完整的PCA类
class SimplePCA:
"""简化版PCA实现"""
def __init__(self, n_components):
self.n_components = n_components
self.components_ = None
self.explained_variance_ = None
self.explained_variance_ratio_ = None
self.mean_ = None
def fit(self, X):
"""训练PCA模型"""
# 保存均值用于中心化
self.mean_ = np.mean(X, axis=0)
# 数据中心化
X_centered = X - self.mean_
# 计算协方差矩阵
cov_matrix = np.cov(X_centered.T)
# 特征值分解
eigenvalues, eigenvectors = np.linalg.eigh(cov_matrix)
# 按特征值大小排序
idx = np.argsort(eigenvalues)[::-1]
eigenvalues = eigenvalues[idx]
eigenvectors = eigenvectors[:, idx]
# 保存结果
self.components_ = eigenvectors[:, :self.n_components].T
self.explained_variance_ = eigenvalues[:self.n_components]
self.explained_variance_ratio_ = self.explained_variance_ / np.sum(eigenvalues)
return self
def transform(self, X):
"""将数据变换到主成分空间"""
X_centered = X - self.mean_
return X_centered @ self.components_.T
def fit_transform(self, X):
"""训练并变换数据"""
return self.fit(X).transform(X)
def inverse_transform(self, X_pca):
"""从主成分空间变换回原始空间"""
return X_pca @ self.components_ + self.mean_
# 测试自制PCA
simple_pca = SimplePCA(n_components=2)
X_simple = simple_pca.fit_transform(X)
print("自制PCA的解释方差比例:", simple_pca.explained_variance_ratio_)实际应用示例
1. 图像压缩
from sklearn.datasets import fetch_olivetti_faces
def pca_image_compression():
"""使用PCA进行图像压缩"""
# 加载人脸数据集
faces = fetch_olivetti_faces()
X_faces = faces.data # 每行是一个64x64=4096维的人脸图像
print(f"原始数据形状: {X_faces.shape}")
print(f"每个图像的维度: {X_faces.shape[1]}")
# 应用不同压缩比的PCA
compression_ratios = [50, 100, 200, 400]
fig, axes = plt.subplots(2, len(compression_ratios) + 1, figsize=(15, 6))
# 显示原始图像
original_image = X_faces[0].reshape(64, 64)
axes[0, 0].imshow(original_image, cmap='gray')
axes[0, 0].set_title('原始图像\n(4096维)')
axes[0, 0].axis('off')
# 显示不同压缩比的结果
for i, n_components in enumerate(compression_ratios):
# PCA压缩
pca = PCA(n_components=n_components)
X_compressed = pca.fit_transform(X_faces)
X_reconstructed = pca.inverse_transform(X_compressed)
# 重构图像
reconstructed_image = X_reconstructed[0].reshape(64, 64)
# 计算压缩比和误差
compression_ratio = X_faces.shape[1] / n_components
mse = np.mean((X_faces[0] - X_reconstructed[0])**2)
explained_variance = np.sum(pca.explained_variance_ratio_)
# 显示结果
axes[0, i+1].imshow(reconstructed_image, cmap='gray')
axes[0, i+1].set_title(f'{n_components}维\n压缩比: {compression_ratio:.1f}:1')
axes[0, i+1].axis('off')
axes[1, i+1].bar(['保留信息', '丢失信息'],
[explained_variance, 1-explained_variance])
axes[1, i+1].set_title(f'信息保留: {explained_variance:.1%}')
axes[1, i+1].set_ylim(0, 1)
axes[1, 0].axis('off')
plt.tight_layout()
plt.show()
# pca_image_compression() # 取消注释运行2. 数据可视化
def pca_data_visualization():
"""使用PCA进行高维数据可视化"""
# 创建高维数据
from sklearn.datasets import make_classification
X_high, y = make_classification(n_samples=300, n_features=20,
n_informative=10, n_redundant=10,
n_clusters_per_class=1, random_state=42)
print(f"原始数据维度: {X_high.shape}")
# 降维到2D进行可视化
pca_2d = PCA(n_components=2)
X_2d = pca_2d.fit_transform(X_high)
# 降维到3D
pca_3d = PCA(n_components=3)
X_3d = pca_3d.fit_transform(X_high)
# 可视化
fig = plt.figure(figsize=(15, 5))
# 2D可视化
ax1 = fig.add_subplot(131)
scatter = ax1.scatter(X_2d[:, 0], X_2d[:, 1], c=y, cmap='viridis')
ax1.set_title('PCA 2D可视化')
ax1.set_xlabel('第一主成分')
ax1.set_ylabel('第二主成分')
plt.colorbar(scatter, ax=ax1)
# 3D可视化
ax2 = fig.add_subplot(132, projection='3d')
scatter3d = ax2.scatter(X_3d[:, 0], X_3d[:, 1], X_3d[:, 2], c=y, cmap='viridis')
ax2.set_title('PCA 3D可视化')
ax2.set_xlabel('第一主成分')
ax2.set_ylabel('第二主成分')
ax2.set_zlabel('第三主成分')
# 解释方差比例
ax3 = fig.add_subplot(133)
n_components = min(10, X_high.shape[1])
pca_full = PCA(n_components=n_components)
pca_full.fit(X_high)
ax3.bar(range(1, n_components+1), pca_full.explained_variance_ratio_)
ax3.set_title('各主成分解释的方差比例')
ax3.set_xlabel('主成分编号')
ax3.set_ylabel('解释方差比例')
plt.tight_layout()
plt.show()
print(f"前2个主成分解释的方差: {np.sum(pca_2d.explained_variance_ratio_):.2%}")
print(f"前3个主成分解释的方差: {np.sum(pca_3d.explained_variance_ratio_):.2%}")
pca_data_visualization()如何选择主成分数量？
1. 解释方差比例法
def choose_n_components_variance():
"""基于解释方差比例选择主成分数量"""
# 创建示例数据
from sklearn.datasets import load_digits
X, y = load_digits(return_X_y=True)
# 计算所有主成分
pca_full = PCA()
pca_full.fit(X)
# 计算累计解释方差比例
cumsum_variance = np.cumsum(pca_full.explained_variance_ratio_)
# 绘制解释方差图
plt.figure(figsize=(12, 5))
plt.subplot(1, 2, 1)
plt.plot(range(1, len(cumsum_variance)+1), cumsum_variance, 'bo-')
plt.axhline(y=0.95, color='r', linestyle='--', label='95%阈值')
plt.axhline(y=0.90, color='g', linestyle='--', label='90%阈值')
plt.xlabel('主成分数量')
plt.ylabel('累计解释方差比例')
plt.title('累计解释方差比例')
plt.legend()
plt.grid(True)
plt.subplot(1, 2, 2)
plt.plot(range(1, 21), pca_full.explained_variance_ratio_[:20], 'ro-')
plt.xlabel('主成分编号')
plt.ylabel('单个主成分解释方差比例')
plt.title('前20个主成分的贡献')
plt.grid(True)
plt.tight_layout()
plt.show()
# 找到达到95%解释方差的主成分数量
n_95 = np.argmax(cumsum_variance >= 0.95) + 1
n_90 = np.argmax(cumsum_variance >= 0.90) + 1
print(f"原始维度: {X.shape[1]}")
print(f"达到90%解释方差需要: {n_90}个主成分")
print(f"达到95%解释方差需要: {n_95}个主成分")
print(f"压缩比(95%): {X.shape[1]/n_95:.1f}:1")
choose_n_components_variance()2. 肘部法则
def elbow_method_pca():
"""使用肘部法则选择主成分数量"""
from sklearn.datasets import load_breast_cancer
X, y = load_breast_cancer(return_X_y=True)
# 标准化数据
from sklearn.preprocessing import StandardScaler
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)
# 计算不同主成分数量下的重构误差
n_components_range = range(1, min(21, X.shape[1]))
reconstruction_errors = []
for n in n_components_range:
pca = PCA(n_components=n)
X_pca = pca.fit_transform(X_scaled)
X_reconstructed = pca.inverse_transform(X_pca)
# 计算重构误差
mse = np.mean((X_scaled - X_reconstructed)**2)
reconstruction_errors.append(mse)
# 绘制肘部图
plt.figure(figsize=(10, 6))
plt.plot(n_components_range, reconstruction_errors, 'bo-')
plt.xlabel('主成分数量')
plt.ylabel('重构误差(MSE)')
plt.title('PCA肘部法则')
plt.grid(True)
# 寻找肘部点
differences = np.diff(reconstruction_errors)
second_differences = np.diff(differences)
elbow_point = np.argmax(second_differences) + 2
plt.axvline(x=elbow_point, color='r', linestyle='--',
label=f'肘部点: {elbow_point}个主成分')
plt.legend()
plt.show()
print(f"建议的主成分数量: {elbow_point}")
elbow_method_pca()PCA的优缺点分析
优点 ✅
- 降维效果好：能显著减少数据维度
- 去除冗余：自动去除特征间的相关性
- 数学基础牢固：基于线性代数理论
- 计算效率高：算法复杂度相对较低
- 可解释性强：主成分有明确的数学含义
缺点 ❌
- 线性假设：只能捕捉线性关系
- 全局方法：需要所有数据才能计算
- 可解释性有限：主成分通常是原特征的复杂组合
- 对缩放敏感：需要预先标准化数据
- 信息丢失：降维必然丢失一些信息
何时使用PCA？
def when_to_use_pca():
"""演示何时适合使用PCA"""
# 案例1: 高度相关的特征
print("案例1: 高度相关的特征")
n_samples = 1000
x1 = np.random.randn(n_samples)
x2 = x1 + 0.1 * np.random.randn(n_samples) # 与x1高度相关
x3 = 2 * x1 + 0.1 * np.random.randn(n_samples) # 与x1线性相关
X_correlated = np.column_stack([x1, x2, x3])
pca_corr = PCA()
pca_corr.fit(X_correlated)
print("相关性特征的解释方差比例:", pca_corr.explained_variance_ratio_)
print("第一主成分解释了", f"{pca_corr.explained_variance_ratio_[0]:.1%}", "的方差\n")
# 案例2: 噪声数据
print("案例2: 含噪声的数据")
# 真实信号
t = np.linspace(0, 1, 100)
signal = np.sin(2 * np.pi * t)
# 添加噪声
noise_level = 0.5
noisy_signal = signal + noise_level * np.random.randn(len(t))
# 创建延迟版本作为额外特征
X_signal = np.column_stack([
noisy_signal,
np.roll(noisy_signal, 1), # 延迟1
np.roll(noisy_signal, 2), # 延迟2
])
pca_signal = PCA()
X_denoised = pca_signal.fit_transform(X_signal)
print("信号数据的解释方差比例:", pca_signal.explained_variance_ratio_)
# 可视化去噪效果
plt.figure(figsize=(12, 4))
plt.subplot(1, 3, 1)
plt.plot(t, signal, 'b-', label='真实信号')
plt.plot(t, noisy_signal, 'r-', alpha=0.7, label='噪声信号')
plt.title('原始数据')
plt.legend()
plt.subplot(1, 3, 2)
# 使用第一主成分重构
X_reconstructed = pca_signal.inverse_transform(X_denoised[:, :1])
plt.plot(t, signal, 'b-', label='真实信号')
plt.plot(t, X_reconstructed[:, 0], 'g-', label='PCA去噪')
plt.title('PCA去噪(仅第一主成分)')
plt.legend()
plt.subplot(1, 3, 3)
plt.bar(range(1, 4), pca_signal.explained_variance_ratio_)
plt.title('各主成分的贡献')
plt.xlabel('主成分')
plt.ylabel('解释方差比例')
plt.tight_layout()
plt.show()
when_to_use_pca()实战技巧和最佳实践
1. 数据预处理
def pca_preprocessing_tips():
"""PCA的数据预处理技巧"""
# 创建不同尺度的数据
X_mixed_scale = np.column_stack([
np.random.randn(100) * 1, # 标准正态分布
np.random.randn(100) * 100, # 大尺度
np.random.randn(100) * 0.01, # 小尺度
])
print("原始数据的标准差:")
print(np.std(X_mixed_scale, axis=0))
# 不标准化的PCA
pca_no_scaling = PCA()
pca_no_scaling.fit(X_mixed_scale)
# 标准化后的PCA
from sklearn.preprocessing import StandardScaler
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X_mixed_scale)
pca_scaled = PCA()
pca_scaled.fit(X_scaled)
print("\n不标准化的PCA解释方差比例:")
print(pca_no_scaling.explained_variance_ratio_)
print("\n标准化后的PCA解释方差比例:")
print(pca_scaled.explained_variance_ratio_)
print(f"\n结论: 标准化让方差分布更均匀!")
pca_preprocessing_tips()2. 处理缺失值
def pca_with_missing_values():
"""处理含有缺失值的PCA"""
# 创建含缺失值的数据
X_complete = np.random.randn(100, 5)
X_missing = X_complete.copy()
# 随机添加缺失值
missing_indices = np.random.choice(X_missing.size, size=int(0.1 * X_missing.size), replace=False)
X_missing.flat[missing_indices] = np.nan
print(f"缺失值比例: {np.isnan(X_missing).sum() / X_missing.size:.1%}")
# 方法1: 删除含缺失值的样本
X_dropna = X_missing[~np.isnan(X_missing).any(axis=1)]
print(f"删除法保留样本数: {len(X_dropna)}/{len(X_missing)}")
# 方法2: 均值填充
from sklearn.impute import SimpleImputer
imputer = SimpleImputer(strategy='mean')
X_imputed = imputer.fit_transform(X_missing)
# 比较PCA结果
pca_complete = PCA().fit(X_complete)
pca_imputed = PCA().fit(X_imputed)
print("\n完整数据PCA解释方差:", pca_complete.explained_variance_ratio_[:3])
print("填充数据PCA解释方差:", pca_imputed.explained_variance_ratio_[:3])
pca_with_missing_values()总结：PCA的核心思想
PCA就像是一个智能的数据摄影师：
- 🎯 找最佳角度：寻找数据变化最大的方向
- 📐 保持垂直：确保各个角度互不干扰（正交）
- 🎭 突出重点：按重要性排序主成分
- ✂️ 精简表达：用最少的维度表达最多的信息
- 🔄 可逆变换：能够从低维恢复到高维（有损）
记忆口诀
- 找方向：找到数据的主要变化方向
- 排顺序：按方差大小排列主成分
- 降维度：选择前几个主要成分
- 保信息：尽可能保留原始信息
PCA不仅是一个强大的降维工具，更是理解数据结构的窗口。它告诉我们：复杂的高维数据往往蕴含着简单的低维结构！
作者: meimeitou