盧闖,胡海棠,覃苑,2,淮賀舉,李存軍*
(1.北京市農林科學院信息技術研究中心,北京 100097;2.遼寧科技大學,遼寧鞍山 114051)
我國東北平原是世界四大黑土區之一,發揮著“糧食穩壓器”的重要作用。近幾十年來,隨著機械化的發展和耕地資源的開墾利用,該地區農業生產逐步實現了集約化,但農田管理愈發單一,廣袤的農田往往采取相同的管理措施,而實際上土壤在空間上具有很大的變異性,黑土區“漫川漫崗”的地形特點進一步增大了土壤的空間變異程度,很大程度限制土地的產出率和投入農資的高效利用[1-2]。因此,提高農田管理的針對性是當地農業可持續發展的重要內容。
農田管理分區是根據各種因素的相似性和差異性,將地塊劃分成不同的子區域來指導變量管理,可顯著提高農田利用潛力。在分區數據源方面,基于土壤養分測試數據的分區方法準確度較高,但大尺度農田的密集采樣工作繁瑣,采集效率低下且分析費用高昂[3-4]。植被指數能夠表征作物生長信息,有學者利用較易獲取的高空遙感影像進行作物長勢監測和分區,張澤等[5]根據Landsat-5 TM歸一化植被指數劃分棉田管理區,分區結果與產量的符合度為75.47%;Breunig等[6]基于PlanetScope影像對黑麥等作物田進行分區,分區間地上生物量和作物產量差異顯著;Amanda等[7]利用2期大豆植被指數將地塊分為2個管理區,分區間土壤有效磷、黏土和粉土含量存在顯著差異。以上研究均證實了使用遙感影像進行分區的有效性和合理性,但高空遙感影像獲取依賴于衛星的運轉周期,獲取過程中存在分辨率粗糙和氣象干擾等影響因素,分區精度有待進一步提高。近年來無人機技術的迅速發展為農田遙感應用研究提供了新的數據源,低空遙感影像在克服天氣因素、提高分辨率、實時測量、控制成本等方面具有顯著優勢,有研究利用無人機影像對作物進行營養診斷,發現近地影像不僅能夠準確反映作物長勢情況,還能間接反映與作物長勢密不可分的土壤狀況[8-9],因此利用關鍵生育期的無人機影像進行精準管理分區具有一定的應用前景。
試驗地點位于黑龍江省北安市趙光農場,試驗地屬于寒溫帶季風氣候,年均氣溫-0.5℃,無霜期120 d,年均降雨量670 mm,年均日照時數2 700 h。選取典型地塊作為研究區(圖1),面積約為30.8 hm2,地塊中心經緯度為48.04°N、126.73°E,試驗地高程變化為296.04—312.27 m,高差16.23 m,2017年種植作物為春玉米,品種選用恒基利馬格蘭種業華美2號,播種時間為2017年5月8日,種植方式為大壟雙行栽培,壟距110 cm,壟上行間距44 cm,株距20 cm,播種密度90 000株·hm-2。底 肥 一 次 性 施 入 氮 肥60 kg N·hm-2,磷 肥80 kg P2O5·hm-2,鉀肥50 kg K2O·hm-2,底肥分別為尿素、磷酸二銨、氯化鉀。2017年6月24日追施尿素,壟上中間開溝施肥。2017年9月25日收獲測產,玉米生育期內其他管理措施與當地大田一致。
1.2.1 無人機影像采集 2017年7月22日,利用無人機獲取多光譜影像,此時正值春玉米吐絲期,葉面積指數達到生育期內最大值,玉米長勢與植被指數相關性達到最大,并很大程度上決定了最終產量[12]。采用大疆S1000+型八旋翼無人機,機身凈質量4.4 kg,最大載物質量11 kg,續航時間15 min,搭載Parrot Sequoia多光譜相機。獲取4個波段的信息:綠光(green,G)波長550 nm,帶寬40 nm;紅光(red,R)波長660 nm,帶寬40 nm;紅邊光(red edge,RE)波長735 nm,帶寬10 nm;近紅外 光(near infrared,NIR)波 長790 nm,帶 寬40 nm。航拍時間為10∶00—11∶00,天氣晴朗無云、風力較小,設置飛行航線為S型,飛行3個架次,共計9條航線,單條航線長度1 140 m,覆蓋度70%×70%,飛行高度100 m,飛行速度5 m·s-1。采集完成后,將影像數據導入Pix4D Mapper進行拼接處理,影像重采樣為2 m分辨率。
1.2.2 土壤數據測定 在室內利用ArcGIS預設115個50 m間隔的采樣點,確定各點經緯度(圖1)。2017年4月23日春玉米播種前,利用GPS定位,使用TDR和電導率速測儀(Spectrum,美國)在田間測定土壤體積含水率和電導率,每個點位3次重復取平均值,同時用3點取樣混合法采集0—20 cm耕層土壤樣品,帶回實驗室測定土壤有 機 質(organic matter,OM)、堿 解 氮(alkali nitrogen,AN)、速效磷(available phosphorus,AP)、速效鉀(available kalium,AK)、pH[12]。

