王 凡，陳龍躍，段丹丹, 4*，曹 瓊，趙 鈺，藍(lán)玩榮

1.國家農(nóng)業(yè)信息化工程技術(shù)研究中心，北京 100097 2.清遠(yuǎn)市智慧農(nóng)業(yè)農(nóng)村研究院，廣東 清遠(yuǎn) 511500 3.農(nóng)芯科技(廣州)有限責(zé)任公司，廣東 廣州 510000 4.湖南農(nóng)業(yè)大學(xué)，湖南 長沙 410125 5.江門市農(nóng)業(yè)技術(shù)服務(wù)中心，廣東 江門 529000

引 言

茶樹是世界上廣泛種植的作物，因其獨(dú)特的滋味以及強(qiáng)大的健康功效而受到人們的喜愛。眾所周知，氮(N)是植物光合作用和生長必需的元素，實(shí)時(shí)進(jìn)行氮素管理, 對于保證茶葉產(chǎn)量, 減少環(huán)境污染[1]具有指導(dǎo)意義。

傳統(tǒng)的作物氮素檢測方法費(fèi)時(shí)耗力，對樣品具有破壞性，需具備較高的操作技能，無法滿足實(shí)時(shí)準(zhǔn)確的測定需求[2]。高光譜作為一種成熟的無損檢測技術(shù)已被廣泛應(yīng)用于茶葉的營養(yǎng)診斷、種類識別、長勢檢測[3-4]。光譜測量植被生化變量中有光譜信息冗余等問題，易受到噪音和背景因素的影響。如何從復(fù)雜重疊的光譜信息中得到微弱的化學(xué)成分變化的有效信息，提高檢測精度，是光譜分析的關(guān)鍵。研究表明，偏最小二乘方法(PLSR)，通過使用主成分分析(PCA)將自變量投影到新空間以概括潛在變量，可以處理大量的噪聲和相關(guān)變量，被證實(shí)同樣適用于少量變量的情況[5]。選擇合理的特征提取方法可以降低模型的復(fù)雜程度并提高預(yù)測準(zhǔn)確性和魯棒性，且少量變量即可提供與原始全光譜波段相媲美的信息。已有的研究報(bào)道了許多方法，如SPA，UVE，PCA，RF以及CARS等[6-7]，這些方法都是選擇信息變量非常有效的策略，然而其變量選擇方法大多是利用蒙特卡洛采樣法在樣本變量空間進(jìn)行重復(fù)取樣后一次性建模，卻沒有考慮到變量之間的組合效應(yīng)。此外，在整個(gè)光譜區(qū)間進(jìn)行變量搜索，這些方法僅消除了非信息性變量而未考慮干擾變量的影響，當(dāng)存在干擾變量時(shí)，將無法獲得令人滿意的結(jié)果。連續(xù)小波變換(CWT)用于處理高光譜數(shù)據(jù)，通過縮放和移動(dòng)小波函數(shù)將原始光譜分解為多個(gè)尺度，促進(jìn)對細(xì)微特征信號的識別，從而構(gòu)建對光譜敏感的生物和化學(xué)成分的預(yù)測模型，已廣泛應(yīng)用于植被參數(shù)評估。變量組合集群分析(VCPA)算法，是基于達(dá)爾文進(jìn)化論的“適者生存”以及依據(jù)交叉驗(yàn)證均方根誤差(RMSE)的一種變量選擇方法[8]，從全光譜數(shù)據(jù)中選擇光譜變量的最佳組合，利用變量的隨機(jī)組合考慮變量間的相互作用，進(jìn)行多次迭代的分析，直到剩下的變量沒有無信息變量和干擾變量，在特征提取方面效果較好[9-10]。連續(xù)投影算法(SPA)、競爭性自適應(yīng)重加權(quán)算法(CARS)廣泛用于高光譜信息篩選中以提高模型性能，而關(guān)于VCPA的應(yīng)用研究較少；另外，基于小波分析的品質(zhì)參數(shù)估算，在茶鮮葉監(jiān)測上的應(yīng)用鮮有報(bào)道，目前尚未有開發(fā)CWT與VCPA特征變量選擇算法相結(jié)合的模型優(yōu)化方法的應(yīng)用實(shí)例。

