感染TMV与PVY两类病害的烟草样本由中国农业科学院烟草研究所提供。 供试烟草为分别在室内与室外培养的普通烟草品种K326。 室外种植区位于青岛即墨试验基地(N36° 26'52.54″, E120° 34'55.24″), 海拔高度为54.47 m, 四季变化和季风进退比较明显。 室内培养条件, 如表1所示。 待烟草植株长到7~8片叶期(生根期), 将TMV与PVY两种病毒分别接种于烟草植株上, 烟草接种病毒7 d之后开始测量烟草光谱数据。

表1

Table 1

表1(Table 1)

表1 室内培养条件表 Table 1 Indoor culture conditions 培养条件 参数 光周期 14L∶ 10D 温度/℃ 251± 1 相对湿度/% 65± 5

表1 室内培养条件表 Table 1 Indoor culture conditions

采集感病烟草光谱数据的仪器为美国ASD公司生产的Field Spec4便携式地物光谱仪, 光源采用光谱仪自带的卤素灯。 该套设备可测量的光谱范围是350~2 500 nm, 采样间隔为1 nm。 光谱采集过程中, 每株烟草视为一个样本。 为保证测量的是叶片表征范围, 测量时使光谱仪的探头对准烟草叶片中间, 保持仪器与待测叶片之间距离固定。 为降低测量时的偶然误差, 在测量一个样本时, 分别于烟株上、 中、 下叶位选择一片叶片, 每个叶片按照叶基、 叶中和叶尖部位(避开叶脉)进行光谱测量, 最终将所获光谱取平均得到一个样本数据。

数据采集过程中, 每隔2 d观察一次烟株病情进展, 根据《GBT-23222— 2008 烟草病虫害分级调查方法》, 对感染TMV与PVY的烟草按严重程度分成0, 1, 3, 5, 7和9六个等级, 其中0级表示健康烟草。 数据采集结束后共获11类, 累计1 697条合格的烟株样本, 各类样本数量如表2所示。 训练集主要选择室内数据, 测试集主要选择室外数据。 训练集比例为75%, 训练集与测试集配置如表3所示。

表2

Table 2

表2(Table 2)

表2 样本数量介绍 Table 2 Introduction to the number of samples Health TMV1 TMV3 TMV5 TMV7 TMV9 PVY1 PVY3 PVY5 PVY7 PVY9 室内 151 29 160 156 46 53 12 18 114 153 126 室外 180 126 66 37 29 28 84 85 9 20 15 编号 1 2 3 4 5 6 7 8 9 10 11

表2 样本数量介绍 Table 2 Introduction to the number of samples

表3

Table 3

表3(Table 3)

表3 训练集与测试集配置 Table 3 Configuration of training set and test set Health TMV1 TMV3 TMV5 TMV7 TMV9 PVY1 PVY3 PVY5 PVY7 PVY9 训练集 248 116 169 144 56 60 72 77 92 129 105 测试集 180 126 66 37 29 28 84 85 114 153 126

表3 训练集与测试集配置 Table 3 Configuration of training set and test set

为了减少随机误差, 更好地分析感病烟草光谱, 首先对每类光谱数据进行平均处理, 处理公式如式(1)所示

${{R}_{i,\text{mean}}}=\frac{1}{{{N}_{i}}}\overset{{{N}_{i}}}{\mathop{\underset{j=1}{\mathop \sum }\,}}\,{{R}_{i,\text{ }\!\!~\!\!\text{ }j}}$(1)

式(1)中, i(i=1, …, k)为病害编号, N_i为第i类样本的数量。 R_{i, j}为病害编号为i的第j个样本的高光谱数据。 经过平均处理, 共获得11条光谱曲线, 分别在图1(a— c)中展示。 为分析不同病害之间的差异, 同时对各类病害数据进行平均处理, 相应曲线如图1(d)所示。

图1

Fig.1

Figure Option ViewDownloadNew Window

