基于變量優選和近紅外光譜技術的紅富士蘋果產地溯源

張立欣，楊翠芳，陳杰，張曉果，張楠楠，張曉*

1(塔里木大學 信息工程學院，新疆維吾爾自治區 阿拉爾，843300)2(南京理工大學 理學院，江蘇 南京，210094) 3(河南城建學院 數理學院，河南 平頂山，467036)

中醫認為蘋果具有生津止渴、潤肺除煩、健脾益胃、養心益氣等功效，并且它的味道酸甜適口，營養豐富，因此，成為老幼皆宜的水果之一。由于環境因素、土壤特征等的不同，各產地蘋果的品質存在差異，其口感也不盡相同，而這些差異無法通過肉眼直接準確地辨別，化學鑒定法費時費力，不適合進行大批量操作。

近紅外光譜檢測以其方便、高效、快速、無污染、無需對樣品預處理等優點，被廣泛應用到現代農業檢測分析中。主要應用領域包括：農產品成熟度的鑒別[1]、復溶果品損傷的鑒別[2]、復溶同一農產品不同基因型的鑒別[3]、復溶農產品產地的鑒別[4]、農產品新鮮度[5]、貨架期[6]、霉變程度[7]、摻假的檢測[8]等。在光譜分析中，經常會受到背景等隨機因素的干擾，因此需對光譜數據進行預處理[9]，常用的光譜預處理方法有：一階導數(first derivative，1-DER)[10]、二階導數(second derivative，2-DER)[11]、標準正態變換(standard normal transformation，SNV)[12]，多元散射校正(multivariate scatter correction，MSC)[13]、平滑變換(smooth transformation，SG)[14]、標準化[15]、歸一化(normalization，NOR)[16]、中心化(centralization,CEN)[17]等。為降低模型的復雜度，減少共線性的干擾，需要提取特征波段[18]，常用的方法有：連續投影算法(successive projection algorithm，SPA)[19]、競爭性自適應重加權算法(competitive adaptive reweighted sampling，CARS)[2]、主成分分析(principal component analysis，PCA)[20]、隨機蛙跳算法(random frog，RF)[21]等，也有將幾種方法聯合起來選擇特征變量的。程介虹等[22]提出一種改進聯合區間的RF選擇特征波長，通過聯合區間偏最小二乘法對全譜進行變量初選，此時得到的波長對目標變量變化最為敏感，將其作為RF的初始變量子集，以解決其運行時間較長、效率較低的問題。袁凱等[23]采用3步混合策略，提出了間隔偏最小二乘、區間變量迭代空間收縮法和迭代保留信息變量聯用的特征變量選擇方法，對生鮮雞胸肉的近紅外光譜進行特征波長選擇，建立了雞肉水分R檢測模型。結果表明，建模波長數量經3步選擇后減少為全光譜建模的0.76%，但模型精確度和穩定性逐步提高。FANG等[2]將SPA、CARS、過濾式特征選擇。3種方法選取的特征變量組合起來建模，取得了很好的預測效果。

研究表明，光譜技術在農產品的檢測分析中具有廣泛的應用，但是，對蘋果產地溯源的研究相對較少。目前見刊的有馬永杰等[24]采用多種數據降維方法，建立k近鄰(k-nearest neighbor, KNN)模型，對紅富士蘋果進行產地溯源，建模集和預測集的正確率分別達到97.3%和92.3%，模型的正確率有待進一步提高。本研究在前人研究的基礎上，基于近紅外光譜技術，以新疆阿克蘇、甘肅靜寧、河南靈寶和山東煙臺的紅富士蘋果為研究對象，利用光譜分析的理論和方法，借助于光譜預處理算法，特征波長選擇方法，建立概率神經網絡模型對蘋果的產地進行判別分析，重點研究不同的光譜預處理方法、不同的變量篩選方法對預測模型的影響，為實現紅富士蘋果的產地溯源提供理論參考。

1 材料與方法

1.1 材料與儀器

