基于分数阶微分光谱指数的冬小麦根域土壤含水率估算模型

刘 浩，杨锡震，张 蓓，黄嘉亮，赵 笑，吴雨箫，向友珍，耿宏锁，陈皓锐，陈俊英※

（1.西北农林科技大学水利与建筑工程学院，杨凌 712100；2.西北农林科技大学旱区农业水土工程教育部重点实验室，杨凌 712100；3.西北农林科技大学综合素质教育学院，杨凌 712100；4.中国水利水电科学研究院流域水循环模拟与调控国家重点实验室，北京 100037）

0 引言

在农业生产中，土壤含水率（soil moisture content，SMC）是影响作物生长发育的重要条件[1]，也是指导灌溉的理论依据[2]，而监测作物根系所在位置的土壤含水率更有助于实现作物精准灌溉[3]。然而，土壤水分极易受环境的影响[4]，并且传统的土壤含水率监测的方法由于费时费力而难以满足精准农业管理对土壤含水率的监测需求[5]。随着高光谱遥感技术在农业中的推广，及时地无损监测土壤含水率逐渐成为可能。高光谱遥感技术以其可以提供丰富的光谱信息和大面积非接触监测等优势，在土壤含水率的监测研究中得到重视[6-7]。然而，从高维且明显存在光谱噪声的作物冠层高光谱数据中挖掘出敏感波段仍存在一定难度[8]。因此，尽可能少地损失光谱特征波段、较为彻底地去除光谱噪声是依据高光谱数据建立土壤含水率监测模型过程中的关键步骤。

目前，在高光谱数据的挖掘研究中，将光谱进行数学微分变换和筛选敏感光谱指数成为探寻敏感波段的常见方法。敏感光谱指数是通过光谱指数与对应的反演指标经二维相关分析筛选确定[9]，既考虑了波长间的关系，又考虑反演指标与波长之间的相关关系，可减少无关波长的影响[10]。大量国内外研究表明[11]，根据微分变换后的光谱信息构建的光谱指数可以有效克服人工观测主观性强、效率低的弊端。如GE 等[12]将光谱经一阶（first differential，FD）和二阶微分（second differential，SD）光谱变换后，筛选出敏感波段并构建了SMC 最优监测模型；刘爽等[13]对光谱指数进行FD 处理后建立敏感光谱指数，从而提升了模型定量估算叶片叶绿素的效果。上述研究均体现了整数阶微分联合光谱指数可以有效探寻对反演指标敏感的波段，但可能会造成部分光谱信息丢失，并引入大量噪声峰[14]。与整数阶导数相比，分数阶微分（fractional-order differentiation，FOD）能以更小的阶数间隔变化，细化光谱反射率信息，并有效避免引入噪声峰的问题，实现数据的深度挖掘[15]。如张俊华等[16]利用FOD 技术构建优化光谱指数，结果表明采用FOD技术可大幅提高土壤有机质含量反演模型的性能。这表明分数阶微分技术的优势在于挖掘了光谱中细微的信息，并且光谱保持了光谱原始的特征。此外，采用不同模型建立的土壤含水率反演精度存在差异，机器学习算法作为解决非线性问题的可靠方法，在反演土壤含水率方面有不错的效果，目前偏最小二乘回归（partial least squares regression，PLSR）模型是高光谱建模中最常用的数学模型，随机森林（random forest，RF）和BP 神经网络（back propagation neural network，BPNN）[17]模型在高光谱反演中也有不错的表现。目前，FOD 技术多应用在以土壤为观测目标反演土壤指标如土壤有机质及土壤含盐量等，而利用FOD 技术协同光谱指数反演作物覆盖条件下的土壤含水率的研究仍鲜有报道。

