蚁群算法在土壤速效磷近红外光谱波长选择中的应用

张小鸣　冒智康　李绍稳　金秀　朱正伟

摘要：光谱变量多、光谱信息易重叠、数据存在冗余且存在大量噪声，导致近红外光谱分析建模的时间长、模型的泛化能力不强、预测精度不高等问题。为此，提出了一种基于变量有效性精英蚁群系统（EAS）的波长选择算法，通过蚁群搜索选出特征波长。将变量有效性作为蚁群的初始信息素，在信息素更新过程中引入变量重要性投影方法，加快算法的收敛速度，避免陷入局部最优。为了验证算法的有效性，选取193个沙姜黑土土壤样本的近红外光谱数据为研究对象，采用变量有效性精英蚁群算法选择特征波长，构建土壤速效磷含量偏最小二乘回归预测模型，对全谱偏最小二乘模型、蚁群算法偏最小二乘模型进行比较研究。结果表明，采用变量有效性精英蚁群算法选出的特征波长建立的回归模型精度比普通蚁群算法选出的特征波长和全光谱高，建立回归模型更加简单。

关键词：土壤速效磷;特征变量选择;蚁群算法;变量有效性;精英策略

中图分类号： TP391文献标志码： A

文章编号：1002-1302（2019）19-00 -

收稿日期：2018-07-18

基金项目：农业农村部农业物联网技术集成与应用重点实验室开放基金（编号：2016KL07）。

作者简介：张小鸣（1958—），男，安徽合肥人，博士，教授，主要从事嵌入式系统应用研究。E-mail：mdd2942942@163.com。

近年来，过度耕种、施肥已经严重影响了土壤氮、磷、钾循环及其平衡[1]。快速有效测定土壤氮、磷、钾含量，对指导农业耕种中作物种植和施水施肥有一定的指导意义[2]。光谱检测技术具有分析速度快、可多组分同时测定、非破坏性、低成本和操作简单等特点，已被广泛应用于土壤养分含量检测[3]。但是近红外光谱存在严重的共线性，信噪比低，并且包含大量无关信息，因此特征波长选择成为近红外光谱分析中不可或缺的重要步骤，是光谱分析技术领域的研究热点[4]。近年来，国内外学者提出了大量特征变量选择算法，主要包括基于偏最小二乘回归（PLSR）模型参数的算法，基于连续投影策略的算法，基于波长区间选择的算法，基于群智能优化算法和基于模型集群分析策略等[5]，其中群智能算法成为研究热点[6]。郭志明等使用普通蚁群算法优选预测苹果可溶性固形物含量的近红外光谱特征波长，建立偏最小二乘模型，使得波长变量数减少，模型预测性能提高，决定系数和预测均方根误差分别为 0.97 和0.51[7]。但是这种方法将蚁群初始信息素浓度都设为1，随机性强，容易出现搜索停滞现象。Shamsipur等提出一种蚁群算法进行特征波长选择，但是波长需要人为设定，具有一定的不确定性[8]。基于以上问题，本研究提出基于变量有效性的精英蚁群算法，同时结合前向选择方法优选特征波长，采用精英蚁群搜索，并将变量有效性作为蚁群初始信息素，以避免蚁群陷入搜索停滞，搜索时间过长。最后，通过前向选择实现特征波长自主选择，以期提高土壤速效磷含量分析模型的准确性和稳定性。

1 土壤样本采集与处理

试验所用的193个土壤样本取自安徽省蒙城县、宿州市埇桥区和怀远县，这些地区的主要作物是小麦和玉米，土壤类型为沙姜黑土。根据作物种植密度、施肥浓度梯度的不同选取采样点，以保证样本速效磷含量存在一定差异。将采得的土壤样本封存到实验室进行风干后碾磨，过20目筛。将处理后的土壤样品等分为2份，一份用于标准定量测试，另一份用于光谱数据采集。速效磷含量测定采用碳酸氢钠浸提-钼锑抗分光光度法[9]。使用OFS1700地物光谱仪（海洋光学亚洲公司）和50 W卤钨灯接触式反射探头采集土壤近红外光谱数据，测量时每个样本旋转90°测量6次取平均值作为光谱数据;光谱范围为350～1 655 nm，光谱分辨率为1 nm，共采集波长数1 306个。去除检验值相差较大的异常样本，最终选出188个样本用于建模分析。

2 光谱数据预处理

