基于分数阶微分的土壤重金属高光谱遥感图像反演

丁松滔，张霞，尚坤，李儒，孙伟超

1.中国科学院空天信息创新研究院,北京 100101;2.中国科学院大学,北京 100049;3.自然资源部国土卫星遥感应用中心,北京 100048

1 引言

成像高光谱遥感能够获取图像上连续且精细的土壤反射光谱，具有低成本大范围快速监测土壤重金属含量的潜力（Khosravi 等，2018）。然而目前大多数研究采用实验室理想条件下采集的光谱，基于高光谱图像的重金属反演研究较少（Ding 等，2022），并且由于土壤样本采集受制于地形、交通可达性等因素，且土壤样本的重金属实验室化验分析成本高，导致样本不充足的问题普遍存在，而对于高光谱图像而言，图像上的土壤像元数量与土壤样本数相差悬殊，小样本问题更加突出。光谱特征选择是解决小样本问题的有效途径之一，而服务于光谱特征选择的特征增强方法也是研究热点，研究显示经过特征增强的光谱曲线比原始光谱曲线可以更有效地提取出光谱的特征波段（沈强 等，2019）。

微分处理是常用的增强光谱特征的方法，一阶微分和二阶微分是最常用的光谱微分方式，能够突出土壤光谱中与重金属相关的光谱信息，提高反演精度（郭学飞 等，2020）。然而，传统的整数阶微分（即一阶、二阶微分）缺乏对可能包含土壤重金属有益信息的渐变倾斜度或曲率的敏感性（Hong等，2018）。分数阶微分作为整数阶微分的扩展，对光谱微分的处理更加细致，阶数上更精细的量化意味着能更敏感地捕捉各波段的曲率和倾斜度变化，突出光谱特征。已有学者将FOD应用于实验室光谱（Chen 等，2022），成功地反演了土壤重金属Cr、Pb、Zn的含量；Meng等（2021）将FOD应用于GF-5图像，使用FOD处理后的像元光谱计算光谱指数作为自变量，构建的有机质反演模型精度R2达到0.84，但FOD 应用于高光谱图像上进行土壤重金属反演的有效性还有待验证。

特征选择算法通过选出特征性较强的波段，提高模型的反演性能。现有很多研究都已经证明，使用优选出的特征波段比使用全波段可以获得更好的反演结果（Tan 等，2020；Zhang 等，2021），这表明波段选择的重要性。土壤重金属反演中最为常见的筛选方法是皮尔森相关系数法，但该方法只考虑单个波段与理化性质间的关系（周冰 等，2021），没有考虑多个波段间的共线性。GA 是一种以“适者生存”为原则的优化算法，以最大化精度为目标，在特征空间中进行启发式搜索，可以快速获取近似最优解（柏晗 等，2022），已有将其应用于遥感图像上反演重金属含量的先例（Wang等，2022），然而GA的缺点是会过早收敛，导致结果难以跳出局部解（Pavez-Lazo 和Soto-Cartes，2011）。CARS 同样是遵循“适者生存”的一种筛选算法（Li 等，2009），通过逐步去除冗余和不重要变量来选择信息，使模型计算效率与稳定性提升（Vohland 等，2017），能够从全波段数据中选择最优的波段组合（Vohland 等，2016），Cheng 等（2021）在土壤全氮反演研究中指出CARS 算法在有效变量的选择上优于GA。综上所述，CARS 算法在这3 种算法中有利于筛选出更有效的波段组合，达到全局最优解。

本文以新疆维吾尔自治区哈密市的黄山南矿区为研究区，应用航空高光谱图像开展对土壤重金属Pb、Zn、Ni 含量的反演研究，提出了一种基于FOD 估算高光谱图像上土壤重金属浓度的反演方法，通过扩充样本、FOD 增强光谱特征、CARS筛选波段选取优光谱特征，探索高光谱图像反演中小样本问题的解决途径，分析确定了研究区3种重金属反演的最佳FOD 阶数，并将CARS 与CC、GA 的建模结果进行对比；进一步分析了CARS 的波段筛选结果，得出对3种重金属反演贡献最高的波段范围及其物理意义，分析了重金属含量分布图的可靠性。

2 研究区与研究数据

2.1 研究区概况

