基于作物空間物候差異提取黃淮海夏玉米種植面積*

王雪婷，張 莎，鄧 帆**，張佳華**

基于作物空間物候差異提取黃淮海夏玉米種植面積*

王雪婷1，張 莎2,3，鄧 帆1**，張佳華2,3**

（1.長(zhǎng)江大學(xué)地球科學(xué)學(xué)院，武漢 430100；2.中國(guó)科學(xué)院大學(xué)，北京 100049；3.中國(guó)科學(xué)院遙感與數(shù)字地球研究所，北京 100049）

考慮大區(qū)域內(nèi)不同緯度間玉米物候差異，利用MODIS EVI時(shí)序曲線提取黃淮海夏玉米種植面積?；贚andsat影像和MOD13Q1數(shù)據(jù)集，提取參考區(qū)夏玉米MODIS EVI時(shí)間序列曲線；根據(jù)研究區(qū)內(nèi)農(nóng)業(yè)氣象站夏玉米生育期觀測(cè)數(shù)據(jù)，構(gòu)建夏玉米各物候期與緯度的關(guān)系，以緯度作為參數(shù)修正參考區(qū)夏玉米MODIS EVI時(shí)序曲線，獲取研究區(qū)夏玉米EVI標(biāo)準(zhǔn)時(shí)序曲線，結(jié)合平均絕對(duì)距離（MAD）和p-分位數(shù)法提取黃淮海平原夏玉米面積。結(jié)果表明，利用遙感影像提取的北京、天津、河北、河南以及山東夏玉米面積分別為125.3×103、162.6×103、2231.8×103、2963.6×103和2731.9×103hm2，各省提取精度均達(dá)到80%以上。在市級(jí)尺度上，決定系數(shù)R2為0.82，均方根誤差RMSE為147.8×103hm2；在縣級(jí)尺度上，決定系數(shù)R2為0.62，均方根誤差RMSE為17.7×103hm2。說(shuō)明利用本方法能夠準(zhǔn)確有效地提取大區(qū)域內(nèi)夏玉米種植面積，為其它農(nóng)作物在大范圍內(nèi)估計(jì)種植面積提供新思路。

夏玉米；種植面積提取；物候差異；多時(shí)相MODIS EVI；黃淮海平原

玉米是中國(guó)主要糧食作物，在農(nóng)業(yè)生產(chǎn)中占有重要地位，客觀、及時(shí)、準(zhǔn)確獲取玉米種植面積，對(duì)保障國(guó)家糧食安全和經(jīng)濟(jì)持續(xù)發(fā)展尤為重要[1?2]。黃淮海平原是夏玉米種植的主要區(qū)域，玉米種植面積占全國(guó)玉米面積的30%左右[3]，準(zhǔn)確估計(jì)該區(qū)域夏玉米種植面積對(duì)中國(guó)農(nóng)業(yè)生產(chǎn)具有重要意義。

遙感作為一種地球觀測(cè)技術(shù)，能夠有效地獲取地物的空間分布和光譜信息，具有更新快、易獲取、覆蓋范圍廣等優(yōu)勢(shì)，被廣泛應(yīng)用于農(nóng)作物的識(shí)別與分類(lèi)[4?6]。植被在遙感影像上呈現(xiàn)出不同的光譜特征，多時(shí)相遙感植被指數(shù)反映了植被光譜特征在時(shí)間維上的生長(zhǎng)變化，能夠有效區(qū)分不同植被類(lèi)型，實(shí)現(xiàn)作物面積提取[7?9]。顧曉鶴等[10]基于小波變換的MODIS（moderate resolution imaging spectroradiometer）和TM（thematic mapper）融合方法，構(gòu)建中分辨率尺度的NDVI（normalized difference vegetation index）時(shí)間序列數(shù)據(jù)集，提取河南原陽(yáng)縣玉米種植面積，提取精度達(dá)到88%以上；劉劍鋒等[11]利用MODIS NDVI數(shù)據(jù)集，結(jié)合玉米的物候信息，對(duì)伊洛河流域夏玉米進(jìn)行識(shí)別，識(shí)別精度達(dá)到90%以上；郭昱杉等[12]通過(guò)對(duì)比待分像元MODIS NDVI時(shí)序曲線與參考曲線的相似性，提取黃河三角洲地區(qū)玉米面積，提取精度達(dá)到85%以上。上述研究表明當(dāng)研究區(qū)域范圍較小時(shí)，即使不考慮區(qū)域間玉米的物候差異，基于多時(shí)相遙感影像提取的玉米種植面積也可以達(dá)到較高精度。

當(dāng)研究區(qū)域范圍較大時(shí)，區(qū)域間氣候（溫度、降水、輻射等）不同，玉米物候在不同區(qū)域呈現(xiàn)出顯著差異[13?15]。利用多時(shí)相遙感影像提取大區(qū)域內(nèi)夏玉米種植面積，若不考慮物候差異的影響，玉米面積提取結(jié)果精度較低[16?17]。因此，在大區(qū)域內(nèi)利用MODIS EVI時(shí)序曲線提取夏玉米面積時(shí)，需要考慮區(qū)域間不同氣候?qū)ψ魑锔魃诘挠绊?，以保障作物面積提取的精度。

在黃淮海平原，降水與日照時(shí)長(zhǎng)對(duì)夏玉米的生育進(jìn)程有一定影響，但不明顯；同期平均溫度變化與夏玉米物候之間呈現(xiàn)出顯著的負(fù)相關(guān)關(guān)系[14,18]；同時(shí)，緯度是溫度區(qū)域差異的決定要素之一，溫度呈現(xiàn)出緯向梯度變化，在該區(qū)不同緯度地帶的夏玉米物候差異明顯[19?21]。本研究將緯度表示的作物空間物候差異信息與多時(shí)相MODIS EVI影像結(jié)合，獲取研究區(qū)夏玉米EVI標(biāo)準(zhǔn)時(shí)序曲線方程，并計(jì)算研究區(qū)夏玉米EVI標(biāo)準(zhǔn)時(shí)序曲線與實(shí)際MODIS EVI時(shí)序曲線之間的平均絕對(duì)距離（MAD），利用p-分位數(shù)法，快速準(zhǔn)確實(shí)現(xiàn)黃淮海平原夏玉米種植面積的提取，為大區(qū)域內(nèi)作物面積快速有效提取提供新思路。

