引用本文
邓昀, 王君, 陈守学, 石媛媛. 基于PCA-DBO-SVR的林地土壤有机质高光谱反演模型[J]. 光谱学与光谱分析, 2025,45(2): 569-583.
DENG Yun, WANG Jun, CHEN Shou-xue, SHI Yuan-yuan. Hyperspectral Inversion Model of Forest Soil Organic Matter Based on PCA-DBO-SVR[J]. Spectroscopy and Spectral Analysis, 2025,45(2): 569-583.
Doi:10.3964/j.issn.1000-0593(2025)02-0569-15  
Permissions
Copyright©2025, 《光谱学与光谱分析》期刊社
《光谱学与光谱分析》期刊社 所有
基于PCA-DBO-SVR的林地土壤有机质高光谱反演模型
邓昀1,2, 王君1,2, 陈守学2,*, 石媛媛3
1. 广西嵌入式技术与智能系统重点实验室, 广西 桂林 541004
2. 桂林理工大学计算机科学与工程学院, 广西 桂林 541004
3. 广西壮族自治区林业科学研究院, 广西 南宁 530002
*通讯作者 e-mail: 7666817@qq.com

作者简介: 邓 昀, 1980年生,桂林理工大学计算机科学与工程学院教授 e-mail: 574359451@qq.com

收稿日期: 2023-12-27 接受日期: 2024-05-29
基金: 国家自然科学基金项目(32360374),广西自然科学基金项目(2018GXNSFAA281235)资助
摘要

森林土壤有机碳(SOC)是土壤中的有机物质(SOM)的碳部分, 它对维持森林生态系统的平衡和稳定非常重要。 传统实验通过化学方法分析土壤中有机物质的含量进而计算土壤中的有机碳, 此类化学方法费时费力且产生化学废水污染环境。 高光谱技术可以非接触、 高效率地检测出土壤的养分信息。 针对现有机器学习土壤有机质预测模型的精度和计算效率方面的不足, 以广西国有黄冕林场和国有雅长林场为土壤样品采集点, 基于全光谱数据利用主成分分析算法(PCA)筛选特征波段的最佳波长数量, 并利用比一阶微分处理数据更加精细且能平衡光谱噪声和光谱分辨率之间的关系的分数阶微分为预处理方法之一对光谱数据进行变换处理, 最后采用相对于传统的中心化算法拥有较高鲁棒性和容错能力的蜣螂算法(DBO)对支持向量回归机(SVR)的高斯核函数的参数组合进行优化。 研究结果表明, PCA-DBO-SVR模型可以有效提高土壤有机质预测的决定系数 R2并降低预测均方根误差(RMSE)。 PCA-DBO-SVR在对比预测模型中表现出最佳的泛化性能和准确度, 其验证集 R2为0.942, RMSE为2.989 g·kg-1, 展现了较好的准确性。

关键词: 近红外光谱; 分数阶微分; 蜣螂优化算法; 土壤养分预测; 支持向量回归机
中图分类号:TP399 文献标志码:A
Hyperspectral Inversion Model of Forest Soil Organic Matter Based on PCA-DBO-SVR
DENG Yun1,2, WANG Jun1,2, CHEN Shou-xue2,*, SHI Yuan-yuan3
1. Guangxi Key Laboratory of Embedded Technology and Intelligent System, Guilin University of Technology, Guilin 541004, China
2. School of Information Science and Engineering, Guilin University of Technology, Guilin 541004, China
3. Guangxi Zhuang Autonomous Region Forestry Research Institute, Nanning 530002, China
*Corresponding author
Abstract

Soil Organic Carbon (SOC) is the carbon component of Soil Organic Matter (SOM) and is crucial for maintaining the balance and stability of forest ecosystems. Traditional methods for analyzing the organic matter content in soil involve chemical analysis, which is time-consuming and labor-intensive, and generates chemical wastewater that pollutes the environment. Hyperspectral technology offers a non-contact, efficient means of detecting soil nutrient information. Addressing the limitations in the accuracy and computational efficiency of existing machine learning models for soil organic matter prediction, this study uses soil samples from Guangxi State-owned Huangmian Forest Farm and State-owned Yachang Forest Farm. Using full-spectrum data, Principal Component Analysis (PCA) was employed to select the optimal wavelength number for feature bands. Fractional-order differentiation, which processes data more precisely than first-order differentiation and balances spectral noise and resolution, was used as one of the preprocessing methods to transform the spectral data. Finally, the Dung Beetle Optimizer (DBO), known for its higher robustness and fault tolerance compared to traditional centralized algorithms, was used to optimize the parameter combination of the Gaussian kernel function in Support Vector Regression (SVR). The results indicated that the PCA-DBO-SVR model effectively improved the coefficient of determination (R2) for soil organic matter prediction and reduced the Root Mean Square Error (RMSE). The PCA-DBO-SVR model demonstrated the best generalization performance and accuracy among the compared prediction models, with a validation set R2 of 0.942 and an RMSE of 2.989 g·kg-1, showcasing excellent accuracy.

Keyword: Near-infrared spectroscopy; Fractional derivative; Dung beetle Optimization algorithm; Soil nutrient prediction; Support vector machine regression
引言

