化橘红胎切片购于广东省化州市某大药材公司。 各年份化橘红切片样本如图1所示, 可见切片色泽及大小与贮藏时间有一定联系, 但肉眼无法准确判断其年份。

图1

Fig.1

	Figure Option ViewDownloadNew Window
	图1 各年份化橘红胎切片Fig.1 Slices of citri grandis fructus immaturus aged with different years

利用实验室搭建的高光谱成像检测系统, 包括: 手持智能型高光谱相机SPECIM IQ, 高亮度光源, 白色标定板等。 高光谱相机可探测波段范围为400~1 000 nm, 输出高光谱图像大小为512× 512像素, 成像距离为150 mm至无穷远。 如图2所示, 在黑色载物板上均匀摆放化橘红胎切片, 使两侧光源对准载物板, 调节高光谱相机水平并对焦后, 触发相机快门开始扫描。

图2

Fig.2

	Figure Option ViewDownloadNew Window
	图2 高光谱成像平台Fig.2 Hyperspectral imaging platform

采集到的原始高光谱数据需经式(1)进行黑白校正, 本实验使用的高光谱成像仪内置黑白校正功能。

O=I-H0H1-H0(1)

式(1)中, I为原始光谱图像, O为校正后光谱图像, H₁为扫描白色标定板(反射率约为1)得到的标定图像, H₀为关闭光源扫描黑色载物板(反射率约为0)得到的标定图像。

读取黑白校正后的高光谱图像并分割出切片图像, 获取整片切片截面为感兴趣区域(region of interest, ROI)。 提取每一份切片400.2~997.4 nm波段的平均反射光谱。 剔除各年份的化橘红切片平均反射光谱离群样本, 得到共176个样本的平均光谱反射曲线。 采用Kennard-Stone算法对每个年份的化橘红切片平均反射光谱以3∶ 1的比例提取训练集与测试集, 该算法基于变量间的欧几里得距离挑选样本进入训练集直到数量达到要求, 能涵盖数据集中大多数差异, 确保训练模型更能代表整个数据集。

首先, 对原数据集进行PCA聚类分析, 探究PCA聚类直接区分化橘红切片的可行性。 另外, 对数据集分别采用SG平滑、 MSC、 SNV三种预处理方法, 观察平均反射光谱曲线变化, 分别建立KELM分类判别模型并对比准确率, 选择判别准确率最高的预处理方法。 最后, 将选得的数据集分别经CARS算法和MCUVE算法筛选特征波长, 进行KELM和HHO-KELM建模, 比较分类效果。 数据处理和建模的大致流程, 如图3所示。

图3

Fig.3

	Figure Option ViewDownloadNew Window
	图3 数据分析步骤Fig.3 Data analysis steps

1.4.1 核极限学习机

假设有n个样本, 以x_i表示每个样本的特征向量; y_i表示每个样本期望的响应, 对分类器而言y_i是样本类别标签; 有L个隐藏节点; G(x)表示激励函数。 那么, 单隐层前馈网络的输出可以表示为

fout(x)=∑i=1LβiG(ωixj+bi)=∑i=1Lβih(xj), j=1, 2, …, n(2)

式(2)中, ω _i为连接输入层节点到第i个隐藏节点的权重; β _i为连接第i个隐藏节点到输出节点的权重; b_i为第i个隐藏节点的阈值; f_out(x)为第j个输入特征向量由网路求得的输出值。 以H表示隐藏层输出, β 表示隐藏层到输出层的输出权重, T表示目标输出。

H(ω1, …, ωL, b1, …, bL, x1, …, xn)=g(ω1x1+b1)…g(ωLx1+bL)︙…︙g(ω1xn+b1)…g(ωLxn+bL)n×L(3)

β=β1︙βL(4)

T=y1︙yL(5)

标准的ELM在训练前随机设定且在训练过程中不更新W和B的参数值, 训练目标为输出结果与期望响应差值最小

‖Hβ-T‖2=0(6)

解出输出层权重β 即完成训练, β 取最小二乘解

β=H+T=HT(HHT)-1T=HTIC+HHT-1T(7)

式(7)中, C为惩戒系数; I为单位矩阵; H为隐藏层矩阵; H⁺为H的广义逆矩阵。 经过训练标准ELM模型的输出为

f(x)=h(x)β=h(x)HTIC+HHT-1T(8)

为避免使用随机参数的不稳定隐患, 引入核矩阵Ω _ELM和核函数k(x_i, x_j)。 采用径向基函数作核函数, 使模型更稳定并提高数据区分度。 式中, σ 为径向基函数的宽度系数。

