基于高光譜技術(shù)及SPXY和SPA的玉米毒素檢測模型建立

于慧春，婁 楠，殷 勇*，劉云宏

新鮮玉米在存貯過程中，由于其胚部大、水分含量高、帶菌量多，高溫高濕環(huán)境下極易霉變[1]，不僅給經(jīng)濟(jì)造成重大損失，而且霉變玉米在代謝過程中會產(chǎn)生多種對人體具有極強(qiáng)致病性、致癌性的毒素[2-3]，危害人畜健康。黃曲霉毒素B1（aflatoxin B1，AFB1）和玉米赤霉烯酮（zearalenone，ZEN）是其中2 種比較穩(wěn)定的代謝產(chǎn)物，容易在霉變玉米中積累，導(dǎo)致含量升高，AFB1和ZEN的多少與玉米霉變情況密切相關(guān)[4]。因此，可以通過監(jiān)測玉米中AFB1和ZEN含量變化表征玉米的霉變情況，實(shí)現(xiàn)對玉米霉變程度的準(zhǔn)確評價(jià)[5]。但是常規(guī)的AFB1和ZEN檢測方法存在樣品處理繁瑣、費(fèi)時(shí)、對樣品有破壞性等缺點(diǎn)，難以實(shí)現(xiàn)簡單、快速、無損檢測，無法滿足實(shí)際生產(chǎn)的需要[6-8]。

高光譜技術(shù)通過提取被測對象的圖像和光譜信息特征并與其內(nèi)外品質(zhì)指標(biāo)間建立聯(lián)系，從而實(shí)現(xiàn)對其內(nèi)外綜合品質(zhì)的評價(jià)，因而在果蔬品質(zhì)[9-11]、禽肉品質(zhì)[12-14]、蛋[15]和煙葉[16]等農(nóng)產(chǎn)品檢測領(lǐng)域獲得廣泛應(yīng)用研究。在霉變玉米的高光譜檢測方面，也有學(xué)者進(jìn)行了一些探索[17-19]，但目前研究主要集中在高光譜技術(shù)對不同霉變程度玉米的定性判別方面，而對霉變玉米中AFB1和ZEN等的定量檢測分析方面尚鮮有報(bào)道。袁瑩[20]和褚璇[21]等將不同濃度的AFB1溶液滴在玉米籽粒上，然后利用高光譜技術(shù)進(jìn)行檢測，但其最終的分析結(jié)果仍然只是對表面含有不同濃度AFB1的玉米樣品進(jìn)行分類識別，也未嘗試通過高光譜技術(shù)對測試樣本在特定條件下自身代謝生成的AFB1含量進(jìn)行預(yù)測。此外，上述霉變玉米高光譜技術(shù)檢測研究中，對所使用的校正集樣本質(zhì)量也未進(jìn)行考察，大量的樣本數(shù)據(jù)之間可能存在相同或者差異過小的情況，從而影響分析結(jié)果精度和耗費(fèi)建模時(shí)間。

本研究以玉米霉變過程中產(chǎn)生的AFB1和ZEN 2 種代謝產(chǎn)物含量為玉米霉變程度的表征指標(biāo)，通過測定不同霉變程度玉米樣本的光譜信息和相應(yīng)的AFB1和ZEN含量，建立基于較少校正集樣本和特征波長光譜信息的AFB1和ZEN含量偏最小二乘回歸（partial least squares regression，PLSR）預(yù)測模型，從而實(shí)現(xiàn)對玉米霉變程度的準(zhǔn)確分析，以期為實(shí)現(xiàn)玉米霉變的在線、快速、精確檢測提供借鑒。

1 材料與方法

1.1 材料

玉米來自于2016年10月在洛陽當(dāng)?shù)厥斋@，品種為中單909，購于洛陽市中原農(nóng)貿(mào)城，由實(shí)驗(yàn)室自行培育出不同霉變程度的玉米。

1.2 儀器與設(shè)備

