油茶果自然霉變程度的可見(jiàn)/近紅外與中短波近紅外光譜檢測(cè)

姜洪喆，楊雪松，李興鵬，蔣雪松，周宏平，施明宏

（南京林業(yè)大學(xué)機(jī)械電子工程學(xué)院，江蘇 南京 210037）

油茶（Camellia oleiferaAbel.）是我國(guó)南方地區(qū)及少數(shù)東南亞國(guó)家特有的優(yōu)質(zhì)木本油料作物[1]。油茶果內(nèi)茶籽加工產(chǎn)品山茶油被譽(yù)為“東方橄欖油”，營(yíng)養(yǎng)成分豐富、具有很好的養(yǎng)生保健功效[2]。油茶果成熟期一般在10月底11月初，此時(shí)南方產(chǎn)地多陰雨天，為完成油茶果堆漚后熟過(guò)程貯藏不當(dāng)短時(shí)間內(nèi)極易發(fā)生外部自然霉變，霉變油茶對(duì)山茶油出油率、酸價(jià)、油氣滋味等品質(zhì)有很大影響[3]。

油茶霉變主要由曲霉屬、青霉屬與鐮刀菌屬等霉變微生物引起[4]，此類霉變微生物產(chǎn)出二級(jí)代謝產(chǎn)物，即真菌毒素，具有致癌、致畸和致基因突變的嚴(yán)重影響[5]，其中黃曲霉毒素B1等已被世衛(wèi)組織列為一級(jí)致癌物，常出現(xiàn)在霉變糧油制品中[6]。霉變果與正常果混合脫蒲、加工會(huì)帶來(lái)山茶油品質(zhì)下降、引入生物毒素等后果，因此在油茶籽榨油前對(duì)霉變油茶果“原材料”進(jìn)行檢測(cè)剔除具有重要意義。

常規(guī)的霉變農(nóng)產(chǎn)品檢測(cè)方法有酶聯(lián)免疫吸附法[7]、高效液相色譜法[8]、熒光染色法[9]、薄層色譜法[10]、高通量測(cè)序[11]等，此類方法操作繁瑣、樣品準(zhǔn)備復(fù)雜且具有破壞性。近紅外光譜技術(shù)是一種快速、無(wú)損、成本低、操作簡(jiǎn)便的現(xiàn)代光譜分析技術(shù)，近年已廣泛研究于各類農(nóng)產(chǎn)品品質(zhì)安全檢測(cè)可行性，具有較強(qiáng)的在線檢測(cè)潛力[12-13]。目前，在油茶產(chǎn)業(yè)領(lǐng)域研究中，近紅外光譜已經(jīng)嘗試應(yīng)用于茶油茶籽含油量[14-15]、產(chǎn)地鑒別[16]、脂肪酸含量[17]、摻假[18]等方面的檢測(cè)，但以上研究均以茶籽或茶油實(shí)驗(yàn)室檢測(cè)為對(duì)象，對(duì)于直接產(chǎn)地采后的油茶果檢測(cè)分選并沒(méi)有實(shí)際意義。針對(duì)油茶果霉變的檢測(cè)剔除，確保后續(xù)處理步驟中規(guī)避生物毒素風(fēng)險(xiǎn)、提升產(chǎn)品品質(zhì)安全、增加出口創(chuàng)匯，目前此方面相關(guān)研究鮮見(jiàn)報(bào)道。

針對(duì)以上問(wèn)題，為把好原材料關(guān)，本研究分別探索采用試可見(jiàn)/近紅外與中短波近紅外光譜結(jié)合化學(xué)計(jì)量學(xué)方法對(duì)不同霉變程度油茶果進(jìn)行檢測(cè)判別，提取指示霉變的特征光譜信息，建立準(zhǔn)確、魯棒的判別模型并對(duì)比結(jié)果，研究近紅外光譜檢測(cè)自然霉變油茶果的可行性，旨在為霉變油茶果剔除的產(chǎn)地處理需求提供一種快速、無(wú)損的檢測(cè)方案。

1 材料與方法

1.1 樣品及處理

