王 倩，葛宏義,2,*，蔣玉英,2，張 元,2，秦一菲,2

（1.河南工業(yè)大學信息科學與工程學院，河南 鄭州 450001；2.河南工業(yè)大學糧食信息處理與控制教育部重點實驗室，河南 鄭州 450001）

隨著人們生活品質(zhì)的不斷提高，消費者越來越熱衷于購買地域性特色大米，如：黑龍江五常大米、河南原陽大米、寧夏珍珠糯米和云南遮放貢米等。不同產(chǎn)地的大米不僅具有特殊的品質(zhì)和信譽，且經(jīng)濟價值也比普通大米高。不良商家為獲取不當高額利潤，通過對低價大米染色或工業(yè)加工等方式偽造特色大米，牟取暴利，造成大米市場紊亂[1]，嚴重損害消費者和生產(chǎn)者權(quán)益。因此，亟需相應的檢測方法識別大米品種，以規(guī)范大米市場。

針對大米品種識別方法較多，傳統(tǒng)檢測方法主要有感官評價法[2]和理化指標技術(shù)法[3]等，具有分析簡單、方法直觀等優(yōu)勢，但也存在一些不足，如感官評價法主觀性強，對于采用化學劑染色大米或工業(yè)加工研磨大米，僅從外觀上難以判別真?zhèn)?；理化指標技術(shù)法操作繁瑣、檢測速度慢等。目前常見識別技術(shù)主要有電子舌技術(shù)[4]、紅外光譜技術(shù)[5]和高光譜圖像技術(shù)[6]，諸多學者以此為基礎結(jié)合化學計量方法對大米識別做出研究，如惠延波等[7]利用電子舌技術(shù)結(jié)合PCA、DFA 兩種模式對不同品種粉碎大米進行識別研究，區(qū)分指數(shù)D1 為95；王靖會等[8]以高光譜成像技術(shù)為基礎結(jié)合光譜和紋理特性實現(xiàn)不同大米品種的識別，識別率為96.57%；王朝輝等[9]利用近紅外光譜技術(shù)結(jié)合SG 平滑預處理方法和偏最小二乘判別分析方法對松原不同品種大米識別，識別率為100%；劉亞超等[10]利用近紅外二維相關光譜對摻和大米識別，識別率最高為100%；雖然以上技術(shù)能有效解決傳統(tǒng)技術(shù)檢測精度不高、主觀性強等缺點，但也存在一些不足，如電子舌技術(shù)檢測結(jié)果易受傳感器靈敏度和環(huán)境因素影響[11]、設備集成度不高、價格昂貴、使用周期短等[12]；紅外光譜技術(shù)靈敏度相對較低[13]、檢測結(jié)果易受各種環(huán)境影響[14]，紅外光譜解釋性差、化學計算依賴性高[15]；高光譜圖像技術(shù)圖像處理復雜性高、檢測速度慢等[16]。太赫茲（Terahertz,THz）[17]是介于微波和紅外之間的電磁波，頻率在0.1～10 THz 之間，具有穿透性強、電離輻射小、分子指紋等特點，還具有檢測分析速度快、多成分同步分析、無損檢測等優(yōu)勢，能有效彌補傳統(tǒng)檢測方法和目前常見技術(shù)的不足。其原理是利用飛秒脈沖產(chǎn)生THz 電磁波，通過探測設備獲得待測樣品的光譜信息，再通過傅里葉變換獲得測物品吸收和折射光譜信息。由于大分子振動和轉(zhuǎn)動能級大多處于THz 波段，并表現(xiàn)出明顯吸收特性，因而可以通過光譜的指紋特性對大米品種進行識別。太赫茲時域光譜技術(shù)已被廣泛應用到安全檢查[18]、生物醫(yī)學[19]、無線通信[20]等領域，目前在農(nóng)產(chǎn)品品質(zhì)檢測方面的應用也在不斷增多，如農(nóng)產(chǎn)品摻假[21]、轉(zhuǎn)基因作物檢測[22]和農(nóng)產(chǎn)品主要成分含量檢測[23-24]等。太赫茲光譜技術(shù)的廣泛應用，為大米品種識別提供了新思路。