1 資料與方法

1.1 研究區(qū)和參考區(qū)

黃淮海平原（圖1a）位于31?41°N，110?123°E，包括河北（不含承德和張家口）、北京、天津、河南和山東。該區(qū)域地勢(shì)平坦，土壤資源豐富，是中國(guó)夏玉米種植的主要區(qū)域。區(qū)內(nèi)屬溫帶季風(fēng)性氣候，四季分明，全年降水量500～900mm，降水多集中于夏季，占全年降水量的70%，年平均溫度8～15℃，不同緯度間夏玉米物候差異明顯[19?21]。

圖1 研究區(qū)（a）和參考區(qū)（b）的空間分布

本研究選取種植結(jié)構(gòu)相對(duì)單一的河南省汝州市作為參考區(qū)（圖1b），該市位于河南省中部（33°56′? 34°20′N(xiāo)，112°31′?113°07′E），主要秋收作物為夏玉米，種植面積約占秋收作物的70%，選取該市作為參考區(qū)能夠減少遙感信息“同物異譜”與“同譜異物”的幾率。

1.2 數(shù)據(jù)來(lái)源及處理

1.2.1 Landsat7 ETM+影像

Landsat7 ETM+影像（https://glovis.usgs.gov）的空間分辨率為30m×30m，適合小區(qū)域高精度作物面積提取。根據(jù)參考區(qū)的天氣狀況、影像質(zhì)量以及夏玉米的物候特征，分別選取2012年6月9日（夏玉米播種期）、9月13日（乳熟期）和10月31日（收獲期）3期ETM+影像，用以提取參考區(qū)夏玉米種植面積。分類(lèi)前對(duì)3幅影像進(jìn)行去條帶、輻射定標(biāo)、基于FLAASH模型的大氣校正以及裁剪等預(yù)處理。

1.2.2 MODIS EVI數(shù)據(jù)

MODIS EVI數(shù)據(jù)來(lái)源于NASA的MOD13Q1數(shù)據(jù)集（https://search.earthdata.nasa.gov）。該數(shù)據(jù)集空間分辨率為250m×250m，時(shí)間分辨率為16d。夏玉米生長(zhǎng)季（6月1日?10月31日）共包含9個(gè)時(shí)相的數(shù)據(jù)，數(shù)據(jù)范圍覆蓋h26v04、h26v05、h27v04、h27v05和h28v05。先將影像進(jìn)行拼接，再對(duì)拼接后的影像進(jìn)行投影轉(zhuǎn)換處理，投影坐標(biāo)系為UTM_Zone_50N，地理坐標(biāo)系為CGS_WGS_1984，最后利用黃淮海和汝州的邊界數(shù)據(jù)分別裁剪拼接后的影像。

1.2.3 土地利用數(shù)據(jù)

利用2010年30m×30m空間分辨率的土地利用覆蓋數(shù)據(jù)（Globeland30，http://www.globallandcover. com）對(duì)所提取的參考區(qū)夏玉米面積進(jìn)行修正。提取土地利用類(lèi)型中的耕地（編碼為10）數(shù)據(jù)集，將耕地?cái)?shù)據(jù)集與參考區(qū)提取的夏玉米種植面積進(jìn)行掩膜處理，以排除其它土地利用類(lèi)型對(duì)參考區(qū)夏玉米提取結(jié)果的影響。

1.2.4 農(nóng)業(yè)氣象站數(shù)據(jù)

2012年夏玉米物候信息來(lái)源于中國(guó)氣象數(shù)據(jù)網(wǎng)（http://data.cma.cn/）的28個(gè)農(nóng)業(yè)氣象站（站點(diǎn)分布情況見(jiàn)圖1a）生育期觀測(cè)數(shù)據(jù)，包括各站點(diǎn)經(jīng)緯度、生育期名稱(chēng)和生育日期數(shù)據(jù)，將日期轉(zhuǎn)換為日序DOY（Ordinal day of a year），詳細(xì)信息見(jiàn)表1。

1.2.5 統(tǒng)計(jì)數(shù)據(jù)

收集2012年玉米種植面積統(tǒng)計(jì)數(shù)據(jù)[22?26]，包括北京市、天津市、河北省7個(gè)市、河南省15個(gè)市、山東省17個(gè)市共41個(gè)市級(jí)統(tǒng)計(jì)數(shù)據(jù)；北京市6個(gè)縣、天津市6個(gè)縣、河北省47個(gè)縣、河南省90個(gè)縣、山東省71個(gè)縣共220個(gè)縣級(jí)統(tǒng)計(jì)數(shù)據(jù)。

1.3 研究方法

黃淮海平原夏玉米種植面積的提取流程如圖2所示。由圖可見(jiàn)，整個(gè)提取方法可分為三大步驟。

第一步，在研究區(qū)內(nèi)選定種植作物相對(duì)單一的參考區(qū)，提取并重構(gòu)參考區(qū)夏玉米MODIS EVI時(shí)間序列曲線。選取2012年6月9日、9月13日和10月31日3期Landsat7影像，利用支持向量機(jī)方法[27?28]（Supported Vector Machine，SVM）提取參考區(qū)夏玉米種植分布區(qū)域，將提取結(jié)果與MODIS時(shí)序影像進(jìn)行掩膜處理，利用多時(shí)相MODIS EVI影像計(jì)算得到參考區(qū)夏玉米MODIS EVI時(shí)間序列曲線，進(jìn)而利用S–G[29?31]（Savitzky–Golay）濾波對(duì)參考區(qū)夏玉米MODIS EVI時(shí)間序列曲線進(jìn)行平滑重構(gòu)，得到參考夏玉米MODIS EVI時(shí)序曲線方程。

第二步，考慮作物物候空間差異性構(gòu)建研究區(qū)夏玉米EVI標(biāo)準(zhǔn)時(shí)序曲線的模擬方程。根據(jù)研究區(qū)內(nèi)28個(gè)農(nóng)業(yè)氣象站的夏玉米物候數(shù)據(jù)，分析夏玉米各生長(zhǎng)階段DOY與緯度之間的關(guān)系，得到夏玉米各生長(zhǎng)階段DOY隨緯度增加的延遲值（Lag），將緯度作為參數(shù)修正參考區(qū)夏玉米MODIS EVI時(shí)間序列曲線以得到研究區(qū)夏玉米EVI標(biāo)準(zhǔn)時(shí)序曲線方程。

