基于不同建模方法的土壤pH值无人机遥感影像反演

王 羿， 丁鸣鸣， 何 菁， 邵光成

(1.河海大学农业科学与工程学院，江苏南京 210098； 2.南京市水务局，江苏南京 210036；3.南京市水利建筑工程检测中心有限公司，江苏南京 210036)

土壤酸碱性是判断作物生长和预测作物产量的重要因素，土壤酸性或碱性物质的输入会导致土壤理化性质的改变，从而影响土壤肥力，对作物生长发育产生抑制作用。同时土壤pH值与土壤养分的转化息息相关，土壤过酸或过碱都会导致其物理化学性质变差，造成土壤和植物的抗逆性减弱，抵御旱涝等自然灾害能力下降，生产能力降低。因此，及时、有效、准确、大面积的土壤pH值估算对于精准农业的实施有着很重要的意义。传统的土壤酸碱度一般采用化学分析法测定，这种方法虽然可以准确反映农田中土壤酸碱度的状况，但其时效性差，由于样本量限制，很难大尺度评估土壤酸碱度问题。搭载多光谱传感器的无人机不仅能够实现大规模、大数据采集，且其低成本，高时效、非接触的优势为获取土壤理化指标指导农业生产提出有效解决对策。在搭建好光谱分析模型后，无人机只需采集光谱，即可对当地土壤pH值进行定性定量分析预测。目前，通过光谱反演土壤理化指标已成为一种可靠并有效的方法，并取得了一定的进展，但不同预测算法反演模型精度和适用条件也各不相同。国内外学者已从反演模型的选择、光谱数据预处理构建光谱指数等多个角度对土壤理化指标开展大量研究，以期能够提高模型的预测精度和适用性。基于多光谱反演的方法主要有统计分析和机器学习2种类型，统计分析是利用土壤理化指标与光谱数据的相关关系建立回归模型，主要方法有多元线性回归(MLR)、偏最小二乘法回归(PLSR)等，机器学习是通过训练样本建立光谱数据与土壤指标关系模型，包括BP神经网络(BPNN)、支持向量机(SVR)、随机森林(RF)等算法。陈思明等采用多元逐步线性回归(MLSR)、BP神经网络(BPNN)和支持向量机(SVM)3种方法对湿地土壤有机质含量进行多光谱反演，通过引入最佳指数(OIF指数)法提取敏感波段，最终采用支持向量机(SVM)建模方法将确定系数提高0.124，降低0.106，模型精度得到提高。王凯龙等比较主成分回归(PCA)、偏最小二乘回归(PLSR)以及BP神经网络模型(BPNN)对土壤pH值的反演精度，发现采用PCA与PLSR方法获得的主成分作为BP神经网络输入变量的建模方法能够明显提高模型的精度与预测能力。赵静等分别采用普通克里格法、多元线性回归克里格法、遥感反演方法以及地理加权回归克里格法预测宁夏回族自治区石嘴山市城市土壤有机碳的空间分布，将地理要素加入到预测算法中，为多光谱遥感反演的普适性提供依据。上述研究通过不同反演方法或改进算法均提高了不同指标的模型预测精度，但对于最适宜土壤pH值的反演模型研究尚少，且对于光谱变量之间共线性问题并未讨论。因此，以江苏省南京市江宁区淳化街道农田的60个样本作为研究对象，分别采用多元逐步线性回归(MLSR)、BP神经网络模型(BPNN)、LASSO回归对于土壤pH值进行多光谱反演，对3种模型精度及适用条件进行评价，以期寻找出最适宜土壤pH值的反演模型，为快速准确土壤pH值指标反演提供理论依据。

1 材料与方法

1.1 研究区概况

研究区位于江苏省南京市江宁区淳化街道前桥头村农田(中心坐标为31°54′N、119°4′E)，研究区为亚热带季风气候，全年平均气温15 ℃，降水量1 069 mm，无霜期233 d。土壤质地为黏棕土，土壤pH值总体呈中性。研究区面积约为6 000 m，地形南高北低、西高东低。

1.2 土壤样本的采集和处理