研究区位置如图1所示，位于新疆维吾尔自治区哈密市黄山南铜镍矿区，是重要的铜、镍、铁、铅、锌等大型矿床集中区（王京彬 等，2006），该矿区于2007年开始勘探，至今仍在开采中。地处新疆维吾尔自治区东部，属于典型的温带大陆性干旱气候，干燥少雨，年降水量33.8 mm，年蒸发量3300 mm，夏季酷热、蒸发强。研究区内主要为空旷戈壁，无植被覆盖，土壤表面多为裸露状态。

图1 研究区位置示意图Fig.1 The geographical location of the study area

2.2 数据获取及预处理

航空高光谱图像获取于2021 年8 月23 号，传感器为中国科学院上海技术物理研究所研制的全谱段多模态成像光谱仪，本次航飞航高3 km，波长范围覆盖350—2500 nm，350—1000 nm 光谱区间的空间分辨率为0.5 m，有251 个波段，1000—2500 nm 波长范围空间分辨率为1 m，有508 个波段。研究区由两条航带所覆盖，图像经过辐射定标、大气校正后得到地表反射率数据，对两航带图像进行配准、拼接等处理后，裁剪出研究区，经纬度范围为：42.197°N—42.216°N，94.625°E—94.684°E。

成像光谱仪在获取光谱数据时，大气中的水汽对1400 nm 和1900 nm 波长附近的辐射能量存在强吸收，且1900—2500nm 的像元光谱存在较严重噪声。为去除水汽等噪声的影响，并保留尽可能多的波段，剔除1350—1450 nm 和1800 nm 之后的波段。对剩余的光谱进行SG 滤波去噪处理，由于不同区间的光谱噪声不同，对其采用不同的窗口大小进行滤波，350—1000 nm 区间的噪声较弱，采用窗口大小为9 的二次多项式，1000—1800 nm区间为中度噪声，采用窗口大小为13 的二次多项式。

与航空飞行准同步，2021年7月22日—7月30日开展了土壤样本地面获取实验，由于研究区属于干旱区，且该时段内土壤未受降雨影响，两者获取数据时土壤状况可认为一致。土壤样本采集主要针对采矿区及周边，沿道路两旁选择土壤的颜色和粒径有明显差异的样点，采集0—20 cm 表层土，共采集土壤样本72 个。图像及土壤样本点位置见图2，根据土壤样本的经纬度信息，提取对应的图像像元反射光谱数据。

图2 航空遥感图像及采样点位置Fig.2 Aerial remote sensing image and sampling locations

将土样在实验室风干研磨后过100目筛制成标准样，用火焰原子吸收分光光度法测定重金属含量。图3 所示为样本中3 种重金属含量的直方图，可见重金属Pb 和Zn 的含量分布都趋向于正态分布，而重金属Ni 含量呈现出明显的偏分布，本文采用对数变换对其进行校正，转换后Ni 含量分布如图4所示，偏分布现象得到明显改善。

图3 土壤样本的重金属含量直方图Fig.3 Histogram of the heavy metal content of soil samples

图4 对数变化后Ni含量直方图Fig.4 Histogram of the Ni content after log change

3 实验方法

图5为本文研究技术路线图，首先根据土壤样本的经纬度提取对应像元和相似的邻近像元光谱，对像元光谱数据进行SG 滤波及分数阶微分处理，处理后的光谱采用CARS 算法进行特征波段优选，选出的波段组合用于偏最小二乘回归PLSR（Partial Least Squares Regression）建立反演模型，最终将构建的最优模型应用于航空高光谱图像上反演重金属含量，获得重金属含量分布图。

图5 技术路线图Fig.5 Flowchart of the proposed method

3.1 样本扩充

由于研究区的地形较平坦，地表一致性较好，可以认为所采集土壤样本的重金属含量可以代表该样本小范围内的土壤重金属含量，又因为本研究所获取的高光谱图像空间分辨率为1 m，可认为邻近像元所对应的土壤中重金属含量差异极小。因此，本研究根据样本经纬度信息定位到图像上的像元后，以该像元为中心像元，提取中心像元以及其八邻域的像元光谱，依次计算八个邻近像元光谱与中心像元光谱的欧氏距离，距离越近则认为其与中心像元的相似度越高，选出欧氏距离最小的两个像元，将3个像元（中心像元和两个相似邻近像元）的光谱共同作为该土壤样本的反射率光谱信息，达到将土壤样本的光谱数据扩大3 倍的效果。本研究采集的72 个土壤样本，用此方法进行样本扩充后，获得216条像元反射率光谱数据。