油茶果樣品采集自湖北省隨州市隨縣洪山鎮(zhèn)油茶園，采摘于2020年10月23日霜降節(jié)氣，采后即寄運(yùn)至實(shí)驗(yàn)室。油茶果等分為3 組，每組至少120 個(gè)樣品，為獲得實(shí)際貯藏環(huán)境真實(shí)自然霉變樣本，在果農(nóng)指導(dǎo)下模擬不同貯藏條件，分別在常溫常濕（（25±2）℃、相對(duì)濕度（40±10）%）、常溫高濕度（（25±2）℃、相對(duì)濕度（70±20）%）、通風(fēng)、常溫高濕度（（25±2）℃、相對(duì)濕度（70±20）%）不通風(fēng)環(huán)境下貯藏2 周致其自然霉變[19-20]，最終獲得輕度霉變、中度霉變和重度霉變3 組類別的外部霉變油茶果。輕度外部霉變表觀上不光滑顏色較鮮果略暗，中度外部霉變顏色較深、可聞到異味，重度外部霉變可以明顯觀測(cè)到生長(zhǎng)的菌絲體、能聞到腐敗刺鼻氣味。

1.2 近紅外光譜采集

光譜采集均在室溫（（25±2）℃）下進(jìn)行，近紅外光譜采集使用儀器為Ocean Optics USB2000＋（350～1000 nm）與Ocean Optics Flame-NIR（950～1650 nm）兩臺(tái)光譜儀，光譜儀應(yīng)用前預(yù)熱30 min并進(jìn)行黑白校正，應(yīng)用OceanView軟件同時(shí)采集兩波段光譜，圖1為整套油茶果近紅外光譜采集系統(tǒng)示意圖，采用光纖漫反射模式進(jìn)行光譜采集。本研究光譜儀配套光源使用海洋光學(xué)公司HL-2000鹵鎢燈光源，該光源適用于360～2000 nm，輸出功率為8.8 mW，色溫3000 K。光源應(yīng)用SMA905光纖連接器與光纖進(jìn)行連接，光纖采用海洋光學(xué)6繞1光纖束設(shè)計(jì)的FOP-UV型號(hào)光纖，入射光纖和接收光纖分布于同一探頭內(nèi)，借助光纖探頭支架將光纖探頭與樣品表面充分貼合，經(jīng)參數(shù)優(yōu)化掃描次數(shù)設(shè)定為16 次，積分時(shí)間設(shè)定為100 ms。油茶果霉變不同于水果霉變從某一點(diǎn)展開(kāi)，油茶果樣品表面各點(diǎn)霉變狀態(tài)基本一致，因此在油茶果赤道處采集3 測(cè)點(diǎn)光譜進(jìn)行平均代表該樣品，減少采集光譜誤差，每測(cè)點(diǎn)間隔大致120°，確保采集到陽(yáng)面、陰面以及接合面處。最終輕度、中度和重度霉變各采集180 個(gè)樣品光譜，其中每組隨機(jī)選取120 個(gè)樣品光譜進(jìn)入建模訓(xùn)練集（120×3=360），每組其余60 個(gè)放入預(yù)測(cè)集（60×3=180）。

圖1 近紅外光譜采集系統(tǒng)示意圖Fig.1 Schematic diagram of near-infrared spectral collection system

1.3 光譜數(shù)據(jù)預(yù)處理

主成分分析（principal component analysis，PCA）是挖掘光譜內(nèi)部隱藏信息的有效方法，將存在線性相關(guān)的變量轉(zhuǎn)換為新的不相關(guān)變量，即PC，并基于累計(jì)的方差貢獻(xiàn)選取前幾個(gè)有效PC。PC得分圖可以用于反映不同樣本間的差異性和相似性，而PC載荷大小則可用來(lái)選取有效霉變程度判別波長(zhǎng)。

為消除近紅外光譜基線漂移、噪聲及散射效應(yīng)等影響，應(yīng)用標(biāo)準(zhǔn)正態(tài)變換（standard normal variate，SNV）、SNV結(jié)合去趨勢(shì)、一階導(dǎo)數(shù)和二階導(dǎo)數(shù)等預(yù)處理方法。其中，一階或二階導(dǎo)數(shù)預(yù)處理應(yīng)用的是2 次5點(diǎn)Savitzky-Golay平滑求導(dǎo)。

1.4 模型的建立與評(píng)價(jià)

偏最小二乘判別分析（partial least squaresdiscriminant analysis，PLS-DA）是基于偏最小二乘回歸模型預(yù)測(cè)值對(duì)比閾值進(jìn)行分類的判別分析方法[21]，適用于變量較多且存在多重共線性的情況。為提高建立模型魯棒性和通用性，本研究中PLS-DA模型采用“留一法”百葉窗式交叉驗(yàn)證確定最優(yōu)潛在變量（latent variable，LV）數(shù)量建立PLS-DA模型。

