0 前言

　　化橘紅又名化皮、化州橘紅，為蕓香科植物化州柚的未成熟果實的外層果皮。前者習稱“毛橘紅”，后者習稱“光七爪”、“光五爪”。化橘紅不僅具有治咳化啖、健胃行氣、醒酒功能，而還是人體美容的原料，有廣闊的市場前景。研究表明，揮發油，黃酮類化合物，多糖以及香豆素類化合物等是化橘紅的主要有效成分。吳宋夏證實了化橘紅中的黃酮苷具有鎮痛、抗炎、解熱的作用；香豆素類化合物具有抗氧化、抗菌抗病毒之功效。不同的品種有效成分的含量不一樣，功效不一樣，并且在價格上也相差較大，以正品皮的效果佳。因此市場上存在許多用化橘紅的正品果、偽品果、偽品皮冒充正品皮，損害了消費者利益，也沖擊了種植優良品種的農民們的利益。

　　目前對正品皮常用的鑒別方法主要有性狀鑒定、顯微鑒定、高效液相色譜法。這些方法雖然各有優勢，但是存在不同程度上的主觀性強、需要預處理、實驗過程復雜等缺點，不能滿足市場快速、可靠檢測的需要。本研究分別利用400-1000nm、1000-2500nm兩款成像高光譜相機獲取正品皮、偽品皮、正品果、偽品果四種化橘紅成分的高光譜信息，利用光譜指數(DVI、NDVI)、偏最小二乘判別分析(PLS-DA)和極限學習法(ELM)三種方法分別構建四種不同成分的鑒別模型，并用獨立樣本數據對不同的模型進行驗證。

　　1 材料與方法

　　1.1 材料

　　試驗用的化橘紅四種不同成分正品皮、偽品皮、正品果、偽品果由廣東省食品藥品職業技術學校提供，其中正品皮樣本32個、正品果10個、偽品果11個，偽品皮7個。樣本經粉碎均勻后，各取5g放置于培養皿上，標號，用于高光譜相機的光譜采集。

　　1.2 高光譜圖像采集

　　高光譜圖像數據采集采用江蘇雙利合譜科技有限公司的 GaiaSorter高光譜分選儀系統（V10E、N25E-SWIR）。該系統主要由高光譜成像儀、面陣列相機、鹵素燈光源、暗箱、計算機組成，如圖1。實驗儀器參數設置如表1。

圖1 GaiaSorter 高光譜分選儀

表1   GaiaSorter 高光譜分選儀系統參數

　　在進行高光譜圖像采集時，需要設置相機曝光時間，平臺移動速度以及物鏡之間的距離。這 3 個參數相互影響，圖像調節的目的是使采集的圖像大小合適，清晰，不變形失真。經過反復嘗試，物鏡高度設置為 31 cm，曝光時間設置為10ms，平臺移動速度分別設置為 6.0 mm/s（400-1000 nm）、16mm/s（1000-2500 nm）。圖像采集軟件采用江蘇雙利合譜科技有限公司提供的高光譜成像系統采集軟件完成。圖像處理采用 ENVI5.3 軟件進行處理。在進行圖像處理之前，先要對采集的光譜圖像進行圖像校正，圖像校正公式如下：

  (1)

　　式中，R_ref 是校正過的圖像，DN_raw 是原始圖像，DN_white為白板校正圖像，DN_dark 是黑板校正圖像。高光譜圖像的光譜與圖像之間有著對應的關系，據此，在正品皮、偽品皮、正品果、偽品果四種樣本上選取 25×25 像素點的感興趣區域，以感興趣區域所有像素的光譜平均值作為該樣本的平均光譜。

　　1.3  光譜噪聲去除

　　試驗得到光譜含有由儀器和試驗條件等引起的噪聲，對這些噪聲的處理有助于減少噪聲對光譜分析的影響，突出光譜的有效信息。SG 平滑算法可以有效消減光譜數據中的隨機噪聲，消噪效果受平滑點數的影響，本文中選擇SG 二次多項式 7 點平滑對光譜數據進行處理（何勇，2013）。