188个土壤样本在波长为350～1 655 nm处的原始近红外光谱反射率如图1-a所示，可以发现，光谱在400～500 nm处反射率较低，在500～700 nm处明显上升;在 930 nm 处有明显的波动，可能是土壤中铁氧化物（针铁矿、赤铁矿）对光谱吸收引起的[10];在1 400 nm处有明显吸收峰，可能是由于土壤样品中残留水分的影响[11]。光谱图像出现较多毛刺，伴有较多噪声，可能是由于土壤颗粒大小、均匀度不同;光谱仪受高频噪声干扰、基线漂移的影响等，因此需要对光谱进行预处理。Savitzky-Golay（S-G）平滑算法可消除基线漂移、倾斜等引起的噪声。标准正态变量变换（SNV）和多元散射校正（MSC）可减少固体颗粒大小不同、光散射造成的影响。矢量归一化可使数据保持相同标度，减少数值差异对建模结果的影響[11]。根据比尔定律，土壤某养分的近红外光谱吸光度与该养分的浓度成正比，因此将原始反射率光谱图转换成吸光度光谱图，具体方式为将原始近红外光谱反射率数据取倒数再取对数转换为吸光度数据，然后采用S-G平滑算法+SNV对吸光度数据进行预处理，得到188个土壤样本在波长为350～1 655 nm近红外光谱波段内的吸光度光谱图（图1-b）。

3 光谱特征变量选择与建模方法

3.1 校正集与验证集的划分

采用浓度梯度法将样本集按照速效磷含量的参照值进行排序，以3 ∶1的比例划分为校正集与验证集。由表1可知，校正集与验证集有相似的数据分布特征结构，说明选出的校正集与验证集可以代表整体数据集的分布特征。

3.2 精英策略蚁群算法

3.2.1 蚁群算法波长选择原理

蚂蚁在觅食过程中，会在途经线路上分泌信息素，后面的蚂蚁会根据信息素判断觅食方向。蚁群算法通过模拟蚂蚁觅食行为，根据各路径上的信息素沉积，蚁群最终收敛到最优解，此最优解即为特征波长变量。

3.2.2 精英策略蚁群算法

波长变量有效性可以表征变量对于模型预测指标的重要性，其值越大，重要性越强。变量重要性计算方法如下：

采用蒙特卡罗法（MCS）采样，从n个样本中随机抽取部分样本建立偏最小二乘回归（PLSR）模型，得到相应的回归系数向量（βi）、决定系数（R2），经过m次采样，第i个波长变量的有效性计算公式为

从公式中可以看出，回归模型的决定系数越大、均方根误差（RMSE）越小，该变量的有效性越大，重要性越强。

在模拟蚁群过程中，若没有信息素，蚂蚁就会随机选择一个方向，而且不会向着同一方向运动。引入精英蚁群思想并将变量有效性作为蚂蚁初始信息素的浓度选择依据，为蚁群提供初始方向。在信息素更新时引入变量重要性投影（VIP）系数，加快算法的收敛，最后通过前向选择方法选出最终特征波长数量。

传统蚁群算法计算时间长，且容易出现所有蚂蚁求出的解完全一致的现象[12]。基于精英策略的蚁群信息素更新机制是根据重要性分配不同权重，将蚂蚁按权重排序，避免算法陷入局部最优，通过迭代实现全局寻优[13]（图2）。每次迭代将蚂蚁所选波长变量代入偏最小二乘回归模型，计算RMSE，RMSE越小，说明所选的波长变量精英度越高。

基于精英策略的蚁群信息素更新计算公式为

τi（t+1）=（1-ρ）×τi（t）+ρ×Δτ。（2）

式中：τi（t+1）为第i个蚂蚁（第i个波长变量）在t+1时刻更新后的信息素含量;τi（t）为第i个蚂蚁在t时刻的信息素含量;ρ为信息素衰减率常数，取值在0～1之间;Δτ为蚁群信息素增量。

蚁群信息素增量Δτ的计算公式为

Δτ=∑[DD（]mi=2[DD）]ωi×τi+λ×F。（3）

式中：ωi为第i个蚂蚁信息素更新所占的权重;τi（t）为第i个蚂蚁在t时刻的信息素含量;λ是信息素更新系数常量，取值为0～1;F为回归模型的RMSE极小值的倒数。

F和ωi计算公式如下：

F=Q1+minRMSEωi=（1-λ）i。（4）

式中：Q是一个显著因子常量，取值为0～1。若RMSE越小，则F值越大，Δτ值也越大，表示模型预测能力越好，蚂蚁精英度越高。若ωi值越大，则Δτ值也越大，代表精英蚂蚁分泌的信息素越多。Δτ值越大，表示蚁群算法更新速度越快。

