高光谱成像（三）主成分分析 PCA

我们已经知道：高光谱图像的每一个像素，本质上是一个高维光谱向量，其波段数量可能达到上百维。
这固然提供了更多的信息，但于此同时也极大增加了计算量。
并且，结合实际生活，你会发现一个新的问题：

大量特征间的重要性并不等价。

举个例子：假设我们要区分三种动物：猫、狗和老虎。
现在我们为每只动物记录一系列特征，如体重、体长、毛发颜色、是否有胡须、是否会游泳、心跳频率、是否有条纹、爪子长度、是否食肉等等。
仔细想想就会发现，这些特征在这个分类任务中并非同等重要，而高光谱图像也是如此，举例如下：
动物分类高光谱图像“是否有条纹”对区分猫和老虎非常关键某些波段区分度高，信息量丰富“是否食肉”对三者几乎没有区分度某些波段对分类贡献有限“是否会游泳”可能完全无关某些波段可能属于无关信息“体重”和“体长”高度相关，信息重复部分波段之间高度相关，存在冗余个别测量误差可能影响判断部分波段可能包含噪声也就是说，虽然我们通过采集高光谱图像获得了每个像素的更多特征，但是同时也会发现这种情况：

光谱向量中有的特征信息量大，有的特征几乎是噪声，还有的特征彼此冗余。

如果我们把这些特征全部平等对待，不仅计算复杂，而且可能让真正重要的差异被淹没。
那么，一个自然的问题便出现了，那就是能不能自动找出“最重要”的那些方向？
这其实本质上是计算机领域的经典的优化问题：

如何在尽可能保留信息的前提下，降低维度来减少计算成本？

于是便有了非常经典的方法： 主成分分析（Principal Component Analysis，PCA）。
主成分分析是一种非常通用的线性降维方法，它广泛应用于几乎所有涉及“高维数据”的领域。高光谱成像只是其中一个典型应用场景。
1. 什么是 PCA ？

PCA 的起源极早，在1901 年，论文 On Lines and Planes of Closest Fit to Systems of Points in Space 就提出了一种用于寻找多维数据的最优拟合方向的数学方法，这就是 PCA 的几何雏形。
而后，1933 年的论文 Analysis of a Complex of Statistical Variables into Principal Components （链接仅为 part1 ）对 PCA 进行了系统推广，延续至今，PCA 仍是现代统计学与机器学习中的经典降维方法。
在开始阐述前，我们要强调的是：PCA 的降维逻辑不是在原数据中找到最有区分度的几个特征后删去其他特征，而是用所有原始特征线性组合生成区别度高的新特征。
而把高维空间的信息压缩到低维空间，“投影”是最自然的压缩方式。
现在，我们举一个最简单的例子：把二维信息降维到一维，来说明 PCA 的核心思想，来看：

如果不计代价，我们很容易就能做到，就像这样：

我们以一种十分简单的方式完成了降维，但就像图中所示，如果我们只取样本在原坐标轴上的投影，因为坐标轴间的正交关系，就会完全损失样本在其他维度上的信息。
因此，如果想在降维的同时保留更多信息，我们就要跳出原本的坐标轴，或者说变换坐标轴，让二维信息综合在新的一维坐标轴上。 就像这样：

这样，我们就实现了在降维数据来减少计算量的同时，在一定程度上保留了信息。
但新的问题马上就来了：

什么方向上的坐标轴，可以让信息得到最大程度的保留呢？

显然，不同方向的坐标轴，其保留信息的能力是不同的，我们继续刚刚的例子：

由此，我们可以总结出一个关键点：合理的降维方法应该在降维后仍能体现出原数据在向量空间里的相对距离关系。
那么，我们要如何量化这个“合理”的程度呢？PCA 给出的答案是：信息主要存在于数据变化最大的方向上。 来看看为什么：

由此，我们的目标就确定了：找到那些能够保留最大方差的方向。这些方向能够最大程度地保持数据的变化结构，从而尽可能保留原始数据中的信息。
这就是 PCA 的核心思想，下面就来分步阐述其具体过程：
2. PCA 的过程