样本采集于2021年4月29日，经前期采用无人机调查后，选取1块未种植作物的田块作为研究区，无人机采集时间为11：00—14：00，并在研究区开展同步样本采集，随机采集60份土壤样本，样本分布均匀覆盖整个研究区域，采集过程中全程同步采用GARMIN品牌的eTrex 221x的GPS记录样本位置信息，为提高定位精度，同一点位置信息记录 3～5次，取平均值作为最终结果，每个点位采集土壤表层(0～30 cm)土壤样本约50 g。收集的土壤样本采用密封袋封装，带回实验室后剔除杂物和石块，将样本自然风干，过1 mm筛，取10 g土壤样本按照水土比为2.5 ∶1配制土壤浸提液，振荡 30 min 后静置180 min，采用电极法测定土壤pH值，将电极插入试样悬浊液中，待数值稳定后读取记录土壤溶液pH值，样本pH值分布基本情况见表1。

表1 土壤pH值含量统计特征

1.3 多光谱数据获取

飞行前规划航线飞行任务，拍摄图像为正射影像(相机垂直于地面)，试验当天天气晴朗无云。本试验采用大疆精灵4多光谱无人机，为提高定位精度，采用网络实时差分定位(RTK)辅助定位，可将定位精度提升到水平方向±0.1 m、垂直方向±0.1 m(表2)。

表2 飞行平台及传感器主要性能参数

1.4 数据预处理

多光谱影像数据导出后采用Agisoft PhotoScan软件对6个波段图像进行拼接校正，地理参考坐标系统采用WGS84，得到研究区拼接图像后，由于各波段之间位置信息存在误差，需要进行地理配准消除各波段地理误差，将波段图像地理对齐。采用ENVI5.3对各波段影像进行地理配准后进行图层叠加，对叠加好的影像进行辐射定标，计算出环境入射光照度，后采用FLASSH大气校正模型对各波段光谱进行大气校正，消除大气和光照对土壤反射率的影响，将各波传感器记录的DN值转换为反射率信息。根据土壤取样点的GPS位置信息构建感兴趣区，以感兴趣区范围内地物的光谱平均反射率作为该取样点土壤表层光谱反射率，得到60个土壤样本对应的5个波段光谱信息。采用SPSS 26、Matlab 17a进行建模分析，Origin 2018进行图形绘制，ENVI 5.3进行反演制图。

2 模型的构建与评价

2.1 LASSO回归算法

套索(LASSO)算法由Robert Tibshirani于1996年首次提出，是一种从最小二乘法原理出发的收缩估计方法。它通过构造1个惩罚函数得到1个较精炼的模型，使得它压缩一些回归系数，即通过限制回归系数绝对值之和小于某个常数的约束条件，使残差和平方最小化，从而提高模型精度，因此保留了子集收缩的优点，是一种处理具有复共线性数据的有偏估计。使用惩罚值越大，进一步估计会使缩小值越趋近于0。通过调整惩罚值，可将系数矩阵整体变小，起到特征选择的作用。LASSO估计为

式中：表示惩罚值；表示响应变量(因变量)；表示数据样本自变量；为回归系数。

2.2 BP神经网络模型(BPNN)

BP神经网络模型通过模拟人脑学习的过程，能够学习和存储大量的输入-输出模式映射关系，训练过程是由信号的正向传播与误差的反向传播共同组成的，最终使输出模型的误差平方和最小，建立土壤pH值与多光谱反射率动态模型。

2.3 多元逐步线性回归(MLSR)

多元逐步线性回归是反演模型中经常使用的一种预测模型的方法，是将多个波段光谱反射率作为自变量，多波段进行最优组合共同预测或估计土壤属性的线性建模方法。即建立土壤pH值与多个波段光谱反射率之间的回归模型，多元线性回归为

=+++…+-1-1+。

式中：表示土壤pH值；表示对应波段光谱反射率，=1，2，3，…，-1，本研究-1=5；表示回归系数；表示随机误差项。

本研究分别将土壤pH值作为因变量，将5个多光谱波段反射率作为自变量进行多元逐步线性回归分析，即逐步回归事先设定光谱波段反射率进入、剔除的标准，按照自变量对因变量作用的显著程度，由大到小逐个引入模型，对于因变量不显著的变量也可能不引入回归模型，每向模型引入1个新的变量，同时考察原来在模型中的自变量是否还存在统计学意义，直至回归模型中只存在对因变量影响显著的自变量，回归结束，最后所得模型即为最优回归模型。