第三步，根據(jù)研究區(qū)夏玉米EVI標(biāo)準(zhǔn)時(shí)序曲線方程提取夏玉米種植面積。計(jì)算研究區(qū)模擬夏玉米EVI標(biāo)準(zhǔn)時(shí)間序列曲線與實(shí)際MODIS EVI時(shí)間序列曲線的平均絕對(duì)距離[32]（MAD），并根據(jù)p–分位數(shù)法[33?34]設(shè)置各省的MAD閾值，將MAD值小于各省對(duì)應(yīng)的MAD閾值的柵格區(qū)域確定為夏玉米種植分布區(qū)域，利用這一規(guī)則，提取黃淮海夏玉米種植面積。

表1 2012年黃淮海區(qū)域農(nóng)業(yè)氣象站點(diǎn)夏玉米物候信息

Table 1 Phenological information of summer maize in Huanghuaihai agro-meteorological sites in 2012

注：“–”表示無(wú)觀測(cè)數(shù)據(jù)。

Note: “–”is no observation data.

2 結(jié)果與分析

2.1 參考區(qū)夏玉米MODIS EVI時(shí)間序列曲線的提取和重構(gòu)

2.1.1 參考區(qū)夏玉米種植面積提取

分析參考區(qū)秋收作物的物候特征（圖3）可知，參考區(qū)夏玉米6月上中旬處于播種期，在遙感影像上表現(xiàn)為裸地信息；花生和棉花的播種期在4月下旬?5月中旬之間，6月上中旬，在遙感影像上表現(xiàn)為植被信息，此時(shí)是利用遙感影像區(qū)分夏玉米與其它秋收作物種植區(qū)域的最佳時(shí)期。夏玉米在7?9月處于三葉–成熟期，在遙感影像上表現(xiàn)為植被信息；10月下旬夏玉米完成收割，冬小麥處于剛剛播種或未播種時(shí)期，夏玉米種植區(qū)域在遙感影像上表現(xiàn)為裸地信息。提取參考區(qū)遙感影像6月9日裸地、9月13日植被、10月31日裸地信息，做空間交集運(yùn)算獲取參考區(qū)夏玉米種植面積。具體步驟為：

（1）采用目視解譯方法從3幅Landsat7遙感影像上分別選取道路、植被、水體、居民地以及裸地的訓(xùn)練樣本，然后利用SVM方法[27?28]對(duì)6月9日、9月13日和10月31日3期影像進(jìn)行分類(lèi)。

（2）分別提取影像中6月9日裸地、9月13日植被、10月31日裸地信息，并做空間交集運(yùn)算，得到參考區(qū)夏玉米種植分布。

（3）利用30m×30m分辨率土地利用覆蓋數(shù)據(jù)中的耕地?cái)?shù)據(jù)集與獲取的參考區(qū)夏玉米分布進(jìn)行疊加分析，修正提取的參考區(qū)夏玉米種植面積。

圖2 夏玉米種植面積提取流程

圖3 參考區(qū)主要作物物候歷

Table基于Landsat-ETM+影像數(shù)據(jù)獲取參考區(qū)夏玉米種植面積（38.5×103hm2）空間分布，如圖4所示。根據(jù)統(tǒng)計(jì)年鑒[24]中汝州市夏玉米面積（43.4×103hm2）進(jìn)行驗(yàn)證，結(jié)果表明，提取精度為88%。表明在種植結(jié)構(gòu)相對(duì)單一的區(qū)域，利用夏玉米關(guān)鍵期不同光譜特征信息進(jìn)行交集運(yùn)算方法提取的夏玉米面積具有較高精度。

圖4 參考區(qū)夏玉米種植面積分布

2.1.2 參考區(qū)夏玉米MODIS EVI時(shí)間序列曲線重構(gòu)

將參考區(qū)提取夏玉米分布區(qū)域與時(shí)序MODIS EVI影像做掩膜處理，得到參考區(qū)夏玉米時(shí)序MODIS EVI空間分布，分別計(jì)算每一期夏玉米MODIS EVI平均值，將平均值作為參考區(qū)夏玉米的EVI值，得到參考區(qū)夏玉米MODIS EVI原始時(shí)間序列曲線（時(shí)序EVI）。利用S–G濾波[29?31]方法對(duì)時(shí)序EVI進(jìn)行平滑，得到參考區(qū)夏玉米平滑的MODIS EVI時(shí)間序列曲線（圖5）。采用高斯函數(shù)對(duì)參考區(qū)夏玉米平滑的MODIS EVI時(shí)間序列曲線進(jìn)行函數(shù)擬合，得到參考區(qū)模擬的夏玉米MODIS EVI時(shí)序曲線方程，即

式中，i代表DOY，EVIi代表在第i天參考區(qū)夏玉米模擬EVI值，模擬EVI函數(shù)各參數(shù)值見(jiàn)表2。

由圖5可見(jiàn)，參考區(qū)夏玉米MODIS EVI值隨時(shí)間推移呈現(xiàn)出先增后減的趨勢(shì)，夏玉米在第157?225天，處于播種–抽穗期，是夏玉米快速生長(zhǎng)階段，EVI逐漸升高；第225天前后進(jìn)入抽穗期，是夏玉米生長(zhǎng)最為旺盛時(shí)期，EVI值也達(dá)到峰值；第225?273天為抽穗–成熟期，夏玉米逐漸成熟，EVI值不斷下降。

表2 參考區(qū)夏玉米模擬EVI函數(shù)參數(shù)值

Table 2 Simulated EVI function parameters of summer maize in reference area

2.2 研究區(qū)夏玉米EVI標(biāo)準(zhǔn)時(shí)間序列曲線構(gòu)建

2.2.1 研究區(qū)夏玉米物候空間延遲特性分析

