基于細粒度校正的育種小區小麥株高無人機測量方法

吳婷婷 劉昕哲 聶睿琪 劉 佳 武 璐 李 濤

(1.西北農林科技大學機械與電子工程學院, 陜西楊凌 712100; 2.西北農林科技大學農學院, 陜西楊凌 712100)

0 引言

株高作為小麥的重要表型性狀[1],是生物量估算模型的重要輸入參數[2-4],還可用于小麥產量預測[5]、氮素營養評價[6]、倒伏識別[7]等,對選育高產、穩產且具有抗逆性的小麥品種具有重要意義[8-9]。現有田間小麥株高主要依賴人工標尺測量,效率低下,難以適應高通量表型組學的發展需要[10]。為解決此問題,可以通過效率高、無主觀誤差的間接法測量小麥田間株高性狀。

現有的間接測量方式有利用雷達點云、深度相機、無人機(Unmanned aerial vehicle,UAV)遙感等手段獲取田間育種小麥株高[11-13]。苗艷龍等[14]提出一種基于地基雷達生成三維點云的株高測量方法,與人工測量結果對比,均方根誤差為21.5 mm。孫翔[12]提出一種基于標尺圖像識別的株高測量方法,相對誤差小于3.55%。此外,無人機遙感因成本低、靈活性高等優勢,在大田作物表型中得到廣泛研究。樊意廣等[15]基于UAV生成試驗區域的數字正射影像(Digital orthophoto map,DOM)和數字表面模型(Digital surface model,DSM)提取冠層光譜特征和株高,實現了植株氮含量的估算。文獻[16-17]對基于UAV的多種空間輔助信息組合進行作物株高估算,結果表明最佳株高估算結果中決定系數R2為0.736,RMSE為5.3 cm。OEHME等[18]基于無人機DSM信息對不同生育階段玉米株高進行估測,測量誤差范圍為12.07%～19.62%;JAMIL等[9]發現基于UAV的小麥株高估測和地面實測數據相關系數僅為0.49,RMSE高達22.04 cm;HAN等[19]提出一種利用多層地面控制點改進無人機作物株高估算的方法,RMSE降低約20%。

低成本無人機測量作物株高普遍精度不夠理想,特別是對小麥這樣的低矮作物[20],在進入特定生育期后株高增長緩慢,高度落差不明顯,直接采用上述方法會存在較大誤差。作者通過全生育期高頻次飛行反演育種小區粒度(面積通常是幾平方米)的株高性狀,以平衡單次測量誤差,然而,多次反演出株高數據逆增長現象,顯然不符合作物生長機理,且發現這種逆生長差距隨時序增長被進一步放大。若隨意拋棄逆增長的株高反演結果,則生育期數據空缺易造成整體模型的不穩定,為此,在缺失地面控制點高程信息的條件下,考慮將時序變化中遇到的實際問題轉換到空間變化中來分析,本文提出一種在育種小區尺度下基于空間已有株高真實數據的近鄰校正法以及一種基于多光譜+RGB數據融合的半經驗光譜指數校正法,以實現對田間育種小區小麥株高的無人機遙感精確測量。

1 試驗材料與數據采集

1.1 試驗設計

試驗研究區位于陜西省楊凌農業高新技術產業示范區實驗基地(34°17′27.05″N,108°04′31.59″E),海拔435～563 m、年降水量635.1～663.9 mm,年平均氣溫12.9℃。選取565份(含5個對照品種)具有豐富遺傳多樣性的小麥種質作為本次研究的試驗材料,包括當前國內的主栽品種、高代品系、核心種質、農家種、骨干親本、不同時期培育的代表性品種以及國外引進的關鍵材料等。大田田間試驗設計采用隨機增廣試驗,共種植640個樣本小區,分16個區組,每區組選取西農511、百農207、偃展4110、濟麥22和周麥185等5個品種作為對照品種,每個小區播種6行,行長5 m,寬1.2 m,播種密度每平方米270株基本苗,播種期為2021年10月18日,具體研究區規劃如圖1所示。

圖1 研究區規劃圖及UAV平臺Fig.1 Research area map and UAV platform

1.2 無人機數據采集