　　1.4  特征波長選擇

　　光譜信息之間存在大量的冗余和共線性信息特征，對光譜有效信息的提取產生了較大的干擾，且大量光譜數據造成模型復雜、計算量大的問題。本文采用連續投影算法（successive projections algorithm，SPA）進行特征波長的選擇，以減少信息冗余和共線性的影響，簡化模型，減少計算量。

　　SPA 是一種特征變量前向選擇算法，在光譜特征波長中取得了廣泛的應用。本文采用 SPA 算法對去噪處理后的光譜進行特征波長選擇。進行SPA 計算時，以建模集的光譜數據和類別賦值為輸入，設置選擇特征波長數的范圍為 5～30。

　　1.5 光譜指數

　　光譜指數的產生來源于植被指數，植被指數是指利用衛星不同波段探測數據組合而成的，能反映植物生長狀況的指數。植物葉面在可見光紅光波段有很強的吸收特性，在近紅外波段有很強的反射特性，這是植被遙感監測利用衛星不同波段探測數據組合而成的，能反映植物生長狀況的指數。植物葉面在可見光紅光波段有很強的吸收特性，在近紅外波段有很強的反射特性，這是植被遙感監測的物理基礎，通過這兩個波段測值的不同組合可得到不同的植被指數。光譜指數是通過任意兩波段組合或三波段組合成各種光譜指數，如NDVI、DVI等，探尋波段組合用于各個領域的模型構建等。

歸一化植被指數(Normalized difference vegetable index)

(2)

差值植被指數(Difference vegetable index)

(3)

其中，λ1和λ2代表任意波長的反射率，波段范圍為400-1000 nm與1000-2500 nm。

　　1.6 判別分析方法

　　偏最小二乘法判別分析(PLS-DA ,Partial least squares discrimination analysis)是一種用于判別分析的多變量統計分析方法。判別分析是一種根據觀察或測量到的若干變量值，來判斷研究對象如何分類的常用統計分析方法。其原理是對不同處理樣本（如觀測樣本、對照樣本）的特性分別進行訓練，產生訓練集，并檢驗訓練集的可信度（Luna et al., 2013）。本文分別基于全光譜、特別波長光譜建立 PLS-DA 判別分析模型，通過建立光譜數據與類別特征之間的回歸模型，進行判別分析。

　　1.7 極限學習機

　　極限學習機（extreme learning machine）ELM是一種簡單易用、有效的單隱層前饋神經網絡SLFNs學習算法。2004年由南洋理工大學黃廣斌副教授提出。傳統的神經網絡學習算法（如BP算法）需要人為設置大量的網絡訓練參數，并且很容易產生局部優解。極限學習機只需要設置網絡的隱層節點個數，在算法執行過程中不需要調整網絡的輸入權值以及隱元的偏置，并且產生優解，因此具有學習速度快且泛化性能好的優點。本文中隱含層神經元個數從 1 到 520（288）以步長 1 進行尋優，以最小訓練誤差下的神經元個數為 ELM 模型隱含層神經元個數。

　　1.8  評價指標

　　回歸模型得到的樣本的預測值不是整數，需要設置閾值以判斷樣本的歸屬。本文中閾值設置為 0.5 ，預測值小數點大于或等于0.5則加1歸整，小于0.5則減1歸整。總體識別精度是指正確識別的個數除以總數，正品皮識別精度是指正品皮正確識別的個數除以正品皮的總數，正品皮識別錯誤率指數被錯誤分為正品皮的個數除以正品皮的總數。

　　2  結果與分析

　　2.1  化橘紅不同成分的原始光譜曲線

　　本試驗采用的V10E 相機獲取的是400-1000 nm波長范圍共520個波段的可見/近紅外光譜數據，N25E-SWIR相機獲取的是1000-2500 nm波長范圍共288個波段的近紅外光譜數據，全部樣本的原始光譜圖如圖 2所示，正品皮、偽品皮、正品果、偽品果的光譜比較圖如圖3所示。