夏玉米物候期包括8個(gè)主要階段，分別是播種、出苗、三葉、七葉、拔節(jié)、抽穗、乳熟和成熟。根據(jù)28個(gè)農(nóng)業(yè)氣象站夏玉米物候數(shù)據(jù)，確定夏玉米物候歷與緯度之間的關(guān)系（圖6）。由圖6可以看出，夏玉米物候期從南到北呈現(xiàn)出不斷延遲的規(guī)律性。研究區(qū)內(nèi)夏玉米延遲系數(shù)（即緯度每增加1°，夏玉米物候期推遲的天數(shù)）在播種期為1.5868d·度?1，出苗期為1.6655d·度?1，三葉期為1.8416d·度?1，七葉期為1.632d·度?1，拔節(jié)期為1.314d·度?1，抽穗期為1.8124d·度?1，乳熟期為2.1881d·度?1，成熟期為2.4244d·度?1。為減少計(jì)算冗余，對(duì)夏玉米各物候期的延遲系數(shù)進(jìn)行整合，將整個(gè)物候期的延遲系數(shù)定為1.5d·度?1，并以參考區(qū)汝州市中心緯度（34.17°N）為基準(zhǔn)，計(jì)算黃淮海平原夏玉米物候的延遲值（Lag），即

2.2.2 研究區(qū)夏玉米EVI標(biāo)準(zhǔn)時(shí)間序列曲線的構(gòu)建

同一夏玉米生育期在黃淮海平原南北區(qū)域相差15d左右，參考區(qū)夏玉米MODIS EVI時(shí)間序列曲線不能代表整個(gè)研究區(qū)夏玉米MODIS EVI時(shí)間序列曲線，需要考慮不同緯度間所存在的物候差異。根據(jù)式（2）將黃淮海平原夏玉米物候延遲值（Lag）作為中間參數(shù)重新構(gòu)建參考區(qū)模擬的夏玉米MODIS EVI時(shí)序曲線方程，得到研究區(qū)夏玉米EVI標(biāo)準(zhǔn)時(shí)序曲線方程（式3），使緯度成為其中的一個(gè)參數(shù)，即

式中，i為DOY，Lag為延遲值，SEVIi為在第i天的夏玉米標(biāo)準(zhǔn)EVI值，SEVI函數(shù)各參數(shù)值見(jiàn)表2。

2.2.3 模擬EVI標(biāo)準(zhǔn)時(shí)間序列曲線的驗(yàn)證

在研究區(qū)的五省市隨機(jī)選取10個(gè)區(qū)縣（保證驗(yàn)證站點(diǎn)數(shù)占用于擬合函數(shù)站點(diǎn)數(shù)的30%[30]），對(duì)模擬EVI標(biāo)準(zhǔn)時(shí)間序列曲線進(jìn)行驗(yàn)證。根據(jù)2.1獲取參考區(qū)夏玉米MODIS EVI時(shí)間序列曲線的方法，計(jì)算10個(gè)驗(yàn)證區(qū)域?qū)嶋H夏玉米MODIS EVI時(shí)間序列曲線，并將其與模擬夏玉米EVI標(biāo)準(zhǔn)時(shí)間序列曲線進(jìn)行相關(guān)性分析，結(jié)果如表3所示。由表可見(jiàn)，在99%的置信水平下，實(shí)際夏玉米MODIS EVI時(shí)間序列曲線與模擬夏玉米EVI標(biāo)準(zhǔn)時(shí)序曲線的相關(guān)系數(shù)都大于0.9，表明根據(jù)式（3）模擬研究區(qū)夏玉米標(biāo)準(zhǔn)EVI時(shí)序曲線與實(shí)際夏玉米MODIS EVI時(shí)序曲線保持一致。

2.3 研究區(qū)玉米種植面積提取

2.3.1 計(jì)算夏玉米標(biāo)準(zhǔn)EVI與實(shí)際MODIS EVI間平均絕對(duì)距離

平均絕對(duì)距離[35]（MAD）的計(jì)算式為

式中，為平均絕對(duì)距離MAD值；mk為模板圖像第k時(shí)序像元個(gè)數(shù)，將標(biāo)準(zhǔn)EVI的每個(gè)柵格的像元作為模板，即mk為1；xik和xjk分別為第k時(shí)序?qū)?yīng)像元的實(shí)際MODIS EVI值和標(biāo)準(zhǔn)EVI值；n為時(shí)序總數(shù)（n=9）。MAD值越小，說(shuō)明夏玉米實(shí)際MODIS EVI值與模擬標(biāo)準(zhǔn)EVI值的相似性越大，為夏玉米的可能性越大；反之，為夏玉米的可能性越小[33]。計(jì)算得到研究區(qū)平均絕對(duì)距離（MAD）圖（圖7），由圖可見(jiàn)，2012年MAD的最小值為0，最大值為6。

表3 實(shí)際MODIS EVI曲線和模擬EVI標(biāo)準(zhǔn)曲線的相關(guān)分析

Table 3 Correlation analysis between the actual MOSIS EVI curve and the simulation standard EVI curve

注：**表示通過(guò)0.01水平的顯著性檢驗(yàn)。

Note:**is P＜0.01.

圖7 2012年黃淮海平原平均絕對(duì)距離（MAD）圖

2.3.2 利用p–分位數(shù)法確定各省MAD閾值

p–分位數(shù)法[33?34]（也稱(chēng)p–tile法）是Doyle于1962年提出的一種經(jīng)典閾值選取方法。該方法將研究對(duì)象分為目標(biāo)模式Po與背景模式Pb，用二者比值得到分位數(shù)P值，令其等于先驗(yàn)概率。根據(jù)先驗(yàn)概率找到合適的閾值T，使研究對(duì)象小于閾值T的為目標(biāo)模式，大于閾值T的為背景模式，具體公式為

P = Po/ Pb（5）

f(x, y)≤T （6）

f(x, y)＞T （7）

式中，Po為目標(biāo)模式；Pb為背景模式；f(x, y)為研究對(duì)象（MAD圖）；T為閾值。

計(jì)算各省統(tǒng)計(jì)玉米面積與總面積比值，將其作為分位數(shù)P值。根據(jù)分位數(shù)P值，對(duì)各省MAD圖的直方圖進(jìn)行累積，直到累積值的比例大于或等于P值，此時(shí)的MAD值即為最佳閾值T。根據(jù)p-分位數(shù)法，最終確定北京的最佳閾值T為0.85，天津?yàn)?.82，河北為0.66，河南為0.62，山東為0.65。當(dāng)北京、天津、河北、河南和山東的MAD值小于各省對(duì)應(yīng)的最佳閾值T時(shí)，認(rèn)為該柵格所占區(qū)域就是夏玉米種植區(qū)域，利用這一規(guī)則對(duì)黃淮海夏玉米種植分布進(jìn)行提取。

