向友珍 安嘉琪 趙 笑 金 琳 李志軍 張富倉

(1.西北農林科技大學旱區(qū)農業(yè)水土工程教育部重點實驗室, 陜西楊凌 712100;2.西北農林科技大學水利與建筑工程學院, 陜西楊凌 712100; 3.西北農林科技大學綜合素質教育學院, 陜西楊凌 712100)

0 引言

大豆是人類和動物飲食中顯著的蛋白質和脂肪酸來源[1],對于眾多亞洲國家,大豆對其飲食文化的影響有著十分重要的作用[2]。因此,確保大豆的高質量和高效率生產具有重要意義。

葉面積指數(shù)(Leaf area index, LAI)和地上部生物量是作物冠層結構的重要參數(shù),不僅可為作物生長提供動態(tài)信息,還決定作物生物物理過程[3-4],同時作為重要的輸入?yún)?shù),在作物生長模型與決策支持系統(tǒng)中發(fā)揮著十分重要的作用[5]。通過對LAI和地上部生物量實測值的分析,能夠進行作物生長監(jiān)測和產量預測,也可以進行田間水肥管理[6-8]。因此,快速、準確地測定地上部分的LAI和地上部生物量具有相當重要的意義。當前獲取作物LAI和地上部生物量的方法有兩種,分別是直接測量法和間接測量法。其中前者的精確度較后者更高,但是直接測量法需要對作物進行破壞性采樣,不僅費時費力并且所取樣本不一定具有代表性,因此該方法有一定的局限性;與之不同的是,間接測量法是將作物光譜信息與田間實測數(shù)據(jù)相結合并通過模型對作物生長指標進行估算的方法,且具有快速、高效等特點。無人機多光譜遙感操作簡易、靈活性高并且成本較低,這些特點都使其成為近年來進行精準農業(yè)研究的一項重要工具以及獲取作物冠層參數(shù)的重要手段[9-10]。邵國敏等[11]通過配備Red Edge多光譜相機的六旋翼無人機,對不同灌溉條件下的玉米LAI進行了估算,發(fā)現(xiàn)在不同的灌溉條件下不同的模型對于估算玉米LAI的效果各不相同,其中在充分灌溉條件下,效果最好的是植被指數(shù)的多元線性回歸模型,而在水分脅迫條件下,表現(xiàn)更好的則是植被指數(shù)的隨機森林(Random forest,RF)模型。劉濤等[12]通過得到不同飛行高度的多光譜數(shù)據(jù)來估計小麥冠層LAI,得到采用偏最小二乘回歸技術反演小麥LAI的精度最高。黃林生等[13]通過使用配備MicaSense Red EdgeTM3光譜相機的DJI M600 Pro六旋翼高性能無人機,對水稻LAI進行了估算,提出結合光譜特征、紋理指數(shù)和植物覆蓋度等多個指標的估算模型具有最佳精度。

作物產量與人民生活水平和國家糧食安全息息相關[14],及時準確地預測糧食產量對農業(yè)管理、糧食政策制定、農業(yè)保險和農業(yè)用水效率評估具有重要的意義。目前,主要通過經驗模型、半經驗模型和物理模型3種方法來進行遙感估產[15-17]。而無人機遙感可獲得比衛(wèi)星和地面遙感更高的時間分辨率、空間分辨率和光譜分辨率影像[18-20],已成為農業(yè)遙感數(shù)據(jù)的重要途徑[21]。無人機遙感的應用對于發(fā)展智慧農業(yè)、精準農業(yè)具有重要的影響。目前,對作物產量的估算已有部分研究。余坤勇等[22]根據(jù)光譜數(shù)據(jù),利用不同的光譜指數(shù)進行玉米產量的預測;STEPANOV等[23]在利用無人機遙感數(shù)據(jù)進行估產的研究中對不同的光譜指數(shù)進行討論,發(fā)現(xiàn)歸一化差異植被指數(shù)可以更好地用于產量估算。王鵬新等[24]在估算小麥產量時用最小二乘法結合無人機高光譜數(shù)據(jù)與不同植被指數(shù),發(fā)現(xiàn)最優(yōu)估算指數(shù)為EVI2(無藍帶增強植被指數(shù))。