　　從圖2和圖3可以看到，總體而言，無論是400-1000 nm或1000-2500 nm波長范圍內，正品皮的光譜反射率值低于其他三種成分的光譜曲線，從曲線變化趨勢來看四種不同成分并沒有十分明顯的差異。本研究按照Kennard-Stone 算法將樣本分成建模集和預測集，其中建模集 38 個樣本，預測集32個樣本。正品皮、偽品皮、正品果、偽品果分別賦值為 1、2、3、4（表 2），不同化橘紅成分建模集和預測集樣本的劃分如表 2所示。

圖2 全部化橘紅樣本的原始反射光譜圖

圖3 化橘紅不同成分反射光譜曲線圖

表2  化橘紅不同成分類別賦值與建模集合檢驗集樣本劃分

　　2.2  主成分分析定性分析

　　對化橘紅的四種不同成分的光譜數據進行主成分分析（principal component analysis，PCA） 可知，400-1000 nm范圍內，第一主成分（principle component 1, PC1）的貢獻率為88.36%，PC2 的貢獻率為7.24%，PC1 和 PC2 累計貢獻率為95.6%，PC1 和 PC2 能夠解釋絕大部分的變量；在1000-2500 nm范圍內，第一主成分（principle component 1, PC1）的貢獻率為93.27%，PC2 的貢獻率為3.72%，PC1 和 PC2 累計貢獻率為97.0%，PC1 和 PC2 也能夠解釋絕大部分的變量；圖4分別為400-1000 nm與1000-2500 nm范圍內，第一主成分與第二主成分的散點分布圖。

圖4  400-1000nm（左）、1000-2500nm（右）第一主成分與第二主成分的得分散點分布圖

圖5 主成分組合識別正品皮（左：400-1000 nm，右：1000-2500 nm）

　　由圖 4 可知，在400-1000 nm與1000-2500 nm范圍內，各不同成分之間較難直接區分出來，得分圖中均有重合的地方。從圖5亦可以看出無論是400-1000 nm或1000-2500 nm，主成分的成分組合均未能從正品皮、偽品皮、正品果、偽品果四種樣本中識別出正品皮，從圖5可知，識別出的正品皮亦有假的正品皮。因此需要對光譜數據進行進一步的分析和處理，以鑒別化橘紅不同成分。

　　2.3  光譜指數

　　本研究使用的成像高光譜儀V10E范圍是400-1000nm，N25E-SWIR范圍是1000-2500 nm，通過不同波段的重新組合，形成不同形式的光譜指數(歸一化光譜指數，NDVI；差值光譜指數，DSI)。

　　運用Matlab軟件編程，將兩兩組合的所有波段構建的光譜指數與各類別賦值，計算相應的決定系數(R²)，繪制決定系數圖。圖中顏色從藍色到紅色變化，圖像顏色越紅，表示決定系數越大，圖像顏色越藍，表示決定系數越小。圖6為化橘紅各成分分類賦值與NDVI、DVI決定系數（R2）二維圖。從圖6可知，在400-1000 nm范圍內， DVI模型決定系數最高的兩波段在綠光范圍內，分別是538.91 nm和543.75nm；NDVI模型絕對系數最高的兩波段989.76 nm和670.39 nm。在1000-2500 nm范圍內，NDVI和DVI兩個光譜指數模型決定系數最高的兩波段組合均為是1820.72 nm和1787.12 nm。綜合400-1000 nm與1000-2500 nm

圖 6化橘紅各成分分類賦值與NDVI、DVI決定系數（R²）二維圖

　　兩組數據構建NDVI和DVI模型時，NDVI和DVI兩個光譜指數模型決定系數最高的兩波段組合與1000-2500 nm范圍內相同。比較400-1000 nm和1000-2500 nm的決定系數圖可知，在1000-2500 nm范圍內，NDVI與DVI構建的模型決定系數較高。

