張傲林,李紳,鄒穎芳,王繼芬*,張震,孫一健
1(中國人民公安大學 偵查學院,北京,100038)2(武漢體育學院 運動訓練學院,湖北 武漢,430079)
咖啡作為一種常見的飲品,在國內外具有廣闊的消費市場。近年來摻假咖啡在市面上流行,以星空咖啡為例,摻假咖啡的成分除咖啡因、蛋白質、礦物質等主要成分還包含西布曲明、苯二氮卓等非法藥物。由于其組合藥理學機制為食欲抑制、飽腹感以及可能誘導產(chǎn)熱[1-3],廣泛流售于減肥人群中。微商、電商、物流寄遞等是其常見的銷售流通渠道,為使服用者減肥效果明顯,不法分子通常在咖啡中大量摻入此類食欲抑制類精神藥物,長期食用此類咖啡嚴重危害消費者的身體健康。由于售賣數(shù)量龐大,包裝形式多樣,檢測手法繁冗,針對國內消費市場形式,有必要建立一種快速無損檢驗方法來開展寄遞渠道大量咖啡中摻假非法添加成分的工作。
對于涉及咖啡及其添加成分的檢測,現(xiàn)有研究中的檢驗方法主要是分子光譜技術和液相色譜等方法,邵金良等[4]選擇高效液相色譜-紫外雙波長方法(high performance liquid chromatography-double UV wavelengths,HPLC-DUVW)同時測定咖啡制品中的生物堿成分以及多酚類化合物。陳華舟等[5]分別采用紫外(ultraviolet,UV)和傅里葉紅外(Fourier-transform infrared spectroscopy,FTIR)光譜法測定咖啡因和有機物的含量,在FTIR的指紋區(qū)域的特征吸收峰反映出咖啡中的無機物和多糖物質。陳秀明等[6]使用近紅外光譜并建立Adulterant Screen算法對咖啡中摻假物大麥等進行了識別。近年來,模式識別方法在食品的檢測中逐漸成熟,王艷艷等[7]建立了一種基于小波變換的反向傳播神經(jīng)網(wǎng)絡,對不同品牌咖啡的紅外光譜數(shù)據(jù)實現(xiàn)了100%分類。色譜檢驗方法預處理工序復雜,檢測效率低,對樣本具有破壞性;傳統(tǒng)分子光譜的穿透性較差,對樣本的檢測受環(huán)境噪聲干擾較大,上述方法不適用于大批量樣本的快速、無損檢測,無法滿足打擊食品犯罪的實際需求。
太赫茲時域光譜(THz time-domain spectroscopy,THz-TDS)是一種新型光譜技術,由于0.1~10 THz電磁頻域區(qū)間具有量子能量低、穿透性能優(yōu)異且許多生物大分子的振動和旋轉頻率都處于THz波段,所以利用THz-TDS可以獲得豐富的生物及其材料信息,且不會使檢材由于量子作用產(chǎn)生爆炸[8-10]。此外,THz輻射能以很小的衰減穿透如碳板、布料、塑料等常見物質,因此THz-TDS非常適合用于存在外包裝的檢材無損檢測。楊少壯等[11]選擇不同肉類在0.6~1.4 THz波段的太赫茲光譜進行主成分分析-支持向量機(principal component analysis-support vector machine,PCA-SVM)建模。結果表明,吸收系數(shù)等光譜參數(shù)的建模分類準確率均能達100%,可以實現(xiàn)對不同肉類的無損檢測。王倩等[12]為實現(xiàn)不同品種的大米鑒別,利用標準差和區(qū)間偏最小二乘方法提取0.53~1.21 THz波段的太赫茲吸收光譜信息進行建模,對比多種模型發(fā)現(xiàn)決策樹模型分類準確率為95%,能夠為食品無損檢測提供一種新思路。目前在國內在使用THz-TDS檢測咖啡方面的研究報道較少[13-14],基于維護涉及咖啡類食品安全的需要,本實驗收集7種不同品牌的70份摻假有西布曲明成分的咖啡樣本,通過THz-TDS獲取光譜數(shù)據(jù),比較導數(shù)(derivative)、巴特沃斯濾波器(Butterworth filter)、特征融合(fusion feature)3種不同類型預處理的建模效果,建模方法選擇隨機森林(random forest,RF)、支持向量機(support vector machine,SVM)、貝葉斯判別分析(bayes discriminant analysis,BDA),以確定最優(yōu)的檢測模型,并進一步對產(chǎn)地進行識別。
實驗中70份減肥咖啡樣本來自于海關實際案件中繳獲,根據(jù)品牌可以劃分為星空咖啡、真業(yè)主代、微韻咖啡、我叫小幾、趙氏思思小樂。其中星空咖啡按來源地劃分可以分為星空聯(lián)創(chuàng)、星空專柜、星空總店,為便于實驗將7類樣本按C1、C2、C3、C4、C5、C6、C7依次編號,具體見表1。
表1 實驗樣本及信息Table 1 Experimental samples and information
THz-TDS系統(tǒng)的原理是基于相干探測技術的太赫茲產(chǎn)生與探測系統(tǒng),可以同時提取太赫茲脈沖的振幅信息、相位信息,通過對時間波形傅立葉變換的處理,能直接得到樣品的吸收系數(shù)。本實驗選用QT-TS2000型號太赫茲時域光譜系統(tǒng)(青島青源峰達太赫茲科技有限公司)對樣本進行快速分析,具體儀器參數(shù)見表2。
表2 儀器參數(shù)Table 2 Instrument parameter
1.3.1 隨機森林
RF是基于Bagging框架和決策樹參數(shù)的集成學習算法,其最大特點可以處理高維特征的輸入樣本,并且能夠評估重要特征[15-16]。集成的樹集合是通過替換訓練樣本的子集進行創(chuàng)建的,通過選擇合適樣本特征進行選擇劃分,可以進行交叉驗證。其主要參數(shù)如下:決策樹最大深度(max depth),內部節(jié)點再劃分所需最小樣本數(shù)(min samples split),葉子節(jié)點最少樣本數(shù)(min samples leaf),最大特征數(shù)(max features),袋外(out-of-bag, OOB)誤差。其中OOB被主要用于模型整體誤差評估。
在統(tǒng)計學上其參數(shù)優(yōu)化組合取決于決策樹最大深度(max depth) , 內部節(jié)點再劃分所需最小樣本數(shù)(min samples split), 葉子節(jié)點最少樣本數(shù)(min samples leaf), 最大特征數(shù)(max features)。由于RF模型可以對數(shù)據(jù)進行自動降維,故當數(shù)據(jù)維度消減到一定時通過對參數(shù)的調整不再能夠影響模型變化。
1.3.2 支持向量機
SVM算法是傳統(tǒng)監(jiān)督學習中性能較為優(yōu)異的分類算法[17-18],主要用于解決線性不可分的數(shù)據(jù)集的分類問題,被定義為特征空間上間隔最大的線性分類器模型。主要通過RBF核函數(shù)、Sigmoid核函數(shù)、線性核函數(shù)、Polynomial核函數(shù)等對向量進行投影劃分最優(yōu)超平面。本實驗選用多分類問題效果優(yōu)異的RBF核函數(shù),其表達式和超平面表達式如公式(1)、公式(2)所示:
(1)
ωx+b=0
(2)
1.3.3 貝葉斯判別分析
BDA在統(tǒng)計學上可以根據(jù)全集樣本的先驗概率得到概率密度函數(shù)[19],然后利用貝葉斯判別公式[公式(3)]計算不同樣本判定為不同子域的概率即后驗概率P(Ci|x),最后根據(jù)后驗概率的大小進行判別,得到整體分類準確率。
(3)
式中:P(x)為x樣本發(fā)生的概率。
在實驗中,首先對未經(jīng)預處理的原始數(shù)據(jù)進行建模研究,選擇RF和SVM模型對70份咖啡樣本的光譜數(shù)據(jù)進行分類研究,由圖1-a和圖1-b中結果可知,RF模型最優(yōu)參數(shù)分類準確率為79.8%,SVM模型最優(yōu)參數(shù)分類準確率為90.8%。
a-RF;b-SVM圖1 RF、SVM參數(shù)影響分類準確率Fig.1 RF and SVM parameters affect classification accuracy
RF模型作為集成學習算法,在模型中可以通過樹子集進行交叉驗證提高。通過輸入70組樣本光譜特征數(shù)據(jù)進行建模,結合交叉驗證能夠得出RF模型對此類咖啡的分類結果。如圖1-a所示,實驗中在max depth≤8,min samples leaf≤20的范圍內,通過調節(jié)參數(shù)可以發(fā)現(xiàn)減肥咖啡的分類準確率在27.5%~79.7%呈離散型變化。當控制最大節(jié)點數(shù)=10 000,max depth=8時,隨min samples leaf不斷減小,模型分類準確率逐漸上升,當min samples leaf=2時準確率不再變化達最高79.7%。
在監(jiān)督學習算法中SVM模型的“魯棒性”較為優(yōu)異,可以根據(jù)核函數(shù)的調整適應不同的分類問題。進而在不同的分類中,通過相應核函數(shù)構建最優(yōu)超平面可以達到最佳訓練效果。本實驗選擇適合用于多分類問題的RBF核SVM,回歸精度設置為0.1,通過調整規(guī)則化參數(shù)R和RBFγ對模型進行優(yōu)化,當R在5~30區(qū)間變化且RBFγ在0.1~0.5區(qū)間變化時準確率明顯提升,當R=30,RBFγ=0.5時達到最高分類準確率90.8%。
BDA模型建模首先計算函數(shù)的最大后驗概率,建立6組BDA判別函數(shù)對樣本進行逐一分類。通過選擇Sig值最明顯的3個函數(shù)BDA1、BDA2、BDA3建立三維方向上的特征,取BDA模型計算判別得分作為方向上的變量,建立樣本三維空間分布圖。圖2中各類品牌咖啡樣本在空間中的聚斂程度不一,部分品牌的樣本聚斂較為分散,容易與其他樣本混淆,總體分類效果較差。
圖2 BDA模型樣本空間分布圖Fig.2 Spatial distribution of BDA model samples
在上述3類模型建模效果不理想的基礎上,進一步選擇預處理方法進行比較,對3種模型進行優(yōu)化,并確定最優(yōu)預處理建模模型。
如圖3-a所示,70組樣本在0~2.5 THz波段的太赫茲時域光譜中的吸收曲線形,吸收曲線走向以及特征峰位置基本一致,在0~0.15 THz及1.8~2.5 THz波段內具有8處明顯特征峰,分別在0.11、0.08、1.8、1.98、1.96、2.22、2.33、2.37 THz八個位置。
a-原圖;b-一階導數(shù);c-二階導數(shù);d-低通巴特沃斯;e-高通巴特沃斯;f-帶通巴特沃斯圖3 光譜預處理效果圖Fig.3 Spectral pretreatment rendering
通過進一步觀察樣本在0~2.5 THz內的譜圖,發(fā)現(xiàn)較大區(qū)間范圍內的光譜信號存在嚴重的噪聲和冗余信息,在特征峰頻率區(qū)間內,不同樣本的特征峰會有較多的重合現(xiàn)象。同時,實驗環(huán)境不能完全將濕度、光線等干擾噪聲消除,需要進一步通過光譜預處理和特征提取的方法消除冗余成分,篩選出主要信息,以提高建模的準確率,達到對不同品牌中咖啡摻假的識別分類。圖3-b和圖3-c是分別通過一階導數(shù)(first derivative,FD)、二階導數(shù)(second derivative,SD)處理后的咖啡太赫茲光譜圖。通過導數(shù)對光譜進行預處理可以有效的消除實驗環(huán)境導致的常數(shù)基線漂移,且能夠將不完全重合的譜峰進行分離。對比圖3-b和圖3-c,隨導數(shù)階的增加,咖啡的太赫茲光譜特征峰分離效果趨于明顯,并能夠有效提高光譜的分辨率、信噪比以及模型的檢測效果。但增加導數(shù)的階,會放大光譜的高頻信號,圖3中二階導數(shù)的高頻信號明顯強于一階導數(shù),因此本研究只選擇一階導數(shù)和二階導數(shù)的方法,避免高階導數(shù)的處理造成信噪比降低。
在信號處理領域,巴特沃斯濾波器因通頻帶內的頻率響應曲線最大限度平坦,沒有頻率上的起伏,而在阻頻帶則逐漸趨近于0。根據(jù)常見的頻率可以分為低通、高通、帶通、帶阻四種濾波器。觀察圖3-d、圖3-e、圖3-f中通過低通、高通、帶通3種10階巴特沃斯濾波器對太赫茲光譜圖進行校正的結果,可以有效地篩選光譜的主要特征信息,并且分離重合的特征峰,是一種良好的光譜信號處理方法。
通過振幅和頻率之間的關系函數(shù)表示巴特沃斯濾波器如公式(4)所示:
(4)
式中:n為濾波器的階數(shù);ωp為截止頻率;ωc為通頻帶邊緣頻率。
2.2.1 導數(shù)預處理方法
由圖3可知,選擇常用的一階導數(shù)、二階導數(shù)對咖啡樣本的太赫茲光譜進行預處理,能夠有效的區(qū)分開在0~0.3 THz區(qū)間和2.2~2.5 THz區(qū)間的特征峰,且對于無特征意義的頻率范圍可以進行平滑處理,防止模型計算冗余。何欣龍等[20]在進行車禍現(xiàn)場中保險杠的識別研究中發(fā)現(xiàn),一階導數(shù)和二階導數(shù)處理的BDA模型優(yōu)于原始數(shù)據(jù)和三階導數(shù),這是由于三階導數(shù)會大幅降低光譜的信噪比和識別靈敏度?;谏鲜鲅芯?本實驗選擇FD-BDA、FD-SVM、FD-RF、SD-BDA、SD-SVM、SD-RF六種導數(shù)融合識別模型與原始數(shù)據(jù)建模進行對比,如圖4所示,未經(jīng)導數(shù)預處理的BDA模型準確率僅為88.6%,經(jīng)一階導數(shù)處理后,BDA模型分類準確率均明顯提高,其中FD-BDA模型分類準確率最高為98.6%。FD-RF模型分類準確率僅為63.8%,與原始數(shù)據(jù)建模結果相同,這可能是因為在導數(shù)計算的過程中,一階導數(shù)只能消除常數(shù)項的噪聲,RF模型本身存在特征篩選的作用,所以一階導數(shù)的處理效果對于RF模型訓練影響不大,但經(jīng)二階導數(shù)處理后可以進一步對光譜信息進行降噪,并且對較高頻率的光譜基線進行校正,因此RF模型的分類準確率明顯提高。SVM模型的變化幅度較小,且SD-SVM 圖4 導數(shù)預處理模型分類結果Fig.4 Classification results of derivative preprocessing model a-太赫茲時域光譜;b-紅外光譜圖5 太赫茲時域光譜和中紅外光譜的特征選擇圖Fig.5 Feature selection of THZ time-domain spectra and mid-infrared spectra 2.2.2 光譜信號校正方法 Butterworth濾波器是一種全極點濾波器,函數(shù)零點均在無窮遠點,其主要特點是最大平坦性和幅頻單調下降,可以有效降低頻率以外的干擾噪聲,優(yōu)化特征變量的建模效果。蔣強[21]通過Butterworth濾波器對太赫茲時域光譜信息進行處理,能夠有效降低管道厚度測量模型的誤差效果?;诖?本實驗選擇采樣頻率/FS為48 000,濾波器階為10的低通、高通、帶通3種Butterworth濾波器進行預處理建模比較。其中帶通的帶通頻率和帶阻頻率于低通和高通不同。具體參數(shù)及建模結果見表3。 表3 巴特沃斯濾波器預處理的分類結果Table 3 Classification results of Butterworth filter preprocessing 由表3可知,RF模型和SVM模型在3種不同的Butterworth濾波器下分類準確率提升至94.2%。在兩類模型中,高通濾波器的分類準確率優(yōu)于低通和帶通。高通濾波器能夠抑制低頻分量并選擇高通分量通過,從而降低噪聲,篩選出合適的光譜信號。相較于直接對模型的參數(shù)進行優(yōu)化,通過High-Pass Butterworth濾波器的方法進行預處理,能夠有效提高模型準確率。進一步選擇Band-Pass Butterworth-BDA模型,對每一類樣本的分類準確率進行預測,結果表明70組樣本的總體識別準確率達98.6%,可以有效地對摻假有西布曲明的咖啡進行識別,但尚未達到100%準確分類。 2.2.3 特征選擇融合方法 Pearson相關系數(shù)是選擇重要特征變量的常用方法[22],可以用來檢測兩個連續(xù)型變量之間線性相關的程度,取值范圍為正負1的區(qū)間,在區(qū)間內大于0表示正相關,小于0表示負相關,絕對值越大表示線性相關程度越高。本實驗中通過Pearson方法對咖啡樣本的太赫茲光譜圖進行特征選擇,并選擇對應咖啡樣本中紅外光譜圖的指紋區(qū)域(650~1 300 cm-1)進行特征融合建模。邱薇綸等[23]通過對180份植物油的拉曼光譜和中紅外光譜進行融合建模研究,結果表明特征層光譜融合效果優(yōu)于數(shù)據(jù)層光譜融合,通過BDA和MLP模型實現(xiàn)分類準確率100%和97%。特征層的光譜融合主要是通過相關性計算的篩選主要特征變量。本實驗通過Pearson方法篩選出太赫茲時域光譜的主要特征頻率范圍,并且通過篩選發(fā)現(xiàn)咖啡的中紅外光譜數(shù)據(jù)在指紋區(qū)域內的特征變量均具有較強的關聯(lián)性。表4中給出了實驗樣本的太赫茲光譜在0~2.5 THz區(qū)間內的主要變量和中紅外光譜選擇的指紋區(qū)域。 表4 咖啡樣本的光譜特征頻率選擇Table 4 Selection of spectral characteristic frequencies of coffee samples 通過太赫茲時域光譜的特征頻率和中紅外光譜的指紋區(qū)域進行融合建模,模型準確率有效提高,SVM>BDA>RF其中SVM模型在特征融合下的準確率提高至100%,分類效果最好。能夠滿足快速、無損對咖啡中摻假非法藥物成分的準確鑒別。在確立最優(yōu)預處理方法建模的基礎上,選擇星空咖啡的三類不同產(chǎn)地的樣本進行識別研究,如圖6所示,能夠準確地對品牌-產(chǎn)地實現(xiàn)二維特征識別,為進一步確定摻假咖啡的來源提供偵查線索。 圖6 星空品牌樣本產(chǎn)地分類Fig.6 Classification of Starry Star brand samples 通過太赫茲時域光譜結合機器學習對70組含西布曲明的減肥咖啡的樣本建立一種“品牌分類-產(chǎn)地識別”的快速分析方法。通過一/二階導數(shù)、巴特沃斯濾波器、融合特征提取3種方法對250維THz-TDS光譜數(shù)據(jù)進行預處理,預處理方法能夠有效提取特征信息,消除部分實驗噪聲。結果表明,Band-Pass Butterworth-BDA模型、FD-BDA模型分類準確率均達98.6%,這表明BDA模型的分類效果較好,且2種預處理方法都能夠有效提升模型精度,但均未達100%。High-Pass-RF模型準確率為94.2%,較FD-RF、SD-RF的準確率高,這表明通過高頻分量的濾波器篩選光譜信息相較單一的特征峰校正和導數(shù)降噪可以有效地提升RF模型的準確率。最終,確立基于THz-TDs-FTIR兩種光譜特征變量融合的SVM模型對不同咖啡樣本的分類準確率最高100%,并進一步選擇此方法開展對產(chǎn)地的識別,結果表明能夠準確區(qū)分來自3個產(chǎn)地的同一品牌摻假咖啡。3 結論