森林土壤有机碳(soil organic carbon, SOC)是土壤中的有机物质(soil organic matter, SOM)的碳部分, 它对维持森林生态系统的平衡和稳定非常重要。 SOC在改善土壤物理性质和土壤缓冲方面起着决定性作用[1]。 精准农业[2]和土地可持续区域开发需要快速监测土壤的SOM信息[3]。 传统的化学方法即在实验室内利用重铬酸钾氧化加热等方法测定土壤中的SOM, 这种测量方法的人力物力成本较高, 并且会产生废水。 物质的近红外光谱是其分子结构的反映, 是监测土壤特性的传统方法的快速且无损的补充[4]。 SOM的独特光谱特征源于光谱区域中C— H, O— H和N— H化学键的强吸收模式[5]。 故采集光谱数据进而利用算法精准且快速地测定土壤的SOM对林地土壤的肥力调控和作业管理有实际意义。 目前, 线性算法偏最小二乘法被广泛用于定性和定量分析[6]。 而PLSR在光谱分析领域应用极为广泛, 因为它可以最大化预测光谱数据与预测变量之间的多重协方差来进行回归预测[7, 8]。 然而在处理一些更复杂的样本集时, 非线性回归方法具有更独特的优势。 机器学习的支持向量回归机具有小样本回归学习能力强的特点, 能够解决小样本非线性且高维度的数据难题, 是小样本集非线性回归问题的不错选择。 另一方面, 近红外光谱会受到诸多方面因素的影响, 譬如环境噪声等等, 这将影响到数据样本的信息特征继而影响预测模型的精准度[9, 10]。 在近红外光谱中会存在一些由于光与土壤不同化学成分之间相互反应而导致的重叠波段, 使红外光谱中的特征波长与SOM之间的相关性可能会被影响。 解决这些问题的最佳办法是使用相关预处理策略来尽可能地消除或是减少环境噪声, 突出需要关注的频谱特征, 并为有用频谱信息的提取做好准备[11, 12]。 对于SVR模型, 我们还可以使用差分进化算法对模型的核函数进行超参数优化, 这将大幅度地提高预测模型的泛化能力和精度。 沈强等[13]分析了土壤养分与光谱反射率之间的相关性, 并在此基础上对支持向量机和多元回归进行建模, 但在支持向量机的超参数选择方面未进行算法优化。 郭李娜等[14]利用网格优化算法和交叉验证对支持向量机的高斯函数超参数进行优化, 成功预测了地表土壤的容重。 然而, 该方法只能搜索有限的离散参数空间, 离散化要求高, 计算量大。 群体智能优化算法不仅可以搜索连续空间, 而且计算速度快, 可以对高维、 非线性和复杂问题进行优化。 谭琨等[15]利用粒子群算法优化支持向量机模型对复垦矿区土壤的有机质含量进行预测。 结果表明, 该模型相较于其他模型, 模型精度和稳定性更加优越。 朱成杰等[16]研究了基于土壤湿度低于30%的高光谱数据, 对土壤湿度和温度预测样本进行不同的预处理和相关性分析, 建立了支持向量机和遗传算法优化向量机(GA-SVR)等模型。 实验表明, 经过连续投影法和遗传算法优化的支持向量机模型能较为精准地预测土壤的湿度和温度。

为提高土壤养分预测模型的拟合精度和稳定性能, 针对土壤光谱数据的波段维度高、 数据关系复杂多变等问题, 提出利用支持向量回归机对土壤养分数据进行非线性回归, 并针对大多数波长可能与所研究的目标无关, 因此对全光谱数据进行特征选择, 以滤除与目标无关的波长。 最终构建了基于PCA-DBO-SVR的土壤养分预测, 并进行实际案例的应用, 结果证明该算法模型有效提升了土壤养分预测的精度和效果。

(1)确定土壤光谱预处理方法

采用0.55~1.65阶以0.1阶为间隔的分数阶微分、 一阶微分分别结合多元散射校正组合出14种预处理方法对光谱数据进行预处理, 同时观察数据预处理前后光谱与有机质之间的相关性变化。 并基于不同预处理方法利用线性模型PLSR和非线性模型BPNN分别进行了建模, 最终通过分析对比获得合适的预处理方法。

(2)建立土壤有机质含量高光谱预测模型

针对小样本、 非线性和高维度的预测模型效果不理想等问题, 引入PCA算法对特征变量数量进行选择, 采用非线性回归具有独特优势的支持向量回归机算法, 同时利用相对于传统的中心化算法拥有较高鲁棒性和容错能力的蜣螂优化参数算法对支持向量回归机的高斯核函数参数组合进行优化, 以有效提高模型精度和性能。 最终以广西国有黄冕林场和国有雅长林场为实际案例, 对土壤中有机质养分含量进行预测, 并与SVR、 GWO-SVR、 PSO-SVR和GA-SVR以及DBO-SVR模型的预测结果进行比较, 证明本文提出的PCA-DBO-SVR土壤有机质反演模型的有效性。

1 实验部分
1.1 研究区域和样本采集

研究区域位于广西国有黄冕林场(109° 43'46″E, 24° 37'25″N)和国有雅长林场(106° 16'30″E, 24° 49'30″N), 两者都属于亚热带气候。 区域内气候温和, 年平均温度19 ℃和16.8 ℃。 雨量充沛, 年均降雨量分别为1 750和1 057 mm。 土地利用类型为典型的桉树人工林, 土壤类型以酸性沉积岩发育的红壤为主。 根据中国土壤分类, 红壤起源于泥盆系砂岩和砂岩页岩, 主要用于种植桉树、 杉木和松树。

应用ArcGis10.1网格化分布点方式(图1)在广西壮族自治区柳州市鹿寨县、 百色市田林县的采样区域进行均匀采样, 使用便携式GPS(Garmin GPS 10, 精度< 2 016 m)记录样品点的位置。 共采集206个土壤样本(0~20 cm), 收集后, 立即将土壤样品密封在塑料袋中并带到实验室。 将收集的土壤样品风干并轻轻研磨过滤至小于0.2 mm的粒径。 样品分为两部分: 一部分用于光谱数据扫描, 另一部分用于化学分析。 所有样品的SOM含量均通过重铬酸钾氧化-外加热法测定[13]。

图1 研究区域采样点分布Fig.1 Sampling point distribution

1.2 土壤样品养分含量的测定

实验室内土壤样本经过0.2 mm的筛子过滤后, 再通过重铬酸钾氧化加热以测定土壤中的SOM, 样本SOM的描述性统计分析如表1所示。 峰度描述数据分布形态陡缓程度, 而偏度则是描述数据分布的偏斜程序, 变异系数用于衡量一组数据的离散程度, 变异系数越小, 说明数据的离散程度越低, 其均值的代表性就越好。 在表1中总样本集的SOM为4.26~80.04 g· kg-1, 平均值为25.80 g· kg-1, 标准差为13.74。 表1中的数据表明了数据在分布上有较低的峰值。 同时数据的偏度均为正数, 即说明采集的数据呈现正偏态分布。 变异系数为53.76%, 表明数据属于中等程度变异。

表1 研究区域样本SOM统计特征 Table 1 SOM statistical characteristics of samples
1.3 数据处理方法

1.3.1 SPXY数据划分(sample set partitioning based on joint x-y distance, SPXY)

