基于模型集群分析的激光诱导击穿光谱变量选择方法

摘要

本发明公开了一种基于模型集群分析思想适用于支持向量机的变量选择方法，其通过蒙特卡洛采样从全光谱数据矩阵中获取子数据集，针对每个子数据集建立一个SVM子模型并预测分类，然后使用Mann-WhitneyU检验对所有子模型的预测正确率进行统计分析，挑选出对模型预测能力有显著作用的有用变量。该方法不以一次性建模结果为依据，而是通过有放回的重采样最大限度地有效利用数据信息，充分考察数据集中各变量间的内在关系，对不同结果的统计分布进行分析，因此具有更好的普遍性和稳定性。

说明书

技术领域

本发明涉及一种基于模型集群分析的激光诱导击穿光谱变量选择方法，属于光谱分析技术领域。

背景技术

激光诱导击穿光谱(laser-induced breakdown spectroscopy, LIBS)是一种基于原子发射光谱的检测物质组分与含量的分析技术。强激光脉冲聚焦在样品上形成等离子体，在等离子体冷却过程中，样品中处于激发态的原子和离子向低能级或基态跃迁产生特定频率的特征发射谱线。LIBS分析简便、快速，不需要样品预处理并能同时进行多元素测定，因此在多个领域被广泛应用。由于特征发射谱线和被测样品的元素组成成分之间有着一一对应的关系，因此根据不同种类物质LIBS光谱的特征谱线，借助化学计量学方法能够判别物质所属类别从而实现对物质的分类，即化学模式识别。

支持向量机（Support vector machine, SVM）是一种基于统计学习理论的机器学习算法。它以训练误差作为优化问题的约束条件，以置信范围值最小化作为优化目标，即SVM是一种基于结构风险最小化准则的学习方法，其推广能力明显优于基于经验风险最小化的传统机器学习方法。支持向量机通过引入核函数将输入样本空间非线性映射到新的高维特征空间，在高维空间中进行相应的线性操作，从而实现非线性关系向线性关系的转化，在解决小样本、非线性及高维模式识别中表现出许多特有的优势。因此将支持向量机和激光诱导击穿光谱技术相结合能够很好地对物质进行判别分类。

通常LIBS光谱都是高维数据，而且在光谱中存在大量的无用变量。这些无用变量主要是噪音和测试中的波动信号，与元素成分组成没有关系，对判别分析也没有帮助。因此使用全光谱会影响到分类效果，同时过大的数据量对计算成本也提出了很高的要求，不利于现场快速分析。另一方面，全光谱中包含测试环境中的大量信息，这些信息也会对建模分类产生影响。当测试环境改变时，环境信息也会变化从而严重影响分类模型的稳定性。因此使用一种专属于支持向量机的变量选择方法将光谱中真正能够代表物质成分特征并对分类有重要作用的变量挑选出来就显得很有必要了。

但目前还没有一种普遍能被接受和认可并适用于支持向量机的变量选择方法，传统的遗传算法、粒子群算法等优化算法虽然也能用于提取变量，但这些算法计算繁琐、耗时较长，还需要设置大量初始参数，而且容易出现局部最优解。本发明提出的基于模型集群分析的变量选择方法通过重采样对各种变量组合情况统计分析，直接分析各变量对预测结果的影响，不会陷入局部最优解，且原理简单，计算速度快，结果稳定可靠。

发明内容

本发明的目的是使用一种基于模型集群分析的激光诱导击穿光谱变量选择方法来挑选真正能够代表物质成分特征并对分类效果有促进作用的有用变量，避免噪音、环境信息、冗余数据等无用变量对分类的影响，同时大大减少数据量，显著降低计算成本。

为实现上述目的本发明采用的技术方案是：

一种基于模型集群分析的激光诱导击穿光谱变量选择方法，包括以下步骤：

(1) 利用激光诱导击穿光谱系统对不同牌号的钢材样品分别在不同的测量位点进行光谱数据采集获得全光谱数据；

(2) 对全光谱数据矩阵进行蒙特卡洛采样

将全光谱波段划分为若干变量，每个变量为包含一定波长点数的波段，对于全光谱数据X_m×n，其中m为光谱个数，n为变量数，每次采样无放回地从X_m×n中抽取s个光谱和q个变量组成一个子数据矩阵R_s×q，从X_m×n余下的数据中抽取同样的q个变量组成T_（m-s）×q，重复进行N次采样可得到N个子数据矩阵；

