基于中红外光谱的掺伪牛奶非靶向检测方法研究

刘伯扬,高安平*,杨 戬,高永亮,白 鹏,特日格乐,马利军,赵三军,李雪晶,张慧萍,康俊巍,李 慧,王 慧,杨 斯,李晨曦,刘 蓉

1. 内蒙古蒙牛乳业(集团)股份有限公司,内蒙古 呼和浩特 011500 2. 天津大学精密仪器与光电子工程学院,天津 300072

引 言

牛奶含有多种人体必需的营养成分,在居民日常膳食中占有重要地位[1]。随着经济发展和人民生活水平提高,牛奶及其衍生产品的消费量不断增长,消费者对乳制品质量安全的要求也逐渐提高。原料奶是整个乳制品产业链中的最上游,其质量直接影响乳制品的生产与消费安全。近年来,原料奶掺假事件时有发生,一些非法添加剂,如三聚氰胺、蛋白精、革皮水解物等,是造成重大食品安全事故的重要因素,严重危害消费者健康。现代乳制品生产环节较多,生产设备卫生安全不达标,或者生产环节出现泄漏,也容易造成二次污染及经济损失。为了保证牛奶及乳制品生产消费安全,我国制定了严密的质量标准体系,并且要求在乳制品生产及消费各个环节对各种非法添加物进行严格检测。

牛奶主要成分,如蛋白质等的检测常用凯氏定氮法、酶-比色法等化学方法,步骤繁琐,检测灵敏度低,时效性差。对于三聚氰胺、植物蛋白、糖类等掺假物质,则主要以高效液相色谱、气相色谱-质谱联用、光谱仪为代表的仪器分析方法为主[2-4]。仪器分析方法检测灵敏度较高,可以实现微量成分检测,但是需要根据不同掺假物种类设计预处理及分析流程,检测效率较低,无法满足实时在线检测需要[5-6]。现代乳制品工业要求在生产的各个环节对牛奶质量及可能存在的掺假物质进行快速检测。目前的在线检测方法及设备仍比较落后,难以实现标准化快速检测,影响生产效率及产品质量。

光谱检测技术具有快速、无损等优势,广泛应用于各种在线检测及过程分析[7-9]。红外光谱检测方法可以同时获取多种分子键及官能团信息,根据光谱吸收峰位置、强度对多种物质成分进行定性定量分析,尤其适用复杂样品分析。国内外研究小组针对牛奶质量及掺假物快速检测需要,研究了多种基于红外光谱分析技术的判别方法与定量分析模型。Brand[10]等采用衰减全反射傅里叶变换红外光谱技结合偏最小二乘回归模型法,建立牛奶中掺假物的回归模型,对尿素预测相关系数达到0.96以上。Zaalberg[11]等在牛奶近红外光谱测量基础上,结合主成分分析与回归判别方法,建立不同理化指标及新鲜度判别模型,取得了较好的准确度。Mohamed[12]等采用近红外光谱技术检测掺有植物水解蛋白的牛奶样品,采用三氯乙酸对掺假样品进行预处理后,检测结果更准确。随着光谱检测技术发展,机器学习及深度学习方法也为光谱数据处理及模型研究奠定了基础,扩展了光谱方法应用领域[13-15]。

目前多数掺假牛奶鉴别模型需要根据掺假物的光谱特征,建立用于不同类型掺假物的靶向鉴别模型,以提高识别精度,而实际生产中掺假物质种类多、成分复杂,建立靶向识别模型的复杂度和成本较高。为了满足生产在线快速检测需要,研究非靶向异常牛奶样品鉴别方法成为提高在线检测效率的关键。由于近红外波段主要是O—H、N—H、C—H等含氢基团的倍频、合频信息,吸收特征往往不明显,并且信号较弱,不利于多种掺假物同时鉴别。中红外光谱表征了分子基频吸收,天然牛奶与掺假物质在中红外波段具有较为明显的吸收特征,丰富的谱峰信息为非靶向判别模型提供了基础。然而,由于中红外光谱覆盖范围大,需要鉴别的掺假物种类较多,非靶向鉴别模型面临的主要问题是如何筛选有效波段与建模变量组合,结合不同分类算法,提高鉴别模型准确度及泛化能力[16]。