SPXY算法是基于KS数据集划分的改进算法。 它可以在不排序的情况下对数据集进行分区, 使其更方便处理非连续变量[17]。 基于联合X-Y距离(SPXY)的样本集划分, 将剩余样本分为训练集和预测集。 其中, 训练集用于训练模型, 预测集用于评估模型性能。 SPXY方法不仅可以考虑样品之间的距离, 还可以平衡光谱变量和内容变量对除法结果的影响。 也就是说, 它可以同时确保差异和平衡。 其基本原理和应用方法可参考相关文献[18]。

1.3.2 分数阶微分(fractional order derivative, FOD)

用FOD算法获得非整数阶微分。 研究表明, 对得到的原始光谱进行分数阶微分变换可以增强光谱差异, 减少光谱的损失[19]。 根据定义方式, 常用的分数阶微分变换算法包括Grunwald-Letnikov(G-L)、 Liouville和Caputo。 G-L算法作为一种离散定义方法, 使用加权和表示函数, 与其他两种算法相比, 该函数更易于计算, 因此选择它进行分数阶微分变换的计算。 在区间[a, b]中的f(θ )的o阶导数公式为

dof(θ)=limh→01ho∑m=0(b-a)/h(-1)mΓ(o+1)m!Γ(o-m-1)f(θ-mh)(1)

式(1)中, h是步进间隔, 将h设置为0.1, (b-a)/h是(b-a)/h的整数部分。 式中Γ 的公式为

Γ(x)=∫0∞exp(-u)ux-1du=(x-1)!(2)

分数阶微分变换的阶数o以0.1的间隔从0.55到1.65变化。 o=0、 1, 分别表示原始未处理的反射光谱和一阶微分变换处理后的反射光谱。

1.4 蜣螂优化算法

蜣螂优化算法(dung beetle optimizer, DBO)是2022年由沈波等提出的一种全新全局优化群体智能算法[20]。 该算法的灵感来自于蜣螂在日常生活中的行为, 例如滚动、 跳舞、 觅食、 偷窃和繁殖。 并相应地设计了5种不同的更新规则。 DBO算法包括滚球、 繁殖、 觅食和偷窃这四个蜣螂的日常活动, 主要思想是将每只蜣螂视为给定搜索空间中可能的可行解。 根据设计的不同变化规律, 朝着更大的适应函数值和实时更新位置的趋势不断迭代, 最后输出最佳位置。 方法的实现总结如下:

(1)对于给定规模的蜣螂种群, 将所有蜣螂按分工划分为滚动蜣螂、 跳舞蜣螂、 觅食蜣螂、 繁殖蜣螂和偷窃蜣螂, 并在搜索空间中随机初始化每只蜣螂的位置和相应的适应值。

(2)每次迭代后, 不同的蜣螂根据各自的更新规则更新自己的位置, 并判断目标位置参数是否超出边界, 同时比较所有蜣螂的适应值, 以记录当前迭代时间下最优蜣螂的相关信息。

(3)重复上述步骤, 直至满足终止标准, 输出全局最优蜣螂个体信息。

此外, 假定光源强度也会影响蜣螂路径, 蜣螂在滚动过程中的位置更新公式如式(3)和式(4)

xi(t+1)=xi(t)+n×m×xi(t-1)+b×X-(3)

X-=|xi(t)-Xc|(4)

式(3)和式(4)中, t为当前迭代次数, xi(t)为第i个蜣螂在第t次的迭代中的位置信息, m表示一个方向偏转系数, 其取值范围为[0, 0.2], b为一个取值为[0, 1]之间的常数项, n为一个取值为-1和1的常数, Xc为全局最差的位置, X-为模拟光亮度的强度变化。

在式(4)中, 如何选取参数m和参数b对整个优化算法的结果是非常重要的; 根据概率方法来确定参数n, 这个参数模拟了自然因素导致蜣螂偏移原始方向, -1代表偏离了既定方向, 1则代表没有偏离方向。 同样地, 参数所模拟的光源强度越大则表示光源越弱, 此时在优化过程中便能较大概率地探索整个搜索空间, 避免陷入局部最优的情况出现。

在蜣螂滚球过程中受到阻碍无法前进时, 会通过跳舞行为重新定位, 并选择前进的路径, 在式(5)中采用切线函数来计算前进角度, 若蜣螂重新锁定前进路径后就会继续进行滚球行为。 跳舞行为的更新公式如式(5)

xi(t+1)=xi(t)+tan(θ∈[0, π])×|xi(t)-xi(t-1)|(5)

在式(5)中|xi(t)-xi(t-1)|表示了第i个蜣螂在第t次和第t-1次迭代之间的位置差异, 因此当前的历史信息对蜣螂位置的更新影响较大。

为提供给后代安全的生活环境, 蜣螂会将粪球推到安全的洞穴进行产卵。 此时为了模拟寻找产卵位置, 需要定义产卵的上下界变量。

Lb* =max(X* ×(1-R), Lb)(6)

Ub* =max(X* ×(1+R), Ub)(7)

式(6)和(7)中, X* 为当前局部最优位置, Lb* 和Ub* 分别为产卵的上界和下界, R=1-(t/maxiter), maxiter表示最大迭代次数。 Lb和Ub分别表示优化的上界和下界。

在寻得产卵位置后蜣螂就会进行繁殖行为, 在每次迭代中每个蜣螂都会产出一颗蛋, 且产卵位置的边界也会由R的变动而不断变化。 而此时的粪球也就是幼雏球的位置更新公式如式(8)

xi(t+1)=X* +b1×(xi(t)-Lb* )+b2×(xi(t)-Ub* )(8)

式(8)中, xi(t)为第i个幼雏球在第t次迭代时的位置信息, b1和b2为两个独立的大小为1× dim维的随机向量, dim为搜索空间的维度, 并且产卵的位置会被严格限制在一定范围内。

一些成年蜣螂会从地下出来寻找食物, 这需要建立起最优的蜣螂区域来对蜣螂觅食进行引导, 最佳的位置区域的边界的上下界变量如式(9)和式(10)

