基于原位光谱的滨海滩涂土壤含水量预测模型

宋敬茹，满卫东，2，3，4*，高均海，张永彬，刘明月，2，3，4，郝玉峰，郑 浩，杨晓芜

（1.华北理工大学 矿业工程学院，河北 唐山 063210；2.唐山市资源与环境遥感重点实验室，河北 唐山 063210；3.河北省矿区生态修复产业技术研究院，河北 唐山 063210；4.河北省矿业开发与安全技术重点实验室，河北 唐山 063210；5.中煤科工生态环境科技有限公司唐山分公司，河北 唐山 063012）

【研究意义】土壤是生态系统的重要组成部分，是植物生长的必需因素。土壤含水量是土壤中必不可少的成分之一，是水文、林业和农业领域关注的热点[1]。植物依托土壤水分获取土壤盐分，土壤含水量直接影响土壤-植被生态系统中的水分和养分循环[2]。土壤含水量是评价土壤水分的重要指标，可以调节土壤所受的侵袭程度影响土壤持水能力。常规土壤含水量获取是一项需消耗大量的人力、物力和资金的工程，且大多是地面测定实验，难以实现大尺度的快速监测[3]。滩涂位于海洋与陆地交会处，因其独特的地理位置，使土壤样品采样难度增加，对滩涂土壤含水量进行区域尺度监测耗费更大。如何实现滨海滩涂土壤含水量预测，是目前海岸带生态环境管理亟待解决的重要问题，构建滨海滩涂土壤含水量预测模型对快速、高效地监测其土壤含水量具有重要意义。【前人研究进展】高光谱技术能够直接获取连续的土壤反射率信息，其对土壤含水量的细微变化较为敏感，可为滩涂土壤含水量定量反演提供技术支持[4]。众多学者通过寻求土壤含水量高光谱最佳波段，进行基于各种反射率数据的土壤含水量定量预测研究。Li 等[5]开发了一种利用手持式光谱辐射计估算土壤含水量水平的方法，发现当土壤含水量低于19% db时，土壤水分含量与对应高光谱反射率呈线性相关，得到了利用地面遥感获得的高光谱反射率来估算土壤含水量是可行的结论。Wu 等[6]采用流行的线性回归算法，选取695，711，736，747，767，778，796 nm 波段构建灌溉农田土壤含水量高精度反演模型。王静等[7]采集175个土样及相应含水量，并结合对应土样实验室光谱数据，建立了在1 423，1 524，1 746 nm 光谱波段下土壤含水量的预测回归模型。刘伟东等[8]利用相对反射率法、一阶微分法、反射率对数一阶微分、差分法等方法分别对10 种土壤含水量进行预测，结果表明：反射率对数一阶微分、差分法对土壤含水量具有很好的预测效果。以上研究多基于室内的高光谱数据，通过算法或光谱反射率变换研究发现反映土壤含水量信息的主要波段集中在1 400，1 900，2 200 nm波段周围，而关于野外原位土壤光谱对土壤含水量的响应以及二者之间的预测模型研究较少。自20世纪90年代起，土壤属性预测模型经历了从线性回归方法到机器学习方法的演变过程[9]，提高模型预测精度一直是相关领域关注的焦点。Moore 等[10]采用线性回归与判别分析结合的方法对土壤属性进行预测。Hong 等[11]利用支持向量机（SVM）估算土壤有机质（SOM）含量，比较了分数阶导数（FOD）与光谱变量选择技术（CARS）、弹性网（ENET）和遗传算法（GA）不同组合优化土壤近红外光谱模型，并取得了较好的效果。Tomislav 等[12]研究发现RF 模型相较于线性回归模型对土壤属性的预测误差明显下降，说明RF模型对于非洲土壤属性制图方面具有良好的适用性。Curtis等[13]使用统计学习方法与机器学习算法分别对玉米种植区土壤氮含量进行预测，对比发现基于机器学习算法构建的模型预测精度更高。Zhang等[14]使用多元线性回归模型和随机森林模型预测复垦土壤的有机碳含量，发现随机森林模型可以获得更高的R2和更低的误差指数。Forkuor等[15]利用遥感变量采用机器学习和多元线性回归对土壤性质进行了对比研究。多数研究表明机器学习模型较传统线性模型能更好地预测土壤的属性，但仍需基于不同土壤类型和区域特性选择适用的建模方法。

