原料: 82个松脂样本由广西林业科学院梁忠云工程师负责采集鉴定, 单个样本采集量为50 g, 采样时间为2016年6月— 7月, 采集地点为广西壮族自治区内的武鸣、 防城、 富川、 梧州、 百色、 乐业, 松脂样本包括马尾松松脂31个和湿地松松脂51个。 松脂样本采集后去除树皮、 枝叶等杂质导入样品瓶, 密封, 存于冰箱冷藏待测定。

仪器: i-spec型近红外光谱仪(美国必达泰克公司), 附件: 积分球采集模块。 所有化学计量学算法均通过matlab2010b编程实现, 部分算法通过化学计量学软件ChemSolv实现。

移取适量松脂装入3 cm× 4 cm的透明密封袋, 压平松脂, 制备得到厚度约5 mm, 面积约2 cm× 2 cm的饼状松脂样本, 呈粘稠固体形态。 近红外光谱仪预热30 min后开始测量, 首先关闭光源测得仪器暗电流, 以白色底板测得参比光谱, 再通过积分球以漫反射方式直接采集松脂样本的近红外光谱信息。 实验参数: 扫描波段900~1 700 nm, 积分时间4 000 ms, 扫描次数10次, 单个样本采集光谱数据点511个。 建模过程包括光谱预处理(平滑去噪、 奇异值筛选、 特征选择等)、 模型训练、 模型校验等步骤, 待测样本引入模型中并进行种类识别, 流程如图1所示。

图1

Fig.1

	Figure Option ViewDownloadNew Window
	图1 分类模型建模流程图Fig.1 Flow chart of modeling

按照1.2的实验方法测得松脂样本的近红外光谱。 由于测量过程中存在的光谱噪声以及因样本物理性状差异所引起的散射影响, 建模之前需要对近红外光谱进行光谱预处理。 考虑到松脂样本兼具固体和浆状物的特点, 采用S-G平滑和多元散射校正(multiplicative scatter correction, MSC)相结合的预处理技术, S-G平滑能够有效减小光谱的噪声, 而MSC是基于统计的方法, 能够修正光谱因散射所发生的线性变化。 预处理后的光谱如图2所示。

图2

Fig.2

	Figure Option ViewDownloadNew Window
	图2 预处理后的松脂近红外光谱图Fig.2 NIR spectrum after pre-treatment

主成分分析是一种经典的数据降维技术, 通过将高维数据在低维空间可视化, 从而了解样本分布的基本情况。 从图3可以看出, 82个样本在第一主成分和第二主成分的投影上具有比较明显的两类趋势, 这说明利用近红外光谱作为松脂种类识别具有理论上的可行性。 但是, 马尾松松脂和湿地松松脂依然存在部分重叠, 因此用无监督模式的分类方法很难对两种松脂的种类进行准确识别。

图3

Fig.3

	Figure Option ViewDownloadNew Window
	图3 样本空间的主成分投影Fig.3 Projection of principal component for sample space

奇异值影响模型的性能, 光谱矩阵和响应矩阵两个方面的离群现象均有可能导致奇异值的产生, 也被称为x奇异值和y奇异值, 其中x奇异值主要是由于样本的组成与整体差异过大、 光谱测量误差等因素所导致, 而y奇异值通常是因为样本的指标值测量或识别错误所致。 为了消除奇异值对建模的影响, 采用基于蒙特卡洛采样的奇异值诊断技术对校正集中的奇异值进行筛选(方法原理见文献[7])。 本实验中均值取0.80, 方差取0.181, 超出阈值的12号、 15号、 33号、 35号、 77号样本为奇异值, 见图4, 所有奇异值均不参与建模。 数据集中剩余的样本数为77个(湿地松48个, 马尾松29个)。

去除奇异值的样本被随机划分为训练集与测试集, 训练集样本数50个(18个马尾松松脂和32个湿地松松脂), 测试集样本数27个(11个马尾松松脂和16个湿地松松脂)。 利用训练集的50个样本建模, 测试集中的样本不参与建模的任何过程, 仅用于检验模型的预测能力。

图4

Fig.4

	Figure Option ViewDownloadNew Window
	图4 基于蒙特卡洛采样的异常样本诊断结果示意图Fig.4 Outlier diagnosis based on Monte Carlo sampling

偏最小二乘判别分析(partial least squares discriminant analysis, PLSDA)是一种广泛应用的有监督模式识别方法, 即用已知的样本进行训练, 让计算机从这些样本的数据中“ 学习” , 从而得到一个对同类样本具有预测功能的判别模型, 该方法特别适合于样本数较少、 变量数较多的数据集。 本研究中, 设训练集的光谱矩阵为X_{50× 511}, 即有50个松脂样本, 每个样本有511个光谱数据点; 训练集的响应矩阵为y_{50× 1}, 其中马尾松松脂由“ +1” 表示, 湿地松松脂由“ -1” 表示。