ΩELM=HHT(9)

ΩELMi, j=h(xi)h(xj)=k(xi, xj)=exp-‖xi-xj‖22σ2(10)

基于径向基函数的核极限学习机的输出函数可表示为

f(x)=k(x, x1)︙k(x, xN)TIC+ΩELM-1T(11)

对于KELM, 调节惩戒系数C和核函数参数σ 即影响模型训练效果。

1.4.2 哈里斯鹰算法优化参数选择

HHO是由Heidari等在2019年提出的一种模拟哈里斯鹰捕食行为实现的优化算法, 具有全局搜索范围广, 收敛速度快等优点。 HHO算法以猎物位置描述待寻优参数, 以猎物逃逸能量描述优化目标函数, 算法寻优过程模仿哈里斯鹰捕食行为, 分为追踪、 围攻、 捕捉等几个阶段进行迭代, 逐步包围猎物, 使目标函数逐步减小。 当迭代次数达到要求时, 哈里斯鹰种群将猎物逼迫在包围圈中, 目标函数达到迭代过程的最小值^[15]。

KELM分类模型的分类效果受惩戒域子C、 核函数宽度系数σ 影响, 因而可采用HHO算法对C、 σ 进行优化。 直接以KELM的分类错误率作为优化的目标函数容易陷入局部最优或出现过拟合现象, 笔者对KELM训练过程加入K-Fold交叉验证, 并以交叉验证得到的均方根值为优化目标函数(适应度值)。 K-Fold交叉验证是常用的模型评估技术, 它将数据集分成K个相等的子集, 然后进行K次模型训练和测试, 每次使用不同的子集作为测试集, 其余的子集作为训练集。 该技术可以减少模型的方差, 提高模型的泛化能力, 因此可以被用于优化算法的参数选择^[16]。 采用K值为10进行交叉验证。

HHO-KELM算法步骤如下:

step1: 初始化HHO种群大小、 迭代次数;

step2: 训练集作输入, 交叉验证训练并测试KELM, 得到均方根值;

step3: 未达到迭代次数, 则采用HHO对C、 σ 进行寻优, 更新KELM;

step4: 达到迭代次数退出HHO, 输出最终适应值和寻优参数, 更新KELM;

step5: 测试集输入得到的HHO* KELM测试训练效果。

采集波长范围为397~1 003 nm的204个可见近红外光谱数据, 因首尾部噪声明显, 截取400.2~997.4 nm的201个波长的光谱数据。 所得全样本原始光谱曲线和按年份分组平均光谱曲线如图4所示。

图4

Fig.4

	Figure Option ViewDownloadNew Window
	图4 化橘红胎切片原始光谱图(a)及平均光谱图(b)Fig.4 Original reflection spectra (a) and mean spectra (b) of samples

由图4可知, 贮藏时间越长的化橘红胎切片可见近红外波段反射率越小, 各组化橘红胎切片的反射光谱曲线线性相近, 有所区别。 但观察图4(a), 除2021年化橘红胎切片的光谱曲线有较为明显的区分外, 其他三个年份化橘红胎切片得到的反射光谱曲线在整个波段范围内都有不同程度的交叠。 因此, 无法直接由反射光谱数据鉴别化橘红胎切片的生产年份。

对化橘红胎切片原始反射光谱进行PCA聚类分析, 所得主成分贡献率信息如表1所示, PC1和PC2累计贡献率达到98.25%, 足以描述化橘红胎切片原始光谱特征。 采用PC1和PC2作得分图, 由图5可以看出, 各年份样本点有明显聚类, 2021年产化橘红胎切片区别明显, 而2018年、 2015年、 2012年产化橘红胎切片样品点有重合, 难以区分。 根据对样品外表观察, 推测原因是2021年产化橘红胎切片未完全脱水, 含水量偏高, 表面未有挥发油类物质析出, 而另外三组化橘红胎切片干制较完全且化学成分变化不明显。 所以, 无法通过主成分分析结果区分生产年份三年以上的化橘红胎切片样本, 仍需要进一步研究。

表1

Table 1

表1(Table 1)

表1 PCA主成分分析 Table 1 Principal component analysis PC summary Eigenvalue Proportion of variance/% Cumulative proportion of variance/% PC1 191.7 95.37 95.37 PC2 5.785 2.88 98.25 PC3 3.099 1.54 99.79 PC4 0.306 4 0.15 99.94

