張靈枝， 黃 艷， 于英杰， 林 剛， 孫威江

（1. 福建農(nóng)林大學園藝學院，福建 福州 350002；2. 福建農(nóng)林大學安溪茶學院，福建 泉州 362400；3. 中國茶葉流通協(xié)會，北京 100801；4. 福建融韻通生態(tài)科技有限公司，福建 福州 350025；5. 福建農(nóng)林大學福建省茶產(chǎn)業(yè)工程技術研究中心，福建 福州 350002；6. 福建農(nóng)林大學海峽兩岸特色作物安全生產(chǎn)省部共建協(xié)同創(chuàng)新中心，福建 福州 350002）

茶葉富含茶多酚、氨基酸、生物堿等多種對人體健康有益的成分，是世界上最受歡迎的三大無酒精飲料之一。 根據(jù)加工工藝及發(fā)酵程度不同，茶葉可分為綠茶（未發(fā)酵）、白茶（微發(fā)酵）、黃茶（輕發(fā)酵）、烏龍茶（半發(fā)酵）、紅茶（全發(fā)酵）以及黑茶（后發(fā)酵）[1-3]六大類。 六大茶類識別主要通過專業(yè)人士對茶葉外形、風味進行感官審評，結(jié)果易受審評人員的身體、精神狀況等影響，審評主觀性和經(jīng)驗性強。 國際層面，茶葉評價術語尚未實現(xiàn)規(guī)范化、統(tǒng)一化和標準化， 化學分析的國際標準制定嚴重滯后，使得國外非專業(yè)人士對六大茶類的分辨更加困難[4-6]，嚴重制約著茶葉國際貿(mào)易和流通。 因此，不少學者在化學分析技術、計算機視覺技術、光譜技術等領域開展了多茶類判別研究。 Ning 等以兒茶素、咖啡堿等生化成分含量結(jié)合Fisher 判別分析對六大茶類進行逐步判別， 其最佳模型的識別正確率達到88.30%以上[7]；Wu 等以茶多酚、咖啡堿等生化成分含量結(jié)合逐步線性回歸判別分析，建立綠茶、白茶、烏龍茶以及紅茶的判別函數(shù)，識別正確率達97.80%[8]；Jiang 等基于氨基酸、兒茶素等生化成分含量結(jié)合線性判別分析， 對六大茶類可實現(xiàn)100.00%準確鑒別[9]；Zhang 等基于非靶向代謝組學識別綠茶、 白茶以及黃茶，找出了關鍵分類指標，分類識別正確率可達100.00%[10]。上述研究雖然建立了優(yōu)異的判別函數(shù)和模型，但存在特征成分檢測耗時長、成本高、分析復雜等問題，不易進行快速無損鑒別，難以在國際上實現(xiàn)產(chǎn)業(yè)化應用。 Wu 等通過采集紅茶、綠茶、烏龍茶的外形圖片， 利用機器學習和計算機視覺技術，結(jié)合K 近鄰（K-nearest neighbor，KNN）分類器構建判別模型，識別正確率可達94.7%[11]。Ning 等將高光譜成像技術結(jié)合化學計量學建立的Lib-SVM 模型對綠茶、黃茶、白茶、黑茶、烏龍茶的識別正確率可達98.39%[12]。 由于茶葉的年份差異導致的波動性較大、高光譜儀普及率低等問題，上述兩項技術并未實現(xiàn)產(chǎn)業(yè)化應用。