2.4 模型验证与评价

模型检验通过确定系数、均方根误差、值、方差膨胀因子(VIF)进行模型的评价，并讨论各模型的适用性，最终选取最适宜土壤pH值的反演建模方法。其中，反映因变量的全部变异能够通过回归关系被自变量解释的比例，越接近于1，表明回归关系越稳定，模型精度越高。为回归系统的拟合标准差，是衡量预测值与实测值参数，E越小，表明实测值与真实值误差越小，预测效果越好。检验为联合假设检验，值越大，表明模型越显著，拟合程度越好，当置信区间为95%，对应值<0.05时，可以认为组间存在显著差异，假设成立。方差膨胀因子VIF是衡量共线性问题的指标，范立新等认为，≥5或≥10时，可认为自变量间存在严重的共线性，模型自变量间存在重叠部分，VIF过大，会使模型对于普遍性解释不足，实际应用存在问题。

3 结果与分析

3.1 多元逐步线性回归(MLSR)分析

以多光谱反射率作为自变量，土壤pH值为因变量，采用多元逐步线性回归方法进行建模。由表3可知，基于全波段模型的明显高于单波段或组合敏感波段模型的，回归模型的与多光谱反射率输入量成正比，表明土壤pH值与5个波段均存在相关关系，值均小于0.05，表明自变量参数模型独立性较强。由表4可知，各波段VIF均大于10，表明各波段间均存在共线性问题，且0.65、0.73 μm 波段VIF远远大于其他波段，存在严重共线性问题。

表3 MLSR模型中各参数值

表4 MLSR共线性分析结果

多元逐步回归模型为

=-19807+28451+85256-187374+95802+4202。

式中：表示土壤pH值；分别表示0.45、0.56、0.65、0.73、0.84 μm波段光谱反射率。

3.2 BP神经网络模型分析

以全波段反射率为自变量，土壤pH值为因变量，45组数据作为样本，其中80%数据作为训练集，10%作为验证集，验证集目的是衡量网络泛化，决定网络是否终止训练，10%作为测试集，选取2层BP神经网络，设置隐含层节点函数为sigmoid，输出层节点函数为线性，经过反复测试，最终设定5个隐层神经元，选用Levenberg-Marquardt训练算法。训练模型见图1。

经过反复训练得到BPNN模型，输出土壤pH值预测值与实测值， 由图2可知， 经过15次迭代后模型收敛，其中，第9次训练结果精度最高，=0.896，=0.314，训练模型内部达到最优。将训练结果应用到外部验证集时=0.704，=0.351，表明模型应用精度不够理想。

3.3 LASSO回归模型分析

LASSO回归模型通过设定惩罚值，降低自变量维度来解决模型共线性问题，惩罚值越接近0，模型结果接近多元线性回归，采用ERIC准则进行验证，LASSO回归结果表明，惩罚值设定越小，回归系数越大，模型越大，越小，模型精度越高，但惩罚值设定太小不能解决自变量共线性问题，需要通过调整惩罚值大小来平衡模型精度与共线性问题，压缩模型回归系数，使参数更符合实际情况。

本研究分别设置惩罚值为0.01、0.015、0.019进行建模，建模结果见图3、表5。使模型既满足精度又解决共线性问题，增大惩罚值将回归系数不断压缩，尽可能减少自变量共线问题，并保证、参数在精度范围。经过测试，最终选择设定惩罚值为0.019时的结果作为LASSO回归最终模型，该模型将0.65 μm红光波段回归系数压缩到0，其他回归系数进一步压缩，可以很好地解决自变量间的共线性问题。

表5 LASSO模型汇总

3.4 结果与分析

