主成分分析

主成分分析（Principal Component Analysis, PCA）是数学建模中最常用的无监督降维方法：把 $p$ 个彼此相关的特征通过正交线性变换压缩成 $k$ （ $k \le p$ ）个互不相关的新变量（主成分），在尽量保留原始信息（方差）的前提下实现降维、去相关、可视化与去噪。它在竞赛中常出现在"高维指标降维后再建模""多指标综合评分""高维数据可视化"等环节，也是处理多重共线性（主成分回归）的经典手段。本文从原理、适用场景、评价指标、可视化诊断到可运行代码，完整梳理 PCA 的竞赛实战用法。

一、算法含义

1.1 通俗理解

假如题目给了你 50 个彼此高度相关的指标（比如 50 道题的考试成绩、50 个城市经济指标），直接建模会遇到三个麻烦：

维数灾难：50 维数据无法画图、难以理解，模型计算量大且容易过拟合；
多重共线性：指标两两相关，回归系数估计不稳定（方差膨胀，见 02_多元线性回归）；
噪声冗余：每个指标都含测量噪声，信息被稀释在大量相关变量中。

PCA 的想法非常朴素：找一个新的坐标系，让数据在新坐标轴上的投影尽量"散开"。"散开"的程度用方差衡量——投影方向上的方差越大，说明这个方向携带的原始信息越多。于是：

第一主成分（PC1）是方差最大的投影方向；
第二主成分（PC2）是在与 PC1 垂直（从而不相关）的所有方向中方差最大者；
依此类推，直到第 $p$ 个主成分， $p$ 个新变量重新分配了全部方差；
每个主成分都是原特征的线性组合： $\text{PC}_j = v_{1j}x_1 + v_{2j}x_2 + \cdots + v_{pj}x_p$ ，权重 $v_{ij}$ 称为载荷（loading）；
主成分之间两两不相关（新坐标轴两两正交）；
方差被"集中"到前面几个主成分上，后面的大多是噪声方向，可以直接丢弃 → 实现降维。

举个直观例子：200 名学生考了 6 门课（数学、物理、化学、语文、英语、历史）。6 门成绩彼此相关，PCA 会发现"背后其实只有 2 个潜在因子"——一个"综合能力因子"（解释约 56% 方差）、一个"文理倾向因子"（解释约 35%），二者合计约 91%。于是用 2 维替代 6 维，只损失约 9% 的信息。这正是第六节合成数据的构造思路，也是 7.1 节实测输出的结论。

1.2 数学模型：方差最大的投影方向

设数据矩阵 $X \in \mathbb{R}^{n \times p}$ （ $n$ 个样本、 $p$ 个特征），第一步中心化（把坐标原点平移到数据中心）：

$X_c = X - \bar{x}, \qquad \bar{x} = \frac{1}{n}\sum_{i=1}^{n} x_i$

样本协方差矩阵（对称半正定）：

$S = \frac{1}{n-1} X_c^\top X_c \in \mathbb{R}^{p \times p}, \qquad S_{ij} = \frac{1}{n-1}\sum_{t=1}^{n}(x_{ti} - \bar{x}_i)(x_{tj} - \bar{x}_j)$

若先对每个特征做 z-score 标准化，则 $S$ 退化为相关矩阵（对角线全为 1，见 2.2 节"必须先标准化"）。

第一主成分的推导：寻找单位向量 $w$ （ $\|w\| = 1$ ），使投影得分 $z = X_c w$ 的方差最大：

$\max_{w} \; \mathrm{Var}(X_c w) = \max_{w} \; w^\top S w, \qquad \text{s.t.} \; w^\top w = 1$

用拉格朗日乘子法：构造 $L(w, \lambda) = w^\top S w - \lambda (w^\top w - 1)$ ，令 $\dfrac{\partial L}{\partial w} = 2Sw - 2\lambda w = 0$ ，得到

$S w = \lambda w$

即 $w$ 必须是 $S$ 的特征向量，而此时目标函数值 $w^\top S w = \lambda$ 恰为对应的特征值。因此：

取最大特征值 $\lambda_1$ 对应的特征向量 $v_1$ 作为 PC1 方向，且 $\mathrm{Var}(z_1) = \lambda_1$ ；
PC2 在与 $v_1$ 正交的方向中找方差最大者。由于 $S$ 是对称矩阵，不同特征值对应的特征向量自动正交，故 PC2 就是第二大特征值 $\lambda_2$ 对应的特征向量 $v_2$ ， $\mathrm{Var}(z_2) = \lambda_2$ ；
依此类推。特征分解把"总方差最大化"问题彻底解决： $S$ 的特征值与特征向量按大小排序就是全部主成分的方向与方差。

1.3 算法步骤（五步）

中心化： $X_c = X - \bar{x}$ ，消除均值影响；
协方差矩阵： $S = \dfrac{1}{n-1} X_c^\top X_c$ （量纲不一致时必须先标准化，此时 $S$ 就是相关矩阵 $R$ ）；
特征分解：解 $S v = \lambda v$ ，得特征值 $\lambda_1 \ge \lambda_2 \ge \dots \ge \lambda_p \ge 0$ 与对应的单位正交特征向量 $v_1, v_2, \dots, v_p$ ；
取前 k 个特征向量：组成投影矩阵 $V_k = [v_1, v_2, \dots, v_k] \in \mathbb{R}^{p \times k}$ ；
投影： $Z = X_c V_k \in \mathbb{R}^{n \times k}$ ， $Z$ 的第 $j$ 列就是第 $j$ 个主成分的得分（score）。

主成分 = 线性组合：第 $j$ 个主成分是原特征（中心化后）的加权和

$z_j = X_c v_j = v_{1j} x_{c1} + v_{2j} x_{c2} + \cdots + v_{pj} x_{cp}$

权重 $v_{ij}$ 的正负与大小刻画了每个原特征对该主成分的贡献方向与强度。主成分的两条核心性质：

方差集中： $\mathrm{Var}(z_j) = \lambda_j$ ，且 $\sum_{j=1}^{p} \lambda_j = \mathrm{tr}(S)$ ，即总方差不变、只是重新分配；
两两不相关： $\mathrm{Cov}(z_i, z_j) = v_i^\top S v_j = \lambda_j v_i^\top v_j = 0 \; (i \ne j)$ ，去相关是 PCA 处理共线性的根基。

方差贡献率与累计贡献率：第 $j$ 个主成分的方差贡献率为 $\gamma_j = \dfrac{\lambda_j}{\sum_{l=1}^{p} \lambda_l}$ ，前 $k$ 个的累计方差贡献率为 $\Gamma_k = \dfrac{\sum_{j=1}^{k} \lambda_j}{\sum_{l=1}^{p} \lambda_l}$ 。 $\Gamma_k$ 就是"降到 $k$ 维保留了百分之多少的信息"，是最常用的选 $k$ 依据（见 3.3 节）。

1.4 几何直觉