近紅外光譜作為綠色分析技術[13]，具有高效、便捷、準確性高等優(yōu)點，已在食品[14-16]、煉油[17-18]、藥物[19-20]等領域廣泛應用。本研究中以NIRS 結(jié)合機器學習算法，探索六大茶類鑒定的可行性，包括以下4個步驟：1）獲取六大茶類在3 600~12 500 cm-1波段的近紅外光譜；2） 應用最小最大歸一化（minmax scaler，Minmax）、連續(xù)小波變換（continuous wavelet transform，CWT）、 標準正態(tài)變換 （standard normal variate，SNV）及多元散射校正（multiplicative scatter correction，MSC）4 種預處理算法， 建模分析他們對OS 的去噪和散射校正性能；3） 比較主成分分析（principal component analysis，PCA）、線性判別分析（linear discriminant analysis，LDA） 和連續(xù)投影算法（successive projections algorithm，SPA）3 種方法提取茶葉光譜特征的能力；4）建立基于RF 和SVM 分類器的茶類鑒別模型，實現(xiàn)對六大茶類的快速、無損識別，為近紅外光譜技術在茶類識別的產(chǎn)業(yè)應用上奠定理論基礎和科學依據(jù)。

1 材料與方法

1.1 材料與儀器

共收集370 份六大茶類樣品，包括122 種烏龍茶、110 種綠茶、55 種紅茶、34 種黑茶、29 種白茶、20 種黃茶，來自中國福建、廣東、臺灣等地區(qū)，以及日本、斯里蘭卡等國家。 所有茶葉樣品均具有該茶類正常的商品外形及特有的色、香、味，無異味、無劣變、無污染、無非茶類夾雜物、無任何添加劑，滿足實驗材料要求，具體信息詳見表1。

表1 各茶類樣品信息Table 1 Sample information of each tea category

MPA 型傅里葉變換近紅外光譜儀：德國布魯克光譜儀器公司；高速粉碎機：上海鼎廣機械設備有限公司；CFJ-II 茶葉篩分機： 杭州大吉光電儀器有限公司。

1.2 實驗方法

1.2.1 樣品制備茶葉粉碎后，置于茶葉篩分機中過80 目篩，取篩下茶粉5 g，密封編號，放于4 ℃冰箱中備用。

1.2.2 光譜采集使用MPA 型傅里葉變換近紅外光譜儀采集樣品光譜信息。 儀器工作時溫度控制在25 ℃，相對濕度＜70%。 光譜采集工作流參數(shù)：波數(shù)為3 600~12 500 cm-1，光譜掃描次數(shù)為64 次，分辨率為8.0 cm-1。 為確保近紅外光譜檢測數(shù)據(jù)的可靠性，每個樣品掃描3 次，取平均光譜作為原始光譜數(shù)據(jù)進行后續(xù)分析。

1.2.3 光譜預處理采用Minmax 算法增強數(shù)據(jù)；選用CWT 算法校正基線漂移并消除高頻噪聲；使用SNV 及MSC 算法校正散射， 消除因茶粉粒徑的不均勻、光程不恒定等因素所帶來的影響。

1.2.4 特征提取為提升模型性能、 運算效率，選用PCA、LDA 及SPA 方法進行光譜數(shù)據(jù)的特征提取，降低數(shù)據(jù)維數(shù)。

1.2.5 模型構建與評價數(shù)據(jù)挖掘分類器廣泛應用于NIRS 數(shù)據(jù)的分析與利用， 不存在始終保持最優(yōu)效果的分類器，因此使用多種分類器建模更利于優(yōu)質(zhì)模型的構建[21]。 本文中基于RF、SVM 兩種分類器，結(jié)合不同預處理、特征提取方法，優(yōu)化模型參數(shù)，探究六大茶類最佳識別模型構建流程。

為確保所建模型的適用性， 將數(shù)據(jù)按照3∶1 的比例劃分為訓練集和驗證集兩個子集，樣本數(shù)分別為277 和93 個，其中訓練集用于模型訓練，驗證集用于測試模型穩(wěn)健性。使用RA、AUC 以及混淆矩陣作為模型精度及性能的評價指標。

1.2.6 數(shù)據(jù)處理軟件數(shù)據(jù)處理軟件包括MATLAB 2016a、Origin 2019b、Excel 及Python。

2 結(jié)果與分析

2.1 光譜分析

370 份茶葉樣品在3 600~12 500 cm-1波段的近紅外光譜如圖1（a）所示，各波段吸光度變化趨勢趨于一致。 隨著波數(shù)的增加，吸光度總體呈現(xiàn)下降趨勢，變化范圍處于0.249～2.196。

