烟叶样本来源于中国八个省份, 包括福建省、 广西壮族自治区、 贵州省、 河南省、 湖北省、 湖南省、 四川省和云南省, 共824份烟叶样本, 由浙江中烟工业有限责任公司提供。 实验采用Antaris Ⅱ FT-NIR(Thermo Fisher Scientific)分析仪, 波数范围为10 000~3 800 cm^-1, 分辨率为8 cm^-1, 每条光谱代表64次扫描的平均值。 由于八个烟叶产区地理位置比较接近, 因此产地来源分辨具有一定难度。

为了消除基线漂移, 降低噪声干扰, 同时又能保证校正后不会过多地降低频谱的信噪比。 采用Savitzky-Golay(SG)卷积平滑结合一阶导数(first derivative, 1^st D)对光谱进行预处理。

为了获得更稳定的近红外判别分析模型, 剔除光谱中的无关变量^[9], 分别采用竞争自适应重加权采样(competitive adaptive reweighted sampling, CARS)和随机青蛙(random frog, RF)算法对全波长的特征变量进行提取以降低数据维度, 从而可以输入较少变量以提高分析速度。

1.4.1 支持向量机(SVM)

SVM是Vapnik等在1995年基于统计学习理论(statistic learning theory, SLT)所提出的新型模式识别方法。 SVM可以很好地解决实际中线性不可分问题, 它以最大间隔超平面为原则, 首先将非线性原始数据映射到一个新的高维空间, 然后在新的维度空间利用核函数从训练数据中学习分类模型, 找到恰好能将数据分开的超平面, 这样就将问题成功转化为高维线性可分问题^[10]。

1.4.2 粒子群优化(PSO)

PSO算法又称为鸟群觅食算法。 在PSO中, 粒子是随机生成的, 这些粒子会根据当前的特征行为和位置动态变化随时调整搜索位置, 每组粒子会跟踪先前的最佳位置和适应度进行最优解的寻找, 而每一个粒子找到的最优解会引导其他粒子收敛到最佳方案, 利用粒子群之间的竞争和协作完成每一次迭代, 在粒子不断演化的过程中, 根据当前最优值搜索到适应性最大的解, 则成为全局最优解^[11]。

1.4.3 遗传算法(GA)

遗传算法是基于生物在自然界中“ 适者生存, 优胜劣汰” 的进化机制所演化出的一种寻找全局最优解的搜索算法, 通过模拟自然界中生物进化、 遗传变异来实现全局寻优^[12]。 在初始种群中, 某些适应度大的个体被选为父代, 父代间进行交叉和突变操作, 以此产生适应度更大的子代, 组成新种群, 再评估每一个体的适应度值, 选择值得进行下一代遗传的个体, 进行下一次迭代。 为了找到合适的遗传个体, 每一个个体都会经过预先设定的适应度函数进行选择, 以使得新一代的染色体比上一代有所提升, 在这样反复的遗传演化过程, 好的基因被遗传, 呈现劣势的基因被淘汰。

1.4.4 灰狼算法(GWO)

GWO的提出是受到了来自于自然界中灰狼生存等级和捕猎的启发^[13], GWO模拟了灰狼族群中严格的统治阶级, 并通过灰狼种群捕猎阶段中跟踪、 追逐、 包围、 攻击等行为寻找最优解, 以此实现对支持向量机等智能算法的优化。 采用GWO优化函数进行建模时, α 表示最合适的解决方案, 其次分别是β 和δ , β 和δ 主要负责辅助α 进行决策, 并且听从α 的指令进行猎捕, 狼群最低等级是ω , 它们必须跟随α 、 β 和δ 的领导对猎物进行围攻, 以此寻找全局最优解。

在狩猎过程中, 狼群的包围行为由式(1)和式(2)表示

D⃗=|C⃗X⃗p(t)-X⃗(t)|(1)

X⃗t+1=X⃗pt-A⃗D⃗(2)

式中, D⃗为当前灰狼与猎物间的距离; A⃗和 C⃗为系数向量; t为当前迭代次数; X⃗p为猎物位置向量。

A⃗=2a⃗r⃗1-a⃗(3)

C⃗=2r⃗2(4)

式中, a⃗会随着迭代次数从2线性减小到0; r⃗1和 r⃗2是[0, 1]之间的随机向量。

