基于Sentinel-2A影像的縣域冬小麥種植面積遙感監測

王 蓉，馮美臣，楊武德，張美俊

（1.山西農業大學資源環境學院，山西 太谷 030801；2.山西農業大學農學院，山西 太谷 030801）

遙感是獲取作物播種信息的關鍵技術，在農情監測領域具有對地表信息獲取的覆蓋面廣、信息量大、周期短、受地面條件限制少、調查成本相對較低等明顯優勢[1]。在農作物遙感估產中，農作物種植面積的遙感估算是農作物產量預測的基礎和主要內容；準確而及時的農作物類型及其空間分布更新信息是優化種植結構的基本依據，為制定合理、有效的糧食宏觀調控措施、保證國家糧食安全提供科學支撐[2]。

農作物遙感識別分類與面積提取的研究主要考慮遙感數據源、特征變量和分類算法[3]3 個方面。到目前為止，適合開展大范圍高空間分辨率農作物種植區域識別與提取的遙感數據源主要有GF-1，Sentinel-2 和Landsat-8。田海峰等[4]使用多景GF-1全色多光譜（Pan Multiple Spectral，PMS）影像，對河南省滑縣冬小麥種植分布信息進行了提取研究，并取得良好效果。BELGIU 等[5]采用Sentinel-2 多光譜傳感器（Multi Spectral Sensor，MSI）數據，對時間加權動態時間規整（Time-Weighted Dynamic Time Warping，TWDTW）算法在農作物分類的應用效果進行研究。ZHANG 等[6]提出利用 Landsat-8 OLI 時間序列和物候參數的融合數據，在多云雨地區提取水稻種植面積。與MODIS 數據相比，GF-1，Sentinel-2 和Landsat-8 等影像數據具有高空間分辨率和高時間分辨率的雙重優勢，并在一定程度上減少了混合像元對作物分類結果的干擾，提高了農作物分類制圖精度與效率。

遙感光譜數據及其衍生的各種指數信息、紋理特征、時間序列閾值和變化信息常被國內外學者用作研究農作物遙感分類的特征變量。其中，光譜及其衍生指數的信息相對容易獲得也是最常用的[3]。SONOBE 等[7]在研究Sentinel-2A 多光譜傳感器數據計算植被指數的潛力中發現，植被指數在確定特定作物類型方面貢獻最大。史飛飛等[8]采用多源數據集進行了農作物分類信息提取研究，發現在NDVI 時間序列數據中融入高光譜數據可以提高分類精度。紋理信息作為區分特征，對高分辨率（米級或亞米級以下）影像數據的分類效果較好[9]。CHUANG等[10]將主成分、灰度共生矩陣紋理信息作為分類器的輸入變量，解譯出臺灣中部地區茶葉種植區域。

另外，選擇不同的算法對農作物面積的提取也有著較大的影響。隨機森林（Random Forest，RF）方法由于2 個隨機性（隨機且有放回地抽樣訓練集、隨機地從樣本特征維度中選取特征子集）的引入，在分類中不易陷入過擬合，具有分類精度高、泛化能力強、數據挖掘能力優異、抗噪聲能力良好等優點，對處理樣本數據分布不平衡等方面表現優良，成為當前機器學習領域的研究熱點[2]。GAO 等[11]研究探討了新型高光譜快照馬賽克相機在雜草和玉米分類方面的潛力，發現RF 模型表現總體上優于K-最近鄰（K-Nearest Neighbor，K-NN）模型。

以上研究都不同程度提高了作物識別精度，但綜合數字高程模型、遙感影像數據及其衍生的植被指數信息、紋理特征，利用單時相遙感影像進行雨養區冬小麥與灌溉區冬小麥空間分布信息提取方面的研究少有報道。

筆者在前人研究的基礎上，利用冬小麥生長關鍵期的Sentinel-2A 遙感影像數據，選擇紅邊歸一化植被指數、主成分分量與幾何紋理特征作為特征變量，研究基于隨機森林算法的冬小麥種植面積提取方法，并結合DEM 數據實現雨養區和灌溉區冬小麥的分類，以期為不同灌溉類型冬小麥長勢遙感監測以及產量估測提供科學依據。

1 材料和方法

1.1 研究區概況

研究區位于山西省運城市北部的聞喜縣，總面積 1 167.1 km2，地理坐標為 110°59′33″～111°37′29″E，35°9′38″～35°34′11″N。聞喜縣地處運城盆地與臨汾盆地的交界處，地形多樣，河谷、塬地、丘陵、山地共存。氣候屬暖溫帶大陸性季風氣候，晝夜溫差大，四季分明。農作物以冬小麥種植為主，屬南部冬麥區。