基于此，本研究以关中平原种植的冬小麦为研究对象，建立典型光谱指数的冬小麦SMC 反演模型，并将光谱反射率经FOD 处理后筛选敏感光谱指数构建出SMC最优变量子集，以此尽可能地挖掘光谱细节、较为彻底地去除噪声，同时对无信息变量进行有效去除，构建基于最优变量集的SMC 监测模型，以期为精准无损地监测冬小麦根域土壤含水率及指导当地精准农业灌溉和农业生产提供科学支撑和参考。

1 材料与方法

1.1 研究地概况

本试验于2020 年10 月至2022 年6 月在西北农林科技大学旱区农业水土工程教育部重点实验室节水灌溉试验站进行（108°24′E，34°20′N）。该地区是典型的半湿润易旱区，气候为暖温带季风半湿润气候。试验站海拔524.7 m，多年平均气温12.9 ℃，年平均降雨量561.1 mm，土壤质地为中壤土，田间持水率为26.0%（质量分数），凋萎系数为8.6%（质量分数），平均土壤干容积密度为1.44 g/cm3，0～30 cm 土层的土壤pH 值为8.14，有机质含量为12.02 g/kg，全氮含量为0.89 g/kg。

1.2 试验设计

试验设置了4 个灌溉水平，分别为充分灌溉（T1）、轻度水分胁迫（T2）、中度水分胁迫（T3）、严重水分胁迫（T4），土壤含水率分别控制在田间持水量的40%～50%、＞50%～65%、65%～80%以及＞80%～100%。试验地块的布局完全随机。每个灌溉水平有3 个重复，总共12 个地块，地块大小为4 m×4 m，行距为0.25 m，种植密度为100 万株/ha²。灌溉采用滴灌方式进行，滴灌器流速为1.2 L/h，间距为0.1 m。播种前施氮肥和磷肥（分别为240 kg/hm²）作为基肥。冬小麦品种为小偃22 号。

1.3 数据采集和处理

光谱采集仪器为美国ASD 公司生产的Field spec-3 高光谱仪。测定时间选在11：00—13：00，测定选择天气晴朗无风、视野良好的天气状况，以尽量减少太阳角度和树冠阴影的变化。光谱仪的采集范围为350～2 500 nm，光谱采样间隔为1 nm。测量时光谱仪探头距冬小麦冠层约15 cm，始终与地面保持90°，视场角为25°。光谱测定前，采用反射率为1 的漫反射参考板对仪器进行校正，工作过程中每5 min 优化1 次，且每隔5 min 采集1 次暗电流对仪器优化。每个小区以五点梅花状采样方式对土壤样点周边小麦进行测量，取其平均值作为监测点最终光谱值。

本研究分别于2021 年1 月9 日、2 月2 日、3 月13 日、3 月23 日、4 月13 日、4 月30 日和2022 年1 月12 日、2 月26 日、3 月28 日、4 月3 日、4 月11 日、5 月11 日测取土壤含水率数据，每次试验土壤含水率数据收集日期与高光谱反射率数据测量日期相同，每次试验12 个数据样本，试验阶段内共获得了144 个数据样本。

为保证建模样本和验证样本范围一致且分布均匀。样品集的划分采用了梯度划分的方式，根据冬小麦土壤样本中体积含水率的数值大小，将样本进行排序，根据排序的顺序，随机选用2/3 的样本用于建模，剩下1/3 的样本用于验证。在本研究中，建模集选用样本数为96 个，验证集选用样本数为48 个，各数据集之间的样本分布基本一致。

1.4 光谱数据的处理和转化

本试验中的144 个样品通过View Spec Pro 软件获得原始光谱，光谱波段为350～2 500 nm。由图1a 的原始光谱全波段图看出在1 350～1 450、1 800～2 050、2 300～2 500 nm 波段范围内出现水汽吸收峰，噪声明显较大，不适用于建模[18]。因此对光谱进行删选，结果如图1b 所示。对光谱数据进行一定的处理可以一定程度上消除测试环境和其他干扰因素的影响，提高数据的信噪比。本文对光谱数据进行SG 平滑预处理之后进行了分数阶微分（FOD）处理。