在搜索空间, α , β 和δ 这三个最优解的位置起着关键的作用, 由这三者预估猎物的大致位置, 而ω 会根据这三者更新自己在猎物周围的位置, 以此实现对猎物的最终包围狩猎。

D⃗α=|C⃗1X⃗α(t)-X⃗(t)|D⃗β=|C⃗2X⃗β(t)-X⃗(t)|D⃗δ=|C⃗3X⃗δ(t)-X⃗(t)|(5)

X⃗1=X⃗α-A⃗1D⃗αX⃗2=X⃗β-A⃗2D⃗βX⃗3=X⃗δ-A⃗3D⃗δ(6)

X⃗(t+1)=X1+X2+X33(7)

式(5)中, X⃗α, X⃗β和 X⃗δ分别表示α , β 和δ 狼当前的位置向量, D⃗α, D⃗β和 D⃗δ分别是三只狼与猎物之间的距离。 由式(6)可知, ω 狼的位置在α , β 和δ 狼指导下不断调整以更靠近猎物, X⃗(t+1)表示ω 狼最终的位置。

1.4.5 GWO-SVM

采用选择高斯核函数作为支持向量机的核函数, 此时SVM受到惩罚因子C和高斯核参数σ 的影响。 在此, 使用GWO算法来对SVM进行参数优化, 通过不断迭代优化, 寻找最佳的参数C和σ , 提高模型的准确率, 降低识错率^[14]。 GWO-SVM算法过程如下:

(1) 将烟叶样本数据划分为训练集和测试集, 并对光谱数据进行归一化处理。

(2) 初始化SVM的参数, 并设置SVM中C和σ 作为GWO优化值。

(3) 根据C和σ 初始化狼群的位置, 使用训练样本对SVM模型进行训练, 并将初始位置输入模型, 以计算灰狼个体的适应度值。

(4) 保留具有最优适应度值的α , β , δ 三只狼, 并根据当前最优适应度值, 更新每只狼的位置。

(5) 若已达到最大迭代次数, 则终止迭代循环, 输出由最优位置得出的全局最优参数C和σ , 若未到达最大迭代次数, 则返回步骤4继续进行迭代。

(6) 输出最优参数C和σ , 并对训练数据重新进行SVM训练。

(7) 通过GWO-SVM模型对测试样本进行预测, 并对模型性能进行评价。

以上数据处理均在Windows 10操作系统下, 使用MATLAB R2018 b进行。

图1用不同颜色的曲线展示了824份不同产区烟叶样品的原始近红外漫反射光谱, 八个省份样本的原始近红外光谱并无太大差距, 可见烟叶产地的鉴别有一定难度。 使用基于x-y距离样本集划分(set partitioning based on joint x-y distance, SPXY)算法按照7∶ 3的比例进行样本划分, 将824个烟叶样本划分为617个校正集样本和207个验证集样本。

图1

Fig.1

	Figure Option ViewDownloadNew Window
	图1 烟叶样品的原始近红外光谱图Fig.1 Raw NIR spectra of tobacco leaf samples

表1

Table 1

表1(Table 1)

表1 不同产地烟叶样本划分情况 Table 1 Training and testing sample sizes of tobacco from different regions 分类编号 产地 样品编号 校正集 验证集 总数 1 福建 98 40 138 2 广西 52 8 60 3 贵州 68 26 94 4 河南 79 33 112 5 湖北 55 44 99 6 湖南 96 16 112 7 四川 86 17 103 8 云南 83 23 106

表1 不同产地烟叶样本划分情况 Table 1 Training and testing sample sizes of tobacco from different regions

2.2.1 CARS法

在CARS算法中, 将蒙特卡洛采样运行次数设置为100次, 通过五折交叉验证选择具有最小均方根误差(root mean square error of cross validation, RMSECV)模型所对应的波长变量。 图2(a)显示了蒙特卡洛采样次数和采样变量数之间的关系, 前20次采样中, 采样变量数快速下降, 主要反映的是消除无信息变量的过程, 下降至一定程度后趋于平缓。 从图2(b)中可以观察到当采样次数达到37时, 对应的RMSECV最小为1.72。 图2(c)显示的是随着采样次数不同而变化的变量系数, 图中用“ * ” 描述了最小RMSECV所对应的最佳特征波长。 结果, CARS选择了141个变量。

图2