LHS-HC-100恒溫恒濕培養(yǎng)箱 上海資一儀器設(shè)備有限公司；IST 50-3810高光譜成像光譜儀 德國Inno-spec公司；康標(biāo)達(dá)鏡頭 日本Computar公司；RK90000420108線性鹵素?zé)簦?00 W） 德國Esylux公司。

高光譜圖像采集系統(tǒng)如圖1所示，該系統(tǒng)主要包括高光譜成像光譜儀、康標(biāo)達(dá)鏡頭、4 個(gè)500 W線性鹵素?zé)簟?shí)驗(yàn)室自制傳送裝置和計(jì)算機(jī)等主要部件。其中高光譜成像光譜儀由攝像機(jī)和光譜儀2 部分組成，攝像機(jī)為CCD相機(jī)，光譜攝像儀為可見-近紅外光譜儀，光譜波長范圍為371.05～1 023.82 nm，可采集到的光譜數(shù)據(jù)點(diǎn)數(shù)為1 288 個(gè)，光譜分辨率為2.8 nm，傳送帶移動速率1.25 mm/s，曝光時(shí)間100 ms，樣本與鏡頭的距離為285 mm。高光譜成像儀通過USB 2.0接口數(shù)據(jù)線連接計(jì)算機(jī)，通過SICap-STVR V1.0.x軟件平臺驅(qū)動控制成像儀，及時(shí)記錄和存貯高光譜數(shù)據(jù)。

圖1 高光譜圖像采集系統(tǒng)Fig. 1 Schematic presentation of the hyperspectral image acquisition system

1.3 方法

1.3.1 玉米樣本制備

將新鮮玉米放入恒溫恒濕培養(yǎng)箱，在培養(yǎng)溫度30 ℃、相對濕度90%[22]的條件下進(jìn)行培養(yǎng)。培養(yǎng)時(shí)間不同，玉米霉變程度也不同，本實(shí)驗(yàn)用培養(yǎng)初始（新鮮玉米）、第2天、第4天、第6天、第8天標(biāo)記玉米的不同霉變程度，5 個(gè)霉變程度的玉米看作5 個(gè)霉變等級，每個(gè)等級樣品各取50 個(gè)測試樣本，每個(gè)樣本含60 g（±0.1 g）玉米。

1.3.2 AFB1和ZEN含量的測定

為了減少AFB1和ZEN的損失，制備的樣本在完成光譜數(shù)據(jù)采集后立即按照GB/T 18979—2003《食品中黃曲霉毒素B1的測定》[23]和GB/T 5009.209—2008《食品中玉米赤霉烯酮的測定》[24]分別測定AFB1和ZEN含量，并將其測定結(jié)果作為建模的標(biāo)準(zhǔn)參考值。

1.3.3 高光譜圖像數(shù)據(jù)采集與校正

在進(jìn)行高光譜圖像采集時(shí)，為了滿足實(shí)際玉米霉變檢測的需要，將玉米樣本散亂均勻平鋪在規(guī)格為10 cm×10 cm的白色載物盒中，然后將裝有玉米的載物盒放置在移動載物平臺上，依次對玉米樣本進(jìn)行高光譜圖像采集，設(shè)定高光譜的掃描圖像大小為550×600像素。本實(shí)驗(yàn)共采集到250 幅玉米樣本高光譜圖像。

同時(shí)，為了消除光源強(qiáng)度分布不均以及光譜暗電流噪聲的影響，對每個(gè)樣本按校正公式進(jìn)行黑白校正，計(jì)算公式如下：

式中：R為校正后的樣本高光譜圖像；R0為樣本原始高光譜圖像；B為全黑標(biāo)定圖像；W為全白標(biāo)定圖像[25]。樣本的高光譜圖像采集在SICap-STVR V1.0.x軟件平臺上完成，實(shí)驗(yàn)數(shù)據(jù)分析與處理在MATLAB 2014a（美國The Math Works公司）和遙感圖像處理平臺ENVI5.1（美國Boulder公司）兩個(gè)軟件上完成。

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