图1 烟草光谱曲线图 (a): 不同种感病烟草光谱曲线图; (b): 感染TMV烟草的不同严重度等级光谱曲线图; (c): 感染PVY烟草的不同严重度等级光谱曲线图; (d): 健康及感染TMV、 PVY病害的烟草光谱曲线图Fig.1 Spectrum curve of tobacco (a): Spectral curves of different kinds of susceptible tobacco; (b): spectral curves of different severity levels of TMV infected tobacco; (c): spectral curves of different severity levels of PVY infected tobacco; (d): spectral curves of healthy and TMV and PVY infected tobacco

由图1(a)可知, 感病烟草光谱曲线整体上保持近红外波段高反射、 红波段强吸收的光谱特性; 在部分光谱范围内, 不同感病烟草光谱曲线相似性高, 不易区分。 由图1(b)可知, 不同严重度的TMV烟草在1 600~1 800 nm范围内差异较大。 由图1(c)可知, 不同严重度的PVY烟草在500~700 nm范围内差异较大。 综上, 感病烟草光谱数据存在较高的冗余, 仅在少数波段存在较强的区分性, 因此必须进行数据降维。

提出的FNG-NMFW算法流程如图2所示, 一共包含三个步骤: 采用FNGBS^[14]对烟草高光谱数据进行精细分组, 之后根据NMFW^[15]确定的权重系数选择各波段组特征波段, 最后采用SVM进行特征分类, 实现感染TMV与PVY烟草的精准检测。

图2

Fig.2

	Figure Option ViewDownloadNew Window
	图2 烟草病害检测方法流程图Fig.2 Flowchart of tobacco disease detection method

为了将FNGBS应用到烟草病害检测, 首先构建离散烟草高光谱数据X={x₁, x₂, …, x_i, …, x_L}^T∈ R^{L× N}, 其中x_i为各样本数据在第i个波段处的反射率, N为样本的数量, L为波段数。

首先设置波段分组数M, 即输入SVM的特征波段数。 高光谱分组共分为两步: 粗分组与精细分组。 粗分组的原则是使各波段组包含的波段数尽可能相等, 因此各组别之间的分割波段可由式(2)表示

G(m)=1, m=1(L-mod(L, M))×mM, 2≤m≤M-1L, m=M(2)