图1 原始及去除水汽带的光谱反射率曲线Fig.1 Spectral reflectance curves of original and removed water vapor belts

1.4.1 Savitzky-Golay（SG）平滑处理

采用SG 平滑处理可以在一定程度上提高降低噪音的干扰，提高光谱曲线的平滑性[19]。本研究中SG 平滑处理在The Unscrambler X 10.4 中实现。

1.4.2 分数阶微分处理

根据本研究的适应性选择了分数阶微分的Grünwald-Letnikov（G-L）定义。利用式（1）在MATLAB R2016b对光谱数据进行处理，在0～2 阶区间内进行了10 次分数阶变换，步长为0.2。

式中 λ为自变量，υ 为 阶数，Γ为Gamma 函数[20]，n为微分上下限差。当υ=0.0 时则表示原始光谱，υ=1.0 时表示一阶微分，υ=2.0 时表示二阶微分。分数阶微分处理在Matlab R2016a 中编写代码实现。

1.5 光谱指数

1.5.1 典型光谱指数

已有研究表明[21-26]，表1 中各指数与植被覆盖下的土壤含水率有着较强的相关性。通过反射率计算的典型光谱指数的具体计算式见表1。

表1 土壤含水率的代表性典型光谱指数Table 1 Representative typical spectral index of soil moisture (SMC)

1.5.2 分数阶微分光谱指数

基于经分数阶微分处理后的光谱数据，构建比值指数[27]（Ratio index,RI）、归一化指数[28]（Normalized difference index,NDI），其计算式分别为

式中Ri和Rj代表i和j波长位置的原始光谱反射率，i≠j。

1.6 机器学习建模

本研究使用PLSR、RF、BPNN 这3 种模型来建立高光谱遥感反演冬小麦根域土壤含水率模型[29]。通过The Unscrambler X 10.4 软件实现PLSR 算法，通过Matlab R2016a 软件来实现RF 和BPNN 模型算法。根据模型的最优反演结果，PLSR 模型采用的主因子数为均为6；RF 模型所取的决策树数目和最小叶子数均为100和30，迭代次数为2 000 次；BPNN 模型所采用的隐藏层神经元个数为5 个，迭代次数为8 000 次。

1.7 模型评价指标

模型对比分析利用以下指标：建模集、验证集决定系数（coefficient of determination，R2）、均方根误差（root mean square error，RMSE）、相对分析误差（relative percentage diff erence，RPD）。R2用来判定模型的稳定程度，越接近1 说明模型的稳定性越好；RMSE 用于表征模型的准确性，其值越小，表明模型的精度越高[30]；RPD 越大表明模型的预测能力越好，反演精度越高。其中RPD 一般被分为6 个等级当RPD≥2.5和2.0≤RPD＜2.5 时，分别表明模型有极好和很好的预测能力；当1.8≤RPD＜2.0时，表明模型具有定量预测能力；当1.4≤RPD＜1.8 时，表明模型具有一般的定量预测能力；1.0≤RPD＜1.4 时，表明模型只有区别高低值的能力；RPD＜1.0 时表明模型没有预测能力[31]。

2 结果与分析

2.1 典型光谱指数相关性分析

将原始反射率数据分别经SG 平滑处理后建立典型光谱指数，并分析典型光谱指数与土壤含水率的相关性，由分析结果可知，各典型光谱指数与土壤含水率之间的相关性差异较大。其中OSARI、NSRI、NDVI、WI 等光谱指数与土壤含水率显著相关（P＜0.001），相关系数分别为0.291、-0.603、-0.281、-0.463，而DWI、MSI、SRWI、NDII、NDred edge 等光谱指数与土壤含水率相关性较差，相关性系数分别为-0.255、0.192、-0.237、-0.248、-0.200。