1.4.1 原始光譜數(shù)據(jù)進(jìn)行預(yù)處理

為了減少背景噪聲、雜散光等無用信息對原始光譜數(shù)據(jù)的干擾，需要對原始光譜數(shù)據(jù)進(jìn)行預(yù)處理。本實(shí)驗(yàn)對比標(biāo)準(zhǔn)化、變量標(biāo)準(zhǔn)化、多元散射校正（multiplicative scatter correction，MSC）、標(biāo)準(zhǔn)正態(tài)變量校正（standard normal variate，SNV）和卷積平滑（Savitzky-Golay smoothing，S-G）5 種預(yù)處理方法。為保證模型校正集最大程度地表征樣本均勻分布，提高模型穩(wěn)定性，本實(shí)驗(yàn)采用光譜-理化值共生距離（sample set partitioning based on joint x-y distance，SPXY）算法[26]劃分模型的初始校正集和預(yù)測集。在此基礎(chǔ)上，為降低或者消除校正集樣本間的共線性[27]，簡化模型運(yùn)算量，采用SPXY算法結(jié)合PLSR法分析不同校正集樣本子集預(yù)測AFB1和ZEN含量的差異，從而進(jìn)一步對劃分的初始校正集樣本進(jìn)行優(yōu)選。在采用均勻光譜間隔（uniform spectral spacing，USS）法對原始光譜變量進(jìn)行初步篩選的基礎(chǔ)上對比連續(xù)投影算法（successive projections algorithm，SPA）[28]和競爭性自適應(yīng)重加權(quán)算（competitive adaptive reweighted sampling，CARS）法[29]2 種特征波長提取方法，最大程度地剔除原始光譜矩陣中的冗余信息，提高模型精度，減少模型運(yùn)算量。

1.4.2 模型建立與評價(jià)

霉變玉米的光譜信息與其AFB1和ZEN含量之間的關(guān)系屬于非確定性問題，采用回歸分析構(gòu)造變量間的數(shù)理統(tǒng)計(jì)模型可用于該類問題的研究分析[30]。PLSR法集主成分分析、普通多元線性回歸和相關(guān)分析于一體，它在描述光譜矩陣X變量的同時(shí)也描述了指標(biāo)矩陣Y變量，較好地解決了自變量的多重共線性問題，在光譜分析領(lǐng)域得到了廣泛應(yīng)用。

2 結(jié)果與分析

2.1 不同霉變等級的玉米AFB1和ZEN含量

在進(jìn)行AFB1和ZEN含量檢測時(shí)，對每類玉米樣品的50 個(gè)平行樣本進(jìn)行測定，以獨(dú)立測定結(jié)果的絕對值差不超過算術(shù)平均值的10%為準(zhǔn)，并將其平均值作為此類玉米樣品的實(shí)際指標(biāo)值，結(jié)果如表1所示。由表1可知，第4、5個(gè)等級的玉米中，這2 類毒素含量已超過國家標(biāo)準(zhǔn)（AFB1≤20 μg/kg，ZEN≤60 μg/kg）。

表1 5 個(gè)霉變等級玉米樣品中AFB1和ZEN含量Table 1 AFB1 and ZEN values in 5 grades of moldy maize samples

2.2 光譜反射值曲線

提取每個(gè)玉米樣本在371.05～1 023.83 nm波長范圍內(nèi)的高光譜圖像的平均光譜反射值，結(jié)果如圖2所示，371.05～480.55 nm和999.66～1 023.83 nm兩個(gè)波段受噪聲的影響較大，因此剔除了這2 個(gè)波段，最終采用的波段范圍為481.06～999.15 nm。由圖3可以看出，5 個(gè)等級玉米樣品光譜曲線變化趨勢基本相似，但不同霉變等級玉米的光譜值有差別，總體看來，5 個(gè)等級玉米樣本具有一定的可分性，這為預(yù)測建模提供了可能。

圖2 250 個(gè)玉米樣本的光譜反射值曲線Fig. 2 Spectral reflectance curves of 250 maize samples