【本研究切入点】沧州市滨海区域是典型的淤泥质海滩，该区域内滩涂广泛分布，植物种类丰富，构建滨海滩涂土壤含水量预测模型对土壤-植被生态系统水分循环研究具有重要意义。【拟解决的关键问题】基于实测滨海滩涂土壤原位高光谱和对应土壤含水量数据，采用逐步线性回归（MSR）和支持向量机回归（SVR）方法构建滨海滩涂土壤含水量预测模型，比较不同变换形式的土壤反射率下表层滩涂土壤含水量预测模型精度，分析出最优预测模型，以期更加快速精准地获取滨海滩涂土壤含水量，解决土壤含水量获取复杂的问题，推动沧州沿海区域滩涂生态系统科学管理。

1 材料与方法

1.1 研究区概况

沧州（37°29′～38°57′N，115°42′～117°50′E）地处河北省东南部，北依京津，南接山东，总面积约为14 304.26 km2（图1），属于暖温带大陆性季风气候，年均降水量达600 mm 左右。四季分明，河流众多，素有“九河下梢”之称。沧州海岸带的近岸滩涂资源丰富，滨海滩涂因渤海潮汐堆积形成，以沙土、黏土和淤泥为主要土壤类型。滨海滩涂表层土壤（0～20 cm）粒径变化约为2.15～55.73 µm，属粉砂质土壤。

图1 研究区概况与采样点分布Fig.1 Survey of the study area and distribution of sampling sites

1.2 样品采集与处理

1.2.1 滨海滩涂土壤样品采集与土壤含水量测定 2019 年8 月，沿沧州海岸线从北到南依次确定14 个滨海滩涂土壤采样点（图1），每点采集表层（0～20 cm）土壤样品2 个并密封保存，用于测定土壤含水量。采用烘干法测定土壤含水量。烘干前将采集的滨海滩涂土壤样品称重记为m1，在105 ℃的烘箱中烘干土壤样品，直至前后两次称重恒定不变，称重记为m2，土壤水分计算公式如下：

其中，ω为土壤含水量，m1为烘干前土样的质量（g），m2为烘干后土样的质量（g）。

1.2.2 滨海滩涂土壤原位光谱采集与处理 利用FieldSpec 4 便携式光谱仪（350～2 500 nm）实测滩涂土壤原位光谱。高光谱数据采集时，天气晴朗无风，在每天10:00—14:00 获取滨海滩涂土壤高光谱数据。采样前先去除暗电流影响，并利用白板进行定标。FieldSpec 4便携式光谱仪传感器探头垂直于地面，并与地面相距1 m，每个采样点采集5条土壤反射率光谱曲线。

用ViewSpecPro 软件对每个采样点的5条实测光谱数据取平均值作为原始光谱反射率[16]。将取均值后的土壤光谱数据进行Savitzky-Golay 卷积平滑处理（窗口大小：50，多项式阶数：2），以去除噪声影响[17]。在对原始光谱反射率进行土壤含水量特征分析的基础上，采用倒数（1/R）、反射率倒数的对数（log（1/R））、一阶微分（R′）以及去包络线（CR）等变换方法[18]，获取土壤含水量预测模型构建的变量。

1.3 研究方法

1.3.1 相关性分析 利用Pearson 相关系数分析原位光谱和4 种光谱变换形式下的光谱数据与滨海滩涂土壤含水量的相关性[62]。根据土壤含水量与光谱数据的相关程度，快速准确地选取了土壤对应的特征光谱波长。这为进一步构建滨海滩涂含水量预测模型提供了便利。相关性分析具有很强的目的性，精度高，可靠性好。相关系数高于0.7 表示关系非常紧密；0.4～0.7 表明关系密切；从0.2～0.4 表明关系正常。

式中：xi表示光谱反射率值，表示光谱反射率值的平均测量值，yi表示土壤含水量的测量值表示土壤含水量的平均测量值，n表示样品的数量，i=1，2，3…，n，r表示相关程度。

1.3.2 多元逐步回归模型（MSR）基于多个变换形式下光谱反射率，采用多元逐步回归方法分别建立滩涂表层土壤含水量的高光谱预测模型，多元逐步回归模型的结构形式为：

式中：x1、x2、x3…xn为光谱反射率值，为表层滩涂土壤含水量预测值，α1、α2、α3…αk为自变量回归系数，b0为常数项。

1.3.3 支持向量回归模型（SVR）SVR 基于支持向量机（SVM）来拟合曲线，其是由SVM 扩展而来的[19]。SVR 与SVM 类似，SVR 使用条带拟合数据，其优点是可以高精度近似复杂的非线性连续函数[20]。优化SVR模型的目的是通过调整超参数来实现的，以找到超平面满足所有数据的最小距离的位置。