3.2 分数阶微分（FOD）

分数阶微分是数学中的重要分支，它将经典的整数阶微分推广为任意阶，能更敏感地捕捉光谱反射率细节的变化。到目前为止，还没有一个统一的公式来定义分数阶微分。数学家们从不同的角度分析，得出了分数阶微积分的不同定义。目前，分数阶微分常见的表达形式主要包括Riemann-Liouville（R-L）、Grunwald-Letnikov（G-L）和Caputo。其中最常用的形式为G-L（Wang等，2017），本研究采用G-L形式进行微分，其定义如下：

式中，α为任意阶数；h为微分步长；t与α分别为微分的上、下限；Γ(α)为Gamma 函数，Γ(α)满足：

令h=1，式（1）能够推导出函数f(x)分数阶微分的表达式为

3.3 竞争自适应重加权采样法（CARS）

CARS 是一种以“适者生存”为原则的筛选算法，基于与模型性能相关的统计数据，通过逐步去除冗余和不重要的变量来选择信息量大的变量，将其应用于高光谱数据进行筛选时，其本质是使用自适应加权采样技术来保留回归中绝对系数较大的光谱波段（Tan等，2021）。

图6 为CARS 的方法流程图，CARS 通过N次蒙特卡罗抽样迭代地选择N个波段子集，最终目的是选择一个能使交叉验证的均方根误差（RMSEcv）最低的最佳变量子集。在高光谱数据上使用CARS筛选的具体步骤如下：

图6 CARS算法流程图Fig.6 Flowchart of CARS

（1）使用蒙特卡罗抽样方法从数据集中随机选择具有固定比例的样本，然后用这些样本建立PLSR 模型。PLSR 模型可用公式表达为y=bX+e，e为常数项，b为回归系数向量，b=［b1，b2，…，bn］，b中第i个元素的绝对值|b|（1≤i≤n）表示第i个波段对y的贡献，波段对目标变量的贡献越大，该波段就越重要。为评价每个波长的重要性，定义权重Wi：

通过CARS 算法去掉的变量，其权重Wi均设为0；

（2）利用指数递减函数EDF（Exponentially Decreasing Function）和自适应重加权采样ARS（Adaptive Reweighted Sampling）分别强行和竞争性地消除权重低的波段；

（3）重复步骤（1）—（2），直到达到N次采样运行，最后选择RMSEcv 最低的波段子集作为最优波段组合。

由于蒙特卡洛抽样每次迭代抽取固定比例的样本参与运算，不使用所有样本，因此选出的波段组合具有更好的适应性；步骤（2）中使用指数递减函数去除波段，能在迭代前期筛选掉大量重要性低的波段，使筛选过程的运算量显著下降。根据多次实验结果，将蒙特卡罗抽样次数设置为50次，每次迭代抽取90%的样本用于运算。

3.4 模型构建与评价指标

PLSR 是土壤光谱分析中最为常用的反演模型，在处理共线性强、计算复杂度高的问题上具有天然优势（Rossel 和Behrens，2010），能够应对自变量多于样本个数的情况，因此其非常适合处理高光谱数据，并且已被广泛的应用于土壤重金属反演。本文采用CARS+PLSR 的模型构建方法，先通过CARS筛选出特征波段，再将特征波段输入PLSR 建立反演模型。在PLSR 建模中，通过留一交叉验证的最小均方根误差确定PLSR 的潜变量个数。由于CARS的蒙特卡罗抽样具有随机性，为了避免实验的偶然性影响，本文对每组反演实验均进行5 次，以5 次反演实验中的精度最优值作为该模型精度将精度最优的模型应用于高光谱图像，得到土壤重金属的含量分布图。

精度评定采用预测均方根误差RMSEP（Root Mean Square Error of Prediction）、相对分析误差RPD（Ratio of Prediction to Deviation）和决定系数R2（Coefficient of Determination）3 个评价指标。3 个评价指标，RMSEP 值越小，RPD 值越大，R2值越接近1，说明反演模型的精度越高；反之，RMSEP值越大，RPD 值越小，R2值越小，说明反演模型的精度越低，本文模型优劣参考现有的土壤属性含量高光谱反演的评价标准（Wang 等，2014）：出色模型，R2>0.9；良好模型，0.9>R2>0.8；近似模型，0.8>R2>0.65；具有一定反演能力，0.65>R2>0.50；不具备反演能力，0.50>R2。