圖1 研究區域Fig.1Research area
1.2.3 玉米生長性狀及產量 2017年7月23日獲取115個點位的玉米生長數據。每個點位選取3株玉米,測定株高、莖粗,帶回實驗室烘干后測定地上生物量;采用手持式SPAD-502型葉綠素計(美能達,日本)測定上位展開葉的葉綠素相對含量SPAD值,每個點位選取6株取平均值;采用LAI-2200C冠層分析儀(美國,LI-COR)測定冠層葉面積指數。
2017年9月25日玉米成熟后獲取產量數據,每個點位收獲6顆玉米果穗,曬干脫粒稱重,折算成單元格產量。
無人機影像與土壤屬性數據類型不同,為了充分利用數據所包含的信息,分別選擇適用于無人機影像和土壤屬性的算法進行分區。
1.3.1 基于無人機影像的管理分區 從春玉米吐絲期多光譜影像中選取歸一化植被指數(normal differential vegetation index,NDVI)、優化土壤調節植被指數(optimized soil adjusted vegetation index,OSAVI)、歸一化差異紅色邊緣指數(normalized difference red edge index,NDREI)、比值植被指數(ratio vegetation index,RVI)、差 值 植 被 指 數(difference vegetation index,DVI)、土壤調節植被指數(soil-adjusted vegetation index,SAVI)、重歸一化植被指數(renormalized difference vegetation index,RDVI)共7種植被指數[13],利用ENVI 5.1提取115個點位周圍5 m圓形區域的植被指數均值,分析其與春玉米生長的相關關系。ECognition是基于目標信息的提取軟件,采用面向對象方法可提高數據的自動識別精度,為保證空間的緊致度和邊界的光滑度,在軟件中將形狀和光譜的權重設置為1∶1,將分割尺度范圍定在10~100內,以5為步長執行分割操作。通過計算每個分割尺度的平均分割評價指數,得出其最大值,確定為此參數下的最優分割尺度,計算公式如下。
由心電圖結果發現,FHCM親屬中G+/P-組Ⅰ、aVR、aVF導聯QRS時限及V1、V2導聯R波振幅均顯著大于對照組(G-/P-組),差異有統計學意義(P<0.05);V4~V6導聯R波振幅均小于對照組(G-/P-組),差異有統計學意義(P<0.05)。見圖1。T波改變(T波低平、T波倒置、T波雙向)發生率在兩組中差異有統計學意義(P<0.05),左室肥厚、異常Q波、碎裂QRS波發生率兩組比較,差異均無統計學意義(P>0.05)。兩組的心率、PR間期、QRS電軸、QTc間期比較,差異均無統計學意義(P>0.05)。見表1。


