冷杉葉片氮含量高光譜反演技術(shù)研究

宋雪蓮,王志偉,張文，2，張威，丁磊磊，柳嘉佳，阮璽睿,王普昶

(1.貴州省農(nóng)業(yè)科學院草業(yè)研究所，貴州 貴陽 550006；2.貴州陽光草業(yè)科技有限責任公司，貴州 貴陽 550006；3.貴州省水利水電勘測設(shè)計研究院有限公司,貴州 貴陽 550002；4.貴州省農(nóng)業(yè)科學院,貴州 貴陽 550006)

氮是植物生長的基本營養(yǎng)元素，是陸地生態(tài)系統(tǒng)植物生長的主要限制因子，與植物的光合作用和細胞生長分類等重要生理活動相關(guān)[1]。缺氮會影響植物的光合作用能力，導致葉片生理及形態(tài)結(jié)構(gòu)的相應變化，因而引起葉片光譜反射特性的變化[2]。傳統(tǒng)的化學方法能夠精確地監(jiān)測出植物氮含量，但具有高損耗、復雜、時滯性等缺點，不能進行全面、大范圍、快速的營養(yǎng)診斷。高光譜光譜分辨率高，能夠獲得連續(xù)的光譜信息，已成為監(jiān)測植物氮含量的有效手段。

目前有不少學者采用敏感波段篩選和植被指數(shù)進行高光譜植物氮含量檢測的相關(guān)研究。朱艷等[4]研究表明，單波段光譜在610 nm和680 nm處的水稻冠層反射率與葉片氮含量具有較好的相關(guān)性，提出采用回歸系數(shù)來提高水稻葉片氮含量估測的準確性。劉冰峰等[5]指出720 nm處的反射光譜以及550，720 nm處的一階光譜導數(shù)是夏玉米氮含量的敏感參數(shù)。Tarpley等[6]研究認為，利用紅邊位置和短波近紅外波段光譜反射率比值可預測棉花葉片氮濃度。Pacheco-Labrador等[7]采用歸一化指數(shù)NDIS和三波段指數(shù)TBIs有效估算Holm oak葉片氮含量；Abdel-Rahman等[8]采用一階導數(shù)的SR指數(shù)(波段743，1 316 nm)(波段743，1 317 nm)(波段741，1 323 nm)估算的葉片氮含量的決定系數(shù)在0.75左右；Yao等[9]采用歸一化指數(shù)NDSI和比值指數(shù)RSI一階導數(shù)對葉片氮的累積量進行了估算，其決定系數(shù)高達0.81；Ullah 等[10]利用MERIS數(shù)據(jù)采用植被指數(shù)(NDVI，NBDI，SAVI，TSAVI，REIP，MTCI)以及波段深度參數(shù)對草地生物量及氮含量進行了估算和對比。Sanches等[11]采用反射率，吸光度及其衍生形式和偏最小二乘方法對植物氮磷鉀進行估算，結(jié)果表明一階導數(shù)形式的估算精度最高。另有學者著重從算法角度進行了植物氮含量高光譜反演技術(shù)研究，Zhang等[12]利用單變量線性回歸、逐步多元線性回歸、偏最小二乘對濕地植物蘆葦?shù)窟M行估算，并采用留一交叉驗證比較3種模型的精度，認為氮的敏感波段集中在紅波段和綠波段，單變量線性回歸對665 nm和680 nm處的歸一化植被指數(shù)效果最好，3種方法中逐步線性回歸的精度最高。Yao等[13]采用原始反射率、一階導數(shù)、植被指數(shù)等變量利用SMLR、PLSR、ANNS、SVMs等反演了冠層葉片氮含量，結(jié)果表明一階導數(shù)的支持向量機方法精度更高。Wang 等[3]驗證了多核支持向量機在估算小麥葉片氮含量的有效性，并與多元線性回歸，偏最小二乘，人工神經(jīng)網(wǎng)絡(luò)，單核支持向量機進行了對比?？傮w來說，現(xiàn)有的關(guān)于高光譜植物氮含量監(jiān)測的大部分研究集中在植被指數(shù)反演氮含量以及敏感波段篩選方面，關(guān)于建模算法的研究較少[3]。且已有的研究或集中于建模之前敏感特征的篩選或集中于建模算法的比較，本文嘗試將不同特征提取方法與建模方法相結(jié)合，通過不同方法間的組合，篩選出一套完整的能有效估算葉片氮含量的流程。

1 材料和方法

1.1 光譜數(shù)據(jù)與氮含量數(shù)據(jù)獲取