本文采用太赫茲時域光譜技術(shù)結(jié)合標準差分析、區(qū)間偏最小二乘和決策樹分類模型對不同大米品種進行識別。分析了大米在光譜技術(shù)作用下的特點，建立了一種準確識別大米品種的分析模型，為大米品種識別提供了一種精準檢測方法。

1 材料與方法

1.1 材料與儀器

實驗樣品分別是河南紅米、珍珠糯米、黑米和富硒大米 均采購于永輝超市，如圖1 所示。大米顆粒利用研磨機研磨成粉末狀，再利用電子稱每次取0.2 g 的粉末放入模具，并對模具使用10 MPa 的壓力施壓3 min，使粉末壓制成直徑約為13 mm、厚度約為1.1 mm、表面均勻的圓形薄片。每種大米制備一組樣品，每組樣品制備14 份，四組樣品共56 份。

圖1 實驗樣品Fig.1 Experimental samples

太赫茲時域光譜儀 由河南工業(yè)大學教育部重點實驗室提供，激光器波長800 nm，光譜范圍0～3.5 THz，分辨率0.03 THz，重復頻率80 MHz，脈沖波長100 fs，信噪比5000:1，具體系統(tǒng)結(jié)構(gòu)如圖2所示。

圖2 太赫茲時域光譜系統(tǒng)透射式（a）和反射式（b）Fig.2 Terahertz time domain spectroscopy system transmission type (a) and reflection type (b)

時域光譜系統(tǒng)工作原理如下：飛秒激光器產(chǎn)生的光束作為輸入光源，經(jīng)分束鏡分成能量較強的泵浦光和能量相對弱的探測光。泵浦光先傳輸至時間延遲控制系統(tǒng)，再入射到光導天線或半導體晶體上激發(fā)出太赫茲脈沖，最后經(jīng)兩組拋物面鏡聚焦到待測物體上。探測光經(jīng)多次反射后和太赫茲脈沖共線觸發(fā)太赫茲探測器，探測器通過對脈沖偏振狀態(tài)的檢測獲得待測樣品的時域波形。在整個測試過程中，由于水分對太赫茲脈沖有強吸收性，因此為減少空氣中水分含量對實驗結(jié)果的影響，需要不斷對時域光譜系統(tǒng)補充氮氣，以保持實驗環(huán)境干燥。

1.2 實驗步驟

實驗分為三個階段：第一階段主要是對儀器參數(shù)和實驗環(huán)境進行調(diào)整，光譜系統(tǒng)設置為透射模式，測試溫度保持在19 ℃左右，濕度保持在2.8%左右，并在測試正式開始前，儀器持續(xù)運行2 min，以確保各項參數(shù)指標達到穩(wěn)定狀態(tài)；第二階段是對樣品進行檢測，把制備好的樣品放置在檢測臺旁，按標簽順序依次放入太赫茲時域光譜系統(tǒng)中進行檢測，重復測量三次，并取平均值作為該樣品時域光譜信息。為減少環(huán)境因素對實驗結(jié)果的影響，取樣間隔時間不超過5 s；第三階段是利用Origin 和Python 軟件結(jié)合分類模型算法對獲得的太赫茲光譜數(shù)據(jù)信息進行數(shù)據(jù)預處理和數(shù)據(jù)分類。

1.3 模型與建模方法

1.3.1 標準差 標準差（SD）值反映的是數(shù)據(jù)集的波動程度，即離散程度[25]。標準差值大表示數(shù)據(jù)和平均值之間差異較大，標準差值小表示數(shù)據(jù)和平均值之間差異小，在光譜信息中標準差值大小表示光譜波動范圍，且兩者為正比關系。標準差的計算公式如下：

式中：xi為樣品數(shù)據(jù)；u 為算術(shù)平均值；N 為數(shù)值個數(shù)。

1.3.2 區(qū)間偏最小二乘 區(qū)間偏最小二乘（iPLS）主要用于光譜波段的選擇，其特點是采用偏最小二乘實現(xiàn)光譜波段的局部回歸分析，并根據(jù)局部回歸分析的最小均方根誤差確定光譜特征波段的選擇。偏最小二乘的原理是將全光譜劃分為多個具有同樣寬度的區(qū)間，然后在各個子區(qū)間內(nèi)使用偏最小二乘回歸，根據(jù)各個子區(qū)間均方根誤差的比較，選擇均方根誤差最小的子區(qū)間對應的光譜波段為光譜特征波段[26]。