2.2 典型光谱指数建模

以原始光谱反射率经SG 光滑变换后建立的9 个典型光谱指数为自变量，以土壤含水率作为因变量分别建立PLSR、RF、BPNN 模型。

由图2 可知，在典型光谱指数建模结果中，PLSR模型具有不错的反演精度（R²均为0.55），PLSR 模型建模集和验证集数据较为均匀地分布在1∶1 线附近。

图2 典型光谱指数反演模型反演SMC 的散点图Fig.2 Scatter diagram of SMC inversion by typical spectral index inversion model

2.3 分数阶微分图像特征

为了研究分数阶微分对光谱指数的影响，更好地探究分数阶微分处理后的光谱指数特征相应情况，在本研究中，将原始光谱反射率进行了以0.2 阶为步长的0～2.0 阶的分数阶微分处理。为了更清晰地观察到微分阶数对光谱曲线的特征影响，取各分数阶微分下的均值光谱绘制如图3。

图3 均值FOD 光谱曲线图Fig.3 Mean FOD (fractional-order differentiation) spectral graph

从图中可以看出，随着微分阶数的增加，光谱的微分反射率值逐渐向0 靠近，光谱曲线在0～1.0 阶时还能基本保持原有的曲线形状，在1.2～1.6 阶变化过程中曲线特征变弱，在＞1.6～2.0 阶时已经趋近于一条y=0 的直线。光谱曲线上的峰谷，在低阶时随着微分阶数的增加逐渐明显，且愈发尖锐，光谱的基线漂移基本已经被消除，表明分数阶微分可以显著突出光谱的曲线特征，帮助凸显对土壤含水率响应敏感的特征波长区间。

2.4 土壤含水率与不同阶数RI、NDI 相关性分析

光谱指数可以将反射率通过比值、差值、线性组合等多种组合来增强或揭示隐含的植物信息。由典型光谱指数与土壤含水率的相关性分析可知，与SMC 相关性较好的光谱指数多为光谱反射率经比值、归一化计算得到的指数。

因此为了精准地提取出土壤含水率与光谱指数相关性最高的波长组合，本研究选取了RI 和NDI 两种光谱指数，在删选后的波长范围内，分别将不同阶数的反射率的所有波长两两组合计算RI 和NDI 两种光谱指数，再利用相关矩阵法分别对不同阶数的RI、NDI 和土壤含水率间关系进行相关性分析，并绘制相关性矩阵图，如图4 所示。以最大相关系数所在i和j波长位置作为最优的波长组合。

图4 不同阶数光谱指数与土壤含水率的相关矩阵Fig.4 Correlation matrix of spectral indices of different orders and soil moisture content

由图4 提取20 个敏感光谱指数与土壤含水率的相关系数最大值均高于0.63，均高于典型光谱指数中相关性表现最好的NSRI，其中与SMC 相关系数最高的NDI 值为0.2 阶，相关系数为0.76，波长组合为816 和752 nm；与SMC 相关系数最高的是RI 值为0.6 阶，相关系数为0.77，波长组合为817 和737 nm。由图4 可知，敏感光谱指数与SMC的相关系数变化大致呈阶梯状变化，敏感RI 和NDI 与SMC 的相关性在各阶表现基本接近，在0.2～1.0 阶敏感光谱指数与SMC 的相关系数稳定于0.77附近，而在1.2 阶开始相关系数开始逐步降低并在1.6 阶稳定在0.65 附近。

2.5 改进光谱指数与波长组合的提取

图5a、5b 分别表示各阶RI 和NDI 与SMC 相关系数最大时的波长位置，由图可知，8 个指数（0.2、0.4、0.6、0.8 阶RI 和NDI）优选出的波长组合均在815 和750 nm 附近，7 个指数（1.4 阶RI、1.6、1.8、2.0 阶RI和NDI）优选出来的波长组合均在700 和600 nm 附近，而1.0、1.2 阶处的敏感光谱指数的波长组合较为混乱，这可能是由于光谱反射率曲线经1.0、1.2 阶分数阶微分处理后曲线特征开始变弱所导致。