(3) 子数据矩阵建模分类

将每个经蒙特卡洛采样得到的子数据矩阵R_s×q作为训练集，对应的T_（m-s）×q作为测试集，用一种改进的支持向量机分类方法——组合模型——对训练集建立子模型，然后用测试集预测分类并计算预测正确率；最后N次采样共得到N个子模型和N个预测正确率；

（4）挑选候选有用变量

对每一个变量依次进行统计分析检验其是否是有用变量，对于变量i，将N个子模型对应的预测正确率分为两组：包含变量i的子模型对应的预测正确率分为A组，不包含变量i的子模型对应的预测正确率分为B组；两组预测正确率的平均值之差Dmean_i=Mean_i,A-Mean_i,B, 如果Dmean_i>0，则将变量i作为候选有用变量，否则变量i被认为是无用变量；

（5）通过无参数检验确定真实有用变量

如果变量i是候选有用变量，利用Mann-Whitney U检验变量i的A组预测正确率和B组预测正确率进行分析并计算p值，如果p值小于预设阈值，则认为变量i是真正有用变量，否则判定变量i为无用变量。

上述步骤（2）中进行蒙特卡洛采样时需要确定三个参数，蒙特卡洛采样次数N, 每次采样抽取光谱个数s和每次采样抽取变量个数q，为了充分获取各种变量组合的情况，N越大越好，但考虑到计算成本，N通常为5000—10000，s值取光谱总数m的1/2—2/3，q值往往会对预测结果造成影响，因此设置一系列取值范围，通过对比选择最适合本组数据的q值，通常q值的取值范围为5-500。

上述步骤（3）中，所述的组合模型先将数据进行一对多建模模糊分类筛选出候选类别，然后在候选类别范围内进行一对一精细分类，通过循环迭代，逐步缩小候选类别，最终确定预测类别。

上述在步骤（5）中，利用Mann-Whitney U检验对A组预测正确率和B组预测正确率进行统计分析，其中p值是评估两组预测正确率是否有显著性差异的参数，当p值小于预设阈值，说明变量i的加入对预测能力有显著提高，是真正有用变量，否则认为变量i对模型预测没有影响或者有消极影响，是无用变量，其中阈值通常为0.01—0.05。

本发明的优点与积极效果：

  (1) 本发明基于模型集群分析思想，通过蒙特卡洛采样获得子模型，对子模型进行建模，然后统计分析模型预测正确率的分布得出结果，通过随机重采样计算得到稳定结果比一次性建模分析的结果更具有可靠性和普遍性；

(2)使用一种改进的支持向量机分类方法：组合模型。先将数据进行一对多建模模糊分类筛选出候选类别，避免无用类别的干扰。然后在候选类别范围内进行一对一精细分类，由于没有无用类别的干扰，而且一对一分类考虑到所有候选类别间的差异，通过循环迭代，逐步缩小候选类别，最终确定预测结果，因此模型预测能力被显著提高；

(3) 通过Mann-Whitney U检验统计分析各种变量组合对模型预测的影响，挑选出能够显著性提高模型预测能力的变量，这样得到的结果相比于传统寻优算法是具有统计学意义，且不会陷入局部最优解的情况。

附图说明

图1是本发明操作流程图；

图2是有用变量和无用变量预测正确率分布图；

图3是不同q值时本专利选择变量的预测正确率；

图4是实例中钢铁样品20#的LIBS全光谱图；

图5是只包含有用变量的20#LIBS光谱图。

具体实施方式

如图1所示，本发明基于模型集群分析的激光诱导击穿光谱变量选择方法包括以下步骤：

(1) 利用激光诱导击穿光谱系统对不同牌号的钢材样品分别在不同的测量位点进行光谱数据采集获得全光谱数据。

(2) 对全光谱数据矩阵进行蒙特卡洛采样

将全光谱波段划分为若干变量，每个变量为包含一定波长点数的波段，对于全光谱数据  X_m×n，其中m为光谱个数，n为变量数，每次采样无放回地从X_m×n中抽取s个光谱和q个变量组成一个子数据矩阵R_s×q，从X_m×n余下的数据中抽取同样的q个变量组成T_（m-s）×q，重复进行N次采样可得到N个子数据矩阵。此步骤需要确定三个参数：蒙特卡洛采样次数N, 每次采样抽取光谱个数s和每次采样抽取变量个数q。为了充分获取各种变量组合的情况，N越大越好，但考虑到计算成本，N通常为5000—10000，s值取光谱总数m的1/2—2/3，q值往往会对预测结果造成影响，因此设置一系列取值范围，通过对比选择最适合本组数据的q值，通常q值的取值范围为5-500。