以阿克蘇、靜寧、靈寶和煙臺的紅富士蘋果為實驗對象，在實驗中所使用的蘋果均為2020年10月份在各蘋果產地郵寄所得。挑選表面沒有缺陷、直徑范圍為65～85 cm、大小均勻的蘋果樣品256個，去除表面的污垢，放置在冰柜內保存，溫度控制在4 ℃，實驗前分批拿出，待其恢復到室溫(20～25 ℃)環境后開始實驗。

實驗中所用的推掃式高光譜分選系統(Hyperspectral Sorting System)為北京卓立漢光公司生產。光譜測定的范圍為900～1 700 nm(實際可測量到1 750 nm)，光譜分辨率5 nm，光譜采樣點4 nm，選取果身中心前后左右4個方位，提取大小為20像素×15像素，4個面均進行提取，共1 200像素點，選取平均值為該樣本反射率。通過自帶的ENVI5.3軟件提取ROI的光譜值，最后導出為Excel文件。

1.2 光譜數據的采集和校正

為了得到清晰的圖像，在采集光譜數據前需要多次的相機聚焦和移動平臺的速度測試，反復嘗試后，確定平臺的移動速度為0.35 cm/s，相機曝光時間為0.09 s。為了減少光照不均勻和暗電流對實驗的影響，需要對采集到的光譜數據進行黑白校正，校正公式如公式(1)所示：

R=(I-B)/(W-B)

(1)

式中：R是校正后的光譜數據，I是原始光譜數據，W為對準白板采集到的數據，B是蓋上相機鏡頭采集到的數據。

1.3 數據處理

1.3.1 光譜數據的預處理

在光譜檢測的過程中，會受到樣品背景等隨機因素的影響，導致光譜數據中含有噪聲，為提高模型的準確性和穩健性，需要對數據進行預處理，采用的方法有NOR、CEN、1-DER、2-DER、SNV、MSC、小波變換(wave transformation，WT)、SG、傅里葉變換(Fourier transformation，FT)。

1.3.2 特征波長提取

光譜能夠體現物質所含成分及含量，但同時包含大量的冗余信息，為降低模型的復雜性，減少共線性的影響，因此，需要提取特征波長。采用的方法有PCA、SPA、CARS和RF。

1.3.3 判別分析模型

概率神經網絡(probabilistic neural network，PNN)是徑向基網絡的一個分支，屬于前饋網絡的一種。它具有學習過程簡單、訓練速度快、分類更準確、容錯性好等優點。

1.3.4 模型驗證

將樣本以2∶1的比例間隔分為訓練集和測試集，依靠訓練集建立上述判別分析模型，測試集將通過已經建立好的PNN模型進行驗證，以正確率為標準來評判各類方法的優劣。

2 結果與分析

采集的紅富士蘋果光譜數據中，剔除異常值后，共得到阿克蘇、靜寧、靈寶、煙臺4個產地的蘋果樣本各60個，其原始光譜曲線如圖1所示。

圖1 原始光譜圖Fig.1 Original spectrogram

不同產地的紅富士蘋果樣本的光譜曲線變化趨勢大致相同，反映不同產地的蘋果之間也有著極大的相似性。但是在一些波峰和波谷處存在明顯的偏離，這是由于不同產地的蘋果內部成分含量的多少存在差異。在1 060 nm處也有1個明顯的波峰，這與C—O—C基團有關；在1 440 nm處有1個明顯的波谷，這與H2O的二倍頻吸收帶有關。成分含量上的差異會導致光譜曲線的差異，這為基于光譜技術對蘋果產地的識別研究提供了信息。

阿克蘇、靜寧、靈寶、煙臺4個產地的蘋果樣本均以2∶1的比例間隔產生訓練集和測試集，得訓練集中4個產地的蘋果樣本各為40個，即160個蘋果樣本組成訓練集，測試集中4個產地的蘋果樣本各為20個，即80個蘋果樣本組成測試集。

2.1 光譜數據的預處理

