應用光譜指數和機器學習反演紫丁香葉片的葉綠素面密度1)

杜菲菲 安慧君 李賀新

(內蒙古農業大學,呼和浩特,010020)(內蒙古第二林業和草原監測規劃院)

紫丁香(SyringaoblataLindl.)是木犀科丁香屬植物,為落葉灌木或小喬木。因其獨特的觀賞價值和較強的環境適應性,被廣泛應用于我國北方各省區的城市園林綠化[1]。作為城市綠化應用最普遍的花木之一,對其生長狀態的監測尤為重要。研究表明[2],植物體內葉綠素含量能夠反映植物的營養狀況和生理代謝水平。傳統的葉綠素含量測定方法費時費力且難以大面積開展,20世紀末出現的高光譜技術能對植被進行微弱光譜差異的定量分析,具有信息獲取量大、實時、無損等優勢,彌補了傳統方法的局限性。

高光譜技術在葉綠素含量預測方面的研究,大多建立在一維層面[3-4],對原始或預處理光譜提取敏感波段或全波段進行建模,但由于輸入參數過多導致模型較為復雜。隨著研究的深入[5-8],二維光譜指數得到了發展,因其考慮到了光譜間的重疊吸收和相互影響,提高了對葉綠素的敏感程度,且經篩選出的最優光譜指數能夠有效降低模型復雜程度。有研究表明構建和篩選全波段三維光譜指數在土壤有機質含量[9]、土壤含水率[10]預測方面比二維更有優勢,但在葉綠素含量預測方面尚未有研究。近年來,機器學習算法在解析光譜參數與葉片生化組分之間的非線性問題時表現能力較一元和多元線性回歸好,如反向傳播神經網絡(BPNN)、支持向量回歸機模型(SVR)等。其中,2007年由Jayadeva提出的孿生支持向量回歸機(TSVR)模型[11]有良好的泛化性能,其運行效率是SVR的4倍,在金融[12]和工業[13-14]領域能很好地反映特征指標的非線性響應,但在植物生理含量預測方面尚未有研究。

本研究以內蒙古農業大學東校區內紫丁香為研究對象,測定葉片反射率和葉綠素面密度,在一維和二維光譜指數的基礎上,引入三維光譜指數,基于原始光譜(R)及其對應的一階微分(RFD)、二階微分(RSD)光譜構建全波段不同維度的光譜指數,通過皮爾遜相關系數法(PCC)進行最優光譜指數篩選,構建基于最優光譜指數和孿生支持向量回歸機(TSVR)的紫丁香葉片葉綠素面密度反演模型。為進一步加強模型的穩定性和泛化能力,采用海洋捕食者[15]仿生優化算法(MPA)對TSVR進行優化[16],提出MPA-TSVR融合模型,并與PLSR、SVR、TSVR模型進行對比,得到最優模型。旨在量化不同維度的光譜指數對紫丁香葉片葉綠素面密度的響應,拓寬優化算法融合機器學習算法在植物生理預測方向的研究,為城市綠化植物葉片葉綠素面密度反演提供新的思路與理論依據。

1 材料與方法

1.1 數據采集

試驗在7月下旬開展,在內蒙古農業大學東校區內無遮擋開闊綠地上選擇5叢生長正常且冠幅約2 m的紫丁香。每叢按樹冠中心分為8個方位,每一方位在表層選4枚生長正常的葉片,每叢采集32枚,共采集160枚;葉片采集完成后放入已編號的塑封袋中當即帶回實驗室,用蒸餾水反復清洗后晾干。用ASD公司生產的FieldSpec4 Hi-Res波譜儀采集葉片反射光譜(350～2 500 nm);每枚葉選取3個樣點,每個樣點用植物探頭測量3次;儀器每隔5 s記錄1次數據,每隔10 min進行1次白板校正。紫丁香葉片葉綠素面密度測定采用分光光度法[15-16],與光譜采集同步進行。

1.2 光譜數據預處理

為消除原始光譜數據(R)受到雜散光、儀器噪聲、基線漂移等因素的干擾,本研究對原始光譜數據進行光譜曲線的一階微分(RFD)及二階微分(RSD)預處理,以降低基線漂移,突出光譜特征[6]。除此之外,剔除位于兩端(350～399 nm和2 401～2 500 nm)存在高頻噪聲的數據[7],并將得到的光譜數據(400～2 400 nm),按5 nm的間隔重采樣,以降低光譜數據的冗余[9-10]。

