尼格拉·吐爾遜， 蘇磊·乃比， 高 健， 沈江龍， 鄭江華*， 余丹林

1. 新疆大學(xué)資源與環(huán)境科學(xué)學(xué)院， 新疆 烏魯木齊 830046 2. 新疆大學(xué)數(shù)學(xué)與系統(tǒng)科學(xué)學(xué)院， 新疆 烏魯木齊 830046 3. 新疆林業(yè)科學(xué)院現(xiàn)代林業(yè)研究所， 新疆 烏魯木齊 830063 4. Department of Earth and Environmental Studies, Montclair State University, New Jersey 07043， USA

引 言

葉綠素含量是棗樹光合作用能力、 生長狀況、 營養(yǎng)狀況的指示劑[1]， 通常采用便攜式葉綠素計(SPAD-502)測定植物葉片SPAD值來直接表征植物葉綠素含量的相對大小， 但使用過程中需要將葉片反復(fù)插入測量， 難以用于大范圍的葉綠素檢測， 研究表明SPAD值能與無損、 無污染、 價格低的高光譜遙感數(shù)據(jù)準(zhǔn)確對應(yīng)， 近年來成為葉綠素含量估算的強(qiáng)有力工具[2]。

20世紀(jì)90年代， Pinar[3]和 Blackburn[4]等研究得到葉綠素與高光譜波段之間的相關(guān)關(guān)系。 隨后， 許多學(xué)者在高光譜估算葉綠素模型方面開展了大量的研究， 杜華強(qiáng)基于高斯核函數(shù)變換的偏最小二成回歸模型建立了馬尾松針葉葉綠素含量與光譜反射率及9個特征參數(shù)之間的預(yù)測模型， 其精度遠(yuǎn)大于傳統(tǒng)線性回歸模型[5]。 劉京等用實例證實了支持向量機(jī)具有更好的SPAD值反演效果[6]。 馮海寬等基于特征光譜參數(shù)， 利用隨機(jī)森林模型較好的估算了蘋果葉片葉綠素[7]。 李曉麗等證實了最小二乘支持向量機(jī)(least sqares support vector regression, LS-SVR)在植物參數(shù)估算方面具有較好效果[8]。

上述研究常選用相關(guān)系數(shù)較高的波段或者植被指數(shù)建模使得變量選擇隨機(jī)、 單一、 缺乏定量化， 模型估算能力低下。 本文通過CP統(tǒng)計量在預(yù)測角度選擇重要性較高的自變量， 篩選重要程度高的特征波段(characteristic band， CB)。 其次， 以往的高光譜估算應(yīng)用廣泛的多元線性回歸(multiple linear regression, MLR)、 支持向量機(jī)(support vector regression, SVR)、 LS-SVR模型較多， 并沒有考慮到地理位置可能對葉片SPAD值產(chǎn)生的影響。 2017年Hwang和Shim對于LSSVM模型加入地理位置影響， 提出了地理加權(quán)最小二乘支持向量機(jī)模型[9](GWLS-SVM)， 證實了其估計精度顯著高于傳統(tǒng)的GWR、 LS-SVR模型。 本研究對于GWLS-SVR模型是否適用于葉綠素含量估算， 能否在紅棗樹葉片葉綠素估算中得到較好效果還需要進(jìn)一步的驗證。 用CP統(tǒng)計量選擇特征波段， 計算若羌紅棗樹葉片SPAD值的全局莫蘭指數(shù)， 分析紅棗樹葉片SPAD值分布是否與空間位置有關(guān)， 再運(yùn)用GWLS-SVR模型， 將建模結(jié)果與傳統(tǒng)模型進(jìn)行對比分析， 檢驗并比較模型的擬合效果。

1 實驗部分

1.1 研究區(qū)概況及采樣點布設(shè)