在光譜檢測的過程中，會受到樣品背景等隨機因素的影響，導致光譜數據中含有噪聲，為提高模型的準確性和穩健性，需要對數據進行預處理。分別采用NOR、CEN、1-DER、2-DER、SNV、MSC、WT、SG、FT等9種預處理方法對原始光譜進行預處理，以消除各類隨機因素對模型性能的干擾。分別以沒有經過預處理的原始光譜(no pretreat, NO)和經過預處理之后的光譜數據作為輸入自變量，分別建立線性的偏最小二乘判別分析(partial least squares discrimination analysis，PLSDA)模型和非線性的PNN模型，對蘋果的產地進行判別，其結果如表1和表2所示。

表1 不種預處理方法下PLSDA模型的準確率 單位：%

由表1可知，不同的預處理方式下，PLSDA模型的總準確率為60%～85%，以總準確率最高為標準，最優預處理方法為MSC，此時阿克蘇、靜寧、靈寶、煙臺4個產地的蘋果識別率分別為100%、95%、85%、60%。由表2可知，PNN模型總準確率為86.25%～97.50%，最優預處理方法亦為MSC，此時阿克蘇、靜寧、靈寶、煙臺4個產地的蘋果識別率分別為100%、100%、90%、100%。

表2 不種預處理方法下PNN模型的準確率 單位：%

比較表1和表2可以發現，在相同的預處理方式下，PNN模型的總準確率高于PLSDA模型的總準確率，這是因為蘋果內部的結構復雜，光在蘋果內部的傳輸是一種復雜的結構，擁有非線性判別能力的PNN優于線性的PLSDA。為提高模型的識別率，后續的判別分析模型中，均采用MSC預處理之后的光譜數據，建立PNN模型。

2.2 特征波長的選取

對MSC預處理之后的光譜數據，分別采用PCA、SPA、CARS提取光譜數據中的重要變量作為建模輸入自變量。

2.2.1 PCA選取特征變量

采用PCA算法提取主成分，各主成分解釋的總方差如表3所示。

表3 解釋的總方差Table 3 Eexplained total variance

前3個主成分累計貢獻率達到95%以上，因此，選取前3個主成分，作為下一步判別分析模型的輸入自變量。不同產地的蘋果光譜數據前3個主成分的得分分別如圖2所示。

圖2 主成分圖Fig.2 Principal component diagram

由圖2可知，阿克蘇紅富士蘋果的前3個主成分可以和其他產地的分開，而其他3個產地的蘋果光譜數據之間重疊比較多，這將影響模型的識別效果。

2.2.2 SPA選取特征變量

SPA進行波長變量選擇，指定波長數為1～20，采用均方根誤差最小來確定最終參與建模的波長變量個數，選取過程如圖3所示。

圖3 SPA選取變量過程Fig.3 Variable selection process

由圖3可知，隨著所選的特征波長數的增加，參與建模的波長變量數增加，而均方誤差整體呈現遞減的趨勢。當選取特征波長變量數為14時，均方根誤差為0.318 4，之后均方誤差逐漸趨于平緩，此時，若再多選取變量，過多的波長變量參與建模會增加模型的復雜性，而均方根誤差并不會有太大的變化。因此，最終選取14個波長變量參與建模，選取的特征波長如圖4所示。

圖4 SPA選取的特征變量Fig.4 Selected variable

最終選取的14個特征波長為圖4中的小方塊對應的橫坐標，即對應波長分別為911.06、923.53、929.78、980.02、1 092.04、1 175.88、1 368.46、1 402.29、1 453.03、1 529.06、1 693.72、1 704.37、1 711.48、1 715.03 nm。

2.2.3 CARS選取特征變量

采用CARS方法進行特征波長的選取，蒙特卡羅方法抽樣迭代200次，抽取過程如圖5所示。

a-變量優化過程；b-均方根誤差變化趨勢；c-回歸系數變化圖5 CARS 選取變量結果Fig.5 Variable selection results