SVR 模型有两个非常重要的超参数cost 与gamma，cost 是惩罚系数，即对误差的宽容度，cost 值越大则越不能容忍出现误差，容易过拟合；cost值越小，容易欠拟合[21]。RBF 函数作为核函数（kernel），其自带超参数为gamma，其隐含地决定了数据映射到新的特征空间后的分布，gamma 值越大，支持向量个数越少；反之gamma 值越小，支持向量个数越多[22]。而支持向量的个数影响训练与预测的速度，故调整参数cost和gamma可以达到优化SVR的作用。

1.3.4 模型精度评定 滩涂土壤预测模型建立时，采用留一交叉验证法（leave-one-out cross-validation，简称LOO-CV）对模型进行验证[23]。该方法是一种K 折叠交叉验证法，使K 等于数据集中的数据量。每次只使用一个测试集，其余的都作为训练集[24]。这种方法得到的结果最接近于训练整个测试集的预期值，适合于基于本研究的小样本数据。

模型的检验精度采用预测值与实测值的调整型决定系数（Adjusted-R2）、均方根误差（RMSE）和相对分析误差（RPD）来评价[25]。Adjusted-R2越大、RMSE 越小，表明模型估算精度越高。RPD 值可以用来解释模型的预测能力，当RPD＜1.4 时，表明模型无法进行准确预测；当1.4≤RPD＜2.0 时，表明模型的预测能力一般；RPD≥2.0时，表明模型具有较好的预测能力[26]，公式如下：

式中，yi为土壤含水量实测值，̂为土壤含水量预测值，为土壤含水量实测值的平均值，n为样本数量[27]，用i=1，2，3，…，n表示，SD为实测值标准偏差。

2 结果与分析

2.1 滨海滩涂土壤原位光谱特征分析

滨海滩涂土壤含水量会对土壤光谱反射率造成较强的非线性干扰，比较不同滨海滩涂土壤含水量反射光谱，发现随着滨海滩涂土壤含水量增加，在原位光谱区间可见近红外（VIS-NIR）（350～2 500 nm）的光谱反射率有明显下降趋势（图2），与其他研究者得到规律一致[28]。不同滨海滩涂土壤样品的光谱走势基本相同，波峰与波谷所在波段一致，均在500～1 500 nm和2 000～2 200 nm光谱区间滨海滩涂土壤光谱反射率随波长增加而增加，1 600～1 800 nm 光谱区间滨海滩涂土壤光谱反射率波动较小，2 200～2 500 nm光谱区间滨海滩涂土壤光谱反射率随波长增加反射率降低。野外实测原位光谱受大气水汽的影响，在1 900 nm光谱处存在较为明显的反射峰。1 400～1 600 nm和2 000～2 200 nm光谱区间处分别存在较为明显的土壤含水量吸收谷。

图2 滨海滩涂土壤光谱曲线Fig.2 Spectral curve of tidal flat soils

2.2 滨海滩涂土壤含水量与光谱数据相关性分析

滨海滩涂土壤含水量与原位光谱相关性通过Pearson 相关系数进行分析，快速、准确筛选出用于构建滨海滩涂土壤含水量预测模型的光谱区间。原始光谱反射率（R）、倒数（1/R）、反射率倒数的对数（log（1/R））、一阶微分（R′）以及去包络线（CR）与滨海滩涂土壤含水量的相关程度在0.4以上的光谱波段如图3 所示。研究发现原位光谱R 和R′、1/R、log（1/R）、CR 均在1 000～1 300 nm、1 400～1 600 nm 和1 900～2 400 nm 维持与滨海滩涂土壤含水量相关性紧密，而光谱R′与滨海滩涂土壤含水量相关性高低分布较为明显。结合光谱特征区间选择用于建模的波段为1 400～1 600 nm 和1 900～2 400 nm。

图3 不同滨海滩涂土壤光谱反射率变换形式对应的滨海滩涂土壤含水量相关系数及波段Fig.3 Correlation coefficient and band of soil water content in tidal flats corresponding to different tidal flats soil spectral reflectance transformation form

2.3 滨海滩涂土壤含水量原位光谱预测模型对比分析

不同滨海滩涂土壤含水量预测模型精度显示：基于MSR构建的滨海滩涂土壤含水量预测模型中，精度最高的是基于R构建的滨海滩涂土壤含水量预测模型，精度最低的为基于CR 构建的滨海滩涂土壤含水量预测模型；利用MSR 方法构建的滨海滩涂土壤含水量预测模型由于Adjusted-R2值均小于0.4，RPD均低于1.4且RMSE较高，导致模型均不具有良好的预测能力。

