基于無人機(jī)多光譜遙感的春玉米葉面積指數(shù)和地上部生物量估算模型比較研究

樊鴻葉, 李姚姚, 盧憲菊, 顧生浩, 郭新宇, 劉玉華

(1.河北農(nóng)業(yè)大學(xué)農(nóng)學(xué)院， 河北 保定 071000； 2.北京農(nóng)業(yè)信息技術(shù)研究中心， 數(shù)字植物北京重點(diǎn)實(shí)驗(yàn)室， 北京 100097； 3.國(guó)家農(nóng)業(yè)信息化工程技術(shù)研究中心， 北京 100097； 4.中國(guó)農(nóng)業(yè)科學(xué)院作物科學(xué)研究所， 北京 100081)

無人機(jī)遙感具有時(shí)空分辨率高、成本低、靈活性和可重復(fù)性強(qiáng)等優(yōu)勢(shì)[1]，能夠通過搭載的傳感器獲取作物冠層反射的電磁波信息，進(jìn)而提取與作物長(zhǎng)勢(shì)相關(guān)的參量，為獲取田塊尺度上即時(shí)、無損、可靠的作物長(zhǎng)勢(shì)信息提供了一種重要手段[2]。地上部生物量與葉面積指數(shù)(leaf area index，LAI)是反映作物生長(zhǎng)狀況的重要農(nóng)藝參數(shù)[3-4]，且二者間存在密切的聯(lián)系[5]。在作物不同的生長(zhǎng)發(fā)育階段，由于作物冠層結(jié)構(gòu)、葉片形態(tài)以及生理生態(tài)特征的變化，作物冠層在圖像中表現(xiàn)出的光譜特征亦有不同[6]。植被指數(shù)(vegetation index，VI)是兩個(gè)或多個(gè)波段的反射率經(jīng)過線性或非線性組合運(yùn)算來增強(qiáng)植被信息，以削弱環(huán)境背景對(duì)植被光譜特征的干擾[7]，如歸一化植被指數(shù)(normalized difference vegetation index，NDVI)[8]、增強(qiáng)植被指數(shù)(enhanced vegetation index，EVI)[9-10]和差值植被指數(shù)(difference vegetation index，DVI)[11]等。

目前，光譜植被指數(shù)法已經(jīng)作為一種快速獲取地表信息的手段，被國(guó)內(nèi)外學(xué)者用于作物長(zhǎng)勢(shì)監(jiān)測(cè)研究。楊貴軍等[12]通過無人機(jī)搭載多傳感器，實(shí)現(xiàn)了對(duì)作物產(chǎn)量、LAI、冠層溫度等多種信息的監(jiān)測(cè)。王亞杰[13]基于無人機(jī)獲取的多光譜數(shù)據(jù)，系統(tǒng)比較了玉米葉面積指數(shù)的不同監(jiān)測(cè)方法，結(jié)果表明，在不同水分處理下基于EVI構(gòu)建的一元線性模型能夠較好的預(yù)測(cè)玉米灌漿期和成熟期的LAI。高林等[14]基于RVI、NDVI和DVI等采用經(jīng)驗(yàn)?zāi)Ｐ头?gòu)建了大豆LAI反演模型，結(jié)果表明，在鼓粒期基于NDVI構(gòu)建的大豆LAI效果最佳(R2=0.829，RMSE=0.301)。孫詩(shī)睿等[15]基于無人機(jī)多光譜影像，通過對(duì)比不同植被指數(shù)建模方法分析冬小麥葉面積指數(shù)的反演精度，發(fā)現(xiàn)基于多植被指數(shù)的隨機(jī)森林預(yù)測(cè)具有更好的擬合效果(R2=0.882，RMSE=1.218)。鄒楠等[16]研究結(jié)果表明，種植密度在6、12和18萬(wàn)株·hm-2條件下，春玉米R(shí)VI、DVI和DVI分別與LAI的相關(guān)性最好，R2分別為0.762、0.691和0.648，可用于監(jiān)測(cè)春玉米LAI。肖武等[17]利用玉米灌漿期地上部生物量與各種植被指數(shù)建立一元曲線回歸模型，結(jié)果表明，植被指數(shù)NDVI和GNDVI效果較好(R2>0.68)。莊東英等[18]研究表明，RVI比NDVI對(duì)冬小麥生物量的估算具有明顯優(yōu)勢(shì)。