對(duì)于定性模型建立后的評(píng)價(jià)，應(yīng)用分類準(zhǔn)確率，混淆矩陣的靈敏度、特異度以及Kappa系數(shù)進(jìn)行模型評(píng)價(jià)。各參數(shù)越大說(shuō)明模型分類效果越好，其中靈敏度和特異度分別代表真陽(yáng)性率和真陰性率，Kappa系數(shù)用于評(píng)價(jià)分類結(jié)果一致性和可信度[22]。

1.5 特征波長(zhǎng)選取方法

連續(xù)投影法（successive projections algorithm，SPA）的大體思路為在光譜矩陣中找尋含有最低限度的冗余信息的變量組，使得變量間的共線性達(dá)到最小，比較投影向量的大小，以投影向量最大的波長(zhǎng)作為待選波長(zhǎng)，然后基于模型選擇最終的特征波長(zhǎng)[23]。一般SPA選取的波長(zhǎng)數(shù)量不會(huì)太多，循環(huán)一次后進(jìn)行多元線性回歸分析，得到預(yù)測(cè)的均方根誤差，最小的誤差值對(duì)應(yīng)的波長(zhǎng)即為所選。

1.6 數(shù)據(jù)分析

本研究中光譜數(shù)據(jù)分析以及繪圖主要應(yīng)用了Matlab 2017b（The MathWorks Inc.,MA,USA）、Unscrambler X10.1（CAMO,Trondheim,Norway）以及OriginPro 2017（OriginLab Inc.,MA,USA）軟件。

2 結(jié)果與分析

2.1 不同霉變程度油茶果光譜特征

不同霉變程度油茶果的400～1000 nm 和900～1700 nm兩波段漫反射光譜曲線以及各自對(duì)應(yīng)光譜標(biāo)準(zhǔn)差如圖2所示。不管是400～1000 nm還是900～1700 nm區(qū)間范圍，不同霉變程度油茶果展現(xiàn)出類似的光譜輪廓，且光譜均隨著霉變程度加深可以發(fā)現(xiàn)反射強(qiáng)度逐漸降低。此種相似的輪廓主要是因?yàn)橛筒韫麅?nèi)共性的成分結(jié)構(gòu)，而差異主要是霉變程度導(dǎo)致的內(nèi)部有機(jī)物化學(xué)成分變化。最左側(cè)456 nm左右的反射谷與邵氏帶有關(guān)）[24]，676 nm與植被固有的“紅邊”有關(guān)[25]，948 nm和967 nm附近的波谷主要與水的O—H拉伸態(tài)的二階泛音有關(guān)[26]，1215 nm與甲基（—CH3）的C—H鍵二階泛音有關(guān)[27]，1488 nm則與水的O—H鍵拉伸態(tài)的一階泛音[26]。

圖2 不同霉變程度油茶果帶標(biāo)準(zhǔn)差近紅外平均光譜圖Fig.2 Average infrared spectra with standard deviations of C.oleifera fruit with different mildew degrees

通過(guò)對(duì)比發(fā)現(xiàn)不同霉變程度組光譜有一定的差別，這是可以進(jìn)行判別的信息基礎(chǔ)。但還存在較多區(qū)域的模糊重疊，仍需后續(xù)從光譜數(shù)據(jù)角度結(jié)合化學(xué)計(jì)量學(xué)進(jìn)行深度解析，作出準(zhǔn)確判別。

2.2 PC解釋性分析

由于光譜變量的高維度性，不可避免存在冗余信息，本研究分別對(duì)可見(jiàn)/近紅外和中短波近紅外光譜進(jìn)行全光譜PCA，對(duì)數(shù)據(jù)降維、觀測(cè)不同組間光譜異同。其中可見(jiàn)/近紅外光譜的PC1（96.13%）和PC2（3.38%）累計(jì)占有99.51%方差貢獻(xiàn)率，而中短波近紅外光譜PCA變換前兩個(gè)PC（PC1與PC2）分別占有98.20%、1.52%的方差貢獻(xiàn)率，前兩個(gè)PC合計(jì)占有99.72%方差貢獻(xiàn)率，這說(shuō)明前兩個(gè)PC具有判別自然霉變程度的潛力。PC得分圖和載荷圖如圖3所示，在PC空間中可以初步看出不同自然霉變程度樣品具有一定聚類效果，圖3a、b中用橢圓標(biāo)出3 組樣品聚類的大致覆蓋區(qū)域，PC得分圖顯示了不同霉變程度油茶果的可分性[28]。

圖3 PC得分圖和載荷圖Fig.3 PC score plots and PC loading plots of PC1 versus PC2