以茶鮮葉為研究對象，利用小波分解對茶鮮葉光譜展開分析，確定了CWT在茶葉N素含量鑒定中的性能，結(jié)合實(shí)驗(yàn)室理化檢測茶鮮葉N素含量，采用SPA，CARS和VCPA特征波段選擇算法提取有益信息，建立茶鮮葉N素含量預(yù)測模型，并探究在小波算法不同分解尺度下，不同特征波段選擇算法對茶鮮葉N素估算精度的影響。

1 實(shí)驗(yàn)部分

1.1 研究區(qū)概況

實(shí)驗(yàn)區(qū)位于廣東省清遠(yuǎn)市英德市積慶里、創(chuàng)美、德高信、八百秀才和石門山5個(gè)茶園(圖1)，該區(qū)域位于北緯24.50°，東經(jīng)113.55°，平均海拔為181 m，地處南亞熱帶向中亞熱帶的過渡地區(qū)，屬亞熱帶季風(fēng)氣候。茶鮮葉光譜數(shù)據(jù)采集在2020年5月、7月和9月開展，覆蓋3個(gè)品種(包括主栽品種英紅九號，鴻雁十二號，金萱)，采樣點(diǎn)均勻分布于實(shí)驗(yàn)區(qū)內(nèi)，共選取了151個(gè)樣本點(diǎn)，每個(gè)樣點(diǎn)采集茶鮮葉樣品約200 g，進(jìn)行化驗(yàn)。

圖1 研究區(qū)域空間分布

1.2 數(shù)據(jù)采集

葉片的光譜反射率采用美國ASD Fieldspec Pro3光譜儀測量，其波譜覆蓋范圍為350～2 500 nm，光譜重采樣間隔為1 nm。利用葉片夾持器(自帶光源)選取葉片的中部較寬的部分進(jìn)行光譜測定，手柄中的光源部分直接夾緊茶葉葉片所測部位，葉片夾持器測定使葉片處于密閉條件下，且所測面積相同, 保證光譜數(shù)據(jù)獲取的精確性和穩(wěn)定性。取10片標(biāo)準(zhǔn)葉片在相同樣點(diǎn)采集光譜數(shù)據(jù)，10個(gè)光譜反射率的平均值作為該樣點(diǎn)葉片的光譜反射率。

采集光譜數(shù)據(jù)后的茶鮮葉放于120 ℃烘箱中殺青固樣30 min，然后置于80～90 ℃下烘至足干。干燥后的樣品磨碎后參照GB5009.5—2016，通過凱氏定氮法測量其全氮量。

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

1.3.1 常規(guī)光譜預(yù)處理方法

采用卷積平滑(SG)、去趨勢(Detrending)、一階導(dǎo)數(shù)(1st)、多元散射校正(MSC)和標(biāo)準(zhǔn)正態(tài)變量變換(SNV)五種預(yù)處理方法，對原始光譜數(shù)據(jù)進(jìn)行預(yù)處理作為參考。

1.3.2 小波變換

連續(xù)小波變換(CWT)起源于傅里葉變換信號處理技術(shù)，同時(shí)開展頻率和時(shí)間分析，從復(fù)雜的原始光譜信號中分辨細(xì)微的變換，捕捉有益信息，表征光譜局部特征信息。小波變換分析廣泛應(yīng)用于高光譜數(shù)據(jù)的降維壓縮、噪聲去除、特征選擇以及弱信號提取等方面。小波分析包括連續(xù)分解與離散分解兩種類型的小波變換，Wang等[11]研究表明離散小波隨著分解尺度的增加，光譜維度減半，不能對光譜數(shù)據(jù)進(jìn)行較好的分析。連續(xù)小波變換是一種連續(xù)分解信號的方法，其輸出光譜可以與原始光譜直接進(jìn)行對比，比離散小波變換后獲得的光譜更易于解釋。通過以式(1)縮放和移動(dòng)小波母函數(shù)