式中,σL為斑塊內像元值的標準差,n為所有像元的個數,CLi表示斑塊內像元i的像元值,根據光譜指數確定,CL表示斑塊內的像元均值;L為對象邊界長度,Li為與第i個相鄰對象的公共邊長度,C’Li為第i個相鄰斑塊的像元均值,N為與當前對象鄰接的斑塊個數;SEI為分割評價指數,ASEI為平均分割評價指數,A為整個目標田塊的斑塊總面積,Ai為第i個斑塊的面積,m為斑塊的總數量。
SEI通過對象間的同質性(公式1)和異質性(公式2)來判斷對象與鄰域間的異同,理想的結果應該與其標準差成反比,而與平均差分的絕對值成正比。SEIi表示第i個對象的分割評價指數。為了避免小面積對象對評價引起的不穩定性,使對象面積的大小對評價有著不同的貢獻,在引入了對象面積后計算研究區域所有對象的SEI平均值而得到ASEI。當ASEI指數達到最大值ASEImax時,對應的分割尺度為最優分割尺度。
1.3.2 基于土壤數據的管理分區 模糊均值聚類具有無監督、軟分類、信息豐富等優勢,在土壤領域被廣泛使用[16],MZA 1.0軟件采用模糊無監督分類方法,可對多個輸入變量進行分類[17],輸出的模糊性能指數(fuzzy performance index,FPI)和歸一化分類熵(normalized classification entropy,NCE)可幫助確定最佳分區數,當2個指標同時達到最小值時所對應的分類數即為最佳分類數。本研究參數設置:最大迭代次數=300,收斂標準=0.000 1,最小分區數=2,最大分區數=8,模糊指數=1.5。
1.3.3 分區結果驗證 使用軟件SPSS 18.0分別對2種分區方法下的產量、土壤屬性進行統計分析和交叉驗證,比較分區間差異,同時比較分區內變異系數,評價分區效果。
在2種方法下,固定的田間樣點會被劃入相同或不同的分區級別,分別統計2種分區方法的落點數,根據公式(5)計算樣點重合率,借此反映2種分區的空間重合關系。

式中,Mi為基于無人機影像的第i級分區的落點數,Si為基于土壤屬性的第i級分區的落點數。
描述性統計結果如表1所示,除pH外的其他屬性變異系數均大于10%,屬于中等程度變異,具有管理分區的必要,從均值來看土壤養分含量豐富,該地塊具有較大的節肥潛力,管理分區可提高土壤養分的利用率。

表1 土壤屬性描述性統計Table 1 Descriptive statistics of soil property
2.2.1 植被指數與春玉米生長性狀相關性分析如表2所示,歸一化紅邊植被指數NDRE與春玉米葉面積指數、地上生物量、株高均具有極顯著的相關性(P<0.01),適合用來表征玉米長勢的差異并作為數據源進行后續處理。基于NDRE進行多尺度分割(multiresolution segmentation),通過比較不同的分割尺度從而在最適宜的尺度層中提取相似的地物類別,分割的原則為分割對象內部同質性高,對象之間則有較高的異質性[14-15]。

表2 植被指數與春玉米生長性狀的相關性Table 2 Correlation between vegetation index and spring maize growth indicators
2.2.2 分割尺度分析 從分割評價指數曲線(圖2A)可以看出,隨著分割尺度的增大,斑塊數量逐漸減少,評價指數總體表現為先增大后降低,后趨于平穩的趨勢。當分割尺度為45 m時,評價指數達到最大值0.23,為該田塊的最優分割尺度,此時的分割斑塊數為38個,分割后的斑塊分布如圖2B所示。