無人機數據獲取平臺如圖1所示,RGB圖像獲取平臺采用小型多旋翼高精度Phantom 4 RTK型航測無人機(DJI大疆創新公司,深圳),搭載RGB相機對大田農情信息進行采集,相機分辨率2 000萬像素,搭配8.8 mm定焦鏡頭(等效焦距24 mm),多光譜圖像獲取平臺采用大疆Inspire2型無人機搭載MicaSense RedEdge-MX型多光譜相機對農田進行采集,多光譜相機每次拍攝的原始圖像尺寸為1 280像素×960像素,相機焦距為5.5 mm,視場角為47.2°,相機配有標準反射板,每次采集數據需要對相機進行校正,多光譜相機各通道中心波長與半波寬如表1所示。

表1 RedEdge-MX型多光譜相機的通道參數Tab.1 Channel of RedEdge-MX multispectral camera nm

2021年10月31日至2022年5月31日,間隔采集試驗田的RGB與多光譜影像共31次。為保證數據質量,選取晴天拍攝,拍攝時間為11:00—15:00之間,設定無人機自動巡航的飛行區域為150 m×100 m,飛行高度25 m,飛行速度2 m/s,航向重疊率85%,旁向重疊率80%,相機垂直拍攝,拍照間隔時間1 s。

1.3 地面數據采集

為了獲取作物高度,以評估株高估算精度和對應無人機圖像質量,在冬小麥的拔節到開花期選取試驗田中的107個小區進行了6次完備株高標尺測量,采樣日期為:2022年3月7、11、28日和2022年4月5、13、27日;采樣時間為10:00—11:00。測量土壤至小麥冠層的高度(不包括麥芒部分)。

2 基于無人機遙感的小麥株高精確測量

2.1 數據預處理

分別將多光譜圖像、RGB圖像導入到Pix4Dmapper軟件中,利用多光譜相機廠家提供的標定系數依次對小麥多光譜單幅影像進行標定,依據對應的位置與姿態數據生成稀疏點云,采用多視點立體匹配算法生成密集點云數據,最后導出為TIFF格式的DOM和DSM數據。生成的多光譜DOM影像空間分辨率約為每像素0.8 cm、可見光DOM影像空間分辨率約為每像素0.56 cm,其投影坐標系設置為WGS 1984坐標系下墨卡托投影分度帶的北半區49°帶。由于攜帶多光譜相機的無人機不具有RTK,故需選擇標準日期,利用地面控制點(GCP)將其余日期以標準日期為參考進行校準。為避免田間試驗破壞性采樣和雜草對試驗結果精確度造成影響,剔除田間破壞性采樣區域保留感興趣區域,將DOM及DSM(DTM)分割為單獨小區,同時利用閾值分割法去除多光譜圖像中的雜草信息。

2.2 高頻次UAV株高測量的誤差機理分析

采用UAV測量小麥株高時,通常將小麥的自然株高(小麥最頂端不包括麥芒部分到土壤的距離)作為度量株高變化的評價指標[16];依據研究側重點選取日長、積溫、生長度日等參數作為時序指標[21-23]。為更好觀測細粒度尺度下株高變化,本文選取日歷日期作為衡量小麥時序變化的指標。

表2為利用傳統方法[24](Traditional method,TM)估測的小麥株高和地面真值,在2022年3月11日、4月27日誤差較大,這是由于無人機飛行高度定位誤差、小麥冠層高程提取誤差而導致。如圖2所示,隨著時間的增長,小麥冠層發育逐漸封壟使得拍攝的RGB圖像特征點變少,同時由于進入多風多雨季節,采樣時由于風的影響使得冠層產生擾動,從而導致圖像在拼接生成DSM時誤差被進一步放大。此外由于無人機處于25 m的較高高度進行群體測量,麥穗間間隙導致小區冠層頂面非水平面,以及部分高桿小麥自然下垂但測量人員會將其扶正后測量,均會導致測量結果比真值偏小,且小麥封壟將進一步導致后期株高監測更為復雜,測量誤差也更大。

表2 基于傳統方法的無人機測量與地面測量株高結果對比Tab.2 Comparison of plant height based on traditional UAV measurement and ground measurement cm

圖2 冠層點云分布示意圖Fig.2 Distribution of point cloud in canopy