(1)

式(1)中，λ為茶葉葉片高光譜數(shù)據(jù)的波段數(shù)，a和b分別代表伸縮因子(即寬度)和平移因子(即位置)，這兩個(gè)變量都為正實(shí)數(shù)。CWT的輸出公式如式(2)

(2)

式(2)中，f(λ)為高光譜反射率，CWT系數(shù)(Wf(ai，bj)，i=1，2，…，m;j=1,2，…，n)包含兩維(即m×n矩陣)，其中一維為波長(1, 2, …,n)，另一維是分解尺度(1, 2, …,m)，因此，小波系數(shù)是一個(gè)行為尺度，列為波長的矩陣。

由于Mexh小波在高光譜分解方面表現(xiàn)出色[12]，且光譜吸收特性與高斯函數(shù)相似，因此選取Mexh作為小波母函數(shù)，使用連續(xù)小波分析(continuous wavelet transformation，CWT)以執(zhí)行尺度21～210分解反射光譜為10個(gè)尺度，并獲取不同分解尺度下的不同波長的一系列小波系數(shù)，進(jìn)行了小波系數(shù)與茶葉N素含量的相關(guān)性分析，并用于構(gòu)建校準(zhǔn)模型。

1.3.3 特征變量選擇方法

連續(xù)投影(SPA)是一種前向變量選擇方法，使用簡單的投影操作來獲得具有最小共線性的變量子集，并可以根據(jù)交叉驗(yàn)證的最小均方根誤差(RMSE)選擇最終的最佳波長。

競爭性自適應(yīng)重加權(quán)采樣(CARS)結(jié)合蒙特卡洛采樣和PLSR，基于達(dá)爾文進(jìn)化論，遵循“適者生存”的原則，根據(jù)每個(gè)變量的重要性，以迭代和競爭的方式從N個(gè)蒙特卡洛采樣順序中選擇N個(gè)變量子集。

變量組合集群分析(variable combination population analysis, VCPA)是一種新的有效的變量選擇方法，它由Yun等建立[13]。VCPA的實(shí)現(xiàn)過程分為兩個(gè)基本步驟，首先基于“適者生存”原理的指數(shù)遞減函數(shù)(exponentially decreasing function，EDF)選擇有用的變量并縮小變量空間，其次通過二進(jìn)制矩陣采樣(binary matrix sampling，BMS)選擇最佳光譜波長子集，通過生成隨機(jī)變量組合，確保每個(gè)變量都有機(jī)會被選擇并產(chǎn)生各種變量子集，在每個(gè)EDF和BMS運(yùn)行的周期中都采用模型集群分析(model population analysis，MPA)根據(jù)較低的RMSECV來確定最優(yōu)變量子集。分別利用SPA，CARS和VCPA對連續(xù)小波變換處理后的光譜數(shù)據(jù)進(jìn)行進(jìn)一步的特征波段提取。

1.4 模型評價(jià)標(biāo)準(zhǔn)

經(jīng)過連續(xù)小波變換后的光譜，通過SPA、CARS和VCPA方法進(jìn)行變量篩選，利用偏最小二乘法(partial least squaresreg ression，PLSR)構(gòu)建茶鮮葉N含量的估測模型。采用Kennard-Stone(KS)方法以近似3∶1比例選取113個(gè)樣本數(shù)據(jù)作為建模集，用于構(gòu)建估測模型，剩余的38個(gè)樣本數(shù)據(jù)作為驗(yàn)證集，用于預(yù)測模型精度和穩(wěn)定性的驗(yàn)證。模型精度采用模型評價(jià)參數(shù)決定系數(shù)(R2)與均方根誤差(RMSE)來綜合評價(jià)模型表現(xiàn)，其計(jì)算方法如式(3)和式(4)

(3)

(4)

1.5 分析軟件

使用Matlab R2018b(The Math Works, Natick, USA)、Unscrambler X10.4(CAMO AS，Oslo, Norway)、SPSS 26(IBM, USA)、Origin 2019b(OriginLab, USA)等軟件。

