基于無人機多時相植被指數的冬小麥產量估測

程 千 徐洪剛 曹引波 段福義 陳 震

(中國農業科學院農田灌溉研究所，新鄉 453002)

0 引言

收獲前，及時準確地估測農作物產量有助于相關糧食政策的調整及制定。遙感技術是監測作物生長狀況的有效手段，遙感估產則是將遙感信息作為輸入變量，通過建立遙感估產模型來估測作物產量。隨著科技的發展，無人機能夠以低成本和更實用的方式獲取高時空分辨率的遙感數據，并將其應用于作物生長監測[1-6]和產量估測中[7-9]，彌補了衛星遙感運行周期長、空間分辨率低和受氣象條件影響大等不足。

多光譜傳感器包含的近紅外波段能夠較好地捕捉作物整個生育期內葉片結構變化產生的光譜反射特征，其與其他波段構建的植被指數包含更多的作物生長信息，有助于產量估測。通過建立植被指數與產量間的統計關系進行遙感估產的方法，簡單易行，目前已有很多學者利用無人機平臺搭載多光譜傳感器對不同作物進行產量估測[6,10-13]。

在冬小麥生長階段，冠層反射率在不斷發生變化，冬小麥生殖生長階段的光譜指數與產量具有更高的相關性，可更好地用來估測產量，如開花期[14]或灌漿初期[15]的光譜植被指數等。隨著生長過程的延續，產量與單一時期光譜指數之間的相關性達到非常顯著的水平，但是無法在單個時期內確定產量形成的過程。相反，由多個生育期生長信息構成的數據集能更好地反映作物生長變化過程，可以提供更多有用的信息，有利于提高算法的估測精度。文獻[16-17]利用多個生育期的植被指數的累積值進行作物產量估測，提高了估測精度。本文嘗試將新生育期的植被指數作為樣本數據新的屬性值添加到上一個生育期的數據集中，形成新的多時相植被指數數據集，從而對冬小麥產量進行估測。

本文采用無人機搭載多光譜傳感器，對冬小麥多個生育期進行遙感監測，提取多種光譜植被指數，并利用多時相的植被指數構建數據集，訓練機器學習算法，進行冬小麥產量的估測分析。

1 材料與方法

1.1 試驗區概況

試驗區位于河南省新鄉市七里營鎮中國農業科學院新鄉綜合試驗基地(113°45′38″N，35°8′10″E)，該試驗區地形為黃河北部沖積扇平原，海拔73 m，土壤類型為粘壤土。氣候類型屬于溫帶季風氣候，四季分明，雨熱同期，多年年平均氣溫14.5℃，年平均降水量560 mm左右，主要集中在6—9月，平均無霜期199 d，試驗研究區位置見圖1。本次研究共有3種灌溉水量處理(Irrigation treatment，IT)：240 mm(IT1)、190 mm(IT2)和145 mm(IT3)，每種灌水處理包含60個小區；試驗選取了30個小麥品種，每個小區1個品種，隨機分布種植，則每個灌水處理對應30個品種，重復兩次，試驗設計如圖1所示。

1.2 無人機平臺與遙感數據獲取

采用大疆經緯M210型四旋翼無人機平臺，無人機凈質量4.69 kg，最大承載質量1.45 kg,續航時間20 min左右。無人機搭載MicaSense RedEdge-MX型多光譜相機拍攝多光譜影像。多光譜相機質量232 g,圖像尺寸1 280像素×960像素,具有藍(Blue，B，(475±20)nm)、綠(Green，G，(560±20)nm)、紅(Red，R，(668±10)nm)、紅邊(Red-edge，RE，(717±10)nm)和近紅外(Near infrared，NIR，(840±40)nm)共5個通道。