1.3.3 決策樹 決策樹（DT）模型是一種類似于二叉樹的網(wǎng)絡結(jié)構(gòu)，具有可讀性、分類速度快等優(yōu)點，在分類、回歸等方面有廣泛的應用。決策樹中每一個非葉子節(jié)點都是一個決策點，即判斷條件，滿足條件的放在節(jié)點的右側(cè)，不滿足的放在節(jié)點的左側(cè)。如何確定決策點是構(gòu)造決策樹的重點，根據(jù)決策點特征選擇的不同，常見的方法有ID3 算法、C4.5 算法和CART 算法。ID3 算法的核心是在決策樹各個節(jié)點上應用信息增益準則選擇特征，選擇信息增益最大的特征作為決策點。但ID3 算法沒有考慮到數(shù)據(jù)的連續(xù)性和過擬合問題，C4.5 算法在彌補ID3 算法的不足的同時也提出了根據(jù)信息增益率作為特征選擇的準則。C4.5 和ID3 算法都能有效的解決分類問題，但不能處理回歸問題，而根據(jù)CART 算法構(gòu)建的決策樹不僅能處理分類問題也能做回歸分析，在分類問題時選擇基尼系數(shù)最小的特征作為決策點，在回歸分析中CART 算法使用平方誤差最小值對應的特征作為決策點[27]。

1.3.4 SD-iPLS-DT 由于原始數(shù)據(jù)集維度高、信噪比低，直接使用決策樹模型進行分類處理，易造成數(shù)據(jù)處理速度慢，模型分類準確率低等問題。因此采用SD-iPLS 與DT 聯(lián)用的方法對測試樣品進行分類，即先采用SD 和iPlS 選擇合適的太赫茲波段作為模型的輸入數(shù)據(jù)，再利用DT 模型進行分類識別。本次實驗中原始數(shù)據(jù)維度為220，通過SD 選擇穩(wěn)定性較好的光譜波段，同時結(jié)合iPLS 選擇均方根誤差最小的光譜區(qū)間，根據(jù)兩者最優(yōu)結(jié)果，選擇其中的38 維吸收光譜數(shù)據(jù)作為分類模型的輸入數(shù)據(jù)，最后使用DT 算法進行大米品種識別研究。

1.4 數(shù)據(jù)處理

由于光譜信息中包含實驗樣品的振幅和相位信息，因此可以將時域光譜信號利用快速傅里葉變換得到實驗樣品的頻域信息，并結(jié)合Dorney[28]和Dubillaret[29]等提出的數(shù)據(jù)處理模型處理光譜數(shù)據(jù)，得到大米樣品的折射率和吸收系數(shù)。計算公式如下：

式中：n(ω) 表 示光譜折射率，n； ω代表角速度，rad/s；φs、 φr分別表示樣品和參考信號的相位信息，p；c 表示光譜傳播速度，m/s；d 表示樣品的厚度，mm。

式中： α (ω)表 示光譜吸收系數(shù)，cm-1； Ar、 AS分別表示參考信號和樣品的振幅值，a.u。

2 結(jié)果與分析

2.1 光譜分析

2.1.1 時域和頻域 為減少實驗誤差，對每份實驗樣品的時域信號分別測量3 次，并取平均值作為樣品的時域信號。圖3（a）是樣品的時域光譜，時域光譜信息經(jīng)傅里葉變換得到樣品的頻譜信息，如圖3（b）所示?？梢钥闯鰳悠沸盘栂鄬τ趨⒖夹盘栍幸欢ǖ臅r延、衰減和重疊現(xiàn)象，產(chǎn)生時延可能是由于光譜穿透樣用時較長，衰減可能是由樣品表面的大顆粒物質(zhì)對太赫茲光譜反射和內(nèi)部樣品的吸收造成的，而光譜之間區(qū)分度不明顯可能是由于各種大米中的主要成分都是碳水化合物、蛋白質(zhì)和脂肪。

圖3 時域光譜圖（a）和頻域光譜圖（b）Fig.3 Time domain spectroscopy (a) and frequency domain spectroscopy (b)