2.3.3 提取黃淮海平原夏玉米面積

研究區(qū)夏玉米種植面積提取結(jié)果見(jiàn)圖8。由圖可知，夏玉米在北京市、天津市以及河南省西部和南部分布較少，在黃淮海平原的其它區(qū)域分布較均勻。這是由于北京市和天津市地區(qū)多為市區(qū)，農(nóng)業(yè)種植區(qū)域較少；河南南部以水稻和小麥為主要種植作物[36]；河北、山東以及河南的中部和北部地勢(shì)平坦，光溫條件良好，是夏玉米的主要種植區(qū)域。與王紅營(yíng)等[37]提取的2012年華北地區(qū)夏玉米面積相比，本研究提取夏玉米分布區(qū)域與其在河北北部和南部、山東西北以及河南北部提取結(jié)果基本一致。

圖8 2012年黃淮海平原夏玉米種植面積分布圖

2.3.4 驗(yàn)證黃淮海平原夏玉米面積提取結(jié)果

利用MODIS EVI時(shí)序曲線提取的北京、天津、河北、河南以及山東夏玉米種植面積（表4）分別為125.3×103、162.6×103、2231.8×103、2963.6×103和2731.9×103hm2，河南省夏玉米種植面積最大，北京市最小，與統(tǒng)計(jì)年鑒數(shù)據(jù)相符。利用統(tǒng)計(jì)值與遙感提取值分析夏玉米面積提取的高估和低估情況[17]，提取面積大于統(tǒng)計(jì)面積表現(xiàn)為正即為高估，為負(fù)即為低估，由表4可知，各省市均表現(xiàn)為低估，其中河南省提取精度最高，達(dá)95.60%，河北省提取精度最低，為83.14%，研究區(qū)五省市提取結(jié)果精度均達(dá)到80%以上。

將夏玉米提取面積分別以市和縣為單位進(jìn)行匯總，并利用統(tǒng)計(jì)年鑒[22?26]中黃淮海41個(gè)市級(jí)玉米統(tǒng)計(jì)面積和220個(gè)縣級(jí)玉米統(tǒng)計(jì)面積對(duì)提取結(jié)果進(jìn)行驗(yàn)證。根據(jù)圖9可知，在市級(jí)尺度上，決定系數(shù)R2為0.82，均方根誤差RMSE為147.8×103hm2；在縣級(jí)尺度上，決定系數(shù)R2為0.62，均方根誤差RMSE為17.7×103hm2。與劉珺等[38]利用多時(shí)相MODIS EVI提取黃淮海平原夏玉米面積結(jié)果對(duì)比，本研究提取精度顯著提高。

表4 各省夏玉米提取面積與統(tǒng)計(jì)數(shù)據(jù)對(duì)比

Table 4 Comparison of extracted area and statistical data of summer maize for each province

注：高估或低估（%）=（提取面積?統(tǒng)計(jì)面積）/統(tǒng)計(jì)面積×100，100-Abs（高估或低估）即為精度（%）。下同。

Note: High(+) or low(?)(%)=(statistical data?extracted area)/ statistical data×100, and accuracy(%) is 100-Abs (high or low). The same as below.

圖9 黃淮海平原夏玉米遙感提取面積與統(tǒng)計(jì)面積相關(guān)性

隨機(jī)選取北京、天津、河北、河南以及山東的5個(gè)縣級(jí)區(qū)域（通州、薊縣、武安、沈丘、高密），作為面積提取抽樣驗(yàn)證區(qū)域，利用2.1.1提取參考區(qū)夏玉米種植區(qū)域方法，提取5個(gè)縣級(jí)驗(yàn)證區(qū)域的夏玉米面積。將基于30m×30m空間分辨率Landsat影像提取的夏玉米種植分布區(qū)域與基于250m×250m空間分辨率的MODIS EVI影像利用本研究方法獲取夏玉米種植分布區(qū)域進(jìn)行疊加分析，得到兩種遙感影像提取夏玉米種植分布重疊部分面積，計(jì)算重疊面積相對(duì)統(tǒng)計(jì)面積的高估和低估情況[17]。結(jié)果顯示（表5），除武安外，其它驗(yàn)證區(qū)域面積提取精度均達(dá)到80%以上，精度均值82.90%，滿足精度要求。

表5 Landsat和MODIS影像提取的夏玉米種植分布重疊部分面積與統(tǒng)計(jì)面積的對(duì)比

Table 5 Comparing the overlapping area of summer maize planting distribution extracted from Landsat and MODIS images with the statistical area

2.3.5 不同年份的提取結(jié)果

將MODIS EVI時(shí)序曲線與夏玉米空間物候差異結(jié)合提取2013年黃淮海夏玉米種植面積分布，北京、天津、河北、河南以及山東的夏玉米種植面積分別為102.2×103、155.6×103、2389.09×103、2975.72×103和3191.09×103hm2，與統(tǒng)計(jì)面積對(duì)比表明，五省市夏玉米提取面積精度均達(dá)到80%以上。以市級(jí)尺度為單位匯總2013年夏玉米種植面積提取結(jié)果（共41個(gè)市），并與統(tǒng)計(jì)面積[22?26]進(jìn)行對(duì)比驗(yàn)證，計(jì)算得到R2為0.79，RMSE為120.17×103hm2。2013年夏玉米的空間分布特征與2012年基本保持一致。可見(jiàn)，在不同年份間，利用夏玉米EVI標(biāo)準(zhǔn)時(shí)序曲線的方法提取大區(qū)域內(nèi)夏玉米種植面積的精度較高，表明該方法具有一定的普適性。

3 結(jié)論與討論

3.1 結(jié)論

（1）根據(jù)28個(gè)農(nóng)業(yè)氣象站生育期觀測(cè)數(shù)據(jù)，分析夏玉米各生長(zhǎng)階段DOY與緯度之間的關(guān)系，結(jié)果表明，在黃淮海平原，緯度每增加1°，夏玉米從出苗–成熟期各生長(zhǎng)階段延遲1.5d左右。

