基于GEE的寶雞峽灌區(qū)耕地灌溉面積遙感監(jiān)測(cè)方法

徐超，呂婧妤，劉昱君，降亞楠,2*

(1. 西北農(nóng)林科技大學(xué)水利與建筑工程學(xué)院，陜西 楊凌 712100；2. 西北農(nóng)林科技大學(xué)旱區(qū)農(nóng)業(yè)水土工程教育部重點(diǎn)實(shí)驗(yàn)室，陜西 楊凌 712100)

耕地是農(nóng)業(yè)生產(chǎn)最基本的生產(chǎn)要素，耕地及其灌溉面積的變化會(huì)造成糧食生產(chǎn)波動(dòng)，進(jìn)而影響到糧食有效供給及糧食的安全水平[1].準(zhǔn)確識(shí)別耕地灌溉區(qū)域，對(duì)于糧食安全和農(nóng)業(yè)水資源的有效管理至關(guān)重要[2].然而，目前對(duì)于中國(guó)灌區(qū)實(shí)際耕地灌溉區(qū)域識(shí)別的研究卻非常有限，現(xiàn)有的灌溉遙感產(chǎn)品如MIRCA2000產(chǎn)品[3]和國(guó)際水資源管理協(xié)會(huì)(IWMI)的灌溉區(qū)域地圖[4]等，由于空間分辨率較低，統(tǒng)計(jì)出來(lái)的灌溉面積往往偏大[5]，難以滿足實(shí)際需要.

近年來(lái)，基于遙感的耕地識(shí)別與提取技術(shù)與方法已經(jīng)比較成熟[6-7]，但是根據(jù)衛(wèi)星觀測(cè)數(shù)據(jù)對(duì)灌溉區(qū)域進(jìn)行識(shí)別的技術(shù)還不夠成熟和完善[8]，是一項(xiàng)更加復(fù)雜且富有挑戰(zhàn)性的工作.目前國(guó)內(nèi)外的研究中，對(duì)于耕地灌溉區(qū)域的識(shí)別主要有決策樹(shù)算法[9]、光譜匹配法[4,10]、支持向量機(jī)算法[11]、隨機(jī)森林算法[12]、半自動(dòng)訓(xùn)練法[13]等方法.其中隨機(jī)森林算法是以決策樹(shù)為基礎(chǔ)分類(lèi)器的集成分類(lèi)器，與傳統(tǒng)的決策樹(shù)相比，有更強(qiáng)的泛化能力和更好的分類(lèi)效果，是目前最有效的機(jī)器學(xué)習(xí)算法之一，近年來(lái)在遙感影像分類(lèi)領(lǐng)域得到了廣泛的應(yīng)用[14-15].

現(xiàn)階段耕地灌溉面積遙感監(jiān)測(cè)研究主要存在兩方面的困難：一是樣本點(diǎn)的獲取，由于灌溉區(qū)域識(shí)別的特殊性，樣本點(diǎn)往往只能通過(guò)實(shí)地采樣調(diào)查獲得，而歷史樣本點(diǎn)的獲取則更加困難；二是目前研究中使用的指數(shù)或指標(biāo)[8]大多是只對(duì)土壤水分或者植物綠度比較敏感，例如歸一化植被指數(shù)(NDVI)[16-17].但是有研究表明NDVI對(duì)于作物的灌溉情況不夠敏感[17-18]，因此DEINES等[12,19]提出了土壤水分和植物綠度結(jié)合的新型灌溉指數(shù)：水調(diào)節(jié)綠色指數(shù)(WGI)，這種綜合性指數(shù)往往對(duì)于作物的灌溉情況更加敏感.

文中將綜合考慮以上問(wèn)題，結(jié)合實(shí)地采樣樣本點(diǎn)，并根據(jù)寶雞峽灌區(qū)實(shí)際，提出一種獲取歷史樣本點(diǎn)的思路和方法，基于GEE云平臺(tái)提供的隨機(jī)森林(RF)算法，結(jié)合水調(diào)節(jié)綠色指數(shù)WGI等，在中國(guó)西北干旱半干旱地區(qū)的典型灌區(qū)(陜西省寶雞峽灌區(qū))開(kāi)展耕地灌溉面積的遙感監(jiān)測(cè)研究，并分析2010—2020年間寶雞峽灌區(qū)耕地面積及其灌溉面積的變化.

1 研究區(qū)概況