1.2 數據來源

1.2.1 Sentinel-2A 影像數據 Sentinel-2A 的多光譜成像儀含有13 個通道，主要應用于地表覆蓋變化監測、資源調查、近海域污染檢測、災害監測等，波段信息如表1所示。結合運城地區冬小麥與其他作物主要生育期，為有效提取冬小麥種植面積[12]，選擇使用2018年3月27日獲取的上層大氣反射的L1C 級別的Sentinel-2A 數據。

表1 Sentinel-2A 波段信息

1.2.2 地面調查 為獲取山西省運城市聞喜縣冬小麥的分布情況，2018年3月進行了實地調查。野外調查時采用手持GPS 獲取冬小麥的經緯度信息，調查路線覆蓋聞喜縣絕大部分鄉鎮（石門鄉除外）。在野外調查GPS 數據基礎上，結合目視解譯方法，確定訓練樣本數據集和驗證樣本數據集。

1.3 研究方法

1.3.1 數據預處理 數據下載自歐洲航空局的數據共享網站（http://scihub.copernicus.du/s2/#/home）。原始數據已進行過幾何校正處理，由于處理前各波段分辨率不同，因此，使用三次卷積內插法對各波段（不包括 Band1，Band8，Band9，Band10）進行重采樣，經處理后的各波段分辨率為20 m。

1.3.2 隨機森林分類算法 其通過對數據及特征變量進行隨機重采樣，構建多個CART 決策樹（不剪枝），最終采用多決策樹投票的方式確定數據的類別歸屬。它能夠處理具有高維特征的輸入樣本，又對過度擬合不敏感[13]，在生成的過程中就可以對誤差建立一個無偏估計。因而，對于遙感影像農作物面積提取具有很好的抗噪性能，在農田制圖研究方面取得了有效的分類結果[14-15]。

1.3.3 特征變量的選擇

1.3.3.1 主成分分析 主成分分析（Principal component analysis，PCA）又稱 K-L 變換，利用波段間的相互關系，在盡可能不丟失信息的條件下，用幾個綜合性波段代表多波段的原圖像，使處理的數據量減少。

1.3.3.2 植被指數 本研究選擇窄帶綠度植被指數（Narrowband greenness）中的紅邊歸一化植被指數（Red Edge Normalized Difference Vegetation Index，RENDVI）[16]。其計算公式如下。

式中，ρ705和 ρ750分別對應表2中的 Band5 和Band6 的反射率，RENDVI 值的范圍是 -1～1，一般綠色植被區的范圍是0.2～0.9。

1.3.3.3 紋理特征 紋理是指圖像色調作為等級函數在空間上的變化[17]。由HARALICK 提出的灰度共生矩陣（GLCM）方法是公認的紋理特征提取的有效方法，具有較強的適應能力和魯棒性[18]。該方法首先計算圖像的GLCM，然后由GLCM 導出描述紋理的二階統計特征[19]。在合適的空間分辨率遙感影像中，農作物區別于其他地物具有更為鮮明的空間紋理特征[20]。因此，加入紋理特征，有利于冬小麥分布區域的準確提取。

1.3.4 基于DEM 的地形分析 數字地面模型（Digiatl Terrain Model，DTM）中地形屬性為高程時稱為數字高程模型（Digital Elveatoin Model，DEM）。數字高程模型是指描述地球表面形態多種信息空間分布的有序數值陣列。根據文獻[21]，本研究采用三次卷積方法對研究區DEM 數據進行重新采樣，以減少精度損失。

1.3.5 精度評價方法 基于地面樣點數據驗證是精度驗證的主要手段之一，也是說明分類結果準確程度的指標之一。本研究以Kappa 系數、總體精度、制圖精度、用戶精度等4 項指標表述基于地面樣點數據精度驗證結果[22]。此外，由于區域冬小麥提取面積在產量估產中具有重要作用，因此，同樣需要評價冬小麥面積提取值與實際冬小麥面積值之間的比值。

2 結果與分析

2.1 特征變量的選擇

2.1.1 波段數據的統計特征分析 從表2可以看出，B8A 波段標準差（1 247.058 207）為所有波段中最大，其次是B11 波段、B7 波段，標準差分別為1 235.138 932 和1 155.616 279。標準差最大，說明該波段內地物的亮度取值距均值的離散程度最大，即地物間的差異可能表現最明顯，信息量最豐富。因此，B8A 波段包含的信息量最豐富，同時也進一步表明原始影像的B8A 波段在植被類型分類中具有顯著作用。

表2 原始影像波段基本信息量