1.3 光譜指數構建與篩選

選擇3個二維光譜指數(IDI、ISRI、INDI)和3個三維光譜指數(ITBI1、ITBI2及ITBI3),對R、RFD和RSD分別構建全波段光譜指數[8-10],計算如式(1～6)。用皮爾遜相關系數法(PCC)分析光譜指數與葉綠素面密度間的相關系數(r),并以最大相關系數(rmax)為標準對構建的光譜指數進行變量篩選,減少冗余變量,消除共線性[20]。

IDI(Rλ1,Rλ2)=Rλ1-Rλ2,

(1)

ISRI(Rλ1,Rλ2)=Rλ1/Rλ2,

(2)

INDI(Rλ1,Rλ2)=(Rλ1-Rλ2)/(Rλ1+Rλ2),

(3)

ITBI1(Rλ1,Rλ2,Rλ3)=(Rλ1-Rλ2)/(Rλ2+Rλ3),

(4)

ITBI2(Rλ1,Rλ2,Rλ3)=(Rλ1-Rλ2)/(Rλ1+Rλ3),

(5)

ITBI3(Rλ1,Rλ2,Rλ3)=Rλ1/(Rλ2+Rλ3)。

(6)

式中：Rλ1、Rλ2和Rλ3為400～2 400 nm中任意波長下的反射率。

1.4 模型的建立

海洋捕食者算法(MPA)是一種仿生智能算法[14-15],與其他優化算法相比,MPA能有效的提升優化效果,能高效地獲取更優的解。MPA在初始化種群后根據捕食者和獵物的不同速度比進行優化,過程分為3個階段。

(1)勘探階段。此階段I∈(0,Imax/3),獵物移動速度與捕食者移動速度之比較高,該階段為全局搜索階段,模型為：

(7)

式中：Ei代表捕食者種群;Pi為獵物種群;捕食者和獵物種群中個體數均為n;Si為運動的步長;NR表示布朗運動,服從正態分布;D為常數項,通常取D=0.5;N∈[0,1],服從均勻分布;I為當前迭代次數;Imax為最大迭代次數;?表示逐項相乘。

(2)勘探轉向開發的過度階段。此階段I∈(Imax/3,2Imax/3),獵物移動速度與捕食者移動速度相近。捕食者通過布朗運動搜索獵物,獵物通過萊維運動更新自身的位置。此階段兼顧全局搜索與局部搜索尋優,使用公式(8)進行局部尋優,使用公式(9)進行全局搜索。

局部尋優模型：

(8)

式中：NL為服從萊維分布的隨機數的向量,表示萊維運動。

全局搜索模型：

(9)

式中：CF表示步長的自適應參數,計算公式為：

(10)

(3)開發階段。該階段I∈(2Imax/3,Imax),獵物移動速度與捕食者移動速度之比較低,捕食者主要進行局部尋優,數學模型為：

(11)

此外,在每次迭代結束后,MPA通過利用魚群聚集裝置效應(FADs)使捕食者進行更長的跳躍,以避免陷入局部最優。數學模型：

(12)

f1(x)=(w1,x)+b1,

(13)

f2(x)=(w2,x)+b2,

(14)

式中：x∈Rn,w1和w2表示權重,b1和b2表示偏置,(.,.)表示內積。

通過求解如下問題即可得到最優決策函數：

(15)

(16)

式中：‖·‖表示L2范數,C1、C2是懲罰參數,ε1、ε2>0是常數,ξ、θ是松弛向量,e是m×1維單位列向量。

通過引入拉格朗日乘子δ、η,根據Karush-Kuhn-Tucker(KKT)條件得到如下對偶形式：

(17)

(18)

式中：G=[Ae],h1=y-ε1e和h2=y+ε2e,式(5)和(6)的最優解分別為：

(19)

(20)

求解出f1(x)和f2(x)后,決策函數如下：

(21)

在非線性情況下,通過引入核函數映射到高維空間進行求解,即輸入矩陣A∈Rm×n變為K(A,AT),式(1)和(2)變為：