图5 不同阶敏感光谱指数RI、NDI 与SMC 相关性最大时的波长位置Fig.5 Wavelength position of different order sensitive spectral indices RI,NDI and SMC with the greatest correlation

根据上述研究提取出的20 组最优波长组合，将其分别代入对应的光谱指数计算式中进行计算，得到对应改进光谱指数值，再建立基于改进光谱指数的冬小麦根域土壤含水率反演模型。

本研究为了对比基于不同改进光谱指数组合的建模效果，将20 个最优光谱指数划分为3 组作为模型输入变量：第一组变量为由比值指数（RI）计算的10 个改进光谱指数，分别为0.2～2.0 阶（步长0.2），称为组合1；第二组变量为由归一化指数（NDI）计算的10 个改进光谱指数，分别为0.2～2.0 阶，称为组合2；第三组变量为所有敏感光谱指数中与SMC 的相关系数由高到低排序的前10 位，分别为0.2～1.0 阶（步长0.2）的RI 和0.2～1.0 阶（步长0.2）的NDI，称为组合3。

2.6 基于FOD 的改进光谱指数的土壤含水率反演模型

不同输入变量和建模方法组合的冬小麦土壤含水率反演模型及验证集预测结果分别如图6 所示。由图6 可以看出以3 组改进光谱指数为自变量结合3 种机器学习模型反演土壤含水率的精度均高于典型光谱指数的反演模型，研究结果说明分数阶微分结合敏感光谱指数的建模决策改进了光谱指数的反演土壤含水率的精度。其中3 种自变量组合建立同一模型的精度相差不大，但3 个自变量组合的RF 模型精度都优于较PLSR 和BPNN，说明同一自变量组合所建不同模型精度存在较大差异。由图6 可以看出RF 模型所得的预测值较PLSR 和BPNN更加接近于实测值。3 个自变量组合建立的RF 模型的建模集和验证集的决定系数均高于0.71。其中精度最高的是组合2-RF，其建模集和验证集的决定系均达到0.75，均方根误差分别为0.021 和0.024 g/g，表明模型具有最好的定量反演能力。而组合1 和组合3 变量建模集和验证集的决定系数最高分别可以达到0.71（组合1-RF）、0.73（组合3-RF），也具有较高的反演精度。

图6 不同输入变量和建模方法组合的SMC 反演的散点图Fig.6 Scatter diagram of SMC inversion with different input variables and modeling methods

3 讨论

光谱技术在快速有效监测土壤属性方面具有较大优势，其中高光谱以光谱信息丰富、波段连续性强而常用于监测土壤含水率[32]，但同时存在数据冗余、共线性强等特点导致反演模型适用性较差。而光谱指数能够充分利用高光谱丰富的信息，减少其他土壤信息对土壤含水率的影响，去除冗余信息，降低模型的复杂性。如REN等[23]利用光谱指数和红边参数对冬小麦的SMC 进行了反演，结果表明基于优化的NSRI 与SMC 有较高的相关性，决定系数R²达到0.52（P＜0.01），其中NSRI 经对数和一阶导数处理后的模型表现最佳，其决定系数达到0.81～0.92。研究在分析典型光谱指数与SMC 的相关性时发现NSRI 表现最好，这也符合REN 等[23]研究的结果规律。赵英时[27]指出，在光谱指数采用的波段中，绿色植物的（叶绿素引起的）强吸收的可见光波段（0.6～0.7 nm）及对绿色植物（叶内组织引起的）高反射和高透射的近红外波段（＞0.7～1.1 μm）,这两个波段反射率随着叶冠结构、植物覆盖度而变化。而SMC 与上覆植被叶片的含水率之间具有一定的相关性，SMC 的高与低将不同程度对叶片结构造成影响[33]。本文将光谱数据经分数阶微分处理后筛选的敏感光谱指数的波长组合多在815 和700、700 和600 nm 附近，由此表明本文筛选的敏感光谱指数的波段组合可以用于冬小麦根域SMC 的反演。