综上述分析可知，单从和等2个角度对模型精度进行评价不能满足模型实际应用的要求，实际应用存在问题。采用BPNN方法构建的模型在验证时，由于模型解释变异性不足，导致模型应用中变化过大(图4-d)。MLSR模型与均表现较平稳，但由于0.65、0.73 μm波段的回归系数远大于其他变量回归系数，会导致模型稳定性不足，易受到2个波段影响，且在实际应用中出现不符合实际情况的表现，因此不选用MLSR与BPNN方法构建的模型。LASSO回归方法构建模型可以消除MLSR模型回归系数相差较大以及BPNN模型训练集不足的问题，最终选用LASSO回归模型进行反演及制图，采用决策树对土壤pH值分布进行分类(3种建模方法比较结果见表6及图4，经决策树分类统计后，分类统计结果见图5，土壤pH值反演及决策树分类图见图6)。夏子书等认为，土壤养分有效性的最高值大多出现在土壤pH值6.5～7.5之间，经过分类统计，研究区92.5%土壤的pH值分布在6.5～7.5，0.064%土壤的pH值分布在7.5～8.0之间，其余土壤pH值分布均小于0.05%(图5)，说明研究区土壤pH值较均匀，适宜作物生长。反演效果符合实际情况，未出现异常点，表明基于该方法进行土壤pH值的反演准确可行。

表6 建模方法比较

4 讨论与结论

研究结果表明，土壤pH值与0.45 μm蓝光波段、0.73 μm红边波段呈负相关，且与0.73 μm红边波段相关性最大，其余波段与土壤pH值为正相关关系。3种建模方法的均在0.6之上，表明土壤pH值与0.45、0.56、0.65、0.73、0.84 μm等5个波段存在相关关系，这与王凯龙等的研究结果相同。

采用MLSR、BPNN、LASSO回归3种建模方法对土壤pH值进行预测， 发现多元逐步线性回归模型=0.776，=0.432，从数学模型角度达到不错精度，但由于多光谱自变量之间存在严重的共线性问题，0.65、0.73 μm波段的VIF分别为822.821、1 213.727，远远大于其他波段值，因此存在严重共线性问题，波段间存在数据重叠，回归模型的估计值极不稳定，对自变量的变化非常敏感，这是由于研究区下垫面为农田，导致下垫面属性较单一，建模集差异性较小，过于均匀的土壤对作物生长有利，但对模型构建不利，这也是模型对土壤pH值差异较大的区域解释性不足等问题的原因，因此采用MLSR构建的模型无法应用到实际中。神经网络模型的=0.896，=0.314，训练模型的精度是所述模型中最好的，但将模型应用到训练集以外的15组数据中进行验证时，模型的=0.704，该值的降幅为0.192，大于其他2个模型，说明训练集数据量不足以解释该模型的变异，也表明神经网络对局部的拟合能达到最优，这与齐琳的研究结果一致，如需将模型应用到实际中，还需对8模型进行大量的数据训练，提高模型对变异的解释能力。LASSO回归算法通过调整惩罚值将共线性强自变量回归系数进行压缩或去除以解决模型的共线性问题，最终模型=0.618、=0.430，将0.65 μm波段回归系数压缩到0，0.73 μm波段回归系数从-187.374压缩到-48.058，并保证模型的预测精度。综上，对应惩罚值为0.019下的LASSO回归模型为该区域最适宜土壤pH值反演的模型。

土壤pH值与0.45、0.73 μm波段光谱反射率呈负相关，与0.56、0.65、0.84 μm波段光谱反射率呈正相关。5个波段MLSR建模方法=0.776，=0.432，但土壤光谱反射率0.45～0.84 μm波段间VIF均大于10，存在多元共线性问题，0.65、0.73 μm 波段VIF分别为822.821、1 213.727，与其他波段存在严重的共线性问题，该模型对土壤pH值解释性不足。BPNN模型建模集预测精度最高，(建模集)=0.896，(建模集)=0.314，但将模型应用到验证集数据时(验证集)=0.704，相较于其他2种建模方法，降幅偏大，BPNN建模方法能达到局部最优，但由于样本量限制，模型解释变异性不足，还需经过大量数据训练。LASSO回归方法=0.618，=0.430，回归系数将MLSR模型回归系数整体压缩，将0.65 μm波段回归系数压缩为0，0.73 μm波段回归系数从-187.374压缩到-48.058，压缩系数同时保证建模精度，LASSO回归模型为该区域土壤pH值与多光谱最适宜的建模方法。经过预测，研究区92.6%的土壤pH值分布在6.5～7.5之间，土壤pH值较均匀，适宜作物生长。