研究區(qū)位于中國新疆若羌縣， 范圍在東經(jīng)87°00′—89°0′、 北緯38°40′—39°30′之間， 屬暖溫帶大陸性荒漠干旱氣候， 是新疆名牌產(chǎn)品“若羌紅棗”種植區(qū)[10]。 于若羌紅棗果實成熟期2019年9月28日—10月2日采樣， 為了保證實驗結(jié)果的全面性和精確性， 在去除野外數(shù)據(jù)異常值后最終保留均勻覆蓋若羌縣的67個棗林樣點， 在預(yù)先設(shè)計的棗林內(nèi)確定代表性棗樹1～3棵進(jìn)行數(shù)據(jù)采集， 再通過手持GPS記錄地理位置信息， 共采集219條紅棗樹葉片高光譜數(shù)據(jù)和219個棗樹葉片SPAD值數(shù)據(jù)， 研究區(qū)位置和采樣樣點地理位置分布情況如圖1。

圖1 研究區(qū)位置和采樣點分布圖

1.2 紅棗樹葉片光譜測定與處理

紅棗樹的葉片光譜反射率在晴朗無風(fēng)無云條件下于北京時間11:00—17:00使用PSR-3500便攜式地物光譜儀在野外測定， 波段范圍是350～2 500 nm， 每隔1 nm輸出一個數(shù)據(jù)， 一共2 151個光譜通道。 在選擇的代表性1～3顆棗樹上、 中、 下層各隨機(jī)采集3片葉片。 為減少誤差， 每次光譜測定之前均進(jìn)行白板標(biāo)定， 同時用干燥紙巾去除葉片表面浮塵， 測量時將葉片鋪平放置在反射率近似為零的黑板上， 將光纖探頭垂直固定于葉片上方約5 mm， 每個葉片樣本避開葉脈重復(fù)測量3次， 取光譜曲線的算術(shù)平均值作為該樣點的原始葉片光譜反射率。 為減少噪聲影響， 剔除1 050～2 500 nm噪聲較大波段， 并利用Origin軟件平滑去噪[11]。 另外， 導(dǎo)數(shù)光譜可以反映植被中生化物質(zhì)的吸收引起的波形變化還能夠揭示光譜峰值的內(nèi)在特征進(jìn)而估算植被內(nèi)部葉綠素含量信息[2]。 因此， 對原始光譜反射率(raw reflectance， RR)求光譜一階導(dǎo)數(shù)(first derivative of reflectance， FD)。

1.3 SPAD值測定

使用葉綠素計(SPAD-502Plus, Konica Minoita, Japan)對現(xiàn)場采集的多個棗樹葉片SPAD值進(jìn)行測定， 測量時避開葉脈部分， 從葉柄至葉尖分段隨機(jī)測量3次， 將多個葉片測定結(jié)果取算術(shù)平均值作為該樣點SPAD值。 SPAD值測定時間與葉片光譜測定同步進(jìn)行， 測定位置與葉片光譜保持一致。

1.4 變量選擇

本工作采用CP統(tǒng)計量進(jìn)行變量選擇。CP統(tǒng)計量可以通過預(yù)測的角度選擇重要性較高的自變量。 其原理為由部分變量預(yù)測的均方誤差可能比利用所有變量進(jìn)行預(yù)測的均方誤差更小， 故可以去除重要程度不是很高的變量。 其計算方式如式(1)

(1)

1.5 地理加權(quán)最小二乘支持向量回歸模型(GWLS-SVR)

f(xi，Ui)=ωTφ(xi)+bi

(2)

設(shè)給定x與Ui下的權(quán)重矩陣為Wi， 則可以將回歸模型轉(zhuǎn)化成如式(3)優(yōu)化問題

(3)

其中， C>0為懲罰參數(shù)， wij為用于表示Ui和Uj之間的距離的權(quán)重函數(shù)。

(4)

2 結(jié)果與討論

2.1 不同范圍SPAD值紅棗樹葉片平均光譜反射特征

