基于模型集群的马铃薯叶绿素检测光谱变量筛选讨论

刘 宁， 邢子正， 乔 浪， 李民赞， 孙 红*， Qin Zhang

1. 中国农业大学现代精细农业系统集成研究教育部重点实验室， 北京 100083 2. Center for Precision & Automated Agricultural System, Washington State University, Pullman WA 99350, USA

引 言

叶绿素含量是评价马铃薯作物光合作用能力与营养水平的重要指标之一[1-2]。 在可见光-近红外区域， 分析含氢基团(O—H， N—H， C—H)振动合频和各级倍频的特性， 是开展作物叶绿素、 氮素、 水分等参数光谱学检测的理论基础， 取得了重要进展[3]。

作物叶绿素光谱学检测中， 常通过筛选特征波长来达到解析光谱变量、 剔除冗余信息、 压缩计算量、 提高诊断模型精度与鲁棒性等目的[4]。 因为相关分析筛选变量存在高度自相关导致的多重共线性问题， 在主成分分析的基础上， 连续投影算法(successive projection algorithm， SPA)、 无信息变量消除法(uninformative variables elimination， UVE)、 间隔最小二乘波长选择方法(interval partial least square， iPLS)、 变量投影重要程度系数法(variable importance in the projection， VIP)等算法被用于筛选特征波长并建立诊断模型[5-6]。

上述一次性建模筛选特征波长的方法， 数据处理易受样本个数的影响[7]。 针对此问题Li等提出基于模型集群思想的蒙特卡洛无信息变量消除(Monte Carlo uninformative variables elimination， MC-UVE)[8]、 随机蛙跳(random frog， RF)[9]、 竞争自适应重加权采样(competitive adaptive reweighted sampling， CARS)[10]等变量筛选算法。 有报道应用CARS算法设置迭代次数为50， 选取10个波长建立南瓜叶绿素检测模型， 精度为0.846。 郑涛等[11]采用MC-UVE算法迭代次数为500， 选出12个马铃薯叶绿素特征波长。 程萌等[12]基于RF算法筛选小麦叶绿素特征波长， 迭代次数为10 000， 选出8个最优波长。

此类研究中尚有如下问题需要深入讨论， 一方面应用不同算法选取变量是否存在差异， 建立的模型是否最优且稳健； 另一方面， MC-UVE， RF和CARS等算法中初始参数迭代次数普遍采用固定值， 修改迭代次数与其他约束是否对变量筛选结果有影响， 需要开展比较和分析。

因而， 在马铃薯作物叶绿素光谱学检测中， 分别应用MC-UVE， RF和CARS算法， 讨论迭代次数(number of iteration，N)参数和特征变量个数(latent variable， LV)对特征波长筛选结果的影响。 通过建立PLS模型， 阐明特征波长分布与叶绿素含量的解析能力， 以模型验证集精度为评价标准， 明确参数最优组合， 以期为马铃薯叶绿素光谱降维与高鲁棒性诊断建模奠定基础， 也为同类研究提供参考。

1 实验部分

1.1 材料

2018年在北京市昌平区小汤山国家精准农业示范基地开展实验， 马铃薯品种为“大西洋”。 30 m×40 m范围内设80个采样区， 在发棵期(M1)、 块茎形成期(M2)、 块茎膨大期(M3)和淀粉积累期(M4)4个生长期跟踪采集马铃薯冠层光谱并进行理化测试。

1.2 田间光谱数据采集与叶绿素含量测定

采用ASD FieldSpec HandHeld2 便携式地物光谱仪测定325～1 075 nm内751个波长处作物冠层光谱反射率， 采样间隔1 nm， 每点重复采集3次取平均值。 同步随机采集叶片经浸提后， 利用紫外分光光度计测定叶绿素含量， 测定方法参考相关文献。 每个生长期采集80组数据， 其中M1因植被覆盖度较低导致无效数据， 保留74组有效数据后， 全生长期共获取314组数据。 数据采集预处理总体流程如图1所示。 其中， 采用标准正态变量(standard normal variate， SNV)方法， 对原始光谱曲线进行预处理来消除环境噪声的干扰。 光谱与处理、 特征波长筛选以及PLSR建模均在matlab2014.a环境中完成。

图1 数据处理总体流程图

1.3 光谱特征波长筛选方法

基于模型集群分析的思想， 比较MC-UVE， RF和CARS 3种变量筛选算法， 在matlab2014.a libpls软件中实现。

(1)MC-UVE算法

MC-UVE算法基于偏最小二乘回归(partial least squares regression， PLSR)提出， 从训练集中取出一定数目(M个)样本构建PLS子集， 重复M次计算PLS回归系数矩阵， 引入变量稳定指数为筛选标准， 计算得到每个变量稳定指数值， 并从高到低排序筛选变量[6]。 其中， 保留的LV数量决定着模型的预测能力和模型的稳定性。