圖3c、d詳細(xì)的展示了前兩個(gè)PC的載荷曲線，其中背離水平線的峰和谷對(duì)判別自然霉變程度有效，針對(duì)可見(jiàn)/近紅外光譜選取4 個(gè)特征波長(zhǎng)（461、650、674 nm和945 nm），中短波近紅外光譜選取5 個(gè)特征波長(zhǎng)（973、1066、1146、1482 nm和1599 nm），這些波長(zhǎng)也嘗試應(yīng)用在后續(xù)簡(jiǎn)化模型建立中。

2.3 基于PLS-DA的全光譜判別模型

針對(duì)原始以及4 種不同預(yù)處理后的全光譜，建立的最優(yōu)PLS-DA模型參數(shù)及判別結(jié)果如表1所示。整體上，兩個(gè)光譜集全光譜建立的PLS-DA模型性能差異不大。對(duì)于可見(jiàn)/近紅外光譜，SNV預(yù)處理全光譜預(yù)測(cè)集結(jié)果最優(yōu)，分類準(zhǔn)確率達(dá)85.0%。原始光譜建模預(yù)測(cè)集準(zhǔn)確率達(dá)83.3%，而其他預(yù)處理如二階導(dǎo)數(shù)則降低了原始光譜建模性能。但原始全光譜建立的模型訓(xùn)練集、交叉驗(yàn)證集以及預(yù)測(cè)集3 者之間分類準(zhǔn)確率差值最小，因此原始光譜建立模型最為魯棒。針對(duì)中短波近紅外光譜，原始光譜建立模型最優(yōu)LV為10，訓(xùn)練集、交叉驗(yàn)證集與預(yù)測(cè)集分類準(zhǔn)確率分別達(dá)90.0%、86.7%和85.0%。結(jié)果發(fā)現(xiàn)，無(wú)論400～1000 nm還是900～1700 nm，各種預(yù)處理方法并沒(méi)有對(duì)原始光譜起到明顯的模型改善作用，原始光譜具有足夠的有效信息，可以不需要對(duì)原始光譜數(shù)據(jù)進(jìn)行額外的預(yù)處理[29]，后續(xù)選取原始光譜建立的PLS-DA模型進(jìn)行分析。

表1 不同預(yù)處理方式全光譜PLS-DA模型性能比較Table 1 Performance comparison of full-spectral PLS-DA models with different spectral pre-processing methods

2.4 波長(zhǎng)篩選

全光譜模型的使用會(huì)由于信息冗余帶來(lái)過(guò)擬合、噪聲等風(fēng)險(xiǎn)，導(dǎo)致模型較為復(fù)雜，不適合實(shí)際應(yīng)用。從全光譜中進(jìn)行有效波長(zhǎng)的選取，找到攜帶重要信息的波長(zhǎng)建立簡(jiǎn)化模型將很有實(shí)際意義。圖4的SPA中，根據(jù)以往文獻(xiàn)設(shè)定所選變量最小數(shù)量設(shè)定為5，為了實(shí)際應(yīng)用的需要最大數(shù)量設(shè)定為15[30]。對(duì)于兩組光譜集，分別在變量選取8 個(gè)和7 個(gè)時(shí)均方根誤差達(dá)到最小，最終分別選取了8 個(gè)和7 個(gè)波長(zhǎng)，如圖4所示。綜合分析，不同方法選取的波長(zhǎng)見(jiàn)表2。

圖4 應(yīng)用SPA方法選取波長(zhǎng)Fig.4 Wavelength selection using SPA

表2 根據(jù)不同方法選取波長(zhǎng)匯總Table 2 Summary of wavelengths selected using different methods

2.5 簡(jiǎn)化模型

表3中列出了不同波長(zhǎng)選擇方法選取的特征波長(zhǎng)建立簡(jiǎn)化PLS-DA模型性能對(duì)比，在400～1000 nm波段，PC-PLS-DA模型表現(xiàn)略差，訓(xùn)練集、交叉驗(yàn)證集以及預(yù)測(cè)集3 者判別準(zhǔn)確率均沒(méi)有超過(guò)80%。SPA-PLS-DA模型獲得了最好的分類準(zhǔn)確率，預(yù)測(cè)集達(dá)到了84.4%。對(duì)于900～1700 nm波段，同樣是SPA-PLS-DA模型表現(xiàn)最好，在訓(xùn)練集、交叉驗(yàn)證集和預(yù)測(cè)集分別獲得84.7%、84.2%和84.4%的分類準(zhǔn)確率。因此SPA選擇的波長(zhǎng)可以兼顧降低變量個(gè)數(shù)以及保持模型性能，波長(zhǎng)數(shù)分別由1506 個(gè)和128 個(gè)下降到8 個(gè)和7 個(gè)，而PC選取的波長(zhǎng)有效性不足的原因可能是二者選取波長(zhǎng)過(guò)程僅針對(duì)光譜變量（X變量）沒(méi)有考慮霉變組別（Y變量），選取的有效波長(zhǎng)較少，難以覆蓋識(shí)別中的所有有效信息[31]。