我们仍以高光谱图像来演示这个过程，假定采集 \(n\) 个波段，那么一个像素就可以表示为：

\[\mathbf{x} = (I_1, I_2, \dots, I_n)  \]
再推广到整幅图像，如果有 \(m\) 个像素，则形成数据矩阵：

\[X =  \begin{bmatrix}  I_{11} & I_{12} & \dots & I_{1n} \\  I_{21} & I_{22} & \dots & I_{2n} \\  \vdots & \vdots & \ddots & \vdots \\  I_{m1} & I_{m2} & \dots & I_{mn}  \end{bmatrix}_{m \times n}  \]
其中，每行代表一个像素，每列代表一个波段。
因此，数据矩阵的每一行就是一个样本，每一列就是一个特征。
2.1 数据中心化

对数据矩阵的第一步处理并不复杂，一句话概括：对每个特征（列）进行零均值处理，在数学和机器学习里，这个操作通常叫做中心化（mean centering）。
用三个像素点来进行进行简单举例演示，每个像素包含 100 个波段：
像素编号波段1波段2…波段100Pixel 10.120.35…0.48Pixel 20.100.37…0.50Pixel 30.110.36…0.47对每一列数据进行零均值处理，即：

\[\tilde{I}_{ij} = I_{ij} - \mu_j  \]
其中，\(\mu_j\) 为该列的均值，\(\tilde{I}_{ij}\) 即为该位置中心化后的数据。
已知第 1 波段三个像素的反射率为 0.12、0.10、0.11，那么：$$\mu_1 = (0.12+0.10+0.11)/3 = 0.11$$当列数据与该均值分别做差，即可得到中心化后的数据，表格如下：
像素编号波段1波段2波段3…波段100Pixel 10.01-0.01……-0.003Pixel 2-0.010.01……0.017Pixel 30.000.00……-0.013中心化消除了波段间的均值差异，让每个波段的数据均值为 0。
这样做减去了数据的偏移，在更准确地反映数据变化方向的同时，也方便了下一步协方差矩阵计算。
2.2 计算协方差矩阵

计算协方差矩阵之前，先简单摆一下协方差的定义，对于两个随机变量 \(X\) 和 \(Y\)，协方差定义为：

\[\text{Cov}(X,Y) = \frac{1}{m-1} \sum_{i=1}^{m} (X_i - \bar{X})(Y_i - \bar{Y})  \]
其中：

• \(X_i, Y_i\) 是第 \(i\) 个样本的值。
  
• \(\bar{X}, \bar{Y}\) 是样本均值。
  
• \(m\) 是样本数。
  

说人话，就是两个变量中心化后对应样本的乘积平均值。
现在，再来看协方差矩阵，我们已经有了中心化后的数据矩阵：

\[\tilde X =  \begin{bmatrix}  \tilde{I}_{11} & \tilde{I}_{12} & \dots & \tilde{I}_{1n} \\  \tilde{I}_{21} & \tilde{I}_{22} & \dots & \tilde{I}_{2n} \\  \vdots & \vdots & \ddots & \vdots \\  \tilde{I}_{m1} & \tilde{I}_{m2} & \dots & \tilde{I}_{mn}  \end{bmatrix}_{m \times n}  \]
那么协方差矩阵公式就是：

\[\mathbf{C} = \frac{1}{m-1} \tilde X^\top \tilde X  \]
显然，得到的协方差矩阵 \(C\) 是 \(n \times n\) 的矩阵,其中的元素 \(C_{ij}\) 就是波段 \(i\) 与波段 \(j\) 的协方差。
还没完，我们再来专门看看协方差矩阵的对角线元素 \(C_{ii}\)，将其代入公式后，协方差公式就变成了：

\[ \frac{1}{m-1} \sum_{i=1}^{m} (X_i - \bar{X})^2  \]
显然，这就是方差公式，也就是说：协方差矩阵的对角线元素是每个特征自身的方差。而非对角的元素就可以看作”两个特征联合的波动程度“，也就是协方差。
到这里，我们就完成了 PCA 降维思想中，对方差的准备部分。
2.3 特征值分解（Eigenvalue Decomposition，EVD）