在空間分布上,試驗區小麥種質資源豐富,不同品系之間小麥固有發育速度有較大差異,最慢品種和最快品種的日株高差可達146 cm,這導致無人機每次采樣的株高數據之間存在較大梯度。以表2中2022年3月11日為例,雖然相對均方根誤差(RMSE100)為59.44%,但相關系數仍然處于較高的顯著性水平,進一步說明時序上出現株高逆增長現象,是設備與算法的固有誤差被相對放大而導致。

2.3 小麥株高精確測量方法

基于前述分析,本文在缺失地面控制點高程信息的條件下,基于空間分布上已有部分隨機采樣真值,提出一種在育種小區粒度下的近鄰校正法和一種基于多光譜+RGB數據融合的半經驗光譜指數校正法。試驗區整體平整,南北海拔差異為0.8 m、東西海拔差異為0.3 m,近鄰小區南北海拔落差在0.05 m以內、東西海拔落差在0.03 m以內,與各小區群體后期生長發育不齊帶來的株高落差相比,可以忽略不計,這為近鄰校正法提供現實基礎。此外,基于不同波段光譜反射率組合的植被指數(VIs)能夠有效反映植物生長狀況[25],以及無人機遙感能夠提供高分辨率的VIs測量,為光譜指數校正法提供了理論依據,故而可以在育種小區粒度下選取合適的VIs建立回歸模型預測小麥株高。圖3展示了數據預處理流程與兩種模型的建模思路。

圖3 總體技術流程圖Fig.3 Overall technical roadmap

2.3.1基于高分辨率RGB圖像的近鄰校正法

本方法在出苗之前采集試驗田裸地數據的DTM信息,并采用下述步驟計算小麥株高。

(1)小區群體高程獲取

不同采樣時間下,由于RTK信號會有輕微偏移,因此將不同日期各小區的遙感影像進行經緯度配準。依據地面控制點進行配準,配準方法采用一階多項式,其又稱為仿射變換,除了平移、縮放和旋轉之外,具有糾正地圖變形的功能,保證每個日期配準殘差小于1個像素。將獲取到的每兩期圖像做差,得到作物株高數據(Wheat plant height, WPH),計算式為

(1)

式中Sday——數字高程數據

T——裸土數字高程數據

(2)結合地埂高程校正

為解決由于設備與算法的固有誤差被相對放大而導致株高時序上出現逆增長問題,利用遙感影像中育種小區之間裸露地埂的土壤基底信息進行修正。由于小區分割時感興趣區域可能包含部分未種植的裸地信息,將每個生長階段DSM值的95%百分位數視為小麥頂部,將DSM值的5%視為裸土高程DEM,得到WPH的估計值,計算式為

(2)

其中

式中Eday——數字高程數據

α——算法調節因子

(3)依據真值滑動校正

在小麥生長發育后期,葉片會將小區間的地埂覆蓋,對無法通過式(1)、(2)修正的小區,采用與其鄰近的地面株高真值修正土壤基底信息。具體地,UAV測量的冠層海拔減去地面真值,得到各小區土壤海拔,并將這一海拔共享給周圍n個近鄰,n個近鄰的最終株高則為UAV測量的冠層海拔扣除上述共享海拔,將n+1個小區作為一個滑動窗口依次計算所有小區的株高(圖4)。

圖4 近鄰小區海拔共享Fig.4 Altitude sharing of adjacent cells

(3)

其中

2.3.2基于多光譜和RGB數據融合的光譜指數法

根據采集的五通道多光譜數據,計算植被指數及其與株高的相關系數,優選合適的植被指數與DSM數據參與建模得出株高回歸模型。

(1)VIs計算

多光譜遙感影像生成的植被指數地圖中像素點灰度變化能夠反映株高生長特征的差異,故對各育種小區發育過程中的植被指數求取均值。計算與作物高度相關程度高的15種植被指數[26-28]:歸一化植被指數(NDVI)、藍色歸一化差異植被指數(BNDVI)、藍色歸一化差異紅邊植被指數(BNDREI)、歸一化綠藍差異指數(NGBDI)、標準藍色指數(NormB)、綠葉指數(GLI)、綠/紅比率指數(GR)、過綠指數(ExG)、修正簡單比率(MSR)、紅邊葉綠素指數(CIRE)、改良葉綠素吸收率指數(MCARI)、轉化葉綠素吸收反射指數(TCARI)、三角植被指數(TVI)、根號歸一化植被指數(RNVI)、非線性植被指數(NLI)等。