表3 不同波長(zhǎng)選擇方法選擇變量建模結(jié)果對(duì)比Table 3 Comparison of modeling results using different wavelength selection methods

為進(jìn)一步直觀解析最優(yōu)的簡(jiǎn)化SPA-PLS-DA模型，表4中列出了預(yù)測(cè)集混淆矩陣以及靈敏度、特異度、Kappa系數(shù)計(jì)算結(jié)果。兩組簡(jiǎn)化模型的Kappa系數(shù)均為0.7667，說(shuō)明在一致性分類精度上兩簡(jiǎn)化模型能力相似。兩模型可以較為準(zhǔn)確地判別出輕度霉變和重度霉變果，可見(jiàn)/近紅外光譜均判別正確60中的57 個(gè)（靈敏度≥0.95），中短波近紅外光譜此兩組靈敏度≥0.88，特異性兩模型均超過(guò)0.84。但兩組簡(jiǎn)化模型對(duì)中度霉變油茶果判別準(zhǔn)確率均不足，主要原因是該組樣品既有可能誤判為輕度霉變和重度霉變，而此方面中短波近紅外光譜模型（靈敏度0.72）優(yōu)于可見(jiàn)/近紅外光譜（靈敏度0.63）。

表4 最優(yōu)簡(jiǎn)化模型預(yù)測(cè)集混淆矩陣Table 4 Confusion matrices for prediction set using optimal simplified models

3 討論與結(jié)論

采集540 個(gè)不同自然霉變程度油茶果近紅外光譜，包括可見(jiàn)/近紅外以及中短波近紅外兩個(gè)波段分別進(jìn)行判別可行性分析。獲取的光譜PCA發(fā)現(xiàn)前兩個(gè)PC得分圖上具有一定的聚類效果，說(shuō)明光譜具有判別不同霉變程度油茶果的能力。建模中發(fā)現(xiàn)采用SNV、SNV＋去趨勢(shì)、一階導(dǎo)數(shù)、二階導(dǎo)數(shù)等預(yù)處理后全光譜建立的PLS-DA模型效果不如原始光譜，說(shuō)明原始全光譜具有足夠的有效信息，不需要再進(jìn)行額外的預(yù)處理。另外兩個(gè)光譜集建模能力相當(dāng)，最優(yōu)全光譜模型的預(yù)測(cè)集判別準(zhǔn)確率均在82%以上。為了進(jìn)一步適用于實(shí)際，應(yīng)用PC載荷和SPA選取特征波長(zhǎng)，建立的簡(jiǎn)化模型結(jié)果顯示無(wú)論400～1000 nm還是900～1700 nm波段均為SPA-PLS-DA模型性能最佳，預(yù)測(cè)集判別準(zhǔn)確率相比于全光譜沒(méi)有顯著下降（84.4%），說(shuō)明SPA選取的特征波長(zhǎng)因考慮了目標(biāo)組別信息在判別霉變程度中更為有效。

結(jié)合混淆矩陣對(duì)最優(yōu)簡(jiǎn)化模型進(jìn)行深入分析，兩模型Kappa系數(shù)均為0.7667，靈敏度和特異度相差無(wú)幾，900～1700 nm對(duì)中度霉變樣品分辨力更強(qiáng)，但是從儀器制造成本上400～1000 nm更便宜更實(shí)用。綜合分析，本研究表明近紅外光譜為油茶果產(chǎn)地處理提供了一種快速、原位、無(wú)損的現(xiàn)場(chǎng)檢測(cè)方法，克服了以往常規(guī)化學(xué)測(cè)定方法費(fèi)時(shí)費(fèi)力、試劑污染以及儀器昂貴等缺點(diǎn)，也為木本糧油林果霉變的快速無(wú)損檢測(cè)提供了有效借鑒。未來(lái)為了提高簡(jiǎn)化模型的實(shí)際應(yīng)用性和魯棒性，還需要加大近紅外光譜采集樣品量不斷進(jìn)行模型的修正。