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

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

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

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

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í)驗(yàn)室化學(xué)分析數(shù)據(jù)、光譜數(shù)據(jù)、小型冠層試驗(yàn)的化學(xué)分析數(shù)據(jù)及冠層建模數(shù)據(jù)。本研究選擇ACCP中91組道格拉斯冷杉幼苗期新鮮葉片光譜數(shù)據(jù)及相應(yīng)的實(shí)驗(yàn)室化學(xué)分析數(shù)據(jù)，采集于1992年11月。光譜數(shù)據(jù)為400～2 498 nm處的吸光度[log(1/反射率)]，光譜間隔為2 nm，分辨率為10 nm。在采集光譜數(shù)據(jù)后的6 h，采集葉片測(cè)量葉片氮及葉綠素含量。采用常規(guī)的實(shí)驗(yàn)室方法，對(duì)每棵樹(shù)幼苗期的葉片樣本進(jìn)行干燥、研磨和全氮、全葉綠素(a+b)分析。樣品經(jīng)硫酸-氧化汞催化劑(Perstorp分析)在塊狀消化器中消化后，用Alpkem連續(xù)流自動(dòng)分析儀測(cè)定總氮,氮含量的單位為單位干重百分比。

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

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

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

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

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

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

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

1.3 建模方法與精度評(píng)價(jià)

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

2 結(jié)果與分析

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

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

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

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

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

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

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

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

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

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

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

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

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

2.3 LASSO算法

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

表3 LASSO篩選的波段組合

2.4 模型預(yù)測(cè)結(jié)果與對(duì)比

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

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

3種變量篩選方法和3種建模方法的組合中，RMSE在0.19～0.38，R2在0.6～0.89，其中，隨機(jī)森林算法對(duì)相關(guān)系數(shù)法和LASS0算法篩選出的Log(1./R)′形式變量的反演誤差最小，RMSE在0.19～0.20，能夠解釋氮含量的變化，但該算法對(duì)篩選出的其他形式的變量反演誤差較大，表現(xiàn)并不穩(wěn)定;不論采用何種建模方法，對(duì)于相關(guān)系數(shù)法，其log(1/R)′形式的變量反演效果較好；對(duì)于連續(xù)投影法，其篩選出的R′變量，反演誤差較??；對(duì)于LASSO算法，其篩選出的R'和log(1/R)′形式的變量反演效果要明顯優(yōu)于其余兩種變量形式。同時(shí)，分別采用偏最小二乘，隨機(jī)森林、支持向量機(jī)進(jìn)行全波段的建模，反演結(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可以看出兩種導(dǎo)數(shù)形式估算出的氮含量與實(shí)測(cè)氮含量1:1的對(duì)應(yīng)關(guān)系更好，散點(diǎn)基本上沿著y=x線(xiàn)分布，具有較好的線(xiàn)性關(guān)系。對(duì)于R及l(fā)og(1/R)形式的反演，三種建模方法的結(jié)果與實(shí)測(cè)值相比普遍偏高，大部分點(diǎn)分布在y=x線(xiàn)以上。兩種導(dǎo)數(shù)形式的反演結(jié)果較為均勻緊致地分布在y=x線(xiàn)兩側(cè)，反演結(jié)果與實(shí)測(cè)數(shù)據(jù)更接近。從實(shí)測(cè)數(shù)據(jù)與反演結(jié)果的擬合圖可以看出，相關(guān)系數(shù)法&log(1/R)′以及LASSO&log(1/R)′中實(shí)測(cè)數(shù)據(jù)與反演結(jié)果的擬合函數(shù)具有較大的斜率，較小的截距，且R2較大，RMSE較小，具有較高的反演精度(圖3)。

圖3 實(shí)測(cè)氮含量與不同方法反演氮含量散點(diǎn)圖

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

4 結(jié)論

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

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

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

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

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