　　根據圖6可知，在400-1000 nm范圍內，DVI_{(538.91/543.75)}與NDVI_{(989.76/670.39)}與類別賦值構建的模型決定系數最高，分別為0.643和0.640。在1000-2500 nm范圍內，DVI_{(1820.72 /1787.12)}與NDVI_{(1820.72 /1787.12)}與類別賦值構建的模型決定系數最高，分別為0.861和0.834。圖7分別為DVI_{(538.91/543.75)}、NDVI_{(989.76/670.39)} 、DVI_{(1820.72 /1787.12)}、NDVI_{(1820.72 /1787.12)}與類別賦值的線性擬合散點分布圖，從圖中可知，在1000-2500 nm范圍內構建的NDVI、DVI光譜指數與類別賦值擬合度最高，且變化曲線顯著。

 

 

圖7  DVI、NDVI光譜指數與類別賦值的散點分布圖

　　運用獨立的數據，分別對圖7中DVI_{(538.91/543.75)}、NDVI_{(989.76/670.39)} 、DVI_{(1820.72 /1787.12)}、NDVI_{(1820.72 /1787.12)}與類別賦值構建的模型進行檢驗，檢驗結果如圖8和表3所示。根據查表可知，1-15與28-32為正品皮，15-25號為16-20正品果、21-24為偽品果，25-27為偽品皮。從圖8可知，正品皮、正品果、偽品果、偽品皮之間均存在不同程度的錯誤識別。表3基于DVI、NDVI光譜指數模型檢驗化橘紅樣本的精度評價表。表中分別統計了總體識別精度、正品皮識別精度和正品皮識別錯誤率。從表3中可知，總體識別精度、正品皮識別精度最高的均為DVI_{(1820.72 /1787.12)}構建的模型，分別是66%和75%，正品皮識別錯誤率低的則為DVI_{(1820.72 /1787.12)}、NDVI_{(1820.72 /1787.12)} 構建的模型，均為5%。

  

圖 8 基于NDVI、DVI檢驗集預測結果對比圖

表2 基于光譜指數模型檢驗化橘紅樣本的精度評價

	總體識別精度%	正品皮識別精度%	正品皮識別錯誤率%
400-1000 (DVI)	53	65	15
400-1000 (NDVI)	50	70	20
1000-2500 (DVI)	66	75	5
1000-2500 (NDVI)	59	65	5

 

　　2.4  特征波長選擇

　　本文以建模集樣本的光譜數據和類別賦值為輸入，利用連續投影算法SPA選擇特征波長。選出的特征波長的個數如表 3所示。從表 3 可知，400-1000 nm范圍內所選擇的特征波段為15個，1000-2500 nm范圍內所選擇的特征波段較少，為5個。

表 3在400-1000 nm與1000-2500 nm范圍內SPA 算法選擇的特征波長個數

范圍	波段位置/nm
400-1000 nm	395.23, 396.40, 399.90, 401.06, 422.13, 501.59, 670.39, 735.79, 743.39, 872.13, 940.31, 962.16, 974.58, 981.45, 992.73
1000-2500 nm	1461.59, 1714.23, 2038.97, 2329.25, 2574.24

 

　　2.5  偏最小二乘判別分析

　　分別將建模集合檢驗集的全波段光譜及特征波段光譜作為偏最小二乘法判別分析輸入變量，從而獲取預測樣本類別的賦值。圖9為檢驗集與預測值的類別賦值圖，并針對圖9進行總體識別精度、正品皮識別精度、正品皮識別錯誤率的統計，如表4所示。從表4中可知，總體識別精度、正品皮識別精度最高的均為在1000-2500 nm范圍內的PLS-DA構建的模型，分別是78%和90%，正品皮識別錯誤率低的則為1000-2500 nm范圍內的PLS-DA、PLS-DA 構建的模型，均為5%。從表4可以看出，無論是400-1000 nm或1000-2500 nm，基于全波段的PLS-DA模型總體識別率和正品皮識別率均高于基于特征波段的PLS-DA模型，而正品皮的錯誤識別率，無論是400-1000 nm或1000-2500 nm范圍內，基于全波段的PLS-DA模型與基于特征波長的PLS-DA模型的錯誤識別率相同，400- 1000 nm范圍錯誤識別率均為10%，1000-2500范圍錯誤識別率則為5%。