在得到协方差矩阵 \(C\) 后，下一步就是对其进行特征值分解，以提取数据中的主成分方向，这是 PCA 的最核心步骤。
首先，在多维向量空间中，我们会定义单位向量 \(\mathbf{w}\)，其长度固定为 1，目的是保证我们只选择方向，而不会因为向量长度不同而改变投影尺度。

\[|\mathbf{w}| = 1  \]
现在，先不提上一步的协方差矩阵，现在假定我们只有一个样本 \(\mathbf{x}\) ，现在，我们要计算 \(\mathbf{x}\)
在某个方向上的投影，其公式就是：

\[z = \mathbf{w}^T \mathbf{x}  \]
得到的 \(z\) 即为 \(\mathbf{x}\) 在该方向上的投影数值，是一个标量。
但我们要计算的是 投影后的方差，一个数怎么计算方差呢？这便涉及到线代里的结论：

\[{Var(z)} = \mathbf{w}^T C \mathbf{w}  \]
其中，\(C\) 就是我们上一步得到的协方差矩阵。
继续下一步就是我们的目标：找到投影后最大的方差 即：

\[\max_{|\mathbf{w}|=1} \mathbf{w}^T C \mathbf{w}  \]
现在，问题就转变成了一个求最值，也就是优化问题，在数学上，其解决方法是用拉格朗日法构造函数：

\[L(\mathbf{w},\lambda) = \mathbf{w}^T C \mathbf{w} - \lambda (\mathbf{w}^T \mathbf{w} - 1)  \]
然后对 \(w\) 求导并让其为 0 得到：

\[C \mathbf{w} = \lambda \mathbf{w}  \]
这个式子就是特征值方程，简单补充一下：
给定一个方阵 A(行列数相等的矩阵)，假设我们有一个非零的向量 \(\mathbf{v}\)，如果这个向量通过矩阵变换后只是按比例放大或缩小，即：

\[A \mathbf{v} = \lambda \mathbf{v}  \]

• \(\mathbf{v}\) 就是矩阵的特征向量。
  
• \(\lambda\) 就是对应的特征值。
  

需要强调的是，这里虽然使用特征值方程，但特征值方程本身并没有“最大方差方向”的语义，只是我们的数学推导恰好得到了这个形式。
继续，把特征值方程代回方差公式，就得到了：

\[{Var(z)} = \mathbf{w}^T C \mathbf{w} = \mathbf{w}^T (\lambda \mathbf{w})=\lambda \mathbf{w}^T \mathbf{w}= \lambda   \]
现在，我们就得出了最终结论：特征值 = 这个方向上的方差，而与特征值对应的特征向量的方向，就是我们要找的方差最大的方向。
自此，我们就完成了 PCA 最核心的内容，要说明的是，PCA 的理论基础是协方差矩阵的特征值分解，但现代机器学习库通常使用另一种叫奇异值分解来完成计算，因为其省略了显式计算协方差矩阵的过程，但二者的结果是相同的。
2.4 选择主成分

根据线代的基本知识，我们知道一个 \(n \times n\) 的方阵，其特征值方程最多有 \(n\) 个解。
而协方差矩阵是对称矩阵的同时也是半正定矩阵，因此其特征值方程一定有 \(n\) 个正交的解。
也就是说，在上一步，我们最终可以得到 \(n\) 组特征值和对应的特征向量。
到这里，我们就称这些特征向量所代表的方向为“成分”。
特征值 的大小决定每个主成分的重要性。特征值越大，意味着该主成分方向承载的数据方差越多，信息量也越大。
最大特征值对应的特征向量，就是第一主成分方向，第二大的特征值对应第二主成分方向，以此类推。
而我们在这一步的工作如下：

• 计算所有特征值之和 = 数据的总方差
  
• 计算每个主成分的方差占比并排序
  
• 一般挑选满足方差占比之和大于 95% 的前 \(k\) 个主成分，最终降维到 \(k\) 维。
  