2 結(jié)果與討論

2.1 常規(guī)光譜預(yù)處理結(jié)果

圖2展示了5種預(yù)處理后的高光譜。與原始光譜(R)相比，預(yù)處理1st和Detrending方法明顯改變了光譜特性，在具體波段上凸顯了光譜信息，而SG，MSC和SNV預(yù)處理方法沒有明顯改變光譜吸收特性，光譜特征變化幅度小，與原始光譜相比更加平滑與集中。上述5種預(yù)處理算法有效地平滑光譜曲線，且可以保留光譜信息的有用信息。

圖2 不同預(yù)處理方法下光譜反射率

2.2 小波分析結(jié)果

2.2.1 基于小波變換的光譜分析

原始光譜經(jīng)CWT處理后的結(jié)果如圖3所示，D1—D10分別代表了以21—210為分解尺度經(jīng)過CWT分解得到的10層小波能量系數(shù)模擬光譜分辨率。可以看出，1—7尺度上的光譜分辨率，隨著分解層數(shù)的增加，光譜吸收、反射波譜特征均加強(qiáng)，光譜波動(dòng)范圍加大，8—10尺度的光譜反射率起伏變化較大?？傮w來看，隨著CWT分解尺度的不斷增加，原始光譜反射率的波峰特征(761，1 070，1 217和1 645 nm)與波谷特征(670，902，1 447和1 931 nm)被不斷放大，并且均有一定的紅移與藍(lán)移現(xiàn)象。小尺度分解可以較好的保持原始光譜特性，隨著分解層數(shù)的增加，光譜特征在不斷凸顯的過程中也不斷被拉伸，在大尺度分解中又可以凸顯具體波段的細(xì)微特征信息。這是因?yàn)镃WT是基于光譜曲線的形狀特征，捕獲光譜形狀的差異并凸顯特征[14]，更有利于捕捉原始光譜反射率的微弱變化，有利于進(jìn)行N含量對茶葉鮮光譜響應(yīng)特征的分析。

圖3 小波變換光譜特性

2.2.2 相關(guān)性分析

小波系數(shù)能量與N素實(shí)測數(shù)據(jù)進(jìn)行相關(guān)性分析，如圖4所示。可看出各尺度小波系數(shù)與茶鮮葉N含量的相關(guān)系數(shù)差異較大，且隨著分解尺度的增加呈逐漸集中分布的趨勢。由圖可知，相關(guān)性較高的波段表現(xiàn)為更明顯的紅色區(qū)域，相關(guān)系數(shù)較大的波長主要與N素、葉綠素以及類胡蘿卜素等物質(zhì)的復(fù)雜性有關(guān)，對于選定的493，987，1 193，1 195和1 825 nm，這些波段大多為N素的重要光譜區(qū)域[15]，另外可以看出，小波變換后的光譜在小尺度分解下與茶鮮葉N素存在良好的相關(guān)性，選取1—6尺度作進(jìn)一步分析。

圖4 小波變換相關(guān)性分析

2.3 特征波長選擇

為了比較不同變量篩選方法對茶鮮葉N含量高定量模型的影響，將經(jīng)過CWT處理的光譜數(shù)據(jù)分別利用SPA，CARS和VCPA進(jìn)行變量篩選。SPA設(shè)置選擇的最小和最大變量數(shù)的參數(shù)分別為1和30，分別對CWT處理后的小波系數(shù)進(jìn)行光譜特征變量的選取，結(jié)果如圖5(a)所示，光譜變量分別減少到15，16，17，11，14和18個(gè)。CARS設(shè)置蒙特卡洛抽樣運(yùn)行次數(shù)為50，十倍交叉驗(yàn)證用于評估每個(gè)子集的效用，分別對CWT處理后的小波系數(shù)進(jìn)行光譜特征變量的選取。結(jié)果如圖5(b)所示，光譜變量分別減少到16，19，17，16，15和17個(gè)。VCPA設(shè)置EDF剩余變量數(shù)為100，迭代次數(shù)為50，最優(yōu)子模型的比例為10%，BMS運(yùn)行次數(shù)為1 000，VCPA步驟中每個(gè)BMS采樣的平均數(shù)設(shè)置為20，并采用5折交叉驗(yàn)證分別對CWT處理后的小波系數(shù)進(jìn)行光譜特征變量的選取，結(jié)果如圖5(c)所示。可以看出基于CWT小波分解后的光譜數(shù)據(jù)通過VCPA的五次迭代將全波段光譜變量分別減少到12，11，10，11，8和8個(gè)。相比SPA和CARS，VCPA選擇的變量更少，這表明該方法在數(shù)據(jù)降維方面有著巨大的優(yōu)勢。