4 实验结果与分析

反演中将216个样本按照2∶1的比例划分为训练集和测试集。由于样本扩充出的3个样本重金属含量相同，如果将样本重金属含量排序后采用分层抽样方法，按2∶1进行划分，会使训练集和测试集的样本含量分布完全一样，从而导致反演精度偏高。为了解决该问题，选出具有代表性的训练集，采用Kennard-Stone（KS）算法（Zhang 等，2017）来划分样本集，具体的实施步骤如下：

（1）首先计算两两样本之间的光谱欧氏距离，选择光谱距离最大的两个样本；

（2）然后分别计算剩余样本与已选两样本之间的光谱欧氏距离；

（3）对于每个剩余样本而言，计算其与已选各样本之间的最短光谱距离，选择这些最短光谱距离中相对最大的距离所对应的样本，作为新入选的样本；

（4）重复步骤（3），直至所选样本的个数等于事先设定的数目为止。

通过KS 算法选出144 个样本作为训练集，剩下的72 个样本作为测试集，样本集划分结果通过直方图展示于图7，3 种重金属的训练集、测试集和总样本集的含量直方图分布一致且有一定差异性，表明该方法划分出的训练集和测试集具有较好的代表性。

图7 各样本集的重金属含量直方图Fig.7 Histogram of heavy metal content for each sample set

4.1 分数阶微分的最佳阶数确定

由于分数阶微分需使用相邻的波段进行计算，理想条件下的土壤反射率光谱应是一条连续的光滑曲线，相邻波段的反射率不应出现较大变化。本研究使用的航空高光谱图像采集350—1000 nm和1000—2500 nm 波长范围的传感器不同，导致光谱在1000 nm 波长处的反射率曲线不够平滑。此外，由于1350—1450 nm 的水汽吸收带被剔除后，导致该光谱区间两端波长相差100 nm 的两个波段变为相邻波段。为了确保微分时相邻波段反射率的差异在合理范围内，将光谱分为350—1000 nm、1000—1350 nm 和1450—1800 nm 这3 个区间来分别进行分数阶微分处理。

从图8中可以看出，随着阶数的增加，光谱的特征变的越来越明显，峰谷差异越发增大；0—1阶的微分结果显示出，阶数增大的过程中，平缓区间的波段反射率值越发趋近于0，而波峰波谷区间的波段反射率值被逐渐放大；1—2 阶的微分结果中，反射率曲线的波形已较为相似（由于1—2阶内的微分曲线较为相似，图8中仅展示1阶、1.5阶和2 阶的微分结果），均与原始光谱曲线差异较大，特征放大的效果明显，随着阶数的增加，变化主要体现光谱曲线的极值不断增大。

图8 不同阶数微分后的像元光谱反射曲线Fig.8 Image spectral curves with different orders of differential

在不同阶数的微分处理后，反射率曲线中的光谱信息有明显差异，为了选出适合于反演各重金属的最佳微分阶数，在0—2阶的范围内，以0.1阶为间隔，将不同阶数微分后的光谱输入CARS+PLSR构建反演模型，以测试集精度为选择依据。为了避免实验的偶然性影响，确保选出的阶数能使反演精度达到最优，每组阶数下均进行了5次反演建模。

图9显示了在不同阶数微分条件下，CARS+PLSR反演3种重金属的测试集R2最大值和平均值。对于Pb，阶数为1.2 时测试集精度的最大值和平均值都达到了峰值（分别为0.7974 和0.7431）；对于Zn，当阶数为0.8 时，R2最大值曲线达到峰值0.8690，此时平均值曲线的值为0.8096，与平均值曲线的峰值0.8126 相差极小，因此将阶数0.8 认为是反演重金属Zn 的最佳微分阶数；对于Ni，当阶数为0.3 时最大值曲线和平均值曲线均达到峰值（分别为0.8303 和0.7681）。同时考虑精度的最优值和平均值，保证选出的最佳阶数能实现反演精度最优且模型稳定性好。

图9 不同阶数下反演3种重金属的测试集精度Fig.9 Estimation accuracy of three heavy metals at different orders

