群体小麦条锈病发病动态无人机遥感监测方法

苏宝峰，刘昱麟，黄彦川，蔚 睿，曹晓峰，韩德俊

（1. 西北农林科技大学机械与电子工程学院，杨凌 712100； 2. 农业农村部农业物联网重点实验室，杨凌 712100；3. 陕西省农业信息感知与智能服务重点实验室；4. 西北农林科技大学农学院，杨凌 712100；5. 西北农林科技大学旱区作物逆境生物学国家重点实验室，杨凌 712100）

0 引 言

无人机成像多光谱遥感在作物病害检测方面已经取得了一定的成果。Liu 等[13]通过无人机搭载传感器采集感染赤霉病的小麦灌浆期高光谱图像，利用误差反向传播（Back Propagation，BP）算法建立了赤霉病监测模型，总体准确率达98%。Rodríguez 等[14]对比了5 种机器学习算法在马铃薯晚疫病检测和评估中的性能表现，证明了可以通过无人机多光谱技术采集马铃薯冠层反射率，并建立线性支持向量机（Liner Support Vector Machine，LSVM）和随机森林分类模型来检测和评估马铃薯晚疫病发病状况。Guo 等[15]通过无人机搭载光谱成像传感器采集小麦条锈病冠层反射率，提取其植被指数和纹理特征，建立了不同感染期的偏最小二乘回归的小麦条锈病监测模型，病害指数与所提取特征的决定系数可达0.82。Su等[16]利用无人机采集了接种条锈病的小麦生长季的冠层多光谱数据，利用不同的冠层植被指数在不同发病期识别条锈病发病区域。Bhandari 等[17]通过无人机获取小麦冠层RGB 图像计算植被指数，证明了其与条锈病侵染系数之间具有较好的相关性。然而上述针对作物病害识别监测的研究中，只对少数品种进行分析，所得出的结论在多基因型育种群体中的应用还有待进一步研究。Chivasa 等[18]基于无人机多光谱技术检测玉米条纹病，建立不同品种抗病等级的随机森林分类模型，表明无人机成像光谱能够提高作物表型鉴定效率。上述研究表明无人机低空遥感技术不但能够提高对作物病害监测效率，还能为育种群体的表型分析提供支持[19-21]，然而仅通过单一时间点的冠层数据无法实现对作物发病情况的动态分析[22]。

因此，本文以田间感染条锈病的育种小麦自然群体为研究对象，获取其无人机成像多光谱数据与条锈病严重度人工调查数据，提取冠层光谱植被指数，并利用特征选择和分类算法建立群体小麦田间条锈病发病阶段和严重度的分类模型，并筛选出对发病阶段和病害严重度分类敏感的特征，以探究通过多时相的光谱植被指数进行群体小麦田间条锈病动态分析的可行性；并通过指数响应时间序列，量化并分析田间条锈病发病状况，以期为群体小麦条锈病的发病动态分析提供一种客观、高效的高通量方法。

1 材料和方法

1.1 试验区概况及试验材料

试验区域位于陕西省咸阳市杨陵区曹新庄试验农场（108°5′34″E，34°8′17″N），海拔约480 m，属于温带季风性气候，年均气温约13 ℃，年均降水量635 mm。试验期间的日平均空气温湿度情况如图1 所示，2021 年4月19—26 日，由于降水，日最大空气湿度保持在较高的水平，且空气温度较为适宜，利于条锈病夏孢子繁殖与传播，使试验区域内大部分小麦品种受到感染，产生条锈病夏孢子堆，出现春季流行症状。

试验材料为中国小麦联合攻关和小麦产业技术体系各单位选送参加区试的510 份小麦新品种（系）与重复种植15 次的6 个参考品种（济麦22、百农207、周麦18、西农511、偃展4110、周麦36）构成的自然群体，按育种试验增广设计实行小区种植，共600 个小区，小区大小为1 m2。试验区域的正射影像如图2 所示。

1.2 数据采集

1.2.1 无人机成像多光谱数据获取

本研究中无人机成像多光谱数据采集系统采用深圳大疆创新科技有限公司的经纬Matrice 100 四旋翼无人机作为平台，搭载美国Micasense公司的RedEdge成像系统，如图3 所示。

设B层中与Aj相关的因素成对比较判断矩阵经过了一致性检验，求得单排序一致性指标为CI(j)，(j=1，…，m)，相应的，平均随机一致性指标为RI(j)，CI(j)、RI(j)已经在层次单排序时求得，那么B层总排序随机一致性比例为：