圖1 370 份茶葉樣品近紅外光譜和六大茶類平均光譜圖Fig. 1 Near-infrared spectra of 370 tea samples and average spectra of six major tea types

六大茶類的平均光譜在大多波段范圍內(nèi)趨勢一致（見圖1（b）），吸收光譜趨于平行，各茶類對應譜圖于3 750 cm-1處初步分開，于3 750~9 000 cm-1處波動最為明顯。受C—H+C—H 組合頻伸縮、變形振動影響，4 000 cm-1處出現(xiàn)明顯波峰；參照相關文獻[22-25]，4 500~4 950 cm-1的譜帶變化可歸因于N—H+O—H 的組合頻；5 000~5 300 cm-1的譜帶變化主要受O—H+O—H 的組合頻振動影響；5 600~6 300 cm-1的譜帶變化受C—H 和S—H 的一級倍頻振動影響；6 700~7 400 cm-1的譜帶變化主要受N —H 的一級倍頻振動影響；8000~9000cm-1中出現(xiàn)的波峰可能與CHCH及C—H的二級倍頻振動有關。

3 750~9 000 cm-1中各平均光譜能基本分開，說明六大茶類樣品吸光度在該波段增減性不同，即光譜信息與茶類間具有相關性。 首段譜帶中各平均光譜間交叉重疊現(xiàn)象頻發(fā)，說明該波段體現(xiàn)的與產(chǎn)地相關的有效信息較少，信噪比低；末端譜帶波動同原始光譜一致，趨于平緩，無明顯波峰波谷，特征信息不顯，因此結(jié)合Meng 等結(jié)論[23]，本文中的模型構建使用3 750~9 000 cm-1波段的光譜數(shù)據(jù)。

2.2 光譜預處理

為校正光譜采集過程中因環(huán)境、光程不恒定及樣品粒徑差異等因素所帶來的誤差， 從數(shù)據(jù)增強、基線校正、 散射校正等角度， 使用Minmax、CWT、SNV 及MSC 算法對OS 進行預處理。 CWT 中小波參數(shù)（小波基、分解尺度）的選擇至關重要，直接決定了后續(xù)模型的優(yōu)劣，經(jīng)比對分析后，本研究中選擇應用最為廣泛的db（daubechies）族中的db4 小波基，分解尺度定為100[26]。

從光譜變化情況（見圖2）可知，4 種處理方式都使光譜形態(tài)發(fā)生了較大改變。Minmax（見圖2（a））將光譜吸光度凝練到-1.0~1.0，增強了數(shù)據(jù)，消除了數(shù)據(jù)量綱及取值范圍的影響，后續(xù)可使所建模型收斂速度加快，提高模型性能。 應用CWT 進行光譜預處理如圖2（b）所示，其形態(tài)變換程度為4 種預處理方式中最大，基線漂移、背景干擾、噪音現(xiàn)象等得到明顯消除，譜峰更清晰，差異信息段更明顯。 受茶粉顆粒大小不均、產(chǎn)生的散射影響，采用SNV 與MSC進行預處理（見圖2（c）和圖2（d）），處理后光譜中的散射干擾得到明顯消除，特征信息更加突出。 相較于OS，預處理可有效消除光譜中因光散射、基線漂移等造成的信號干擾，但處理后的光譜圖仍無法直觀分辨茶類間的差異，這可能是由于不同茶類在內(nèi)含物的組成與含量上具有較多的相似性。 為了進一步評估各處理對建模結(jié)果的影響，將各處理所得的光譜吸光度分別作為模型的輸入變量，后續(xù)基于模型評價指標判斷預處理效果。

圖2 光譜預處理效果Fig. 2 Effect of spectral pretreatment

2.3 光譜數(shù)據(jù)降維

近紅外光譜的連續(xù)波數(shù)中存在大量的冗余信息，其與特征信息之間存在很強的相關性[27]。 通過選取特征向量或波數(shù)降低數(shù)據(jù)維數(shù)，可保留原始數(shù)據(jù)中的主要特征信息， 減少后續(xù)處理的計算任務。本研究中采用PCA、LDA 和SPA 方法進行數(shù)據(jù)降維。