圖2 基于無人機影像的多尺度分割結果Fig.2 Calculated Segmentation evaluation index and patchesbased on UAV images
2.2.3 分割結果聚類 在形狀要素的影響下,斑塊間的指數會具有一定的同質性,因此有必要對斑塊進行聚類以合并同質性較高的區域。對所有斑塊的NDRE均值進行模糊聚類,如圖3所示,FPI和NCE在聚類數為4時最小(圖3A),將38個斑塊合并為4個管理分區(圖3B),經合并后空間連續性更強,方便農田管理。將4個分區分別命名為M1、M2、M3和M4,M1主要分布于地塊北部的洼地以及地塊南部的坡背區域,M2分布于地塊西北部及中部的坡肩位置,M3和M4分別分布于地塊西部和東北部的坡趾區域。

圖3 基于無人機影像的管理分區結果Fig.3 Optimal management zones based on UAV images
2.3.1 主成分分析 為反映不同土壤屬性之間的關系及其對分區的影響,利用主成分分析對重疊信息進行壓縮,對7個土壤屬性進行主成分分析(表3)表明,前3個主成分特征值大于1,累計貢獻率達到88.28%,可作為后續分析的主成分。從因子載荷(表4)可以看出,PC1主要反映土壤有機質和堿解氮變化的信息,PC2中土壤有效磷和速效鉀的載荷占據主導地位,PC3中土壤水分載荷最大,反映了研究區土壤水分的信息。

表3 土壤屬性主成分分析Table 3 Principal Component analysis of soil variables and factor loadings

表4 主成分因子載荷Table 4 Principal component factor loads
2.3.2 模糊聚類結果分析 將前3個主成分作為輸入變量利用MZA1.0軟件進行模糊聚類,得到不同分區數下的FPI和NCE,如圖4A所示,當分區數為4時FPI和NCE值最小,因此最優分區數為4個,分別命名為S1、S2、S3和S4。S1主要位于東南部地勢較高的坡頂以及坡肩區域,S2分區主要位于地塊西北部的洼地,S3和S4分別分布于地塊西部和東北部(圖4B)。

圖4 基于土壤屬性的管理分區結果Fig.4 Optimal management zones based on soil properties
2.4.1 不同分區方法間樣點重合率分析 對比基于土壤和無人機影像的分區結果可以看出,2種結果不僅都為四級分區,在空間分布上也具有一定的相似性。經統計,2種結果的M1和S1分區包含田間樣點數分別為14、28,其中12個樣點同時分布于M1和S1中,重合率為40.00%(表5),重合區域主要為地塊南部的坡背區域;M2和S2分區的重合率為46.51%,M3和S3分區重合率為57.45%,農田西部的坡趾區域是M3和S3的一致區域;M4和S4的重合率為59.38%,重合區域主要分布于農田東北部。從整個農田地塊來看,基于無人機影像、土壤屬性2種分區結果的重合率為51.32%。