为了进一步加快蚁群算法收敛速度，减少计算时间，在基于精英策略的蚁群信息素基础上，再引入VIP，计算公式如下：

VIPj=p∑w2k（SSYa/SSYt）SSY=∑[DD（]ni=1[DD）]（yi-y）2。（5）

式中：VIPj表示第j个波长变量的VIP得分;p表示波长变量数量;wk表示对应的权重系数;SSYa表示用a个潜变量建模时的SSY值;SSYt表示用全部潜变量建模时的SSY值。SSY表示每个土壤样本速效磷含量预测值（yi）与所有土壤样本速效磷含量预测值的平均值 y 之差的平方和;蚁群信息素改进计算公式如下：

τ1=τ×VIP。（6）

式中：τ1为算法改进后的蚁群信息素浓度;τ是蚁群原始信息素浓度;VIP表征波长变量对预测的重要性[14]。当VIP>1时，τ1>τ，信息素被放大，波长变量被选中的概率变大。反之，信息素被缩小，波长变量被选中的概率变小。

最后，将通过精英策略蚁群算法选出的特征波长变量按选择次数从大到小排列，并将排序后的特征波长变量依次代入偏最小二乘交互检验模型计算RMSE，根据采样前向选择（FS）规则进行特征波长变量选择，如果RMSE变小，则选择此变量作为特征波长变量，否则删除。

4 试验结果与分析

4.1 普通蚁群算法+PLSR

使用MATLAB R2014a进行光谱预处理和PLS回归建模分析。以1 306个（3 50～1 655 nm）全波长变量作为选择对象，根据经验值与优选策略，设置普通蚁群算法的初始种群数量为100，单个蚂蚁寻找最大变量数为40，显著性因子为0.01，衰减系数为0.65。对每次迭代选出的最优波长变量组中波长变量的相应贡献加1，迭代完成后，根据波长变量的贡献率，按照前向选择规则选出最佳特征波长变量组合。将普通蚁群算法（ACO）的迭代次数分别设置为30、50、100次，对产生的特征波长变量数及其建立的PLSR模型与全波PLSR模型的评价参数进行比较，如表2所示。

从表2可以看出，使用普通蚁群选出的特征波长可以提高模型精度。随着迭代次数的增加，模型的校正集均方根误差（RMSECV）逐渐降低，这是因为算法中将目标函数设置为模型的均方根误差。但验证集均方根误差（RMSEP）有先变小后增大的趋势，综合各项指标以及算法时间复杂度，选择迭代次数为50次。运行MATLAB R2014a软件，得到基本蚁群算法运行后波长变量累计贡献率图（图3）和特征波长点筛选图（图4）。图3中纵向直线代表每个波长被选择的次数即贡献率，波长被选择次数越多说明该波长贡献越大，代表该波长被选为特征波长的概率越大。以被选中的波长贡献率为依据，通过前向选择选出最终的特征波长;图4中圆圈代表对应波长被选择为特征波长，可以清楚地看出最终选择出的特征波长分布情况，同时可以看出，该方法选出的特征波长点数较多。

4.2 精英策略蚁群算法+PLSR

在变量有效性精英策略蚁群算法（EEAS）的波长变量选择方法研究中，参数设置与普通蚁群算法相同，初始种群大小为100，最大变量数为40，显著性因子为0.01，衰减系数为

0.65，迭代次数分别为30、50、100次。引入前向选择规则，根据波长点贡献率大小进行再次筛选，最终得到最佳特征波长，图5为通过MATLAB R2014a画出的基于变量有效性精英蚁群算法的所有波长点累计贡献率图，图中纵向直线表示精英蚁群算法运行过程中波长被选择的次数，即贡献，可以清楚地看出哪些波长点被选中的频次高，即该波长点被选为特征波长的概率大。图6为通过MATLAB R2014a画出的基于变量有效性精英蚁群算法的特征波长点筛选图，图中圆圈代表最终选择出的特征波长，可以清楚了解最终选择出的特征波长分布，与图4相比，本算法所选出的特征波长明显减少。其对应偏最小二乘回归（PLSR）模型评价参数如表3所示。与普通蚁群算法相比，迭代次数为30、50、100次时，模型的校正集均方根误差均有降低，且模型的决定系数提高，说明基于变量有效性的精英蚁群算法可以更准确选出与建模相关的特征变量。另外，从所选择的波长点图（图6）可以看出，基于变量有效性的精英蚁群算法选择出的特征波长明显减少，波长点分布在 400～480 nm、550～750 nm之間以及950、1 400、1 600 nm 附近，与文献[15-16]对土壤速效磷含量预测有重要作用的波长大致相同，证明了本算法的正确性。