從抽取結果來看，到第82次迭代時，均方根誤差達到最小為0.290 1，此時選出35個波長變量，對應的波長為：914.17、920.41、929.78、932.90、967.41、1 127.08、1 133.56、1 136.81、1 146.55、1 375.21、1 378.59、1 405.68、1 409.08、1 412.47、1 415.87、1 419.27、1 456.80、1 484.23、1 487.67、1 491.11、1 494.55、1 497.99、1 501.43、1 504.88、1 508.33、1 515.23、1 522.14、1 529.06、1 581.17、1 584.66、1 588.15、1 591.64、1 697.27、1 700.82、1 736.41 nm。

2.2.4 RF選取特征變量

采用RF算法提取特征波長變量，迭代1 000次，每個波長變量被選中的頻率如圖6所示。

圖6 每個波長變量被選中的概率圖Fig.6 Probability diagram of each wavelength variable being selected

由圖6可知，波長變量被選中的概率范圍為0～0.941 0，優先選擇率較大的波長變量，采用交叉驗證的均方根誤差最小來確定最終選取的變量個數，均方根誤差和所選變量個數的關系如圖7所示。

圖7 RF選取變量過程Fig.7 Process of selecting variables by RF method

由圖7可知，隨著所選變量個數的增加，均方根誤差迅速下降，當所選變量個數為82時，交叉驗證的均方根誤差最小，為0.272 2。所選取的變量為圖6中概率大于0.088 0的波長，即圖6中概率在水平線上方的波長變量，所對應的波長分別為：911.06、914.17、917.23、920.41、923.53、929.78、932.90、945.43、957.98、967.41、970.56、976.86、980.02、989.49、998.98、1 002.15、1 008.48、1 014.83、1 049.87、1 065.87、1 078.70、1 101.23、1 114.14、1 130.32、1 136.81、1 143.30、1 146.55、1 149.80、1 153.06、1 156.31、1 159.57、1 166.09、1 169.35、1 225.09、1 234.98、1 238.28、1 241.59、1 244.89、1 278.04、1 308.03、1 344.89、1 348.25、1 351.61、1 365.09、1 371.83、1 375.21、1 385.35、1 398.90、1 405.68、1 412.47、1 415.87、1 419.27、1 484.23、1 491.11、1 494.55、1 497.99、1 501.43、1 504.88、1 508.33、1 511.78、1 515.23、1 522.14、1 529.06、1 532.52、1 549.85、1 570.71、1 591.64、1 605.63、1 658.34、1 665.41、1 672.47、1 676.01、1 679.55、1 683.09、1 686.63、1 690.17、1 693.72、1 697.27、1 700.82、1 711.48、1 715.03、1 736.41 nm。

由于CARS和RF算法選出的 特征變量的個數仍然比較多，變量之間可能還存在共線性，為簡化模型，將CARS和RF選出的特征變量再進一步優選。

2.2.5 CARS-SPA選取特征變量

對CARS方法選取的特征波長變量，采用SPA算法進一步優選，變量的選取過程如圖8所示。

圖8 CARS-SPA優選變量過程Fig.8 Process of optimizing variables by CARS-SPA method

隨著所選變量個數的增加，均方根誤差呈現遞減的趨勢，當選取變量個數為16時，均方根誤差達到最小0.261 9。此時選取的特征波長變量如圖9中小方塊所對應的橫坐標。

圖9 CARS-SPA優選的變量Fig.9 Preferred variables by CARS-SPA method

最后選出16個特征波長變量，對應的波長為：920.41、929.78、967.41、1 127.08、1 133.56、1 146.55、1 375.21、1 378.59、1 405.68、1 412.47、1 415.87、1 504.88、1 581.17、1 591.64、1 700.82、1 736.41 nm。

2.2.6 RF-SPA選取特征變量

對RF方法選取的特征波長變量，采用SPA算法進一步優選，變量的選取過程如圖10所示。

圖10 RF-SPA優選變量過程Fig.10 Process of optimizing variables by RF-SPA method