（2）將緯度作為參數(shù)構(gòu)建研究區(qū)夏玉米EVI標(biāo)準(zhǔn)時(shí)間序列曲線方程，利用10個(gè)驗(yàn)證點(diǎn)對(duì)擬合曲線方程進(jìn)行驗(yàn)證，驗(yàn)證結(jié)果表明，在99%的置信水平下，擬合的夏玉米EVI標(biāo)準(zhǔn)時(shí)序曲線與實(shí)際MODIS EVI時(shí)序曲線相關(guān)，說(shuō)明擬合曲線是可靠的。

（3）基于遙感影像提取的北京、天津、河北、河南和山東的夏玉米種植面積分別為125.3×103、162.6×103、2231.8×103、2963.6×103和2731.9×103hm2。將統(tǒng)計(jì)面積與提取的夏玉米面積對(duì)比，顯示在市級(jí)尺度上R2（RMSE）為0.82（147.8×103hm2），縣級(jí)尺度上R2（RMSE）為0.62（17.7×103hm2），說(shuō)明提取方法具有較好的精度。

3.2 討論

在大尺度上利用遙感影像時(shí)序曲線提取作物面積的過(guò)程中，由于區(qū)域間存在作物物候差異，會(huì)對(duì)面積提取結(jié)果產(chǎn)生一定影響[17,19,33,38]。本研究分析夏玉米各生育期（DOY）與緯度之間的關(guān)系，將MODIS EVI時(shí)序曲線和緯度表示的作物空間物候差異相結(jié)合，建立研究區(qū)夏玉米EVI標(biāo)準(zhǔn)時(shí)序曲線，使其與研究區(qū)內(nèi)各地夏玉米實(shí)際MODIS EVI時(shí)序曲線保持一致，并利用夏玉米EVI標(biāo)準(zhǔn)時(shí)序曲線提取大區(qū)域內(nèi)夏玉米種植面積。

由驗(yàn)證結(jié)果可以看出，在大區(qū)域上利用緯度修正的夏玉米標(biāo)準(zhǔn)EVI時(shí)序曲線提取夏玉米種植面積精度明顯提高，各省夏玉米面積提取精度均達(dá)到80%以上；隨機(jī)選取5個(gè)縣級(jí)區(qū)域提取的夏玉米種植面積，將其與Landsat影像提取面積的重疊部分與統(tǒng)計(jì)面積進(jìn)行對(duì)比，提取面積的平均精度值為83.26%，基本滿足研究精度要求，表明在大區(qū)域內(nèi)該方法能夠較準(zhǔn)確地提取出夏玉米種植面積；同時(shí)，基于2013年多時(shí)相MODIS EVI影像與夏玉米生育期觀測(cè)數(shù)據(jù)，構(gòu)建夏玉米標(biāo)準(zhǔn)EVI時(shí)序曲線提取黃淮海夏玉米種植面積，各省提取精度均達(dá)到80%以上，說(shuō)明該方法具有一定的普適性，可為大區(qū)域上其它農(nóng)作物面積提取提供一定技術(shù)參考。

本研究是在地勢(shì)平坦，光照條件良好的黃淮海平原進(jìn)行，只考慮了在該地區(qū)不同緯度間溫度差異造成的玉米物候差異，若在地形和氣候條件更復(fù)雜的區(qū)域提取夏玉米種植分布或者獲取更高精度的面積提取結(jié)果時(shí)，還需考慮海拔、降水、光照時(shí)長(zhǎng)以及種植習(xí)慣等對(duì)玉米各物候期的影響，以快速準(zhǔn)確提取大區(qū)域內(nèi)夏玉米種植面積。

[1]Chen S,Kung J K-S.Of maize and men:the effect of a new world crop on population and economic growth in China[J].Journal of Economic Growth,2016,21(1):71-99.

[2]趙久然,王榮煥.中國(guó)玉米生產(chǎn)發(fā)展歷程?存在問(wèn)題及對(duì)策[J].中國(guó)農(nóng)業(yè)科技導(dǎo)報(bào),2013,15(3):1-6.

Zhao J R,Wang R H.Development process,problem and countermeasure of maize production in China[J].Journal of Agricultural Science and Technology,2013,15(3):1-6.(in Chinese)

[3]中華人民共和國(guó)國(guó)家統(tǒng)計(jì)局.中國(guó)統(tǒng)計(jì)年鑒[M].北京:中國(guó)統(tǒng)計(jì)出版社,2013.

National Bureau of Statistics.China statistical yearbook[M].Beijing:China Statistics Press,2013.(in Chinese)

[4]Carr?o H,Gon?alves P,Caetano M.Contribution of multispectral and multitemporal information from MODIS images to land cover classification[J].Remote Sensing of Environment,2008,112(3):986-997.

[5]Mccullough I M,Loftin C S,Sader S A.High-frequency remote monitoring of large lakes with MODIS 500m imagery[J].Remote Sensing of Environment,2012,124:234-241.

[6]張莎,張佳華,白雲(yún),等.基于MODIS-EVI及物候差異免閾值提取黃淮海平原冬小麥面積[J].農(nóng)業(yè)工程學(xué)報(bào),2018,34(11):150-158.

Zhang S,Zhang J H,Bai Y,et al.Extracting winter wheat area in Huanghuaihai Plain using MODIS-EVIdata and phenology difference avoiding threshold[J].Transactions of the CSAE,2018,34(11):150-158.(in Chinese)

[7]Zhang M,Zhou Q,Chen Z,et al.Crop discrimination in Northern China with double cropping systems using Fourier analysis of time-series MODIS data[J].International Journal of Applied Earth Observation & Geoinformation,2008,10(4):476-485.

[8]郭建茂,施俊怡.基于HJ-CCD時(shí)序數(shù)據(jù)提取鎮(zhèn)江市水稻種植分布[J].中國(guó)農(nóng)業(yè)氣象,2014,35(2):228-234.

Guo J M,Shi J Y. Mapping rice paddy distribution in Zhenjiang city based on Multi-temporal HJ-CCD imagery[J]. Chinese Journal of Agrometeorology,2014,35(2): 228-234. (in Chinese)

[9]曹丹,景海濤,孫金珂,等.多時(shí)相MODIS指數(shù)提取東北三省水稻種植面積[J].測(cè)繪科學(xué),2018,43(7):54-60.

Cao D,Jing H T,Sun J K,et al.Extraction of rice area in Northeast China based on MODIS timing indices[J].Science of Surveying and Mapping,2018,43(7):54-60.(in Chinese)