Lbb=max(Xb×(1-R), Lb(9)

Ubb=max(Xb×(1+R), Ub)(10)

式(9)和式(10)中, R为当前全局最优的觅食位置, Lbb和Ubb分别为最优觅食区的下界和上界。 蜣螂觅食的位置更新公式如式(11)

xi(t+1)=xi(t)+c1×(xi(t)-Lbb)+c2×(xi(t)-Ubb)(11)

式(11)中, xi(t)为第t次迭代时第i个觅食蜣螂的位置, c1则表示一个服从正态分布的随机数字, c2表示一个范围在[0, 1]的随机向量。

在蜣螂群体中会有一些蜣螂偷窃其他蜣螂的粪球, 偷窃蜣螂在优化的过程中会不断地更新位置, 最后输出最优的位置Xb。 此时偷窃蜣螂的位置更新公式如式(12)

xi(t+1)=Xb+S×g×(|xi(t)-X* |)+|xi(t)-Xb|)(12)

式(12)中, S表示恒定值, g是一个大小为1× dim维的随机向量, 服从正态分布规律。

本研究试图改进现有SVR模型, 利用蜣螂优化算法(DBO)进一步提高模型精度和稳定性。 目前DBO-SVR已经应用于多个方面, 但在土壤高光谱建模研究中尚未见报道。

1.5 模型评估

以决定系数(R2)、 均方根误差(RMSE)和相对分析误差(RPD)来评价模型的精度。 决定系数表示土壤有机质含量实测值与模型预测值之间的相关性, 其值越接近1, 则表示模型的精度越好。 公式为

R2=1-∑i=1n(ti-pi)2∑i=1n(pi-t-)2(13)

RMSE表示实际土壤有机质含量和预测值的偏差, 其偏差值越低则表示预测效果越好, 模型的精度越高。 公式如式(14)

RMSE=∑i=1n(ti-pi)2n(14)

相对分析误差(RPD)表示模型预测的能力, 用以评价模型精度和可靠性, 其值越大表示模型预测性能越好, 一般的标准: 当RPD< 1.4, 认为模型不可靠; 当1.4< RPD< 2, 认为模型较可靠; 2< RPD< 3认为模型具备较高可靠性, 可用于模型分析; 3< RPD, 认为模型非常可靠, 可用于实际应用。 公式如式(15)

RPD=∑i=1n(ti-t-)2n-1RMSE(15)

式(13)— 式(15)中, ti为土壤有机质含量的实际测量值, pi为土壤有机质含量的模型预测值, t-为土壤有机质含量的实际测量值的平均值。

1.6 实验室光谱与预处理

土壤样品的光谱数据由ASD FieldSpec4 Hi-Res地面物体光谱仪(Analytical Spectral Devices, Inc., Boulder, CO, USA)测量。 为了避免外部光源的干扰, 在暗室中进行测量[15]。

光谱可测范围包含可见光和近红外区域(350~2 500 nm), 分辨率为1 nm, 探头的场角为15° , 50 W卤素灯作为光源, 入射角为45° 。 样本被装入直径为7 cm深度为1.5 cm的样本杯中, 用尺子将样品表面压平压紧, 探头与土壤表面之间的距离设置为1.5 cm。 每10次测量后, 使用标准化白板(25 cm× 25 cm, 100%反射率)校准光谱反射率。 每个样品测量10次; 去除异常光谱, 取算术平均值作为样品的光谱反射率。

考虑到近红外光谱数据在光谱两端(350~399和2 401~2 500 nm)产生较多的频率噪声, 将这些数据范围进行删除。 此外, 对400~2 400 nm的光谱数据以10 nm间隔进行下采样以减小维数和光谱混叠, 提高数据运行速度。 通过重采样获得的201个波段被定义为原始反射光谱。

如图2所示, 林场采集重采样后的206个土壤样品的光谱反射率在0至0.9之间。 可以明显看出, 400~2 400 nm波段范围的土壤样本光谱反射率曲线在总体上的变化形态几乎类似, 大体的变化趋势比较平缓。 而在不同的波长范围内, 光谱反射率的大小随着波段的变化而变化。 在400~700 nm范围内, 所有样品的光谱曲线相对集中, 反射率普遍较低。 700~2 400 nm光谱反射率开始呈现发散趋势, 在经历700~900 nm的反射率快速增长后, 随后的光谱反射率出现了缓慢增长甚至局部下降的现象。 在1 400、 1 900和2 300 nm波段附近有三个明显的吸收谷。

图2 土壤样品光谱曲线图Fig.2 Spectra of soil samples

以SOM每10 g· kg-1进行分组, 取其光谱反射率平均值进行绘图, 观察图3中不同SOM含量对反射率的影响, 我们发现这些光谱曲线的形态特征基本相同。 在400~1 320 nm波段范围内当SOM较低(小于30 g· kg-1)时, 随着SOM的增加, 反射率下降更加明显。 在可见光波段400~780 nm有较大的曲线斜率。 在900 nm附近存在与氧化铁相关的吸收谷, 在近红外波段的780~2 400 nm, 不同SOM水平在900、 1 430、 1 950和2 250 nm处均有明显的吸收峰。 在480~1 860 nm波段范围内当SOM为30~90 g· kg-1时, 其包含的6组曲线出现了交叉现象, 且在这波段内有部分样本的反射率高于1 860~2 400 nm波段内SOM为60~70和70~80 g· kg-1的光谱反射率。 故我们可以得出初步结论, SOM与光谱反射率大小有关系, SOM越大, 光谱反射率越小。 SOM越小, 光谱反射率越大。

图3 不同SOM土壤样本的光谱曲线Fig.3 Spectral curves of soil samples with different SOM values

采用比一阶微分处理数据更加精细且能平衡光谱噪声和光谱分辨率之间的关系的FOD来降低噪声, 提高信号的峰值信噪比。 同时结合MSC来消除滤波后对材料表面的散射影响, 最终提高数据的精度和准确性, 使得实验结果更加可靠和可信。

