简介面向高光谱图像分析与机器学习应用的MATLAB工具包集成主成分分析(PCA)、K近邻(KNN)与卷积神经网络(CNN)三种典型算法适用于遥感监测、农业勘察、地质探测等场景下的光谱分类、特征提取与异常检测任务。压缩包共376个文件、约15.8MB以MATLAB脚本(.m/.p)和函数为主体辅以JPG示例图像、FIG图形界面文件、HTML/HTM说明文档以及少量MAT数据与PDF结构完整可直接在MATLAB环境中加载运行并对照学习。目前已有946人学习下载适合具备一定MATLAB基础、希望快速上手高光谱处理链路的科研人员和学生也适合课程设计及算法对比实验参考。资源内含PCA降维预处理、KNN像素级光谱分类、CNN网络构建与训练等关键实现并配有图形界面及相关fig图表能够帮助读者理解从高维光谱数据到分类结果的完整处理流程也可为毕业设计、科研实验提供可扩展的代码基础。1. 为什么说高光谱图像处理工具包绕不开 PCA、CNN 和 KNN 的组合高光谱图像处理最让人头疼的不是数据量大而是几百个波段高度相关、标注样本稀少模型在验证集上一个样推到真实场景又一个样。这个工具包把 PCA 降维、CNN 特征提取和 KNN 分类串成一条流水线先用 PCA 砍掉波段间的冗余信息再用 CNN 在空谱联合窗口里抽特征最后用 KNN 做分类兜底。适合手里有高光谱立方体、想在 Matlab 里跑通“降维 分类”全流程的从业者也适合正在纠结 CNN 最后一层该用 Softmax 还是 KNN 的人。下面每一步都按可复现的标准来写参数也按工程经验给好。2. 光谱数据读入与预处理从 ENVI 文件到能喂给算法的三维张量2.1 读取高光谱数据的三条路径与 multibandread 参数先说读取。高光谱数据在工程里最常见的不是一份现成的 .mat 变量而是 ENVI 标准格式一个 .dat 二进制数据文件配一个 .hdr 文本头。.hdr 里记录了 lines、samples、bands、data type、byte order 和 band names。Matlab 读它不需要额外工具箱multibandread一个函数就能完成关键是头信息不能抄错。% 从 .hdr 文件里确认的尺寸与类型信息 lines 145; samples 145; bands 200; data multibandread(indian_pines.dat, ... [lines samples bands], ... uint16, 0, ... bsq, ieee-le);multibandread的五个核心参数从左到右依次是二进制文件名、[行 列 波段] 维度向量、数据类型、头文件字节偏移ENVI 标准格式一般为 0、波段排列方式bsq / bil / bip 三选一。这里 145×145×200 是 Indian Pines 的经典尺寸换成你自己的数据时把 lines、samples、bands 和 .hdr 对好即可少一个维度都会直接读崩。波段排列方式要单独说一句这是第一个坑。bsq 是逐波段整块存储bil 是逐行存所有波段bip 是逐像元存光谱曲线。如果原始文件是 bip而你按 bsq 读Matlab 不会报错但波段维度会错位后面 PCA 的协方差结构全乱。我一般会让脚本自动读 .hdr 里的interleave字段来定这个参数而不是手写死。第二常见的数据来源是 .mat 文件通常已经是 lines×samples×bands 的三维数组直接load就能用。第三类是无人机或卫星厂商导出的逐波段 TIFF 堆栈常见做法是循环读入再沿第三维拼接起来顺序不能乱很多国产无人机导出的波段顺序和波长列表是反的读完之后一定要跟波长相核对。2.2 光谱预处理四件套除吸水波段、归一化、降噪和样本筛选拿到原始张量后不能直接进 PCA先做四件事剔除吸水波段、归一化、降噪、样本筛选。吸水波段指 1350–1450 nm 和 1800–1950 nm 附近大气水汽吸收导致信噪比极低这些波段的光谱曲线基本是噪声留着只会干扰 PCA 的主成分方向。% 根据波长列表剔除 1350-1450 nm 与 1800-1950 nm 两个吸水区间 wavelength hdr.wavelength; % hdr 中读出的中心波长单位 nm keep ~((wavelength 1350 wavelength 1450) | ... (wavelength 1800 wavelength 1950)); data data(:, :, keep);这里用的是逻辑索引keep是跟波段数等长的逻辑向量保留位置为 true。如果你手里没有波长列表还有一种常见做法是直接观察各波段均值曲线找出剧烈下陷又回升的区间手动删掉。删除之后记得同步更新波长度量否则后面画特征曲线时横轴会错位。归一化的方式会直接影响 PCA 结果。% 逐波段做 min-max 归一化压缩到 [0,1] for b 1:size(data, 3) band data(:, :, b); bmin min(band(:)); bmax max(band(:)); data(:, :, b) (band - bmin) / (bmax - bmin eps); endeps是为了防止某个波段全为常数时分母为 0。注意这里不能把整个三维矩阵拉平后做一次全局 min-max因为各波段的动态范围差异很大全局缩放会抹掉低反射率波段的有效细节。有的传感器导出的 DN 值范围是 0–65535反射率要除以 10000 或 65535这个换算要在 min-max 之前完成。降噪方面轻度使用中值滤波或高斯滤波可以提升信噪比但空间滤波会改变光谱曲线形态做分类可以接受做光谱分析要慎重。样本筛选指的是剔除背景像元比如反射率全 0 或全 1 的像素、水体阴影区常见做法是用一个简单的阈值掩膜过滤掉否则这些非目标像元会把 PCA 的主成分方向带偏。3. PCA 降维用累计贡献率决定保留多少波段才不算浪费3.1 高光谱为什么必须降维波段相关性与 Hughes 效应高光谱数据几百个波段中相邻波段往往只差几个纳米相关性动辄 0.9 以上。把这么多相关特征直接喂给分类器不仅计算量大还会触发 Hughes 效应当训练样本数固定时特征维数越高分类器精度反而越低。这不是玄学而是因为高维空间中样本分布极其稀疏距离度量稳定性变差模型把小噪声也当成判别信号学进去了。PCA 的作用是把原始波段线性组合成一组正交的主成分第一主成分方向是数据方差最大的方向第二主成分在剩余方向上取方差最大以此类推。这样既去除了波段间的线性相关又把信息集中到前面少数几个分量上。结合 ENVI 主成分分析 PCA 的经验很多人习惯性看前三个分量做 RGB 合成图但做分类任务时这么选往往会丢掉大量可判别信息需要按累计贡献率来定。3.2 用 pca 函数做主成分分析的正确姿势Matlab 自带的pca函数可以直接用但它默认只做中心化、不做标准化而高光谱各波段的量纲和方差差异很大必须先自己做 z-score 标准化否则方差大的波段会主导主成分方向。% 高光谱立方体 145x145x200 - 样本矩阵 21025x200 [H, W, B] size(data); X reshape(data, H*W, B); % 用掩膜去掉背景像元只保留有效样本 valid all(X 0 X 1, 2); Xv X(valid, :); % 逐波段 z-score 标准化 Xz (Xv - mean(Xv)) ./ std(Xv); % PCA 主函数返回系数、得分、方差贡献率 [coeff, score, ~, ~, explained] pca(Xz); % 累计贡献率达到 99% 所需的主成分个数 cumvar cumsum(explained); nComp find(cumvar 99, 1);pca返回的coeff是 B×B 的系数矩阵每列是一个主成分在原始波段上的权重向量score是样本在主成分空间中的坐标尺寸是样本数×主成分数explained是每个主成分解释方差占总方差的百分比。find(cumvar 99, 1)取达到阈值的最小索引这就是后续喂给 CNN 的波段数。Indian Pines 这类数据通常 99% 贡献率落在 20–40 个主成分之间具体数值取决于原始数据的信噪比和预处理质量。这里有一个容易踩的细节pca函数默认对输入去做中心化所以 Xz 其实不需要再手动减去均值。但标准化这一步必须自己做因为pca不提供缩放功能。如果数据已经做过分波段归一化这里可以跳过 z-score但需要确认各波段方差在同一个量级。3.3 累计贡献率选 95% 还是 99%目视分析和分类任务的取舍在 ENVI 里做主成分分析很多人习惯看前三个分量合成假彩色图突出的是空间纹理和主要地物差异这对目视判读够用但直接拿来做分类通常精度不够。Matlab 里做分类层面的降维常见做法是把累计贡献率阈值设到 99%因为高光谱的判别信息往往藏在尾部贡献率只有 0.1% 的分量里。贡献率阈值与后续分类器有关。CNN 对噪声有一定鲁棒性95% 也能跑但 KNN 对噪声敏感尾部主成分里的噪声会被当作有效距离计算建议保留到 99%。我做对比实验时发现同一批数据95% 取到 15 个主成分、99% 取到 32 个主成分KNN 的总体精度差了 4–6 个百分点这个差距在矿物分类这种目标差异小的场景里尤其明显。主成分数量最好固定成一个参数放在脚本顶部而不是写死在后续代码里这样调起来会高效很多。另外降维后的score矩阵要 reshape 回三维立方体维度是 H×W×nComp后面 CNN 提取 patch 时直接在这个张量上切窗口。4. pcacnn 与 KNN 组合先让卷积网络提特征再用近邻分类兜底4.1 pcacnn 处理流程的搭建逻辑pcacnn 的核心思路不是把 PCA 和 CNN 简单拼在一起而是让 PCA 先做光谱去相关和降维再把每个像元连同空间邻域一起送入 CNN提取空谱联合特征。高光谱单像元光谱曲线包含物质判别信息但忽略了空间上下文CNN 的卷积窗口能同时看到光谱和空间这也是高光谱图像处理工具包里把 pcacnn 作为标准流程的原因。流程上分四步PCA 降维、构造 patch 数据集、训练 CNN 特征提取器、用 KNN 分类。注意这里 CNN 不是必须输出分类结果它可以只做特征提取器把最后一个全连接层前的激活值当成向量特征再交给 KNN 决策。这种设计的好处是当训练样本很少时CNN 的 softmax 层容易过拟合而 KNN 不做训练直接在特征空间中找最近邻泛化更稳。4.2 在 Matlab 里定义 CNN 特征提取网络先把 PCA 降维后的数据 reshape 回三维立方体再以每个已标注像元为中心取邻域 patch。% 将 PCA 降维结果还原为三维张量 data_pca reshape(score(:, 1:nComp), H, W, nComp); % 以每个训练样本为中心取 11x11 邻域 patchSize 11; numTrain length(trainLabels); patches zeros(patchSize, patchSize, nComp, numTrain); for i 1:numTrain r trainRows(i); c trainCols(i); r1 r - floor(patchSize/2); r2 r floor(patchSize/2); c1 c - floor(patchSize/2); c2 c floor(patchSize/2); patches(:, :, :, i) data_pca(r1:r2, c1:c2, :); end % 定义轻量 CNN 特征提取网络 layers [ imageInputLayer([patchSize patchSize nComp], Name, input) convolution2dLayer(3, 16, Padding, same, Name, conv1) batchNormalizationLayer(Name, bn1) reluLayer(Name, relu1) maxPooling2dLayer(2, Stride, 2, Name, pool1) convolution2dLayer(3, 32, Padding, same, Name, conv2) batchNormalizationLayer(Name, bn2) reluLayer(Name, relu2) fullyConnectedLayer(64, Name, fc_feat) softmaxLayer(Name, softmax) classificationLayer(Name, output) ]; options trainingOptions(adam, ... InitialLearnRate, 1e-3, ... MaxEpochs, 30, ... MiniBatchSize, 128, ... Verbose, false); net trainNetwork(patches, categorical(trainLabels), layers, options);网络设计上有几个参数值得展开。输入层第三个维度必须是 nComp即 PCA 保留的主成分数不是原始波段数。第一个卷积层用 3×3 卷积核、16 个滤波器是为了在小 patch 上提取局部光谱和空间纹理后面加 batchNormalization 是因为高光谱特征尺度不稳定BN 能保证每层输入分布相对稳定。第二个卷积层把滤波器加到 32最后用 64 维全连接层fc_feat作为输出特征维度。patchSize取 11 是个折中太小看不到空间纹理太大会把周围不同类别的像元卷进来反而干扰特征。对空间分辨率 3–10 m 的高光谱影像5–15 都有人在用常见做法是先试 11再根据混淆矩阵调。4.3 用 KNN 替换 Softmax 做最终分类训练完成后把fc_feat层的输出提取出来作为每个 patch 的特征向量。% 从训练好的网络提取 64 维特征向量 featTrain activations(net, patches, fc_feat); featTrain squeeze(featTrain); % 转成 样本数 x 64 % 训练 KNN 分类器 mdl fitcknn(featTrain, double(trainLabels), ... NumNeighbors, 5, ... Standardize, true, ... Distance, euclidean); % 测试集同样提取特征后预测 featTest activations(net, testPatches, fc_feat); featTest squeeze(featTest); pred predict(mdl, featTest); acc sum(pred testLabels) / numel(testLabels);activations的第二个参数是网络层名称第三个参数是要提取特征的层名返回的维度是 1×1×64×N所以要squeeze后转置。fitcknn里NumNeighbors控制近邻数量默认 5 在小样本场景表现稳定Standardize必须打开因为 CNN 的 64 维特征各维度尺度差异可能很大不标准化时距离会被少数维度主导。为什么最后不用 CNN 自带的 softmax 层直接预测小样本场景下CNN 的 softmax 输出容易对训练集过拟合而 KNN 是惰性分类器不学习分类边界直接在特征空间中度量距离。特征本身由 CNN 的深度结构保证判别力KNN 只负责决策两者互补。实际对比中训练样本每类少于 50 个时KNN 方案比 softmax 方案高 2–5 个百分点样本充足时两者接近KNN 的优势就不明显了。5. 高光谱处理避坑指南样本泄漏、波段尺度和 CNN 输入尺寸5.1 现象PCA 之后分类精度反而下降有人在跑 PCA 时先对整幅影像做了标准化再把标准化后的结果划分训练集和测试集。这样做出来的主成分方向其实已经“偷看”了测试集信息训练精度看着不错推到新影像上立刻翻车。原因在于 PCA 属于无监督变换它的系数应当只从训练集统计量中求得。正确做法是先划分训练集和测试集只对训练集计算均值和标准差、求 PCA 系数再用同一套系数变换测试集。% 只对训练集做 PCA测试集直接乘系数 mu_train mean(Xtrain); std_train std(Xtrain); Xtrain_z (Xtrain - mu_train) ./ std_train; [coeff, score_train] pca(Xtrain_z); Xtest_z (Xtest - mu_train) ./ std_train; score_test Xtest_z * coeff;5.2 现象验证集精度虚高推到真实场景就翻车patch 窗口存在重叠时训练像元和测试像元的邻域都包含同一片空间空间相关性被当成可判别信号学进去了。这种现象在验证集精度上表现得很“好看”但在独立场景验证时精度暴跌。解决思路是阻止 patch 重叠。常见做法是按区域划分把影像分成互不重叠的若干块整块分配给训练或测试或者保证所有标注像元两两之间距离大于 patchSize这样邻域窗口不会交叉。Indian Pines 这类小场景里样本本来就少若能找到另一张同时相的影像做独立验证比任何数据划分都可信。5.3 现象KNN 距离被个别波段主导打开Standardize后仍然出现某类样本全部被分错多半是特征空间里个别维度方差过大且该维度携带的是噪声而非区分信息。用 CNN 提取的 64 维特征通常比较均衡但如果是从 PCA 得分直接喂 KNN前几个主成分方差占比极大距离计算被它们主导。解决方法是先对特征矩阵再做一次 z-score或者用马氏距离替代欧氏距离。fitcknn里把Distance设为mahalanobis也能缓解尺度问题但计算量大不少。工程上我一般保留欧氏距离手动做一次标准化。5.4 现象CNN 输入 patch 越界报错维度不匹配边界像元取邻域时r1或c1小于 1Matlab 直接报索引越界。最常见做法是边缘复制把边界外的像元用最近的像元值填充。% 用复制边缘方式处理边界越界 padData padarray(data_pca, [floor(patchSize/2) floor(patchSize/2)], replicate, both); % 后续从 padData 中取 patch索引统一加 padSize注意padarray的填充方式还有symmetric和circular高光谱场景下replicate最稳因为边界外的像元大概率与最近像元同类。填充之后再取 patch 就不会越界代价是影像外围的区域分类结果不参与评估。5.5 现象吸水波段没剔除特征图出现条带噪声如果跳过吸水波段剔除直接归一化吸水波段的反射率本来就接近 0min-max 归一化会把微小的传感器噪声放大成整个动态范围PCA 里这些波段会形成独立的噪声主成分CNN 特征图对应位置出现横向条带。解决方法是把吸水波段删除放在归一化和 PCA 之前顺序不能颠倒。如果数据已经走完归一化才发现这问题唯一靠谱的后悔药是回到原始数据重新处理没有别的捷径。参数速查表环节推荐参数说明patch 大小11×11小场景 5–7大场景可到 15PCA 累计贡献率99%分类任务别按 ENVI 目视习惯选 95%第一卷积层3×316 滤波器先轻后重防止小样本过拟合特征维度64过小丢信息过大 KNN 距离不稳KNN K 值5用奇数样本多可试 7–9避免偶数平票6. 验证与调参用混淆矩阵和 K 值敏感性确认结果不是玄学6.1 用混淆矩阵与 Kappa 系数验证分类结果总体精度一个数说明不了问题地物类别不均衡时尤其如此。算一下混淆矩阵和 Kappa 系数才能判断精度高是因为真的分对了还是因为把多数类全猜对了。cm confusionmat(testLabels, pred); accuracy sum(diag(cm)) / sum(cm(:)); % Kappa 系数计算 po accuracy; pe sum(sum(cm, 2) .* sum(cm, 1)) / (sum(cm(:))^2); kappa (po - pe) / (1 - pe);Kappa 大于 0.8 认为分类结果可靠0.6–0.8 属于中等一致低于 0.6 说明流程里大概率存在前面提到的某个坑。这时不要急着调网络先看混淆矩阵对角线之外哪些类别互相混淆通常相邻地物类型比如草地和牧场容易混这时可以回到 PCA 保留数或 K 值上找原因。6.2 K 值敏感性与特征图检查调参时最值得做的一个实验是 K 值敏感性扫描。for k 1:2:15 mdl fitcknn(featTrain, double(trainLabels), ... NumNeighbors, k, Standardize, true); pred predict(mdl, featTest); acc(k) sum(pred testLabels) / numel(testLabels); end正常结果是一条先上升然后平台或缓慢下降的曲线。如果精度随 K 值剧烈抖动说明特征空间里同类样本并不聚集问题通常出在 PCA 保留数太少或 CNN 特征提取不足。我习惯把这一步当成每次实验的标配几分钟跑完能拦住一大批无效调参。特征图检查也不可跳过。用activations提取第一个卷积层的输出可视化前几个通道如果激活值集中在目标地物区域而不是随机散点说明网络学到的是空间散射特征可以继续调参如果全是噪声先回查预处理。这套链路走下来真正花时间的往往不是网络本身而是数据准备和验证设计。我自己的教训是第一次跑通时验证集精度很高兴奋地直接换场景推理结果惨不忍睹最后发现是 patch 重叠泄漏。从那以后每批数据都强制做独立区域划分并把 PCA 系数、K 值、patch 大小这些关键参数写成脚本头部的配置项跑实验清晰很多。希望这些经验能帮你少走一趟弯路。本文还有配套的精品资源点击获取