数据点云在 $p$ 维空间中大致呈一个"椭球"形状。主成分就是这个椭球的主轴：PC1 是椭球最长的轴（数据散得最开的方向），PC2 是与 PC1 垂直的最长轴，……，各主轴两两垂直；特征值 $\lambda_j$ 就是椭球沿第 $j$ 轴方向的"半径平方"。
两个高度正相关的特征（点云呈斜放的"雪茄"形）是最经典的示意：PC1 沿雪茄长轴（几乎包含全部信息），PC2 沿短轴（只剩噪声）。把坐标系旋转到长轴方向后，第一维"几乎等于全部数据"，第二维可以直接扔掉——这就是降维。
与回归的区别（竞赛答辩常被问）：一元线性回归最小化的是点到直线的竖直距离（残差），且区分自变量/因变量；PCA 最小化的是点到直线的垂直距离（正交回归/总体最小二乘），没有自变量因变量之分，完全由方差驱动、对称地对待所有特征。
重构视角： $k$ 维主成分重构 $\hat{X} = Z_k V_k^\top + \bar{x}$ 是所有秩 $k$ 矩阵中对 $X$ 的最佳近似（Eckart–Young 定理），所以 PCA 也是低秩近似/矩阵压缩工具。

1.5 与特征选择（feature selection）的区别

PCA 是"特征提取"（构造新变量），特征选择是"挑选原变量"，两者常被混淆：

对比维度	PCA（特征提取）	特征选择（如方差过滤、LASSO、互信息）
输出	原特征的线性组合（新变量）	原变量的一个子集
是否用标签 $y$	不用（无监督）	可用（过滤式/包裹式/嵌入式常监督）
可解释性	新变量含义需靠载荷推断	保留原变量，含义不变
信息保留	前 $k$ 个主成分保留了大部分方差	保留对任务最有用的变量，未必是方差大的
共线性处理	天然去相关	通常不能消除入选变量间的相关性

竞赛中的选择：目标是"降维后建模/可视化/去相关" → PCA；目标是"找出哪几个原始指标最关键、结论要落到原始指标上" → 特征选择。两者也可以串联（先 PCA 降维，再在低维空间建模）。

1.6 优缺点

优点：

无监督、无需标签，任何数据都能用，原理简单、结论直观；
去相关：主成分两两不相关，直接解决多重共线性（主成分回归 PCR 的基础）；
可压缩：只保留前 $k$ 个主成分就能重构出原数据的绝大部分信息，实现降维、压缩与提速；
去噪：小特征值方向通常对应噪声，丢弃后相当于给数据"滤波"，下游模型（聚类、回归、神经网络）反而更稳；
有闭式解（特征分解），实现简单（numpy 十几行即可，见第六节），无迭代、无局部最优、结果唯一（方向符号除外）；
可视化利器：前两个主成分的散点图是竞赛中最常用的高维数据展示手段。

缺点：

可解释性下降：主成分是"全部特征的线性组合"，业务含义要靠载荷矩阵推断甚至命名，有时解释不清；
线性假设：只捕捉线性相关（Pearson 相关），曲线/流形结构无效；
对量纲极其敏感：不标准化时方差大的特征垄断主成分（2.2 节）；
对离群点敏感：协方差矩阵会被极端值扭曲（平方项放大影响）；
计算 $O(p^3)$ ： $p$ 极大时需用随机化 SVD、增量 PCA 等近似算法；
小主成分可能恰恰携带对 $y$ 最有预测力的信息（PCA 只保方差、不保预测力），PCR 回归可能反而不如直接回归（见 2.1 与 8.5 节）。

二、何时使用（适用场景与条件）

2.1 适用场景

高维数据降维与可视化：特征数十上百，无法画图、难以理解时，先 PCA 投影到前 2~3 个主成分画散点图，观察样本的聚类结构、离群点与分布形态。竞赛中的经典套路是"先 PCA 降维，再 K-means/层次聚类"，既缓解维数灾难又去噪；
多重共线性处理：自变量高度相关导致回归系数不稳定时有两种用法——（a）主成分回归（PCR）：用前 $k$ 个主成分得分 $Z_k$ 代替原自变量做回归，共线性被"去相关"彻底消除，代价是可解释性下降，且被丢弃的小主成分可能对 $y$ 有预测力（此时可用 06_岭回归或偏最小二乘 PLS）；（b）仅作诊断：特征值中出现 $\lambda_j \approx 0$ 说明特征间存在近似线性关系，可据此删减指标；
数据压缩与去噪：保留前 $k$ 个主成分再重构 $\hat{X} = Z_k V_k^\top + \bar{x}$ ，即低秩近似，用于图像压缩、信号去噪、存储与传输提速；竞赛中"先 PCA 去噪再建模"经常能小幅提升模型稳健性；
竞赛题目中的指标降维再建模：评价类题目中把几十个相关指标压缩成少数主成分再做综合评价（第一主成分常被当作综合得分，或按贡献率加权求和）；预测类题目中先降维再进回归/分类/神经网络，降低过拟合风险；
探索性分析：用载荷矩阵发现指标分组与潜在结构（哪些指标"同涨同跌"、哪些"此消彼长"），如第六节的例子中 PCA 从 6 个可观测特征恢复了 2 个潜在因子。

2.2 使用前提（建模前检查清单）

必须先标准化（最高优先级）：PCA 基于协方差矩阵，对量纲极度敏感。若特征单位不同（元 vs 人 vs %），方差大的特征会垄断主成分。标准做法是先做 z-score 标准化

$x_{ij}' = \frac{x_{ij} - \bar{x}_j}{s_j}, \qquad s_j = \sqrt{\frac{1}{n}\sum_{i=1}^{n}(x_{ij} - \bar{x}_j)^2}$

标准化后每个特征方差为 1、协方差矩阵等于相关矩阵，主成分反映的是相关结构而非量纲。唯一例外：全部特征同量纲且希望保留"绝对规模"的含义（如均为"亿元"时可能有意让大规模指标权重更大）——竞赛中几乎总是标准化，7.2 节给出了不标准化的实测反例；

线性结构假设：PCA 只刻画线性相关。若特征间是曲线关系（如 $x_2 = x_1^2$ ）或数据呈流形结构（瑞士卷），PCA 效果差，应换核 PCA、t-SNE 或 UMAP（见 2.4 与第八节）；
连续数值型特征：分类变量需先编码（独热编码后 PCA 可运行但解释性差；以分类变量为主的数据更适合对应分析 MCA）；
样本量充分：经验上建议 $n \ge 5p$ （更严格 $n \ge 10p$ ）。 $n < p$ 时协方差矩阵估计不稳定，需改用稀疏 PCA、正则化方法；
无明显极端离群点：一个极端值就能把主成分方向拽歪，建议先用箱线图或 $3\sigma$ 准则筛查。

2.3 不适用 / 慎用的情形

非线性流形结构：瑞士卷、环形、嵌套球面等数据，线性投影必然把结构"压扁" → 核 PCA、t-SNE、UMAP；
必须保留原变量含义：结论必须落到"哪个原始指标最重要"（如政策建议要指名道姓到指标）→ 特征选择、LASSO；若需要"可命名的公共因子"→ 因子分析（可旋转，见 2.4）；
主成分解释不清的场景：载荷矩阵杂乱、每个主成分都"六不像"，且评审要求逐变量解释 → 慎用 PCA 作为最终结论，可只用它做预处理/可视化；
要求稀疏性：希望每个新变量只用到极少数原特征（便于理解、便于落地）→ 稀疏 PCA。

2.4 与因子分析、t-SNE 的对比与选择