f1(x)=K(xT,AT)w1+b1,

(22)

f2(x)=K(xT,AT)w2+b2。

(23)

相應的決策函數與式(21)類似。

應用MPA-TSVR對紫丁香葉片葉綠素面密度反演模型構建。

(1)模型輸入參數及數據集劃分

把以不同預處理方法和不同維度為標準篩選出的最優光譜指數分別作為輸入參數構建模型。將160個樣本隨機劃分,選取80%的數據作為訓練集(128個樣本)用于模型訓練,其余20%作為測試集(32個樣本)用于模型性能評估,數據集劃分如表1所示。本研究中所有試驗均采用此標準劃分數據集。

表1 紫丁香葉片葉綠素面密度實測數據集劃分統計

(2)模型內置參數優化

TSVR的懲罰因子C1、C2和核參數g1、g2對模型預測精度影響較大。用MPA對TSVR中C1、C2和g1、g2進行優化,以訓練集的均方根誤差(RMSE)為適應度函數,循環迭代計算個體適應度值,直到找到滿足I

圖1 MPA-TSVR反演紫丁香葉片葉綠素面密度流程圖

(3)模型評價指標

選擇決定系數(R2)和均方根誤差(RMSE)作為模型精度評價指標[21],R2越趨近于1,RMSE越趨近于0,模型預測精度越高、穩定性越強、整體效果越好。

2 結果與分析

2.1 葉片光譜預處理結果

紫丁香葉片原始(R)和預處理(RFD、RSD)后的160條光譜曲線如圖2所示。R在700、1 400和1 900 nm處有明顯的吸收谷,其中700 nm處的吸收谷與葉綠素面密度有關。相比R,經RFD處理后的光譜曲線在一定程度上消除了基線漂移,凸顯R的部分肩峰,但在凸顯特征峰的同時噪音也有所增加。經RSD處理后的光譜曲線,反射率值明顯減小,特征峰和噪音較RFD更為明顯。說明利用微分求導預處理光譜,在增強光譜特征的同時也會增強噪音。

圖2 原始光譜和預處理后光譜

2.2 光譜指數PCC分析及最優光譜指數

2.2.1 一維光譜指數與葉綠素面密度的相關性

一維光譜指數與紫丁香葉片葉綠素面密度相關系數矩陣圖如圖3所示,rmax達到0.874 0。不同光譜反射率與葉綠素面密度之間的敏感波長主要集中在500～650 nm以及700～800 nm區間,不同預處理下光譜反射率與葉綠素面密度間的rmax均在0.8以上。經RFD處理之后,葉綠素面密度敏感波段范圍縮小,敏感波長更為精準,葉綠素面密度敏感度不高的范圍被削弱,而RSD處理后的光譜更加突出此現象,只有幾個較強波長的r較高,其余的波長的r都相對較低,甚至接近于零。提取不同預處理下rmax對應的反射率作為最優一維光譜指數用以后續建模,具體見表2。

圖3 一維光譜指數與葉綠素面密度相關矩陣

表2 最優光譜指數

2.2.2 二維光譜指數與葉綠素面密度的相關性

二維光譜指數與紫丁香葉片葉綠素面密度相關系數矩陣圖如圖4所示。不同的預處理方法和不同的二維光譜指數均會對r產生影響。

圖4 二維光譜指數與葉綠素面密度相關矩陣

總體來看,與一維光譜指數相比,二維光譜指數的r有所提升,說明二維光譜指數通過波段組合能夠加強波段間的相互關系,有效提高光譜信息的利用率。二維光譜指數與一維光譜指數的敏感波長范圍基本一致,主要集中在可見光區域,但存在最優光譜指數INDI(RFD744,RFD1399)第2個波長位置在1 399 nm,說明在近紅外區也有與葉綠素相關性的信息。不同光譜指數,R構建的二維光譜指數rmax普遍偏低,均小于0.883 5;RFD構建的二維光譜指數rmax均最高,其中RFD-ISRI(圖4(e))具有最好的效果,rmax達0.891 1,其次為RSD構建的二維光譜指數,說明RFD和RSD均能有效提高光譜指數與葉綠素面密度之間的r。選取各個矩陣中rmax對應的指數作為最優二維光譜指數用以后續建模,具體見表2。