已有研究大都基于光譜參數(shù)對(duì)作物L(fēng)AI和地上部生物量進(jìn)行模擬估算，但鮮有研究對(duì)作物不同品種LAI和地上部生物量的光譜估算模型進(jìn)行探討，本研究以此為切入點(diǎn)，以無人機(jī)搭載多光譜相機(jī)所獲多光譜影像數(shù)據(jù)為基礎(chǔ)，選取7種常用的植被指數(shù)與玉米LAI和地上部生生物量進(jìn)行相關(guān)性分析；綜合4種基于光譜植被指數(shù)的模擬方法(指數(shù)、線性、對(duì)數(shù)、冪)，構(gòu)建2個(gè)玉米品種在吐絲期和灌漿期的LAI和地上部生物量估算最優(yōu)模型，為進(jìn)一步提高春玉米長(zhǎng)勢(shì)監(jiān)測(cè)的精度提供依據(jù)和方法。

1 材料與方法

1.1 試驗(yàn)材料

供試材料為緊湊型鄭單958(Zhengdan 958，ZD958)和半緊湊型先玉335(Xianyu 335，XY335)，由河南豫研種子科技公司和山東登海先鋒種業(yè)有限公司提供，施用化肥為尿素(總氮≥46.4%)、過磷酸鈣(有效磷≥12%，水溶性磷≥7.0%，硫≥8.0%)和氯化鉀(K2O≥60%)，分別由山東潤(rùn)銀生物化工股份有限公司、江蘇省威力磷復(fù)肥有限公司和中化化肥有限公司生產(chǎn)。所用無人機(jī)為大疆M600 Pro，最大載重5 kg，最大續(xù)航時(shí)間30 min。搭載于無人機(jī)上的多光譜傳感器為Rededge-MX(MicaSense，USA)，該相機(jī)有5個(gè)波段，相機(jī)焦距為5.5 mm，視場(chǎng)角為47.2°，圖像分辨率為1 280 pixels×960 pixels，該相機(jī)配備了光強(qiáng)傳感器及標(biāo)定板，其中光強(qiáng)傳感器可校正航拍過程中外界光線的變化對(duì)光譜影像造成的影響，標(biāo)定板具有固定的反射率，可以利用標(biāo)定板進(jìn)行輻射定標(biāo)。

1.2 試驗(yàn)處理設(shè)計(jì)

試驗(yàn)于2018—2019年4—10月在位于吉林省公主嶺市的中國(guó)農(nóng)科院公主嶺試驗(yàn)基地(43°11′—44°9′N、124°02′—125°18′E)進(jìn)行。該地區(qū)屬于溫帶半濕潤(rùn)大陸性季風(fēng)氣候，長(zhǎng)年平均降雨量為595 mm， 土壤類型主要為黑土，土壤耕層pH 6.0，有機(jī)質(zhì)含量29.68 g·kg-1、全氮1.28 g·kg-1、速效磷28.16 mg·kg-1、速效鉀187.74 mg·kg-1。設(shè)置4個(gè)氮肥處理，分別為無肥N0(0 kg·hm-2)、低肥N1(150 kg·hm-2全部一次性基肥)、中肥N2(225 kg·hm-2，基施150 kg·hm-2，拔節(jié)追施75 kg·hm-2)和髙肥N3(300 kg·hm-2，基施150 kg·hm-2，拔節(jié)期、吐絲期各施75 kg·hm-2)，氮、磷、鉀分別為尿素、過磷酸鈣(P2O542.5 kg·hm-2)、氯化鉀(K2O 42.5 kg·hm-2)，種植密度為67 500株·hm-2，行距0.65 m，株距0.22 m。除施肥處理不同外，各處理田間管理方式均保持一致。

1.3 玉米地上部生物量及LAI值獲取

于2018年的7月21日(吐絲期)、8月25日(灌漿期)和2019年的7月27日(吐絲期)、8月27日(灌漿期)獲取田間數(shù)據(jù)，將每個(gè)處理的試驗(yàn)小區(qū)所選取3株具有代表性的玉米植株帶回實(shí)驗(yàn)室后，先放置在烘箱中105 ℃殺青30 min，然后80 ℃烘干48 h以上，直至恒重，最后稱取地上部生物量。

測(cè)量植株所有葉片的長(zhǎng)和最大葉寬，根據(jù)式(1)計(jì)算葉面積指數(shù)。

(1)

式中，0.75表示玉米葉面積校正系數(shù)；ρ表示種植密度，株·m-2；N表示取樣株數(shù)，株；n為總?cè)~片數(shù)；li和wi分別表示第i片葉的葉長(zhǎng)和最大葉寬，m。