RedEdge 传感器具有5 个通道，其名称、中心波长/半峰全宽分别为：蓝（Blue，B）475 nm /20 nm、绿（Green，G）560 nm/20 nm、红（Red，R）668 nm/10nm、红边（Red Edge，RE）717 nm/10nm、近红外（Near Infrared，NIR）840 nm/40 nm，每个通道空间分辨率为1 080×756 像素；成像系统下行光照传感器（Downwelling Light Sensor，DLS）模块用于数据采集时环境光线与太阳角度的校正；平台集成的全球定位系统（Global Positioning system，GPS）用于记录图像的位置信息；反射率校正板用于校正各个通道的反射率。使用该平台在2020 年3 月31 日至5月14 日采集15 次小麦冠层成像多光谱数据，数据采集时间为能见度良好的10:00－15:00，无人机飞行高度为20～21 m，前向与侧向重叠率为80%～90%，飞行速度保持在1.4～1.6 m/s。

1.2.2 小麦条锈病人工鉴定数据采集

按照小麦叶片上条锈病夏孢子所占整个叶片面积的百分比，将严重度划分为13 个等级（0、1%、5%、10%、20%、30%、40%、50%、60%、70%、80%、90%、100%）[23]。在2021 年5 月12 日，按照上述标准对研究区域内600 个小区进行病害严重度鉴定数据采集，将其划分为3 个严重度等级，即病害严重度在0～30%为高抗反应型（Resistance，R），在40%～60%为中抗反应型（Moderate Resistance，MR），在70%～100%的为感病反应型（Susceptible，S）[24]。

依照发病后时间顺序将条锈病发病划分为5 个阶段，对应数据采集日期为：4 月22 日、4 月26 日、5 月1 日、5 月6 日、5 月12 日，分别记录为S0、S1、S2、S3 和S4，其中S0 为田间第一次出现明显夏孢子堆的状态。（Optimized Soil Adjusted Vegetation Index，OSAVI）图像分割小麦冠层与背景土壤，通过统计比较确定每次多光谱数据的阈值生成二值化掩膜图像。最后，通过图像裁剪处理得到试验群体600 个种植小区的反射率指数图像和掩膜，并计算和提取各小区光谱植被指数。本文中所用到的光谱植被指数的定义及计算公式如表1 所示（人工选择的阈值会影响OSAVI 指数的提取结果，故该指数，仅用作阈值分割，不参与后续建模与分析）。

表1 本研究中使用植被指数公式Table 1 Spectral vegetation indices used in this study

1.3 研究方法

首先对无人机多光谱图像进行预处理。由于光谱植被指数拥有明确的物理意义，能够探测病害引起的生理生化过程[25]，故构建由光谱植被指数、条锈病发病阶段和条锈病严重度等级数据组成的数据集；通过特征筛选和机器学习算法建立条锈病发病阶段和严重度的分类模型，同时筛选出对分类结果敏感的特征；并基于筛选出的光谱植被指数响应的时间序列，分析部分试验材料在不同阶段的田间条锈病发病动态。本研究整体技术路线如图4 所示。

1.3.1 无人机图像预处理与植被指数提取

首先，使用Pix4Dmapper（Pix4DInc., Switzerland,http://www.pix4d.com/）软件对无人机采集的多光谱图像进行拼接和辐射校正，最后生成5 个通道的单波段反射率指数地图；为了后续批量处理多次遥感数据，利用QGIS 软件（https://www.qgis.org/）结合地面控制点对拼接完成的单波段反射地图进行图像几何校正；而后，采用阈值法进行图像分割，以优化土壤调节植被指数

1.3.2 光谱响应分析

在判断条锈病严重度时，由于群体小麦不同品种（系）间的初始状态、对病害的敏感程度和染病速率不同，单一时间点的光谱特征无法评估不同品种（系）小麦对条锈病胁迫的响应，为了消除单一时间点的误差，根据公式（1）计算小麦冠层22 个光谱植被指数的平均变化率作为条锈病胁迫下小麦的光谱响应[37]，并根据响应的时间序列作为构建严重度等级分类模型的特征和条锈病的发病动态分析指标。

式中response 为某一植被指数的响应，S0为该植被指数的初始状态，本文中取4 月22 日的各个光谱植被指数作为初始状态；Si为S0后第i个时期的光谱植被指数，本文中取i取1、2、3、4，对应1.2.2 中介绍的S1 ～ S4。

1.3.3 特征筛选与分类算法