图4为在区间为0.5~1.6阶, 间隔为0.05的不同FOD光谱反射率图。 由于图4(a)原始光谱曲线具有吸收特性不明显、 波段之间相关重叠、 光谱信号弱和波峰波谷不明显以及反射率范围宽等特点, 所以通过FOD处理光谱能更好地捕捉丰富的光谱信息。 图4(b)和图4(c)分别为0.5阶和0.55阶的光谱数据曲线图, 曲线没有出现显著的吸收特征变化。 从0.6阶开始, 光谱反射率首次出现了负值, 且随着阶数的增加, 反射率逐渐趋向于0, 这均说明随着微分变换的阶数增加, 光谱曲线的混合重叠峰和基线漂移问题正在逐渐被消除。 从图中可以看出FOD为一阶时, 大多数光谱反射率逐渐趋向于0, 而随着FOD的阶数越来越大[图4(l— x)], 反射率的光谱信息也得到了增强, 同时特定的吸收特征也消失了。 当FOD阶数超过一阶时, 光谱的反射率基本趋向于0, 图像的基本信息几乎丢失, 曲线相对于FOD之前产生了较为明显的上下波动, 局部反射率得到了不同程度的放大, 出现了多个峰值波段。 通过以上对比分析发现分数阶微分变换较一阶微分变换更能捕捉到光谱数据的潜在信息, 从而提高模型的预测精度。

图4 FOD变换后光谱反射率曲线图
(a)原始光谱反射率; (b): 0.5阶; (c): 0.55阶; (d): 0.6阶; (e): 0.65阶; (f): 0.7阶; (g): 0.75阶; (h): 0.8阶; (i): 0.85阶; (j): 0.9阶; (k): 0.95阶; (l): 1阶; (m): 1.05阶; (n): 1.1阶; (o): 1.15阶; (p): 1.2阶; (q): 1.25阶; (r): 1.3阶; (s): 1.35阶; (t): 1.4阶; (u): 1.45阶; (v): 1.5阶; (w): 1.55阶; (x): 1.6阶Fig.4 Plots of spectral reflectance after FOD transformations
(a): Original; (b): 0.5-order; (c): 0.55-order; (d): 0.6-order; (e): 0.65-order; (f): 0.7-order; (g): 0.75-order; (h): 0.8-order; (i): 0.85-order; (j): 0.9-order; (k): 0.95-order; (l): 1-order; (m): 1.05-order; (n): 1.1-order; (o): 1.15-order; (p): 1.2-order; (q): 1.25-order; (r): 1.3-order; (s): 1.35-order; (t): 1.4-order; (u): 1.45-order; (v): 1.5-order; (w): 1.55-order; (x): 1.6-order

在图4中利用微分变换预处理探测并放大特征波段信息, 图5则利用皮尔逊相关系数法分析了土壤中有机质养分在不同阶数的FOD+MSC处理后的光谱相关系数。 如图5(a)所示原始光谱数据的反射率和有机质含量在整体上呈现负相关性。 曲线在530~1 300 nm波段有良好的相关性, 在此波段范围内的相关系数的绝对值基本大于0.65, 在620 nm处达到最大相关系数。 在2 130~2 400 nm的波长范围内, 相关系数呈下降趋势, 在2 310 nm处出现一个反射峰。 如图5(b— o)所示, 在不断增加分数阶微分变换阶数的同时, 光谱数据与有机质之间的相关系数也在不断地显示出来, 通过这些变换曲线可以大概地看出阶次与相关性之间表现出了正负相关的关系。 而相关系数的走势曲线由初始光谱数据时的较为平滑慢慢演变成局部尖锐的情况。 图5中我们发现与预测值有着高相关性的波段在数量上呈现出先增加后减少的情况, 说明微分变换的阶次数并不是越高越好, 阶次越高相关性则出现降低的趋势。 当光谱的波段处于610 nm时, 阶次为0.55微分变换结合多元散射校正所呈现出的相关性系数是最高的, 其值达到了0.764。 而当不断增加微分阶次时, 我们发现光谱反射率会越来越趋近于0, 此时对应的相关性系数也开始降低。 此外, 整数阶微分的一阶微分在结合多元散射校正后的光谱数据只有4个波段的绝对相关性系数大于0.65, 它们分别是-0.678、 -0.697、 -0.689、 -0.657, 对应540、 560、 550和570 nm的波段。 由此我们得出一些结论, 结合多次散射校正的一阶微分变换对原始光谱信息中的细微变化缺乏处理, 这则很有可能漏掉光谱信息中的一些有效特征信息。 而作为对比, 更为细致的分数阶微分则可以以微小的阶次量来进行移动处理, 从而从光谱的斜率曲线中获得尽可能多的有效信息。

图5 不同阶数的FOD+MSC处理后的光谱相关系数
(a): 原始光谱; (b): 0.55-FOD+MSC; (c): 0.65-FOD+MSC; (d): 0.75-FOD+MSC; (e): 0.85-FOD+MSC; (f): 0.95-FOD+MSC; (g): 1-FOD+MSC; (h): 1.05-FOD+MSC; (i): 1.15-FOD+MSC; (j): 1.25-FOD+MSC; (k): 1.35-FOD+MSC; (l): 1.45-FOD+MSC; (m): 1.55-FOD+MSC; (n): 1.65-FOD+MSCFig.5 Spectral phase relationships after FOD+MSC processing with different orders
(a): Original; (b): 0.55-order FOD+MSC; (c): 0.65-order FOD+MSC; (d): 0.75-order FOD+MSC; (e): 0.85-order FOD+MSC; (f): 0.95-order FOD+MSC; (g): 1-order FOD+MSC; (h): 1.05-order FOD+MSC; (i): 1.15-order FOD+MSC; (j): 1.25-order FOD+MSC; (k): 1.35-order FOD+MSC; (l): 1.45-order FOD+MSC; (m): 1.55-order FOD+MSC; (n): 1.65-order FOD+MSC

1.7 软件

所有算法均配备AMD Ryzen 5 5600X CPU、 16 GB RAM、 在Windows 10操作系统下运行PyCharm Community Edition 2022.3.3的计算机上实现。

2 结果与讨论
2.1 样本划分