对比维度	PCA	因子分析（FA）	t-SNE / UMAP
目标	最大化方差，降维压缩	解释变量间的协方差结构（潜变量模型）	保持样本间局部邻近关系，可视化
数学本质	特征分解 $S v = \lambda v$ ，无概率模型	$x = \Lambda f + \varepsilon$ ，噪声单独建模、迭代估计	非线性优化（t-SNE 随机，UMAP 近似确定）
输出	载荷 + 得分，有显式线性变换	因子载荷（可 varimax 旋转）+ 因子得分	仅低维坐标，无可复用变换
假设	线性相关	潜变量线性结构 + 特殊因子独立	几乎无分布假设
可解释性	载荷可解读，无旋转惯例	旋转后更易命名，可解释性更强	低（坐标本身无业务含义）
后续建模	可作新特征进任何模型	因子得分可进模型	不能（仅可视化）
竞赛典型用途	降维建模、去相关、综合评价	探索性结构分析、"公共因子"命名	高维数据可视化、发现簇结构

选择口诀：降维后还要接着建模 → PCA；要给"公共因子"起名字、做结构分析 → 因子分析；只想画张漂亮的二维图看簇/流形 → t-SNE/UMAP（注意其结果不可复用到新样本）。

三、算法指标

3.1 特征值 λ_j

中文名：特征值（第 $j$ 主成分的方差）
公式： $S v_j = \lambda_j v_j$ ，且 $\lambda_j = \mathrm{Var}(z_j)$
含义：第 $j$ 个主成分方向上数据的分散程度。 $\lambda_j$ 越大，该方向携带的信息越多。
解读：① 拐点法——画出 $\lambda_j$ 随 $j$ 的变化（碎石图），在"骤降之后趋于平缓"的位置截断；② Kaiser 准则——只保留 $\lambda_j > 1$ 的主成分（适用于标准化数据： $\lambda_j = 1$ 表示该主成分的信息量恰好等于一个标准化特征的平均信息量， $\lambda_j < 1$ 表示连一个特征都比不上，通常只含噪声）；③ $\lambda_j \approx 0$ 说明特征间存在近似精确的线性关系，可作共线性诊断。

3.2 方差贡献率

中文名：方差贡献率（单个主成分解释的信息比例）
公式： $\gamma_j = \dfrac{\lambda_j}{\sum_{l=1}^{p} \lambda_l} \times 100\%$
含义：第 $j$ 个主成分单独解释了总方差的百分之多少。
解读： $\gamma_j$ 越大说明该主成分越重要。若 $\gamma_1$ 很高（如 > 60%），说明数据高度共线，几乎可以用一个综合指标概括；若前几个 $\gamma_j$ 都很小且差不多，说明数据近似"球形"，各方向方差均匀，没有降维空间（PCA 帮不上忙）。

3.3 累计方差贡献率

中文名：累计方差贡献率（降维后的信息保留比例）
公式： $\Gamma_k = \dfrac{\sum_{j=1}^{k} \lambda_j}{\sum_{l=1}^{p} \lambda_l} \times 100\%$
含义：前 $k$ 个主成分共保留了原数据的百分之多少的方差（信息）。
解读：这是竞赛中最常用的选 $k$ 标准。经验阈值通常取 85%（也有取 70%、75%、90% 的，视题目对精度的要求与"降维力度"的权衡而定），即选最小的 $k$ 使 $\Gamma_k \ge 85\%$ 。注意：85% 是经验惯例而不是定理（见 7.2 节"阈值机械套用"坑）；选 $k$ 时应把拐点法、Kaiser 准则、累计贡献率曲线三者互相印证，并在论文中说明。

3.4 重构误差 MSE(k)

中文名：重构均方误差（信息损失的直接度量）
公式：先用 $k$ 个主成分重构 $\hat{X}(k) = Z_k V_k^\top + \bar{x}$ ，再计算

$\mathrm{MSE}(k) = \frac{1}{np}\sum_{i=1}^{n}\sum_{j=1}^{p}\left(x_{ij} - \hat{x}_{ij}(k)\right)^2$

含义：只用 $k$ 维信息重构全部特征时的平均误差（标准化数据下单位是"标准化方差"）。
解读：① $\mathrm{MSE}(k)$ 随 $k$ 增加单调下降， $k = p$ 时为 0（完全重构）；② 标准化数据下 $\mathrm{MSE}(0) = 1$ （每个特征方差为 1，只用均值重构），且近似有 $\mathrm{MSE}(k) \approx 1 - \Gamma_k$ ——重构误差与累计贡献率是一枚硬币的两面，可互相印证；③ $\mathrm{MSE}(k)$ 的"骤降→平缓"拐点与碎石图拐点一致，可交叉验证 $k$ 的选择。

3.5 载荷矩阵

中文名：载荷矩阵（主成分与原始特征的权重关系）
公式： $L = [v_1, v_2, \dots, v_k] \in \mathbb{R}^{p \times k}$ ，元素 $l_{ij} = v_{ij}$ （特征 $i$ 在主成分 $j$ 上的载荷，即组合系数）。有的教材把载荷定义为 $v_{ij}\sqrt{\lambda_j}$ ——标准化数据下它恰好等于特征 $i$ 与主成分 $j$ 的相关系数 $r(x_i, z_j)$ ，两者只差一个常数倍，解读方式相同；
含义： $|l_{ij}|$ 大 → 特征 $i$ 对主成分 $j$ 贡献大；符号 → 贡献方向。
解读：同一主成分中同号的特征"同涨同跌"、异号的特征"此消彼长"。据此给主成分命名（如"综合规模因子""结构差异因子"），这是把 PCA 结果翻译成业务语言的唯一途径，也是论文出彩点。注意载荷整体变号是允许的（ $v_j$ 与 $-v_j$ 都是合法解），解读时以相对大小与符号关系为准。

3.6 指标汇总表

指标	中文名	公式	含义与解读
$\lambda_j$	特征值	$S v_j = \lambda_j v_j$	第 $j$ 主成分的方差；拐点法、Kaiser 准则（ $\lambda_j>1$ ）选 $k$
$\gamma_j$	方差贡献率	$\lambda_j / \sum \lambda_l$	单个主成分解释的信息比例；判断数据有无降维空间
$\Gamma_k$	累计方差贡献率	$\sum_{j \le k} \lambda_j / \sum \lambda_l$	降维后信息保留比例；经验阈值 85% 选 $k$
$\mathrm{MSE}(k)$	重构均方误差	$\frac{1}{np}\sum (x_{ij}-\hat{x}_{ij})^2$	信息损失； $\approx 1-\Gamma_k$ （标准化数据），与拐点互证
$L = [v_{ij}]$	载荷矩阵	特征向量按列排布	特征对主成分的贡献与方向；用于主成分命名与业务解读

四、可视化图表

4.1 四张图速查表

图名（文件）	用途	关键解读点
① 碎石图 `pca_scree.png`	拐点法选主成分个数	特征值从大到小排列；在"骤降变平缓"处截断；标注 Kaiser 线 $\lambda=1$ ，拐点与 $\lambda_j>1$ 的个数一致时结论最硬
② 累计方差贡献率曲线 `pca_cumulative.png`	阈值法选主成分个数	曲线与 85%（红线）/90%（绿线）阈值线的交点对应 $k$ ；曲线越陡说明降维越有效
③ 前两个主成分散点图 `pca_scores.png`	高维数据可视化	样本在 PC1-PC2 平面的分布；有无簇结构、离群点、梯度方向；坐标轴标注各主成分贡献率
④ 双标图 `pca_biplot.png`	样本 + 载荷联合解读	红色箭头 = 特征载荷：箭头越长贡献越大，方向相近的特征正相关、相反的特征负相关；箭头方向=该特征在样本散点中的变化方向

