高光谱成像基础（三）主成分分析 PCA
<p>我们已经知道：高光谱图像的每一个像素，本质上是一个高维光谱向量，其波段数量可能达到上百维。<br>这固然提供了更多的信息，但于此同时也极大增加了计算量。<br>
并且，结合实际生活，你会发现一个新的问题：</p>
<blockquote>
<p><strong>大量特征间的重要性并不等价。</strong></p>
</blockquote>
<p>举个例子：假设我们要区分三种动物：猫、狗和老虎。<br>
现在我们为每只动物记录一系列特征，如体重、体长、毛发颜色、是否有胡须、是否会游泳、心跳频率、是否有条纹、爪子长度、是否食肉等等。</p>
<p>仔细想想就会发现，这些特征在这个分类任务中并非同等重要，而高光谱图像也是如此，举例如下：</p>
<table>
<thead>
<tr>
<th><strong>动物分类</strong></th>
<th><strong>高光谱图像</strong></th>
</tr>
</thead>
<tbody>
<tr>
<td>“是否有条纹”对区分猫和老虎非常关键</td>
<td>某些波段区分度高，信息量丰富</td>
</tr>
<tr>
<td>“是否食肉”对三者几乎没有区分度</td>
<td>某些波段对分类贡献有限</td>
</tr>
<tr>
<td>“是否会游泳”可能完全无关</td>
<td>某些波段可能属于无关信息</td>
</tr>
<tr>
<td>“体重”和“体长”高度相关，信息重复</td>
<td>部分波段之间高度相关，存在冗余</td>
</tr>
<tr>
<td>个别测量误差可能影响判断</td>
<td>部分波段可能包含噪声</td>
</tr>
</tbody>
</table>
<p>也就是说，虽然我们通过采集高光谱图像获得了每个像素的更多特征，但是同时也会发现这种情况：</p>
<blockquote>
<p><strong>光谱向量中有的特征信息量大，有的特征几乎是噪声，还有的特征彼此冗余。</strong></p>
</blockquote>
<p><strong>如果我们把这些特征全部平等对待，不仅计算复杂，而且可能让真正重要的差异被淹没。</strong></p>
<p>那么，一个自然的问题便出现了，那就是能不能自动找出“最重要”的那些方向？<br>
这其实本质上是计算机领域的经典的优化问题：</p>
<blockquote>
<p><strong>如何在尽可能保留信息的前提下，降低维度来减少计算成本？</strong></p>
</blockquote>
<p>于是便有了非常经典的方法： <strong>主成分分析（Principal Component Analysis，PCA）</strong>。<br>
主成分分析是一种非常通用的线性降维方法，它广泛应用于几乎所有涉及“高维数据”的领域。高光谱成像只是其中一个典型应用场景。</p>
<hr>
<h1 id="1-什么是-pca-">1. 什么是 PCA ？</h1>
<p>PCA 的起源极早，在1901 年，论文 <strong><em>On Lines and Planes of Closest Fit to Systems of Points in Space</em></strong> 就提出了一种用于寻找多维数据的最优拟合方向的数学方法，这就是 PCA 的几何雏形。<br>
而后，1933 年的论文 <strong><em>Analysis of a Complex of Statistical Variables into Principal Components</em></strong> （链接仅为 part1 ）对 PCA 进行了系统推广，延续至今，PCA 仍是现代统计学与机器学习中的经典降维方法。</p>
<p>在开始阐述前，我们要强调的是：PCA 的降维逻辑不是在原数据中找到最有区分度的几个特征后删去其他特征，而是<strong>用所有原始特征线性组合生成区别度高的新特征</strong>。<br>
而把高维空间的信息压缩到低维空间，“<strong>投影</strong>”是最自然的压缩方式。</p>
<p>现在，我们举一个最简单的例子：把二维信息降维到一维，来说明 PCA 的核心思想，来看：<br>
<img src="https://img2024.cnblogs.com/blog/3708248/202603/3708248-20260305155518870-1989614428.png" alt="image.png" loading="lazy"><br>
如果<strong>不计代价</strong>，我们很容易就能做到，就像这样：<br>
<img src="https://img2024.cnblogs.com/blog/3708248/202603/3708248-20260305155518519-1928796848.png" alt="image.png" loading="lazy"><br>
我们以一种十分简单的方式完成了降维，但就像图中所示，如果我们只取样本在原坐标轴上的投影，因为坐标轴间的正交关系，就会完全损失样本在其他维度上的信息。</p>
<p><strong>因此，如果想在降维的同时保留更多信息，我们就要跳出原本的坐标轴，或者说变换坐标轴，让二维信息综合在新的一维坐标轴上。</strong> 就像这样：<br>
<img src="https://img2024.cnblogs.com/blog/3708248/202603/3708248-20260305155518885-1616637769.png" alt="image.png" loading="lazy"><br>
这样，我们就实现了在降维数据来减少计算量的同时，在一定程度上保留了信息。<br>
但新的问题马上就来了：</p>
<blockquote>
<p><strong>什么方向上的坐标轴，可以让信息得到最大程度的保留呢？</strong></p>
</blockquote>
<p>显然，不同方向的坐标轴，其保留信息的能力是不同的，我们继续刚刚的例子：<br>
<img src="https://img2024.cnblogs.com/blog/3708248/202603/3708248-20260305155518534-802585489.png" alt="image.png" loading="lazy"><br>
由此，我们可以总结出一个关键点：<strong>合理的降维方法应该尽可能保留数据原本的变化结构。</strong></p>
<p>那么，我们要如何量化这个“合理”的程度呢？PCA 给出的答案是：<strong>信息主要存在于数据变化最大的方向上。</strong> 来看看为什么：<br>
<img src="https://img2024.cnblogs.com/blog/3708248/202603/3708248-20260305155518486-1206138105.png" alt="image.png" loading="lazy"></p>
<p>由此，我们的目标就确定了：<strong>找到那些能够保留最大方差的方向</strong>。这些方向能够最大程度地保持数据的变化结构，从而尽可能保留原始数据中的信息。<br>
这就是 PCA 的核心思想，下面就来分步阐述其具体过程：</p>
<h1 id="2-pca-的过程">2. PCA 的过程</h1>
<p>我们仍以高光谱图像来演示这个过程，假定采集 <span class="math inline">\(n\)</span> 个波段，那么一个像素就可以表示为：</p>
<p></p><div class="math display">\[\mathbf{x} = (I_1, I_2, \dots, I_n)
\]</div><p></p><p>再推广到整幅图像，如果有 <span class="math inline">\(m\)</span> 个像素，则形成数据矩阵：</p>
<p></p><div class="math display">\[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}
\]</div><p></p><p>其中，每行代表一个像素，每列代表一个波段。<br>
<strong>因此，数据矩阵的每一行就是一个样本，每一列就是一个特征。</strong></p>
<h2 id="21-数据中心化">2.1 数据中心化</h2>
<p>对数据矩阵的第一步处理并不复杂，一句话概括：<strong>对每个特征（列）进行零均值处理，在数学和机器学习里，这个操作通常叫做中心化（mean centering）。</strong></p>
<p>用三个像素点来进行进行简单举例演示，每个像素包含 100 个波段：</p>
<table>
<thead>
<tr>
<th>像素编号</th>
<th>波段1</th>
<th>波段2</th>
<th>…</th>
<th>波段100</th>
</tr>
</thead>
<tbody>
<tr>
<td>Pixel 1</td>
<td>0.12</td>
<td>0.35</td>
<td>…</td>
<td>0.48</td>
</tr>
<tr>
<td>Pixel 2</td>
<td>0.10</td>
<td>0.37</td>
<td>…</td>
<td>0.50</td>
</tr>
<tr>
<td>Pixel 3</td>
<td>0.11</td>
<td>0.36</td>
<td>…</td>
<td>0.47</td>
</tr>
</tbody>
</table>
<p>对每一列数据进行零均值处理，即：</p>
<p></p><div class="math display">\[\tilde{I}_{ij} = I_{ij} - \mu_j
\]</div><p></p><p>其中，<span class="math inline">\(\mu_j\)</span> 为该列的均值，<span class="math inline">\(\tilde{I}_{ij}\)</span> 即为该位置中心化后的数据。</p>
<p>已知第 1 波段三个像素的反射率为 0.12、0.10、0.11，那么：$$\mu_1 = (0.12+0.10+0.11)/3 = 0.11$$当列数据与该均值分别做差，即可得到中心化后的数据，表格如下：</p>
<table>
<thead>
<tr>
<th>像素编号</th>
<th>波段1</th>
<th>波段2</th>
<th>波段3</th>
<th>…</th>
<th>波段100</th>
</tr>
</thead>
<tbody>
<tr>
<td>Pixel 1</td>
<td>0.01</td>
<td>-0.01</td>
<td>…</td>
<td>…</td>
<td>-0.003</td>
</tr>
<tr>
<td>Pixel 2</td>
<td>-0.01</td>
<td>0.01</td>
<td>…</td>
<td>…</td>
<td>0.017</td>
</tr>
<tr>
<td>Pixel 3</td>
<td>0.00</td>
<td>0.00</td>
<td>…</td>
<td>…</td>
<td>-0.013</td>
</tr>
</tbody>
</table>
<p>中心化消除了波段间的均值差异，让每个波段的数据均值为 0。<br>
<strong>这样做减去了数据的偏移，在更准确地反映数据变化方向的同时，也方便了下一步协方差矩阵计算。</strong></p>
<h2 id="22-计算协方差矩阵">2.2 计算协方差矩阵</h2>
<p>计算协方差矩阵之前，先简单摆一下协方差的定义，对于两个随机变量 <span class="math inline">\(X\)</span> 和 <span class="math inline">\(Y\)</span>，<strong>协方差定义为</strong>：</p>
<p></p><div class="math display">\[\text{Cov}(X,Y) = \frac{1}{m-1} \sum_{i=1}^{m} (X_i - \bar{X})(Y_i - \bar{Y})
\]</div><p></p><p>其中：</p>
<ul>
<li><span class="math inline">\(X_i, Y_i\)</span> 是第 <span class="math inline">\(i\)</span> 个样本的值。</li>
<li><span class="math inline">\(\bar{X}, \bar{Y}\)</span> 是样本均值。</li>
<li><span class="math inline">\(m\)</span> 是样本数。</li>
</ul>
<p>说人话，就是<strong>两个变量中心化后对应样本的乘积平均值</strong>。<br>
现在，再来看协方差矩阵，我们已经有了中心化后的数据矩阵：</p>
<p></p><div class="math display">\[\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}
\]</div><p></p><p>那么<strong>协方差矩阵公式</strong>就是：</p>
<p></p><div class="math display">\[\mathbf{C} = \frac{1}{m-1} \tilde X^\top \tilde X
\]</div><p></p><p>显然，得到的协方差矩阵 <span class="math inline">\(C\)</span> 是 <span class="math inline">\(n \times n\)</span> 的矩阵,其中的元素 <span class="math inline">\(C_{ij}\)</span> 就是波段 <span class="math inline">\(i\)</span> 与波段 <span class="math inline">\(j\)</span> 的协方差。<br>
还没完，我们再来专门看看<strong>协方差矩阵的对角线元素 <span class="math inline">\(C_{ii}\)</span></strong>，将其代入公式后，协方差公式就变成了：</p>
<p></p><div class="math display">\[ \frac{1}{m-1} \sum_{i=1}^{m} (X_i - \bar{X})^2
\]</div><p></p><p>显然，这就是方差公式，也就是说：<strong>协方差矩阵的对角线元素是每个特征自身的方差。而非对角的元素就可以看作”两个特征联合的波动程度“，也就是协方差。</strong><br>
到这里，我们就完成了 PCA 降维思想中，对方差的准备部分。</p>
<h2 id="23-特征值分解eigenvalue-decompositionevd">2.3 特征值分解（Eigenvalue Decomposition，EVD）</h2>
<p>在得到协方差矩阵 <span class="math inline">\(C\)</span> 后，下一步就是对其进行特征值分解，以提取数据中的主成分方向，这是 <strong>PCA 的最核心步骤</strong>。<br>
首先，在多维向量空间中，我们会定义单位向量 <span class="math inline">\(\mathbf{w}\)</span>，其长度固定为 1，目的是保证我们只选择<strong>方向</strong>，而不会因为向量长度不同而改变投影尺度。</p>
<p></p><div class="math display">\[|\mathbf{w}| = 1
\]</div><p></p><p>现在，先不提上一步的协方差矩阵，现在假定我们只有一个样本 <span class="math inline">\(\mathbf{x}\)</span> ，现在，我们要<strong>计算 <span class="math inline">\(\mathbf{x}\)</span></strong><br>
<strong>在某个方向上的投影</strong>，其公式就是：</p>
<p></p><div class="math display">\[z = \mathbf{w}^T \mathbf{x}
\]</div><p></p><p>得到的 <span class="math inline">\(z\)</span> 即为 <span class="math inline">\(\mathbf{x}\)</span> 在该方向上的投影数值，是一个标量。<br>
但我们要计算的是 <strong>投影后的方差</strong>，一个数怎么计算方差呢？这便涉及到线代里的结论：</p>
<p></p><div class="math display">\[{Var(z)} = \mathbf{w}^T C \mathbf{w}
\]</div><p></p><p>其中，<span class="math inline">\(C\)</span> 就是我们上一步得到的协方差矩阵。<br>
继续下一步就是我们的目标：<strong>找到投影后最大的方差</strong> 即：</p>
<p></p><div class="math display">\[\max_{|\mathbf{w}|=1} \mathbf{w}^T C \mathbf{w}
\]</div><p></p><p>现在，问题就转变成了一个求最值，也就是优化问题，在数学上，其解决方法是用<strong>拉格朗日法</strong>构造函数：</p>
<p></p><div class="math display">\[L(\mathbf{w},\lambda) = \mathbf{w}^T C \mathbf{w} - \lambda (\mathbf{w}^T \mathbf{w} - 1)
\]</div><p></p><p>然后对 <span class="math inline">\(w\)</span> 求导并让其为 0 得到：</p>
<p></p><div class="math display">\[C \mathbf{w} = \lambda \mathbf{w}
\]</div><p></p><p><strong>这个式子就是特征值方程</strong>，简单补充一下：<br>
给定一个<strong>方阵 A</strong>(行列数相等的矩阵)，假设我们有一个非零的向量 <span class="math inline">\(\mathbf{v}\)</span>，<strong>如果这个向量通过矩阵变换后只是按比例放大或缩小</strong>，即：</p>
<p></p><div class="math display">\[A \mathbf{v} = \lambda \mathbf{v}
\]</div><p></p><ul>
<li><span class="math inline">\(\mathbf{v}\)</span> 就是矩阵的<strong>特征向量</strong>。</li>
<li><span class="math inline">\(\lambda\)</span> 就是对应的<strong>特征值</strong>。</li>
</ul>
<p><strong>需要强调的是，这里虽然使用特征值方程，但特征值方程本身并没有“最大方差方向”的语义，只是我们的数学推导恰好得到了这个形式。</strong></p>
<p>继续，把特征值方程代回方差公式，就得到了：</p>
<p></p><div class="math display">\[{Var(z)} = \mathbf{w}^T C \mathbf{w} = \mathbf{w}^T (\lambda \mathbf{w})=\lambda \mathbf{w}^T \mathbf{w}= \lambda
\]</div><p></p><p>现在，我们就得出了<strong>最终结论：特征值 = 这个方向上的方差，而与特征值对应的特征向量的方向，就是我们要找的方差最大的方向。</strong></p>
<p>自此，我们就完成了 PCA 最核心的内容，要说明的是，PCA 的理论基础是协方差矩阵的特征值分解，但现代机器学习库通常使用另一种叫奇异值分解来完成计算，因为其省略了显式计算协方差矩阵的过程，但二者的结果是相同的。</p>
<h2 id="24-选择主成分">2.4 选择主成分</h2>
<p>根据线代的基本知识，我们知道一个 <span class="math inline">\(n \times n\)</span> 的方阵，其特征值方程最多有 <span class="math inline">\(n\)</span> 个解。<br>
而协方差矩阵是对称矩阵的同时也是半正定矩阵，因此其特征值方程一定有 <span class="math inline">\(n\)</span> 个正交的解。<br>
也就是说，在上一步，我们最终可以得到 <span class="math inline">\(n\)</span> 组特征值和对应的特征向量。<br>
到这里，我们就称这些特征向量所代表的方向为“<strong>成分</strong>”。<br>
<strong>特征值</strong> 的大小决定每个主成分的重要性。特征值越大，意味着该主成分方向承载的数据方差越多，信息量也越大。<br>
最大特征值对应的特征向量，就是第一主成分方向，第二大的特征值对应第二主成分方向，以此类推。</p>
<p>而我们在这一步的工作如下：</p>
<ol>
<li><strong>计算所有特征值之和 = 数据的总方差</strong></li>
<li><strong>计算每个主成分的方差占比并排序</strong></li>
<li><strong>一般挑选满足方差占比之和大于 95% 的前 <span class="math inline">\(k\)</span> 个主成分，最终降维到 <span class="math inline">\(k\)</span> 维。</strong></li>
</ol>
<p>来看一个简单的例子，假设我们对数据做完特征值分解后，得到如下结果：</p>
<table>
<thead>
<tr>
<th>排名</th>
<th>特征值 <span class="math inline">\(\lambda\)</span></th>
<th>解释的方差比例</th>
</tr>
</thead>
<tbody>
<tr>
<td>1</td>
<td>5.0</td>
<td>5.0 / (5.0+0.8+0.2)=83.3%</td>
</tr>
<tr>
<td>2</td>
<td>0.8</td>
<td>0.8 / (5.0+0.8+0.2)=13.3%</td>
</tr>
<tr>
<td>3</td>
<td>0.2</td>
<td>0.2 / (5.0+0.8+0.2)=3.3%</td>
</tr>
</tbody>
</table>
<p>你会发现：前两个主成分已经解释了约 96% 的信息<br>
这时我们就可以说：<strong>数据的主要结构，几乎都集中在前两个方向上。</strong></p>
<p>因此：如果我们把数据从 3 维降到 2 维，仍然可以保留 96% 的信息，丢失的信息只有 4%，这就是 PCA 降维的核心思想。<br>
到此，我们就只剩下最后一步了。</p>
<h2 id="25-构建降维后的数据">2.5 构建降维后的数据</h2>
<p>在得到重要的主成分方向后，最后就是将原始数据投影到这些主成分方向上，得到降维后的数据了。这个过程就是将原始数据映射到新的低维空间中。<br>
<strong>投影公式</strong>为：</p>
<p></p><div class="math display">\[Y = \tilde{X} \mathbf{V}_k
\]</div><p></p><p>其中：</p>
<ul>
<li><span class="math inline">\(Y\)</span> 是降维后的数据矩阵，尺寸为 <span class="math inline">\(m \times k\)</span>。</li>
<li><span class="math inline">\(\mathbf{V}_k\)</span> 是包含前 <span class="math inline">\(k\)</span> 个主成分的特征向量矩阵，尺寸为 <span class="math inline">\(n \times k\)</span>。</li>
<li><span class="math inline">\(Y\)</span> 每行对应一个样本，通过前 <span class="math inline">\(k\)</span> 个主成分表示。</li>
</ul>
<p>我们继续上一步的例子来演示，假定主成分对应的单位特征向量为：</p>
<p></p><div class="math display">\[\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}
\]</div><p></p><p>然后将两个特征向量组合为矩阵：</p>
<p></p><div class="math display">\[\mathbf{V}_k =
\begin{bmatrix}
0.58 & 0.71 \\
0.58 & -0.71 \\
0.58 & 0
\end{bmatrix}
\]</div><p></p><p>现在，假设中心化后的数据为：</p>
<p></p><div class="math display">\[\tilde{X} =
\begin{bmatrix}
2 & 1 & 1 \\
1 & 0 & 1 \\
3 & 1 & 2
\end{bmatrix}
\]</div><p></p><p>代入PCA 的投影公式即可得到：</p>
<p></p><div class="math display">\[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}
\]</div><p></p><p>计算结果为：</p>
<p></p><div class="math display">\[Y =
\begin{bmatrix}
2.32 & 0.71 \\
1.16 & 0.71 \\
3.48 & 1.42
\end{bmatrix}
\]</div><p></p><p>最终，我们就在尽可能保留信息的前提下，将数据从 3 维压缩到 2 维。</p>
<h1 id="3pca-的优缺">3.PCA 的优缺</h1>
<p>天底下没有完美的方法，PCA 自然也有其优劣所在。</p>
<h2 id="31-pca-的优势">3.1 PCA 的优势</h2>
<table>
<thead>
<tr>
<th>优点</th>
<th>说明</th>
</tr>
</thead>
<tbody>
<tr>
<td>计算效率高</td>
<td>PCA 主要依赖协方差矩阵的特征值分解，算法结构简单、无需复杂迭代，计算稳定且易实现</td>
</tr>
<tr>
<td>能去除特征冗余</td>
<td>通过协方差矩阵分析特征之间的相关性，将高度相关的特征压缩到少数几个主成分中</td>
</tr>
<tr>
<td>数据压缩能力强</td>
<td>在尽量保留主要方差信息的前提下降低数据维度，从而减少计算量和存储需求</td>
</tr>
<tr>
<td>主成分互不相关</td>
<td>PCA 得到的主成分两两正交，因此在统计意义上互不相关，有利于后续建模</td>
</tr>
</tbody>
</table>
<h2 id="32-pca-的不足">3.2 PCA 的不足</h2>
<table>
<thead>
<tr>
<th>缺点</th>
<th>说明</th>
</tr>
</thead>
<tbody>
<tr>
<td>只能处理线性结构</td>
<td>PCA 本质是线性投影方法，对于存在非线性流形结构的数据表达能力有限</td>
</tr>
<tr>
<td>对特征尺度敏感</td>
<td>协方差计算会受到特征量纲影响，因此通常需要先进行标准化处理</td>
</tr>
<tr>
<td>主成分可解释性较弱</td>
<td>每个主成分通常是多个原始特征的线性组合，难以直接解释其物理或业务含义</td>
</tr>
<tr>
<td>对异常值敏感</td>
<td>方差对极端值非常敏感，少量异常样本可能显著改变主成分方向</td>
</tr>
</tbody>
</table>
<p>总之，PCA 通过特征值分解和主成分选择，从高维数据中提取最重要的信息，极大地减少了计算复杂度。它不仅用于高光谱图像的降维，也广泛应用于其他数据分析和机器学习任务中。以传统的 PCA 为基础，后续也发展出了Kernel PCA，稀疏 PCA等改进方法。</p><br><br>
来源：https://www.cnblogs.com/Goblinscholar/p/19673621
頁:
[1]