1.4 無人機(jī)多光譜影像數(shù)據(jù)獲取

試驗(yàn)采用大疆M600 Pro六旋翼小型電動(dòng)無人機(jī)，搭載多光譜傳感器Micasense Rededge-MX，上述田間取樣當(dāng)天選擇在晴朗、無風(fēng)或微風(fēng)、無云的天氣狀況下在10:00—14:00采集多光譜影像數(shù)據(jù)，無人機(jī)飛行高度為40 m，航向?yàn)?0%，旁向?yàn)?5%。

1.5 多光譜影像數(shù)據(jù)處理

將所獲取的四個(gè)生育期多光譜影像使用軟件Pix4Dmapper進(jìn)行拼接，獲得覆蓋整個(gè)試驗(yàn)區(qū)的多光譜影像，拼接前先人為篩選剔除航帶外以及無人機(jī)調(diào)頭時(shí)的影像，并輸入每個(gè)波段對(duì)應(yīng)標(biāo)定板的反射率，最終可獲得五張?jiān)囼?yàn)區(qū)單波段的正射圖像。將拼接好的影像導(dǎo)入ENVI 5.1中進(jìn)行幾何校正，以高清數(shù)碼影像為參考影像，均勻選取20個(gè)點(diǎn)校正所獲無人機(jī)多光譜影像，誤差控制在0.5個(gè)像元以內(nèi)，主要是為消除因無人機(jī)飛行高度及姿態(tài)不穩(wěn)定對(duì)影像造成的影響。然后將校正后的影像進(jìn)行波段組合，將單波段的影像疊加到一起，最后進(jìn)行波段運(yùn)算，根據(jù)需求計(jì)算相關(guān)植被指數(shù)。將計(jì)算所得的植被指數(shù)影像在ENVI 5.1軟件中進(jìn)行參數(shù)提取，每個(gè)小區(qū)繪制3個(gè)均勻大小的感興趣區(qū)，繪制時(shí)避開小區(qū)邊界四個(gè)方向各一行玉米，以排除邊行長(zhǎng)勢(shì)不均勻植株，進(jìn)而提取出每個(gè)感興趣區(qū)的平均光譜植被指數(shù)，每個(gè)感興趣區(qū)域作為一次重復(fù)，共三次。

1.6 光譜植被指數(shù)計(jì)算

本研究選取了與作物長(zhǎng)勢(shì)相關(guān)的7個(gè)植被指數(shù)：NDVI[6]、GNDVI[17]、RVI[19]、GOSAVI[20]、VIopt[21]、EVI[8]、DVI[11]，計(jì)算分別公式如下。

NDVI=(RNIR-RRED)/(RNIR+RRED)

(2)

GNDVI=(RNIR-RGREEN)/(RNIR+RGREEN)

(3)

RVI=RNIR/RRED

(4)

GOSAVI=(1+0.16)(RNIR-RGREEN)/(RNIR+RGREEN+0.16)

(5)

VIopt=1.45(RNIR×RNIR+1)/(RRED+0.45)

(6)

EVI=2.5(RNIR-RRED)/(RNIR+6RRED-7.5RRLUE+1)

(7)

DVI=RNIR-RRED

(8)

式中,RNIR、RRED、RGREEN和RBLUE分別表示近紅外波段、紅光波段、綠光波段和藍(lán)光波段的反射率。

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

采用SPSS22.0軟件進(jìn)行試驗(yàn)數(shù)據(jù)統(tǒng)計(jì)分析，用Duncan法進(jìn)行處理間的多重比較(P<0.05)，用Sigmaplot12.5軟件制圖。

1.8 模型精度評(píng)價(jià)指標(biāo)

2019年數(shù)據(jù)用于建模，2018年數(shù)據(jù)用于外部獨(dú)立驗(yàn)證。采用決定系數(shù)(coefficient of determination,R2)、均方根誤差(root mean square error,RMSE)和歸一化均方根誤差(normal root mean square error，NRMSE)評(píng)價(jià)和檢驗(yàn)?zāi)Ｐ汀2表示回歸方程與實(shí)測(cè)值的擬合程度，其值越接近1效果越好；RMSE用來衡量模擬值與實(shí)測(cè)值偏差，該值越小模擬效果越好；NRMSE<10%模擬效果非常好，10%～20%效果較好，20%～30%效果可以接受，>30%效果差[22]。本研究評(píng)價(jià)不同估算模型時(shí)綜合考慮R2、RMSE和NRMSE三個(gè)指標(biāo)[23-24]。