2.3.1 PCAPCA 是一種常用于降低大數(shù)據(jù)集維數(shù)的無監(jiān)督特征提取方式，能從大量數(shù)據(jù)中提取出特征，轉(zhuǎn)換為仍包含絕大部分有效信息卻擁有較小維數(shù)的數(shù)據(jù)集， 最大程度保留原始數(shù)據(jù)信息， 因而PCA 是一種最優(yōu)、最為常用的方法[28]。 使用PCA 方法對3 750~9 000 cm-1中OS 及經(jīng)Minmax、CWT、SNV 和MSC 預處理后的4 種光譜進行降維， 截取前15 個主成分特征值與累積貢獻度。 結(jié)果如表2所示， 以特征值大于1、 累積貢獻度大于80%為原則， 篩選模型輸入主成分個數(shù)。 OS 以及Minmax、CWT、SNV 和MSC 預處理后光譜分別篩出6、11、13、12、12 個主成分，累積貢獻度分別達到99.89%、99.71%、99.67%、99.82%、99.82%，符合原則?；诤Y選的主成分構建模型。

表2 PCA 特征值及累積貢獻度Table 2 PCA feature values and cumulative contribution

2.3.2 LDALDA 是一種有監(jiān)督的特征提取方法[29]，在茶葉領域常作為分類器使用。而利用LDA 進行光譜特征提取、降低數(shù)據(jù)維數(shù)，并結(jié)合分類器建立茶類識別模型的研究，目前尚未見相關報道。 LDA 最多可使數(shù)據(jù)矩陣降至類別數(shù)減1 的維數(shù)，降低維數(shù)的同時不過多丟失原始信息，LDA 降維后所得維數(shù)將被用于后續(xù)建模。

2.3.3 SPASPA 是一種前向循環(huán)特征提取方法，其可以通過連續(xù)投影的方式從原始光譜矩陣中提取有效預測響應變量的信息，最大限度減少光譜變量之間的共線性效應，達到所選響應變量預測能力的最大化。 其主要通過將波數(shù)投影于其他波數(shù)，比較投影向量大小， 波數(shù)間投影向量大者為待選波數(shù)，最終投影向量最大且與特征集內(nèi)波數(shù)共線性最小的波數(shù)選入特征集合[30]。 特征波數(shù)的數(shù)目由校準集內(nèi)部完全交叉驗證的均方根誤差 （root mean square error，RMSE）確定，與最小RMSE 對應的特征波數(shù)數(shù)目和特征波數(shù)為最佳值[31-32]。

由SPA 方法提取的特征向量如圖3 所示。 SPA中4 種預處理后的光譜信息的RMSE 迭代下降曲線分別如圖3（a）、圖3（c）、圖3（e）、圖3（g）和圖3（i）所示。 從圖中可以看出，當選擇特定數(shù)量的波數(shù)時，RMSE 達到最小值； 而后RMSE 雖仍波動下降，但降幅很小且導致所選波數(shù)增加，沒有必要為了追求微小的RMSE 而增加維數(shù)。因此，最終從OS 及經(jīng)Minmax、CWT、SNV 和MSC 預處理光譜中獲得的特征波數(shù)數(shù)目分別為15、16、19、15、7 個 （如圖3 所示，特征波數(shù)具體信息如表3）。

圖3 連續(xù)投影算法（SPA）提取特征波數(shù)Fig. 3 Extraction of feature wavenumbers by successive projections algorithm (SPA)

表3 通過SPA 篩選出的特征波數(shù)Table 3 Feature wavenumbers selected by SPA

2.4 模型參數(shù)優(yōu)化