過往的研究大多只通過單一植被指數(shù)或光譜反射率分析反演估算的LAI、地上部生物量或者產量,通過不同植被指數(shù)的組合作為輸入變量進行建模的較少。此外,在利用無人機多光譜數(shù)據(jù)預測西北地區(qū)大豆產量的研究較少。因此本文以西北地區(qū)大豆為研究對象,通過無人機獲取大豆的田間多光譜數(shù)據(jù),分析各生育期植被指數(shù)、LAI、地上部生物量和產量之間的關系,采用RF、反向傳播神經網(wǎng)絡(Back propagation neural network, BPNN)和支持向量機(Support vector machine, SVM)3種機器學習方法建立定量模型。并討論不同機器學習方法對大豆LAI、地上部生物量和產量估算模型精度的影響以及不同生育期的光譜信息對產量估算模型的精度影響,以期為更加準確、快速地獲取大豆的生長指標和產量預測提供理論參考。

1 材料與方法

1.1 研究區(qū)概況

試驗于2021年6—9月在西北農林科技大學旱區(qū)農業(yè)水土工程教育部重點實驗室節(jié)水灌溉試驗站(34°20′N,108°24′E,海拔521 m)進行(圖1)。該地區(qū)氣候類型為暖溫帶季風半濕潤氣候,是典型的半濕潤偏(易)旱區(qū)。根據(jù)當?shù)貧庀筚Y料顯示降水年際間差異較大(1995—2020年),平均降水量為561 mm,其中全年的降雨主要發(fā)生在7—9月,平均占年降水量的一半以上(約55.5%)。試驗區(qū)土壤為重壤土,0～100 cm土層的田間持水率(質量分數(shù),下同)為23%～25%,凋萎含水率為8.5%,平均土壤干容積密度為1.58 g/cm3。

圖1 研究區(qū)概況Fig.1 Study area

1.2 試驗設計

本試驗設置4種覆蓋類型:平作秸稈覆蓋(SM)、壟覆地膜溝覆秸稈(SFM)、壟覆地膜溝無覆蓋(FM)和平作無覆蓋方式。設置5個施氮量:0 kg/hm2(N0)、30 kg/hm2(N1)、60 kg/hm2(N2)、90 kg/hm2(N3)和120 kg/hm2(N4),共計20個處理。

設置7個種植密度:1.5×105株/hm2(R0)、2.0×105株/hm2(R1)、2.5×105株/hm2(R2)、3.0×105株/hm2(R3)、3.5×105株/hm2(R4)、4.0×105株/hm2(R5)和4.5×105株/hm2(R6),共計7個處理。

設置各處理小區(qū)時采用隨機排列原則,并設置2個重復,共54個小區(qū),小區(qū)面積為2.4 m×6 m=14.4 m2,試驗區(qū)周圍布設2 m保護帶。每個小區(qū)磷肥和鉀肥保持一致,均為30 kg/hm2。在2021年6月18日按行距50 cm、株距約10 cm人工點播大豆,其他田間生產管理(施藥、除草等)均與當?shù)乇３忠恢隆Ｓ?月30日收獲大豆,并回收地膜。大豆全生育期為苗期(2021年6月18—25日)、分枝期(2021年6月26日—7月18日)、開花結莢期(2021年7月19日—8月17日)、鼓粒期(2021年8月18日—9月10日)與成熟期(2021年9月11—30日)。其中分枝期大豆正處于營養(yǎng)生長階段,是決定整個生育期植株健壯與否、分枝與開花多少的關鍵時期[25];結莢期是大豆生長最旺盛的時期,營養(yǎng)生長與生殖生長并進,是生物量積累的關鍵時期[26];鼓粒期是大豆產量和品質形成的重要生育時期,大豆在鼓粒期進行蛋白質和脂肪的儲藏[27]。而苗期和成熟期大豆各項生長參數(shù)變化不大,因此,本研究選擇分枝期、結莢期和鼓粒期3個生育期獲取無人機影像數(shù)據(jù)。