冬小麥抽穗之后，進入生殖生長階段，各器官營養逐步轉化為籽粒的形成，小麥葉片開始衰落，此時雖然冠層光譜反射率不斷降低，但與產量的相關性逐步升高[12]。因此，從冬小麥開花期到灌漿期，選取多次采集的多光譜影像，進行產量預測。具體采集日期為：開花期(S1，2020年4月14日)、灌漿初期(S2，2020年4月23日)、灌漿中期(S3，2020年4月30日)和灌漿末期(S4，2020年5月10日)。

影像采集時使用大疆GSPro地面站對試驗區覆蓋范圍規劃航線自主飛行。為保證影像質量，減少大氣影響，選擇天氣晴好、無云低風速且11:00—14:00時間段飛行。飛行高度30 m，影像重疊率為航向85%和旁向80%,地面空間分辨率為2 cm。

1.3 數據處理

1.3.1影像預處理

采用Pix4Dmapper軟件對無人機的多光譜影像進行拼接處理，并導入每次飛行前后拍攝的校準板的影像進行輻射校正，生成該試驗區地物反射率的正射影像圖。將正射影像圖導入開源軟件QGIS，進行幾何校正、試驗小區矢量文件制作、多光譜植被指數計算、各小區多光譜植被指數統計等，得到各小區多光譜植被指數均值作為產量估測的輸入變量。

1.3.2多光譜植被指數選取

目前已有多種多光譜植被指數用于作物產量估測研究，如表1所示，大致分為以下幾類[18]：①基于紅波段和近紅外波段運算的多光譜指數：差異植被指數(Difference vegetation index，DVI)、增強植被指數(Enhanced vegetation index，EVI)、比率植被指數(Ratio vegetation index，RVI)、改善簡單比率指數(Modified simple ratio，MSR)、歸一化植被指數(Normalized difference vegetation index，NDVI)、植被色素比率(Plant pigment ratio，PPR)、結構不敏感色素指數(Structure insensitive pigment index，SIPI)和改進的非線性植被指數(Modified nonlinear vegetation index，MNVI)。②優化土壤背景的紅波段和近紅外波段運算的多光譜指數：土壤調節植被指數(Soil-adjusted vegetation index，SAVI)和優化土壤調節植被指數(Optimization of soil-adjusted vegetation index，OSAVI)。③利用紅邊波段特征的多光譜指數：紅邊葉綠素指數(Red-edge chlorophyll index，CIRE)、紅邊歸一化植被指數(Red-edge NDVI，NDVIRE)、改良葉綠素吸收率指數(Modified chlorophyll absorption ratio index，MCARI)和轉化葉綠素吸收反射指數(Transformed chlorophyll absorption reflectance index，TCARI)。④基于從綠波段到近紅外波段光譜反射率曲線形狀的指數：三角植被指數(Triangular vegetation index，TVI)和改善三角植被指數2(Modified triangular vegetation index 2，MTVI2)等。

表1 植被指數計算公式Tab.1 Vegetation index calculation formula

1.4 數據分析方法

水分虧缺條件下冬小麥生長特征出現差異，遙感影像的像元灰度變化能夠反映生長特征的變化。通過多光譜植被指數評價水分虧缺條件下，冬小麥生育期的生長差異，并采用皮爾遜相關系數和決定系數評價不同生育期的多光譜植被指數與冬小麥產量的線性相關程度,選擇合適的生育期和植被指數進行產量估測。

分別采用偏最小二乘回歸(Partial least squares regression, PLSR)、支持向量機回歸(Support vector regression, SVR)、隨機森林回歸(Random forest regression, RFR)等統計學習算法進行冬小麥產量估測；并采用決定系數R2、均方根誤差(Root mean squared error, RMSE)和平均絕對誤差(Mean absolute error, MAE)等指標評價不同統計學習算法的估測精度。以上數據分析均使用Python語言中的Scikit-learn[34]機器學習包完成。

2 結果與分析

2.1 植被指數與產量相關性分析