表1 PCA主成分分析 Table 1 Principal component analysis

图5

Fig.5

	Figure Option ViewDownloadNew Window
	图5 PCA主成分得分分布Fig.5 Principal component score distribution

依次采用S-G平滑、 MSC、 SNV方法对化橘红胎切片原始反射光谱曲线进行预处理, 处理后反射光谱曲线如图6所示。 分别将无预处理和预处理得到的共四份高光谱数据训练KELM并分别验证训练集与测试集的预测准确率, 得到结果如表2所示。 经比较, 选取建立KELM模型效果相对最好的化橘红胎切片 SNV 处理反射光谱数据作进一步的分析和处理。

图6

Fig.6

	Figure Option ViewDownloadNew Window
	图6 不同预处理方法化橘红胎切片光谱图 (a): S-G平滑; (b): MSC; (c): SNVFig.6 Pretreated spectra of samples (a): S-G smoothing; (b): MSC; (c): SNV

表2

Table 2

表2(Table 2)

表2 不同预处理方法KELM 判别模型效果 Table 2 Accuracies of KELM discriminant models with different preprocessing methods Pretreatment Accuracy rate Training set/% Test set/% NONE 90.84 80.00 S-G 90.84 80.00 MSC 94.66 91.11 SNV 99.24 95.56

表2 不同预处理方法KELM 判别模型效果 Table 2 Accuracies of KELM discriminant models with different preprocessing methods

为比较CARS与MCUVE筛选特征波长后的创建判别模型的效果, 由CARS筛选得到19个特征波长, 对应地选用MCUVE得到的重要性系数最高的前19个波长作为特征波长。 经CARS和MCUVE所得特征波长信息如表3所示。 由图7, CARS所选特征波长集中于两端, MCUVE所选特征波长分布相对均匀。 两种方法得到的特征波长在500~800 nm波段范围有明显差别, 具体效果需要建立分类判别模型进行实际对比。

图7

Fig.7

	Figure Option ViewDownloadNew Window
	图7 特征波长分布Fig.7 Characteristic wavelength distribution

表3

Table 3

表3(Table 3)

表3 特征波长信息 Table 3 Feature wavelength information MODEL Wavelength selection/nm CARS 449.35 455.16 501.72 513.4 598.6 601.55 604.51 853.33 874.64 889.9 905.18 920.48 938.88 945.02 966.55 981.96 991.22 994.31 997.4 MCUVE 400.2 403.09 648.95 651.92 654.89 726.53 771.61 777.64 780.65 850.29 920.48 938.88 945.02 954.24 975.79 981.96 991.22 994.31 997.4

表3 特征波长信息 Table 3 Feature wavelength information

以全波段、 CARS和MCUVE方法筛选特征波长得到的数据集创建KELM模型, 分别检验训练集和测试集分类判别效果。 初始化HHO算法参数, 种群数量为30, 迭代次数500, 对KELM的惩戒系数C和核函数参数σ 寻优。 待HHO寻优结束后, 采用更新的C和σ 重新训练KELM并测试判别效果。 得到不同波段范围创建KELM和HHO-KELM分类判别模型的测试结果如表4所示。

表4

Table 4

表4(Table 4)

表4 分类判别模型结果 Table 4 Classification discriminant modeling results Wavelength selection method Model Accuracy rate Training set/% Test set/% NONE KELM 99.24 95.56 HHO-KELM 100 100 CARS KELM 95.42 93.33 HHO-KELM 93.89 93.33 MCUVE KELM 100 100 HHO-KELM 100 100

表4 分类判别模型结果 Table 4 Classification discriminant modeling results

全波段信息最为充足, 但信息的冗余会带来很大的运算压力, 使训练时间变长, 且可能使判别模型劣化。 通过波段筛选可以有效减少信息冗余, 加快训练速度, 可能提升判别效果。 对化橘红胎切片样品高光谱采用波段筛选后, 训练速度和判别效果有所提升。 并且, 在保留特征波长数量相同的前提下, MCUVE比CARS效果更好, 无论优化与否都能达到最高的判别准确率, 认为MCUVE筛选的19个特征波长可以描述原光谱数据的全部信息。 可以看出, 采用哈里斯鹰算法优化参数的核极限学习机对训练集检验的准确率有明显的提高, 且优化效果在测试集中也有所体现。 无论是在数据集包含全波段还是少数特征波长的情况下, HHO-KELM都能很好的通过高光谱数据实现化橘红胎切片的年份鉴别。