SVM 是近年來茶葉中應用最廣、效果最好的機器學習方法之一。 它是一種利用核函數(shù)將n 維輸入向量映射到K 維特征空間（K>n），從而通過高維特征空間進行分類的算法[33]。為提高模型質(zhì)量，本文中所有SVM 模型皆基于高斯 （radial basis function，RBF）核函數(shù)，該核函數(shù)可降低訓練過程的計算復雜度，在一般平滑假設下具有良好性能；與此同時，懲罰參數(shù)及gamma 參數(shù)最優(yōu)值的確定也至關重要，SVM 模型精度取決于這兩個參數(shù)的組合。根據(jù)初步試算結(jié)果，將懲罰參數(shù)取值定為1×103、1×104、1×105、1×106，gamma 參數(shù)取1×100、1×10-1、1×10-2、1×10-3。將3 750~9 000 cm-1的OS 數(shù)據(jù)作為輸入量，結(jié)合模型的識別正確率進行參數(shù)優(yōu)化。 結(jié)果如表4 所示，當懲罰參數(shù)為1×106、gamma 參數(shù)為1×10-2時模型具最佳識別正確率，后續(xù)模型皆基于此參數(shù)構建。

表4 SVM 模型參數(shù)優(yōu)化Table 4 Parameter optimization of SVM model

RF 是一種有監(jiān)督的集成分類算法， 主要是為解決單一決策樹可能出現(xiàn)的很大誤差和過擬合的問題，在分類問題中表現(xiàn)優(yōu)異，具有成為各情況下效果最優(yōu)分類器的巨大潛力[34]。 該模型由許多的決策樹組成，但每一棵決策樹之間沒有關聯(lián)，得到森林之后對新樣本進行判斷或預測時，將由森林中的每一棵決策樹分別進行判斷， 分辨該樣本屬于哪類，比對出選擇數(shù)最多的類別，從而對樣本類別做出判斷， 因此該模型中樹木數(shù)量的選擇極為重要。OS 數(shù)據(jù)試算后， 選定樹木數(shù)量為1~100 進行參數(shù)優(yōu)化，結(jié)果如圖4 所示，樹木數(shù)量為70 時，RF 模型具最佳識別正確率，后續(xù)模型皆基于此參數(shù)構建。

圖4 RF 模型參數(shù)優(yōu)化Fig. 4 Parameter optimization of RF model

2.5 模型性能評價

RA 常被用于模型預測能力的評估，AUC 則常被用于模型泛化能力的評價， 取值區(qū)間為0.5～1.0，值的大小與模型質(zhì)量呈正相關。 因此， 采用RA、AUC 聯(lián)合評估模型性能。 此外，為直觀呈現(xiàn)所建模型對各茶類識別性能的優(yōu)劣，引入混淆矩陣對模型預測結(jié)果進行評價。

由表5 可知， 各茶類NIRS 數(shù)據(jù)結(jié)合不同預處理、 特征提取方法及數(shù)據(jù)挖掘分類器最終獲得40個茶類識別模型， 識別正確率介于59.14%~100.00%，AUC 處于0.70~1.00， 大多數(shù)模型識別正確率高于70%，模型性能良好。 OS 結(jié)合RF、SVM 所得模型識別正確率分別為69.89%及92.47%，光譜經(jīng)預處理、特征提取后，建立的RF 模型絕大多數(shù)識別正確率、AUC 顯著提升，模型精度、泛化能力均得到改善， 最佳模型OS-LDA-RF 的識別正確率可達94.62%，AUC 可達0.96；SVM 結(jié)合預處理后光譜建模效果欠佳，多出現(xiàn)識別準度下降的問題，不同特征提取方法中效果最佳的為LDA，結(jié)合不同預處理后的光譜數(shù)據(jù)皆優(yōu)化了模型性能， 基于SVM 建立的茶類識別模型中最佳的是OS-LDA-SVM，識別正確率為100.00%，AUC 為1.00。 總體而言，不同分類器結(jié)合不同預處理、特征提取方法所取得的效果也不盡相同，在茶類識別模型構建中，RF 適合與不同化學計量學方法結(jié)合，多數(shù)預處理、特征提取方法對RF 模型性能優(yōu)化效果顯著；SVM 適合基于去除頭尾信息匱乏波段的原始光譜結(jié)合特征提取方法進行建模，所得模型特征數(shù)可得到簡化，提升運算速度及模型質(zhì)量。