表5 不同分區方法間樣點重合率Table 5 Comparison of the percentage of points classified by different strategies
2.4.2 不同分區方法效果分析 為評價2種分區方法的效果,對各分區間春玉米產量、土壤養分和土壤含水率的差異進行了方差分析,結果如表6所示,基于無人機影像為數據源的分區,M1、M2、M3、M4各分區間產量差異顯著;4個分區土壤水分較為接近,其中M1分區內土壤水分變異系數較大,達到19.68%;土壤有機質、堿解氮、速效鉀、速效磷含量在M1、M2、M3分區間差異顯著,M4和M3無顯著差異,區內土壤養分變異性均有所降低,以土壤有機質為例,M1、M2、M3、M4有機質變異系數分別較農田整體變異性降低31.82%、14.36%、35.50%、30.10%,說明基于無人機影像的分區方法可以使管理分區內土壤養分差異減小,適宜在同一分區內實施相同管理。
基于土壤屬性為數據源的分區,S1、S2、S3分區間產量差異顯著,土壤養分排序均為S4>S3>S2>S1,各分區間差異均達到顯著水平,變異系數降低幅度較大,其中S1、S2、S3、S4分區土壤有機質的變異系數分別較農田整體變異降低43.48%、51.81%、63.58%、41.18%;土壤含水率排序為S2>S4>S3>S1。
分區管理的目的是因地制宜調整施肥等管理措施,優化農業資源利用。針對土壤差異進行分區是目前的常用方法,但因土壤數據不便獲取而應用局限性較大。作物長勢一方面反映土壤本底差異,另一方面也反映了其他環境要素如高程、地形、微氣候差異等,其所表征的綜合信息可為分區提供重要參考[18]。Chang等[19]利用手持式光譜儀Greenseeker采集了煙草的歸一化植被指數(NDVI),得到的5個管理區之間土壤養分和產量均差異顯著,為指導煙田施肥決策提供了依據。有研究表明,紅邊波段更能穿透作物冠層,對LAI變化敏感,歸一化紅邊植被指數(NDREI)更加適用于分區[20-21],本研究同樣表明,NDREI與作物生長具有較好的相關性,作為數據源將農田分為4個管理分區;從產量實測結果來看分區內的產量變異性均有所降低,分區間產量差異顯著(表6),這表明了利用無人機影像判斷作物產量的差異及進行分區的可行性。

表6 分區間春玉米產量及土壤養分和含水率差異Table 6 Variance analysis of spring maize yield and soil properties among management zones
除產量表現出明顯的分區效果外,分區內的土壤養分也向均一化方向發展,分區間土壤養分異質性明顯,表明當土壤樣品難以獲取時基于無人機影像的分區結果對于土壤養分管理同樣具有一定的參考意義,對比基于無人機影像數據和土壤數據劃分管理分區的結果也可以看出,三級分區M3與S3、四級分區M4與S4具有較高的相似性,進一步說明利用無人機影像進行分區和養分管理的可能性。值得注意的是,M1與S1、M2與S2重合度相對較低,土壤水分可能是造成2種方法差別的主要原因:水分過高或過低都會對作物生長產生脅迫,在坡耕地條件下土壤水分還具有較強的時空間動態變化特征[22-23],對于基于土壤數據的分區方法,只將1個時期的土壤水分作為參考要素可能會帶來誤差,本研究中S1和S2分區內土壤水分的變異性較小,而產量變異性相對較大,可見播前土壤水分和產量的變化并不一致。基于影像的分區結果對土壤水分管理也存在一定偏差,M1和M2分區內的產量變異性較小而土壤水分變異性較大,可能是1個時期的影像對水分盈虧狀態的評估不足,考慮到分區精度的提升,還需加強基于多期無人機影像的分區研究以提高分區的精度和指導性。
在實際應用中,基于無人機影像得到的4個分區具有較好的連續性,便于集中管理,也能夠滿足不同幅寬農機具進行作業路徑設計,在不同區域可以采取針對性的管理措施。對于長勢較差的M1分區,可適當增加施肥量以滿足作物對速效養分的需求;此外,對于研究區的旱作雨養農業需額外關注地形特征對作物長勢的影響,根據高程和坡度變化調節降水的再分配過程,M1分區中位于地塊北部的洼地可通過平整土地、提高深松深度等措施以避免土壤淹水,地塊南部的坡背區域水分易流失,可增大秸稈還田量改良土壤的持水性能。M2區域玉米籽粒達到了當季的目標產量,在后茬作物管理中可保持常規的施肥量。長勢較好的M3和M4區域可適當減少施肥以避免奢侈吸收,保證經濟產量并提高產投比,M3和M4分區內地形起伏相對較小,能夠保證適宜的土壤含水率。
本研究的主要目的是在土壤樣品不易獲取時為分區提供低成本、快捷的實施方法,對于分區后的農田管理措施多為定性分析,由于影響產量的多種因子間具有共線性,下一步的研究還需要結合多源數據探尋各因子對產量的影響閾值,并制定詳細完備的分區管理技術參數。