式(2)中, mod(· , · )为取余操作, 按照算法G(m)的原则, 最终粗分组结果{X'_m }m=1m=M可以由式(3)表示

X'm={xi}i=G(m)i=G(m+1), m=1, 2, …, M-1(3)

在进行精细分组之前, 需要引入一个辅助参数p_m, p_m可以根据式(4)获得。

${{p}_{m}}=\left\lfloor m\times \frac{L}{M}-\frac{L}{2M} \right\rfloor , \text{ }m=1, 2, \cdots , M$ (4)

则 xpm可表示第m组的聚类中心。 为加快算法执行速度又保证挖掘相关信息, 对于当前聚类中心, 只考虑{x_j }j=aj=b波段范围内的波段。 其中a, b定义如式(5)所示

(a, b)=[1, pm+1), m=1(pm-1, pm+1), 2≤m≤M-1(pm-1, L], m=M(5)

获得各波段组的${{x}_{{{p}_{m}}}}$之后就可以进行精细分组操作。 首先基于欧氏距离计算x_i与x_j之间相似性, 其中i, j=1, 2, …, L, 并构建相似矩阵D^{L× L}。 判断D(${{x}_{{{p}_{m}}}}$ , x_j)是否小于R(j)(初始值设为无穷大), 其中a≤ j≤ b。 若D(${{x}_{{{p}_{m}}}}$ , x_j)< R(j)则使R(j)=D(${{x}_{{{p}_{m}}}}$ , x_j), 并利用向量T^{1× L}记录对应波段所属波段组m, 即T(j)=m。 将迭代总数设为n, 重复上述步骤。 待迭代次数等于n或向量T不再改变, 停止处理。 最终根据标签向量T, 选择每个波段组中最长波段对应的索引作为新分割波段, 以获得精细分组结果分组流程如图3所示。

图3

Fig.3

	Figure Option ViewDownloadNew Window
	图3 烟草高光谱分组流程图Fig.3 Flowchart of Tobacco Hyperspectral fine grouping

在进行精细分组的同时, 需要利用NMFW计算波段重要性, 以此选择每个波段组中的最佳波段。 在计算NMFW时, 首先将数据矩阵X转置, 此时X=[X₁, X₂, X₃, …, X_d, …, X_N]^T, 1≤ d≤ N, 其中X_d=[x_{d, 1}, x_{d, 2}, …, x_{d, L}]代表第d条样本的光谱数据。

NMFW算法的核心是式(6)

w=κT(K-1(X-m')T)(6)

式(6)中, w=[w1, w2, …, wi, …, wL]T∈RL×N, wi=[wi, 1, wi, 2, …, wi, d, …, wi, N]。 $\text{m}'=\frac{1}{N}\overset{N}{\mathop{\underset{d=1}{\mathop \sum }\,{{X}_{d}}}}\,$是样本的光谱平均值, κ是标准化常数, 计算方法为$\text{ }\!\!\kappa\!\!\text{ }=\frac{1}{\overset{L}{\mathop{\mathop{\sum }_{i=1}}}\,\kappa {{'}_{i}}}$,$\kappa '=\left( X-m' \right){{K}^{-1}}\left( X-m' \right)$,$\kappa '=\left[ \kappa {{'}_{1}},\text{ }\!\!~\!\!\text{ }\kappa {{'}_{2}},\text{ }\!\!~\!\!\text{ }\ldots ,\text{ }\!\!~\!\!\text{ }\kappa {{'}_{i}},\text{ }\!\!~\!\!\text{ }\ldots ,\text{ }\!\!~\!\!\text{ }\kappa {{'}_{L}} \right]\in {{R}^{N\times L}}$,$K=\frac{{{\left( X-m' \right)}^{\text{T}}}\left( X-m' \right)}{N-1}$代表协方差矩阵。

已知滤波器w, 以及光谱数据X={x₁, x₂, x₃, …, x_L}^T∈ R^{L× N}。 光谱数据经过滤波器, 获得输出Y, 可由式(7)表示

Y=wTX=w1x1+w2x2+…+wjxj+…+wLxL(7)

由式(7)可知, Y为L条波段数据的加权和。 其中第j个权重的绝对值|w_j|代表第j个波段对于Y的重要性。

为减少劣质波段对检测结果的影响, 需要对烟草高光谱数据进行标准化。 首先定义一个标准化矩阵A, 如(8)式。

A=1‖X1‖0…00001‖X2‖…000︙︙⋱︙︙︙00…1‖Xd‖…0︙︙⋱︙︙︙00…0…1‖XN‖(8)

式(8)中, A是一个N维对角矩阵, 对角线(d, d)处对应的参数是第d条样本数据的2范数的倒数。 之后利用$\tilde{X}=AX$, $\tilde{m}'=\frac{1}{N}\overset{N}{\mathop{\underset{d=1}{\mathop \sum }\, }}\, {{\tilde{X}}_{d}}$计算参数$\tilde{X}$和$\tilde{m}'$。 基于上述参数可获得标准化常数$\tilde{\kappa }$与标准化权重$\tilde{w}$。

由于共有N条样本, 为分析各波段对整体检测结果的影响, 引入标准化权重的绝对值的平均值$|\tilde{w}{{|}_{avg}}$, 如式(9)

$|\tilde{w}{{|}_{\text{avg}}}=\frac{1}{N}\overset{N}{\mathop{\underset{d=1}{\mathop \sum }\, }}\, |{{\tilde{w}}_{d}}|$ (9)

在完成$|\tilde{w}{{|}_{avg}}$计算之后, 对各个波段进行降序排序, 每个波段获得唯一的序号。 NMFW算法流程如图4所示。 而最终获得的特征光谱即为各波段组中最小序号对应的波段。

图4

Fig.4

	Figure Option ViewDownloadNew Window
	图4 烟草高光谱排序流程图Fig.4 Flowchart of Tobacco Hyperspectral sorting

获取特征波段之后, 为进一步挖掘烟草光谱特征并实现烟草病害分级检测, 采用SVM对特征光谱进行分类。 SVM通过核函数映射将低维空间线性不可分的样本映射到高维空间, 通过在高维空间构造最优分类平面, 实现复杂数据分类与拟合^[16]。 采用RBF核函数, 并利用网格搜索法对SVM的两个关键参数C和σ 进行优化, 两个参数的搜索范围分别为C∈ {10^-1, 10⁰, 10¹, 10², 10³}和σ ∈ {2^-3, 2^-2, 2^-1, 2⁰, 2¹, 2², 2³}。 根据文献[17]及实验观察, 在给定的参数搜索区间内可得到较高的验证精度。

波段分组的数目会影响最终的波段选择数量, 因此, 首先研究不同波段分组数对检测结果的影响。 实验样本共有2 151个波段, 波段分组数变化范围为[2, 800], 不同分组数的检测结果如图5(a— j)所示。

图5

Fig.5

Figure Option ViewDownloadNew Window

图5 波段分组实验结果图 (a), (b), (c), (d)表示不同感病样本的生产者精度分类结果, (e), (f), (g), (h)是用户精度(UA)分类结果, (i), (j)感病烟草总体Kappa系数与OA随波段分组数的变化图Fig.5 Experimental results of band grouping (a), (b), (c) and (d) show the producer accuracy (PA) of tobacco disease detection algorithm; (e), (f), (g) and (h) show the user accuracy (UA) of tobacco disease detection algorithm; (i) and (j) show the kappa and OA oftobacco disease detection algorithm

由图5可知, 随着波段分组数的上升, 各类样本检测精度整体呈上升趋势, 且在波段分组数较少时上升较快。 OA在特征光谱数仅为20时即已突破80%。 除此之外, 当特征光谱数上升到一定数量之后, TMV1(PA), Health(UA)以及PVY1(UA)随特征光谱数上升而下降。 这表明少量光谱即可实现较高精度的病害检测, 光谱数量过多反而不利于烟草病害检测。 因此, 选择特征光谱对于烟草病害精准防治是必要的。

为选择合适的波段分组数, 提取了各类变化曲线的关键节点, 如表4所示。 表4中精度指标的最大值点表示该指标第一次取得最大值时对应的分组数; 平稳点为精度趋于平稳时对应的波段分组数, 当波段分组数选择大于平稳点时, 指标变化小于3%。 最优波段分组数的选择基于以下标准: (1)对应较高的检测精度。 (2)从稳定的曲线中选择。 (3)能够兼顾其他类别的精度。 (4)满足以上条件的最小波段分组数。

表4

Table 4

表4(Table 4)

表4 各参数曲线关键点 Table 4 Key points of each parameter curve 精度指标 最大值点 最大值点 精度/% 平稳点 平稳点 精度/% Health(PA) 674 95.45 330 90.56 TMV1(PA) 126 96.03 50 91.27 TMV3(PA) 111 96.97 40 92.42 TMV5(PA) 92 89.19 75 83.78 TMV7(PA) 54 86.21 352 75.86 TMV9(PA) 97 85.71 608 75.00 PVY1(PA) 163 90.48 113 84.52 PVY3(PA) 118 94.12 92 89.41 PVY5(PA) 283 99.12 116 93.86 PVY7(PA) 26 100 21 97.62 PVY9(PA) 254 100 82 94.77 Health(UA) 283 89.07 94 83.08 TMV1(UA) 127 95.52 76 90.77 TMV3(UA) 105 100.00 200 94.12 TMV5(UA) 785 94.12 无 无 TMV7(UA) 3 100.00 582 95.45 TMV9(UA) 37 100.00 346 95.65 PVY1(UA) 118 93.67 392 89.29 PVY3(UA) 175 95.12 105 91.76 PVY5(UA) 308 100.00 93 98.13 PVY7(UA) 195 100.00 21 94.81 PVY9(UA) 52 100.00 23 95.35 OA 397 93.72 58 87.56 Kappa 397 94.46 58 89.01

表4 各参数曲线关键点 Table 4 Key points of each parameter curve

为使选择最佳波段分组数时各类别均可获得较高精度, 将最大值点和平稳点的最大值作为候选最优波段分组数。 当最大值点小于平稳点时, 最大值点可能具有偶然性, 因此首先排除此类最大值点。 TMV5(UA)变化曲线无平稳点, 因此, 最大平稳点无法代表全部病害检测精度。 OA以及Kappa系数的最大值点相同, 不妨选取OA的最大值点。 综上, 最终将674与397作为候选最佳波段分组数, 基于候选最佳波段分组数获得的结果如表5所示。

表5

Table 5

表5(Table 5)

表5 不同候选波段分组数烟草病害检测精度比较 Table 5 Comparison of tobacco disease detection accuracy of different candidate band group numbers 波段分组数 OA/% AA/% kappa 建模时间/s 397 94.46 91.07 0.937 2 1.14 674 93.97 90.33 0.9316 2.61

表5 不同候选波段分组数烟草病害检测精度比较 Table 5 Comparison of tobacco disease detection accuracy of different candidate band group numbers

由表5可知, OA的最大值点在精度、 稳定性以及建模时间上均获得最优结果。 为实现烟草病害精准识别, 采用OA的最大值点对应的分组数— — 397作为最终波段分组数。

为了分析FNG-NMFW算法的性能, 将其与全波段病害检测算法进行对比。 两个模型均采用SVM分类器, 核函数参数C与σ 由网格搜索法获得, 分别为10²和2^-3。 实验结果如表6所示。

表6

Table 6

表6(Table 6)

表6 不同波段选择方法下的SVM建模结果 Table 6 SVM modeling results under different band selection methods 波段选择 方法 特征 光谱 烟草病害识别精度/% OA/% AA/% kappa 数据处理 时间/s Health TMV1 TMV3 TMV5 TMV7 TMV9 PVY1 PVY3 PVY5 PVY7 PVY9 FNG-NMFW 397 93.9 94.4 95.5 86.5 72.4 78.6 89.3 93.0 98.3 100 100 94.46 91.07 0.937 2 1.1 全波段 2 151 96.1 90.5 92.4 83.8 72.4 67.9 82.1 88.2 98.3 100 100 92.80 88.34 0.918 2 12.9

表6 不同波段选择方法下的SVM建模结果 Table 6 SVM modeling results under different band selection methods

由表6可知, 全波段病害检测算法的OA为92.80%, FNG-NMFW达94.46%。 在运行效率上, 相较于全波段算法, 虽然FNG-NMFW增加了波段选择步骤, 但是由于选择的波段较少(特征光谱数由2 151减为397), 大大提高了运行效率, 数据处理时间由12.9 s降到1.1 s。 综上所述烟草病害检测算法更有助于实现高精度烟草病害检测。

由表6中各类样本的检测精度可知, 虽然全波段分类在Health类别的精度高于FNG-NMFW算法, 但是在实际生产中, 烟草患病种类与患病等级的精准识别对烟草病害的精准防治至关重要。 可以看到, FNG-NMFW算法在不同患病等级的烟草检测上均优于全波段分类算法。 由图1(a)可知, 处于患TMV与PVY不同阶段的烟草的光谱曲线差别较小, 因此具有较高的检测难度。 但是, 通过特征光谱的选择, 本算法能够对不同病害种类进行有效识别, 对早期烟草患病等级(患病等级1, 3), 除PVY1的89.3%外, 其余患病等级检测均高于93%。 因此, FNG-NMFW具有较好的实用价值。

为进一步验证FNG-NMFW的稳定性, 将其与全波段算法在1%~100%的训练集比例下进行实验, 实验结果如图6所示。 当训练集比例较高时, 两算法精度相差较小, 但当训练集比例处于1%~75%时, FNG-NMFW与全波段算法相比具有明显优势。 这说明FNG-NMFW在小样本情况下能够实现精度的有效提升。 在实际生产中, 因训练样本不足而造成烟草病害检测精度较低的问题是普遍存在的, 而FNG-NMFW能够通过选择特征光谱有效缓解这一问题。 因此, FNG-NMFW在烟草病害精准识别领域具有较好的应用前景。

图6

Fig.6

	Figure Option ViewDownloadNew Window
	图6 不同训练集比例下OA与Kappa系数变化图 (a): OA变化图; (b): Kappa系数变化图Fig.6 Variation chart of overall accuracy OA and kappa coefficient under different training set proportion (a): Overall accuracy variation diagram; (b): kappa coefficient variation diagram