1.3 無人機影像數(shù)據(jù)獲取與處理

無人機遙感影像獲取時間為2021年7月25日(四節(jié)期,V4)、2021年8月20日(結莢期,R4)和2021年9月5日(鼓粒期,R6)。采用大疆精靈4多光譜版一體化無人機(Phantom4-M,P4M)獲取大豆冠層多光譜影像數(shù)據(jù)。該設備集成1個可見光傳感器通道和5個多光譜傳感器通道(藍光、綠光、紅光、紅邊和近紅外),每次拍攝可獲得6幅影像。且每幅影像具有2.0×106像素,最高飛行速度14 m/s,最大續(xù)航時間27 min,搭配有TimeSync時間同步系統(tǒng),可獲得厘米級定位精度,此外P4M頂部集成光強傳感器,可捕捉太陽輻照度數(shù)據(jù)用于影像的光照補償,排除環(huán)境光對數(shù)據(jù)的干擾,提高不同時段采集數(shù)據(jù)的準確度與一致性。試驗利用DJITerra軟件進行無人機航線規(guī)劃,航向和旁向重疊率均為80%,影像分辨率分別為1.6 cm/像素。使用Pix4D mapper軟件進行影像拼接,可獲得不同飛行高度下的RGB和單波段正射影像。將5個波段的多光譜正射影像導入ENVI軟件進行波段合成,使用快速輻射校正工具將像素DN值轉換為反射率。

1.4 地面實測數(shù)據(jù)獲取

在無人機多光譜數(shù)據(jù)影像采集的當天進行大豆LAI和地上部生物量的數(shù)據(jù)采集。

LAI測定方法:使用美國LI-COR LAI-2200C型植物冠層分析儀測定。在54個試驗小區(qū)內各隨機進行6次大豆LAI數(shù)據(jù)采集,以其均值作為各試驗小區(qū)的實測值。

地上部生物量測定方法:在每個小區(qū)選取6株能代表該小區(qū)長勢的大豆進行破壞性取樣,在105℃殺青30 min,隨后75℃干燥后使用電子天平稱量大豆總質量。最終將測量結果換算成單位面積的干物質量。

產量測定方法:在大豆成熟后,選取每個小區(qū)的中間4行測量其產量,并進行考種工作(隨機選擇10株),籽粒產量按含水率13%計[28]。

1.5 光譜指數(shù)構建與選擇

根據(jù)植被的光譜吸收特性,將遙感影像不同波段反射率之間進行線性或非線性組合形成的無量綱指標參數(shù)稱為植被光譜指數(shù)。該指數(shù)是對地表植被狀況簡單、有效的經驗度量,同時能夠體現(xiàn)植被在可見光、近紅外波段反射與土壤背景之間的差異[10]。植被的LAI、地上部生物量和產量等指標與差異化植被指數(shù)(DVI)、比值植被指數(shù)(RVI)、歸一化植被指數(shù)(NDVI)、增強型植被指數(shù)(EVI)和土壤調劑植被指數(shù)(SAVI)等經驗植被指數(shù)密切相關。因此,選取15種植被光譜指數(shù)構建模型,表達式如表1所示。

表1 經驗植被指數(shù)及計算公式Tab.1 Empirical vegetation spectral indexes and calculation formulas

1.6 模型構建與精度檢驗

基于田間試驗共獲取162組LAI和地上部生物量樣本及54組產量樣本和三生育期的光譜數(shù)據(jù)樣本(共162組)。對162組LAI和地上部生物量樣本進行從小到大排序,同時也將54組產量樣本進行從小到大排序,在這些樣本中隨機選取2/3的樣本作為建模集,剩余1/3的樣本作為驗證集。表2為建模集和驗證集的樣本數(shù)量及產量的統(tǒng)計特征。

表2 大豆LAI、地上部生物量和產量描述性統(tǒng)計Tab.2 Descriptive statistics for leaf area index, above-ground biomass and yield of soybean

