许敬诚,吕新,林皎,张泽,姚秋双,范向龙,洪延宏
(石河子大学农学院/ 新疆兵团绿洲生态农业重点实验室,新疆石河子832003)
随着转苏云金芽孢杆菌(Bacillus thuringiensis,Bt)基因抗虫棉的种植,棉铃虫种群增长被有效控制, 但棉蚜成为新疆棉田主要害虫之一[1-2]。棉蚜多群居于棉花叶背,受到棉蚜危害的棉叶片失水变形,并出现褶皱,严重时向下卷曲[3];叶背部棉蚜代谢的糖分滴落在叶表面阻碍植株的光合作用和呼吸作用,影响植株的生长发育。 传统棉蚜信息监测需要人工观察受害棉花叶片变化,统计棉蚜数量,费时费力,并且由于统计的延时性造成虫情信息滞后, 不利于棉蚜危害精准防控。 农田传感器的普及使棉蚜信息快速获取成为可能。 机器视觉技术因其识别速度快而被广泛应用于虫情监测,Ghyar 等利用灰度共生矩阵和颜色矩阵分离叶片病变区域,实现水稻的病虫害的机器识别[4]。基于 LAB 颜色空间,Madhuri 等通过提取目标对象的纹理特征构建支持向量机模型完成田间害虫的分类[5]。 高光谱成像技术同时集成了机器视觉技术和高光谱技术特点,能同时反映被测样品外观形态与内部结构等特性,作为一种快速、无损的识别方法在植物虫害的症状检测方面具有很大优势。 李震等通过基于叶绿素敏感波长建立叶绿素含量预测模型检测柑橘红蜘蛛为害[6]。Wu 等基于偏最小二乘判别分析法和反向传播神经网络对菜心叶片高光谱成像影像上幼虫部分进行识别,在近红外波段可以有效提取出幼虫所在区域[7]。 Cao 等通过主成分分析法选择505、659 和955 nm 特征波长,以反向传播神经网络模型实现粮仓害虫快速的识别[8]。 Huang 等采用连续投影法和竞争性自适应重加权抽样法选择特征波长,建立偏最小二乘判别分析和最小二乘支持向量机模型完成桑螟幼虫为害等级的快速诊断[9]。 所以,基于高光谱成像技术,可以多角度对虫害进行检测,指导农作物的植保工作。
现阶段纹理特征获取多基于RGB 图像[10-12]。高光谱成像技术因其光谱分辨率高、图像信息更为丰富,应用于提取目标特征是虫害监测的新方向之一[13-14]。 本研究以受棉蚜为害的棉花叶片为研究对象,采集健康与受棉蚜为害的棉花叶片正面高光谱图像,通过不同光谱信息降维法提取特征波长,采用灰度共生矩阵方法提取特征波长下灰度图像的纹理特征,以纹理特征向量为输入建立受棉蚜为害棉花叶片的判别模型,以期为棉花蚜害快速监测与精准农业靶向施药提供技术支持。
本试验于2019 年7 月在新疆维吾尔自治区石河子市石河子大学农学院试验场进行。 试验田种植的棉花品种为新陆早45 号, 依照新疆传统高密度膜下滴灌栽培模式种植,采用常规水肥管理模式, 不同处理的棉花植株分别于2 顶孔径0.075 mm(200 目)封闭防虫网帐篷中生长。
在棉花盛花期采集叶片,上午以活体形式采集,用果枝剪沿叶柄基部剪下叶片,除去寄生在叶背的棉蚜, 并轻拭棉花叶片正面擦除杂质,随后装入单独的自封袋中编号,迅速放入装有生物冰袋的冰盒中冷藏保存,带回室内获取高光谱影像数据。 本研究共采集142 片棉花上部叶片,其中健康叶片71 片,蚜虫为害叶片71 片。
高光谱影像数据采集采用SOC710VP 可见光- 近红外地物高光谱成像仪 (Surface Optics Corporation,美国),光谱波段范围351~1 044 nm,光谱分辨率可设置。 本研究中,光谱分辨率设定为1.3 nm,光谱波段数128。采集暗室内置带刻度升降平台控制样本与相机之间的距离, 光源为2盏75 W 卤素灯(OSRAM,德国),光线直射暗室内壁漫反射于待测物体上。 在高光谱成像信息采集前通过对相机曝光度、焦距、光圈、扫描速率的调整,确保采集到的图像纹理清晰、不失真、无形变,仪器和样品不会被卤素灯高温损伤。 通过尝试确定高光谱成像仪的参数分别设置如下:样品距离镜头64 cm,光圈1.4,曝光时间19 ms,扫描速率 150~200 帧·s-1, 图像分辨率为 696×520像素, 光源焦点在镜头与拍摄叶片二分之一处,如图1 所示。
将新鲜棉花叶片依次置于底面垫有低反射率黑色背景板的高光谱成像暗室中,同时在距叶片边缘2 cm 处放置标准灰板辅助黑白校正。SOC710VP 可见光-近红外地物高光谱成像仪获取的光谱图像原始数据是像元亮度值(Digital number,DN), 通过自带软件 SRAnal710 进行光谱标定、空间辐射与光谱辐射标定可以将像元亮度值转化成带光谱反射率。 全波段光谱首尾存在噪声,且光谱影像有椒盐噪声。 去除噪声明显的光谱波段,对400~840 nm 光谱采用多项式平滑法(Savizky-Golay smoothing,SG)去噪。避开棉花叶片主叶脉所在区域, 于单片棉花叶片上提取2处面积为50×50 像素纹理清晰的感兴趣区域(Region of interest,ROI)作为样本,共得 284 个样本。 随机划分67%的样本作为建模集,其余部分作为预测集。
图1 高光谱图像采集平台Fig.1 Hyperspectral image acquisition platform
全波段光谱数据的模型存在大量冗余数据和同质性数据影响判别精度,因此为了压缩数据量,降低模型复杂度,本研究采用主成分载荷法(Principal component analysis-Loading, PCA-Load ing)[15]、随机蛙跳算法(Random frog, RF)[16]和连续投影法(Successive projections algorithm,SPA)[17]选择与受棉蚜为害后棉花叶片变化相关性更强的特征波长。
主成分载荷法根据主成分分析的结果,计算出每个波长的载荷,反映主成分对每个变量的相关性。 较大的正或负载荷意味着该变量对于样本主成分具有重要的影响, 对应的变量相对较为重要。
随机蛙跳算法是衡量变量重要性的一种方法,通过建立变量的正态分布空间并预选初始变量集,将初始变量子集中的变量通过不断选择选入候选子集,候选子集中的变量随着迭代而不断变化,统计每次迭代各变量出现的频率,子集被选频率越高说明其越重要,依据样本频率排名选择变量序号对应的波段为特征波长。
连续投影法主要解决共线性问题,能从信息变量中选择冗余信息最少的变量,提取共线性最小的变量,将各个波长组合分别建立子集,每个子集逐一采用多元线性回归法计算均方根误差(Root mean square error, RMSE)值,选择平稳且最小的RMSE 值对应子集作为特征波长子集。
图像视觉特征提取方法主要有颜色、形状和纹理提取3 类。 健康棉花叶片光滑平展,受棉蚜为害的棉花叶片皱缩卷曲,叶表面纹理特征变化明显。 本研究主要采用灰度共生矩阵(GLCM)[18]算法进行纹理特征提取, 在提取纹理特征时,以0°、45°、90°、135°作为常用方位角,通过计算 2 个像素之间的相近关系的概率来获得能体现物体特征的二阶统计量。 灰度共生矩阵计算所获得的无量纲参数数量较大, 不宜直接作为纹理特征,通常基于其构建统计量,结合方位角建立纹理特征向量。 常用的4 个不相关的纹理特征有能量(Energy)、熵(Entropy)、对比度(Contrast)、相关性(Correlation)[19]。 能量反映图像灰度粗细分布的均匀程度,图像的灰度分布越均匀,能量值越大;熵反映图像灰度值散乱度,灰度值分布随机度越高熵值越大; 对比度反映图像强度差异,图像灰度值差异大时,对比度大;相关性反映图像纹理在行或列上一致性,图像横向或纵向纹理均一时,相关性值增大。 本研究中,基于特征波长选择方法选择的特征波长下的感兴趣区域灰度图像进行纹理特征提取。
本试验采用2 种分类方法建立模型:线性判别分析算法偏最小二乘线性判别分析(Partial least-square-lineardiscriminantfunction,PLS-LDA)模型[20]和非线性判别分析算法支持向量机(Support vector machine,SVM)模型[21]。
PLS-LDA 模型是一种线性分类方法,它基于PLS 回归模型预测每个样本的类别数,选取平方和预测误差最小值处的最小潜在变量数(Latent variables,LVs)。本研究采用去中心化法对数据进行归一化处理,通过10 000 次蒙特卡洛验证获取建模集的最佳潜在变量值LVs,通过多次寻优选择最佳阈值来为分类后的样本归类。
构建支持向量机(SVM)模型需要对2 个参数寻优:核函数参数和模型的惩罚系数。 由于选择合适的核函数没有系统的方法,而本试验的样本呈现非线性分布,径向基(Radial basis function,RBF)核函数能够更好处理非线性关系并降低训练过程计算的复杂程度。 本试验选择径向基核函数作为SVM 的核函数构建RBF-SVM 模型,随后通过蒙特卡洛交叉验证对惩罚参数c 和样本影响半径的倒数g 寻优,c、g 参数的变化范围设为-5~5,基于多次交叉验证后的结果选择最优的c 和 g 的组合。
试验获取的可见光-近红外高光谱波长范围为 351~1 044 nm, 共 128 个波段。 由于试验环境、仪器的影响以及暗电流干扰,造成获取光谱前后端有明显噪声, 因此剔除噪声明显的波段,截取波长400~840 nm 范围共86 个波段光谱影像用于后续建模分析,健康棉花叶片和棉蚜为害棉花叶片的平均光谱如图2 所示。
棉花叶片高光谱特征在415 nm 有一个明显的吸收谷,绿光波段开始出现反射峰,在557 nm达到最高值; 红光波段出现吸收谷,677 nm 为最低值。
图2 切除两端噪声后的平均光谱Fig.2 The mean spectra after removing both two ends noise
主成分分析法可以将原可能相关的变量正交变换为更多不相关的变量,显示数据的内部结构。 利用主成分分析法对401~842 nm 波段棉花叶片光谱数据进行分析, 将86 个波段的信号转化为若干个主成分(Principal component,PC),得到前3 个主成分第一主成分PC1(96.16%)、第二主成分 PC2(2.88%)、第三主成分 PC3(0.46%),累计贡献率达到99.5%, 能解释大部分变量。PC1、PC2 和 PC3 三维得分分布图如图3 所示,健康样本和棉蚜为害样本之间存在非常明显的聚类和区分,同时可以发现,2 种样本在三维空间中存在一定的重合,需要进一步对光谱数据进行分析。
图3 主成分聚类分布图Fig.3 Cluster plots based on the PCs
选择特征波长对于去除高维数据中的冗余信息,优化校准模型,获得良好的结果具有重要意义。 本试验利用主成分载荷法(PCA-Loading)、随机蛙跳算法(RF)和连续投影法(SPA)提取特征波长(Characteristic wavelengths, CWs),3 种算法分别筛选 9、10、13 个特征波长, 如图4 和表1表示。
图4 PCA–Loading (a), RF (b), SPA (c)法提取特征波长示意Fig.4 The optimal wavelengths selection by PCA-loading (a), RF(b) and SPA (c)
PCA-Loading 法选取主成分载荷图曲线绝对值大于0.05 的波峰和波谷作为特征波长点,共筛选出9 个特征波长,占总波长数的10.7%。 RF 法设置算法迭代次数N 为10 000 次, 提取数量10个,蛙跳初始种群数目Q 为2 个,以每个光谱被选择的可能性为筛选依据,运行结果为降序排列的被选择可能性, 设定被选择可能性阈值为0.829, 共筛选出10 个特征波长, 占波长总数的12.8%。 SPA 法以全光谱范围的波长依次对应每个变量,设置被提取的最小变量个数为5,最大变量个数为30,以寻找最小共线性变量组,最终共筛选出 13 个特征波长,RMSE 值为 0.001。3 种算法筛选出的特征波长主要集中在绿峰波段和红边波段。
健康的棉花叶片正面表面平展,由于棉蚜从棉花叶片背部吸食汁液,影响养分供给,导致受到棉蚜为害的棉花叶片表面与叶肉组织结构遭到破坏,出现皱缩,分支叶脉向下凹陷形成阴影,与向上突起的叶肉形成明显灰度差。 在纹理特征感兴趣区域选择上, 由于主叶脉灰度均一性较高,与叶肉有明显差异,且阴影面积大,不利于纹理特征的提取,因此对感兴趣区域的提取要避开叶片主叶脉位置。 现从特征波长对应灰度图像中选取50×50 像素无主叶脉区域作为感兴趣区域。 试验中将原始图片256 等级的灰度值压缩到16 级以减小计算量, 以 4 个方向 0°、45°、90°、135°,固定取样距离为1 像素构建4 个方向的灰度共生矩阵,分别统计各个矩阵能量、熵、对比度、相关性,如图5 和表2。
表1 对光谱样本提取的特征波长Table 1 The effective wavelengths selection for hyperspectral sample
图5 不同波长健康样本与为害样本对比Fig.5 The contrast of healthy and aphid cotton leaf’s ROI in different wavelengths
表2 不同样本纹理特征向量参数Table 2 The average feature vector of GLCM for two types of leaves
根据不同波长下灰度图像可知,健康叶片表面平滑, 分支叶脉的纹理在可见光波段不能辨认,在近红外波段能看到少量较粗的部分,单张图像中纹理分布均匀,灰度变化小。 为害叶片可见光波段和近红外波段均能看见清晰的分支叶脉轮廓并由边缘始形成向心突起,纹理呈斜向分布,向心突起处有亮斑,与分支叶脉处形成明显灰度差。
通过4 个角度的灰度共生矩阵计算得到4个二阶统计量,包括能量、熵、对比度、相关性,共得到16 组特征向量。 分别计算不同处理样本的特征向量,得到2 组样本基于这4 个角度的纹理特征平均值和标准差。 健康样本的能量、对比度均超过为害样本,熵低于为害样本,由于为害样本纹理分布散乱, 垂直方向和斜向灰度变化杂乱,导致相关性差异不显著。4 个角度的特征向量及均值、标准差有显著差异,可以解释叶片受棉蚜为害后的变化。 由于纹理特征标准差数据分布接近,因此采用特征波长图像的纹理特征标准差作为输入,建立了棉蚜识别模型。
分别基于3 种特征波长选择方法选择的特征波长下的灰度图像获取的纹理特征, 建立PLS-LDA 和 RBF–SVM 模型, 其判别分析结果如表3 所示。
表3 基于纹理特征的PLS-DA 和RBF-SVM 判别分析模型判别结果Table 3 The discriminant results of PLS-DA and RBF-SVM models using texture features
比较模型分析效果,使用全部纹理特征数据集进行分析,SPA 法提取的特征波长组识别率最高,2 种模型识别率在82.98%以上,RF-PLS-LDA模型效果最好,预测集识别率达到91.49%。 结果显示,使用全部纹理特征数据集建模识别效果较好,但输入数据量较大,模型计算时间效率较低。逐一以纹理特征二阶统计量作为输入建立判别模型,2 类模型预测集的识别精度均在74%以上,可以实现对棉花上棉蚜为害叶片的精确识别。 以能量作为输入的模型对棉蚜为害植物的识别效果最佳,预测集平均识别率达到89%;其次是对比度作为输入的模型, 预测集平均识别率达到84%。由于相关性纹理特征差异不显著,识别率最低。 其中,最优模型为SPA-PLS-LDA 模型,预测集识别率达到 92.55%, 其次是 PCA-Loading-PLS-LDA-Energy 模型(91.49%)和 SPA-RBFSVM-Energy 模型(90.43%)。
基于高光谱技术的作物监测模型主要以作物反射率作为模型输入数据[22-23],但作物在不同胁迫环境生长下形态、 颜色等特征会发生变化,出现“同物异谱”或“异物同谱”现象。 高光谱成像技术可以同时获取光谱信息和图像信息,避免此类现象发生。 因此本研究使用高光谱成像仪进行健康、 棉蚜为害叶片正面高光谱成像数据采集,通过不同方法提取特征图像。 基于灰度共生矩阵提取图像纹理特征, 依此建立蚜害判别模型,证明使用高光谱成像对受棉蚜为害棉花叶片进行分类预测的可行性,为基于多光谱成像的棉蚜监测装置提供技术支持。
棉花受到蚜虫为害后光谱特征上会出现变化[24]。本试验中受蚜虫为害的叶片蓝光、红光所在波段与近红外波段吸收率下降, 而植株光谱红、蓝光波段反射率降低是由于叶绿素含量降低所引起,近红外的反射率下降的原因是植株细胞结构变化所产生,这也从光谱角度证实了蚜虫为害会使叶片叶绿素和相对含水量下降,细胞结构出现变化[25-28]。 纹理特征上,蚜虫为害棉花叶片会导致其出现皱缩、卷曲,叶片局部形成凸起、向下卷曲等症状。 本研究利用灰度共生矩阵计算特征图像灰度值的变化,经由得出的二阶统计量可以直观地反映图像部分特征[29-30]。 纹理特征参数中健康样本的能量对比度超过蚜害样本,熵低于蚜害样本。 由于健康叶片表面平滑,图像灰度分布均匀,因此表面灰度值相近,反映图像像元灰度平方和的能量较高;而蚜虫为害叶片表面有明显由褶皱、蚜蜜造成的阴影和亮斑,因此图像的灰度值差异性较大,反映图像无序性的熵较大。
为精确蚜虫为害后光谱反射率变化的敏感区域,降低数据冗余度,使用主成分载荷法、随机蛙跳算法和连续投影法提取特征波长,从全波长信息中分别优选出 9、10、13 个特征波长。 以特征波长高光谱图像作为分析对象,利用灰度共生矩阵提取图像的纹理信息,并取其二阶统计量的标准差作为模型的输入,最终通过偏最小二乘判别分析法和径向基支持向量机建立分析判别模型。通过对比试验发现,RF-PLS-LDA 模型效果最好,预测集识别率达到91.49%。由于模型输入数据量较大,导致建模效率低,通过减少模型输入数据对模型进行优化。 逐一以纹理特征二阶统计量作为输入建立判别模型, 其中 PCA-Loading-PLS-LDA-Energy 模型性能最好,预测集识别率达到92.55%。该方法为棉花棉蚜虫情识别提供了新的思路,但其预测性能还须进一步改进。 下一步的研究重点为找寻不同蚜害等级下纹理特征的变化特点。
由于棉蚜为害导致棉花叶片叶绿素和相对含水量下降、细胞结构产生变化,因此棉蚜为害样本高光谱曲线蓝光、红光波段以及近红外波段的反射率较健康样本更低。 通过对比试验发现,RF-PLS-LDA 模型对棉蚜为害叶片的识别率最好,预测集识别率达到91.49%。 由于模型输入数据量较大,导致建模效率低,通过减少模型输入数据对模型进行优化。 逐一以纹理特征二阶统计量作为输入建立判别模型, 其中PCA-Loading-PLS-LDA-Energy 模型性能最好,预测集识别率达到92.55%。