每种重金属反演时光谱区间内对不同种类重金属敏感有效的波段不同，即各重金属的特征波段不同，不同阶数的微分处理后对各区间波段的特征化效果有所差异，因此各重金属的最佳微分阶数不尽相同。图9 中当阶数为0 时，代表不做微分处理使用原始像元光谱，当阶数为1 和2 时相当于使用一阶微分和二阶微分，3种重金属的反演结果都显示出使用最佳阶数分数阶微分的精度高于原始像元光谱、一阶微分和二阶微分光谱的反演精度，证明了使用分数阶微分能够有效增强光谱特征，提高重金属反演精度，并且分数阶微分比整数阶微分能更敏感地突出对土壤重金属反演有益的光谱信息。根据本节的实验结果，本文在反演重金属Pb、Zn、Ni的各项实验中，分别选定1.2、0.8和0.3阶的分数阶微分对像元光谱进行处理。

4.2 样本扩充有效性分析

为验证样本扩充的有效性，每种重金属均进行两组反演实验，一组对样本进行扩充，使用216个样本参与建模，另一组不扩充样本，使用72 个样本，反演每种重金属时，都以其最佳阶数对像元光谱进行微分处理后，建立CARS+PLSR模型。

表1 中展示了各重金属的两组对照实验结果，在没有进行样本扩充前，小样本问题引起的过拟合现象明显，3 种重金属的反演模型训练集精度R2都大于0.98，然而只有Zn 的测试集精度为0.8178，Pb和Ni的测试集精度都较低，R2小于0.7；对样本集进行扩充后建立的反演模型，过拟合现象得到了很好的缓解，3种重金属的训练集精度和测试集精度的R2差距小于0.05。样本扩充后测试集的精度都得到了明显的提升，Pb 的测试集R2从0.6128提升到0.7974，Zn 的精度从0.8178 提升到0.8690，Ni 的R2从0.6969 提升到0.8303。由此可见，样本扩充不仅缓和了模型的过拟合现象，还有效地提升了3种重金属的反演精度。

表1 样本扩充前后的3种重金属反演精度（ 和 分别代表训练集和测试集的反演精度R2）Table 1 The estimation accuracies of three heavy metals before and after sample expansion（ and represent the inversion accuracies R2 of the training and test sets，respectively）

图10 展示了样本扩充前后重金属Ni 的反演结果散点图，样本中Ni含量大于2000 mg/kg的高含量样本仅有两个，从图10（a）（b）可以看出，两个高含量样本均被选入训练集，测试集中未包含高含量样本，难以通过测试集精度评估该模型应用于图像后在高含量区域的反演能力。从图10（c）（d）可以看出，样本扩充后训练集和测试集都包括了高含量的样本，因此测试集精度可以代表在高含量区域的反演效果。在使用3.1 节所述的样本扩充条件下，使用KS算法划分样本集能够确保高含量样本存在于训练集和测试集，使本研究方法在反演含量偏分布的重金属时，能够得到可靠性较高的反演模型。

图10 样本扩充前后重金属Ni反演散点图Fig.10 Scatter plot of Ni estimation results before and after sample expansion

4.3 波段选择算法的建模精度对比分析

GA 在PLSR 建模中被认为是一种有效的波段选择算法（Leardi 和González，1998）。此外，相关系数法在土壤重金属反演中被广泛使用，因此本研究将GA+PLSR、CC+PLSR 两种建模方法与CARS+PLSR 进行对比分析。参考已有的GA+PLSR研究（Sun 等，2022），GA 的参数设置为染色体个数20，迭代次数150，代际间隙90%，基因变异概率10%。CC 的相关性筛选阈值通过在0—1 范围内，以0.1 为间隔设置阈值进行3 种重金属的反演实验，以CC+PLSR 的反演精度最优为原则设置阈值的参数。

通过分数阶微分对扩充后的216条像元光谱进行处理后，分别用CC、GA 和CARS 等3 种算法筛选波段，选出的3 组波段组合分别构建PLSR 反演模型。

表2 展示了3 种重金属在使用不同波段选择算法下的反演精度，表中的T（s）代表单次模型构建所花费的秒数，分别代表训练集和测试集的反演精度R2。3种建模方法中，CC+PLSR 构建的模型测试集精度最低，甚至在反演Pb 时R2<0.5，模型没有估算能力，并且在3种重金属反演中均出现过拟合现象，反演Pb 和Zn 时过拟合最为显著，训练集精度远高于测试集精度。

表2 不同波段选择算法的重金属反演精度Table 2 Heavy metal estimation accuracy of different band selection algorithms