利用SPXY算法对数据集按照4∶ 1比例进行划分, 其中165个样本作为训练集, 其余41个样本作为预测集, 如表2所示, 划分后SOM的最大最小值都被划分到了训练样本中, 预测样本的数据值都在测试样本的SOM区间内, 从训练集和预测样本的变异系数来看均属于中等程度变异。 峰度(Kurtosis)是一种统计量, 用来描述与正态分布相比, 分布的尖峰程度或平坦程度。 如果峰度为0, 则表明该分布与正态分布具有相同的峰值水平。 峰度大于0表示该分布比正态分布的峰值更大, 也就是说它的峰值更尖锐。 相反, 峰度小于0表示分布比正态分布更平坦, 因此峰值更平坦。 通常, 大于3的峰度表示峰值比正态分布的峰值更尖锐、 更陡峭。 偏度与峰度类似, 是一种用于描述数据分布形状的统计量。 它衡量数据分布的对称性。 与正态分布相比, 偏度为0表明该分布与正态分布具有相同的对称性。 在表2中训练样本和预测样本的峰度均小于3, 偏度均为正数, 以上表明划分后的数据在分布上具有较低的峰值, 而每个样本的数据分布也符合正偏态分布, 故SPXY处理后的数据集离散程度合理。

表2 训练样本和预测集样本统计特征 Table 2 Statistical characteristics of training and prediction samples
2.2 基于全波段光谱不同预处理方法的的偏最小二乘回归模型反演分析

基于不同阶数的FOD+MSC的全波段偏最小二乘回归模建模结果如表3所示, 验证集的R2的范围是[0.811, 0.899], RMSE的范围是[3.849, 5.781], RPD的范围是[3.1, 4.4]。 从表3中数据的对比可以发现, 基于1.25阶FOD光谱的PLS模型预测精度最佳, 验证集的R2、 RMSE和RPD分别为0.899, 3.849和3.2。 图6为实测值和预测值的散点图, 1.25阶FOD最靠近1∶ 1线, 而此时基于1阶的FOD+MSC的PLSR模型的验证集的R2、 RMSE和RPD分别为0.871, 4.390和3.3。 并且分别基于一阶微分变换和MSC的PLSR模型的验证集效果均没有FOD+MSC组合算法的预测精度高。 综上所述, 在线性模型PLSR中, 分数阶FOD+MSC算法的拟合效果优于一阶FOD+MSC。

表3 基于光谱全波段的PLSR模型构建及验证 Table 3 Construction and validation of PLSR model based on full band spectra

图6 偏最小二乘回归模型的有机质实测值与预测值比较
(a): PLSR (RAW); (b): PLSR (FOD(1)); (c): PLSR (FOD(1)+MSC); (d): PLSR (FOD(1.25)+MSC)Fig.6 Comparison of measured and predicted values of organic matter obtained by partial least squares regression models
(a): PLSR (RAW); (b): PLSR (FOD(1)); (c): PLSR (FOD(1)+MSC); (d): PLSR (FOD(1.25)+MSC)

2.3 基于光谱全波段和不同预处理方法的反向传播人工神经网络模型反演分析

基于不同阶次的FOD+MSC组合预处理法构建反向传播人工神经网络模型的结果如表4所示。 表中显示验证集的决定系数为0.808~0.917, 同时均方根误差为3.351~5.266 g· kg-1, 最后相对分析误差为2.3~3.7。 综合比较后基于1.25阶的FOD+MSC的反向传播人工神经网络模型建模效果最佳, 验证集的R2、 RMSE和RPD分别为0.917, 3.351 g· kg-1和3.6。 而基于一阶FOD+MSC的反向传播人工神经网络模型的验证集的R2、 RMSE和RPD分别为0.902, 3.816 g· kg-1和3.2。 图7用y=x函数来比较预测值和实测值。 结果表明, 在非线性模型BPNN中, 基于分数阶FOD+MSC的BPNN估测模型的精度和稳定性优于一阶的FOD+MSC。

表4 基于光谱全波段的BPNN模型构建及验证 Table 4 Construction and validation of BPNN model based on full band spectra

图7 反向传播人工神经网络模型的有机质实测值与预测值比较
(a): BPNN (RAW); (b): BPNN (FOD(1)); (c): BPNN (FOD(1)+MSC); (d): BPNN (FOD(1.25)+MSC)Fig.7 Comparison of measured and back propagation Ann model prediction values for organic matter
(a): BPNN (RAW); (b): BPNN (FOD(1)); (c): BPNN (FOD(1)+MSC); (d): BPNN (FOD(1.25)+MSC)

2.4 光谱预处理在不同模型反演分析

为探究不同阶数的FOD+MSC预处理方法在线性和非线性模型的建模效果, 选取了线性PLSR和非线性BPNN两种模型进行对比分析。 通过比较不同阶次的光谱微分变换结合多种散射校正预处理方法, 分别采用PLSR和BPNN模型对土壤有机质含量进行预测。 结合图8(a)和(b)可以看出, 不同阶数的光谱预处理组合, 较原始光谱而言, PLS和BPNN两种模型的决定系数分别提高5.7%~24.4%和4.3%~13.4%, 均方根误差分别降低了2.7%~3.7%和3.1%~3.3%, 相对分析误差分别降低12%~38%和提升了8.7%~60.9%。 总体上看, 分数阶的FOD+MSC处理后的模型性能比原始光谱和一阶微分构建的预测模型有着更好的精度和性能效果。 同时通过预测值和实测值散点图对比可以发现, 分数阶微分结合多元散射校正能够更好地体现出光谱数据的有效信息。 如图6和图7所示, 当SOM含量小于30 g· kg-1时, 多数样本的模型预测结果偏高。 当SOM含量大于30 g· kg-1时, 多数样本的模型预测结果偏低。 如图3所示, SOM含量大于70 g· kg-1样品的反射率远高于SOM含量较小的样品, SOM含量的较大差异带来的光谱反射率数值的较大跨度, 可能导致线性模型的预测不准确。 因此, 在实验中, 在不同的光谱预处理方法上, 非线性模型BPNN优于线性模型PLSR。 如图8所示, 当SOM含量小于10 g· kg-1时, 模型对大多数样本的预测结果趋于较高。 当SOM含量小于30 g· kg-1时, 多数样本的模型预测结果偏高。 当SOM含量大于30g· kg-1时, 多数样本的模型预测结果偏低。 结果表明, 1.25阶的FOD+MSC的预处理方法搭配非线性预测模型, 相比线性模型的拟合效果更好。 这也与Hong[21]和Gao[22]等采用的分数阶微分组合预处理方法作为光谱预处理, 导致相关性和系数显著提高的结论相似。

图8 验证集不同预测模型的评价指标对比
(a): 决定系数; (b): 均方根误差; (c): 相对分析误差Fig.8 Comparison of evaluation metrics of different prediction models for the validation set