V=f1(Blue,Green,Red,Red_edge,NIR)

(4)

式中,V為多種植被指數,Blue、Green、Red、Red_edge、NIR為五通道光譜反射率。

(2)植被指數優選

分別在時序上計算植被指數向量與株高向量的皮爾遜相關系數和空間上各育種小區的植被指數向量與株高向量的皮爾遜相關系數,形成三維相關系數矩陣,并在合適的顯著性水平下選出用于下一步建模的VIs[25-27]。

(5)

式中Vday——日期day植被指數向量

Gday——日期day株高向量

(3)VIs-WPH模型構建

將選擇出的VIs與DSM數據進行歸一化處理,分別采用偏最小二乘g1(V,S,b)、支持向量機g2(V,S,b)、隨機森林g3(V,S,b)3種回歸算法,計算校正后的株高,再進行模型的選擇。

Hday=min{RMSE(g(V,S,b))}

(6)

其中

g(V,S,b)=

{g1(V,S,b),g2(V,S,b),g3(V,S,b)}

式中Hday——日期day小麥株高

S——數字高程數據

b——調節系數

RMSE(g(V,S,b))——計算g(V,S,b)RMSE的函數

2.4 準確性評價指標

采用決定系數R2、均方根誤差(RMSE)、相對均方根誤差(RMSE100)、相關系數R評價模型的穩健性和準確性。R2越高擬合效果越好,RMSE越小,測量越精確。107個采樣點隨機分為80%的訓練集和20%的驗證集,采用五折交叉驗證的方法以獲得可靠結果。

3 結果與分析

3.1 植被指數與株高的相關性分析

根據相關性系數顯著性檢驗公式可知,自由度為107(即樣本量)時,0.02、0.01、0.05、0.001顯著水平的臨界值分別為0.224 7、0.248 0、0.269 5、0.313 7,即相關性系數絕對值大于0.248 0時達到極顯著水平。表3展示了15種植被指數在不同時期與小麥株高的皮爾遜相關系數,大部分植被指數與株高的相關系數從2022年3月28日才開始出現絕對值大于0.313 7,然后隨時間推移不斷增加,在2022年4月5日達到最高后開始不斷減小,值得注意的是過綠指數(Excess green,ExG)在整個時期均表現出較高的顯著性水平;在2022年3月28日至2022年4月13日這個時間段大部分植被指數與株高的相關性維持在較高的顯著性水平上,說明在挑旗期、抽穗期的植被指數能很好地反映株高變化,而在其它生長期與株高的相關性較低,難以精確反映小區粒度上的株高變化。

表3 不同時期植被指數與株高的皮爾遜相關系數Tab.3 Correlation coefficient between vegetation index and plant height in different periods

3.2 不同模型下的小麥株高時序變化

表4給出了TM、NNCM以及SICM的株高估算結果,可以看出傳統測量方法在2022年3月11日和2022年4月27日出現較大偏差,NNCM的測量精度在2022年4月5日達到最高,其R為0.94,RMSE為3.72 cm,RMSE100為5.62%,此時大部分育種小區完成拔節,株高處在40～100 cm范圍,動態變化快、空間差異大,有利于數據校正;測量精度在2022年4月13日開始降低,主要原因在于大部分品種的株高已經停止變化,營養用于葉片和生殖器官的發育,導致育種小區間地埂被逐步覆蓋,降低了算法的精度。SICM的測量精度在2022年3月28日達到最高,其R2為0.79,R為0.89,RMSE為3.25 cm,RMSE100為6.05%,此時小麥處于快速生長期,很高的R值說明植被指數能夠較好地反映小麥株高變化。從時序上看,2022年3月28至4月13日期間NNCM和SICM結果相比于其它時期具有更高的測量精度,說明這2種方法在小麥挑旗期與抽穗期取得了最好的應用效果。為了進一步觀測小區粒度株高的反演情況,得到單個小區的模型估測結果與真值的散點圖如圖5所示,不同列表示不同方法的株高反演情況,不同行表示不同時期的株高反演情況。可以看出,在2022年3月11日、2022年4月27日兩個時期無人機傳統測量方法估測的小麥株高明顯低于地面測量真值;相比于TM,NNCM與SICM觀測點擬合直線斜率k在6個時期均接近1,擬合直線截距更接近0,各觀測點小區更靠近坐標軸的45°線,說明NNCM與SICM能夠很好地測量小麥株高。SICM與NNCM相比在2022年4月27日結果較差,這是由于在小麥株高生長后期植被指數并不能夠很好地反映小麥株高的變化,故可以依據不同的實驗時期以及實驗條件,合理選擇SICM模型或NNCM模型。