GA+PLSR 的模型精度较高，3 种重金属反演结果仅略低于CARS+PLSR 模型，反演Zn时的测试集精度R2为0.8119，属于良好模型，反演Pb 和Ni时测试集精度R2也均大于0.78，模型没有显示出过拟合现象。但由于GA 在筛选波段时，需要多次迭代选择最优波段组合，每次迭代中所有波段均参与运算，并且由于GA 算法中有种群机制的设置，在本研究的参数设置下，每次迭代需生成20种波段组合，意味着有20倍的总波段数参与运算，导致整个模型的构建时间较长，3种重金属的GA+PLSR反演模型构建时间都在300 s 以上，耗费时间远远多于CARS和CC两种筛选方法。

CARS+PLSR 在3 种建模方法中获得了最高的反演精度，Zn 和Ni 的反演精度R2都高于0.8，Pb的反演精度R2也达到0.7974，并且模型没有呈现出明显的过拟合现象，训练集和测试集都保持着较高的精度。此外，虽然CARS 和GA 算法都需要多次迭代寻找最优波段组合，但由于CARS 内的EDF 和ARS 算法在迭代初期就快速、强力的去除重要性低的波段，因此迭代中后期仅有少数波段参与运算，所需的建模时间显著缩短。因此CARS是3种方法中最优的波段选择算法，能够更快筛选出最优波段组合，使模型具有更好的反演能力。

4.4 CARS的波段选择结果分析

图11 展示了CARS 选出的最优波段组合的直方图，直方图内柱形图越高代表该波长范围内选择的波段越多，3 种重金属的CARS 筛选结果在不同波段区间选择的波段数量不同，说明各波长范围对不同重金属反演的贡献不同。土壤反射光谱分析表明土壤在VNIR-SWIR 区间的吸收特征主要由土壤有机质、铁氧化物等土壤光谱活性物质引起（Kooistra 等，2003）。由于土壤光谱活性物质在土壤中对重金属的吸附具有主导作用，可据此推算土壤重金属浓度，这是通过反射光谱间接反演重金属的主要机理（Rathod 等，2013）。然而有机质、铁氧化物等土壤光谱活性物质对于不同的重金属吸附强度不同，因此反演不同重金属时起主导作用的土壤组分不同（Covelo 等，2007），现有研究表明土壤铁氧化物对重金属Pb 的吸附作用是土壤反射光谱反演Pb 含量的主要机理，土壤有机质对重金属Ni、Zn 的吸附作用是土壤反射光谱反演Ni、Zn 的主要机理，铁氧化物的吸收特征在500 nm（Wu 等，2007），而600—800 nm 附近的吸收峰被认为是土壤有机质的吸收特征（徐彬彬 等，1991）。已有研究证明500 nm 的铁氧化物特征波段对反演重金属Pb 是有效的（张霞 等，2022），Sun 等（Sun 和Zhang，2017；Sun 等，2018）证明了使用600—800 nm的有机质特征谱段对反演Zn和Ni是有效的。根据图11所示，在反演Pb时，CARS筛选结果的直方图中500 nm 范围内的柱形图最高；在反演Zn 和Ni 时，直方图内最高的柱形图落在600—800 nm区间和1600 nm范围内，Li等（2022）指出有机质的官能团对1600 nm 光谱影响显著，证明1600 nm 同样是对有机质敏感的特征波段。综上所述CARS的波段选择结果与已有反演机理研究保持一致，证明CARS能够筛选出对反演重金属有益且合理的特征波段。

图11 3种重金属反演模型的CARS优选波段结果Fig.11 Results of CARS band selection for three heavy metal estimation models

4.5 重金属含量分布图分析

为了保证模型应用到高光谱图像上时反演结果的可靠性，通过光谱角匹配法（Ramirez-Lopez 等，2013）筛选出图像上与土壤样本像元光谱相似度高的区域，将3种重金属的最优模型分别应用于该区域进行重金属含量反演，图12（b）（c）（d）展示了依据此方法制作出的3 种重金属含量分布图。由于Ni 的含量区间跨度较大，连续的图例难以体现出含量差距，图12（d）的Ni含量分布图中采用了分级图例进行展示，根据2018 年发布的《土壤环境质量 建设用地土壤污染风险管控标准》（GB 36600-2018）相关规定，矿区作为二类用地，Ni 含量的筛选值（可能存在风险）为900 mg/kg，因此设置900—2500 mg/kg 的级别代表超过筛选值，用于表示对人体健康可能存在风险的区域，900 mg/kg 以下采用分位数法进行分级，确保每一级别内所含的像元个数大致相等。