圖2是不同范圍SPAD值的紅棗樹葉片平均光譜反射率曲線圖。 由圖可知， 不同范圍SPAD值的紅棗樹葉片平均反射率曲線變化趨勢基本相同。 總體上， 350～750 nm波段內(nèi)反射率比750～1 050 nm波段低。 在350～675 nm波段內(nèi)隨著SPAD值的升高， 紅棗樹葉片平均光譜反射率降低， 光譜差異較明顯， 其中， 在500～551 nm波段范圍內(nèi)反射率緩慢上升， 551 nm附近出現(xiàn)反射峰， 675 nm附近出現(xiàn)吸收谷； 675～750 nm處平均光譜反射率隨著波長呈現(xiàn)快速上升趨勢， 750～1 050 nm范圍內(nèi)， 隨著SPAD值的升高， 平均光譜反射率升高。 紅棗樹的長勢狀態(tài)直接決定了SPAD值的大小， SPAD值也會影響紅棗樹葉片的反射率。

圖2 不同范圍SPAD值紅棗樹葉片平均光譜反射特征

2.2 紅棗樹葉片光譜和SPAD值相關(guān)性分析

為了明確紅棗樹葉片SPAD值相對應(yīng)的敏感波段， 將紅棗樹葉片SPAD值和原始光譜、 光譜一階導(dǎo)數(shù)反射率波段做皮爾遜相關(guān)性分析。 由圖3可知， 紅棗樹葉片SPAD值和原始光譜反射率及光譜一階導(dǎo)數(shù)反射率緊密相關(guān)， 且都存在著極顯著相關(guān)。 對原始光譜來說， 在570～620及690～700 nm間達(dá)到相關(guān)系數(shù)峰值， 通過了0.01的顯著性水平， 相關(guān)系數(shù)分別達(dá)到-0.578及-0.561， 此波段范圍受葉綠素吸收的影響， 相關(guān)系數(shù)呈負(fù)相關(guān)， 選擇這兩組波段的原始光譜反射率作為估測棗樹葉片SPAD值的敏感波段區(qū)間。 SPAD值與光譜一階導(dǎo)數(shù)呈正負(fù)相關(guān)， 相關(guān)性極顯著的波段分布在400～750 nm區(qū)間內(nèi)， 最高值出現(xiàn)在688 nm處。 與原始光譜相比， 在492～510， 542～543， 642～652， 657～670和682～692 nm區(qū)間內(nèi)的SPAD相關(guān)性有所提高， 且分別達(dá)到-0.655， -0.662， -0.697， 0.709和-0.749， 也說明了這些波段的光譜反射率與棗樹葉片SPAD值相關(guān)性好， 適合用于敏感波段的挑選。 綜上所述， 紅棗樹葉片反射率光譜做一階導(dǎo)數(shù)處理后與SPAD的相關(guān)性有較顯著的提高。

2.3 變量選擇

結(jié)合圖3， 在原始光譜570～620 nm范圍內(nèi)選擇了相關(guān)性高的581， 590， 595和602 nm波段， 690～700 nm波段范圍內(nèi)選擇695和696 nm共6個特征波段進(jìn)行CP統(tǒng)計量的計算； 基于光譜一階導(dǎo)數(shù)與SPAD值相關(guān)性高低， 在492～510， 542～543， 642～652， 657～670和682～692 nm共5個波段內(nèi)分別選擇相關(guān)性達(dá)到區(qū)間內(nèi)最高的495， 543， 649， 664和688 nm共5個特征波段計算出其不同組合統(tǒng)計量， 表1為波段的相關(guān)系數(shù)表。

圖3 SPAD值與光譜反射率之間的相關(guān)性

表1 波段的相關(guān)系數(shù)表