觀察最佳模型混淆矩陣可知 （見圖5），OSLDA-RF 對烏龍茶、 黑茶、 白茶的識別正確率為100.00%，誤判發(fā)生于綠茶、紅茶及黃茶的識別中，三者識別正確率分別為89.28%、92.86%以及80.00%；OS-LDA-SVM 對各茶類的識別正確率皆為100.00%，混淆矩陣中數(shù)值皆處于對角線，表明識別效果優(yōu)異，模型質(zhì)量好。

圖5 最佳模型混淆矩陣Fig. 5 Optimal model confusion matrix

對性能最優(yōu)模型OS-LDA-SVM 進行三維空間可視化（見圖6）。 可發(fā)現(xiàn)，識別正確率達100.00%的茶類識別模型中，綠茶與黃茶的光譜特征于三維空間上的分布極為接近，可能與其加工工藝的高度相似有關，黃茶僅比綠茶多了悶黃工藝；還可能受黃茶樣本數(shù)量較少導致特征信息不顯的影響，黃茶在六大茶類中占比最小，僅在我國四川、湖南等少數(shù)省份小規(guī)模生產(chǎn)，后期將通過逐年增加黃茶樣本數(shù)量的方式，強化黃茶光譜特征，提高模型性能。 其他茶類的光譜特征在三維空間分布差異大，可能與加工工藝、茶樹品種不同有關。

圖6 OS-LDA-SVM 模型三維空間效果Fig. 6 Three-dimensional space effect of OS-LDA-SVM model

3 結(jié) 語

基于茶類的370 個近紅外光譜數(shù)據(jù)，剔除認為不含茶類相關信息的波段后， 使用Minmax、CWT、SNV 及MSC 進行預處理，PCA、LDA 及SPA 進行特征提取，最后基于RF、SVM 構建茶類識別模型。 主要結(jié)論如下：

1）茶類識別模型構建中，不同預處理對不同分類器的效果不盡相同。對RF 模型而言，預處理對模型性能提升效果顯著， 而相同處理于SVM 模型中效果欠佳，模型質(zhì)量多不如原始光譜模型。

2）通過特征提取優(yōu)化模型，效果顯著。 PCA、LDA、SPA 方法皆大幅度降低了數(shù)據(jù)維數(shù)，提高了模型運算效率。 其中LDA 效果最佳，與不同預處理方法、分類器結(jié)合所得模型皆質(zhì)量優(yōu)異。

3）RF、SVM 皆適用于茶類識別模型構建， 相較RF 模型，SVM 模型總體效果略勝一籌。 RF 模型中最佳模型為OS-LDA-RF，RA 為94.62%， 對烏龍茶、 黑茶及白茶的RA 可達100.00%，AUC 為0.96，模型性能優(yōu)異、穩(wěn)定；OS-LDA-SVM 為SVM 模型中的最優(yōu)模型，不同茶類的RA 皆達100.00%，AUC 為1.00，模型質(zhì)量高、泛化能力好。

不同茶類的近紅外光譜數(shù)據(jù)， 經(jīng)適當預處理、特征提取方法的選擇后，結(jié)合RF、SVM 分類器可挖掘出近紅外光譜中茶類識別相關的關鍵信息，構建出高識別正確率的茶類識別模型。 基于此，針對光譜特征接近的茶類，可逐年擴大樣本量，優(yōu)化與驗證茶類識別模型性能，提高模型適用性。 除此之外，本研究后續(xù)將利用互聯(lián)網(wǎng)技術構建六大茶類的近紅外光譜數(shù)據(jù)庫，搭建在線茶類識別平臺，開展遠程茶類識別，以期早日在國際市場上實現(xiàn)產(chǎn)業(yè)化應用，促進我國茶產(chǎn)業(yè)的高質(zhì)量發(fā)展。