(3) 子数据矩阵建模分类

将每个经蒙特卡洛采样得到的子数据矩阵R_s×q作为训练集，对应的T_（m-s）×q作为测试集，用一种改进的支持向量机分类方法——组合模型——对训练集建立子模型，用测试集数据预测分类并计算预测正确率, N次采样共得到N个子模型和N个预测正确率。

组合模型的建模预测具体步骤如下：

A．一对多建模模糊分类：对于第a类，将训练集数据中属于第a类的数据设为正标签，其他所有类别的数据都设为负标签，并建立针对第a类的二元分类器，总共k类数据共需建立k个二元分类器。然后将测试集数据依次带入各分类器预测，综合k个分类器的预测值得出一对多模型的模糊分类预测结果，如果一对多模型判断该数据属于某一类别，则整个预测过程结束；如果预测结果出现多分类情况，则把数据可能属于的类别作为候选类别；

B．一对一建模精细分类：将所有候选类别两两组合，每两个候选类别建立一个二元分类器，对于f个候选类别，则共需建立f(f-1)/2个二元分类器。然后测试数据被所有二元分类器依次预测，每个二元分类器做出判断后按如下方法投票：对于a—b类二元分类器，如果该分类器判断测试数据为第a类，则第a类得票数加1，否则第b类得票数加1。所有分类器均按上述方法判断并投票后，统计各类别总得票数，以得票数最高的类别为最终预测结果。如果最高得票数的类别不止一种，则将最高得票数的所有类别作为新的候选类别，然后重复迭代步骤b中一对一建模分类方法直至最终确定唯一类别，即为最终预测类别；或者连续两次候选类别完全相同，此时判定该数据“无法分类”。

(4) 挑选候选有用变量

对每一个变量依次进行统计分析检验其是否是有用变量。对于变量i，将N个子模型对应的预测正确率分为两组：包含变量i的子模型对应的预测正确率分为A组，不包含变量i的子模型对应的预测正确率分为B组。两组预测正确率的平均值之差Dmean_i=Mean_i,A-Mean_i,B, 如果Dmean_i>0，说明变量i加入模型可能对SVM的预测能力有所提升，因此将其作为候选有用变量。如果Dmean_i≦0，则可认为加入变量i对模型的预测能力没有影响或者会降低预测能力，此时认为变量i是无用变量。

(5) 通过无参数检验确定真实有用变量

如果变量i是候选有用变量，说明加入变量i的模型预测正确率平均值高于没有变量i的模型预测正确率平均值。为了验证这种提高是否具有显著性，利用Mann-Whitney U检验对变量i的A组预测正确率和B组预测正确率进行分析并计算p值，如果p值小于预设阈值，则认为变量i对预测能力有显著提高，是真实有用变量，否则判定变量i为无用变量。

实施例1

以下以对九种不同牌号的圆钢样品的建模分类过程中的变量选择为例，结合附图和实例来进一步说明本发明的操作流程，但本发明不限于此例。

本实例使用的的LIBS系统主要由调Q脉冲Nd:YAG激光器，中阶梯光谱仪（ARYELLE-UV-VIS, LTB150, German)，可移动样品台和计算机等组成。激光能量为61mJ, 基频光波长1064nm，脉宽10为ns，重复频率为10Hz，光谱范围为220nm-800nm，共有29888个波长点。

选择九种不同牌号的圆钢样品：20#(Φ20×900mm), 20Cr(Φ20×900mm), 20CrMnTi(Φ30×900mm),20CrMo(Φ20×900mm),20CrNiMo(Φ20×900mm),35#(Φ20×900mm),35CrMo(Φ20×900mm),40Cr(Φ20×900mm),42CrMo(Φ25×900mm)(西宁特殊钢股份有限公司)。每一类样品在不同位置截取三个6mm高的钢柱，经打磨处理后放置在样品台上，利用激光诱导击穿光谱系统对样品进行测量。

 在每个钢柱的每个截面上随机挑选五十个测量点，在每个测量点经20次连续激光脉冲打击后得到一个测量光谱，每五个测量光谱经平均得到一个分析光谱，最终九类钢材共获得540个分析光谱（一类钢材三个样品，一个样品两个截面，一个截面10个分析光谱）。