(2)RF算法

RF算法类似于可逆跳转马尔可夫链蒙特卡洛。 与PLSR相结合， 通过PLSR结果模拟一条服从稳态分布的马尔可夫链来计算每个变量被选择的概率， 从而进行重要变量的筛选[7]。

(3)CARS算法

CARS算法基于自适应重加权采样和指数衰减函数， 选取在PLSR模型中回归系数绝对值大的变量， 得到一系列波长变量子集； 然后对每个波长子集采用交叉验证建模， 从中挑选出模型均方根误差最小的子集[8]。 因此CARS算法筛选得到的特征变量个数一定。

为了检测作物叶绿素含量， 本研究以马铃薯作物为例， 对CARS算法的迭代次数(N)参数、 RF和MC-UVE算法的迭代次数(N)参数和特征变量数(LV)参数对叶绿素特征波长筛选结果的影响进行讨论。 迭代次数设置6个梯度， 分别为N=50， 100， 500， 1 000， 5 000和10 000； 特征变量数设置4个梯度， 分别为LV=15， 20， 25和30， 分析迭代次数(N)和特征变量数(LV)两个参数的最优组合情况。

1.4 PLSR模型建立与模型评价

2 结果与讨论

2.1 马铃薯作物生长期冠层反射光谱响应分析

SNV校正后的各生长期的马铃薯冠层反射光谱曲线如图2所示， 总体而言， 在可见光波段， 由于色素体对蓝、 红光的强吸收存在400～500与611～710 nm低反射率区， 并在400和680 nm附近出现吸收谷； 520～610 nm体现为色素体的强反射， 550 nm附近为绿色反射峰。 受到叶肉内海绵组织结构内的空腔反射率增强影响， 近红外711～760 nm快速攀升后进入761～1 000 nm高反射平台区， 其中970 nm附近出现水分的微弱吸收谷。 由M1至M4推进， 在400～500和740～880 nm反射率降低； 在530～640和910～960 nm反射率升高， 且M4和M1分别呈现与其他生长期较大的差别。 综上说明作物光谱响应是对植物生长过程中色素体、 水分分子、 结构等的综合表现， 针对叶绿素指标， 挖掘全谱中特征波长十分必要。

图2 SNV处理后生长期冠层平均反射光谱曲线

2.2 相关性分析与样本集划分结果

分析光谱反射率与叶绿素含量的相关性， 结果如图3所示。 在387～509， 519～633和744～844 nm波段， 二者相关系数绝对值(|r|)均高于0.6， 在678 nm达正相关峰值0.411； 在702 nm存在负相关峰值-0.715。 845～917 nm正相关系数逐渐降低， 917 nm之后呈负相关。 此结果与叶绿素吸收可见光蓝、 红光， 反射绿光的物理现象一致， 但相关性曲线显示相邻波长之间的相关系数接近。 若选取相关系数较高者为特征波长， 会存在波长冗余与多重共线性问题。 因此， 利用SPXY算法划分样本集结果如表1所示， 后续建模开展特征波长变量筛选方法讨论， 用建模集筛选特征波长、 建立回归模型， 以验证集的结果评价特征波长筛选结果。

图3 光谱反射率与叶绿素含量相关性曲线

表1 建模集与验证集划分统计

2.3 基于模型集群分析的马铃薯叶绿素特征波长筛选

2.3.1 MC-UVE算法

由于MC-UVE算法对于同一批光谱数据， 设置同样的迭代次数， 运行多次计算变量的稳定指数不一致， 因此分别讨论迭代次数(N)和特征波长数量(LV)的影响。

图4 MC-UVE算法在迭代次数为500时的运行结果

图5 MC-UVE在LV梯度下最佳迭代次数时特征波长位置

表2 基于MC-UVE的叶绿素含量检测PLSR模型验证集结果(RMSEV： mg·L-1)

2.3.2 RF算法

RF算法与MC-UVE算法类似， 首先讨论迭代次数N的影响， 分别设置N为50， 100， 500， 1 000， 5 000和10 000次6个梯度， 运行5次取平均值。 以N=10 000为例的运行结果如图6所示， 纵坐标为每个波长的被选择概率(selection probability)， 被选择概率越高说明波长越重要。 其次讨论波长个数LV的影响， 按照选择概率从大到小设置LV分别为15, 20, 25和30建立马铃薯叶绿素检测PLS模型， 共得到24种模型。

图6 RF算法在迭代次数为10 000时的运行结果

图7 RF在四种LV梯度下最佳迭代次数时特征波长位置

表3 基于RF在不同输入参数下的叶绿素含量检测PLSR 模型验证集结果(RMSEV： mg·L-1)