2.3.1 子窗口随机分析变量选择

近红外光谱的信息复杂且不同波长间的信号存在关联, 其中信息变量、 非信息变量以及干扰变量共存。 如何从大量的光谱数据中选择合适的波段进行分析是近红外光谱建模中重要问题, 早期研究中一般认为偏最小二乘利用潜变量建模能够达到去噪和变量选择的目的, 但是最近研究表明, 变量选择与潜变量建模相结合能够简化模型, 提高模型的稳健性和预测精度^{[8, 9]}。

本实验采用子窗口随机化分析法(sub-window permutation analysis, SPA)进行变量选择, 从511个原始光谱数据点中选择了300个变量进行建模。 子窗口随机化分析变量选择是基于模型集群的思想, 通过不同的训练子集和变量子集反复地建立子模型, 通过对大量子模型的统计分析, 筛选出最佳的变量集。 图5将算法选择出的变量与样本的近红外光谱图进行对照, 可以发现有效变量出现的区域与近红外光谱吸收峰的位置一致。

2.3.2 重复双重交叉检验确定潜变量数

偏最小二乘判别分析是一种利用潜变量建模的方法, 因此, 如何选择合适的潜变量数是影响模型性能的关键。 潜变量数选择过少, 模型欠拟合; 选择过多, 则会引起模型过拟合。 交互检验可以有效地利用样本信息, 在进行建模的同时进行样本校验, 并据此获得最佳模型。 重复双重交叉检验(repeated double cross validation, RDCV)是一种新的交互检验方法, 它将蒙特卡洛交叉校验与双重交叉校验结合到一起, 兼具二者的优势, 能够给出预测误差的分布信息。 本实验中利用重复双重交叉检验建立1000个子模型, 统计最优模型潜变量数的最高频度值最终确定为模型的潜变量数(n=7), 结果如图6所示。

图5

Fig.5

	Figure Option ViewDownloadNew Window
	图5 子窗口随机化分析变量选择Fig.5 Variable selection by SPA

图6

Fig.6

	Figure Option ViewDownloadNew Window
	图6 重复双重交叉校验潜变量数频度图Fig.6 Frequency diagram of RDCV

表1

Table 1

表1(Table 1)

表1 偏最小二乘判别分析模型参数及结果 Table 1 Model parameters and results of PLSDA 样本 编号 建模方法及参数 真实类别 预测类别 准确率 /% 1 马尾松松脂 马尾松松脂 2 马尾松松脂 马尾松松脂 3 马尾松松脂 马尾松松脂 4 马尾松松脂 马尾松松脂 5 马尾松松脂 马尾松松脂 6 马尾松松脂 马尾松松脂 7 Ⅰ . 光谱预处理: 马尾松松脂 马尾松松脂 8 ● S-G平滑 马尾松松脂 马尾松松脂 9 ● 多元散射校正 马尾松松脂 马尾松松脂 10 Ⅱ . 奇异值筛选: 马尾松松脂 马尾松松脂 11 ● 蒙特卡洛采样奇异值诊断 马尾松松脂 马尾松松脂 12 Ⅲ . 波长选择: 湿地松松脂 湿地松松脂 13 ● 子窗口随机化分析法 湿地松松脂 湿地松松脂 96.30 14 ● 光谱数据点300个 湿地松松脂 湿地松松脂 15 Ⅳ . 交互检验: 湿地松松脂 湿地松松脂 16 ● 重复双重交叉检验法 湿地松松脂 湿地松松脂 17 ● 潜变量数目n=7 湿地松松脂 马尾松松脂 18 Ⅴ . 校正集分类准确率: 湿地松松脂 湿地松松脂 19 ● 96.88% 湿地松松脂 湿地松松脂 20 Ⅵ . 交互检验分类准确率: 湿地松松脂 湿地松松脂 21 ● 96.10% 湿地松松脂 湿地松松脂 22 湿地松松脂 湿地松松脂 23 湿地松松脂 湿地松松脂 24 湿地松松脂 湿地松松脂 25 湿地松松脂 湿地松松脂 26 湿地松松脂 湿地松松脂 27 湿地松松脂 湿地松松脂

表1 偏最小二乘判别分析模型参数及结果 Table 1 Model parameters and results of PLSDA

2.3.3 分类模型建立及评价

经过光谱预处理、 奇异值筛选、 波长选择和确定潜变量数等步骤后, 利用偏最小二乘判别分析法建立模型, 并通过测试集中的未知样本对模型的预测能力进行评估。 表1可见, 测试集中的27个样本, 除17号湿地松松脂被错判为马尾松松脂, 其余26个样本均能正确识别, 准确率达到96.30%。