圖3 5 類玉米樣本的平均光譜反射值曲線Fig. 3 Average spectral reflectance curves of five grades of maize samples

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

隨機(jī)選取200 個(gè)樣本作為校正集，剩余的50 個(gè)樣本作為預(yù)測集。采用不同預(yù)處理方法對原始光譜數(shù)據(jù)進(jìn)行預(yù)處理，并基于各預(yù)處理后的光譜數(shù)據(jù)建立PLSR模型，結(jié)果如表2所示。由表2可知，與基于原始光譜數(shù)據(jù)的預(yù)測結(jié)果相比，5 種預(yù)處理方法中，基于SNV建立的PLSR模型對這2 種毒素含量的預(yù)測效果最好，對應(yīng)AFB1和ZEN含量的預(yù)測集相關(guān)系數(shù)和均方根誤差（R2pre，RMSEP）分別為（0.994 4，0.984 6）和（0.991 6，2.320 9），因此確定SNV為預(yù)處理方法。

表2 基于不同預(yù)處理方法的PLSR建模結(jié)果Table 2 Results of PLSR models based on different pretreatments

2.4 校正集樣本的優(yōu)選

2.4.1 SPXY算法劃分樣本集數(shù)據(jù)

校正集樣本的劃分在一定程度上能夠決定所建模型的預(yù)測性能，本實(shí)驗(yàn)采用SPXY算法對樣本集進(jìn)行劃分。首先設(shè)定初始校正集樣本個(gè)數(shù)為200 個(gè)，剩余的50 個(gè)樣本為預(yù)測集樣本，SPXY算法的詳細(xì)步驟參考文獻(xiàn)[26]，劃分結(jié)果如表3所示，不同霉變等級的樣品間，校正集和預(yù)測集的樣本個(gè)數(shù)存在較明顯的差別。

表3 SPXY劃分樣本集結(jié)果Table 3 Calibration and prediction set partitioned by SPXY

2.4.2 SPXY算法優(yōu)選校正集樣本

針對表3中初分的初始校正集樣本，采用SPXY算法進(jìn)行進(jìn)一步優(yōu)選，以盡量有效地降低樣本間的共線性。確定樣本數(shù)N的范圍為80～200，步長為10，校正集樣本優(yōu)選過程中，分別基于優(yōu)選出的校正集樣本子集建立PLSR模型，其預(yù)測集的相關(guān)系數(shù)R2pre和RMSEP的變化如圖4所示。

圖4 R2pre和RMSEP隨校正集樣本數(shù)量變化曲線Fig. 4 Plots of R2pre and RMSEP values versus number of calibration set samples

由圖4可以看出，R2pre隨樣本數(shù)N的增加呈遞增趨勢，但總體變化不明顯（0.995～0.999），RMSEP曲線呈遞減趨勢，總體變化較明顯。圖4a中，對于RMSEP變化曲線，N在80～130范圍內(nèi)取值時(shí)，RMSEP變化差異明顯，N在130～200范圍內(nèi)取值時(shí)，RMSEP的變化趨于平緩，N取值為130 為該曲線的拐點(diǎn)，該點(diǎn)對應(yīng)的R2pre和RMSEP為0.997 4和0.672 0，當(dāng)N取值為190時(shí)，R2pre達(dá)到極大值0.997 6，RMSEP降為極小值0.641 4，與N取130 時(shí)對應(yīng)的R2pre和RMSEP相比，R2pre僅增加0.000 2，RMSEP僅減少0.030 6，從數(shù)值上看，兩者差異較小，但校正集樣本數(shù)增加了60 個(gè)，因此，為簡化模型運(yùn)算量，綜合圖4a中R2pre和RMSEP曲線的變化趨勢，對于AFB1，最終從初始校正集中優(yōu)選出130 個(gè)樣本組成模型校正集。同樣根據(jù)圖4b，對于ZEN，最終從初始校正集中優(yōu)選出140 個(gè)樣本組成模型校正集，對應(yīng)的R2pre和RMSEP為0.998 7和0.862 1。

