葉 強, 楊鳳海*, 劉煥軍,2, 孟祥添, 官海翔, 崔 楊
(1.東北農業(yè)大學公共管理與法學院,哈爾濱 150030;2.中國科學院東北地理與農業(yè)生態(tài)研究所,長春 130012)
玉米是中國三大主要糧食作物之一,種植面積及年產量達到中國糧食總產量的33.6%和36.1%。玉米作物長勢動態(tài)監(jiān)測及產量的準確預測,為規(guī)?;r業(yè)經營管理提供科技支撐,是精準農業(yè)發(fā)展的迫切需求,對國家糧食宏觀價格調控,三農振興以及糧食外貿等具有重要意義。
玉米產量的預測方法多采用傳統(tǒng)估產法、農學估產法及氣象統(tǒng)計法等。其中,農學預報法主要根據產量器官參數如有效穗數、百粒重、平均穗粒數等估算產量,受作物種類、品種、耕作制度等限制[1];氣象統(tǒng)計法以天氣要素預測產量,發(fā)展較為成熟,但存在單一氣象因子估產精度不足且各要素空間插值方法不同,難以劃定氣象站數據范圍的問題[2];傳統(tǒng)估產法中地面采樣工作采集樣本數據規(guī)模較小、時間周期長、氣象條件限制,樣本實驗室理化性質檢測結果存在滯后性的問題,難以滿足及時、準確的作物估產和田間作物生態(tài)保護的管理需求[3]。遙感技術突破地面采樣的局限,以快速、無損、獲取成本低、多尺度、多時間序列的優(yōu)點成為精準農業(yè)研究熱點[4]。Chen等[5]采用作物物候和葉面積指數引入MCWLA-Wheat同化模型提升相對估產精度16%。朱婉雪等[6]結合少量產量樣本點與植被指數建立線性模型,決定系數R2最高達到0.7;Muoz等[7]應用產量和遙感參量建立指數模型和Richards非線性回歸模型,偏差信息準則(deviance information criterion,DIC)最低為62.8。但是,對田塊尺度作物遙感估產,測產多采用對角線、十字、“S”形,同級重復采樣以保證單產精度。測產過程耗時耗力,破壞采樣區(qū)間作物樣本,存在主觀因素影響,不同規(guī)模測產工作造成大量糧食損耗。
谷物聯(lián)合收割機搭載測產系統(tǒng)可以實現同步收獲、同步測產,從而節(jié)省人工,提高測產效率,降低作物糧食損耗量。該系統(tǒng)裝載有差分校正的全球定位系統(tǒng)(global positioning system,GPS)能夠收集有空間地理信息的產量數據,測產精度可達97%以上[8]。為提高遙感估產精度,以遙感數據作為主要信息源的統(tǒng)計模型不斷優(yōu)化,結合作物生長模型估產,或加入遙感輔助信息,如土地分類圖、作物氮含量、作物株高[9]、葉面積指數、雷達數據、傾斜測量構建3D點云等[10],以提高統(tǒng)計經驗模型的可靠性和可延展性。但作物生理參數以及土壤化學性質隨作物生長周期變化,多時期采集樣本費時耗力。高精度地形數據空間分辨率高,地形演變在一定時期內穩(wěn)定性強,影響作物水肥運移,間接作用于土壤理化性質和植被生長,造成作物產量空間異質性。地形因子作為輔助信息引入遙感估產模型具有合理性。大尺度區(qū)域研究中高程影響稻田空間格局和生產資源,如光照、土壤類型等[11];李宇宇等[12]研究了赤峰草場發(fā)現低濕區(qū)牧草鮮重分別高于坡地和平地67.94%,45.7%;土壤屬性與甘蔗產量在凹地較高,在平坦地形處空間變異性較弱[13];坡度約為2°區(qū)域大豆長勢最佳[14];自然降雨條件下坡耕地坡度為15°設立高反坡階玉米產量最高[15]。目前,田塊尺度結合地形因子提高變量測產數據遙感反演估產精度的研究尚未完善。
以黑龍江省北安市規(guī)?;N植的春玉米為研究對象,利用作物關鍵生長期SPOT 6影像,提取7種典型植被指數,結合變量測產數據建立經驗統(tǒng)計模型,確定適合春玉米遙感估產的最佳生長時期及最優(yōu)植被指數;基于無人機獲取的地形數據提取6種地形因子,應用多元逐步回歸評價引入地形因子的春玉米遙感估產模型;應用聚類分析和空間統(tǒng)計分析探索產量空間分布格局。研究旨在篩選春玉米監(jiān)測最優(yōu)植被指數和最佳時期,建立以遙感信息結合地形因子的估產模型,確定玉米最佳生長區(qū)域,為農場精細生產管理、耕地高效利用、玉米單產精細尺度制圖提供科學依據。
研究區(qū)位于黑龍江省北安市,地處小興安嶺南麓邊緣,屬于典型黑土區(qū)壟作耕地,年平均氣溫0.5 ℃,年平均降雨量為570 mm。變量測產數據采集區(qū)域位于該市趙光農場,農場經營管理條件一致,具有相同的機械化耕作方式與種植制度。主要種植作物為春玉米,春播晚熟,一年一熟。播種時間為每年4月底,成熟收獲時間為同年9月底10月初。實驗區(qū)田塊面積5 hm2,中心經緯度為48°2′42.45″ N、126°45′19.45″ E,平均海拔為307.41 m,平均壟長為553.09 m,最大高程差為9.37 m,采用玉米大壟雙行種植模式,壟寬為1.1 m。
編程訂購2018年7月27日、8月19日、9月13日的SPOT 6多光譜衛(wèi)星影像,分別為玉米生命周期中吐絲期(S1)、灌漿期(S2)和成熟期(S3)。衛(wèi)星影像空間分辨率為6 m,包含藍、綠、紅和近紅4個波段。使用ENVI 5.1對影像進行輻射定標、大氣校正、正射校正等預處理。高精度地形數據獲取,采用DJI Phantom 4 Pro RTK無人機飛行平臺,搭載1 inch(1 inch=2.54 cm)2 000萬像素CMOS傳感器高清影像相機獲取照片,飛行高度110 m,旁向重疊率為80%,航向重疊率為80%,氣象條件晴朗無風,應用Pix4D Mapper軟件對含有經緯度和高程信息的照片數據進行三維重建,經過軟件自帶的濾波功能處理生產空間分辨率為0.074 m柵格數據,基于衛(wèi)星影像空間分辨率重采樣為6 m。變量測產數據獲取時間為2018年9月30日,利用搭載美國Raven Envizio Pro智能測產系統(tǒng)的聯(lián)合收割機測定含空間地理信息的產量數據,收割機割幅6 600 mm,采集量為6大壟玉米,采樣頻率為1 Hz,平均步進距離1.27 m。因為田間復雜情況聯(lián)合收割機存在軌跡彎曲、停滯不前、信號不穩(wěn)定等問題,造成產量數據點缺失及分布不均。選擇產量數據點相對完整,軌跡較直,行進速度相對穩(wěn)定的田塊變量測產數據。產量校正基于三倍標準差法將超出值域范圍的樣本數據視為異常點刪除?;谟跋駭祿嬎銝鸥裣裨a量均值,最終保留578個產量點。作物產量屬性特征如表1所示。
1.3.1 植被指數與地形因子
為確定最適合春玉米估產的植被指數,根據國內外作物產量預測的相關文獻和玉米自身光譜特性,選擇7種光譜植被指數,構建作物估產模型,即歸一化植被指數(normalized difference vegetation index,NDVI)、差值植被指數(difference vegetation index,DVI)、比值植被指數(ratio vegetation index,RVI)、無藍色波段增強型植被指數(enhanced vegetation index without a blue band,EVI2EVI2)、次生修正土壤調節(jié)植被指數(modified secondary soil adjusted vegetation index,MSAVI2)、綠度歸一化植被指數(green normalized difference vegetation index,GNDVI)、修正比值植被指數(modified simple ratio,MSR)[6]。地形因子中高程、坡度、坡向、曲率等影響作物水肥運移方向和速度,使用ArcGIS 10.3提取6種地形因子高程(elevation)、坡度(slope)、坡向(aspect)、曲率(curvature)、平面曲率(plane curvature)、剖面曲率(profile curvature)。
1.3.2 春玉米遙感估產模型
為快速、準確地估算玉米產量,本研究按照2:1的比例將數據集分為訓練集和驗證集,各含390個和188個樣本點,采用線性回歸模型,利用最小二乘法,構建基于植被指數與產量的經驗統(tǒng)計模型,如式(1)所示。應用多元逐步回歸分析,評價引入地形因子的遙感反演估產模型。精度評價指標選擇決定系數R2、均方根誤差(root mean squared error,RMSE),檢驗遙感估產模型的擬合效果,其計算公式分別為
Yield=aVIt+b
(1)
(2)
(3)
1.3.3 空間格局分析
地理空間相關性在自然界中廣泛存在,距離越近的事物之間空間相關性越強,地統(tǒng)計學分析基于此特性在地理空間研究中應用廣泛?;谧兞繙y產數據,使用ArcGIS 10.3中Geostatistics模塊分析其趨勢面特征及相關地形因子趨勢面特征,應用K-mean聚類和空間統(tǒng)計分析探索春玉米產量的空間分布格局及地形因子對實際產量的影響。IBM SPSS Statistics 22中K-mean聚類分析是空間類型劃分的主要工具,該算法快速、簡單具有可伸縮性,能夠高效處理大數據集,實現多維度的“自然”聚類,是應用范圍最廣泛的聚類算法。其算法過程是在數據集中隨機篩選k個樣本,每個數據作為簇中心,計算其于其他數據中心距離,根據距離最近原則劃分簇類。所有數據劃分后,計算每個簇均值,并以均值重新聚類直至滿足收斂函數的要求,可表示為
(4)
為進一步分析田塊尺度中基于地形因子的實測產量空間分布規(guī)律,本研究采用自然斷點法,將不同地形因子分為高、中、低3類,使用分區(qū)中產量的均值和標準差評價產量組間變化。自然斷點法具有組內間距最小,組間間距最大的特點,應用其對地理空間柵格數據分類,可以更加直觀反映產量的空間分布格局。
探究不同生育期春玉米衛(wèi)星遙感估產精度,分別將春玉米吐絲期、灌漿期和成熟期的植被指數與實測單產數據建立經驗統(tǒng)計模型。如表2所示,春玉米3個關鍵生長期的R2和RMSE中,灌漿期的R2普遍在0.6以上,最高達到0.73。估產模型精度由高到低的順序是灌漿期、吐絲期、成熟期。灌漿期模型的RMSE比吐絲期和成熟期低,即灌漿期的產量估算與實際產量的偏差最小。主要是由于作物出苗、拔節(jié)至吐絲期前,營養(yǎng)物質主要集中在根莖葉等生長部位。吐絲期至結穗灌漿期作物光合作用形成的淀粉、蛋白質等大量營養(yǎng)物質轉移至產量器官,是決定玉米產量的重要時期。而成熟期田間玉米葉片水分、葉綠素下降,存在倒伏現象,估產精度較差。綜上,春玉米灌漿期決定產量器官營養(yǎng)物質積累,是遙感反演估產的最佳時期。
表2 春玉米3個生長期植被指數估產評價模型
不同時期、不同植被指數的產量預測精度存在顯著差異。在吐絲期植被指數EVI2的R2最高,可達0.24。至灌漿期,RVI的預測精度達到0.73,其余植被指數依次為MSR、EVI2、NDVI、MSAVI2、GNDVI,R2均大于0.61,對應RMSE也比DVI的估產誤差小。成熟期DVI的R2達到0.25,比吐絲期大,其余植被指數R2普遍較低,約為0.19。綜上,估產模型R2最高可達0.73,RMSE為388.87 kg/hm2,最優(yōu)植被指數為灌漿期的RVI,遙感預測模型為
Yield=1 482.261 RVI-4 771.111
(5)
分析不同地形因子與春玉米產量的敏感性作用。如表3所示,農作物實測產量與高程、曲率、剖面曲率、坡度和坡向均顯著相關(P<0.05),但不同地形因子對實測產量敏感性存在顯著差異,其中高程和坡度與產量的相關系數達到0.25以上,P值均小于0.01,坡度與產量的相關系數比高程大,說明高程和坡度是決定水肥運移,影響作物產量的重要因素。
表3 不同植被指數和地形因子與農作物實測產量的相關系數
如表4所示,春玉米灌漿期植被指數結合地形因子的產量估算模型及對應的R2和RMSE中,RVI的R2最高達到0.79,其余植被指數依次為MSR、EVI2、NDVI、MSAVI2、GNDVI,R2均大于0.7,相應RMSE比DVI的誤差小。對比灌漿期植被指數與實測產量建立的經驗統(tǒng)計模型,R2提升最大的為DVI估產模型,達到10.17%,其余植被指數相對提升精度均大于5.5%以上,說明地形因子作為遙感輔助信息能夠減小春玉米預測產量與實際產量之間的偏差。
表4 引入地形因子的玉米預測模型結果
如圖1所示,產量與地形因子趨勢面分析圖中x軸表示地理空間中東西方向,y軸表示地理空間中南北方向,z軸表示采集樣本數據屬性值在空間分布中的大小。圖1(a)表示產量在東西方向和南北方向均呈現高-低-高的顯著變化趨勢,與圖1(b)坡度因子趨勢面的規(guī)律相反,而圖1(c)高程趨勢面在南北方向和東西方向上均逐漸降低。說明地形對農作物實測產量存在影響。使用聚類分析工具,設置3個聚類種類,加入實測產量、高程和坡度因子。K-mean聚類分析結果顯示隨高程增大坡度減小,春玉米產量逐漸增高,產量最佳的聚類中心為高程307.93 m與坡度1.39°。
圖1 產量與地形趨勢面分布
為進一步分析產量空間分布規(guī)律,基于自然斷點法將高程與坡度因子分為3類生成9種不同地形因子組合。圖2為地形因子對產量影響的空間分布結果。結果顯示在高海拔-高坡度區(qū)域春玉米實測產量均值最低,為7 502.64 kg/hm2,標準差為802.62 kg/hm2,該區(qū)域海拔和坡度大,是水肥流失最嚴重的區(qū)域,難以滿足玉米生長需求。其他地形組合除了低海拔-低坡度區(qū)域沒有對應產量數據外,春玉米產量平均值均大于8 000 kg/hm2,在中海拔-低坡度區(qū)域春玉米實測產量均值最高,達到9 157.63 kg/hm2,標準差為357.48 kg/hm2(表5)。原因是崗地易出現缺肥而洼地富營養(yǎng)化,中海拔且坡度較小的地方水肥運移速度較慢,易于與田間施肥方案確定的肥量保持一致,滿足作物生長所需的營養(yǎng)。同時在中海拔-低坡度、中海拔-中坡度與高海拔-低坡度區(qū)域春玉米實測產量均值大于田塊產量平均水平,面積占比達到46.88%。高坡度區(qū)域產量均小于田塊平均水平,面積占比達到23.22%,說明不同坡度通過水肥運移過程造成的作物長勢差異顯著。
圖2 實測產量的空間分布
表5 不同地形因子對產量空間分布的影響
田塊尺度的春玉米產量估算,對于農場規(guī)?;洜I管理,糧食產量供給側調整,儲備糧更新,大型合作社生產等均具有實際需求。在制定科學管理決策時,當年農作物實際產量及其空間分布是重要參考依據,因此需要不斷提高產量的預測精度與分析其空間變異性。本研究中使用變量測產設備采集的春玉米單產數據,盡管按照測產標準設置了智能測產系統(tǒng)的參數,但在實際應用中仍會存在產量數據缺失,因此后期需對產量數據做精細校正。遙感影像數據空間分辨率為6 m,基于地理空間相關性原則,相同像元內作物長勢較為均勻,去除異常采樣點數據后,按照影像像元大小計算實測產量算數平均值的方法具有一定的合理性。
地物空間異質性隨農作物生命周期而改變,基于不同時期遙感監(jiān)測的估產模型精度存在顯著差異。玉米從播種至吐絲期主要生長部位為根莖葉,反映植株的干生物量積累,估產精度較低。灌漿期氮素吸收率高,至乳熟期是玉米一生中光合作用最強的時期。該時期由光合作用產生的蛋白質、脂肪、淀粉等有機質轉移至玉米籽粒中,故此時估產精度高。成熟期冠層水分與葉綠素含量下降,葉片發(fā)黃導致冠層光譜植被指數的估產精度降低。本研究中灌漿期估產精度R2大于0.6,最高達到0.73,是春玉米估產的最佳監(jiān)測時期。RVI對綠色作物的植被覆蓋度反應敏感,從出苗至灌漿期葉面積指數增加,估產精度高。成熟期逐漸下降,且易出現玉米倒伏現象,此時估產精度最差。當前研究衛(wèi)星影像空間分辨率較大難以識別倒伏陰影,在今后的研究中可以考慮結合無人機影像濾波去除作物陰影對估產光譜指數的影響。
作物估產引入不同輔助因子優(yōu)化模型能夠提高精度。地形獲取途徑多樣,在短時間內獲得的多尺度空間分辨率影像提升估產精度。高程將田塊尺度內耕地劃分為不同坡位,決定水肥運移的方向,影響農作物生長時期的太陽輻射量和土壤溫度,陽坡面日照時間長,太陽輻射量高,土壤溫度相對較高,而陰坡面相反;洼地匯聚田塊降水,土壤水分含量高,崗地相對較低。坡度決定農作物生長過程中水肥運移速度。水土流失現象一般出現在坡度較大的區(qū)域,基于田間降水徑流方向呈條帶狀分布,土壤養(yǎng)分隨耕地表層土流失造成玉米產量下降,而低洼處易富營養(yǎng)化,造成土壤板結。坡面平緩的區(qū)域保持適合作物生長的土壤養(yǎng)分和水分,區(qū)域產量高。本研究中玉米產量與高程呈正相關,與坡度為負相關。今后研究中可以結合實測土壤理化性質數據構建作物生長模型探究地形對農作物生長要素的影響。
地理空間統(tǒng)計能夠反映春玉米產量數據在三維空間不同方向的變化趨勢,難以定量分析田塊尺度中春玉米產量的空間格局變化。K-mean聚類分析確定實際產量數據與地形因子的空間類型變化,通過地形因子重分類組合分析田塊尺度中不同區(qū)域春玉米產量的均值與標準差,進一步明確田塊尺度中精確的空間格局變化。實際生產管理中,田間地形主要通過水肥運移與太陽輻射量和地溫等因素間接影響農作物生長,在農作物整個生命周期中,它對作物生長的作用具有連續(xù)性、長期性,時間跨度越大,地形對作物產量的影響更加顯著。高海拔-高坡度區(qū)域春玉米產量均值最低,為7 502.64 kg/hm2,標準差為802.62 kg/hm2。在中海拔-低坡度區(qū)域產量最高,達到9 157.63 kg/hm2,標準差為357.48 kg/hm2。其次為中海拔-中坡度區(qū)域與高海拔-低坡度區(qū)域,產量均值大于田塊尺度中產量均值,面積占比總計達到46.88%。說明坡度平緩且海拔較高處春玉米產量相對較高,與K-mean聚類分析結果基本一致。今后研究中可以進一步探究多年間相同地形條件下產量時空分異規(guī)律。
結合田塊尺度地形特征、作物不同生長期SPOT 6影像與變量測產數據,使用經驗統(tǒng)計模型確定春玉米估產最佳生長期與最優(yōu)植被指數,評價引入地形輔助信息的估產模型,分析田塊尺度中產量空間格局分布。得出如下結論。
(1)春玉米遙感估產的最佳時期是灌漿期,其次是吐絲期和成熟期,RVI是最優(yōu)植被指數,R2達到0.73,RMSE為388.87 kg/hm2。
(2)最佳估產模型引入地形輔助信息后精度提升5.6%,R2達到0.79,RMSE為346.03 kg/hm2。
(3)空間格局分析顯示高海拔且高坡度區(qū)域產量最低,中海拔與低坡度區(qū)域產量最高。
(4)下一步將研究無人機影像濾除植被陰影估產與分析作物時空分異規(guī)律。本文研究建立了一種高精度的春玉米估產方法,確定玉米最佳生長區(qū)域的地形特征,為規(guī)?;r業(yè)精細管理、耕地整理與作物種植結構調整提供了科學依據。