(9)

(10)

(11)

2 結(jié)果與分析

2.1 不同施氮量對(duì)玉米冠層植被指數(shù)的影響

由圖1可見，兩品種在N1、N2和N3處理下各植被指數(shù)均顯著高于N0處理。ZD958在吐絲期和灌漿期各植被指數(shù)N2處理顯著高于N1處理，說明N2 處理下的玉米長(zhǎng)勢(shì)比N1處理好，在吐絲期N3處理下的GOSAVI、DVI和VIopt顯著高于N1，而在灌漿期卻顯著低于N1(圖1A、C)。XY335在吐絲期各植被指數(shù)除DVI、VIopt和RVI外在各施肥處理間無顯著性差異(圖1B)，灌漿期的GOSAVI、EVI和RVI在N2和N3處理下顯著高于N1處理(圖1D)。

2.2 植被指數(shù)與LAI和地上部生物量間的相關(guān)性分析

LAI和地上部生物量與植被指數(shù)間的相關(guān)性分析結(jié)果(表2)表明，在吐絲期，兩品種各植被指數(shù)與LAI和地上部生物量均呈極顯著相關(guān)(P<0.01)；在灌漿期，DVI與ZD958的LAI和地上部生物量相關(guān)性不顯著(P>0.05)，而與XY335的LAI和地上部生物量顯著相關(guān)(P<0.05)；灌漿期VIopt與ZD958的LAI和地上部生物量顯著相關(guān)(P<0.05)，與XY335的LAI和地上部生物量極顯著相關(guān)；其余植被指數(shù)與兩品種在吐絲期和灌漿期的LAI和地上部生物量均呈極顯著相關(guān)(P<0.01)，選擇與LAI和地上部生物量分別達(dá)極顯著水平的植被指數(shù)來構(gòu)建估算模型。

表2 玉米LAI模型精度對(duì)比Table 2 Comparison of the accuracy of different models for maize LAI (n =12)

2.3 模型比較分析

2.3.1玉米LAI模型比較分析 從表3中，選取R2最大，RMSE相對(duì)較小的植被指數(shù)和擬合方程(線性、對(duì)數(shù)、指數(shù)和冪)，如果出現(xiàn)R2相近則選擇RMSE最小的擬合方程和植被指數(shù)。在吐絲期，GOSAVI與ZD958的LAI的擬合冪函數(shù)效果最好(R2=0.946，RMSE=0.252)；RVI與XY335的LAI擬合的指數(shù)函數(shù)有較高的R2和最小RMSE(R2=0.947，RMSE=0.179)。在灌漿期，NDVI與ZD958的LAI擬合的冪函數(shù)R2相對(duì)較大，RMSE最小(R2=0.879，RMSE=0.419)；GNDVI與XY335的LAI擬合的冪函數(shù)效果相對(duì)較好(R2=0.861，RMSE=0.261)。

表3 玉米地上部生物量模型精度對(duì)比Table 3 Comparison of accuracy of different models for maize above-ground biomass (n =12)

2.3.2玉米地上部生物量模型比較分析 吐絲期冪函數(shù)對(duì)兩品種的地上部生物量估算效果較好，在灌漿期指數(shù)函數(shù)對(duì)兩品種的地上部生物量估算效果好(表4)。在吐絲期，GNDVI與ZD958地上部生物量擬合的冪函數(shù)效果最好，R2和RMSE分別為0.867和0.488 t·hm-2，GOSAVI與XY335地上部生物量擬合的冪函數(shù)估算效果最好(R2=0.791，RMSE=0.463 t·hm-2)。在灌漿期，EVI與ZD958地上部生物量擬合的指數(shù)函數(shù)建模效果較好(R2=0.772，RMSE=1.477 t·hm-2)，GNDVI與XY335地上部生物量擬合的指數(shù)效果最好(R2=0.852，RMSE=2.0132 t·hm-2)。

表4 LAI和地上部生物量模型驗(yàn)證結(jié)果Table 4 Model validation for LAI and above-ground biomass (n =12)

2.4 LAI和地上部生物量模型驗(yàn)證