經(jīng)SPXY算法劃分的預(yù)測集樣本和初始校正集樣本以及經(jīng)SPXY算法優(yōu)選的校正集樣本的AFB1和ZEN含量變化情況如表4所示。由表4可知，無論是對AFB1還是對ZEN，SPXY法優(yōu)選出的校正集樣本的各指標(biāo)值含量變化范圍在初始校正集范圍內(nèi)，且標(biāo)準(zhǔn)差相近，證明SPXY算法優(yōu)選后的樣本具有一定的代表性。

表4 不同樣本集玉米AFB1和ZEN含量分布Table 4 Distribution profiles of AFB1 and ZEN values in different sample sets

2.5 特征光譜的選擇

2.5.1 光譜數(shù)據(jù)的降維

高光譜圖像具有較高的光譜分辨率，容易對被測物進(jìn)行分辨，但也會導(dǎo)致數(shù)據(jù)量的劇增和數(shù)據(jù)冗余。相鄰波段間差異很小，直接進(jìn)行特征波段的提取可能會漏掉某些有用信息，因此，本實(shí)驗(yàn)首先對光譜變量進(jìn)行初降維。USS法通過控制步長選出間隔均勻、相關(guān)性低、信息量大的少量波段，從而可對光譜數(shù)據(jù)進(jìn)行有效降維。在481.06～999.15 nm范圍內(nèi)，分別以正整數(shù)2、3、4、5、6……20為間隔，從1 023 個(gè)經(jīng)預(yù)處理后的光譜波段中抽取光譜組成新的光譜矩陣并建立PLSR模型，根據(jù)模型的和RMSEP，確定最佳的步長間隔。表5為利用不同步長間隔下抽取的光譜建立的PLSR模型對AFB1和ZEN含量的預(yù)測結(jié)果。

從表5可以看出，當(dāng)間隔為8，波段為128 個(gè)時(shí)，模型對AFB1含量的預(yù)測效果較好，R2pre和RMSEP為0.997 4和0.673 4；當(dāng)間隔為7，波段為147 個(gè)時(shí)，模型對ZEN含量的預(yù)測效果較好，R2pre和RMSEP為0.998 8和0.835 4。對比只采用SPXY算法優(yōu)選校正集樣本而未對光譜變量進(jìn)行篩選所建的PLSR模型結(jié)果可知，通過USS法進(jìn)行變量篩選后所建模型對AFB1含量的預(yù)測精度變化不大，對ZEN含量的預(yù)測精度有一定的提高，但參與建模的變量數(shù)能顯著減少。

表5 基于不同間隔光譜數(shù)據(jù)建立PLSR模型預(yù)測結(jié)果Table 5 Predictive results of PLSR models with different numbers of spectral intervals

2.5.2 光譜特征波段的選取

圖5 基于SPA的變量篩選Fig. 5 Wavelengths selected by SPA

為進(jìn)一步實(shí)現(xiàn)霉變玉米在線檢測的檢測時(shí)間短、準(zhǔn)確率高的要求，本實(shí)驗(yàn)在USS法篩選波長變量的基礎(chǔ)上，對比了SPA和CARS 2 種特征波長選擇方法。

在SPA中，最小的交互驗(yàn)證RMSE對應(yīng)的波長變量個(gè)數(shù)即為最終的選擇結(jié)果，指定波段數(shù)N的范圍為2～25，SPA的具體運(yùn)算步驟可參考文獻(xiàn)[28]。由圖5a可以看出，從128 個(gè)波長變量中提取17 個(gè)波長時(shí)RMSE達(dá)到最小；由圖5b可以看出，從147 個(gè)波長變量中提取17 個(gè)波長時(shí)RMSE達(dá)到最小。表6為基于SPA提取的特征波長建立的AFB1和ZEN含量的PLSR模型預(yù)測結(jié)果。