[10]顧曉鶴,韓立建,王紀(jì)華,等.中低分辨率小波融合的玉米種植面積遙感估算[J].農(nóng)業(yè)工程學(xué)報(bào),2012,28(3):203-209.

Gu X H,Han L J,Wang J H,et al.Estimation of maize planting area based on wavelet fusion of multi-resolution images[J].Transactions of the CSAE,2012,28(3):203-209.(in Chinese)

[11]劉劍鋒,張喜旺.基于光譜和時(shí)相特征的夏玉米遙感識(shí)別[J].遙感技術(shù)與應(yīng)用,2016,31(6):1131-1139.

Liu J F,Zhang X W.Research on remote sensing identification method of maize based on seasonal rhythms and spectral characteristics[J].Remote Sensing Technology and Application,2016,31(6):1131-1139.(in Chinese)

[12]郭昱杉,劉慶生,劉高煥,等.基于MODIS時(shí)序NDVI主要農(nóng)作物種植信息提取研究[J].自然資源學(xué)報(bào),2017,32(10):1808-1818.

Guo Y S,Liu Q S,Liu G H,et al.Extraction of main crops in yellow river delta based on MODIS NDVI time series[J].Journal of Natural Resources,2017,32(10):1808-1818.(in Chinese)

[13]Hou P,Liu Y,Xie R,et al.Temporal and spatial variation in accumulated temperature requirements of maize[J].Field Crops Research,2014,158:55-64.

[14]Wang Z,Chen J,Li Y,et al.Effects of climate change and cultivar on summer maize phenology[J].International Journal of Plant Production,2016,10(4):509-525.

[15]李秀芬,馬樹(shù)慶,趙慧穎,等.基于WOFOST模型的內(nèi)蒙古河套灌區(qū)玉米低溫冷害評(píng)價(jià)[J].中國(guó)農(nóng)業(yè)氣象,2016,37(3):352-360.

Li X F,Ma S Q,Zhao H Y,et al.Evaluation on maize chilling injury based on WOFOST Model in Hetao irrigation region in Inner Mongol[J].Chinese Journal of Agrometeorology,2016,37(3):352-360.(in Chinese)

[16]劉珺,田慶久,黃彥,等.黃淮海夏玉米物候期遙感監(jiān)測(cè)研究[J].遙感信息,2013,28(3):85-90.

Liu J,Tian Q J,Huang Y,et al.Monitoring phenology of summer maize in Huang-Huai-Hai using remote sensing technology[J].Remote Sensing Information,2013,28(3):85-90.(in Chinese)

[17]姜亞珍,張瑜潔,孫琛,等.基于MODIS-EVI黃淮海平原冬小麥種植面積分帶提取[J].資源科學(xué),2015,37(2):417-424.

Jiang Y Z,Zhang Y J,Sun C,et al.MODIS-EVI based winter wheat planting information extraction of zoning on Huanghuaihai Plain[J].Resources Science,2015,37(2):417-424.(in Chinese)

[18]孟林,劉新建,鄔定榮,等.華北平原夏玉米主要生育期對(duì)氣候變化的響應(yīng)[J].中國(guó)農(nóng)業(yè)氣象,2015,36(4):375-382.

Meng L,Liu X J,Wu D R,et al.Responses of summer maize main phenology to climate change in the North China Plain[J].Chinese Journal of Agrometeorology,2015,36(4):375-382.(in Chinese)

[19]Liu Y,Xie R,Hou P,et al.Phenological responses of maize to changes in environment when grown at different latitudes in China[J].Field Crops Research,2013,144:192-199.

[20]Yan H,Zhen Y,Dong B W,et al.Effect of latitude on growth period and quality of maize(L.)[J].Chinese Agricultural Science Bulletin,2009,27(12):38-41.

[21]盧愛(ài)剛,龐德謙,何元慶,等.全球升溫對(duì)中國(guó)區(qū)域溫度緯向梯度的影響[J].地理科學(xué),2006,26(3):345-350.

Lu A G,Pang D Q,He Y Q,et al.Impact of global warming on latitudinal temperature gradients in China[J].Scientia Geographica Sinica,2006,26(3):345-350.(in Chinese)

[22]北京市統(tǒng)計(jì)局,國(guó)家統(tǒng)計(jì)局北京調(diào)查總隊(duì).北京統(tǒng)計(jì)年鑒[M].北京:中國(guó)統(tǒng)計(jì)出版社,2013.

Beijing Municipal Bureau of Statistics,Survey Office of the National Bureau of Statistics in Beijing. Beijing statistical yearbook [M]. Beijing: China Statistics Press,2013. (in Chinese)

[23]河北省人民政府辦公廳,河北省統(tǒng)計(jì)局.河北農(nóng)村統(tǒng)計(jì)年鑒[M].北京:中國(guó)統(tǒng)計(jì)出版社,2013.

The People's Government of Hebei Province Hebei Provincial Bureau of Statistics.HeBei rural statistical yearbook[M].Beijing:China Statistics Press,2013.(in Chinese)

[24]河南省統(tǒng)計(jì)局.河南統(tǒng)計(jì)年鑒[M].北京:中國(guó)統(tǒng)計(jì)出版社,2013.

Henan Province Bureau of Statistics. Henan statistical yearbook[M]. Beijing: China Statistics Press,2013. (in Chinese)

[25]山東省統(tǒng)計(jì)局,國(guó)家統(tǒng)計(jì)局山東調(diào)查總隊(duì).山東統(tǒng)計(jì)年鑒[M].北京:中國(guó)統(tǒng)計(jì)出版社,2013.

Shandong Provincial Bureau of Statistics,Survey Office of the National Bureau of Statistics in Shandong.Shandong statistical yearbook[M].Beijing:China Statistics Press,2013.(in Chinese)

[26]天津市統(tǒng)計(jì)局.天津統(tǒng)計(jì)年鑒[M].北京:中國(guó)統(tǒng)計(jì)出版社,2013.

Tianjin Municipal Bureau of Statistics.Tianjin statistical yearbook[M].Beijing:China Statistics Press,2013.(in Chinese)

[27]Mitra P,Shankar B U,Pal S K.Segmentation of multispectral remote sensing images using active support vector machines[J].Pattern Recognition Letters,2004,25(9):1067-1074.

[28]Mountrakis G,Im J,Ogole C.Support vector machines in remote sensing:a review[J].ISPRS Journal of Photogrammetry and Remote Sensing,2011,66(3):247-259.