来看一个简单的例子，假设我们对数据做完特征值分解后，得到如下结果：
排名特征值 \(\lambda\)解释的方差比例15.05.0 / (5.0+0.8+0.2)=83.3%20.80.8 / (5.0+0.8+0.2)=13.3%30.20.2 / (5.0+0.8+0.2)=3.3%你会发现：前两个主成分已经解释了约 96% 的信息
这时我们就可以说：数据的主要结构，几乎都集中在前两个方向上。
因此：如果我们把数据从 3 维降到 2 维，仍然可以保留 96% 的信息，丢失的信息只有 4%，这就是 PCA 降维的核心思想。
到此，我们就只剩下最后一步了。
2.5 构建降维后的数据

在得到重要的主成分方向后，最后就是将原始数据投影到这些主成分方向上，得到降维后的数据了。这个过程就是将原始数据映射到新的低维空间中。
投影公式为：

\[Y = \tilde{X} \mathbf{V}_k  \]
其中：

• \(Y\) 是降维后的数据矩阵，尺寸为 \(m \times k\)。
  
• \(\mathbf{V}_k\) 是包含前 \(k\) 个主成分的特征向量矩阵，尺寸为 \(n \times k\)。
  
• \(Y\) 每行对应一个样本，通过前 \(k\) 个主成分表示。
  

我们继续上一步的例子来演示，假定主成分对应的单位特征向量为：

\[\mathbf{v}_1 =  \begin{bmatrix}  0.58 \\  0.58 \\  0.58  \end{bmatrix}  \quad  \mathbf{v}_2 =  \begin{bmatrix}  0.71 \\  -0.71 \\  0  \end{bmatrix}  \]
然后将两个特征向量组合为矩阵：

\[\mathbf{V}_k =  \begin{bmatrix}  0.58 & 0.71 \\  0.58 & -0.71 \\  0.58 & 0  \end{bmatrix}  \]
现在，假设中心化后的数据为：

\[\tilde{X} =  \begin{bmatrix}  2 & 1 & 1 \\  1 & 0 & 1 \\  3 & 1 & 2  \end{bmatrix}  \]
代入PCA 的投影公式即可得到：

\[Y =  \begin{bmatrix}  2 & 1 & 1 \\  1 & 0 & 1 \\  3 & 1 & 2  \end{bmatrix}  \begin{bmatrix}  0.58 & 0.71 \\  0.58 & -0.71 \\  0.58 & 0  \end{bmatrix}  \]
计算结果为：

\[Y =  \begin{bmatrix}  2.32 & 0.71 \\  1.16 & 0.71 \\  3.48 & 1.42  \end{bmatrix}  \]
最终，我们就在尽可能保留信息的前提下，将数据从 3 维压缩到 2 维。
3.PCA 的优缺

天底下没有完美的方法，PCA 自然也有其优劣所在。
3.1 PCA 的优势

优点说明计算效率高PCA 主要依赖协方差矩阵的特征值分解，算法结构简单、无需复杂迭代，计算稳定且易实现能去除特征冗余通过协方差矩阵分析特征之间的相关性，将高度相关的特征压缩到少数几个主成分中数据压缩能力强在尽量保留主要方差信息的前提下降低数据维度，从而减少计算量和存储需求主成分互不相关PCA 得到的主成分两两正交，因此在统计意义上互不相关，有利于后续建模3.2 PCA 的不足

缺点说明只能处理线性结构PCA 本质是线性投影方法，对于存在非线性流形结构的数据表达能力有限对特征尺度敏感协方差计算会受到特征量纲影响，因此通常需要先进行标准化处理主成分可解释性较弱每个主成分通常是多个原始特征的线性组合，难以直接解释其物理或业务含义对异常值敏感方差对极端值非常敏感，少量异常样本可能显著改变主成分方向总之，PCA 通过特征值分解和主成分选择，从高维数据中提取最重要的信息，极大地减少了计算复杂度。它不仅用于高光谱图像的降维，也广泛应用于其他数据分析和机器学习任务中。以传统的 PCA 为基础，后续也出现出了Kernel PCA，稀疏 PCA等改进方法。

来源：程序园用户自行投稿发布，如果侵权，请联系站长删除
免责声明：如果侵犯了您的权益，请联系站长，我们会及时删除侵权内容，谢谢合作！