不同生育期多光譜植被指數與冬小麥產量的相關系數如圖2所示，可以看出，進入開花期后，光譜指數與產量出現相關性，隨著冬小麥生長，相關性逐漸增強，灌漿末期多個植被指數與產量相關系數達到0.7。灌漿期后大部分植被指數與產量呈正相關關系，且植被指數之間相關性較強，而TCARI和SIPI與產量及其他植被指數呈較強的負相關關系。

不同生育期多光譜植被指數與冬小麥產量的線性回歸決定系數如表2所示，可以看出，隨著冬小麥生殖階段不斷發展，決定系數不斷增大，在灌漿末期達到最大。除PPR外，其他植被指數均能在一定程度上反映產量的變化特征。可見本文選取的多光譜植被指數與冬小麥產量有較強的相關性，能夠用來估測冬小麥產量。

表2 不同生育期植被指數與產量決定系數Tab.2 Coefficient of determination between vegetation index and yield in different growth stages

2.2 多時相植被指數估測產量

選取不同生育期多光譜植被指數，采用PLSR、SVR和RFR等統計學習算法估測冬小麥產量，結果如圖3和表3所示。單時相時(如S1或S2等)，將文中選取的16種植被指數分別作為自變量的1個維度，形成自變量向量，即每個樣本數據有16個屬性值。當構建多時相植被指數數據集時(如S1+S2等)，將新生育期(S2)的16個植被指數合并到前1個生育期(S1)的自變量向量中，形成新的多時相植被指數數據集，每增加1個時相，自變量向量增加16個維度，即每個樣本數據增加16個屬性值。

在訓練統計學習算法時，為綜合提高算法精度和泛化能力，常采用交叉驗證的方法，由于試驗有3種灌水處理，為減少灌水處理對算法的影響，保證數據分布的均勻性，本文采用10次5折交叉驗證，取均值作為每種算法的評測指標值，圖3和表3中各評測指標值針對測試集的數據，更能代表模型的估測能力。

由圖3和表3可知，當使用單個生育期的植被指數估測冬小麥產量時，由開花期到灌漿初期、灌漿中期和灌漿末期，不同算法的預測精度不斷提高，這與植被指數與產量相關性變化規律一致。

表3 不同回歸模型估測產量統計分析Tab.3 Validation statistics of different regression models for grain yield prediction

使用生育期(S1+S2)的多時相植被指數比單個時相(S2)的植被指數估測產量的精度高，PLSR模型估測精度R2提高約0.021；SVR模型R2提高約0.015；RFR模型R2提高約0.051。同理，生育期(S1+S2+S3)的多時相植被指數比單個時相(S1或S2或S3)的植被指數估測產量的精度高，說明在灌漿末期(S4)到來前，這種數據集構建方式能夠提高產量估測精度，對產量作出較好的提前估測。灌漿末期(S4)的植被指數與產量相關性更好，估測精度更高，此時構建的多時相植被指數，不同模型均有較高的估測精度。

對比3種統計學習算法，針對本文的數據集，當數據屬性值增多時，PLSR估測能力低于SVR和RFR。PLSR模型針對生育期(S1+S2+S4)數據測試集的估測精度最高，平均R2為0.459，RMSE為1 822.746 kg/hm2；SVR模型針對生育期S4數據測試集的估測精度最高，平均R2為0.540，RMSE為1 676.520 kg/hm2；RFR模型針對生育期(S3+S4)數據測試集的估測精度最高，平均R2為0.560，RMSE為1 633.896 kg/hm2。而且，RFR模型針對生育期S4后不同時相植被指數組合的數據測試集，均表現出較高的估測精度。

圖4為不同灌水處理時，不同品種冬小麥產量估測值與實測值對比，圖中估測值針對某一生育期的全部數據樣本，且估測模型分別選取估測精度最高的訓練模型，PLSR模型針對生育期(S1+S2+S4)，SVR模型針對生育期S4，RFR模型針對生育期(S3+S4)。可以看出，3種模型的估測值均反映不同灌溉處理條件下產量變化規律，并且RFR模型的估測值與實測值更為接近。