2.1.2 主成分分析 從表3可以看出，第1 主成分分量（PC1）的信息占9 個波段總信息量的97.3%，第2 主成分占總信息量的2.23%，第3 主成分占總信息量的0.38%，前3 個主成分共占據了99.91%的信息量，其他成分對于影像信息只是噪音，所以，前3 個成分可以代表原始影像的數據信息。在構成第1 主成分的向量中，B8A 最高（為 0.44），占各波段特征向量和的15.3%。這說明第1 主成分中，B8A波段的貢獻最大，其次是B11（14.6%） 和B7（14.2%），說明在原始影像的9 個波段中，B8A 波段包含的地物信息量最豐富。

2.1.3 波段間相關性分析 從表4可以看出，B8A與B7 波段間相關系數最高（r＝0.999 4），其次是與B6 波段（r＝0.997 4） 和 B5 波段（r＝0.959 5），與B11 波段的相關系數為0.958 2。說明B8A 波段與這些波段數據彼此重疊較多。

2.1.4 特征變量選擇結果 選用B8A 波段作為專門提取紋理特征的波段，經多次試驗對比分析，選用3 多3 大小的移動窗口，利用灰度共生矩陣計算該波段的8 種紋理特征：均值、方差、同質性、對比度、差異性、熵（Entropy）、二階矩、相關性。

綜上所述，本研究提取3 類特征波段：紅邊歸一化植被指數（RENDVI）、主成分變換前3 個分量以及原始影像B8A 波段的紋理特征，共計12 個變量。各變量具體信息列于表5。

表3 研究區Sentinel-2A 各波段主成分分析結果

表4 研究區Sentinel-2A 各波段數據相關分析結果

表5 用于分類的特征變量信息

2.2 基于特征波段與隨機森林的冬小麥種植面積提取結果

將 2018年3月27日的 Sentinel-2A 影像輸入到隨機森林分類器[23]中，將前 3 個主成分（PC1，PC2，PC3）、紅邊歸一化植被指數（RENDVI）和 B8A 波段紋理特征組合進行分類提取，獲得聞喜縣冬小麥種植面積空間分布（圖1）。由圖1可知，冬小麥主要分布在中部地區，北部有少量種植分布，原因為聞喜縣地處運城盆地與臨汾盆地的交界夾槽處，三面環山，地勢西北、東南高，中間低，冬小麥適宜生長在地形平坦地區，林地等非冬小麥植被主要分布在東南山地；另外，由于丘陵塬地遍布縣境，冬小麥在西北部分布較為零散，這與實地調查情況相符。

2.2.1 特征波段對分類精度的影響 由表6可知，未加入特征波段，采用原始影像波段作為隨機森林分類器的輸入數據，分類結果的總體精度為94.72%，Kappa 系數為0.92；采用3 項特征波段后，分類結果的總體精度增加了3.39 百分點，Kappa 系數增加了0.05，這與分類結果（圖1）表現一致，說明綜合主成分分析、紅邊歸一化植被指數以及紋理特征后能提高遙感影像分類總體精度，同時冬小麥的分類用戶精度由87.91%提升到92.22%，說明引入特征波段有助于提高冬小麥種植面積提取的可信度。

表6 隨機森林算法條件下原始影像波段與特征波段分類精度對比

2.2.2 決策樹和特征變量的數量對精度的影響考慮到在隨機森林分類過程中，決策樹和特征變量的數量可能影響最終分類的精度，需研究決策樹和特征變量的數量對分類精度的影響。特征變量的數量對分類精度的影響甚微可忽略不計，即分類精度對特征變量數量的設置并不敏感[24]。因此，本研究僅分析決策樹的數量對分類精度的影響。由圖2可知，當保持特征變量數量不變時，將決策樹的數量分別設定為20，40，60，80，100，120，140，160，180進行試驗，得到的分類精度隨決策樹數量變化而改變。由圖2可知，在20～100 的范圍內，當決策樹數量增加，分類精度隨之增加，并且在決策樹數量為100 時分類精度達到最大值；當決策樹數量超過100 時，分類精度隨決策樹數量的增加而出現降低趨勢，且隨著決策樹數量增加，模型訓練時間也在增加。

2.3 不同分類方法分類結果對比

為進一步驗證本研究所采用的隨機森林分類方法的有效性，利用實地調查點信息將隨機森林算法分類結果與支持向量機算法、最大似然算法進行對比分析（圖3）。