寶雞峽灌區(qū)設(shè)計(jì)灌溉面積為197 707 hm2，是陜西省最大的國(guó)有灌區(qū)，位于關(guān)中平原西部地區(qū)；西起寶雞峽，東至涇河，南界渭河，北抵渭北高原腹地.灌區(qū)屬溫帶大陸性季風(fēng)氣候，降水時(shí)空分布不均，易發(fā)生干旱，作物生長(zhǎng)對(duì)灌溉需求較大.灌區(qū)主要灌溉水源為地下水、寶雞峽和馮家山水庫(kù).灌區(qū)主要作物為冬小麥、夏玉米、果樹(shù)(包括獼猴桃、蘋(píng)果、梨、桃等)和少量蔬菜等經(jīng)濟(jì)作物.寶雞峽灌區(qū)田塊面積較小，且種植結(jié)構(gòu)復(fù)雜.因此，文中將灌區(qū)的農(nóng)用地分為耕地和果樹(shù)2類(lèi)，并對(duì)耕地灌溉情況進(jìn)行識(shí)別.

2 數(shù)據(jù)源及預(yù)處理

2.1 衛(wèi)星遙感數(shù)據(jù)

首先對(duì)GEE云平臺(tái)提供的2010年Landsat5 TM和2020年的Landsat8 TM衛(wèi)星圖像進(jìn)行輻射定標(biāo)、大氣校正、影像去云、影像拼接等預(yù)處理.然后調(diào)用GEE云平臺(tái)內(nèi)SRTM DEM數(shù)據(jù)輔助訓(xùn)練，再調(diào)用平臺(tái)函數(shù)計(jì)算歸一化植被指數(shù)(NDVI)、歸一化水體指數(shù)(NDWI)、歸一化建筑指數(shù)(NDBI)以提高對(duì)植被、水體及建筑物的識(shí)別精度.除此之外，為了能更加精確地識(shí)別耕地灌溉區(qū)域，添加了增強(qiáng)植被指數(shù)(EVI)、綠色葉綠素植被指數(shù)(GCVI)、水調(diào)節(jié)綠色指數(shù)(WGI)3個(gè)指數(shù)，其中GCVI和WGI是對(duì)于作物灌溉特別敏感的指標(biāo)[8,12,19].GCVI，WGI計(jì)算公式為

(1)

WGI=NDWI·GCVI，

(2)

式中：NIR，Green分別為近紅外波段、綠波段的反射率.

2.2 樣本點(diǎn)數(shù)據(jù)

2020年實(shí)地采樣樣本點(diǎn)以耕地、果樹(shù)為主，包括典型的灌溉區(qū)與非灌溉區(qū)，共計(jì)427個(gè)，樣本點(diǎn)分布如圖1所示.由于2010年的灌溉情況已無(wú)法進(jìn)行實(shí)地考察，文中基于寶雞峽灌區(qū)實(shí)際，提出了一種獲取歷史樣本點(diǎn)的思路和方法：對(duì)2020年實(shí)地樣本點(diǎn)在2010年的光譜曲線和指數(shù)特征進(jìn)行分析，確定其灌溉情況，然后隨機(jī)添加目視解譯點(diǎn)，依據(jù)光譜特征篩選作為樣本補(bǔ)充，最后確定樣本個(gè)數(shù)為395個(gè)，樣本點(diǎn)分布如圖2所示.土地利用分類(lèi)樣本點(diǎn)可以在Google Earth Pro中對(duì)2010年和2020年的高分辨率影像進(jìn)行目視解譯獲得.

圖1 2020年樣本點(diǎn)分布示意圖Fig.1 Sample point distribution diagram in 2020

圖2 2010年樣本點(diǎn)分布示意圖Fig.2 Sample point distribution diagram in 2010

3 研究方法

3.1 隨機(jī)森林算法

隨機(jī)森林是BREIMAN[20]于2001年提出的一種以決策樹(shù)為基礎(chǔ)，結(jié)合Bootstrap aggregating集成學(xué)習(xí)理論與特征隨機(jī)選取的機(jī)器學(xué)習(xí)算法.隨機(jī)森林的核心即為隨機(jī)性，具有樣本和特征2個(gè)隨機(jī)性.有了這2個(gè)隨機(jī)性，在處理復(fù)雜且高維的數(shù)據(jù)時(shí)，隨機(jī)森林就不會(huì)出現(xiàn)過(guò)擬合現(xiàn)象，抗噪能力很強(qiáng)，而且隨機(jī)森林幾乎不需要人工干預(yù)，處理速度較快.

