PCA 主成分分析
PCA 是"换一组坐标轴看数据"：不改变数据本身，只找最能拉开差异的方向。 这一课讲清它在脑电里的两类用法、成分数怎么选，以及它和 ICA 到底该选谁。
🎯 本节学习目标
- 能用"换坐标轴"的直觉解释 PCA 在做什么，说清特征向量与特征值各自的角色
- 区分脑电中空间 PCA 与时间 PCA 的输入、输出与用途
- 能运行手写 eig 版 PCA，并与内置 pca() 的结果对照验证
- 掌握三种成分数选择方法：碎石图、Kaiser 准则、平行分析
- 会用 PCA vs ICA 对比表决定何时选谁，并理解 varimax 轮换的价值
一、PCA 直觉：换一组坐标轴看数据
给一群高矮胖瘦各不相同的人拍照：从正面拍，大家的高矮拉开了；从头顶拍， 高矮全被压扁，胖瘦倒是分明了。PCA 做的就是自动寻找"把数据拍得最分散" 的那几个拍摄角度——数据本身一个字都不改，只是换一组坐标轴来描述它。
新坐标系的一条轴。第一主成分是数据投影方差最大的方向，第二主成分在与它垂直的方向上方差最大，依此类推。前几个主成分往往就抓住了数据的"主要故事"。
对脑电而言这意味着什么？64 个通道的数据表面上活在 64 维空间里， 但真实的头皮分布模式可能只有少数几种——PCA 用几个主成分就能近似还原 绝大部分信息，降维与去噪的价值正在于此。
二、数学骨架（直觉版）：协方差矩阵的特征分解
整个 PCA 的数学只用到一句线代，不做推导，只要记住对应关系：
- 数据去均值后算协方差矩阵 C：对通道数据而言，它是一个"通道 × 通道"的表，记录谁和谁一起涨落；
- 对 C 做特征分解 C = EΛE′：每个特征向量是一根新坐标轴（一种空间模式），对应的特征值是数据在这根轴上的方差（这种模式的"势力"）；
- 把特征值从大到小排序：第一大的特征向量就是第一主成分。成分的"解释方差比例" = 该特征值 / 所有特征值之和。
就这么多。PCA 不做统计假设、结果唯一（除符号与顺序），这是它比 ICA"老实"的地方。
三、脑电中的两类用法：空间 PCA 与时间 PCA
| 类型 | 输入矩阵 | 得到什么 | 典型用途 |
|---|---|---|---|
| 空间 PCA | 通道 × 时间点 | 成分的空间模式（地形图）+ 各成分的时间过程 | 降维、去噪、提取主要头皮分布 |
| 时间 PCA | 被试（或条件）× 时间点 | 成分的时间波形 + 每个被试的载荷 | 分解重叠的 ERP 成分（N1、P3 混在一条波形里时） |
时间 PCA 是 ERP 研究的经典操作：不同被试的 P3 潜伏期略有差异， 直接平均会把波形抹胖；先用 PCA 提取"潜伏期可伸缩"的成分波形， 再比较各被试在成分上的载荷，就能把"成分强度"与"潜伏期差异"分开。 两者还可以组合成时空 PCA，进阶时再学。
四、Matlab 实现：手写 eig 与内置 pca() 对照
% ---- 手写空间 PCA：对"通道 x 时间点"数据找主要空间模式 ----
X = EEG.data; % 通道 x 时间点
Xc = X - mean(X, 2); % 每个通道沿时间去均值（标准前处理）
C = (Xc * Xc') / (size(Xc, 2) - 1); % 通道 x 通道协方差：谁和谁一起涨落
[Evec, Eval] = eig(C); % 特征分解
ev = diag(Eval);
[ev, order] = sort(ev, 'descend');
Evec = Evec(:, order); % Evec 第 k 列 = 第 k 主成分的空间模式
score = Evec' * Xc; % 各成分的时间过程（投影得分）
explVar = ev / sum(ev); % 各成分解释方差比例
fprintf('前 3 个成分解释方差: %s\n', mat2str(explVar(1:3)*100, 3));% ---- 用内置 pca() 对照（需要 Statistics Toolbox）----
% 注意：pca 按"行 = 观测、列 = 变量"，空间 PCA 要传入转置
[coeff, score2, latent] = pca(Xc');
% coeff(:,k)：第 k 主成分的载荷（空间 PCA 下即空间模式）
% 与手写版只差一个正负号——特征向量方向本来就 ±等价
sgn = sign(coeff(1,1) * Evec(1,1));
disp(max(abs(coeff(:,1) - Evec(:,1) * sgn))) % 应接近 0
% 碎石图：解释方差在哪开始"变平"（肘部法则）
figure; plot(1:numel(latent), 100*latent/sum(latent), 'o-');
xlabel('主成分序号'); ylabel('解释方差 (%)'); title('碎石图 Scree Plot');空间 PCA 把时间点当观测、通道当变量；时间 PCA 把被试当观测、时间点当变量。传给 pca() 前先问自己一句"我到底想压缩哪个维度"，转置与否全由它决定。变量量纲不一时（时间 PCA 中各时间点已经同量纲，一般无需处理），先标准化再分析。
五、成分数怎么选
- 碎石图肘部：解释方差曲线由陡变平的拐点，最直观也最主观；
- Kaiser 准则：保留特征值大于 1 的成分（前提：变量已标准化），偏宽松，容易多留；
- 平行分析（最推荐）：把数据随机打乱（或用随机数）重复做 PCA，只有"真实特征值超过随机水平"的成分才保留——给主观判断配了一个客观参照。
成分数不是数学问题，是"信噪比 + 研究目的"的问题。降维去噪可以宽松些，解释成分则要保守些；论文里报告你用了哪种方法、以及换一种方法结论是否稳定。
六、PCA vs ICA：何时选谁
| 维度 | PCA | ICA |
|---|---|---|
| 优化目标 | 最大化投影方差 | 最大化统计独立性 |
| 用到的统计量 | 二阶（协方差） | 高阶（非高斯性） |
| 成分之间 | 正交（不相关） | 不要求正交 |
| 结果稳定性 | 唯一（除符号与顺序） | 受初始化与数据量影响 |
| 脑电典型用途 | 降维、去噪、ERP 成分提取 | 伪迹分离（眼电 / 心电）、源分解 |
| 何时选它 | 想压缩维度、找主要模式 | 想把混合在一起的"源"拆开 |
PCA 只保证成分二阶不相关；ICA 追求完全统计独立。眼电与脑电可能互不相关却仍然混合在一起——这正是模块 2 里 ICA 能拆开眼电、而 PCA 不能的原因。两个工具没有高下，只有分工。
七、轮换：让成分更好解释
轮换前的前几个主成分常是"过渡混合体"：PC1 一半像 N1、一半像 P3， 没法命名。varimax（正交轮换）与 promax（斜交轮换） 通过旋转坐标轴让载荷趋于"简单结构"——每个成分只在少数时间点（或通道）上 有高载荷，其余接近零，于是成分变得像一个个可命名的 ERP 成分。 轮换不改变总解释方差，只是把它在前 k 个成分间重新分配。
% ---- 时间 PCA + varimax：让成分波形"像 ERP 成分" ----
% 场景：Xz 为"被试 x 时间点"的 ERP 矩阵（已按时间点标准化）
[coeffT, ~, ~] = pca(Xz); % 时间 PCA：主成分是时间波形
coeffR = rotatefactors(coeffT(:, 1:4), 'method', 'varimax');
% 轮换只作用于保留的前 k 个成分；旋转后载荷更接近简单结构
figure; plot(times, coeffR);
legend('PC1', 'PC2', 'PC3', 'PC4');
xlabel('时间 (ms)'); title('varimax 轮换后的成分波形');
% 轮换后各波形应在不同时间窗集中，便于对应 N1 / P2 / P3 等成分📌 本节小结
- PCA = 换一组坐标轴：特征向量是新轴（模式），特征值是该轴上的方差（势力）
- 空间 PCA 压缩通道维度，时间 PCA 分解重叠的 ERP 成分——先想清楚压缩谁
- 手写 eig 与内置 pca() 结果一致（只差符号），两边对照是最好的理解方式
- 成分数三法：碎石图看肘部、Kaiser 偏宽松、平行分析最推荐
- 不相关 ≠ 独立：PCA 去相关，ICA 求独立——降维用前者，拆源用后者
- varimax 轮换把"过渡混合"的主成分变成可命名的简单结构
✏️ 课后练习
- 对 64 通道静息态数据做空间 PCA，用
topoplot画出前 4 个成分的地形图与时间过程。 - 画碎石图并用平行分析确定保留几个成分，比较 Kaiser 准则给出的答案。
- 构造两个潜伏期不同、部分重叠的 ERP 波形（加噪声混合成"被试 x 时间"矩阵），用时间 PCA + varimax 尝试还原两个成分。