表1 不同滨海滩涂土壤含水量预测模型的验证与对比Tab.1 Verification and comparison of soil water content prediction models in different tidal flats

与MSR 对比，同种光谱变换形式下，采用SVR 方法构建的滨海滩涂土壤含水量预测模型结果更优，Adjusted-R2和RPD 值均有提高，且RMSE 值也显著降低。其中，在R′变换形式下，采用SVR 方法构建的滨海滩涂土壤含水量预测模型的Adjusted-R2最高为0.81，RMSE 较小为2.56，且满足RPD≥2.0，模型具有很好的预测能力。比较基于以上10 个模型的滨海滩涂土壤含水量预测值与实测值（图4），分析滨海滩涂土壤含水量预测值与真实值的95%置信带和预测带，发现在R′变换形式下采用MSR 和SVR 方法得到结果中95%置信带和预测带会更加接近真实值和预测值的1∶1线，并结合对模型的精度评价，说明在R′变换形式下采用SVR方法构建的滨海滩涂土壤含水量预测模型为最优模型。

图4 五种光谱反射率处理模式的回归模型验证Fig.4 Regression model validation diagram of five spectral reflectance processing modes

3 讨论

全波段建模既会引起信息冗余，又容易产生严重的多重共线性[29]，对滨海滩涂土壤原位光谱区间进行选取是克服上述问题实现滨海滩涂土壤含水量预测模型构建的重要方法。通过分析原位光谱及其变换后的光谱数据与滨海滩涂土壤含水量相关性，筛选出用于构建滨海滩涂土壤含水量预测模型的光谱（1 400～1 600 nm和1 900～2 400 nm），此研究结果与Han等[30]研究的土壤含水率引起的光谱吸收谷（1 450 nm 和1 940 nm）相一致。产生的微小范围差异可能与滨海滩涂土壤类型、环境因素、植被生长等研究区现状有关[31-33]。

同种原位光谱变换形式下，采用SVR 方法构建的滨海滩涂土壤含水量预测模型精度高于MSR 方法。当土壤含水量与对应的VIS-NIR 光谱之间存在线性关系时，MSR 方法具有结构简单、计算方便、便于理解等优势[34]。当滨海滩涂土壤含水量与原位光谱较难呈现理想化的线性关系，利用SVR 方法处理其非线性关系，构建的滨海滩涂土壤含水量预测结果具有明显的精度优势。相应结论在Forkuor 等[15]的研究中得到了证实：比较多元线性回归（MLR）、随机森林回归（RFR）、支持向量机（SVM）、随机梯度提升（SGB）4 种方法，发现RFR 在大多数情况下提供了最高的准确性，而MLR 无法处理因变量和自变量之间的非线性关系。滨海滩涂土壤含水量与光谱数据可能呈现的线性关系与非线性关系，未来可采用决策树、随机森林和神经网络等多种机器学习方法建模，增加模型的可靠性，以探求更适用于滨海滩涂土壤的最优建模方法。

4 结论

利用FieldSpec 4便携式光谱仪高效获取原位滨海滩涂土壤VIS-NIR光谱数据（350～2 500 nm），对光谱反射率进行变换处理（R、1/R、log（1/R）、R′、CR），采用MSR 和SVR 两种建模方法基于土壤含水量高相关性光谱波段探究滨海滩涂土壤含水量预测模型。主要研究结论如下：

1）滨海滩涂土壤原位光谱反射率随土壤含水量增加而下降，在1 400～1 600 nm、1 900～2 400 nm 光谱存在较为明显的反射峰和吸收谷，光谱反射率经过变换处理后与滨海滩涂土壤含水量密切相关。

2）对比两种建模方法发现无论在何种变换形势下，机器算法SVR比普通线性回归方法MSR得到预测模型更加准确，Adjusted-R2更高，RMSE更小。其中在R′变换形式下利用SVR方法构建的滨海滩涂含水量模型建模效果最好，即为最优预测模型，其能够为在短时间内获取滨海滩涂土壤含水量提供数据支持。

致谢：唐山市科技计划重点研发项目（19150231E）、唐山市科技研发平台培养计划（2020TS003b）、中煤科工生态环境科技有限公司生产力转化基金项目（0206KGST005）同时对研究给予了资助，张阔、李京对土壤样品采集给予了帮助，谨致谢意！