图12 重金属含量分布图Fig.12 Heavy metal content distribution map

研究区属于镍铜矿山区域，可以看出Ni 在整个研究区范围内的含量偏高，但几乎没有超过筛选值的情况，大部分区域的Ni 含量明显高于Pb 和Zn 的含量值，彭再华和蒋素芳（2018）研究也表明在哈密黄山南矿山中最主要的有价元素为铜和镍，其他金属元素含量较低。重金属Pb 和Zn 具有相似的分布趋势，李玲等（2020）在新疆矿区的研究也指出Pb、Zn 含量表现出显著相关性。3 种重金属分布的共同点为：在北部以及南部区域的重金属含量较高，在靠近居住区的东部区域含量较低。

图13展示了研究区的高程分布，可以看出研究区的地势呈现出东北高，西南低。经实地调研得知图12（a）中的A 区域用于堆放采矿渣，图12（a）中的B区域是采矿区，采矿被认为是重金属污染的最重要因素之一，重金属会伴随尾矿渣，废水等进入土壤（王海洋 等，2022），该采矿区于2007年开始勘探，已在此地进行了多年的采矿活动，由于研究区西南地势较低，重金属在重力作用下不断发生迁移，长时间作用下导致研究区南部、西南部区域受到污染，因此3 种重金属在南部的含量较高。

图13 研究区高程图Fig.13 Elevation map of the study area

研究区北部含量较高的区域对应图13 中高程出现显著变化的部分，图12（a）中可看出该区域图像与周围有明显差异，根据实地调研得知该区域为山体斜坡，坡上的土壤粒径极小呈细砂状，坡面平滑。在坡面上采集的4个土壤样本重金属含量均不低，其中3 个样本的Pb、Zn 含量高于平均值，2 个样本的Ni含量远高于平均值，且其中1 个达1220 mg/kg，表明该区域土壤的重金属含量偏高，反演结果与采样分析结果一致。

5 结论

目前针对高光谱图像的重金属反演研究较少，本文以3 种重金属Pb、Zn、Ni 为例，开展面向高光谱图像的反演研究。由于土壤样本点与图像土壤像元数差距悬殊，通过提取采样点邻近像元的方法扩充样本数据，同时增加样本的光谱多样性；采用FOD 增强对重金属反演有益的光谱信息，提高反演精度；由于高光谱数据的波段众多，使用CARS 算法选出特征波段以建立反演模型，从降维角度减弱小样本的影响。以新疆哈密黄山南铜镍矿区为研究区，分析了该反演方法的有效性。

研究结果表明，样本扩充有效地提升了3种重金属的测试集反演精度，缓解了模型的过拟合问题，扩充前3种重金属的训练集精度R2均在0.98以上，远高于测试集精度，模型明显过拟合，而样本扩充后，3 种重金属的测试集反演精度都得到了提升，同时，训练集R2也在0.8 以上，与测试集精度相近，过拟合问题得到显著改善。

3 种重金属在最佳阶数分数阶微分下的模型反演精度均高于使用整数阶微分的反演精度，说明分数阶微分能够更加有效地突出对土壤重金属反演有益的光谱信息。在最佳阶数微分的基础上，将CARS与常用的CC和GA两种波段选择算法进行了对比，CARS+PLSR 模型的精度最高，并且建模耗时短，本研究认为CARS 是3 种算法中最优的波段选择策略。

对CARS的波段选择结果进行分析，反演Pb时，选出的波段大多位于500 nm 范围内，该区间是对Pb 敏感有效的铁氧化物的特征波段，反演Zn、Ni时，选出的波段大多位于600—800 nm 和1600 nm附近，该波长范围属于对Zn、Ni 敏感的有机质特征谱段，与现有反演机理研究相一致，证明CARS算法能够有效地筛选出光谱中对反演重金属有益的特征波段。

本研究方法应用于航空高光谱图像上反演土壤重金属Pb、Zn、Ni 的含量，反演精度较高，并且反演出的3种重金属分布图与实际相符，表明该方法有很好的鲁棒性，具有反演多种土壤重金属的能力。但该方法尚需要在其它矿区、农作区验证其适用性，且土壤类型的影响也需要进一步探讨。