2.2.3 三維光譜指數與葉綠素面密度的相關性

三維光譜指數與紫丁香葉片葉綠素面密度相關系數矩陣圖如圖5～7所示。不同的預處理方法和不同的三維光譜指數均會對r產生影響。三維光譜指數的rmax均高于0.882 2,比同一預處理下一維和二維光譜指數高,說明將組成指數的波長擴展至3個能夠將更多的光譜信息提取出來,對提升葉綠素面密度的敏感度具有一定作用。R處理下的不同三維光譜指數rmax均最低,不同于最優預處理為RFD的一維和二維光譜指數,RSD是構建三維光譜指數的最優預處理方法,其中RSD-ITBI3(圖7中的(g)、(h)、(i))的rmax最高,達0.901 5,說明三維光譜指數有利于降低RSD預處理在凸顯特征時所產生的噪音,能夠有效提高特征篩選的靈敏度。

圖5 三維光譜指數ITBI1與葉綠素面密度相關矩陣

圖6 三維光譜指數ITBI2與葉綠素面密度相關矩陣

圖7 三維光譜指數ITBI3與葉綠素面密度相關矩陣

研究發現,三維光譜指數ITBI1(圖5)構建的最優指數均在可見光(700～750 nm)范圍內;ITBI2和ITBI3中,除了R-ITBI2和RSD-ITBI3,每個指數中至少有1個波長靠近1 450 nm,結合二維光譜指數來看,說明近紅外區存在與葉綠素面密度相關的信息。選取各個矩陣中rmax對應的指數作為最優三維光譜指數用以后續建模,具體見表2。

2.3 模型預測結果

經PCC篩選出的最優光譜指數如表2所示。將不同預處理及不同維度下的最優光譜指數作為輸入參數,分別用TSVR和MPA-TSVR對紫丁香葉片葉綠素面密度進行反演。

不同預處理下最優光譜指數分別建立TSVR和MPA-TSVR紫丁香葉片葉綠素面密度反演模型,結果見表3。MPA-TSVR構建的紫丁香葉片葉綠素面密度反演模型精度均高于TSVR。同一建模方法,不同預處理下構建模型的精度由高到低順序為RSD、RFD、R,說明對R進行預處理是有必要的,RSD是最佳預處理方法,RSD-MPA-TSVR精度最高,測試集R2為0.906 0,RMSE為3.882 7。使用PCC篩選得到的RFD最優光譜指數比RSD最優光譜指數rmax高,而建模結果顯示RSD是最佳預處理方法,說明非線性模型的r和R2之間沒有必然聯系。因此在構建非線性模型時,輸入參數的好壞不能用r進行評價,而要對預處理方法和所構建的光譜指數進行排列組合,對比模型預測精度才能確定最佳輸入參數。

表3 不同預處理下最優光譜指數構建TSVR和MPA-TSVR模型反演葉綠素面密度結果

不同維度下最優光譜指數分別建立TSVR和MPA-TSVR紫丁香葉片葉綠素面密度反演模型,結果見表4。不同輸入參數,同樣是MPA-TSVR模型精度高于TSVR,可以看出不論輸入參數如何變化,MPA-TSVR模型精度均高于TSVR,說明MPA優化算法能夠找到TSVR中懲罰參數和核參數的最優值,避免TSVR在運算過程中陷入局部最優。同一建模方法,不同維度構建模型精度由高到低順序為：三維、二維、一維,說明光譜指數維度增加,其對葉綠素面密度的響應度也會增加。二維光譜指數中,INDI-MPA-TSVR精度最高,測試集R2為0.882 0,RMSE為4.350 2;三維光譜指數中,ITBI3-MPA-TSVR精度最高,測試集R2為0.911 0,RMSE為3.776 3。

表4 不同維度下最優光譜指數構建TSVR和MPA-TSVR反演葉綠素面密度結果

綜合所有模型來看,ITBI3為最優輸入參數,其構建的TSVR和MPA-TSVR模型參數對比如表5所示。

表5 TSVR初始參數與MPA-TSVR優化參數