數(shù)據(jù)集采用Accelerated canopy chemistry program(ACCP)(https://doi.org/10.3334/ORNLDAAC)[14]，該數(shù)據(jù)集旨在研究不同生態(tài)系統(tǒng)中植被冠層氮含量及木質(zhì)素含量遙感反演的理論基礎(chǔ)，包含野外樣品的實驗室化學分析數(shù)據(jù)、光譜數(shù)據(jù)、小型冠層試驗的化學分析數(shù)據(jù)及冠層建模數(shù)據(jù)。本研究選擇ACCP中91組道格拉斯冷杉幼苗期新鮮葉片光譜數(shù)據(jù)及相應的實驗室化學分析數(shù)據(jù)，采集于1992年11月。光譜數(shù)據(jù)為400～2 498 nm處的吸光度[log(1/反射率)]，光譜間隔為2 nm，分辨率為10 nm。在采集光譜數(shù)據(jù)后的6 h，采集葉片測量葉片氮及葉綠素含量。采用常規(guī)的實驗室方法，對每棵樹幼苗期的葉片樣本進行干燥、研磨和全氮、全葉綠素(a+b)分析。樣品經(jīng)硫酸-氧化汞催化劑(Perstorp分析)在塊狀消化器中消化后，用Alpkem連續(xù)流自動分析儀測定總氮,氮含量的單位為單位干重百分比。

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

利用Matlab對試驗數(shù)據(jù)進行了異常數(shù)據(jù)剔除、計算了原始反射率、一階導數(shù)R′、反射率倒數(shù)的對數(shù)log(1/R)及其一階導數(shù)log(1/R)′，以及這些參數(shù)與氮含量的相關(guān)系數(shù)，采用連續(xù)投影算法、LASSO、相關(guān)系數(shù)篩選以上4種參數(shù)形式的特征波段，并采用隨機森林、偏最小二乘以及支持向量機建立氮含量反演模型。其中，相關(guān)系數(shù)、連續(xù)投影算法篩選波段組合，偏最小二乘建模在Matlab中完成，LASSO算法篩選變量組合、隨機森林和支持向量機建模在R中完成。

1.2.1 相關(guān)系數(shù)法 相關(guān)系數(shù)法是最簡單直接的能夠幫助理解特征與相應變量之間關(guān)系的方法，衡量的是變量之間的線性相關(guān)性，其值絕對值越大表示相關(guān)性越強。很多研究通過選擇最大相關(guān)性的變量作為特征變量進行建模研究。λ1-λ2，λ1/λ2表示兩變量間的差異性，其與因變量間的相關(guān)系數(shù)表示兩變量組合對因變量的解釋量，選取相關(guān)系數(shù)最大的兩個波長即為特征波長。

1.2.2 連續(xù)投影算法 連續(xù)投影算法(SPA)是一種矢量空間共線性最小化的前向變量選擇算法，能夠有效從全波段中提取特征波段，消除原始光譜矩陣中的冗余信息，降低模型的復雜度，在波長選取中取得了較好的效果[15]。連續(xù)投影算法是一種前向循環(huán)選擇算法，設(shè)樣本集M和波段數(shù)K組成一個光譜矩陣XMXK，分別記Xk(0)和N為初始的迭代向量和需要提取的波段個數(shù)。從一個波長開始，每次循環(huán)計算它在未選入的波長上的投影，將投影向量最大的波長引入波長組合，直到循環(huán)N次，每一個新入選的波段，都與前一個線性關(guān)系最小。通過循環(huán)會得到N×K對波段組合，對每一對XK(0)和N所決定的組合分別建立多元回歸模型，并用預測均方根誤差RMSE來決定所建模型的優(yōu)劣，最小的RMSE對應的Xk(0)為最佳的波段組合。

1.2.4 隨機森林 隨機森林是一種并行式集成學習法，以決策樹為基礎(chǔ)，并在決策樹訓練過程中引入隨機屬性選擇。它通過自助法重采樣技術(shù)，從原始訓練樣本集N中有放回地重復隨機抽取K個樣本生成新的訓練樣本集合，然后根據(jù)自助樣本集生成K個分類樹組成隨機森林[17-18]。

1.2.5 支持向量機 支持向量機是一種基于統(tǒng)計學習理論的機器學習算法，它建立在VC(Vapnik-Chervonenkis Dimension)維理論和結(jié)構(gòu)風險最小原理基礎(chǔ)上，能較好地解決小樣本、非線性、高維數(shù)和局部極小點等實際問題，同時能獲得較好的泛化能力[19]。

1.2.6 偏最小二乘 偏最小二乘回歸是一種新型的多變量回歸分析方法，可以實現(xiàn)回歸建模、簡化數(shù)據(jù)結(jié)果和分析兩組變量見得相關(guān)性，給多元數(shù)據(jù)統(tǒng)計分析帶來極大便利，它能在回歸建模過程中采用數(shù)據(jù)將為、信息綜合和篩選技術(shù)，提取對系統(tǒng)最佳解釋能力的新綜合成分[20]。

1.3 建模方法與精度評價

分別采用偏最小二乘，隨機森林與支持向量機建立葉片氮含量反演模型，并采用決定系數(shù)R2、均方根誤差來評價各個模型的精度。R2表示相關(guān)密切程度，RMSE用來衡量估測值與真實值之間的偏差程度。R2越大，RMSE的絕對值越小，表示模擬結(jié)果精度更高。

2 結(jié)果與分析

2.1 不同光譜形式與氮含量的相關(guān)性分析

原始反射率R與氮含量呈負相關(guān)關(guān)系，在500～730 nm波段的相關(guān)系數(shù)超過-0.6，在550、702 nm左右相關(guān)性達到最大，且差異顯著(P<0.05)，在550 nm左右相關(guān)性最大是由于葉綠素對綠光的反射作用，702 nm位于紅邊區(qū)域，有相關(guān)研究發(fā)現(xiàn)紅邊波段與植物葉綠素含量、生物量等參數(shù)存在顯著相關(guān)性；氮含量與原始光譜的相關(guān)性在近紅波段較小(表1)。

一階導數(shù)R′能夠有效消除葉面積變化的影響，反射率的一階微分與氮含量的相關(guān)性變化較為復雜，在490～544、556～710和736～800 nm超過0.6，其中R′在530、574、626，694 nm達到最大，近紅波段的相關(guān)性較小，但在1 200、1 672、2 164、2 310 nm相關(guān)性分別達到峰值0.6以上(圖1)。

根據(jù)工程的復雜程度，對設(shè)計資質(zhì)提出具體要求，如對提灌站、壓力輸水管道、100 m以上的深井、規(guī)模較大的建筑物等技術(shù)含量較高的工程必須聘請有相應資質(zhì)的單位進行設(shè)計，從源頭抓好工程質(zhì)量。

圖1 R，R′，log(1/R)，log(1/R)′與氮含量的相關(guān)性

Log(1/R)相關(guān)性曲線與原始反射率的相關(guān)性曲線線型一致，但相關(guān)性相反。log(1/R)′與氮含量的相關(guān)性變化也較為復雜，在490～546，620～670，672～694，702～804 nm相關(guān)性在0.6以上，在紅外波段1 204，1 672，2 158 nm的相關(guān)性形成峰值，達到0.6以上。

應用λ1-λ2，λ1/λ2與氮含量的相關(guān)系數(shù)篩選出的波段組合如表1所示，R與log(1/R)篩選出的波段集中在可見光波段，兩種導數(shù)形式篩選出的波段集中在近紅波段；原始波段反射率及其變換形式篩選出的波段組合都包含522，526 nm；導數(shù)形式篩選出的波段組合均包含2 158，2 044和2 070 nm 3個波段，雖然這幾個單獨的波段與氮含量的相關(guān)性并不強，但其波段組合與氮含量的相關(guān)性達到0.86以上，說明該波段組合包含了氮含量的大部分信息。

表1 相關(guān)系數(shù)法篩選的波段組合

2.2 連續(xù)投影算法篩選敏感波段特征篩選

ACCP葉片數(shù)據(jù)集共91組鮮葉片數(shù)據(jù)，隨機選取其中65組數(shù)據(jù)進行連續(xù)投影算法。連續(xù)投影方法將65組數(shù)據(jù)分為訓練集和校正集，以校正集的預測均方根誤差來確定最佳的光譜變量總數(shù)，均方根誤差和決定系數(shù)來評價最佳的波段變量組合。對以上4種形式的光譜參數(shù)分別采用連續(xù)投影算法，確定氮含量反演的最佳光譜變量(表2)。

表2 連續(xù)投影法篩選的波段組合

原始反射率所選的波段在波段范圍內(nèi)分布較為均勻(圖2)，位于可見光的波段有6個，2個波段位于紅邊區(qū)域，位于波峰或波谷處的波段有6個。一階導數(shù)只選擇了3個變量，其中一個位于可見光波段。Log(1/R)投影算法選擇了兩個紅光波段，其余均為近紅波段。log(1/R)′僅選擇了3個紅外波段。校正數(shù)據(jù)集RMSE最低的是Log(1/R)形式，RMSE最高的是一階導數(shù)形式。

圖2 連續(xù)投影算法篩選出的波段位置

2.3 LASSO算法

采用Lasso算法對4種光譜形式進行變量篩選，各光譜形式篩選出的變量數(shù)均比連續(xù)投影篩選出的變量數(shù)多。原始光譜反射率篩選出22個波段，其中有12個波段位于可見光區(qū)。Log(1/R)篩選出的變量數(shù)為16個，其中9個可見光波段。兩種導數(shù)形式均篩選出23個變量數(shù)，一階導數(shù)篩選出的可見光波段8個，Log(1/R)′篩選出的可見光波段6個(表3)。

表3 LASSO篩選的波段組合

2.4 模型預測結(jié)果與對比

分別將相關(guān)系數(shù)篩選法、連續(xù)投影算法、LASSO算法篩選出的波段作為偏最小二乘、隨機森林、支持向量機的輸入來估算葉片氮含量，以65組數(shù)據(jù)為建模數(shù)據(jù)，26組數(shù)據(jù)為驗證數(shù)據(jù)。不同方法組合的反演結(jié)果如表4及圖3所示。

表4 不同方法組合的反演結(jié)果

3種變量篩選方法和3種建模方法的組合中，RMSE在0.19～0.38，R2在0.6～0.89，其中，隨機森林算法對相關(guān)系數(shù)法和LASS0算法篩選出的Log(1./R)′形式變量的反演誤差最小，RMSE在0.19～0.20，能夠解釋氮含量的變化，但該算法對篩選出的其他形式的變量反演誤差較大，表現(xiàn)并不穩(wěn)定;不論采用何種建模方法，對于相關(guān)系數(shù)法，其log(1/R)′形式的變量反演效果較好；對于連續(xù)投影法，其篩選出的R′變量，反演誤差較小；對于LASSO算法，其篩選出的R'和log(1/R)′形式的變量反演效果要明顯優(yōu)于其余兩種變量形式。同時，分別采用偏最小二乘，隨機森林、支持向量機進行全波段的建模，反演結(jié)果分別為R2=0.605、RMSE=0.4542，R2=0.659、RMSE=0.3352，R2=0.730 4、RMSE=0.348 4。表明，相關(guān)系數(shù)&Log(1./R)′，LASSO&Log(1./R)′，SPA&R′能夠有效篩選出冷杉葉片氮含量的敏感波段組合，且前兩種方法組合效果更優(yōu)。

3 討論

3.1 不同變量篩選方法和建模方法的組合效果比較

從圖3可以看出兩種導數(shù)形式估算出的氮含量與實測氮含量1:1的對應關(guān)系更好，散點基本上沿著y=x線分布，具有較好的線性關(guān)系。對于R及l(fā)og(1/R)形式的反演，三種建模方法的結(jié)果與實測值相比普遍偏高，大部分點分布在y=x線以上。兩種導數(shù)形式的反演結(jié)果較為均勻緊致地分布在y=x線兩側(cè)，反演結(jié)果與實測數(shù)據(jù)更接近。從實測數(shù)據(jù)與反演結(jié)果的擬合圖可以看出，相關(guān)系數(shù)法&log(1/R)′以及LASSO&log(1/R)′中實測數(shù)據(jù)與反演結(jié)果的擬合函數(shù)具有較大的斜率，較小的截距，且R2較大，RMSE較小，具有較高的反演精度(圖3)。

圖3 實測氮含量與不同方法反演氮含量散點圖

3.2 3種篩選方法的敏感波段分析

4 結(jié)論

采用相關(guān)系數(shù)法、連續(xù)投影法、LASSO算法對4種形式的光譜變量R，R′，log(1/R)，log(1/R)′進行敏感波段篩選，對篩選后的敏感波段分別采用偏最小二乘、隨機森林、支持向量機對氮含量進行建模反演，得出以下結(jié)論：

1)兩種導數(shù)形式變量的反演誤差最小。

2)相關(guān)系數(shù)&Log(1/R)′，LASSO&Log(1/R)′能夠有效篩選出冷杉葉片氮含量的敏感波段組合，無論采取何種建模方法，其估算效果在幾種篩選方法組合中最好，R2>0.84，RMSE在0.19～0.24，估算效果明顯優(yōu)于采用3種建模方法進行全波段建模。

3)隨機森林算法對相關(guān)系數(shù)法和LASS0算法篩選出的Log(1/R)′形式變量的反演誤差最小，但對其他形式篩選出的變量反演結(jié)果的誤差變化范圍較大，表現(xiàn)并不穩(wěn)定。

4)3種變量篩選方法篩選出的R形式的變量與前人研究相符，連續(xù)投影算法能篩選出更多與葉片其他化學含量相關(guān)的波段，證明了3種變量篩選方法的有效性。