基于計算得到的15種光譜指數(shù),分析全生育期光譜指數(shù)與LAI、地上部生物量的相關性,篩選出相關系數(shù)(R)最高的5個植被指數(shù),作為模型的輸入變量。采用SVM、RF、BPNN進行建模,對大豆的LAI和地上部生物量進行回歸預測。同理,將各生育期的光譜指數(shù)與產量數(shù)據(jù)進行相關性分析,在各生育期篩選出與產量相關系數(shù)(R)最高的5個植被指數(shù),作為模型的輸入變量。采用SVM、RF、BPNN進行建模,對大豆的產量進行回歸預測。

1.7 模型方法和評價指標

支持向量機模型的核函數(shù)類型設置為“poly”,根據(jù)最小交叉驗證誤差原則分別確定參數(shù)懲罰系數(shù)C和γ為20和0.02;在RF模型的構建中,經過參數(shù)優(yōu)化和多次訓練后,LAI模型和地上生物量模型中的決策樹數(shù)均設置為600。反向神經網(wǎng)絡模型隱蔽層傳遞函數(shù)設置為“TANSIG”,基于數(shù)值優(yōu)化理論的Levenbeger-Marquardt(Train-LM)算法被用作網(wǎng)絡訓練函數(shù),經過多次訓練,確定中間層的神經元數(shù)量為15個[39]。

采用均方根誤差(RMSE)、決定系數(shù)(R2)與平均相對誤差(MRE)對模型擬合結果進行評價[37],模型的預測精度與決定系數(shù)(R2)成正相關;RMSE與MRE越小,說明模型的性能更穩(wěn)定,預測結果也更為集中。

2 結果與分析

2.1 植被指數(shù)與LAI和地上部生物量相關性分析

利用數(shù)據(jù)集中實測的LAI值分別與15個植被指數(shù)進行相關性分析,計算結果如表3所示。其中,RVI的相關系數(shù)最高,為0.706;由大到小排序分別篩選出5個植被指數(shù)作為LAI估算模型的輸入,分別為RVI、NDVI、GNDVI、NDRE和CIRE。地上部生物量分別與15個植被指數(shù)的相關性分析結果如表3所示,OSAVI與地上部生物量的相關系數(shù)最高,為0.676,同理,選擇OSAVI、SAVI、GOSAVI、EVI和MTCI作為地上部生物量估算模型的輸入。

表3 植被指數(shù)與LAI和地上部生物量相關系數(shù)計算結果Tab.3 Calculation of correlation coefficients between vegetation index and leaf area index and above-ground biomass

2.2 不同生育期植被指數(shù)與產量相關性分析

各生育期15組植被指數(shù)與產量的相關性分析結果如表4所示,其中VIopt在V4期與產量相關性最佳,相關系數(shù)高達0.813;而在R4和R6生育期與產量相關性最佳的植被指數(shù)為CIRE,相關系數(shù)分別為0.849和0.827。與LAI和地上部生物量篩選植被指數(shù)的方法同理,在V4選取植被指數(shù)VIopt、RVI、GCI、TVI和MTVI作為該生育期的產量估算模型的輸入;選取植被指數(shù)CIRE、NDRE、GCI、MTCI和GNDVI作為R4的產量估算模型的輸入;選取植被指數(shù)NDRE、CIRE、GCI、MTCI和GNDVI作為R6的產量估算模型的輸入。

表4 植被指數(shù)與產量相關系數(shù)Tab.4 Correlation coefficient between vegetation index and yield

2.3 大豆LAI和地上部生物量估算模型構建

以2.1節(jié)篩選出的植被指數(shù)組合作為自變量,以大豆LAI作為響應變量,分別采用SVM、RF、BPNN構建大豆全生育期LAI估算模型,從R2、RMSE、MRE 3方面綜合評定模型精度,不同建模方法對于大豆葉面積的預測結果如圖2所示。結果表明:基于SVM構建的大豆LAI估算模型建模集和驗證集的R2分別為0.506和0.593,RMSE分別為1.054 m2/m2和1.006 m2/m2,MRE分別為26.454%和23.563%; 基于RF構建的大豆LAI估算模型建模集和驗證集的R2分別為0.783和0.801,RMSE分別為0.713 m2/m2和0.675 m2/m2,MRE分別為17.693%和18.684%; 基于BPNN構建的大豆LAI估算模型建模集和驗證集的R2分別為0.588和0.708,RMSE分別為0.961 m2/m2和0.756 m2/m2,MRE分別為23.564%和19.741%。