用最優輸入參數ITBI3分別構建PLSR、BPNN和SVR模型。其中,PLSR主成分數為2,ITBI3-PLSR測試集R2為0.897 1,RMSE為4.061 8;BPNN隱含層節點數設置為7,ITBI3-BPNN測試集R2為0.899 6,RMSE為4.012 9;SVR選用徑向基核函數,懲罰參數C和核函數g分別為210和2-6,ITBI3-SVR測試集R2為0.901 7,RMSE為3.970 8,精度均低于ITBI3-MPA-TSVR。不同模型反演紫丁香葉片葉綠素面密度預測值與實測值的對比如圖8所示。ITBI3-MPA-TSVR預測值與實測值的散點較其他模型更靠近Y-預測值=Y這條直線,預測精度更高。因此,ITBI3-MPA-TSVR對葉綠素面密度的預測能力、擬合程度和準確性均優于本研究中所有模型。

圖8 不同模型反演紫丁香葉片葉綠素面密度測試集預測結果

3 結論

本研究以內蒙古農業大學東校區內紫丁香為研究對象,采集160個樣品的光譜數據和對應的葉片葉綠素面密度,在一維和二維光譜指數的基礎上,首次在葉綠素面密度預測方面引入三維光譜指數,構建了基于R、RFD和RSD光譜400～2 400 nm范圍內任意波長下不同維度的光譜指數。用PCC分析了光譜指數與葉綠素面密度間的相關性,并篩選出最優光譜指數。首次將MPA-TSVR模型引入植物生理指標預測研究,構建了基于不同預處理和不同維度下的最優指數和MPA-TSVR紫丁香葉片葉綠素面密度預測模型,最優輸入參數分別建立PLSR、BPNN和SVR預測模型,優選出最佳紫丁香葉片葉綠素面密度預測模型,得出以下結論：

(1)最優光譜指數中,一維光譜指數RFD744、二維光譜指數ISRI(RFD704,RFD738)和三維光譜指數ITBI3(RSD714,RSD745,RSD700)的r最大,分別為0.874 0和0.891 1為0.901 5,不同預處理,基于RFD構建的一維和二維光譜指數與葉綠素面密度相關性最高,RSD構建的三維光譜指數與葉綠素面密度相關性最高;不同維度,三維光譜指數相比于一維和二維光譜指數對葉綠素面密度更加敏感。

(2)同一建模方法,不同預處理下最優光譜指數建模精度由高到低的順序為RSD、RFD、R。基于RSD-MPA-TSVR構建的葉綠素面密度反演模型具有最佳預測效果,其R2可達0.906 0,RMSE僅有3.882 7;同一預處理方法,不同維度最優光譜指數建模精度由高到低的順序為：三維、二維、一維。二維光譜指數中,INDI-MPA-TSVR反演效果最好,測試集R2為0.882 0,RMSE為4.350 2;三維光譜指數中,ITBI3-MPA-TSVR反演效果最好,測試集R2為0.911 0,RMSE為3.776 3。不論輸入參數如何,MPA-TSVR的建模精度均高于TSVR。

(3)對于紫丁香葉片葉綠素面密度反演模型而言,MPA-TSVR相比于PLSR、BPNN和SVR模型更穩定。ITBI3-MPA-TSVR是紫丁香葉片葉綠素面密度最優反演模型。本研究光譜采集未在活體上進行,是為同步葉綠素面密度測定,減少誤差。后續光譜數據可在活體采集,通過構建光譜指數輸入ITBI3-MPA-TSVR模型即可實現葉綠素面密度的無損檢測。

通過ITBI3-MPA-TSVR實現了紫丁香葉片葉綠素面密度的反演,三維光譜指數和MPA-TSVR模型的應用對城市綠化植物其他生化組分的定量預測具有一定的參考價值。全波段構建三維光譜指數既減少了光譜無關信息,又因其在二維光譜指數的基礎上增加第3個波段,提升了信息的包含度,能夠很大程度上突出特征信息。但是由于光譜信息十分復雜,預處理方法和光譜指數的選擇在很大程度上影響光譜指數對葉綠素面密度的響應,最優光譜指數的篩選方法直接對模型精度產生影響。因此,對預處理方法、光譜指數以及最優光譜指數篩選方法的選擇有待進一步研究。