基于主成分分析和梯度提升决策树的油菜蜜真伪鉴别方法与流程

1.本发明涉及蜂蜜真伪鉴别技术，尤其涉及一种基于主成分分析和梯度提升决策树的油菜蜜真伪鉴别方法


背景技术：

2.蜂蜜是蜜蜂从开花植物的花中采得花蜜在蜂巢中经过充分酿造而成的天然甜物质，其气味清香浓郁，味道纯真甜美。蜂蜜是一种糖类为主的天然食品，葡萄糖和果糖为其主要成分，可以被人体直接吸收，而不需要酶的分解，也是最常用的滋补品之一，深受消费者喜爱。
3.但是，蜂蜜的成分复杂，且组分含量与蜜源植物种类，花期，气候以及存贮加工工艺等相关，目前国内对蜂蜜的真伪判别，主要采用食品安全国家标准gb/t18932.1-2002和gb/t 18932.2-2002，使用靶向检测技术，如质谱或色谱技术，对蜂蜜的掺入物的特征组分进行判定。然而，造假者可以对标准中的检测项目，有组织地进行规避，故目前的国家标准已不能满足现阶段蜂蜜掺假或混配的检测，急需开发和建立基于机器学习的蜂蜜真伪鉴别的有效方法。
4.核磁共振氢谱检测技术作为非靶向检测技术，对蜂蜜样品有机组分中的氢元素进行检测，采集蜂蜜样品有机组分中氢元素的整体信号，可以反映出各有机组分的含量。外源物质掺入的蜂蜜样品，其核磁共振氢信号谱图的整体形貌会发生变化，比如糖浆的掺入会引起糖区信号的相对增加，与真实蜂蜜的谱图产生差异，借助人工智能机器学习方法去学习真蜂蜜和掺假蜂蜜的谱图特征，进而实现真假蜂蜜的判别，可以有效地检测蜂蜜中是否掺入糖浆。


技术实现要素：