圖2 基于SVM、 RF和BPNN的大豆LAI估算模型Fig.2 Estimation model of soybean LAI based on reverse SVM, RF and BPNN

同理,不同建模方法對于大豆地上部生物量的預測結果如圖3所示。結果表明:基于SVM構建的大豆地上部生物量估算模型建模集和驗證集的R2分別為0.481和0.451,RMSE分別為2 296.584 kg/hm2和 2 213.534 kg/hm2,MRE分別為28.880%和26.649%;基于RF構建的大豆地上部生物量估算模型建模集和驗證集的R2分別為0.765和0.744,RMSE分別為1 550.173 kg/hm2和1 548.140 kg/hm2,MRE分別為21.673%和18.769%; 基于BPNN構建的大豆地上部生物量估算模型建模集和驗證集的R2分別為0.557和0.536,RMSE分別為2 096.711 kg/hm2和1 973.511 kg/hm2,MRE分別為27.127%和25.642%。

圖3 基于SVM、RF和BPNN的大豆地上部生物量估算模型Fig.3 Estimation models of soybean above-ground biomass based on reverse SVM, RF and BPNN

2.4 大豆產量估算模型構建

以2.2節(jié)篩選出的各生育期植被指數(shù)組合作為輸入變量, 以大豆產量作為響應變量,分別采用SVM、RF、BPNN構建大豆各生育期產量估算模型,從R2、RMSE、MRE 3方面綜合評定模型精度,不同建模方法對于大豆產量的預測結果如圖4所示。

圖4 基于SVM、RF和BPNN的大豆產量估算模型Fig.4 Estimation models of soybean grain yield based on reverse SVM, RF and BPNN

結果表明:在V4生育期,基于SVM構建的大豆地上部生物量估算模型建模集和驗證集的R2分別為0.634和0.615,RMSE分別為289.474、388.397 kg/hm2,MRE分別為6.894%和9.735%;基于RF構建的大豆地上部生物量估算模型建模集和驗證集R2分別為0.781和0.741,RMSE分別為225.671、332.876 kg/hm2,MRE分別為5.433%和7.731%;基于BPNN構建的大豆地上部生物量估算模型建模集和驗證集的R2分別為0.752和0.709,RMSE分別為238.264 kg/hm2和338.495 kg/hm2,MRE分別為5.529%和7.990%。

在R4生育期,基于SVM構建的大豆地上部生物量估算模型建模集和驗證集的R2分別為0.704和0.679,RMSE分別為268.616、341.065 kg/hm2,MRE分別為4.224%和7.246%; 基于RF構建的大豆地上部生物量估算模型建模集和驗證集的R2分別為0.847和0.818,RMSE分別為188.617、287.539 kg/hm2,MRE分別為3.707%和7.128%; 基于BPNN構建的大豆地上部生物量估算模型建模集和驗證集的R2分別為0.805和0.792,RMSE分別為271.815、296.111 kg/hm2,MRE分別為4.937%和7.399%。

在R6生育期,基于SVM構建的大豆地上部生物量估算模型建模集和驗證集的R2分別為0.687和0.636,RMSE分別為281.901、372.816 kg/hm2,MRE分別為6.191%和8.432%;基于RF構建的大豆地上部生物量估算模型建模集和驗證集R2分別為0.802和0.776,RMSE分別為214.021、313.524 kg/hm2,MRE分別為4.808%和8.983%;基于BPNN構建的大豆地上部生物量估算模型建模集和驗證集的R2分別為0.773和0.741,RMSE分別為227.772、315.486 kg/hm2,MRE分別為5.617%和8.916%。

3 討論