3.2 耕地灌溉區(qū)域識(shí)別

在進(jìn)行耕地灌溉區(qū)域識(shí)別前，首先進(jìn)行了土地利用分類(lèi)，以便將耕地與其他土地利用類(lèi)型區(qū)分開(kāi).而為了將灌溉作物和雨養(yǎng)作物分開(kāi)，需要準(zhǔn)確捕捉由于灌溉導(dǎo)致的作物光譜特性和指標(biāo)的差異.考慮到農(nóng)民進(jìn)行灌溉的時(shí)間通常是作物生長(zhǎng)的關(guān)鍵時(shí)期，可以根據(jù)作物生長(zhǎng)關(guān)鍵期的作物生長(zhǎng)狀況，以有效地識(shí)別與提取耕地灌溉區(qū)域.在寶雞峽灌區(qū)，小麥生長(zhǎng)的關(guān)鍵期3—5月份、玉米生長(zhǎng)的關(guān)鍵期7—9月份，均存在云量低于15%的影像，滿足使用要求.

基于GEE云平臺(tái)對(duì)實(shí)測(cè)耕地灌溉與非灌溉樣本點(diǎn)的光譜特征及時(shí)相特征進(jìn)行分析，主要采用樣本點(diǎn)的NDVI，NDWI，EVI，GCVI與WGI時(shí)間序列曲線.圖3為2個(gè)典型灌溉區(qū)域樣本點(diǎn)與2個(gè)典型非灌溉區(qū)域小麥種植樣本點(diǎn)的NDVI與WGI時(shí)間序列曲線，圖中t為1 a中的天數(shù).可以明顯看到，灌溉區(qū)域樣本點(diǎn)與非灌溉區(qū)域樣本點(diǎn)的曲線波動(dòng)趨勢(shì)基本一致，但是非灌溉區(qū)域的波動(dòng)幅度明顯小于灌溉區(qū)域，而且灌溉區(qū)域樣本相同時(shí)間的參數(shù)值均大于非灌溉區(qū)域樣本點(diǎn)，這說(shuō)明在干旱半干旱地區(qū)，灌溉對(duì)于作物的長(zhǎng)勢(shì)影響非常大.通過(guò)對(duì)比NDVI與WGI時(shí)間序列圖，可以發(fā)現(xiàn)并確認(rèn)DEINES等[12,19]提出的新型灌溉指數(shù)WGI放大了作物生長(zhǎng)關(guān)鍵期灌溉與非灌溉樣本點(diǎn)參數(shù)值的差距，使灌溉與非灌溉區(qū)域更易區(qū)分，可以有效描述作物的灌溉情況.

圖3 小麥種植期NDVI，WGI時(shí)間序列Fig.3 Time series of NDVI, WGI in wheat planting period

4 結(jié)果分析

4.1 土地利用分類(lèi)結(jié)果及精度驗(yàn)證

通過(guò)在GEE API平臺(tái)中進(jìn)行多次可視化試驗(yàn)對(duì)比，選取Landsat影像的B2，B3，B4，B5，B6，B7輔以NDVI，NDWI，NDBI，EVI以及DEM共11個(gè)波段作為RF算法的測(cè)試變量，以建設(shè)用地、水域、耕地、果樹(shù)4種土地利用類(lèi)型作為分類(lèi)體系，最終應(yīng)用于寶雞峽灌區(qū)的Landsat 影像分類(lèi)，分類(lèi)結(jié)果如圖4所示.在GEE云平臺(tái)通過(guò)Java代碼編程實(shí)現(xiàn)面積統(tǒng)計(jì)和精度評(píng)價(jià)，其中2010和2020年耕地面積分別為121 689，90 328 hm2.2010年分類(lèi)結(jié)果總體精度為92.86%，Kappa指數(shù)為0.900；2020年分類(lèi)結(jié)果總體精度為92.86%，Kappa指數(shù)為0.902，說(shuō)明分類(lèi)精度較高，分類(lèi)效果良好.

圖4 2010年和2020年寶雞峽灌區(qū)土地利用分類(lèi)結(jié)果Fig. 4 Land use classification results of Baojixia Irrigation District in 2010 and 2020

4.2 耕地灌溉區(qū)域識(shí)別結(jié)果及精度驗(yàn)證