表4 不同日期不同模型的評價指標計算結果Tab.4 Results of different models on different dates

圖5 模型估測結果精度及時序變化分析Fig.5 Accuracy of model estimation results and time series change analysis

3.3 不同生育期的小麥株高空間分布

圖6為小區細粒度下采用NNCM繪制的小麥株高時序變化空間分布,不同顏色代表不同株高,2022年3月7日多數小區處于橘黃色,將圖中不同小區的顏色與右側參考值對比,結合圖5可以看出在2022年3月7日小麥株高普遍在20～30 cm。此外,原圖在QGIS軟件中可以查看每個像素點的具體株高。

圖6 NNCM測量小麥株高時序變化的空間分布Fig.6 Spatial distribution of wheat plant height

在2022年3月28日試驗田北部株高在35～50 cm,南部株高在45～60 cm,表明南部種植的小麥品種相較于試驗田北部小麥品種有更快的株高生長速率;在2022年4月13日試驗田,南部種植的小麥株高比北部種植的小麥品種高20～40 cm,說明NNCM能夠精確反演小區粒度尺度下的株高及不同小區株高在空間中的變化,為后續全基因組關聯分析提供精確度高、通量大的小麥田間株高表型數據。

3.4 不同模型的性能對比分析

圖7給出了不同時期株高真值GM及使用兩種校正方法后的分布結果,對比發現,NNCM與SICM核密度分布曲線與地面真值都很接近,中位數、四分位數、最大值、最小值偏差均不超過0.5%,說明這兩種方法解決了細粒度株高精確測量的問題,有效校正了固有誤差帶來的數據逆增長現象。圖7小提琴圖的輪廓即核密度曲線,反映了株高性狀的分布,可以看出,本試驗區株高性狀服從正態分布,種質資源豐富,群體間株高性狀具有較大的差異性。此外還可以看出,在2022年4月13日株高增長較快的小麥已停止增長,較2022年4月5日采樣結果GM中位數增長19 cm、TM模型中位數增長10 cm、NNCM模型中位數增長17 cm、SICM模型中位數增長17 cm,相比于TM模型,NNCM和SICM能夠很好地同步跟隨GM觀測到的株高變化。

圖7 實測株高與模型校正后的株高分布Fig.7 Measured plant height and plant height distribution after model correction

4 討論

本文基于無人機遙感,通過TM、NNCM以及SICM提取小麥冠層株高,并與6個時期的實測株高對比,表明無人機遙感在小麥育種小區粒度下,基于人工測量真值的NNCM和基于多光譜+RGB數據融合的SICM能夠精確獲取小麥株高。NNCM獲取精度略高于SICM,在具有地面隨機采樣真值的條件下,獲取成本低于SICM;在不具有地面隨機采樣真值的條件下應當選用SICM, NNCM在實際農田的效率可能不如SICM,應當依據實際情況合適選擇測量方法。此外,無人機數據采集時,天氣變化、風的有無會對數據質量造成較大影響,所以要在天氣晴朗、無風或微風的條件下采集數據。下一步將在不同試驗田、不同采集設備條件下探究本文所提模型的泛化能力,本文所提思想還可以為小麥其他表型參數在細粒度尺度下精確測量提供借鑒。

5 結論

(1)驗證了基于細粒度校正的育種小區小麥株高無人機測量是可行的,本文提出的兩種株高精確獲取方法能夠在育種小區細粒度下精確獲取小麥株高,兩種方法的RMSE100可達6.00%左右,其中NNCM的RMSE100最低可達5.70%;SICM的RMSE100最低可達6.05%。

(2)NNCM和SICM對真實株高動態變化均具有良好的跟隨性,可以精確測量育種小麥小區粒度的全生育期株高,避免了UAV高頻次飛行傳統測量所出現的株高逆增長現象,在實際生產實踐與輔助育種中有較好的應用前景,可為細粒度下小麥株高無人機遙感精確獲取提供參考。