2.1.2 折射率和吸收系數(shù) 太赫茲時域光譜數(shù)據(jù)中包含豐富的振幅和相位信息，經(jīng)公式（2）和（3）計算得到樣品的折射率和吸收系數(shù)，如圖4 所示?？梢钥闯鏊姆N樣品的折射率曲線和吸收光譜曲線在低頻波段均出現(xiàn)重疊現(xiàn)象，在高頻波段折射率光譜曲線趨于一致，吸收光譜曲線區(qū)別明顯，因此本文選擇吸收光譜作為模型的輸入數(shù)據(jù)。

圖4 折射光譜（a）和吸收光譜（b）Fig.4 Refraction spectrum (a) and absorption spectrum (b)

2.2 光譜處理

由于吸收光譜存在信噪比低和光譜重疊現(xiàn)象，為更好識別大米品種，本文首先對吸收光譜數(shù)據(jù)進行標準化處理解決光譜重疊問題，再利用標準差分析和區(qū)間偏最小二乘（iPlS）選取光譜穩(wěn)定性好、信噪比高的太赫茲波段作為分類模型的輸入數(shù)據(jù)。

2.2.1 吸收光譜預處理和穩(wěn)定性分析 為增大吸收光譜曲線之間的差異性，本次實驗采用標準化預處理方法對吸收光譜數(shù)據(jù)進行預處理操作。圖5 為預處理后四種樣品的吸收光譜數(shù)據(jù)平均值，通過觀察光譜曲線可以看出，四種樣品的吸收光譜曲線差異明顯。同時，利用光譜角度[30]對光譜曲線之間的差異性進行評價，以河南紅米吸收光譜為參照，結(jié)果如表1 所示，可以看出在經(jīng)過預處理之后，光譜角度的數(shù)值明顯增大，即光譜曲線之間的差異性增大。但標準化預處理后的光譜仍存在信噪比低、光譜波動大等問題，因此為減少噪聲和光譜波動性對模型準確率的影響，需選擇合適的光譜波段作為模型分類識別的輸入數(shù)據(jù)。

圖5 吸收數(shù)據(jù)標準化預處理Fig.5 Standardized preprocessing of absorption data

表1 不同樣品光譜角度Table 1 Different sample spectral angles

根據(jù)不同光譜波段的標準差對光譜的穩(wěn)定性進行分析[31]，實驗過程中把吸收全光譜劃分為7 個不同的區(qū)間，并對每個區(qū)間進行標準差分析，結(jié)果如表2、圖6 所示，可以看出0.53～1.21 THz 之間四種大米樣品標準差分別為0.06、0.07、0.08、0.06，平均值為0.07，是8 組光譜數(shù)據(jù)中標準差值最小的一組，即光譜波動最小，穩(wěn)定性最好的一組；0～0.53 THz 之間的四種大米樣品標準差分別為0.65、0.53、0.36、0.12，平均值為0.42，是8 組光譜數(shù)據(jù)中標準差值最大的一組，即光譜波動最大，穩(wěn)定性最差的一組。

圖6 四種樣品預處理后不同波段標準差Fig.6 Standard deviation of different bands after pretreatment of four samples

表2 四種樣品預處理后不同波段標準差Table 2 Standard deviation of different bands after pretreatment of four samples

2.2.2 吸收光譜特征譜區(qū)間選擇 實驗過程中采用SPXY-iPLS 算法選擇吸收光譜特征區(qū)間。首先利用SPXY 算法[32]進行數(shù)據(jù)集的劃分，再利用iPLS算法計算每個子區(qū)間的均方根誤差。在子區(qū)間劃分中，把預處理后的全光譜分為2、3、4、5、6、7、8、10、15 和22 個子區(qū)間，分別對各子區(qū)間建立PLS 回歸模型[33]，進而比較模型的均方根誤差，確定均方根誤差最小的回歸模型對應的子區(qū)間。由表3、圖7可以看出當劃分區(qū)間個數(shù)為2 時，最佳區(qū)間對應的均方根0.84，是10 組區(qū)間中均方根誤差最大的區(qū)間；當劃分區(qū)間個數(shù)為3 時，最佳區(qū)間1 區(qū)間對應的均方根誤差為0.50，是10 組區(qū)間中均方根誤差最小的區(qū)間，即模型效果最好，對應的THz 光譜波段為0～1.32 THz。