目前無人機多光譜遙感技術已經較為成熟,且在作物生理生態(tài)及生化參數(shù)反演方面的應用也越來越廣泛。當前使用無人機多光譜遙感技術對作物的生長生理指標的估算研究主要是通過計算植被指數(shù)并以此構建統(tǒng)計模型或機器學習模型,不同研究地域、植被種類和估算方法得到的研究結果基本相同[40-43],但精度存在一定的差異。

本文通過對植被指數(shù)與目標參數(shù)的相關性分析,發(fā)現(xiàn)植被指數(shù)GCI與LAI、地上部生物量和產量的相關性均較好,這可能是因為GCI對植物中葉綠素含量變化的響應比較敏感[12]。葉綠素含量可以反映植被的生理狀態(tài),葉片的生長及生物量的積累需要葉綠素的驅動進而進行光合作用[44],這直接導致葉綠素與LAI、地上部生物量和產量的相關性均較高,所以通過GCI的變化能很好地實時監(jiān)測大豆的長勢。本研究還發(fā)現(xiàn),隨著生育期的推進,大部分植被指數(shù)與產量的相關性呈先增大后減小的趨勢,同時,在輸入變量和機器學習模型相同的情況下,R4的產量預測模型精度更高,這可能是因為R4生長時期是作物生長最旺盛的時期,其LAI、葉綠素含量和光合速率都達到峰值,而在V4,由于大豆還處于營養(yǎng)生長階段,作物體內葉綠素含量水平較低,光譜鏡頭對該信息的捕捉效果受限;同理,在R6隨著葉片衰老和脫落同樣導致作物體內葉綠素含量變低,這與陶惠林等[45]的研究結果相似。

不同的建模方法顯著影響監(jiān)測模型的預測精度[46]。本文通過篩選15種植被指數(shù),同時使用3種機器學習方法構建大豆LAI、地上部生物量和產量估測模型。經對比分析發(fā)現(xiàn),基于RF模型的估測精度均高于SVM模型和BPNN模型,說明RF在估算大豆LAI、地上部生物量和產量方面更有優(yōu)勢,這與陳曉凱等[47]的研究結果相類似,主要原因是RF模型對自變量共線性問題敏感度較差,具有強抗干擾性和抗過擬合能力,對噪聲異常值容忍性高,且無需檢驗變量的正態(tài)性和獨立性,可以有效地提高模型的精度[48];而SVM模型的核心問題是確定核函數(shù)與懲罰因子,由于其抗干擾性較差,核函數(shù)與懲罰因子等參數(shù)選擇受限,導致預測效果受限[49];BPNN模型的估計精度較RF模型低,可能是因為本研究中樣本數(shù)較少導致泛化能力低造成的[50]。因此,RF模型可以作為大豆LAI、地上生物量和產量建模的首選方法。

此外,本研究結果還發(fā)現(xiàn)在進行大豆LAI與地上部生物量估算時,在運用同一種建模方法時,前者的估算模型精度均高于后者,這一結果與文獻[51-52]的研究結果一致。造成這一現(xiàn)象的原因可能是LAI與生物量相比更為直觀,是作物重要的形態(tài)參數(shù)。目前,基于經驗植被指數(shù)建立作物生長估算模型和產量估算研究中仍有一些問題有待解決,如不同區(qū)域、相同區(qū)域不同作物、甚至相同區(qū)域相同作物但不同時期的不同植被指數(shù)的相關系數(shù)變化及模型反演精度可能有所不同。除此之外,本研究后期也可以嘗試使用波段更多的高光譜儀進行光譜信息的采集,以避免多光譜受部分光譜信息缺失的不利影響。

4 結論

(1)光譜指數(shù)GCI與大豆LAI、地上部生物量和產量的相關性均較好,通過GCI的變化能實時監(jiān)測大豆的長勢。

(2)基于RF模型構建的大豆LAI和地上部生物量預測模型的精度顯著高于SVM與BPNN模型。

(3)在R4生育期,基于RF模型構建的大豆產量預測模型的精度最高,故在R4生育期使用RF模型能最準確地估算產量。