Fig.2

	Figure Option ViewDownloadNew Window
	图2 CARS变量筛选优化结果 (a): SG+1st D处理后的波长数随采样运行次数而变化; (b): SG+1st D处理后的RMSECV随采样运行次数而变化; (c): SG+1st D处理后变量系数随采样运行次数变化Fig.2 CARS variable screening results (a): Number of wavelengths changed with sampling runs; (b): RMSECV changed with sampling runs; (c): Variable coefficient changed with sampling runs

2.2.2 RF法

随机青蛙算法通过计算每个建模变量被选择的概率, 从中选出被选概率高的变量组合成为特征波长; 由于随机青蛙采用蒙特卡洛随机采样, 因此为了降低随机性, 将所有波长变量运行1 000次, 计算平均筛选概率^[15], 设置0.04为筛选阈值, 最终从1 609个变量中选出534个变量作为特征变量, 结果如图3所示。

图3

Fig.3

	Figure Option ViewDownloadNew Window
	图3 RF变量筛选结果Fig.3 RF variable screening results

为了准确对比各优化算法的分类效果, 考虑了不同的控制参数, 三种优化算法均设置相同的参数: 种群规模为20, 最大进化次数k为200, C和σ 的搜索范围为(0.01, 100)。

分类模型的性能使用准确度和混淆矩阵来评估, 混淆矩阵主要用于比较客观结果和实际类别。 下面分别计算了PSO-SVM, GA-SVM和GWO-SVM模型的分类准确度和混淆矩阵, 结果见表2和图4。 混淆矩阵用横坐标表示真实类别, 纵坐标表示预测结果, 以此描述模型对单个产地的判断分类性能, 颜色越深, 代表正确分类的数量越多, 图4中对角线的颜色越深, 代表准确率越接近100%, 分类性能越好。

表2

Table 2

表2(Table 2)

表2 不同模型分类效果 Table 2 Classification effects of different models 建模方法 波长选择方法 变量数 最佳C值 最佳σ 值 运行时间/min 校正集正确率/% 验证集正确率/% PSO-SVM CARS 141 37.07 0.01 47 94.00 87.92 RF 534 33.00 0.01 178 99.68 95.17 GA-SVM CARS 141 4.61 0.37 39 100.00 91.30 RF 534 15.58 0.03 153 100.00 94.20 GWO-SVM CARS 141 6.07 0.11 10 100.00 93.24 RF 534 34.48 0.10 22 100.00 96.62

表2 不同模型分类效果 Table 2 Classification effects of different models

图4

Fig.4

	Figure Option ViewDownloadNew Window
	图4 不同分类模型的混淆矩阵 (a): RF-PSO-SVM; (b): RF-GA-SVM; (c): RF-GWO-SVMFig.4 Confusion matrices of different classification models (a): RF-PSO-SVM; (b): RF-GA-SVM; (c): RF-GWO-SVM

由表2的结果统计可以看出, 在建立烟叶产地分类的SVM模型过程中, 采用RF筛选出的波长组合, 三种模型的判别准确度更高, 其中RF-GWO-SVM的分类正确率最高, 达到96.62%, 说明RF算法具有良好的关键信息筛选能力; 从运行时间来看, GWO-SVM相比于其他两种模型, 寻优时间大大缩短。 这很好的说明, 在相同的测试样本、 相同搜索条件下, GWO算法对SVM的两个参数寻优具有更快速、 更优化的能力, GWO-SVM模型对于烟叶产地的识别效果明显优于PSO-SVM和GA-SVM。

根据混淆矩阵结果, GWO-SVM模型识别正确率达到100%的有五个省份, 分别是福建省、 河南省、 湖北省、 湖南省和云南省, 河南省的分类正确率比其他模型提升了3.03%, 湖南省的分类正确率比其他模型提升了12.5%~18.75%; 但对于来自广西壮族自治区的样本识别仍然不理想, 原因可能是广西的烟叶样本总数仅60个, 分到验证集仅8个, 面对类别不平衡样本时, SVM模型的分类精度受到较大影响。 另外, 将贵州烟叶错分至湖北和湖南, 可能是因为几个省份地理位置接近, 经纬度相似, 烟叶样本中化学成分与物理特征差异不够明显, 从而造成误分。 除此之外, 各个地区的烟叶产地鉴别都得到了很高的精度, RF-GWO-SVM的总体正确率达到了96.62%, 分别比RF-PSO-SVM和RF-GA-SVM提高了1.45%和2.42%。 从预测结果可以得出结论, GWO-SVM的识别效果明显优于其他两个模型。