由图7可以看出，基于变量有效性的精英蚁群算法收敛速度明显加快，而且建立的回归模型的RMSECV更小。在相同的参数下通过EEAS特征波长算法建立的PLSR模型效果更好，选出的波长点少，可降低模型的复杂度。

从图8可以发现，EEAS-PLS模型的点集中在y=x附近，其次是ACO-PLS模型，而全谱模型点比较分散。说明采用特征波长建立的速效磷含量预测模型优于全谱PLS模型，EEAS选出特征波长建立模型效果优于ACO选出特征波长建立的模型。所以使用基于变量有效性的精英蚁群算法选择特征波长，可以有效去除近红外光谱中的光谱信息重叠、冗余，降低模型复杂度、节省时间，同时提高模型预测精度。

5 结论

采用变量有效性精英蚁群算法选择的特征波长变量数少于普通蚁群算法。试验结果表明，利用特征波长变量建立的回归模型预测效果优于全谱波长变量。基于变量有效性精英蚁群算法的PLS模型可以实现土壤速效磷含量的较高精度预测。

参考文献：

[1]杨爱霞，丁建丽，李艳红，等. 基于可见-近红外光谱变量选择的荒漠土壤全磷含量估测研究[J]. 光谱学与光谱分析，2016，36（3）：691-696.

[2]马慧敏，陈 亮，潘康兵，等. 防磷固定剂对土壤有效性的影响[J]. 西南农业学报，2017，30（11）：2533-2536.

[3]王朴杰，王世东，张合兵，等. 基于高光谱的复垦农田土壤有机质含量估测[J]. 土壤，2018，50（3）：558-565.

[4]林志丹，汪玉冰，王儒敬，等. 波长优选对土壤有机质含量可见光/近红外光谱模型的优化[J]. 发光学报，2016，37（11）：1428-1435.

[5]宋相中，唐 果，张录达，等. 近红外光谱分析中的变量选择算法研究进展[J]. 光谱学与光谱分析，2017，37（4）：1048-1052.

[6]宾 俊，范 伟，周冀衡，等. 智能优化算法应用于近红外光谱波长选择的比较研究[J]. 光谱学与光谱分析，2017，37（1）：95-102.

[7]郭志明，黄文倩，彭彦昆，等. 自适应蚁群优化算法的近红外光谱特征波长选择方法[J]. 分析化学，2014，303（40）：513-518.

[8]Shamsipur M，Zare-Shahabadi V，Hemmateenejad B，et al. Ant colony optimisation：a powerful tool for wavelength selection[J]. Journal of Chemometrics，2006，20：146-157.

[9]Bray R H，Kurtz L T. Determination of total，organic，and available forms of phosphorus in soils[J]. Soil Science，1945，59（1）：39-45.

[10]Viscarra Rossel R A，Bui E N，de Caritat P，et al. Mapping iron oxides and the color of Australian soil using visible-near-infrared reflectance spectra[J]. Journal of Geophysical Research：Earth Surface，2010，115：F4.

[11]Rossel R A V，Behrens T. Using data mining to model and interpret soil diffuse reflectance spectra[J]. Geoderma，2010，158（1/2）：46-54.

[12]江明珠，韓邦兴，颜 晖，等. 便携式近红外光谱仪快速无损鉴别霍山石斛枫斗和河南石斛枫斗[J]. 林产化学与工业，2017，37（5）：102-106.

[13]吴碧霞. 一种基于精英策略的蚁群优化QoS单播路由算法[J]. 莆田学院学报，2011，18（5）：73-76.

[14]邸 亿，龙 飞，李卓越，等. 一种基于改进蚁群算法的多目标跟踪数据关联方法[J]. 计算机应用与软件，2013，30（4）：306-309.

[15]Mouazen A M，Maleki M R，De Baerdemaeker J，et al.On-line measurement of some selected soil properties using a VIS-NIR sensor[J]. Soil and Tillage Research，2007，93（1）：13-27.

[16]Daniel K W，Tripathi N K，Honda K. Artificial neural network analysis of laboratory and in situ spectra for the estimation of macronutrients in soils of Lop Buri （Thailand）[J]. Australian Journal of Soil Research，2003，41（1）：47-59.