4.2 每张图的详细解读要点

图① 碎石图（Scree Plot）：横轴主成分序号、纵轴特征值（柱状），右轴叠加累计贡献率曲线。判读三步：先找特征值骤降处（"悬崖"脚下），本例 $\lambda$ 从 2.1004 骤降到 0.2466，拐点即 $k=2$ ；再看 Kaiser 线 $\lambda=1$ ——位于线上方的柱（2 个）就是要保留的个数；最后与右轴累计曲线对照（ $k=2$ 时约 90.9%）。三种信息画在一张图上，选 $k$ 的论证一目了然。

图② 累计方差贡献率曲线：横轴 $k$ 、纵轴 $\Gamma_k$ （%）。曲线与 85%、90% 两条阈值线的交点就是候选 $k$ 。判读要点：① 曲线前段越陡，说明方差越集中在少数主成分上，降维越划算；② 若曲线缓慢爬升（如 $k=4$ 才到 70%），说明特征之间相关弱，PCA 降维收益小，应重新考虑是否用 PCA；③ 阈值线只在论文中作为论证依据之一，不能替代业务合理性（7.2 节坑 3）。

图③ 前两个主成分散点图：每个点一个样本，横纵坐标是 PC1、PC2 得分。判读要点：① 是否自然分成若干簇（为后续聚类提供依据）；② 有无远离主体的离群点；③ 坐标轴标注贡献率，说明"这张图浓缩了多少信息"；④ 若样本有已知标签（类别、时间、数值型协变量），可用颜色/形状标注观察主成分与它的关系——本例按潜在因子 f1 着色，颜色沿 PC1 方向渐变，直观证明 PC1 捕获了 f1。

图④ 双标图（Biplot）：在同一坐标系里画样本得分（灰点）与特征载荷（红色箭头，长度按显示比例缩放）。判读要点：① 箭头长度代表该特征对这两个主成分的贡献大小，短箭头（如 X5 在 PC1 方向几乎无投影）说明该特征与这两个主成分关系弱；② 两箭头夹角小（方向相近）→ 两特征正相关，夹角接近 180° → 负相关，接近 90° → 不相关；③ 样本沿某箭头方向延伸 → 该方向上样本主要被这个特征区分；④ 载荷整体变号不影响结论（箭头全部翻转后解释不变）。

五、符号说明

符号	含义	备注
$X$	数据矩阵， $n \times p$	$n$ 个样本、 $p$ 个特征
$\bar{x}$	各特征的样本均值向量	长度 $p$
$X_c$	中心化后的数据	$X_c = X - \bar{x}$
$S$	样本协方差矩阵， $p \times p$	$S = X_c^\top X_c / (n-1)$ ；标准化后即相关矩阵 $R$
$\lambda_j$	第 $j$ 大特征值	$= \mathrm{Var}(z_j)$ ，第 $j$ 主成分的方差
$v_j$	第 $j$ 个特征向量（载荷向量）	单位向量、两两正交；方向符号可任意（ $-v_j$ 等价）
$V_k$	前 $k$ 个特征向量组成的投影矩阵	$V_k = [v_1, \dots, v_k] \in \mathbb{R}^{p \times k}$
$Z$	主成分得分矩阵， $n \times k$	$Z = X_c V_k$ ；第 $j$ 列即 PC $j$ 的得分
$z_{ij}$	第 $i$ 个样本在第 $j$ 主成分上的得分	常作为降维后的新特征使用
$k$	保留的主成分个数	由拐点法 / Kaiser 准则 / 累计贡献率阈值共同确定
$L$	载荷矩阵， $p \times k$	$L = V_k$ ；有的教材取 $v_{ij}\sqrt{\lambda_j}$ （= 特征与主成分的相关系数）
$\gamma_j$	第 $j$ 主成分的方差贡献率	$\gamma_j = \lambda_j / \sum_l \lambda_l$
$\Gamma_k$	前 $k$ 个主成分的累计方差贡献率	$\Gamma_k = \sum_{j \le k} \lambda_j / \sum_l \lambda_l$
$\hat{X}(k)$	用前 $k$ 个主成分重构的数据	$\hat{X}(k) = Z_k V_k^\top + \bar{x}$
$\mathrm{MSE}(k)$	重构均方误差	$\frac{1}{np}\sum_{i,j}(x_{ij}-\hat{x}_{ij})^2$
$f_1, f_2$	潜在因子（真实结构）	仅合成数据中存在，用于验证 PCA 是否恢复了真实结构

六、可运行程序（完整代码）

环境要求：Python 3.12，依赖 numpy、pandas、scikit-learn、matplotlib（pip install numpy pandas scikit-learn matplotlib）。以下所有代码块按顺序拼接保存为 pca_demo.py，在本文档所在目录运行即可：控制台打印第三节全部指标（特征值、方差贡献率、累计贡献率、选 $k$ 结论、载荷矩阵、重构误差、sklearn 对照、因子恢复验证、不标准化的反例），并在 figures/ 子目录生成 4 张图（文件名以 pca_ 开头）。数据为程序内合成数据（ $n = 200$ ，6 个特征：前 3 个由 2 个潜在因子强线性组合 + 小噪声构成，后 3 个由相同 2 个因子弱线性组合 + 较大噪声构成，np.random.seed(42) 保证可复现），无外部文件依赖。

# -*- coding: utf-8 -*-
"""
============================================================
主成分分析（PCA）完整示例：手写实现 + sklearn 对照
------------------------------------------------------------
数据（合成）：n = 200，6 个特征。前 3 个特征由 2 个潜在因子
（f1、f2）强线性组合 + 小噪声构成；后 3 个特征由相同 2 个
因子弱线性组合 + 较大噪声构成（模拟真实数据中的弱相关指标）。
np.random.seed(42) 保证可复现，无外部文件依赖。
流程：数据生成 → z-score 标准化 → 手写 PCA（中心化、协方差矩阵、
      np.linalg.eigh 特征分解、按特征值降序取前 k、投影、方差
      贡献率与重构误差）→ sklearn PCA 对照（主成分方向符号可能
      相反，方差贡献率必须一致）→ 打印全部指标
      → figures/ 下 4 张图（pca_ 前缀）
依赖：numpy、pandas、scikit-learn、matplotlib（Python 3.12）
============================================================
"""
import os
import numpy as np
import pandas as pd
import matplotlib
import matplotlib.pyplot as plt
from sklearn.preprocessing import StandardScaler
from sklearn.decomposition import PCA as SklearnPCA

# ------------------------- 全局设置 -------------------------
np.random.seed(42)
plt.rcParams["font.sans-serif"] = ["PingFang SC", "Hiragino Sans GB",
                                   "Arial Unicode MS", "Heiti TC"]
plt.rcParams["axes.unicode_minus"] = False   # 解决负号显示为方块的问题

FIG_DIR = "figures"
os.makedirs(FIG_DIR, exist_ok=True)          # 图片保存目录（不存在则创建）