從圖3可以看出，隨機森林算法（RF）、最大似然算法（ML）和支持向量機算法（SVM）分類結果在空間分布上大致趨于一致，但隨機森林算法分出的冬小麥地塊較為完整，且能分出絕大部分冬小麥種植區域，其分類總體精度和Kappa 系數均最高，分別為98.11%和0.97（圖4），其次是支持向量機算法（總體精度和Kappa 系數分別為96.3%和0.95），最大似然算法最差。這是由于最大似然算法分類器基于統計算法，計算給定像元屬于某一訓練樣本的似然度，根據似然度將像元歸并到似然度最大的一類當中，而在3月下旬聞喜縣冬小麥與蔬菜的光譜特征容易混淆，增大了利用似然度劃分地物像元的難度。支持向量機分類是一種建立在統計學習理論基礎上的機器學習方法，可以自動尋找那些對分類有較大區分能力的支持向量，由此構造出分類器，可以將類與類之間的間隔最大化，因而較最大似然法有更高的分類準確率，但是蔬菜與冬小麥的混淆光譜仍對其有一定的影響。而隨機森林算法通過對數據及特征變量進行隨機重采樣，采用多決策樹投票的方式確定數據的類別歸屬，在生成過程中就可以對誤差進行一個無偏估計，因此，其可信度最高，為92.22%。黃健熙等[25]使用MODIS 數據對黑龍江省主要農作物的種植面積提取進行研究，發現SVM算法分類精度優于RF，ML 算法的分類精度，這可能是由于MODIS 數據空間分辨率低，混合像元效應明顯，導致每棵決策樹的分類能力降低，而森林的分類效果（錯誤率）與每棵樹的分類能力有關；本研究中Sentinel-2A 數據空間分辨率明顯優于MODIS 數據，故RF 算法分類效果最佳。

2.4 雨養區與灌溉區冬小麥分類結果

由于不同灌溉條件下冬小麥有著不同的生育進程，因此，在冬小麥長勢監測和估產的研究中應該進行雨養區冬小麥與灌溉區冬小麥分類，以提高監測精度。在晉南地區冬小麥種植區域中，灌溉區冬小麥一般分布在海拔600 m 以下，坡度小于15°的平川區域，此區域有利于土壤水分涵養；反之，則為雨養區冬小麥分布區[26]。

從圖5和表7可以看出，灌溉區冬小麥主要分布在禮元鎮、東鎮鎮、侯村鄉、桐城鎮、郭家莊鎮、裴社鄉、河底鎮以及神柏鄉東南部等各鄉鎮沿涑水河流域的平川區域，面積為12 773.33 hm2，占總面積的34.92%；而雨養區冬小麥主要分布在聞喜縣西北部和東南部各鄉鎮丘陵垣臺區域，在西北部分布較為零散，在后宮鄉西部分布較為集中，面積為23 806.66 hm2，占聞喜縣冬小麥總面積的65.08%。同山西省統計年鑒數據（38 800 hm2）進行對比分析，聞喜縣冬小麥總提取面積的精度可以達到94.28%，雨養區為93.48%，灌溉區為95.8%。本研究得出的雨養區和灌溉區冬小麥提取精度明顯優于文獻[26]的研究結果（86.16%和86.15%），說明Sentinel-2A數據可以明顯降低混合像元效應，更適合區域冬小麥種植面積提取研究。

表7 聞喜縣雨養區、灌溉區冬小麥分類統計結果

3 結論

本研究以2018年冬小麥拔節期內的Sentinel-2A 數據為數據源，以山西省運城市聞喜縣為研究區，通過利用遙感影像主成分變換、波段特征提取技術，將主成分信息、波段空間紋理特征、紅邊歸一化植被指數以及大量調查樣本點輸入到隨機森林分類器中進行決策樹構建，并基于數字高程模型提取的坡度、高程信息，提取了2018年聞喜縣冬小麥的灌溉區與雨養區的空間分布信息。結果表明，Sentinel-2A 遙感數據適合作為縣域尺度冬小麥監測的數據源。Sentinel-2A 數據重訪周期為10 d，在農作物物候時期內更容易獲得質量較好地遙感影像，并且該數據影像空間分辨率較高，獲取成本低，能夠滿足大面積冬小麥監測的要求。采用隨機森林分類器對聞喜縣地物進行識別與分類，其分類總體精度達到98.11%，Kappa 系數達到0.97，優于同等條件下采用支持向量機法與最大似然算法。主成分分析、紋理特征和RENDVI 的引入可以提高單時相遙感影像對縣域冬小麥分類識別能力。2018年聞喜縣冬小麥遙感監測面積為36 580 hm2，提取精度達到94.28%。其中，雨養區冬小麥遙感監測面積為23 806.66 hm2，提取精度達到93.48%，占當年總面積的65.08%；灌溉區冬小麥遙感監測面積為12 773.33 hm2，提取精度為95.8%，占當年總面積的34.92%。