表2為CP統(tǒng)計量計算結(jié)果表， 考慮到所有變量組合方式數(shù)目較大， 且大部分組合方式的CP統(tǒng)計量都遠(yuǎn)高于表2中的幾種組合方式， 只列出CP統(tǒng)計量值靠前的組合。CP統(tǒng)計量越低， 代表該種變量選擇方式重要性程度越高。 且由表2可知， 原始光譜選擇在570～620和690～700 nm范圍內(nèi)分別選擇595和696 nm時CP統(tǒng)計量絕對值最低， 因此將595與696 nm原始光譜重要程度最高的兩個變量作為建模的特征波段。 光譜一階導(dǎo)數(shù)變換后688 nm波段CP統(tǒng)計量絕對值最低， 因此光譜一階導(dǎo)數(shù)的特征波段定為688 nm。 原始光譜特征波段696 nm和光譜一階導(dǎo)數(shù)特征波段688 nm都處于紅邊波段[13]， 說明紅邊與植被的各種理化參數(shù)是緊密相關(guān)的， 是描述植物色素狀態(tài)和健康狀況的重要的指示波段。

表2 特征波段組合及CP統(tǒng)計量計算結(jié)果

不難發(fā)現(xiàn)， 對于同一個區(qū)間的波段組合總有單波段的CP統(tǒng)計量低于多波段組合的CP統(tǒng)計量。 說明相近波段組合建模會使得誤差增大， 這可能是相近波段之間較強(qiáng)共線性造成的， 故每個敏感波段區(qū)間只選取一個波段進(jìn)行建模是合理的。

2.4 模型的構(gòu)建與檢驗

運(yùn)用CP統(tǒng)計量選出的3個特征波段以及實測葉片SPAD值建立MLR， SVR和GWLS-SVR模型。 相比較而言， GWLS-SVR主要的優(yōu)勢是變量系數(shù)隨著地理位置而變化， 具有較強(qiáng)的靈活性。 為了明確紅棗樹葉片SPAD值分布是否與地理位置有關(guān)， 對其進(jìn)行Moran’s Ⅰ的計算結(jié)果為0.125 8(p<0.1)， 呈空間正相關(guān)， 說明棗樹葉片SPAD值的分布有顯著的空間聚集性， 適合運(yùn)用GWLS-SVR模型來建模。

原始光譜(RR)與光譜一階導(dǎo)數(shù)(FD)分別基于MLR， SVR以及GWLS-SVR擬合的MSE與R2如圖4所示。 從建模效果來看， 基于原始光譜建立的三種模型中， MLR與SVR的R2低于0.8， MSE也較高， 說明這兩種模型的穩(wěn)定性較差， 預(yù)測效果不理想； GWLS-SVR的R2為0.915， MSE低至3.679， 表明GWLS-SVR的穩(wěn)定性及估算能力優(yōu)于MLR與SVR模型。 光譜一階導(dǎo)數(shù)變換后的三種模型精度較原始光譜均有所提升， 且MSE整體上都有所降低， 表明數(shù)據(jù)變換后模型的穩(wěn)定性和精度有了一定的提高； 而GWLS-SVR在光譜一階導(dǎo)變換后均顯著優(yōu)于其余兩種模型， 模型的R2提高到了0.975， MLR與SVR的MSE均比GWLS-SVR高約20倍， 綜合上述可得GWLS-SVR模型不僅擬合精度高， 其估計偏差與方差綜合看來均低于其余兩個模型。

圖4 MLR， SVR和GWLS-SVR對實測值與預(yù)測值間的擬合圖

從擬合效果來看， GWLS-SVR在原始光譜與光譜一階導(dǎo)數(shù)的擬合曲線比起其他兩種模型真實值與預(yù)測值均勻分布在1∶1直線周圍， 表明GWLS-SVR的擬合效果較好， 且在光譜一階導(dǎo)數(shù)變換后的擬合效果更佳。

為了檢驗三種模型的擬合效果差異， 對于原始樣本利用Bootstrap再抽樣方法進(jìn)行100次有放回隨機(jī)抽樣， 每次抽取67個樣本。 之后對于隨機(jī)生成的樣本利用上述三種模型分別建模計算100組MSE的均值(mean of MSE)和方差(varionce of MSE)， 并利用配對t檢驗進(jìn)一步比較原始光譜和光譜一階導(dǎo)數(shù)基于GWLS-SVR與其他兩個模型的MSE之間的差異的顯著性。