正常灌溉處理(IT1)時，冬小麥產量平均達到8 700 kg/hm2，水分虧缺處理(IT2和IT3)時，產量逐步降低，且水分虧缺嚴重時，產量下降較多，相比正常灌溉，冬小麥減產將近1/2。對比冬小麥產量估測值與實測值，IT1灌水處理時，估測值偏小；IT2和IT3灌水處理時，估測值與實測值較為接近；但每個灌水處理條件下估測值的分布范圍均比實測值的小。由于試驗區冬小麥品種較多，同一灌水處理條件下，特別是正常灌水處理時，產量相差較大，使得數據分布較為分散。

3 討論

本文將新生育期的植被指數作為新的屬性值，添加到上一個生育期的數據集中，形成新的多時相植被指數數據集，提高了冬小麥產量估測精度，這與文獻[11,13,17]的結論一致。而且文獻[11]研究中同時使用了株高參數，并發現灌漿期株高比植被指數對冬小麥產量有更高的估測精度，說明不同類型傳感器獲取的冠層參數，經過數據融合能進一步提高作物產量估測精度。可見，無論是多時相數據或多源數據，均提高了樣本數據屬性值的維度，有助于產量估測模型的訓練。

采用統計學習模型進行回歸預測時，樣本數據的選擇對模型的訓練效果有較大影響，特別是樣本數據量較小時，訓練數據集的劃分更為重要。交叉驗證法適用于樣本數量不大的情況，而且對模型偏差和方差能夠兼顧。本文采用P次K折交叉驗證法，減小訓練數據隨機選取對模型預測精度帶來的影響。P次為重復K折交叉驗證的次數，通過P次重復減小隨機誤差，本文P選取10。K折即把數據集劃分為K等份，其中K-1份作為訓練集，剩余1份作為測試集，由于樣本總量為180個，單生育期時數據屬性為16個，多時相時數據屬性最多有64個，為保證不同生育期訓練模型的偏差和方差均較低，K值不宜過大或過小，本文K選取5作為數據集劃分標準。

冬小麥產量構成因素較多，與多光譜植被指數之間存在較為復雜的非線性關系，PLSR模型雖然包含了主成分分析過程，可以有效解決自變量之間的多個共線性問題，但當自變量與因變量之間存在復雜的關系時，其仍有較大的局限性。SVR模型和RFR模型對復雜的非線性關系均有一定的回歸能力，并且RFR模型具有更好的抗噪能力，本文冬小麥試驗中有30個小麥品種，數據分布較為離散，針對灌漿末期的不同時相植被指數，RFR模型均表現出較高的估測精度。但是當樣本數據量過大，或者樣本數據特征過多時，傳統統計學習算法的預測能力將下降，需要訓練更優秀的模型，如深度神經網絡模型等。

4 結論

(1)隨著冬小麥的不斷生長，植被指數與產量的相關性不斷增強，灌漿末期多個植被指數與產量相關系數達到0.7，基于此生育期植被指數估測冬小麥產量的精度也較高。

(2)采用多時相植被指數能進一步提高冬小麥產量的估測精度，采用冬小麥開花期和灌漿初期的多時相植被指數進行估產比采用單個生育期的植被指數估測產量的精度高，PLSR模型估測決定系數R2提高約0.021，SVR模型R2提高約0.015，RFR模型R2提高約0.051。灌漿末期的多時相植被指數，3種模型均有較高的估測精度。PLSR模型估測精度最高時R2、RMSE分別為0.459、1 822.746 kg/hm2，SVR模型估測精度最高時R2、RMSE分別為0.540、1 676.520 kg/hm2，RFR模型估測精度最高時R2、RMSE分別為0.560、1 633.896 kg/hm2，相比之下，RFR模型具有更好的穩定性。

(3)長時間的水分虧缺導致冬小麥產量下降，與正常灌溉水量相比，長時間的水分虧缺將造成冬小麥減產約1/2，此時需及時制定灌溉計劃。