表3 iPLS 不同分割區(qū)間下最優(yōu)區(qū)間列表Table 3 List of optimal intervals under different segmentation intervals of interval partial least square

圖7 最優(yōu)區(qū)間均方根誤差Fig.7 Optimal interval root mean square error

結(jié)合光譜穩(wěn)定性和區(qū)間偏最小二乘兩者對應的最優(yōu)波段，最終選擇0.53～1.21 THz 之間的吸收光譜數(shù)據(jù)作為模型的輸入數(shù)據(jù)。

2.3 模型分析

采用SD-iPLS-DT 聯(lián)用模型對樣品進行分類。首先對數(shù)據(jù)集進行劃分，將0.53～1.21 THz 之間的光譜數(shù)據(jù)按照2:1 的比例劃分為訓練集和測試集。為實現(xiàn)模型最優(yōu)分類結(jié)果，結(jié)合網(wǎng)格搜索算法尋找最優(yōu)參數(shù)組合。在決策樹模型中主要的參數(shù)有特征選擇評價標準（criterion）、最大深度（max_depth）。其中在決策樹模型中特征評價標準主要參數(shù)有基尼系數(shù)（gini），信息增益（entropy）。根據(jù)網(wǎng)格搜索算法最終確定模型的最優(yōu)參數(shù)criterion 為gini、最大深度為3。模型對應的分類結(jié)果如表4 所示，分類準確率為95%。

表4 決策樹分類準確率Table 4 Classification accuracy of decision tree

為證明本文模型的性能，將SD-iPLS-DT 模型與常見的分類模型邏輯回歸（LC）和支持向量機[34]（SVM）進行對比。其中在SVM 模型中，分別利用徑向基核函數(shù)、線性核函數(shù)和Sigmoid 核函數(shù)等三種不同核函數(shù)對樣本進行分類。如表5 所示，可以看出在SVM 模型中基于線性核函數(shù)的分類模型效果最好，分類準確率可達88.75%。由表6 可以看出，直接使用邏輯回歸和支持向量機分類準確率分別為80.75%和88.75%，使用標準差和iPLS 結(jié)合選擇合適的光譜波段后，再使用決策樹模型進行分類準確率為95%，因此使用SD-iPLS-DT 聯(lián)用模型效果要優(yōu)于邏輯回歸和支持向量機模型。

表5 不同核函數(shù)分類準確率Table 5 Classification accuracy of different kernel functions

表6 不同模型的分類準確率Table 6 Classification accuracy of different models

3 結(jié)論與討論

本文利用太赫茲時域光譜技術(shù)獲取四種大米的時域和頻域光譜數(shù)據(jù)，并計算得到折射率和吸收系數(shù)。采用標準差對吸收光譜數(shù)據(jù)穩(wěn)定性進行分析，得出0.53～1.21 THz 之間的光譜穩(wěn)定性最好；區(qū)間偏最小二乘選擇吸收光譜特征波段，確定0～1.32 THz 區(qū)間的均方根誤差最小。結(jié)合穩(wěn)定性最好和最小均方根誤差對應的太赫茲波段，選擇0.53～1.21 THz 波段吸收光譜信息作為決策樹模型的輸入數(shù)據(jù)。實驗結(jié)果顯示：所提出的SD-iPLS-DT 的方法分類準確率可達95%。為更好對比模型分類效果，利用邏輯回歸、支持向量機與SD-iPLS-DT 方法進行實驗對比，實驗結(jié)果顯示：邏輯回歸分類準確率為80.75%，支持向量機準確率為88.75%，SD-iPLS-DT 模型分類準確率為95%。識別準確率高于現(xiàn)有利用近紅外光譜技術(shù)進行稻花香大米定性分析研究94%準確率[35]、利用高光譜成像技術(shù)對大米無損檢測92.96%準確率[36]、利用拉曼光譜大米產(chǎn)地識別91.11%準確率[37]等。因此可以得出結(jié)論：時域光譜技術(shù)結(jié)合SD-iPLSDT 方法可以實現(xiàn)不同大米品種的準確識別。同時也為太赫茲時域光譜技術(shù)在農(nóng)產(chǎn)品識別和質(zhì)量安全檢測方面提供了新的檢測方法。