针对牛奶在线快速检测需要,研究了基于中红外光谱的多种类掺假异常牛奶样品非靶向鉴别方法。选择蒙牛公司奶质稳定的六个奶源地收集天然原样品,配制了含有多种掺假物的异常牛奶样品。采集样品中红外吸收光谱,利用预处理方法提高光谱信噪比及一致性。采用蒙特卡洛无信息变量消除(Monte Carlo uninformative variable elimination,MC-UVE)、无变量信息消除-连续投影(uninformative variable elimination-successive projections algorithm,UVE-SPA)与竞争自适应重加权采样(competitive adapative reweighted sampling,CARS)三种方法筛选原始光谱中的特征波长变量,分别建立偏最小二乘判别(partial least squares discriminant analysis,PLS-DA)及支持向量机(support vector machine,SVM)的非靶向识别模型,对比分析不同特征变量组合对于掺假异常牛奶样品的识别效果,为乳制品在线生产及质量检测提供支持。

1 实验部分

1.1 样品光谱测量

采用的牛奶样本来源于蒙牛公司六个奶质稳定的奶源地,共收集1 560例优质生鲜奶,其中470例作为正常牛奶样本集。根据原料奶收储中常见的掺假物,配制了含有植脂末、淀粉、蔗糖、尿素、三聚氰胺、水解蛋白粉、亚硝酸盐、柠檬酸、碳酸氢钠等多种类型的异常牛奶样品集,共计1 090例。每种样品均经过充分混匀,并放置于冰箱冷藏保存,室温环境下完成光谱测量。采用MilkoScan FT1傅里叶变换光谱仪(丹麦Foss公司)测量样品光谱,测量范围4 000～926 cm-1,光谱分辨率4 cm-1,每个样品重复扫描16次取平均光谱,以提高光谱信噪比。

牛奶样品属于高散射介质,透射光谱容易受到脂肪颗粒等多次散射影响,测量过程中的样品变化也会在一定程度上引起光谱基线漂移。在线检测过程中,受到仪器振动、噪声等因素影响,原始光谱信噪比也会下降。针对牛奶样品特点及测量干扰来源,选择平滑滤波、多元散射校正及归一化等预处理方法,提高光谱信噪比及一致性。

1.2 特征波长变量筛选方法

非靶向鉴别模型需要准确识别含有不同种类掺假物的异常牛奶样本,因此测量光谱范围较宽,以保证能够覆盖不同物质的吸收特征。测量得到的样品中红外吸收光谱共包含1 060个波长变量,不同波长信息对于鉴别模型的贡献不同,部分无信息变量还有可能影响分类模型性能。为了提高鉴别模型准确度以及算法效率,首先需要采用变量筛选方法剔除不相关的波长,得到特征波长变量组合作为判别模型输入。本文研究中,根据牛奶样品光谱特征、掺假物种类,并综合信息有效性等指标,选择MC-UVE、UVE-SPA与CARS三种方法对原始光谱中进行变量筛选,得到特征波长变量组合[16]。

MC-UVE 方法根据偏最小二乘回归系数的稳定性来消除无信息变量。首先在样品光谱矩阵中加入相同变量数目的随机噪声,然后通过交叉验证法建立PLS模型,根据得到的回归系数向量的平均值除以标准差得到判定阈值,分别比较波长变量与随机噪声最大阈值的范围,将小于随机噪声阈值的变量视为无信息变量。MC-UVE 方法在交互验证过程中,采用蒙特卡洛采样方法,从训练集中随机选择一定数量的样本作为构建PLS子模型的训练子集,可以在一定程度上提高选择最佳模型的概率。

采用MC-UVE筛选得到的数量往往较多,在一定程度上影响分类算法效率。研究中采用连续投影方法(successive projections algorithm,SPA)进一步筛选,获得最优变量组合。UVE-SPA 方法中,SPA用于选择具有最小冗余信息的变量,从一个变量开始,在每次迭代中合并一个具有最大投影向量的新变量,根据线性判别模型验证集中的最佳准确度筛选得到有效变量。UVE-SPA可以有效减少变量的数量,选出均方根误差最小的变量组合。

CARS方法则基于“适者生存”准则,采用迭代统计方法筛选有效信息变量。算法过程可以分为随机采样,波长选择,自适应重加权采样等过程,通过循环迭代方式筛选出最佳波长变量组合。随机采样过程一般采用蒙特卡洛采样方法,划分数据集进行建模分析,根据原始数据和相应的响应变量计算有效信息变量,将得到的变量按权重进行降序排序,应用偏最小二乘回归来评估所选变量的性能。每个变量的权重基于自适应加权采样计算公式如式(1)

(1)

式(1)中,bi是每个变量的回归系数的绝对值。选择具有最佳精度的变量子集作为有效信息变量集合。

1.3 判别算法