光谱反射率在经微分变换后，能提供较原始光谱精度更高的分辨率和更清晰的光谱轮廓并增强反射率与土壤属性间的相关性[11]，王瑾杰等[34]将光谱进行FOD 处理后发现低阶FOD 处理后的光谱与SMC 的相关系数较高，在0.4 阶处理时相关系数达到0.768，低阶FOD 的相关系数整体上优于高阶FOD，本研究结果与上述研究结果近似相同。本研究发现，当分数阶在低阶增加时，光谱曲线的重叠峰和基线漂移现象逐渐消除，导致敏感光谱指数与SMC 的相关性有略微提升，但随着分数阶增加到高阶（1.2）附近时，由于光谱强度的幅度的下降，相关系数出现显著下滑。本研究以0.2 阶为间隔对冬小麦冠层反射率进行FOD 处理，筛选各阶敏感波长组合，并分组建模，模型的决定系数最高可以达到0.75，较典型土壤含水率光谱指数建立的反演模型精度提升较多（R²提升136%）。说明光谱数据经分数阶微分处理后可以改进光谱指数与土壤含水率的相关性，为获取高相关性敏感指数提供了更大可能，有助于提高冬小麦根域SMC 监测模型的精度。

将机器学习回归模型应用在监测土壤含水率上具有显著优势，然而不同机器学习模型会显著影响监测模型的反演精度[35]。在本研究中，以各组改进光谱指数为自变量建立的RF 反演模型具有最佳的准确性，研究分析这可能是因为RF 是一种具有集成思维的机器学习方法，具有很强的自学能力，对光谱反射率中的噪声和异常值有较强的耐受性，并且不容易出现过拟合等问题[36]。而PLSR 和BPNN 模型的监测精度较低，这可能是由于样本相对较少导致泛化能力较低所致[37]。综上，RF 可作为监测冬小麦SMC 的首选方法。在未来的农业生产中，可以考虑将高光谱遥感作为监测土壤含水率的主要手段之一，从而提升获取土壤水分信息的效率，以更省时省力的方法获取作物根系土壤含水率数据，从而指导农田制定精准灌溉策略。

目前，基于分数阶微分协同光谱指数建立土壤含水率的反演模型研究中仍有一些问题有待解决，如不同区域、相同区域不同作物，甚至相同区域相同作物但不同时期的冠层反射率进行分数阶转化时与SMC 的相关系数变化及模型反演精度可能有所不同。此外，今后的研究还将考虑土壤属性中盐分、pH 值、养分和有机质等因素对通过冠层光谱反演土壤含水率的影响。

4 结论

1）分析现有典型光谱指数与冬小麦根域SMC 的相关性发现，NSRI 与土壤含水率相关性较好，建立的PLSR 模型反演SMC 精度最高。

2）光谱曲线在0.2～1.0 阶能在基本保持曲线形状的基础上增加波峰波谷特征，而在1.6 阶之后接近y=0 的直线。基于分数阶微分改进后的RI 和NDI 与土壤含水率的相关性较典型光谱指数有明显提升，并在0.2～2.0阶上呈阶梯状变化，低阶分数阶微分（0.2～1.0 阶）相较于高阶分数阶微分（1.2～2.0 阶）相关性更好。敏感光谱指数采用的反射率波段主要集中在600～815 nm。

3）基于改进光谱指数的土壤含水率反演模型中表现较好是RF 模型，其中以分数阶微分联合归一化敏感光谱指数为自变量的RF 模型具有最佳反演精度，可为反演冬小麦根域土壤含水率提供依据。