2.5 基于PCA降维的土壤有机质反演模型构建

选择相关性高的敏感波段和最优数量的主要成分对在回归中保留住模型的性能具有重要意义。 选择非最佳数量的因素会导致信息丢失, 而超过最佳数量的因素会也会让模型的预测能力变差。

2.5.1 确定PCA最佳维度数

如图9所示, 在基于1.25阶的FOD+MSC预处理方法下, 不同建模方法的相对分析误差伴随维度数的变化表现出明显的区别。 当维数为61时, DBO-SVR模型反演土壤有机质的RPD值在测试集达到最大值。 因此, 在后续DBO-SVR方法构建土壤SOM反演模型时, 选择61个维数作为PCA降维的最佳维数。

图9 PCA降维过程中的维度数Fig.9 The number of dimensions in PCA dimension reduction

2.5.2 基于最优波长的蜣螂算法优化支持向量回归(PCA-DBO-SVR)

使用DBO算法进行SVR参数搜索的建模过程如图10所示, 为了建立理想的回归模型, 采用DBO算法对SVR模型中的参数进行微调。 DBO算法参数使用以下设置: 最大迭代限制为110, 总体为20, 下限为0.01, 上限为1 000。 此外, 使用1-R2作为适应度函数, 其中R2是决定系数。 适应度函数的值越小, 判定系数R2的值就越大。 有机质含量被作为预测值输出, 算法被用于训练模型。

图10 DBO-SVR流程图Fig.10 Flow chart of DBO-SVR

如图11所示, 以PCA降至61维后的特征波段数据为例, 当进化世代数达到11时迅速收敛, 在迭代至87次达到稳定, 获得模型的最佳适应度和C、 g参数的最佳组合, C和gamma分别为1 000.0和88.927 6。 此时, 此时适应度值为0.358 85, MSE为8.764 3。 经过DBO优化后的SVR模型测试集R2达到0.942, RMSE为2.989, RPD为4.1, 结合表5的分析结果看, 对比SVR模型、 PSO-SVR、 GA-SVR和GWO-SVR以及DBO-SVR等反演模型, PCA-DBO-SVR模型的训练集和测试集的实测值与预测值点更为接近, 如图12所示, 模型具有更为准确的预测效果。

图11 蜣螂算法适应度曲线Fig.11 DBO fitness curve

表5 不同变量筛选方法的土壤SOM反演模型构建与验证 Table 5 Construction and validation of soil SOM estimation models with different variable screening methods

图12 DBO-SVR模型结果Fig.12 DBO-SVR model results

2.6 不同建模方法的结果分析

通过分析不同阶数的FOD+MSC预处理方法对模型精度的影响, 结果证明1.25阶微分变换处理后的多重散射校正是最佳的预处理方法。 以此为基础, 构建了SVR、 GWO-SVR、 PSO-SVR、 GA-SVR以及DBO-SVR模型, 并进行了比较实验, 实验结果如表5所示。

综合上来看, 全波段建模, 土壤SOM反演模型精度排序为PCA-DBO-SVR> DBO-SVR> PSO-SVR> GWO-SVR> GA-SVR, 采用PCA-DBO-SVR模型反演效果最佳(RPD=4.1)。 SVR在训练集和测试集的各项评价指标均不如其他五种优化后的模型, 测试集的模型效果也是六种模型中最差的, 其问题就是支持向量机中的高斯核函数超参数组合没有得到优化处理, 直接影响了模型最终的拟合精度。

从优化原理上看, GA-SVR模型的决定系数和稳定性都低于其他优化后的SVR模型, 而且训练时间也是所有模型中最久的, 这由于遗传优化算法极易陷入局部最优, 导致其优化后SVR模型的泛化性能也表现较差。

PSO-SVR和GWO-SVR作为经典的群体智能优化算法由于具有记忆性, 粒子群的历史最佳位置可以被记住并与其他粒子共享。 此外, 自适应收敛因子和信息反馈机制的存在, 允许在局部优化和全局搜索之间取得平衡。 因此, 提高了模型的预测精度和收敛速度。

PCA-DBO-SVR模型的决定系数最高为0.942, RPD超过了4.0, 其他模型的RPD均小于4.0, PCA-DBO-SVR模型的预测效果和拟合精度较DBO-SVR和其他模型有显著的提升, 充分说明PCA-DBO-SVR能更好地提高土壤预测模型的精度和预测能力。

从训练集和测试集的模型效果可以看出, PCA-DBO-SVR模型的泛化能力较其他模型好。 从测试集的模型效果来看, 其中PCA-DBO-SVR模型的均方根误差最低且相对误差最高, 即模型的拟合效果和可靠性在所有模型中是最优的。 在表5中, 所有模型的预测能力都有较好的表现, 所有优化后模型的决定系数均大于0.92, 均方根误差小于3.5 g· kg-1, 相对分析误差大于3.5。 而与DBO-SVR模型相比PCA-DBO-SVR模型的决定系数R2与相对误差RPD分别升高了1.4%和7.9%, 与SVR模型相比PCA-DBO-SVR模型的决定系数R2与相对误差RPD分别升高了6.2%和36.70%, 说明了PCA-DBO-SVR较DBO-SVR的优越性。

3 结论

近年来, 精准林业不断发展, 对土壤进行测定和分析后可以根据不同树种的养分需求制定合理的施肥方案, 以提高林木生长和产量的研究工作, 因此利用光谱快速且准确地预测土壤有机质含量来计算有机碳至关重要。 通过观察大量相关文献, 构建了基于PCA-DBO-SVR 的土壤养分预测模型并通过广西国有黄冕林场和国雅长林场采集土壤样品进行实例验证。 考虑到土壤环境中有机碳是土壤有机质的碳部分, 因此在试验验证以土壤中的有机质含量作为预测对象。 所有实验过程使用PyCharm编译软件建模仿真。 最终仿真试验结果表明, 相比于其他常用的回归预测模型, 本文所提出的PCA-DBO-SVR的土壤有机质预测模型在精度和泛化性能上表现较好, 具备实际使用价值。 结论如下:

(1)对比了不同阶数的分数阶微分和一阶微分变换分别结合多元散射校正的预处理方法组合基于线性和非线性模型的预测能力。 采用1.25阶微分变换和多元散射校正的预处理组合后, 线性模型 PLSR与非线性模型BPNN两种模型的R2分别提高了19%~24.2%和12.8%~13.5%, 说明分数阶微分变换克服了传统一阶导数对细微变化不敏感的缺点。 这是一种增强土壤养分含量与光谱反射率相关性的有效预处理方法。 此外, 实验结果表明, 非线性模型的整体建模性能优于线性模型。 同时在对比有机质实测值后可以发现非线性模型的样本预测值在0~10 g· kg-1的有机质含量范围内均小于实测值, 在20~30 g· kg-1 范围内均大于实测值, 在超出30 g· kg-1 范围外大部分均小于实测值。

(2)在不同优化算法优化支持向量机模型的对比中, PCA-DBO-SVR模型相较GA-SVR、 GWO-SVR和PSO-SVR以及DBO-SVR模型验证集的决定系数分别提高了2.4%、 2.1%、 1.9%和1.4%; 同时均方根误差分别降低了0.553、 0.393、 0.390和0.368; 相对分析误差分别提高了17.1%、 13.9%、 10.8%和7.9%; 充分说明PCA-DBO-SVR 模型能更好提高土壤有机质反演模型的精度和预测能力。

(3)与其他预测模型进行性能比较后, 仿真结果表明, PCA-DBO-SVR 模型具有更好的预测精度和泛化能力, 验证集的决定系数R2是0.942, 均方根误差是2.989 g· kg-1, 相对分析误差是4.1。 实验证明了本文提出的PCA-DBO-SVR土壤养分预测模型适用于实际土壤有机质养分含量的预测, 拥有较好的实践价值。

在PCA-DBO-SVR模型的实际案例应用中, 利用分数阶微分耦合其他预处理方法在预测模型上显示出了更好的精度和稳定性, 并且支持向量机在小样本数据中表现优秀, 预处理方法又进一步提高了光谱反射率和土壤养分含量的相关性。 虽提出的算法在土壤养分预测试验中有着较好的表现, 但仍需进一步完善, 可以改进的几点如下:

(1)本实验数据集为土壤光谱的一维光谱反射率数据, 其数据信息的丰富度范围具有局限性。 如果可以收集到更高维度的光谱图像类数据, 进一步丰富输入变量的信息量范围, 则能从理论和实践上实现区域内土壤相关养分的精准预测。

(2)为更好地提高光谱数据预处理的效果, 可利用一些相关的数学方法处理数据, 同时搭配相适应的光谱特征选择方法, 进一步组合或改进出拟合效果更好的模型预测效果。

(3)土壤环境组成复杂, 本实验只验证了土壤有机质养分的预测能力, 未来可以考虑采用更多土壤养分数据, 以建立更有实际应用价值的土壤养分预测模型。

参考文献
[1] Mcbratney A, Field D J, Koch A. Geoderma, 2014, 213: 203. [本文引用:1]
[2] Piikki K, Söderström M. Validating and Adapting Regional Soil Maps for Precision Agriculture, Conference on Proximal Sensing SupportingPrecision Agriculture, 2015. [本文引用:1]
[3] Dotto A C, Dalmolin R S D, Caten A T, et al. Geoderma, 2018, 314: 262. [本文引用:1]
[4] O'Rourke S M, Stockmann U, Holden N M, et al. Geoderma, 2016, 279: 31. [本文引用:1]
[5] Guo L, Zhao C, Zhang H, et al. Geoderma, 2017, 285: 280. [本文引用:1]
[6] Chen X, Xu Y, Meng L, et al. Sensors and Actuators B: Chemical, 2020, 311: 127924. [本文引用:1]
[7] Xie Z, Feng X, Chen X. Chemometrics and Intelligent Laboratory Systems, 2022, 221: 104486. [本文引用:1]
[8] Su L, Shi W, Chen X, et al. Food Chemistry, 2021, 338: 127797. [本文引用:1]
[9] Stenberg B, Rossel R A V, Mouazen A M, et al. Advances in Agronomy, 2010, 107: 163. [本文引用:1]
[10] Rossel R A V, Behrens T, Ben-Dor E, et al. Earth-Science Reviews, 2016, 155: 198. [本文引用:1]
[11] Gao Y, Cui L, Lei B, et al. Applied Spectroscopy, 2014, 68(7): doi: DOI:10.1366/13-07031. [本文引用:1]
[12] Vasat R, Kodesova R, Klement A, et al. Geoderma, 2017, 298: 46. [本文引用:1]
[13] SHEN Qiang, ZHANG Shi-wen, XIA Sha-sha, et al(沈强, 张世文, 夏沙沙, 等). Journal of Anhui University of Science and Technology(Natural Science)[安徽理工大学学报(自然科学版)], 2019, 39(4): 39. [本文引用:2]
[14] GUO Li-na, FAN Gui-sheng(郭李娜, 樊贵盛). Chinese Journal of Soil Science(土壤通报), 2018, 49(3): 512. [本文引用:1]
[15] TAN Kun, ZHANG Qian-qian, CAO Qian, et al(谭琨, 张倩倩, 曹茜, 等). Earth Science—Journal of China University of Geosciences(地球科学—中国地质大学学报), 2015, 40(8): 1339. [本文引用:2]
[16] ZHU Cheng-jie, WANG Zheng-quan(朱成杰, 汪正权). China Rural Water and Hydropower(中国农村水利水电), 2024, (1): 30. [本文引用:1]
[17] CHEN Yi-yun, QI Tian-ci, HUANG Ying-jing, et al(陈奕云, 齐天赐, 黄颖菁, 等). Transactions of the Chinese Society of Agricultural Engineering(农业工程学报), 2017, 33(6): 107. [本文引用:1]
[18] Han L, Chen R, Zhu H, et al. Sustainability, 2020, 12(4): 1476. [本文引用:1]
[19] Zhang Z, Ding J, Wang J, et al. Catena, 2019, 185(10): 104257. [本文引用:1]
[20] Xue J, Shen B. The Journal of Supercomputing, 2022, 79(7): 7305. [本文引用:1]
[21] Hong Y, Liu Y, Chen Y, et al. Geoderma, 2019, 337: 758. [本文引用:1]
[22] Gao L, Zhu X, Han Z, et al. Pedosphere, 2019, 29(6): 810. [本文引用:1]