考虑到钢材样品具有整体成分不均一，局部成分均匀的特点，如果将所有光谱随机分为训练集和测试集则容易出现过拟合导致预测准确率虚高。因此，每类样品中随机选择4个截面的光谱数据作为训练集，其余为测试集。

为了避免过拟合，变量选择过程中只使用训练集，待确定真实有用变量后再用测试集数据检验效果。

另外，如果以每一个波长点为一个变量，则计算量过大且没有意义，因此将29888个波长点中每30个波长点做一个变量，最后不足30个波长点的用0补足。

这样整个光谱共有997个变量，每个变量为包含30个波长点的光谱波段。

按图1所示流程进行变量提取，对于全光谱数据X_360×997，每次采样无放回地从X_360×997中抽取240个光谱和q个变量组成一个子数据矩阵R_240×q，其中q值分别为5、10、30、50、100、150、200。

从X_360×997余下的数据中抽取同样的q个变量组成T_120×q。 重复进行5000次采样可得到5000个子数据矩阵。

将每个经蒙特卡洛采样得到的子数据矩阵R_240×q作为训练集，对应的T_120×q作为测试集，用一种改进的支持向量机分类方法——组合模型——对训练集建立子模型。

建模过程中选择多项式核函数，使用训练集数据通过网格法对多项式参数d(1—10)和惩罚因子C(10^-5—10⁵）进行寻优，最终确定最优参数d=1,C=1，然后用测试集数据预测分类并计算预测正确率。5000次采样共得到5000个子模型和5000个预测正确率。 对每一个变量依次进行统计分析检验其是否是有用变量。

对于变量i，将5000个子模型对应的预测正确率分为两组：包含变量i的子模型对应的预测正确率分为A组，不包含变量i的子模型对应的预测正确率分为B组。

两组预测正确率的平均值之差Dmean_i=Mean_i,A-Mean_i,B, 如果Dmean_i>0，则将变量i作为候选有用变量，否则变量i被认为是无用变量。

如果变量i是候选有用变量，利用Mann-Whitney U检验对变量i的A组预测正确率和B组预测正确率进行分析并计算p值。

如果p值小于预设阈值0.05，则认为变量i对预测能力有显著提高，是真实有用变量，否则仍然判定变量i为无用变量。

挑选出真实有用变量后，从原数据矩阵中抽取真实有用变量组成新的矩阵，然后建模对测试集数据进行预测并计算预测正确率。

图2A是典型的有用变量（第303号变量，q=50）对预测正确率分布的影响，其中条纹是A组预测正确率分布图，黑色是B组预测正确率分布图。 从图中可以看出当第303号变量加入模型后，预测正确率的分布明显右移，说明该变量使SVM的预测能力有所提高，而且这种提高是显著的（p值=2.28×10^-21）。

图2B是典型的无用变量（第544号变量，q=50）对预测正确率分布的影响，其中条纹是A组预测正确率分布图，黑色是B组预测正确率分布图。可以看到第544号变量加入模型后，预测正确率并没有明显变化，因此该变量对SVM分类没有帮助，是无用变量。

由于变量选择的过程是基于蒙特卡洛采样，采样过程具有随机性，因此每次运行本专利进行变量选择的结果都会存在一定的差异。

为了检验这种差异对结果稳定性造成的影响，选择一系列q值，在每一个q值下运行本专利20次，用每次挑选出的有用变量组成新的数据矩阵建模预测，然后计算20次预测正确率的平均值和标准偏差。

通过图3可看到，随着q值的增加，预测正确率显著升高，当q=50时，达到最大值98.03%,然后开始下降，而标准偏差的变化趋势和预测正确率的变化趋势相反，q=50时，标准偏差最小。

这表明通过本专利挑选有用变量对提高SVM预测能力有显著帮助，而选择合理的变量采样数q能够降低蒙特卡洛采样随机性造成的偏差，提高本专利的稳定性。

结合表1和图3可确定对于本数据最优q值为50。

对q=50时运行本专利20次挑选的有用信息进行统计，将被选中作为真实有用变量次数最多的前35个变量当做本专利变量提取的最终结果。

图4是20#的全光谱谱图，图5是只包含有用变量的光谱图。

可以看出原光谱中的大量干扰信息和噪音都已被剔除，只利用有用变量对各种牌号的钢材进行分类，可以显著提高分类准确率且计算时间大大减少。