100次Bootstrap再抽樣并基于三種模型建模后的MSE箱線圖如圖5所示。

由表3及圖5可以看出， 整體上GWLS-SVR的100組MSE在原始光譜與光譜一階導(dǎo)變換后均為最低， 且波動也比較小， 說明100次Bootstrap再抽樣后GWLS-SVR相比于傳統(tǒng)的MLR及SVR模型預(yù)測的精度較高且發(fā)揮比較穩(wěn)定， 其中基于光譜一階導(dǎo)數(shù)建立的GWLS-SVR模型的MSE最小， 且波動最小， 說明基于光譜一階導(dǎo)數(shù)建立的GWLS-SVR模型的模型預(yù)測精度最佳且穩(wěn)定。

表3 Bootstrap再抽樣結(jié)果

圖5 100次建模的MSE箱線圖

為了進(jìn)一步檢驗GWLS-SVR的MSE是否顯著小于其他兩個模型的MSE， 以下分別做GWLS-SVR與其他兩個模型的單邊配對t檢驗。 所得T統(tǒng)計量與p值如表4所示。

表4 t檢驗結(jié)果

由表4可見， 4組單邊配對t檢驗的t統(tǒng)計量的絕對值都比較大， 且p值均非常接近0。 所以GWLS-SVR預(yù)測的MSE小于其他兩個模型的MSE這一假設(shè)在統(tǒng)計學(xué)上是高度顯著的。

3 結(jié) 論

利用野外實測67個樣點， 219條紅棗樹葉片高光譜數(shù)據(jù)和棗樹葉片SPAD值數(shù)據(jù)， 對SPAD值與高光譜波段進(jìn)行相關(guān)性分析、CP統(tǒng)計量特征波段選取、 建立基于特征波段的SPAD值估算模型， 結(jié)果表明：

(1)光譜一階導(dǎo)數(shù)起到了對原始光譜數(shù)據(jù)的去噪、 突出高光譜信息的作用， 尤其是在492～510， 542～543， 642～652， 657～670和682～692 nm區(qū)間內(nèi)明顯提高了與SPAD值的相關(guān)性。

(2)根據(jù)統(tǒng)計量計算發(fā)現(xiàn)： 對于同一個敏感波段區(qū)間的波段組合總有單個波段的CP統(tǒng)計量低于多個波段組合的CP統(tǒng)計量， 臨近分布的波段之間的存在的較強(qiáng)共線性可能導(dǎo)致這些波段的組合誤差增大。

(3)基于實地采樣數(shù)據(jù)進(jìn)行地統(tǒng)計分析若羌縣棗樹SPAD值與地理位置的關(guān)聯(lián)性， 發(fā)現(xiàn)若羌縣存在空間聚集性， 全局莫蘭指數(shù)為0.125 8(p<0.1)， 表明地理加權(quán)最小二乘支持向量機(jī)方法適用于估算若羌縣棗樹葉片SPAD值。

(4)基于光譜一階導(dǎo)數(shù)的特征波段建立的GWLS-SVR模型的估算能力(R2為0.975， MSE為1.082)優(yōu)于基于原始光譜特征波段建立的GWLS-SVR模型(R2為0.915， MSE為3.679)， 且由結(jié)合Bootstrap再抽樣方法與t檢驗的結(jié)果來看， 基于光譜一階導(dǎo)數(shù)的加入地理位置信息的GWLS-SVR模型為最優(yōu)的棗樹葉片SPAD值估算模型， 能夠為快速無損的監(jiān)測紅棗樹生長狀況提供參考。

致謝：感謝若羌縣委的一貫支持， 感謝若羌縣委辦、 縣農(nóng)業(yè)農(nóng)村局、 自然資源局、 交通運(yùn)輸局對本項野外調(diào)查工作的具體幫助。 感謝縣委辦戶亮亮同志對本項工作的協(xié)調(diào)和幫助。