# --------------------- 1. 合成数据 ---------------------
n, p = 200, 6
f1 = np.random.normal(0, 1, n)               # 潜在因子 1（不可观测）
f2 = np.random.normal(0, 1, n)               # 潜在因子 2（不可观测）
X1 = 3.0 * f1 + 0.3 * np.random.normal(0, 1, n)            # 强载荷于 f1
X2 = 2.5 * f1 - 1.0 * f2 + 0.3 * np.random.normal(0, 1, n)  # 强载荷于 f1、f2
X3 = -1.5 * f1 + 2.0 * f2 + 0.3 * np.random.normal(0, 1, n)  # 强载荷于 f1、f2
X4 = 0.9 * f1 + 0.5 * np.random.normal(0, 1, n)             # 弱载荷于 f1
X5 = -0.8 * f2 + 0.5 * np.random.normal(0, 1, n)            # 弱载荷于 f2
X6 = 0.7 * f1 + 0.7 * f2 + 0.5 * np.random.normal(0, 1, n)  # 弱载荷于 f1、f2
X_raw = np.column_stack([X1, X2, X3, X4, X5, X6])
feat_names = ["X1", "X2", "X3", "X4", "X5", "X6"]

# 标准化：PCA 前必须做（协方差矩阵对量纲敏感），z = (x - mean) / std
scaler = StandardScaler()
X_std = scaler.fit_transform(X_raw)

print("=" * 72)
print("0. 数据概览（标准化前）")
print("=" * 72)
print(f"样本数 n = {n}，特征数 p = {p}")
print(f"各特征原始标准差: {np.round(X_raw.std(axis=0, ddof=1), 4)}")
print(f"标准化后各特征均值/标准差: 均值≈{np.round(X_std.mean(axis=0), 6)}, "
      f"标准差={np.round(X_std.std(axis=0), 6)}")
print("相关系数矩阵（注意 X1~X3 之间、以及它们与 X4~X6 之间的相关结构）:")
print(pd.DataFrame(np.corrcoef(X_raw.T), index=feat_names,
                   columns=feat_names).round(3))

以上代码块 1 完成数据生成与标准化：两个标准正态的潜在因子 f1、f2 驱动 6 个可观测特征（X1~~X3 强载荷，X4~~X6 弱载荷 + 较大噪声），随后用 StandardScaler 做 z-score 标准化（等价于手写 $x_{ij}' = (x_{ij}-\bar{x}_j)/s_j$ ），并打印原始标准差与相关系数矩阵。

# --------------------- 2. 手写 PCA ---------------------
def pca_manual(X, k=None):
    """
    手写 PCA（五步）：
      ① 中心化：Xc = X - mean(X)
      ② 协方差矩阵：S = Xc^T Xc / (n-1)
      ③ 特征分解：S v = λ v（np.linalg.eigh 按升序返回，需反转为降序）
      ④ 取前 k 个特征值对应的特征向量组成 Vk
      ⑤ 投影：Z = Xc @ Vk
    返回：特征值(降序)、特征向量(按列，降序)、主成分得分 Z、均值向量
    """
    n = X.shape[0]
    mean = X.mean(axis=0)
    Xc = X - mean                              # ① 中心化
    S = (Xc.T @ Xc) / (n - 1)                  # ② 协方差矩阵
    lam, V = np.linalg.eigh(S)                 # ③ 特征分解（升序）
    order = np.argsort(lam)[::-1]              # 降序排列
    lam, V = lam[order], V[:, order]
    if k is None:
        k = X.shape[1]
    Z = Xc @ V[:, :k]                          # ⑤ 投影得到主成分得分
    return lam, V, Z, mean

lam, V, Z, mean = pca_manual(X_std)            # 对标准化数据做 PCA
explained = lam / lam.sum()                    # 方差贡献率
cum = np.cumsum(explained)                     # 累计方差贡献率
k85 = int(np.searchsorted(cum, 0.85) + 1)      # 累计贡献率首次 ≥ 85% 的 k
k90 = int(np.searchsorted(cum, 0.90) + 1)      # 累计贡献率首次 ≥ 90% 的 k
elbow = int(np.argmax(lam[:-1] - lam[1:]) + 1)  # 拐点：特征值骤降最大的位置

print("=" * 72)
print("1. 手写 PCA（基于标准化数据，协方差矩阵 = 相关矩阵）")
print("=" * 72)
print(f"特征值（降序）      : {np.round(lam, 4)}")
print(f"方差贡献率          : {np.round(explained, 4)}")
print(f"累计方差贡献率      : {np.round(cum, 4)}")
print(f"特征值最大骤降位置（拐点）: k = {elbow}"
      f"（λ 从 {lam[elbow-1]:.4f} 骤降至 {lam[elbow]:.4f}）")
print(f"Kaiser 准则（λ>1）建议保留的主成分个数: {int((lam > 1).sum())}")
print(f"累计贡献率首次 ≥85% 的 k = {k85}；首次 ≥90% 的 k = {k90}")
print(f"三种选 k 方法结论一致：拐点法 k={elbow}，Kaiser 准则 "
      f"k={int((lam > 1).sum())}，85% 阈值 k={k85}，90% 阈值 k={k90}")
print("\n载荷矩阵（前 3 个主成分，行=特征，列=主成分）:")
print(pd.DataFrame(V[:, :3], index=feat_names,
                   columns=["PC1", "PC2", "PC3"]).round(4))

以上代码块 2 是手写 PCA 核心：完全对应 1.3 节的五步（中心化 → 协方差 → np.linalg.eigh 特征分解 → 降序取前 $k$ 个特征向量 → 投影）。注意 eigh 返回升序特征值，需 argsort 反转；用 searchsorted 求累计贡献率首次越过 85%/90% 阈值的 $k$ ，用相邻特征值差最大处定位拐点。

# --------------------- 3. 重构误差 ---------------------
print("=" * 72)
print("2. 重构误差（用前 k 个主成分重构标准化数据）")
print("=" * 72)
print(f"{'k':>3} {'MSE(k)':>10} {'累计贡献率':>12}")
for k in range(p + 1):
    if k == 0:
        Xhat = np.zeros_like(X_std)            # k=0：只用均值（=0）重构
    else:
        Xhat = Z[:, :k] @ V[:, :k].T + mean    # 低秩重构
    mse = np.mean((X_std - Xhat) ** 2)
    cum_k = 0.0 if k == 0 else cum[k - 1]
    print(f"{k:>3} {mse:>10.4f} {cum_k:>12.4f}")
print("（标准化数据每个特征的方差为 1，总方差为 6，故 MSE(0)=1）")

# --------------------- 4. sklearn 对照 ---------------------
pca_sk = SklearnPCA()                          # 默认 n_components=p
Z_sk = pca_sk.fit_transform(X_std)

print("=" * 72)
print("3. sklearn PCA 对照")
print("=" * 72)
print(f"sklearn 方差贡献率 : {np.round(pca_sk.explained_variance_ratio_, 4)}")
print(f"手写   方差贡献率  : {np.round(explained, 4)}")
print(f"两者最大绝对差     : {np.abs(pca_sk.explained_variance_ratio_ - explained).max():.3e}")
for j in range(p):
    r = np.corrcoef(Z[:, j], Z_sk[:, j])[0, 1]
    print(f"PC{j+1}: corr(手写得分, sklearn得分) = {r:+.4f}  （±1 即完全一致，"
          f"符号为负说明该主成分方向相反，属正常现象）")