圖5 SPA，CARS和VCPA篩選的特征變量

2.4 不同預(yù)處理方法結(jié)合特征波段篩選的N素含量預(yù)測

2.4.1 基于常規(guī)預(yù)處理方法和特征篩選的N素含量預(yù)測

采用卷積平滑(SG)、去趨勢(Detrending)、一階導(dǎo)數(shù)(1st)、多元散射校正(MSC)和標(biāo)準(zhǔn)正態(tài)變量變換(SNV)5種傳統(tǒng)光譜處理技術(shù)分析茶鮮葉光譜特性，并采用偏最小二乘算法建立茶鮮葉N素含量預(yù)測模型作為參考，結(jié)果如表1所示。

表1 基于不同預(yù)處理方法和變量選擇方法的模型

2.4.2 基于小波變換和特征篩選的N素含量預(yù)測

使用CWT處理后的光譜數(shù)據(jù)篩選后的特征波段作為自變量，以茶鮮葉N含量為因變量，采用偏最小二乘法構(gòu)建茶鮮葉N含量的估測模型(表2)。由表2可知，經(jīng)過CWT處理后光譜構(gòu)建的模型均表現(xiàn)出較好的擬合精度，并且隨著分解尺度的增加模型精度呈現(xiàn)下降趨勢，說明小規(guī)模分解通常能夠改善CWT處理后的小波系數(shù)對分析物的敏感性，這與Liu等的研究結(jié)果一致[16]。其中SPA選擇變量估測茶鮮葉N含量的精度較低，1—6尺度建模R2范圍在0.65～0.83，RMSE范圍在0.25～0.33，預(yù)測R2范圍在0.65～0.83，RMSE范圍在0.25～0.33；以CARS選擇變量為因變量，1—6尺度建模R2范圍在0.73～0.91，RMSE范圍在0.19～0.32，預(yù)測R2范圍在0.64～0.88，RMSE范圍在0.24～0.35；VCPA選擇變量估測模型表現(xiàn)最好，1—6尺度建模R2范圍在0.73～0.95，RMSE范圍在0.16～0.29，預(yù)測R2范圍在0.73～0.90，RMSE范圍在0.23～0.35，其中在1尺度下的模型擬合效果最高。

表2 基于不同分解層數(shù)和變量選擇方法的模型

圖6所示為不同變量選擇方法下最優(yōu)估測模型實(shí)驗(yàn)值與預(yù)測值1∶1對比圖，可以看出CWT-VCPA模型表現(xiàn)最優(yōu)，訓(xùn)練集和驗(yàn)證集的R2分別為0.95和0.90，RMSE分別為0.16和0.23，二者差異較小，相比于CWT-SPA和CWT-CARS，建模及預(yù)測集的預(yù)測精度最高。對比各模型可知，CWT方法結(jié)合特征變量選擇算法明顯提升了模型的精度，且CWT-VCPA篩選變量估算N含量效果最好。

圖6 基于(a)CWT(1尺度)-SPA-PLSR，(b)CWT(1尺度)-CARS-PLSR和(c)CWT(1尺度)-VCPA-PLSR的氮含量最佳模型的擬合1∶1擬合圖