在GEE云平臺(tái)計(jì)算指數(shù)NDVI，NDWI，EVI，GCVI，WGI作為分類(lèi)參數(shù)，應(yīng)用隨機(jī)森林算法實(shí)現(xiàn)分類(lèi)面積統(tǒng)計(jì).

將實(shí)測(cè)樣本點(diǎn)的80%作為訓(xùn)練樣本，20%作為驗(yàn)證樣本，然后根據(jù)篩選補(bǔ)充得到的樣本點(diǎn)使用同樣的分類(lèi)方法對(duì)2010年Landsat衛(wèi)星影像進(jìn)行耕地灌溉區(qū)域識(shí)別.統(tǒng)計(jì)識(shí)別提取得到的2010和2020年寶雞峽灌區(qū)耕地灌溉面積分別為104 211 和90 174 hm2，識(shí)別提取結(jié)果如圖5所示.

對(duì)提取結(jié)果進(jìn)行混淆矩陣精度評(píng)價(jià)，2010年總體精度達(dá)91.67%，Kappa系數(shù)為0.832，生產(chǎn)者精度為93.75%，用戶精度為88.24%；2020年總體精度達(dá)95.35%，Kappa系數(shù)為0.907，生產(chǎn)者精度為97.56%，用戶精度為93.02%.識(shí)別精度較高，耕地灌溉區(qū)域提取效果良好.

圖5 2010年和2020年寶雞峽灌區(qū)耕地灌溉區(qū)域提取結(jié)果Fig.5 Extraction results of cultivated land irrigation area in Baojixia Irrigation District in 2010 and 2020

4.3 耕地及其灌溉面積變化分析

查閱《2010陜西統(tǒng)計(jì)年鑒》[21]得知統(tǒng)計(jì)的寶雞峽灌區(qū)實(shí)際灌溉面積為104 000 hm2，與2010年提取結(jié)果比較接近.根據(jù)土地利用分類(lèi)結(jié)果，2010年耕地面積為121 689 hm2，而2020年耕地面積為90 328 hm2，比2010年少了31 361 hm2.這一方面是由于人口增加，城鎮(zhèn)擴(kuò)張，占據(jù)了部分耕地；另一方面由于獼猴桃、蘋(píng)果等果樹(shù)的推廣，果樹(shù)種植收益大于傳統(tǒng)的小麥、玉米，部分耕地轉(zhuǎn)變?yōu)楣麡?shù)種植區(qū).雖然研究區(qū)耕地面積減少了31 361 hm2，但是耕地灌溉面積僅僅減少了14 037 hm2.2010年研究區(qū)耕地灌溉率為85.6%，而2020年已經(jīng)提高到了99.8%.這與近年來(lái)政府加大投資力度、普及灌溉設(shè)施以及灌區(qū)水資源的合理調(diào)控十分相關(guān)，基本保證了每田皆有水灌；也可能與2020年寶雞峽灌區(qū)降水充沛，導(dǎo)致灌溉面積識(shí)別存在誤差有關(guān).

5 結(jié) 論

1) 基于GEE云平臺(tái)隨機(jī)森林算法，并結(jié)合新型灌溉指數(shù)WGI的耕地灌溉面積遙感監(jiān)測(cè)方法，對(duì)于耕地灌溉區(qū)域的識(shí)別總體精度較高，WGI對(duì)于作物的灌溉狀況十分敏感，可以有效識(shí)別作物的灌溉情況，在耕地灌溉面積遙感監(jiān)測(cè)的研究中可以推廣使用.

2) 基于寶雞峽灌區(qū)實(shí)際，提出了一種獲取歷史樣本點(diǎn)的思路和方法，根據(jù)實(shí)測(cè)樣本點(diǎn)的光譜曲線和指數(shù)特征可以獲取歷史位置的灌溉情況，擴(kuò)大了灌溉區(qū)域識(shí)別研究的時(shí)間尺度.但若要進(jìn)行更大范圍的耕地灌溉面積識(shí)別，則要考慮包括影像融合、氣候、作物種類(lèi)等許多其他因素，難度會(huì)大大增加.

3) 文中研究方法及成果可以為灌區(qū)水資源管理、灌區(qū)灌溉系統(tǒng)運(yùn)轉(zhuǎn)監(jiān)測(cè)、干旱的預(yù)防、灌區(qū)水資源優(yōu)化配置等提供技術(shù)支持和參考.