随机森林是一种由多决策树集成的机器学习算法，其通过重采样方式多次随机抽取原始训练集一部分样本，进行多个决策树的建模，并通过投票方法输出最后的分类结果。随机森林算法在特征筛选上的随机性使其不容易出现过拟合，抗噪能力较好，且能对特征的重要性进行排序[38]。在对条锈病发病阶段的分类模型构建和敏感特征筛选中，利用试验群体600 个小区5 个发病阶段的22 个光谱植被指数（表2 中除OSAVI）构建建模数据集。数据集共3 000 个样本，每个样本具有22 个特征，分为S0、S1、S2、S3、S4 共5 个类别，并按照3:1 的比例划分训练集与测试集，设置随机森林算法的决策树个数为300，建立发病阶段的随机森林分类模型，同时给出22 个特征的重要性排序；在对条锈病病害严重度等级的分类模型构建和敏感特征筛选中，以5 月12 日采集的无人机多光谱数据所提取的22 个植被指数的响应作为特征，构建条锈病严重度等级随机森林分类模型（参数设置和数据集划分与发病阶段分类模型相同），并给出22个特征的重要性排序。

为了对比不同算法在筛选条锈病发病阶段和病害严重度敏感特征的结果，使用随机蛙跳算法对特征进行筛选。随机蛙跳是一种能够利用少量的变量迭代进行建模的高维数据变量选择方法，它能够输出每个变量选择的可能性，从而进行变量的筛选[39]。支持向量机（Support Vector Machine，SVM）算法在解决小样本、非线性及高维模式识别中具有一定的优势[40]，它首先寻找一个最大边际超平面，并将低维数据通过核函数映射到高维空间中，从而使线性不可分的样本变为线性可分的样本，并通过引入模型惩罚因子，提高分类模型的泛化性，达到更好地分类效果。通过随机蛙跳算法（设置族群中青蛙的数量N为10 000，子族群中青蛙的数量Q设为2）实现对条锈病发病阶段和病害严重度敏感特征的筛选；为了验证筛选特征对于分类问题的有效性，将筛选后的被选择概率排在约前1/3 的（前7 个）植被指数作为SVM的输入特征，构建条锈病发病阶段和病害严重度等级分类模型。

2 结果与讨论

2.1 群体发育阶段的动态分析

根据人工调查群体中600 个小区的生育期结果，该群体4 月16 日—22 日基本完成抽穗，4 月28 日—5 月4日期间群体小麦扬花，扬花授粉后灌浆缓慢启动，大约两周内籽粒含水率迅速增加，而干物质的积累很少，群体尚未明显衰老，期间条锈病是影响小麦长势的主导因素。这与利用群体冠层NDVI 指数时间序列分析[41]的结果一致：如图5 所示，NDVI 在拔节期后的饱和效应提前结束。

2.2 条锈病发病阶段分类

使用22 个光谱植被指数作为分类特征构建的群体小麦条锈病发病阶段分类模型在测试集中的混淆矩阵如图6a 所示，筛选出的敏感特征如图7a 所示。通过随机蛙跳对条锈病发病阶段敏感特征筛选的结果如图7b 所示。将被选择概率排名约前1/3 的特征作为输入特征构建的SVM 分类模型在测试集中的混淆矩阵如图6b 所示。

由图6 可知，通过随机蛙跳筛选特征建立的SVM 分类模型和利用22 个植被指数构建随机森林分类模型对条锈病发病阶段的分类都取得了较好的效果，表明可以通过光谱植被指数来描述群体小麦田间条锈病发病阶段，模型的评价结果如表2 所示。SVM 算法的F1 分数为0.985，随机森林算法的F1 分数为0.970，略低于前者，表明随机蛙跳算法和随机森林算法都可以用来筛选对于条锈病发病阶段更为敏感的指数。S2 和S3 时期的精度和召回率相对较低，从混淆矩阵中也可以看出两种分类模型在S2 和S3 时期出现相对较为明显的混淆，可能是由于在两个阶段前后病害表现变化较小，导致光谱植被指数差异较小。

表2 发病阶段分类模型评价Table 2 Accuracy evaluation of classification models of disease stages

DVIRE、NDVIrededge、GVI 指数在随机森林和随机蛙跳算法中具有相似的排名，对条锈病发病阶段较为敏感，表明这些光谱植被指数能够更好的描述群体小麦田间条锈病发病动态。

2.3 条锈病严重度等级分类

使用随机森林算法建立的条锈病病害严重度等级分类模型在测试集中的混淆矩阵如图8a 所示，筛选出的病害严重度敏感特征如图9a 所示。通过随机蛙跳对条锈病病害严重度敏感特征筛选的结果如图9b 所示。将被选择概率排名约前1/3 的特征作为输入特征构建SVM 分类模型在测试集中的混淆矩阵如图8b 所示。