本研究的非靶向鉴别模型需要准确识别添加不同类掺假物的异常牛奶样品,因此对鉴别模型的普适性、无偏性与泛化能力要求更高。根据这一实际要求,分别研究了基于支持向量机(SVM)与偏最小二乘判别(PLS-DA)的鉴别模型,并对比不同波长变量组合的鉴别结果。

偏最小二乘判别在主成分分析基础上,结合了数据降维与回归模型的优势,实现样品准确、高效分类。首先对样本类别进行编码,然后使用偏最小二乘法对样本类别编码数值进行回归,划定阈值用于预测未知样品类别。其优势在于在选择判别主成分时不仅考虑了特征矩阵的信息,同时也兼顾了特征向量方向信息,具有较好的准确度和泛化能力。

支持向量机是一种常用的经典机器学习模型,其核心思想是寻找一个最优分割平面,使得两类样本与分割平面的间隔最大,实现分类与判别。算法核心是对于给定样本数目为m的训练样本集D={(x1,y1),(x2,y2),…,(xm,ym)},y={+1,-1},找到一个超平面(ω,b)使得间隔最大。该目标可表示为式(2)

s.t.yi(ωTxi+b)≥1,i=1,2,…,m

(2)

根据Karush-Kuhn-Tucker条件,只有满足yif(xi)-1=0的样品点才对会对支持向量机分类模型有影响,定义这样的点为支持向量。为了满足高维数据处理需要,支持向量机分类模型可以通过构建核函数的方式,将样品向量从原始空间映射到一个更高维的特征空间,并使其在特征空间内线性可分。利用不同的核函数,可以处理各种非线性分类与回归任务,避免直接计算高维特征空间内及提高算法效率。由于支持向量机最终决策函数只由少数的支持向量确定,避免了“维数灾难”,因此模型具有较好的泛化能力。

1.4 模型训练与评价

模型训练采用蒙特卡洛内部交叉验证方法,按照4∶1的比例随机划分为训练集和测试集,以保证训练后模型准确度及泛化能力。模型评价指标包括准确率(Accuracy),灵敏度(Sensitivity),特异性(Specificity),计算公式如式(3)—式(5)

(3)

(4)

(5)

式(3)—式(5)中,TP表示被正确分类为阳性的样本数目,FP表示被错误分类为阳性的样本数目,TN表示被正确分类为阴性的样本数目,FN表示被错误分类为阴性的样本数目。以上三种评价指标中,准确率是指所有样本中正确分类数量的占比,用于衡量分类器总体性能。但是在实际应用中,正常样品与异常样品数量可能不平衡,在某些情况下需要确定正常样本与异常样品各自分类精度。因此,分别采用灵敏度衡量分类器对正常样品的识别能力,而用特异性衡量分类器对异常样品的识别能力。

2 结果与讨论

2.1 牛奶样品中红外光谱比较

如图1所示,正常牛奶样本与添加不同掺假物的异常牛奶样品的中红外光谱在大部分范围内较为相似,但在某些波段表现出较大差异。其中,1 250、1 550和1 650 cm-1几个波长主要反映了蛋白及含有氨基类物质的吸收,而1 774、2 856和2 928 cm-1附近的吸收峰主要由脂肪引起。正常与异常牛奶光谱差别较大的区域主要分布在1 000～1 100、1 370～1 580和3 030～3130 cm-1三个波长范围。其中,1 000～1 100 cm-1波段主要反映氮杂环化合物的环振动峰,1 370～1 580 cm-1吸收峰主要与含磷成分及多肽类物质有关,3 030～3 130 cm-1主要反映了氨基及多肽等物质的吸收特征。

图1 正常牛奶样品及掺伪样品中混外光谱比较Fig.1 Mid-infrared spectrum of normal and abnormal milk

由于牛奶中87%的成分是水,水分子在中红外波段具有较强的吸收峰,尤其是在1 610～1 650及3 200～3 400 cm-1范围内,由于氢键缔合作用,表现为杂乱的强吸收峰。虽然在这一波长范围内不同样品的光谱差异较大,但是主要是水分子吸收造成的,与掺伪物质无关。因此,需要利用特征波长变量筛选方法剔除这一波段范围内的波长变量,以提高鉴别准确度。

2.2 波长变量筛选结果

鉴别模型的准确度与输入变量密切相关,采用变量筛选方法得到的建模波长变量组合不仅需要包含样品光谱信息,并且需要剔除无关变量以及可能对分类影响较大的波段。分别对比MC-UVE、UVE-SPA及CARS三种变量筛选方法得到波长变量组合的分布情况。