图8 CARS在迭代次数N为100时运行结果

由图7对比相关性分析结果可知， 在LV从15增至30过程中分布愈加广泛， 反映的信息愈加全面。 在LV=15时， 在绿光区域没有筛选到特征波长， 而在LV=20， 25和30时， 筛选到的特征波长在蓝、 绿、 红区域均有分布。 LV=30时， 970 nm附近反映水分弱吸收的波长被选中， 说明该方法筛选波长对含氢基团具有较好的选择性。

2.3.3 CARS算法

CARS算法与RF和MC-UVE不同， 对于同一批数据， 在相同的迭代次数(N)下变量筛选结果唯一， 所以仅考虑设置N为50， 100， 500， 1 000， 5 000和10 000次6个梯度。N=100时的运行结果如图8所示， 图8(a)为筛选过程中变量数随着迭代次数N的变化曲线， 筛选的波长数(LV)随运行次数的增加而减少； 图8(b)为RMSECV随着迭代次数的变化曲线， 在前30次时RMSECV保持不变， 30次后下降， 在迭代61次时RMSECV的值最小为3.928， 之后逐步攀升； 图8(c)为各光谱波长的回归系数的变化趋势， 其中“**”列表示RMSECV最小时所对应的迭代运行次数。 运行后得到的波长变量集采用交叉验证， 根据RMSECV的值来确定最优波长变量子集为21个特征波长。

图9 CARS算法不同迭代筛选到的特征波长位置

表4 基于CARS算法不同迭代次数的叶绿素含量检测PLSR模型验证集结果(RMSEV： mg·L-1)

Table 4 PLSR validation results on the chlorophyll content detection with iteration of CARS(RMSEV： mg·L-1)

迭代次数(N)最佳迭代次数特征波长个数(LV)R2vRMSEV5021670.6454.40810061210.6894.183500249390.6364.4601 000502380.6494.3795 0002960220.6804.29410 0005918220.6724.301

图10 MC-UVE， RF和CARS最优特征波长位置

2.4 三种波长筛选方法建模比较与讨论

对比分析MC-UVE， RF和CARS筛选的最优特征波长， 位置如图10所示， 从特征波长分布角度， 在可见光范围(400～710 nm)， RF算法筛选波长分布均匀； MC-UVE算法对550 nm附近绿光区域不敏感， 而在450 nm附近蓝光区域“波长聚集”现象显著； CARS算法对该区域筛选变量较少。 在近红外区域(711～1 100 nm)， RF算法得到的特征波长分布仍然较为均匀； MC-UVE在800～1 000 nm只筛选到一个特征波长； CARS筛选到的特征波长均聚集在900～1 000 nm内。 综上说明RF算法在可见光和近红外区筛选得到的特征波长对叶绿素光谱吸收和反射等特征具有较为全面的代表性。

从相关性的角度考虑， RF算法筛选得到的特征波长在叶绿素高相关范围(387～509， 519～633， 744～844和845～917 nm)和相关性峰值(702 nm)均有分布。 而MC-UVE算法筛选变量只在387～509和744～844 nm两个范围， CARS算法筛选变量则只有391， 392， 393， 394和896 nm五个波长落入高相关性范围内， 且前四个为相邻波长而存在波长信息冗余。 上述结果在PLSR模型中也得到了验证， RF-PLSR模型的精度最优， MC-UVE-PLSR模型次之， CARS-PLSR模型最差。

表5 MC-UVE-PLSR， RF-PLSR和CARS-PLSR验证集结果

综上表明， 当合理选择N和LV参数时， RF算法对马铃薯叶绿素特征波长筛选能力优于MC-UVE和CARS两种算法， 同时也避免了高相关性区间筛选相邻波长存在的高度自相关导致的多重共线性问题。 所建立的RF-PLSR模型可为马铃薯叶绿素含量诊断提供支持， 而研究讨论的变量筛选方法与参数分析过程， 可为其他同类光谱学检测提供参考。

3 结 论

为了高精度地检测马铃薯作物叶绿素含量， 利用基于模型集群思想的CARS， RF和MC-UVE三种算法筛选叶绿素特征波长， 建立叶绿素含量检测PLS模型。 以PLS模型验证集结果为评价指标， 讨论三种算法的迭代次数(N)和特征变量个数(LV)参数对模型结果的影响， 确定三种算法的最佳输入参数组合， 对比分析MC-UVE， RF和CARS筛选的最优特征波长， 结论如下：

对叶绿素含量和光谱数据做相关性分析， 发现在387～509， 519～633和744～844 nm三个波段内， 叶绿素含量与光谱反射率的相关系数较高， 其相关系数绝对值均高于0.6； 在678和702 nm处存在相关性极值， 相关系数分别为0.411和-0.715。