[29]Savitzky A,Golay M J E.Smoothing and differentiation of data by simplified least squares procedures[J].Analytical Chemistry,1964,36(8):1627-1639.

[30]邊金虎,李?lèi)?ài)農(nóng),宋孟強(qiáng),等.MODIS植被指數(shù)時(shí)間序列Savitzky-Golay濾波算法重構(gòu)[J].遙感學(xué)報(bào),2010,14(4):725-741.

Bian J H,Li A N,Song M Q,et al.Reconstruction of NDVI time-series datasets of MODIS based on Savitzky-Golay filter[J].Journal of Remote Sensing,2010,14(4):725-741.(in Chinese)

[31]陳燕芬,牛振國(guó),胡勝杰,等.基于MODIS時(shí)間序列數(shù)據(jù)的洞庭湖濕地動(dòng)態(tài)監(jiān)測(cè)[J].水利學(xué)報(bào),2016,47(9):1093-1104.

Chen Y F,Niu Z G,Hu S J,et al.Dynamic monitoring of Dongting Lake wetland using time-series MODIS imagery[J].Journal of Hydraulic Engineering,2016,47(9):1093-1104.(in Chinese)

[32]Gu X H,Pan Y Z,He X,et al.Measurement of sown area of winter wheat based on per-field classification and remote sensing imagery[J].Journal of Remote Sensing,2010,14(4):789-805.

[33]Zhang J H,Feng L L,Yao F M.Improved maize cultivated area estimation over a large scale combining MODIS–EVI time series data and crop phenological information[J].ISPRS Journal of Photogrammetry & Remote Sensing,2014,94(94):102-113.

[34]李紅真,楊朝,劉恩海,等. P-tile與直方圖FCM結(jié)合的路面圖像分塊分割[J].計(jì)算機(jī)時(shí)代,2010,(8):32-34.

Li H Z,Yang C,Liu E H,et al.P-tile combined with histogram-based FCM for pavement image partitioning[J].Computer Era,2010,(8):32-34.(in Chinese)

[35]Gu X H,Pan Y C,He X,et al.Estimation of maize planting area through the fusion of multi-source images[J].Computer and Computing Technologies in Agriculture V,2012:470-477.

[36]劉新圣,孫睿,武芳,等.利用MODIS-EVI時(shí)序數(shù)據(jù)對(duì)河南省土地覆蓋進(jìn)行分類(lèi)[J].農(nóng)業(yè)工程學(xué)報(bào),2010,(s1):213-219.

Liu X S,Sun R,Wu F,et al.Land-cover classification for Henan Province with time-series MODIS EVI data[J].Transactions of the CSAE,2010,(s1):213-219.(in Chinese)

[37]王紅營(yíng),潘學(xué)鵬,羅建美,等.基于遙感的華北平原農(nóng)作物時(shí)空分布變化特征分析[J].中國(guó)生態(tài)農(nóng)業(yè)學(xué)報(bào),2015,23(9):1199-1209.

Wang H Y,Pan X P,Luo J M,et al.Using remote sensing to analyze spatiotemporal variations in crop planting in the North China Plain[J].Chinese Journal of Eco-Agriculture,2015,23(9):1199-1209.(in Chinese)

[38]劉珺,田慶久,黃彥,等.2000-2010年黃淮海平原夏玉米種植時(shí)空變化檢測(cè)[J].光譜學(xué)與光譜分析,2012,32(9):2534-2539.

Liu J,Tian Q J,Huang Y,et al.Dynamic monitoring of summer maize planting information for spatial and temporal variations in Huanghuaihai Plain during 2000-2010[J].Spectroscopy and Spectral Analysis,2012,32(9):2534-2539.(in Chinese)

Mapping the Cultivation Areas of Summer Maize Using Spatial Variations of Crop Phenology over Huanghuaihai Plain

WANG Xue-ting1, ZHANG Sha2,3, DENG Fan1, ZHANG Jia-hua2,3

(1. College of Earth Sciences, Yangtze University, Wuhan 430100, China; 2. University of Chinese Academy of Sciences, Beijing 100049; 3. Institute of Remote Sensing and Digital Earth, Chinese Academy of Sciences, Beijing 100049）

Considering the phenological differences of maize along with latitude over a large area, the planting area of summer maize was extracted by MODIS EVI time series curve over Huanghuaihai Plain in this paper. Based on Landsat images and MOD13Q1 dataset, the MODIS EVI time series curve of summer maize was obtained in a reference area. The phenology observation data from the agro-meteorological stations were collected to build the relationship between various summer maize growth stages and latitude. The latitude was used to correct the MODIS EVI time series curve of summer maize in the reference area, and then the standard EVI time series curve of summer maize was obtained. Using the mean absolute distance (MAD) method and p-tile algorithm, the summer maize planting area was estimated in the study area. The results showed that the areas of summer maize extracted from remote sensing images in Beijing, Tianjin, Hebei, Henan and Shandong were 125.3×103, 162.6×103, 2231.8×103, 2963.6×103and 2731.9×103ha, respectively. Meanwhile, the accuracy of extraction result in each province was above 80%. The determinant coefficient (R2) is 0.82 and the root mean square error (RMSE) is 147.8×103ha at city level, and 0.62 and 17.7×103ha at county level, respectively. It indicated that the method in present study has the ability of extracting the summer maize planting area effectively and provides a new idea for estimating the planting area of other crops in a large region.

Summer maize; Planting area extraction; Phenology difference; Multi-temporal MODIS EVI; Huanghuaihai Plain

10.3969/j.issn.1000-6362.2019.10.005

王雪婷,張莎,鄧帆,等.基于作物空間物候差異提取黃淮海夏玉米種植面積[J].中國(guó)農(nóng)業(yè)氣象,2019,40(10):647-659

2019?04?17

。E-mail：dengfan@yangtzeu.edu.cn; zhangjh@radi.ac.cn

國(guó)家重點(diǎn)研發(fā)計(jì)劃（2016YFD0300110）；國(guó)家自然科學(xué)基金（41871253；31671585）；山東自然科學(xué)基金重大基礎(chǔ)研究項(xiàng)目（ZR2017ZB0422）

王雪婷（1989?），女，碩士生，主要從事遙感技術(shù)應(yīng)用研究。E-mail：Wangxt1227@163.com