基于MC-UVE算法的特征波长变量选择结果如图2所示,图2(a)中分界线左侧为测量光谱变量矩阵,右侧为算法添加与测量光谱变量数相同的随机噪声矩阵,根据随机噪声最大阈值决定光谱矩阵中的无用变量,并将其剔除。由于MC-UVE方法根据随机噪声选取阈值,因此筛选得到的变量较多,如图2(b)所示,达到142个。

图2 MC-UVE波长变量筛选过程及结果(a):算法过程;(b):得到波长变量结果及分布Fig.2 Screening process and results of MC-UVE algorithm(a):Stability distribution of each variable;(b):Distributionof the selected variable

图3中,在MC-UVE方法基础上,进一步采用SPA方法进行二次筛选,最终得到变量数量为26个。从筛选得到的变量分布可以看出,大部分波长变量分布在不同物质吸收峰。变量数目减少,一方面可以避免分类模型训练出现过拟合,另一方面也可以提高模型分类识别效率。

图3 UVE-SPA波长变量筛选结果Fig.3 Results of UVE-SPA algorithm

CARS算法提取光谱特征变量,设置采样次数为100,采用5折交叉验证方式,重采样率选择为0.8。图4(a)中展示了筛选过程,随着迭代次数增加,算法选取的特征变量数据逐步减少。在前30次迭代中,特征变量数呈指数衰减趋势,在这之后的迭代过程中,特征变量数缓慢减少并逐渐趋于稳定。当采样数小于30时,RMSECV基本不变,大于48后,RMSECV缓慢增加,综合考虑建模效果,最终根据迭代次数选择最佳变量组合,得到变量数为56个。如图4(b)所示,筛选得到的变量也主要分布在蛋白质、脂肪以及可能存在的掺伪物质的特征吸收波长。

图4 CARS波长变量筛选过程及结果(a):算法过程;(b):得到波长变量结果及分布Fig.4 Screening process and results of CARS algorithm(a):Algorithmic processes;(b):Distribution of the selected variable

2.3 识别模型与结果分析

分别以MC-UVE,UVE-SPA和CARS筛选得到的特征变量组合为分类模型输入变量,建立基于PLS-LDA及SVM分类模型,结果如表1所示。对比不同分类模型效果可知,总体上SVM模型优于PLS-DA模型,特别是在变量较多的情况下,SVM分类模型的核函数映射方法可以避免数据维数过高造成的“灾难”。变量筛选方法可以进一步提高模型分类准确率,其中CARS算法筛选出的变量组织精度最高。MC-UVE与UVE-SPA筛选出的变量数量虽然相差较多,但是最后得到的分类准确度基本相当,并且UVE-SPA方法在预测集上表现更好,也说明采用SPA方法进一步筛选,减少特征变量数,可以在一定程度上减少模型过拟合风险。

表1 不同特征变量组合与分类模型结果比较Table 1 Prediction results by different combinations of characteristic variables with classification model

如图5(a,b)所示,对比不同特征变量组合与分类模型对正常和掺假牛奶误分类结果可知,基于CARS方法筛选得到的特征波长变量组合建立的SVM识别模型可以较为准确识别异常牛奶样品,整体误分类样品较少。而基于MC-UVE筛选特征变量组合与SVM分类方法对于正常牛奶的识别分类精度较高,但是对少部分异常牛奶样本产生误判。在实际应用中,由于异常样本误判可能会对后续生产及产品质量产生比较严重的影响,提高模型对于异常样本的判别准确度是保证生产安全的关键。

图5 异常牛奶鉴别结果混淆矩阵(a):CARS-SVM鉴别模型结果;(b):MC-UVE-SVM鉴别模型结果Fig.5 Confusion matrix detailing the multiclass discrimination results of abnormal milk(a):CARS-SVM model;(b):MC-UVE-SVM model

3 结 论

针对乳制品生产异常牛奶样品检测需要,研究了中红外光谱非靶向鉴别方法,结合变量筛选及机器学习方法,提高模型判别准确度。研究了基于MC-UVE、UVE-SPA与CARS方法的特征波长变量筛选方法,分别建立PLS-DA及SVM非靶向识别模型,鉴别多种掺假成分异常牛奶样品。实验结果表明,采用法筛选得到的特征波长变量组合建立的SVM分类模型效果整体优于PLS-DA模型。采用CARS筛选得到的波长变量组合结合SVM分类模型对训练集与测试集的鉴别准确率分别达到97.84%与94.55%,并且可以准确识别异常牛奶样品。研究表明,基于中红外光谱的非靶向鉴别模型可以实现多种掺假物的异常牛奶样品快速准确识别,为牛奶及乳制品质量与生产过程检测提供了支持。