# 因子恢复验证：前两个主成分应与两个潜在因子高度相关
r_f1 = np.corrcoef(Z[:, 0], f1)[0, 1]
r_f2 = np.corrcoef(Z[:, 1], f2)[0, 1]
print("=" * 72)
print("4. 潜在因子恢复验证")
print("=" * 72)
print(f"|corr(PC1, f1)| = {abs(r_f1):.4f}，|corr(PC2, f2)| = {abs(r_f2):.4f}")
print("（接近 1 说明 PCA 从 6 个可观测特征中成功恢复了 2 个潜在因子）")

# --------------------- 5. 不标准化会怎样（反面教材） ---------------------
lam_raw, V_raw, Z_raw, _ = pca_manual(X_raw)   # 直接对原始数据做 PCA
print("=" * 72)
print("5. 不标准化的后果（对照）")
print("=" * 72)
print(f"原始数据特征值     : {np.round(lam_raw, 4)}")
print(f"原始数据方差贡献率 : {np.round(lam_raw / lam_raw.sum(), 4)}")
print(f"原始数据 PC1 载荷  : {np.round(V_raw[:, 0], 4)}")
print("（方差大的 X1~X3 垄断主成分，弱指标的载荷被压缩，结论被量纲绑架）")

以上代码块 3 完成重构误差、sklearn 对照与反例验证：低秩重构公式 $\hat{X}(k) = Z_k V_k^\top + \bar{x}$ 计算 $\mathrm{MSE}(k)$ ；与 sklearn.decomposition.PCA 对比方差贡献率（应完全一致）与得分（相关应为 ±1，负号 = 主成分方向相反，属正常现象，见 7.2 节）；再把潜在因子 f1、f2 与主成分得分求相关，验证 PCA 是否恢复了真实结构；最后直接对未标准化数据跑一遍 PCA 作为反面教材。

# --------------------- 6. 绘图 ---------------------
# 图1：碎石图（特征值 + 累计贡献率双轴）
fig, ax1 = plt.subplots(figsize=(8, 5))
bars = ax1.bar(range(1, p + 1), lam, color="tab:blue", alpha=0.65,
               label="特征值 λ")
for i, v in enumerate(lam):
    ax1.text(i + 1, v + 0.05, f"{v:.2f}", ha="center", fontsize=9)
ax1.axhline(1.0, ls="--", color="gray", lw=1)
ax1.text(4.35, 0.82, "λ=1（Kaiser 准则线）", va="center", fontsize=9)
ax1.annotate(f"拐点 k={elbow}\n（之后特征值骤降且 < 1）",
             xy=(elbow, lam[elbow - 1]), xytext=(3.3, 2.5),
             arrowprops=dict(arrowstyle="->", color="black"))
ax1.set_xlabel("主成分序号")
ax1.set_ylabel("特征值 λ")
ax1.set_xticks(range(1, p + 1))
ax1.set_title("碎石图（Scree Plot）：拐点法选择主成分个数")
ax2 = ax1.twinx()
ax2.plot(range(1, p + 1), cum * 100, "o-", color="tab:orange", lw=1.8,
         label="累计方差贡献率")
ax2.set_ylabel("累计方差贡献率（%）")
lines1, labels1 = ax1.get_legend_handles_labels()
lines2, labels2 = ax2.get_legend_handles_labels()
ax1.legend(lines1 + lines2, labels1 + labels2, loc="center right")
fig.tight_layout()
fig.savefig(f"{FIG_DIR}/pca_scree.png", dpi=150)

# 图2：累计方差贡献率曲线（标注 85%/90% 阈值线）
fig, ax = plt.subplots(figsize=(8, 5))
ax.plot(range(1, p + 1), cum * 100, "o-", color="tab:blue", lw=2)
for i, v in enumerate(cum):
    ax.text(i + 1, v * 100 + 1.2, f"{v*100:.1f}%", ha="center", fontsize=9)
ax.axhline(85, ls="--", color="tab:red", lw=1.4)
ax.text(1.05, 86.2, "85% 阈值线", color="tab:red", fontsize=10)
ax.axhline(90, ls="--", color="tab:green", lw=1.4)
ax.text(1.05, 91.2, "90% 阈值线", color="tab:green", fontsize=10)
ax.scatter([k85], [cum[k85 - 1] * 100], color="tab:red", zorder=5)
ax.annotate(f"k={k85} 时达 {cum[k85-1]*100:.1f}%",
            xy=(k85, cum[k85 - 1] * 100), xytext=(k85 + 1.0, cum[k85 - 1] * 100 - 7),
            arrowprops=dict(arrowstyle="->", color="tab:red"), fontsize=10)
ax.set_xlabel("主成分个数 k")
ax.set_ylabel("累计方差贡献率（%）")
ax.set_xticks(range(1, p + 1))
ax.set_ylim(0, 105)
ax.set_title("累计方差贡献率曲线：85%/90% 阈值法选 k")
fig.tight_layout()
fig.savefig(f"{FIG_DIR}/pca_cumulative.png", dpi=150)

以上代码块 4 绘制图①碎石图（柱状特征值 + 右轴累计贡献率 + Kaiser 线 $\lambda=1$ + 拐点标注）与图②累计方差贡献率曲线（85%/90% 阈值线与交点标注）。twinx() 双轴技巧让"拐点"与"贡献率"两种信息在同一张图中互相印证。

# 图3：前两个主成分得分散点图（按潜在因子 f1 着色）
fig, ax = plt.subplots(figsize=(8, 6))
sc = ax.scatter(Z[:, 0], Z[:, 1], c=f1, cmap="viridis", s=32,
                alpha=0.85, edgecolors="k", linewidths=0.3)
fig.colorbar(sc, ax=ax, label="潜在因子 f1（真实值，仅用于验证）")
ax.axhline(0, color="gray", lw=0.8)
ax.axvline(0, color="gray", lw=0.8)
ax.set_xlabel(f"PC1（方差贡献率 {explained[0]*100:.1f}%）")
ax.set_ylabel(f"PC2（方差贡献率 {explained[1]*100:.1f}%）")
ax.set_title("前两个主成分得分散点图：颜色沿 PC1 方向渐变\n"
             "（说明 PC1 捕获了潜在因子 f1 的信息，验证见 7.1 节）")
fig.tight_layout()
fig.savefig(f"{FIG_DIR}/pca_scores.png", dpi=150)

# 图4：双标图（biplot）：得分散点 + 载荷箭头
fig, ax = plt.subplots(figsize=(8, 6))
ax.scatter(Z[:, 0], Z[:, 1], s=18, alpha=0.35, color="gray",
           edgecolors="none", label="样本得分")
scale = 0.6 * np.max(np.abs(Z[:, :2])) / np.max(np.abs(V[:, :2]))  # 箭头显示比例
for j, name in enumerate(feat_names):
    ax.annotate("", xy=(V[j, 0] * scale, V[j, 1] * scale), xytext=(0, 0),
                arrowprops=dict(arrowstyle="-|>", color="tab:red", lw=1.6))
    ax.text(V[j, 0] * scale * 1.12, V[j, 1] * scale * 1.12, name,
            color="tab:red", fontsize=10, ha="center", va="center")
ax.axhline(0, color="gray", lw=0.8)
ax.axvline(0, color="gray", lw=0.8)
ax.set_xlabel(f"PC1（{explained[0]*100:.1f}%）")
ax.set_ylabel(f"PC2（{explained[1]*100:.1f}%）")
ax.set_title("双标图（Biplot）：箭头长度=特征对主成分的载荷大小，\n"
             "方向相近的箭头表示特征高度相关")