圖 9 基于PLS-DA檢驗集預測結果對比圖

表4基于PLS-DA模型檢驗化橘紅樣本的精度評價

	總體識別精度%	正品皮識別精度%	正品皮識別錯誤率%
400-1000 (PLS-DA)	72	65	10
400-1000 (PLS-DA-SPA)	63	55	10
1000-2500 (PLS-DA)	78	90	5
1000-2500 (PLS-DA-SPA)	72	70	5

 

　　2.6  極限學習機

　　分別將建模集合檢驗集的全波段光譜及特征波長光譜作為極限學習機的輸入變量，從而獲取預測樣本類別的賦值。圖10為檢驗集與預測值的類別賦值圖，并針對圖10進行總體識別精度、正品皮識別精度、  正品皮識別錯誤率的統計，如表5所示。從表5中可知，總體識別精度、正品皮識別精度高的均為在1000-2500 nm范圍內的ELM與ELM-SPA構建的模型，分別是84%和95%，正品皮識別錯誤率低的則為1000-2500 nm范圍內的ELM與ELM-SPA 構建的模型，均為5%。從表5可知，在400-1000 nm范圍內，基于特征波段光譜的ELM模型總體識別率與正品皮識別率均高于基于全波段的ELM模型，對于正品皮的識別錯誤率，基于全波段與基于特征波段的錯誤識別率相同；在1000-2500 nm范圍內，無論是基于全波段光譜或基于特征波段的ELM模型，其總體識別率、正品皮識別率、正品皮識別錯誤率均相同，分別為84%、95%和5%。

 

  

圖 10 基于ELM檢驗集預測結果對比圖

表4 基于ELM模型檢驗化橘紅樣本的精度評價

	總體識別精度%	正品皮識別精度%	正品皮識別錯誤率%
400-1000 (ELM)	75	75	10
400-1000 (ELM-SPA)	78	85	10
1000-2500 (ELM)	84	95	5
1000-2500 (ELM-SPA)	84	95	5

 

　　2.7  光譜指數模型，PLS-DA 模型，和 ELM 模型的比較

　　綜合對比光譜指數模型，PLS-DA 模型，和 ELM 模型的識別效果可知，無論是光譜指數模型，PLS-DA 模型或ELM 模型，基于1000-2500 nm范圍內構建的模型，其預測值的總體識別率、正品皮識別率均高于400-1000 nm范圍內的模型，且正品皮的識別錯誤率也低于400-1000nm范圍內的模型。在光譜指數模型、PLS-DA 模型和 ELM 模型的模型中，ELM 模型的識別準確性最高，其次是PLS-DA模型，最后是光譜指數模型。基于特征波段光譜的PLS-DA模型其識別準確性低于基于全波段光譜的PLS-DA的模型，但是基于特征波段光譜的ELM模型在400-1000 nm范圍內，其識別準確性高于基于全波段光譜的ELM模型，在1000-25000 nm范圍內，其識別準確性與基于全波段光譜的ELM模型相同。

　　3  結論與討論

　　基于V10E與N25E-SWIR兩款成像高光譜相機，分別獲取正品皮、正品果、偽品皮、偽品果四種化橘紅成分400-1000 nm與1000-2500 nm范圍的光譜反射率，采用 SG 平滑算法對提取出的光譜數據進行去噪處理，同時采用 SPA 算法對去噪后的光譜提取特征波長，并分別基于全波段光譜、特征波段光譜建立 PLS-DA 判別模型和 ELM 模型，同時采用全波段循環，探尋NDVI、DVI兩個光譜指數構建判別模型，用于鑒別正品皮、正品果、偽品皮、偽品果，取得了比較好的識別效果。基于特征波段光譜與全波段光譜建立的 ELM 模型取得了好的效果，總體識別精度、正品皮識別精度、正品皮識別錯誤率分別為84%、95%和5%。在實際運用中，考慮到識別時間與成分，基于SPA算法提取的特征波段構建的ELM模型效果佳。本論文研究結果為高光譜成像技術在藥品真偽等檢測中的應用提供了可行性。