在CARS算法中，設(shè)定蒙特卡羅采樣次數(shù)為50，CARS算法的具體運(yùn)算步驟參考文獻(xiàn)[29]，如圖6所示，基于CARS算法針對AFB1和ZEN含量的變量篩選過程，對這2 種毒素指標(biāo)值，當(dāng)采樣次數(shù)分別為14 次和19 次時(shí)，對應(yīng)的特征波長個(gè)數(shù)分別為18 個(gè)和19 個(gè)。

圖6 基于CARS的變量篩選Fig. 6 Wavelengths selected by CARS

表6 基于SPA特征波長提取方法的PLSR模型結(jié)果Table 6 Parameters of PLSR models with characteristic wavelengths extracted by SPA

表7 基于CARS特征波長提取方法的PLSR模型結(jié)果Table 7 Parameters of PLSR models with characteristic wavelengths extracted by CARS

對比表6和表7可知，基于SPA篩選出的波長所建立的PLSR模型對AFB1和ZEN含量的預(yù)測效果較CARS算法好。

2.6 模型分析

原始光譜數(shù)據(jù)經(jīng)過SNV預(yù)處理后，采用SPXY算法對校正集樣本進(jìn)行劃分與優(yōu)選，運(yùn)用USS法結(jié)合SPA對光譜數(shù)據(jù)進(jìn)行篩選，然后建立基于優(yōu)選后的校正集樣本及特征波長的PLSR模型。如圖7所示，波長優(yōu)選后，PLSR模型的預(yù)測能力依舊表現(xiàn)良好，AFB1含量的預(yù)測精度（，RMSEP）為（0.997 3，0.681 5），與優(yōu)選之前的預(yù)測精度（0.997 4，0.672 0）相近；ZEN含量的預(yù)測精度（，RMSEP）為（0.997 7，1.144 1），略低于優(yōu)選前的精度（0.998 7，0.862 1）。

圖7 基于特征波段建立PLSR模型預(yù)測結(jié)果Fig. 7 Predictive results of PLSR models based on characteristic wavelengths

3 結(jié) 論

本實(shí)驗(yàn)利用不同霉變時(shí)期的玉米光譜實(shí)驗(yàn)數(shù)據(jù)，研究了高光譜技術(shù)用于霉變玉米中AFB1和ZEN含量檢測模型的構(gòu)建，并得到以下結(jié)論：

1）對于AFB1和ZEN含量，經(jīng)SPXY法分別優(yōu)選出的130 個(gè)和140 個(gè)校正集樣本的指標(biāo)值變化范圍在初始校正集范圍內(nèi)，且標(biāo)準(zhǔn)差相近，說明SPXY算法優(yōu)選的樣本具有一定的代表性。

2）基于SPXY法劃分的初始校正集樣本建立的PLSR模型預(yù)測精度明顯高于劃分之前的預(yù)測精度，且通過SPXY法對初始校正集樣本進(jìn)一步優(yōu)選，PLSR模型精度基本不變，但校正集樣本數(shù)減少為初始校正集的65%和70%，在保證模型穩(wěn)健性的前提下，有效降低了校正模型的復(fù)雜性。

3）利用USS法結(jié)合SPA進(jìn)行特征波段選擇，最終從1 023 個(gè)波段中分別篩選出17 個(gè)波段。分別建立2 種毒素值基于各自特征波段的PLSR模型，結(jié)果顯示AFB1含量預(yù)測精度（0.997 3，0.681 5）與優(yōu)選前的預(yù)測精度（0.997 4，0.672 0）相差不大；ZEN含量預(yù)測精度（0.997 7，1.144 1）略低于優(yōu)選之前的預(yù)測精度（0.998 7，0.862 1），R2pre從0.998 7降低到0.997 7，RMSEP從0.862 1上升到1.144 1，但波長減少為17 個(gè)，在保證模型精度的前提下，實(shí)現(xiàn)光譜數(shù)據(jù)的壓縮。

綜合上述研究結(jié)果，基于SPXY算法和SPA建立的高光譜檢測模型，可以實(shí)現(xiàn)對霉變玉米中AFB1和ZEN含量的準(zhǔn)確預(yù)測。