ax.legend(loc="upper right")
fig.tight_layout()
fig.savefig(f"{FIG_DIR}/pca_biplot.png", dpi=150)

print("=" * 72)
print(f"4 张图已保存至 {FIG_DIR}/ 目录：pca_scree.png / pca_cumulative.png / "
      "pca_scores.png / pca_biplot.png")
print("=" * 72)
plt.show()

以上代码块 5 绘制图③前两个主成分散点图（按潜在因子 f1 着色，用于验证）与图④双标图（灰色样本得分 + 红色载荷箭头，箭头长度按显示比例缩放）。末尾 plt.show() 弹窗显示 4 张图，同时全部已保存到 figures/ 目录。

七、结果解读与注意事项

7.1 实例解读（以上一节代码的实际输出为例）

运行第六节代码，控制台的关键输出如下（np.random.seed(42) 固定，结果可完全复现）。

第一步：数据与标准化

6 个特征原始标准差为 2.7688、2.4904、2.2954、0.9428、0.9362、1.1373——X1~~X3 的波动是 X4~~X6 的 2~3 倍，若不标准化，前 3 个特征必然垄断主成分。标准化后每个特征均值≈0、标准差=1，PCA 将在"公平"的尺度上工作。相关系数矩阵印证了合成结构：X1-X2 达 0.902、X2-X3 为 −0.820、X1-X4 达 0.830（共享因子 f1），X5-X6 为 −0.605（共享因子 f2）——6 个特征背后确有两个潜在因子在驱动。

第二步：特征值与拐点（主成分个数）

特征值（降序）为 3.3808, 2.1004, 0.2466, 0.2243, 0.0639, 0.0142。最大骤降发生在 k=2 处（λ 从 2.1004 骤降至 0.2466），且只有前 2 个特征值大于 1（Kaiser 准则），两法一致判定保留 2 个主成分。第 3 个特征值 0.2466 远小于 1，说明第 3 个方向携带的信息连一个标准化特征都不如，主要是噪声——这正是数据由 2 个潜在因子生成的直接证据。

第三步：累计方差贡献率（85% 阈值）

方差贡献率为 56.07%、34.83%、4.09%、3.72%、1.06%、0.24%，累计贡献率为 56.07%、90.90%、94.99%、98.71%、99.76%、100%。累计贡献率在 k=2 时达 90.90%，同时越过 85% 与 90% 两条阈值线——三种选 k 方法（拐点法、Kaiser 准则、85% 阈值）结论完全一致：k=2。用 2 个主成分替代 6 个特征，只损失约 9.1% 的信息，降维收益极高（前 2 个主成分就浓缩了九成信息，说明特征间冗余严重、数据高度共线）。

第四步：载荷矩阵与业务含义（论文的"亮点"所在）

前 2 个主成分的载荷（特征向量）为：

特征	PC1	PC2
X1	−0.5205	+0.1453
X2	−0.5279	−0.1266
X3	+0.3860	+0.4637
X4	−0.4819	+0.1407
X5	−0.0639	−0.6420
X6	−0.2553	+0.5620

PC1 在 X1、X2、X4 上有较大的同号载荷（−0.52、−0.53、−0.48），X3 上为异号载荷（+0.39）——整体变号后（允许的），PC1 ≈ 0.52X1 + 0.53X2 − 0.39X3 + 0.48X4 + 0.06X5 + 0.26X6，可命名为"综合规模因子"（X1、X2、X4 同涨同跌，X3 与之反向）。PC2 在 X3、X6 上同号（+0.46、+0.56）、X5 上异号（−0.64），可命名为"结构差异因子"（X3、X6 与 X5 此消彼长）。这组命名与数据生成的真实结构完全吻合——PC1 对应潜在因子 f1、PC2 对应潜在因子 f2，验证结果：|corr(PC1, f1)| = 0.9573，|corr(PC2, f2)| = 0.9438，PCA 从 6 个可观测特征中近乎完美地恢复了 2 个不可观测的潜在因子。

第五步：重构误差（信息损失的货币单位）

k	MSE(k)	累计贡献率
0	1.0000	0.0000
1	0.4393	0.5607
2	0.0910	0.9090
3	0.0501	0.9499
4	0.0129	0.9871
5	0.0024	0.9976
6	0.0000	1.0000

MSE 从 k=1 的 0.4393 骤降到 k=2 的 0.0910（信息损失锐减 79%），此后每增加一个主成分只能小步改善——与碎石图拐点完全互证。注意每一行都满足 MSE(k) ≈ 1 − 累计贡献率(k)（标准化数据），两种指标本质是同一件事的两个说法：MSE(2)=0.0910 意味着"平均每个特征还有 9.1% 的方差没被重构出来"。

第六步：sklearn 对照（算法正确性自检）

手写实现与 sklearn.decomposition.PCA 的方差贡献率最大绝对差仅 2.2e-16（机器精度），且各主成分得分相关为 ±1——PC1、PC2、PC3、PC5 相关为 −1.0000（方向相反），PC4、PC6 为 +1.0000。符号相反属正常现象：特征向量 $v_j$ 与 $-v_j$ 都是合法解（目标函数 $w^\top S w$ 与符号无关），不同实现（甚至同一实现的不同库版本）可能给出相反方向。解读时只看载荷的相对大小与符号关系，必要时整体反转某一列再解释。

第七步：不标准化的反例

若跳过标准化直接对原始数据做 PCA：特征值为 16.8825、4.5186、0.2571、0.2487、0.1999、0.0894，PC1 的方差贡献率从 56.07% 膨胀到 76.06%，且 PC1 载荷中弱指标 X4 被压缩到 −0.1834（标准化时为 −0.4819）——方差大的 X1~X3 垄断了主成分，量纲绑架了结论。这就是 2.2 节"必须先标准化"的实测依据。

第八步：降维后的数据去哪儿了

得分矩阵 $Z$ （200 × 2）就是降维后的新特征：可替代原 6 维数据进入聚类、回归、分类等任何下游模型；前两个主成分的散点图（pca_scores.png）可直接放进论文展示样本分布。需要注意得分是基于标准化数据的（均值为 0、方差为 λ），若下游模型对量纲敏感需重新说明。

7.2 常见坑