為檢驗(yàn)?zāi)Ｐ偷目煽啃裕?018實(shí)測(cè)數(shù)據(jù)對(duì)所選模型進(jìn)行驗(yàn)證。用RMSE和NRMSE兩個(gè)指標(biāo)來評(píng)價(jià)模型精度，從表5可以看出，在吐絲期，與ZD958擬合的LAI和地上部生物量模型估算精度較高，RMSE和NRMSE分別為1.103、17.91%和1.485 t·hm-2、16.04%，XY335擬合的模型精度相對(duì)較差，RMSE和NRMSE分別為1.369、25.08%和1.730 t·hm-2、20.11%。在灌漿期，ZD958的LAI和地上部生物量估算模型精度較高，RMSE和NRMSE分別為0.619，12.22%和2.482 t·hm-2，15.22%；XY335的LAI估算精度較差(RMSE=1.008，NRMSE=25.74%)，地上部生物量估算模型效果較好(RMSE=2.119 t·hm-2，NRMSE=12.26%)。

XY335的LAI和地上部生物量模型估算效果比ZD958的差，但總體精度可以接受，能夠用來估算LAI和地上部生物量。

3 討論

吐絲期到灌漿期是玉米營(yíng)養(yǎng)吸收和物質(zhì)積累的高峰期[25]，從吐絲期開始玉米植株的養(yǎng)分向籽粒運(yùn)轉(zhuǎn)[26]，及時(shí)、準(zhǔn)確的獲取吐絲期和灌漿期玉米長(zhǎng)勢(shì)參量并據(jù)此進(jìn)行調(diào)優(yōu)栽培對(duì)提高玉米產(chǎn)量具有重要意義。本研究結(jié)果表明，對(duì)于ZD958而言，在吐絲期和灌漿期各植被指數(shù)在N2處理下顯著高于N1處理，說明可能N2處理下玉米葉片葉綠素含量高于N1處理，馮宗會(huì)[27]在冬小麥上的研究表明，孕穗期施氮量270 kg·hm-2處理NDVI顯著高于施氮量150 kg·hm-2處理。經(jīng)相關(guān)性分析發(fā)現(xiàn)，同一植被指數(shù)與不同玉米品種間的LAI和地上部生物量的相關(guān)性差異較大，表明不同玉米品種的冠層光譜特征存在差異，這與張俊華等[28]研究結(jié)果一致，可能是由于品種間株型不同導(dǎo)致冠層光譜差異。多種植被指數(shù)與兩玉米品種的LAI和地上部生物量均呈極顯著或顯著性相關(guān)水平，這與趙巧麗等[29]的研究結(jié)果基本吻合。

本研究結(jié)果表明，吐絲期冪函數(shù)對(duì)ZD958的LAI估算效果最好，指數(shù)函數(shù)對(duì)XY335的LAI估算效果好，灌漿期冪函數(shù)對(duì)兩品種的LAI估算效果最佳 。吐絲期冪函數(shù)對(duì)兩品種的地上部生物量估算效果較好，灌漿期指數(shù)函數(shù)對(duì)兩品種的地上部生物量估算效果好。吐絲期ZD958 的GNDVI對(duì)地上部生物量的估算效果最佳，而灌漿期GNDVI對(duì)XY335的估算效果最優(yōu)，Shanahan等[30]研究也發(fā)現(xiàn)，抽雄期到灌漿期，GNDVI的估產(chǎn)效應(yīng)優(yōu)于NDVI和TSAVI。灌漿期EVI對(duì)ZD958的地上部生物量估算精度較高，韓文霆等[31]發(fā)現(xiàn)，單個(gè)生育期EVI對(duì)產(chǎn)量的估算效果最佳。XY335在灌漿期GNDVI與LAI和地上部生物量估算效果均較好，說明在該時(shí)期GNDVI對(duì)XY335的LAI和地上部生物量的估算可通用。本研究一方面明確了不同玉米品種長(zhǎng)勢(shì)參量的光譜響應(yīng)特征和最優(yōu)預(yù)測(cè)模型，另一方面為下一步開展品種分辨率的光譜響應(yīng)機(jī)理解析奠定了基礎(chǔ)。本研究仍存在一些不足，玉米封壟后覆蓋度較高，易導(dǎo)致植被指數(shù)出現(xiàn)“飽和”的現(xiàn)象，在今后的研究中應(yīng)在多品種氮肥長(zhǎng)期定位試驗(yàn)的基礎(chǔ)上，研究光譜、空間和時(shí)相特征的融合方法，面向田間精確栽培的管理需求，構(gòu)建適應(yīng)于不同品種、不同生育期的玉米LAI和地上部生物量的遙感精確反演模型。

致謝：感謝中國(guó)農(nóng)業(yè)科學(xué)院作物科學(xué)研究所謝瑞芝研究員為本研究的順利開展提供了完備的實(shí)驗(yàn)條件。