以茶鮮葉為研究對象，構(gòu)建并對比分析了常用5種光譜處理與連續(xù)小波分解方法的N含量估測模型。研究結(jié)果表明，與常規(guī)光譜處理方法相比，小波變換方法可明顯提升茶鮮葉光譜對N含量的響應(yīng)，且模型精度與穩(wěn)定性均較高。這是由于在常規(guī)的光譜預(yù)處理方法中，SNV，MSC和Detrending通常用于校正低頻散射信號、1st可以消除信號的基線漂移，而CWT可以同時(shí)消除光譜信號中的高頻和低頻噪聲，并且連續(xù)小波分解本質(zhì)上是一種基于相關(guān)性的方法，可以捕獲對茶多酚含量最敏感的小波特征，這種頻譜相關(guān)性的方法在凸顯光譜細(xì)微變化特征方面效果顯著，更有利于后續(xù)對光譜特征波段的提取。隨著分解尺度的增加，模型預(yù)測精度逐步降低，主要原因是光譜中含有較多的噪聲，造成了光譜的有益信息與噪聲的分布波段區(qū)間有較高的重疊性，且隨著CWT分解尺度的增加，其內(nèi)含的有益信息與無益噪聲信息重新分配，導(dǎo)致引入的噪聲信息較多從而降低模型的預(yù)測精度。在1—6尺度模型精度達(dá)到最低時(shí)相較于傳統(tǒng)光譜變換仍然有一定程度的提升，表明CWT能提取光譜特征中的有益信息，并一定程度上分離無用噪聲，提升光譜信息信噪比。基于三種不同波段篩選方法(SPA，CARS和VCPA)都可以在減少變量的同時(shí)提升模型預(yù)測效果，且相比于SPA和CARS，VCPA選擇的變量更少，預(yù)測效果更好。這歸因于EDF對變量的連續(xù)收縮，在EDF的每次運(yùn)行中，都會消除一些貢獻(xiàn)較小的變量。此外，還與BMS和MPA的使用有關(guān)，BMS是一種采樣策略，它為所有變量提供了相同的機(jī)會進(jìn)行采樣，使用不同的變量組合用于子模型的構(gòu)建。MPA用于提取子模型中包含的最佳變量組合，比較每個(gè)子模型的RMSECV，保留最佳10%子模型中具有最高頻率的變量[17]。因此，VCPA通過最少的選擇變量，保證了光譜數(shù)據(jù)集更好的可預(yù)測性。當(dāng)將最佳波長與全光譜進(jìn)行比較時(shí)，基于CWT分解光譜結(jié)合VCPA篩選特征變量的PLSR模型比其他光譜變換的PLSR模型具有更高的性能(R2=0.95)，并且減少了全光譜數(shù)據(jù)變量的99.34%，極大的簡化了建模過程。綜上可知，CWT-VCPA方法在實(shí)現(xiàn)對茶鮮葉N含量的定量監(jiān)測方面具有強(qiáng)大的潛力。

3 結(jié) 論

采用傳統(tǒng)光譜處理技術(shù)和連續(xù)小波技術(shù)對光譜分別開展處理和分析，并通過SPA，CARS和VCPA方法進(jìn)一步篩選特征波段，建立預(yù)測茶鮮葉N素含量的模型。旨在研究基于連續(xù)小波技術(shù)結(jié)合特征變量選擇方法預(yù)測茶鮮葉N素含量的最優(yōu)模型，結(jié)論如下：

(1)相比于傳統(tǒng)光譜處理方法，經(jīng)過CWT分解光譜結(jié)合變量選擇方法構(gòu)建的茶鮮葉N含量估算模型均表現(xiàn)出令人滿意的結(jié)果，且隨著分解層數(shù)的增加，模型的精度整體有逐步下降的規(guī)律。表明CWT處理可有效提升光譜對N的敏感性，且小規(guī)模分解通常能夠改善這種關(guān)系，提高模型性能。

(2)以SPA，CARS和VCPA篩選出的特征變量為模型輸入因子均在CWT分解尺度為1時(shí)達(dá)到最佳估算效果，CWT-SPA-PLSR建模R2為0.83，RMSE為0.25；CWT-CARS-PLSR建模R2為0.91，RMSE為0.19；CWT-VCPA-PLSR建模R2為0.95，RMSE為0.16。