5.本发明提供一种基于主成分分析和梯度提升决策树的油菜蜜真伪鉴别方法(下文中提到的蜂蜜和样品均为油菜蜜)。该方法主要包括：采集和制备真假蜂蜜样品并生成蜂蜜数据、对蜂蜜数据记录进行真假标注得到油菜蜜数据集、使用主成分分析法(pca)对油菜蜜数据集降维、使用梯度提升决策树算法(gbdt)构建油菜蜜真伪鉴别模型并评估模型的有效性；之后就可以利用该模型对待测蜂蜜进行真伪鉴别。
6.针对上述目的，本发明技术方案是：
7.一种基于主成分分析和梯度提升决策树的油菜蜜真伪鉴别方法，该方法包含下列步骤：
8.a：采集和制备真假蜂蜜样品并生成油菜蜜数据。
9.本方法中蜂蜜数据获取主要通过核磁共振技术，对蜂蜜样品进行noesy 1d检测，得到检测样品的核磁共振氢信号谱图，进而对谱图进行数字化转换，得到反映样品谱图特征的油菜蜜数据。具体方法如下：
10.a1:采集和制备真假蜂蜜样品。直接从蜂场，按照事先制定的程序，采集真实蜂蜜
样品若干；制备掺假蜂蜜样品若干，掺假蜂蜜样品主要是通过在真实蜂蜜样品中添加一定比例的不同来源的糖浆而获得，比例一般在10％～90％之间。
11.a2:采用核磁共振波谱仪对真假蜂蜜样品进行检测，得到蜂蜜样品的核磁共振氢信号谱图。谱图使用核磁波谱仪的采集和处理软件vnmrj 4.2进行相位调整、设定内标物的化学位移原点(内标物tsp)得到。谱图的形状反映了被测样品所含组分及其含量，是鉴别蜂蜜真伪的重要依据。
12.a3：对蜂蜜样品的核磁共振氢信号谱图进行数字化和归一化。使用核磁数据分析处理辅助软件mestronova 11.0对谱图进行等间隔积分并提取核磁数据，选取化学位移在0.53～9.53ppm(不包含4.78～4.98ppm的水峰信号)区间的信号，并对谱图以0.05ppm为间隔进行等间隔分段积分，采用峰面积归一化得到样品的分段积分。分段积分的方法如公式(1)-(3)：
[0013][0014]
p＝{p0,p1,p2,
···
pi,
···
p
n-2
,p
n-1
}
ꢀꢀꢀ
(2)
[0015][0016]
其中，a为化学位移在0.53～9.53ppm区域内总积分面积；ai为该积分区域内第i个分段积分面积；n为分段的个数，本发明具体实施时取值为176；pi为归一化后第i个分段积分值。
[0017]
由此得到代表该样品的一条数据记录，该记录具有176个属性值或特征值，这些特征值的分布反映了该样品核磁共振氢信号谱图的轮廓，即蜂蜜各组分含量分布。对每一个油菜蜜样品进行上述步骤，得到一系列具有176个特征的真假油菜蜜数据记录，即油菜蜜数据。
[0018]
b：对油菜蜜数据记录进行真假标注得到油菜蜜数据集。
[0019]
对步骤a中得到的油菜蜜数据中的每一条记录添加真假标签，真实蜂蜜标注“0”，掺假蜂蜜标注“1”。最后为了在实际测试中更加方便地使用数据，给每一条记录按照一定规则再添加一个样本名称，故此时得到具有176 2＝178维的油菜蜜数据集，记为d。
[0020]
c：使用主成分分析法(pca)对油菜蜜数据集d降维。
[0021]
使用pca对油菜蜜数据集d进行降维(不包括真假标注维和样本名称维)，降维维度按照一定的方差比例确定。首先对油菜蜜数据集d进行规范化，使得该数据中每一特征的取值平均值为0，方差为1；其次将规范化后的数据通过正交变换，得到由多个线性无关的新特征表示的数据；再将新特征按照方差大小进行排序，并将排序后的特征依次称为第一主成分、第二主成分等；计算每个特征的贡献率，当前k个特征的贡献率之和大于等于给定阈值时，前k个特征即为最终的特征集合，得到k维的油菜蜜数据集d
′
。
[0022]
采用主成分分析方法对油菜蜜数据集d进行计算的过程：假设油菜蜜蜜数据集d中有m条数据，n个特征，组成n
×
m蜂蜜数据矩阵x
*
,如式(4)所示:
[0023]
[0024]
其中xj(j＝1,2,
…
,m)为第j条蜂蜜数据。
[0025]
(1)对x
*
的每一行进行零均值化(每一行即各蜂蜜样本在该特征下的取值),即减去该行(该特征)的均值，得到蜂蜜零均值化矩阵x，如式(5)所示。
[0026][0027]
其中μi为第i个特征对应的平均值，计算过程如式(6)：
[0028][0029]
(2)计算蜂蜜零均值化矩阵x的协方差矩阵,得到n
×
n的协方差矩阵c，计算过程如式(7)所示：
[0030][0031]
(3)根据式(8)得到协方差矩阵c的特征值和特征向量：
[0032]
cβi＝λiβi，i＝1，2，...，n
ꢀꢀꢀ
(8)
[0033]
其中，βi为第i个特征向量，λi为第i个特征值且满足：λ1≥λ2≥...≥λn。
[0034]
(4)特征向量单位化，计算过程如式(9)：
[0035][0036]
其中pi为第i个特征向量的单位向量，此时，满足：
[0037][0038]
(5)计算特征值λi的贡献率γi，计算方式如式(10)所示。假定阈值为γ(为了保证保留足够的方差信息，一般γ≧90％)，当前k个特征值的贡献率之和大于等于给定阈值γ时，取前k个特征值对应的单位特征向量得到k
×
n的单位特征向量矩阵q，如式(11)，其中前面所述的k个主成分即为取到的k个特征值。
[0039][0040][0041]
(6)将单位特征向量矩阵q与蜂蜜零均值化矩阵x相乘得到降维后k
×
m的矩阵y，即k维的油菜蜜数据集d
′
，如式(12)。
[0042][0043]
d：使用梯度提升决策树算法(gbdt)构建油菜蜜真伪鉴别模型。
[0044]
对经pca降维后得到的油菜蜜数据集d
′
按一定比例划分为训练集和测试集，进行gbdt模型的训练和参数优化，获得训练好的最优的蜂蜜真伪鉴别模型(记为pca-gbdt)。
[0045]
gbdt模型具体通过基学习器的线性组合，不断减小训练过程产生的残差来达到将数据分类的目的。首先初始化一个弱学习器，计算此时训练集中每个样本损失函数的负梯度值；其次利用负梯度值拟合一棵决策树，并求出最佳的负梯度拟合值；然后使用最佳的负梯度拟合值来更新弱学习器，得到第一个强学习器。将前一轮得到的强学习器作为本轮初始弱学习器，按上述步骤不断迭代更新，当迭代次数达到给定的值时，输出的强学习器为最终的分类器。
[0046]
gbdt模型训练的具体计算过程如下：
[0047]
假设油菜蜜训练集样本s＝(x1,y1),(x2,y2),
…
,(xk,yk),
…
,(xm,ym)，其中xk是第k个样本，yk为第k个样本的真伪值(即0或者1)，最大迭代次数为t，损失函数l，输出是强学习器f(x)。同时，设置训练时的学习率(learning_rate)、训练样本的子采样比例(subsample)、决策树的最大叶子节点数(num_leaves)等参数，上述参数及最大迭代次数t均由参数优化最终确定。
[0048]
(1)初始化弱学习器f0(x)，表示为式(13)。其中，由于本鉴别方法为二元分类，故损失函数使用对数似然损失函数l(y,f(x))＝log(1 e-yf(x)
)，其中y∈{-1,1}，y为真实值，f(x)为预测值，初始时f(x)的值c设置为y(所有油菜蜜训练样本的标签值)的均值。
[0049][0050]
(2)对迭代次数t(t＝1,2,3,
…
,t)有：
[0051]
(2a)：对油菜蜜样本xi(i＝1,2,3,
…
,m),根据式(14)计算负梯度误差r
ti
：
[0052][0053]
(2b):将(2a)中得到的r
ti
作为样本当前的真实值，拟合得到第t棵决策树f
t
(x)，其对应的叶子节点区域为r
tj
，j＝1,2,3,
…
,j，每个叶子节点区域对应的是一个或多个油菜蜜样本。其中j为决策树t的叶子节点数量。
[0054]
(2c)：对每个叶子节点j，计算最佳拟合值的近似值c
tj
，如式(15)：
[0055][0056]
(2d)：更新强学习器，如式(16)：
[0057][0058]
(3)得到最终的强学习器f(x)，如式(17)。对于新的油菜蜜样本，将其数据代入式
(17)，即可得到分类结果:
[0059][0060]
参数优化使用scikit-learn的gridsearchcv方法实现，以对数损失函数(logloss)为模型评判标准获得最优参数。本方法中主要优化的参数包括：
[0061]
(1)决策树的个数或最大的迭代次数(n_estimators)，即式(17)中的t；
[0062]
(2)训练每棵树时数据的采样率，即子采样率(subsample)，以有效减小油菜蜜样本在拟合决策树时的方差和防止过拟合；
[0063]
(3)决策树的最大叶子数(num_leaves)，即步骤(2b)中所述j；
[0064]
(4)学习率(learning_rate),控制每棵决策树(弱分类器)对整体的影响。
[0065]
e:对待鉴别的油菜蜜样品进行noesy 1d检测，得到待测样品的核磁共振氢信号谱图，并对谱图进行数字化转换，得到反映样品谱图特征的油菜蜜数据，再利用训练好的油菜蜜蜜真伪鉴别模型进行真伪鉴别。
[0066]
与现有技术相比，本发明的有益效果是：
[0067]
本发明提出了一种基于主成分分析和梯度提升决策树(gbdt)的油菜蜜真伪鉴别方法。首先采集蜂蜜并制作真实蜂蜜和掺糖浆的假蜂蜜样品；然后使用核磁共振技术获取蜂蜜样品的核磁共振氢信号谱图，并将谱图数字化为176维的特征向量；进而使用pca方法保留一定的方差比例进行降维，以降低gbdt模型训练的复杂度；最后用已标注真假的蜂蜜数据对gbdt模型进行训练，得到结果gbdt模型，用于对新检测的油菜蜜进行真伪判别。与现有方法相比，本发明方法可以有效精准的鉴别油菜蜜的真伪，避免了人工查看谱图进行真伪鉴别的误差。综合使用核磁共振技术、主成分分析和梯度提升决策树gbdt模型，有效提高了油菜蜜真伪鉴别的准确率、均方根误差和auc值，降低了数据特征维数、模型训练时间和模型复杂度，是一种鉴别油菜蜜真伪的有效方法。
附图说明
[0068]
图1是本发明方法的整体流程框图
[0069]
图2(a)为本发明实施例中noesy 1d(0.5-9.5ppm)检测结果谱图。图2(b)为区间在0.5-3ppm，纵坐标放大500倍的放大图；图2(c)为区间在6-9.5ppm，纵坐标放大1000倍的放大图。
[0070]
图3是本发明实施例中对油菜蜜数据集d进行主成分分析之后特征的累计贡献率变化趋势图
具体实施方式
[0071]
下面结合附图，通过实施例进一步描述本发明。但不以任何方式限制本发明的范围。
[0072]
本发明提出了一种基于主成分分析和梯度提升决策树的油菜蜜真伪鉴别方法。首先采集和制备真假蜂蜜样品并生成蜂蜜数据集，其次对蜂蜜数据记录进行真假标注得到油菜蜜数据集，然后使用主成分分析法(pca)对油菜蜜数据集降维，最后使用梯度提升决策树
(gbdt)构建pca-gbdt蜂蜜真伪鉴别模型，基于模型对待测蜂蜜样本进行鉴别。本方法主要针对采集的油菜蜜样品进行测试，此方法综合了核磁共振技术、主成分分析和梯度提升决策树，能够快速高效方便的鉴别油菜蜜样品的真伪性。方法流程如图1所示。
[0073]
本发明中涉及到的仪器设备及试剂如下：
[0074]
(1)仪器设备：安捷伦公司agilent dd2 600mhz核磁共振波谱仪配有5mm二合一宽带onenmr探头，7510自动进样器；默克化工技术有限公司merck产milli-q integral超纯水系统；梅特勒-托利多mettler toledo的xs104天平和fiveeasy plus ph计；北京博宇宝卫实验设备公司thz-82a型恒温水浴振荡器；日本as one trio tm-1n涡旋振荡器；德国艾本德eppendorf公司200μl和1000μl可调节移液器；美国诺雷尔norell公司5mm核磁管。
[0075]
(2)试剂：分析纯磷酸二氢钾和氢氧化钾，购自西陇化工股份有限公司；99.9％氘代的重水和98％氘代的2,2,3,3-氘代三甲基硅烷丙酸钠(tsp)购自美国cambridge isotope laboratories；叠氮钠购自mpm biological technology co.ltd.；分析纯磷酸购自美国西格玛阿德里齐(sigma-aldrich)公司。
[0076]
本发明具体实施例的操作步骤为：
[0077]
a：采集和制备真假蜂蜜样品并生成油菜蜜数据
[0078]
a1:采集和制备真假蜂蜜样品。按照制定的程序从蜂场采集200个真实蜂蜜样品，对其中100个蜂蜜样品添加一定比例的不同来源的糖浆制成糖浆掺假蜂蜜样品，糖浆掺假比例在10％～90％不等。剩余100个蜂蜜样品作为真实蜂蜜样品；
[0079]
a2:对真假蜂蜜样品采用核磁共振波谱仪进行检测，得到蜂蜜样品的核磁共振氢信号谱图。谱图使用核磁波谱仪的采集和处理软件vnmrj 4.2进行相位调整和设定内标物tsp得到，其中变化点数设置为64k，线宽因子设置为0.3hz，设定内标物tsp的峰为0.00ppm进行化学位移定标；
[0080]
a3:对油菜蜜样品的核磁共振氢信号谱图进行数字化和归一化处理。使用核磁数据分析处理辅助软件mestronova 11.0对谱图进行等间隔积分并提取核磁数据，选取化学位移在0.53～9.53ppm(不包含4.78～4.98ppm的水峰信号)区间的信号，以0.05ppm为间隔进行等间隔分段积分，采用峰面积归一化得到样品的分段积分。提取的每个样品数据包含176维属性值或特征值，由此得到代表该油菜蜜样品的数据记录，该记录具有176个特征值，这些特征值的分布反映了该样品核磁共振氢信号谱图的轮廓，即蜂蜜各组分含量分布。
[0081]
对每一个真假油菜蜜样品进行a2，a3步骤，得到一系列具有176个特征的真假油菜蜜数据记录，即油菜蜜数据，从而实现样品检测数据的数字化。在本实施例中共获得100条真实油菜蜜数据记录和100条掺假油菜蜜数据记录组成200条油菜蜜数据。
[0082]
b：对油菜蜜数据记录进行真假标注得到油菜蜜数据集
[0083]
分别对真实蜂蜜和掺假蜂蜜数据记录添加真假标签，对每一条真实蜂蜜数据添加标签，标记为数值“0”，对每一条糖浆掺假蜂蜜数据添加标签，标记为数值“1”，标签特征名记为“result”，即对于真实蜂蜜设定其result＝0，掺假蜂蜜设定其result＝1。其次为了在实际测试中更加方便地使用数据，对每一条数据按照一定规则添加了样本名称(samplename)。此时得到200
×
178的油菜蜜数据集d，如表1所示，该数据是经过核磁共振和分段积分，并已标注真假标签和添加了样本名称的部分数据。其每一行为特定样品的特征(积分段)描述，每一列(除samplename和result)是不同的样品中特定化学位移区间的归一
化后的积分面积。
[0084]
表1油菜蜜数据集d(部分)
[0085]
ꢀꢀ
1234
…
173174175176samplenameresult0.53-0.580.58-0.630.63-0.680.68-0.73
…
9.33-9.389.38-9.439.43-9.489.48-9.53001r00.0033138690.0064941820.0067874190.005728609
…‑
0.000793290.001866126-0.001239190.000608956051r00.002013378-0.001988960.00156360.003106756
…‑
0.000115160.004121909-0.001126460.001833055052r0-0.00581904-0.00317804-0.00872074-0.00195171
…‑
0.00284376-0.00071118-0.00334178-0.00309494053r0-0.00387749-0.00630714-0.00363571-0.00155817
…‑
0.004329-0.00569407-0.00418908-0.00458526054r00.0008951550.0006450080.0005963890.002487297
…‑
0.0022927-0.00469504-0.00328826-0.0012771001r_h04_1010.001661179-0.000329640.0030651260.004724401
…‑
0.00139588-0.00280657-0.00491048-0.00412981051r_h04_101-0.00937495-0.00694828-0.0082971-0.00897548
…‑
0.00245339-0.00383869-0.0033809-0.0031896052r_h04_101-0.0029709-0.00413827-0.00354084-0.00155589
…‑
0.0009939-0.00489534-0.00554901-0.00158525053r_h04_101-0.000572020.0016295970.0014920770.002708676
…‑
0.00468328-0.0024162-0.00509602-0.00472317054r_h04_1010.0019400590.0018739650.0030283170.002668439
…‑
0.00425949-0.00289033-0.00544729-0.00558976
……………………………
[0086]
c：使用主成分分析法pca对油菜蜜数据集d降维。
[0087]
使用pca对具有176维(不包括真假标注维和样本名称维)的油菜蜜数据集进行降维，降维维度按照一定的方差比例确定。
[0088]
表2前9个主成分累计贡献率
[0089][0090]
根据pca计算过程，首先对该矩阵的每一行进行零均值化得到新的矩阵x，然后对新矩阵求协方差矩阵c，求出协方差矩阵c的特征值和特征向量，并对特征向量进行单位化，对特征值进行降序排序，取前k个特征值对应的特征向量，得到降维后的油菜蜜数据。实际计算过程中得到k个主成分，且排序之后前面的主成分所包含的数据信息大于后面主成分所包含的信息。选取k值时根据各个主成分累计贡献的大小选取，也即保留的方差信息比例。在本实施例中保留了99％的方差信息，也即k个主成分的累计贡献率大于等于99％时的k值。对油菜蜜数据集d进行主成分分析之后，当k＝9时，累计贡献率达到99.236％(保留三位有效数字)，得到9维的油菜蜜数据集d
′
，前9个主成分的累计贡献率如表2，累计贡献率变化趋势如图3所示。
[0091]
d：构建油菜蜜真伪鉴别模型。
[0092]
使用梯度提升决策树算法(gbdt)构建油菜蜜真伪鉴别模型。
[0093]
对经pca降维后得到的油菜蜜数据集d
′
按一定比例划分为训练集和测试集，进行gbdt模型的训练，然后对gbdt进行参数调整和优化，获得训练好的最优的蜂蜜真伪鉴别模型(记为pca-gbdt)，其中参数调整和优化使用网格搜索法。
[0094]
gbdt模型具体通过基学习器的线性组合，不断减小训练过程产生的残差来达到将数据分类的目的。首先初始化一个弱学习器，计算此时训练集中每个样本损失函数的负梯度值；其次利用负梯度值拟合一棵决策树，并求出最佳的负梯度拟合值；然后使用最佳的负梯度拟合值来更新弱学习器，得到第一个强学习器。将前一轮得到的强学习器作为本轮初始弱学习器，按上述步骤不断迭代更新，当迭代次数达到给定的值时，输出的强学习器为最终的分类器。
[0095]
gbdt模型训练的具体计算过程如下：
[0096]
对于油菜蜜训练集样本s＝(x1,y1),(x2,y2),
…
,(xk,yk),
…
,(xm,ym)，其中xk是第k个样本，yk为第k个样本的真伪值(即0或者1)，最大迭代次数为t，损失函数l，输出是强学习器f(x)。
[0097]
(1)初始化弱学习器f0(x)，表示为式(13)。其中，由于本鉴别方法为二元分类，故损失函数使用对数似然损失函数l(y,f(x))＝log(1 e-yf(x)
)，其中y∈{-1,1}，y为真实值，f(x)为预测值，初始时f(x)的值c设置为y(所有油菜蜜训练样本的标签值)的均值。
[0098][0099]
(2)对迭代次数t＝1,2,3,
…
,t，有：
[0100]
(2a)：对油菜蜜样本xi(i＝1,2,3,
…
,m),根据式(14)计算负梯度误差rti：
[0101][0102]
(2b):将(2a)中得到的r
ti
作为样本当前的真实值，拟合得到第t棵决策树f
t
(x)，其对应的叶子节点区域为r
tj
，j＝1,2,3,
…
,j，每个叶子节点区域对应的是一个或多个油菜蜜样本。其中j为分类树t的叶子节点数量。
[0103]
(2c)：对每个叶子节点j，计算最佳拟合值的近似值c
tj
，如式(15)：
[0104][0105]
(2d)：更新强学习器，如式(16)：
[0106][0107]
(3)得到最终的强学习器f(x)，如式(17)。对于新的油菜蜜样本，将其数据代入式(17)，即可得到分类结果:
[0108][0109]
将数据集d
′
按9:1的比例划分为训练集和测试集，进行模型的训练和参数优化，使用gridsearchcv进行参数优化。确定分类树的个数(estimators)、样本采样比例(subsample)、分类树的最大叶子树(num_leaves)和学习率(learning_rate)。对上述参数的取值范围设置如表3所示，经过gridsearchcv优化后最佳参数取值如表4所示：
[0110]
表3 gbdt算法参数取值范围
[0111][0112]
表4网格优化最佳参数取值
[0113][0114]
在优化后的模型中使用本实施例中的油菜蜜数据集d
′
准确率达到95％。为了体现pca-gbdt方法的优势，使用支持向量机(svm)和决策树(dt)算法对同样的油菜蜜数据进行
预测，其测试准确率(accuracy)、均方根误差(rmse)和auc值对比如表5所示。从表中可以看到，在油菜蜜测试数据集d上本发明方法pca－gbdt的准确率和auc值均高于svm和dt算法，rmse均低于svm和dt，可以看到pca-gbdt方法优于svm和dt方法。
[0115]
表5测试集实验结果对比
[0116][0117]
对待鉴别的油菜蜜样品进行noesy 1d检测，得到待测样品的核磁共振氢信号谱图，并对谱图进行数字化转换，得到反映样品谱图特征的油菜蜜数据，并使用pca进行降维，再利用训练好的油菜蜜真伪鉴别模型可以进行有效的真伪鉴别。
[0118]
需要注意的是，公布实施例的目的在于帮助进一步理解本发明，但是本领域的技术人员可以理解：在不脱离本发明及所附权利要求的范围内，各种替换和修改都是可能的。因此，本发明不应局限于实施例所公开的内容，本发明要求保护的范围以权利要求书界定的范围为准。