忘记标准化（最高频错误）：直接用原始特征跑 PCA，量纲大的特征垄断主成分，7.1 第七步的实测对照（PC1 贡献率 76.06% vs 56.07%）就是铁证。对策：PCA 之前一律 z-score 标准化（StandardScaler）；论文中写明"先对指标进行 z-score 标准化以消除量纲影响"；
主成分方向符号相反，误以为是错误：手写实现与 sklearn 的 PC1 相关为 −1.0000，得分全部变号、载荷全部变号，这是数学上完全等价的解，不是 bug。对策：不纠结符号，解读以载荷的相对大小与符号关系为准；若论文要展示载荷表，可约定"将 PC1 载荷调为以最大绝对载荷为正"后整体反转；
阈值机械套用：85% 是经验惯例而非定理。本例 85% 与 90% 都指向 k=2 是巧合；真实数据中阈值取 70%、75%、85%、90% 会给出不同的 k。对策：拐点法、Kaiser 准则、累计贡献率曲线三者互证（本例三法一致，论证最硬）；若第 k+1 个主成分有明确业务含义（载荷结构清晰），即使贡献率小也应保留并说明理由；
降维后含义丢失：主成分是全部特征的线性组合，评审老师问"PC1 代表什么"时必须能答。对策：必画载荷矩阵/双标图并给主成分命名（7.1 第四步的"综合规模因子""结构差异因子"模板）；若实在解释不了，改用因子分析（可旋转）或特征选择；
把 PCA 当聚类：PC1-PC2 散点图上的"看起来的分群"只是可视化观察，不能当作聚类结果写结论。对策：先 PCA 降维，再用 K-means/层次聚类/DBSCAN 聚类，聚类结果与散点图互相印证；
忽视线性假设：特征间只有曲线关系时 PCA 降维无效（贡献率分散、没有明显拐点）。对策：先画散点图矩阵检查线性相关；非线性结构改用核 PCA、t-SNE、UMAP；
离群点污染主成分：协方差矩阵对极端值敏感，一个离群点就能把 PC1 方向拽歪。对策：先箱线图/3σ 筛查；必要时用稳健 PCA；
$n < p$ 或 $n$ 远不够大：样本太少时协方差矩阵估计噪声大，主成分方向不稳定（换个样本集拐点可能漂移）。对策：保证 $n \ge 5p$ （最好 $n \ge 10p$ ）；高维小样本用稀疏 PCA、正则化方法；
把得分当"排名"过度解读：PC1 得分常被当作综合排名，但其方向符号可任意翻转，且 PCA 完全无视业务目标。对策：用 PCA 得分做综合评分时，说明"得分是标准化空间中的相对位置"，排序方向需结合载荷方向人工确认，必要时与熵权法、TOPSIS 等方法交叉验证（见 16_熵权法）。

7.3 竞赛论文写作建议（话术模板）

在论文中，PCA 部分通常按"方法一句话 + 选 k 论证 + 结果描述 + 业务解读"四段式展开。可直接套用以下模板（以本示例数据代入）：

由于 6 项指标间存在较强相关性（相关系数最高达 0.902），直接建模存在多重共线性与维数冗余问题。本文先对指标进行 z-score 标准化以消除量纲影响，再采用主成分分析进行降维。特征值依次为 3.3808、2.1004、0.2466、0.2243、0.0639、0.0142，碎石图在 k=2 处出现明显拐点，且仅有前 2 个特征值大于 1（Kaiser 准则）；累计方差贡献率在 k=2 时达 90.90%，超过 85% 经验阈值。三种准则相互印证，故取 k=2，保留原数据约 90.9% 的信息。载荷矩阵显示：第一主成分在 X1、X2、X4 上载荷较大且同号（−0.52、−0.53、−0.48），可解释为"综合规模因子"；第二主成分在 X3、X6 上载荷为正（0.46、0.56）、X5 上为负（−0.64），可解释为"结构差异因子"。两主成分与原数据的重构均方误差仅 0.091，降维损失可控，后续分析将在前 2 个主成分得分上展开。

变体话术（按题型替换划线句）：

降维后再建模："……将前 2 个主成分得分作为新特征输入聚类模型，聚类轮廓系数较直接使用 6 维原始指标提升 0.15，说明 PCA 降维去除了噪声维度，改善了聚类结构。"
多指标综合评价："……按方差贡献率加权（0.5607、0.3483）计算综合得分 $F = 0.5607 \cdot z_1 + 0.3483 \cdot z_2$ ，得分越高代表综合水平越高，排序结果与专家经验一致。"
共线性处理（主成分回归）："……原始指标 VIF 最高达 25（严重共线），改用前 2 个主成分得分对因变量做回归，模型 R² 为 0.xx 且各系数均显著，避免了共线性导致的系数不稳定问题。"
谨慎表态（累计贡献率偏低时）："……前 k 个主成分累计贡献率为 xx%，虽低于 85% 经验阈值，但碎石图拐点明确且主成分业务含义清晰，本文仍取 k，并保留原特征作为对照，主要结论不依赖降维结果。"

八、延伸阅读

核 PCA（Kernel PCA）：先用核技巧把数据映射到高维（甚至无穷维）特征空间再做 PCA，等价于对非线性相似度矩阵做特征分解，能捕捉曲线结构。竞赛中当线性 PCA 的累计贡献率爬升缓慢（特征间非线性相关）时考虑；代价是失去线性可解释性、需选核函数与带宽（sklearn.decomposition.KernelPCA）；
稀疏 PCA（Sparse PCA）：在 PCA 目标上加入 L1 正则，强制每个主成分只用到极少数原特征，兼顾降维与可解释性，适合"既要压缩、又要能说出每个主成分由哪几个指标组成"的场景（sklearn.decomposition.SparsePCA）；
因子分析（Factor Analysis）：把观测变量建模为少数公共因子的线性组合 + 独立特殊噪声（ $x = \Lambda f + \varepsilon$ ），关注解释协方差结构而非最大化方差，支持 varimax 等旋转使因子载荷更"干净"、更易命名。竞赛中做"指标结构分析""公共因子命名"时比 PCA 更合适（sklearn.decomposition.FactorAnalysis）；
t-SNE / UMAP：非线性降维，专长是保持样本间局部邻近关系，二维图上簇结构清晰，是高维数据可视化的事实标准。注意二者不能提供可复用的线性变换（新样本需重新拟合），只用于画图、不用于后续建模；t-SNE 有随机性、超参数敏感（perplexity），UMAP 更快且保留更多全局结构；
主成分回归（PCR, Principal Component Regression）：PCA + 线性回归的组合拳——先对自变量做主成分，再用前 $k$ 个主成分得分回归。彻底消除共线性，但丢弃的小主成分可能恰恰对 $y$ 有预测力（PCA 只保方差、不保预测力），此时偏最小二乘回归（PLS）（同时利用 $X$ 与 $y$ 的结构）通常是更好的替代，可与 06_岭回归对照学习。

本文小结：PCA 把 $p$ 个相关特征通过"方差最大投影"压成 $k$ 个互不相关的主成分，五步实现（中心化 → 协方差 → 特征分解 → 取前 $k$ → 投影），核心指标是特征值、方差贡献率、累计贡献率（85% 经验阈值）与重构误差。三个关键习惯：建模前必标准化；选 $k$ 用拐点法 + Kaiser 准则 + 累计贡献率三者互证；解读必画碎石图、累计曲线、散点图与双标图并用载荷给主成分命名。降维后得分可直接进入聚类、回归等下游模型，让高维数据的建模"先瘦身、再发力"。

一、算法含义​

1.1 通俗理解​

1.2 数学模型：方差最大的投影方向​

1.3 算法步骤（五步）​

1.4 几何直觉​

1.5 与特征选择（feature selection）的区别​

1.6 优缺点​

二、何时使用（适用场景与条件）​

2.1 适用场景​

2.2 使用前提（建模前检查清单）​

2.3 不适用 / 慎用的情形​

2.4 与因子分析、t-SNE 的对比与选择​

三、算法指标​

3.1 特征值 λ_j​

3.2 方差贡献率​

3.3 累计方差贡献率​

3.4 重构误差 MSE(k)​

3.5 载荷矩阵​

3.6 指标汇总表​

四、可视化图表​

4.1 四张图速查表​

4.2 每张图的详细解读要点​

五、符号说明​

六、可运行程序（完整代码）​

七、结果解读与注意事项​

7.1 实例解读（以上一节代码的实际输出为例）​

7.2 常见坑​

7.3 竞赛论文写作建议（话术模板）​

八、延伸阅读​