隨著所選變量個數的增加，均方根誤差整體下滑，綜合考慮均方根誤差和模型的復雜性，當所選變量為17時，均方根誤差為0.260 9，之后，隨著變量個數的增加，均方根誤差并無明顯減少的趨勢，因此，選取17個特征波長變量。所選取的變量波長如圖11中小方塊對應的橫坐標。

圖11 RF-SPA 優選的變量Fig.11 Preferred variables by RF-SPA method

最后選出的17個特征波長變量對應的波長為：911.06、914.17、920.41、929.78、932.90、980.02、1 114.14、1 169.35、1 371.83、1 405.68、1 453.03、1 515.23、1 591.64、1 658.34、1 693.72、1 700.82、1 736.41 nm。

2.3 建模結果分析

以選出的特征波長變量作為輸入變量，建立PNN判模型，對蘋果的產地進行判別分析。模型識別的準確率如表4所示。

表4 不同的變量選取方式下模型的準確率 單位：%

由表4可知，總準確率范圍為71.25%～100%，PCA方法選取的特征變量建模識別率最低，這是因為PCA在選取過程中，只是對樣品的自變量進行重新組合，沒有考慮因變量的影響。SPA算法可以最大程度地降低被選中變量間的共線性問題，但是對于光譜而言，有效變量之間的投影距離不一定最大，因此，篩選出的變量子集中可能包含一些無用信息，甚至是干擾信息，降低模型的泛化能力。CARS和RF算法選出的特征變量個數比較多，無法克服變量之間的共線性的影響，這都將會影響模型的泛化能力。CARS-SPA和RF-SPA在初選特征變量的基礎上再進一步優選，選出的變量個數適中，消除了變量之間的共線性。其中CARS-SPA選出的特征變量建模，阿克蘇、靜寧、靈寶、煙臺4個產地的紅富士蘋果識別率分別達到了100%、100%、95%、100%，總識別率達到了98.75%。綜合考慮識別率和模型的復雜性，最優模型為MSC-CARS-SPA-PNN模型。

2.4 與KNN模型的比較

采用馬永杰等[24]的K最近鄰算法，分別以PCA、SPA、CARS、RF、CARS-SPA、RF-SPA選出的特征波長變量作為輸入變量，建立KNN模型。k的取值范圍為3～9，通過交叉驗證的方式確定最佳的k值，識別結果如表5所示。

表5 KNN的識別結果 單位:%

比較表4和表5可知，在PCA、SPA、CARS、CARS-SPA變量選取方式下，PNN模型的總準確率均高于KNN模型的總準確率；在RF變量選取方式下，PNN模型和KNN模型有相同的總準確率；在RF-SPA變量選取方式下，PNN模型的總準確率略低于KNN模型的總準確率。但是，在CARS-SPA變量選取方式下，PNN模型的總準確率達到最高98.75%。

3 討論與結論

采集阿克蘇、靜寧、靈寶和煙臺的紅富士蘋果近紅外光譜數據，分別進行歸一化、中心化、一階導數、二階導數、標準正態變換、多元散射校正、小波變換、平滑變換、傅里葉變換共9種預處理方法，建立概率神經網絡模型對蘋果產地進行識別，結果表明，多元散射校正預處理的效果最優，4個產地的蘋果樣品識別率分別為100%、100%、90%、100%，總識別率為97.5%。

為了簡化模型提高識別率，將多元散射校正預處理后的光譜數據分別采用主成分法、連續投影算法、競爭性自適應重加權算法、隨機蛙跳算法、競爭性自適應重加權-連續投影算法、隨機蛙跳-連續投影算法選取特征變量建模。研究結果表明，競爭性自適應重加權-連續投影算法效果最佳，MSC-CARS-SPA-PNN模型的總識別率達到了98.75%，4個產地的紅富士蘋果識別率分別達到了100%、100%、95%、100%。因此，近紅外光譜技術能夠快速、無損地識別蘋果的產地，為紅富士蘋果的產地溯源提供理論參考。