如图8 所示，两种分类算法分类精度较之对条锈病发病阶段的分类精度都有所下降。测试集中错误的样本主要集中在严重度相邻的两级之间（“R”与“MR”、“MR”与“S”）。由表3 可以看出，通过随机蛙跳筛选特征建立的SVM 分类模型的F1 分数为0.780，而利用22 个原始特征构建的随机森林模型的F1 分数为0.741。由于数据集中类别“R”的数量较多，因此在两种分类模型中类别“R”的精度和召回率最高；两种分类模型中“MR”的精度和召回率较低，易出现与严重度相邻级别的混淆。由于人工鉴定具有一定的主观性，给出的分类标签可能会与实际条锈病严重度情况有所差异，造成分类结果并不十分准确，但总体上看，人工鉴定数据与机器学习算法构建的分类模型对条锈病严重度等级判断结果基本一致。由此可知，能够通过由无人机多光谱数据提取的冠层光谱植被指数响应量化群体中不同小麦品种（系）条锈病病害严重程度差异。

表3 病害严重度分类模型精度评价Table 3 Accuracy evaluation of classification models of disease severity

如图9 所示，通过随机蛙跳算法筛选出较为重要的光谱植被指数与随机森林算法的分类贡献性结果相似程度不高，其中WI、NPCI、NDVI 在两种敏感特征筛选算法中都具有较高的重要性。

图10 为人工鉴定的条锈病严重度等级与分类效果较好的SVM 分类模型结果的可视化。可以更加直观地看出构建的SVM分类模型基本可以完成对田间小麦条锈病病害严重度的分类，个别混淆主要发生在严重度等级相邻的样本中（“R”与“MR”、“MR”与“S”）。

2.4 条锈病胁迫下光谱响应动态分析

根据两种分类算法对条锈病发病阶段和严重度的分类结果，筛选出6 个对条锈病发病阶段和病害严重度都相对敏感的指数，分别为 DVIRE、GVI、NDVI、NDVIrededge、NPCI 和WI，计算6 个指数相对于4 月22 日的响应，并在群体重复种植的15 组参考品种中随机选择一组，分析其冠层植被指数响应时间序列，如图11 所示。

在6 个参考品种中，济麦22 和百农207 病害严重度等级为“S”，偃展4110、周麦18、西农511 的严重度等级为“MR”，周麦36 的严重度等级为“R”。在5 月12 日，6 个指数基本能够区分3 类病害严重度等级，其中DVIRE、NDVI、NDVIrededge、NPCI 和WI 响应的时间序列中，效果更为明显；然而在上述5 个指数中NPCI与WI 对条锈病发病阶段的分类并不敏感；对发病阶段分类贡献率较大的DVIRE、GVI、NDVI 和NDVIrededge中，GVI 对病害严重度等级也不够敏感；而DVIRE、NDVI和NDVIrededge 同时对发病阶段和病害严重度的分类有较大的贡献。在条锈病胁迫早期，小麦冠层DVIRE 的响应最快，其他植被指数响应较慢，GVI、NPCI 和WI、在条锈病胁迫早期出现了与最终结果相反方向的变化。

如图11a 所示，从DVIRE 的响应时间序列来看，在小麦条锈病胁迫早期，发病速率较快的是偃展4110，而不是最终病害严重度等级最高的百农207；偃展4110 在4 月29 日之后发病速率存在明显的降低，并最终处于“MR”级别，而百农207 的DVIRE 的响应基本呈线性趋势下降；而最终严重度等级最低的周麦18 在4 月26 日之后就已经发生类似于偃展4110 的变化，发病速率明显减缓，并在5 月12 日人工鉴定时呈现出“R”等级。

3 结 论

本文提出了一种群体小麦条锈病发病动态无人机遥感分析方法，并通过所建立的机器学习模型验证了其可靠性，主要结论如下：

1）利用22 个原始特征构建的随机森林分类模型与通过随机蛙跳筛选特征建立的SVM分类模型对群体小麦田间条锈病发病阶段的分类都极为准确，证明了利用无人机遥感采集并提取的光谱植被指数能够描述群体小麦田间条锈病发病随时间的动态变化；筛选出的NDVIrededge、DVIRE 和GVI 对模型分类结果更加敏感，可以更好地实现对群体小麦田间条锈病发病阶段的分类。

2）利用22 个植被指数的响应构建的随机森林分类模型与通过随机蛙跳筛选特征建立的SVM分类模型都能够较好的完成对条锈病严重度等级的分类，证明了无人机采集并提取的光谱植被指数可以描述小麦群体中不同品种（系）之间的田间条锈病严重度差异；筛选出的WI、NPCI、NDVI 在两种特征选择算法中取得相似的重要性，对条锈病病害严重度更为敏感。

3）通过多时相的光谱指数响应的时间序列能够在病害发生的不同阶段，精细、客观地量化并分析育种群体中不同品种（系）的发病状态，实现对群体小麦条锈病表型的动态分析。在全部提取的植被指数中，DVIRE、NDVI 和NDVIrededge 可以同时描述群体小麦条锈病发病阶段和病害严重度，其中小麦冠层DVIRE 可以作为最佳指标精细量化群体小麦条锈病发病动态情况。