張 莎，張佳華※，白 雲(yún) ，姚鳳梅

（1. 中國(guó)科學(xué)院遙感與數(shù)字地球研究所數(shù)字地球重點(diǎn)實(shí)驗(yàn)室，北京 100094；2. 中國(guó)科學(xué)院大學(xué)地球與行星科學(xué)學(xué)院，北京 100049）

0 引 言

冬小麥?zhǔn)侵袊?guó)3大糧食作物之一[1-2]。獲取準(zhǔn)確冬小麥種植面積對(duì)于糧食估產(chǎn)和保障國(guó)家糧食安全有重要意義[1]。黃淮海平原是中國(guó)最大的平原區(qū)，是重要的冬小麥生產(chǎn)基地[3]。因此，準(zhǔn)確估計(jì)黃淮海地區(qū)冬小麥種植面積有重要的意義。

遙感作為一種地球觀測(cè)技術(shù)，在大范圍面積和資源監(jiān)測(cè)上有著比傳統(tǒng)手段更為明顯的優(yōu)勢(shì)[2]，已被廣泛地應(yīng)用于作物面積提取[4-5]。前人在黃淮海區(qū)域冬小麥面積的遙感提取方面已開(kāi)展了大量工作。研究表明，當(dāng)利用中低分辨率影像提取中國(guó)北方冬小麥面積時(shí)，對(duì)綠度指數(shù)進(jìn)行圖像分割就能取得較好的效果[6]。葛廣秀等[7]使用歸一化植被指數(shù) NDVI（normalized difference vegetation index）密度分割方法提取了江蘇省沭陽(yáng)縣的冬小麥面積，精度達(dá)92.37%。根據(jù)不同省市冬小麥的物候特征，使用關(guān)鍵生育期時(shí)相的 NDVI影像并設(shè)置閾值構(gòu)建各省或各市的冬小麥提取模型，都可以達(dá)到較好的精度。黃青等[8]根據(jù)江蘇省冬小麥物候特征，基于MODIS NDVI數(shù)據(jù)構(gòu)建了適合于江蘇省的冬小麥提取模型，總體提取精度達(dá)78%以上；楊小喚等[9]根據(jù)冬小麥的生育期特點(diǎn)，使用關(guān)鍵時(shí)相的 NDVI影像并設(shè)置閾值提取了北京市包含冬小麥在內(nèi)的多種作物種植面積，總體精度達(dá)95%以上。Qiao等[10]也采用 NDVI閾值法較為準(zhǔn)確地提取了河南省冬小麥種植面積。

可以看出，當(dāng)研究區(qū)范圍小于或等于省級(jí)尺度時(shí)，研究者即便使用同一套閾值參數(shù)，即忽略冬小麥在區(qū)域間的物候差異，也可以較為準(zhǔn)確地提取冬小麥面積；而當(dāng)研究范圍緯度跨越較大時(shí)，如同一省份南北方向延伸較大或包含多個(gè)省份時(shí)，物候差異則不能忽略，研究者往往會(huì)選擇分省或分市進(jìn)行閾值設(shè)置，以區(qū)分冬小麥在不同地區(qū)的的物候差異。陳健等[11]基于時(shí)序中分辨率成像光譜儀-增強(qiáng)型植被指數(shù)（MODIS-EVI，moderate resolution imaging spectroradiometer-enhanced vegetation index）數(shù)據(jù)提取河北?。?6°01′N—42°37′N）冬小麥種植面積時(shí)，就采用了分市設(shè)置閾值的方法。張佳華等[2]根據(jù)黃淮海各省冬小麥物候特征，采用EVI閾值模型，分別對(duì)研究區(qū)內(nèi)5個(gè)省市的冬小麥生育期共10個(gè)EVI閾值進(jìn)行設(shè)置，最終得到黃淮海平原典型區(qū)的冬小麥種植面積，在省級(jí)尺度驗(yàn)證整體精度可達(dá) 98%。當(dāng)研究區(qū)包含較多省份時(shí)，這種方法無(wú)疑增加了工作量。還有學(xué)者根據(jù)緯度進(jìn)行水平分帶[3]和傾斜分帶[12]以區(qū)分冬小麥在區(qū)域間的物候差異，根據(jù)各緯度帶冬小麥生育關(guān)鍵期不同而使用不同時(shí)相的影像，各緯度帶的提取模型使用相同的閾值。

采用植被指數(shù)閾值法，無(wú)論采用分區(qū)方案，還是根據(jù)緯度進(jìn)行水平分帶或傾斜分帶方案，都無(wú)法避免對(duì)植被指數(shù)閾值進(jìn)行確定。研究者通常依據(jù)作物在各區(qū)的物候特征根據(jù)經(jīng)驗(yàn)值不斷調(diào)整來(lái)確定閾值，分區(qū)越小提取精度相應(yīng)越高[8,13-14]。再者，不同區(qū)域之間，或同一區(qū)域不同年份之間，作物長(zhǎng)勢(shì)有所差異則閾值也要相應(yīng)調(diào)整。閾值的不確定性增加了冬小麥種植面積提取的不確定性。因此，找到一個(gè)避免閾值設(shè)置的提取特征十分有必要。潘學(xué)鵬等[15]利用MODIS數(shù)據(jù)，以基于復(fù)種指數(shù)的決策樹(shù)分類方法提取了黃淮海平原主體地區(qū)的冬小麥面積，京津冀地區(qū)長(zhǎng)時(shí)間序列的市級(jí)尺度驗(yàn)證結(jié)果 R2為0.885。王學(xué)等[16]利用MODIS數(shù)據(jù)和TM數(shù)據(jù)，根據(jù)小麥生育期內(nèi)波峰和波谷的特征構(gòu)建提取模型，提取了黃淮海平原主體部分冬小麥種植面積。王云峰等[17]使用MODIS NDVI數(shù)據(jù)，根據(jù)冬小麥區(qū)別于同期作物的物候特征，即在收獲期植被指數(shù)突降（此階段其他作物均為生長(zhǎng)季），采用兩景影像提取山西省運(yùn)城市冬小麥面積，并取得了較為合理的結(jié)果。孫振蓉[18]也利用了類似的方法提取了京津冀地區(qū)的冬小麥種植面積。為了減小閾值設(shè)置給提取結(jié)果帶來(lái)的不確定性，同時(shí)在大區(qū)域范圍內(nèi)既能考慮物候差異又能快速提取冬小麥面積，本文將引入波形特征。

鑒于此，本研究以黃淮海平原為研究區(qū)，擬提出一種既能考慮區(qū)域間物候差異又具有普適性的提取方法。首先，采用基本氣候要素（氣溫、降水和輻射）模擬冬小麥的播種期和成熟期；然后，使用濾波后的MODIS EVI數(shù)據(jù)，利用差分法逐像元提取播種期至成熟期內(nèi)的 EVI峰值頻數(shù)；逐像元計(jì)算成熟期前后 EVI的下降速率；最后，利用統(tǒng)一的提取規(guī)則提取研究區(qū)內(nèi)冬小麥種植面積并用統(tǒng)計(jì)數(shù)據(jù)驗(yàn)證。旨在為大區(qū)域范圍內(nèi)快速、準(zhǔn)確地提取冬小麥種植面積提出一套既能滿足作物種植面積監(jiān)測(cè)精度需求又具有普適性、能在大范圍區(qū)域使用的提取模型，服務(wù)于國(guó)家冬小麥面積監(jiān)測(cè)的業(yè)務(wù)化運(yùn)行。

1 材料與方法

1.1 研究區(qū)概況

本文研究區(qū)為黃淮海平原，位于中國(guó)東部地區(qū)（32°~42°N，113°~120°E），覆蓋北京市、天津市、河北省、河南省、山東省、山西省、湖北省、安徽省和江蘇省共計(jì) 9個(gè)省市。黃淮海平原是中國(guó)最大的平原區(qū)，也是中國(guó)冬小麥的主產(chǎn)區(qū)。該區(qū)范圍大，在區(qū)域范圍內(nèi)輻射和氣溫等氣候要素有明顯的區(qū)域變化，在植被上變現(xiàn)出明顯的物候差異。

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

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

所采用的MODIS EVI 數(shù)據(jù)為MOD13Q1的C6數(shù)據(jù)（https://search.earth data.nasa.gov/），空間分辨率為250 m。本文獲取從2009年第241天至2010年第209天每16 d合成一期的影像，共22期。這些影像涵蓋了研究區(qū)（h26v04，h26v05，h27v04，h27v05，h27v06，h28v05，h28v06）冬小麥的全部生育期。MODIS EVI數(shù)據(jù)有效值為-3 000～10 000，比例因子為0.000 1。為去除非植被像元干擾，本文提取EVI大于0的像元以進(jìn)行下一步研究。

1.2.2 土地利用數(shù)據(jù)和旱地?cái)?shù)據(jù)

中國(guó)2010年土地利用數(shù)據(jù)來(lái)源于中國(guó)科學(xué)院資源環(huán)境科學(xué)數(shù)據(jù)中心(http://www.resdc.cn)，空間分辨率為1 km。將該數(shù)據(jù)重采樣為250 m，并用矢量圖裁剪，得到研究區(qū)土地利用（圖1）。除耕地按二級(jí)分類體系顯示外，其余土地利用類型均將編碼合并處理按一級(jí)分類體系顯示。本文提取了土地利用類型中的旱地（編碼為121，122，123和124，分別為山地旱地、丘陵旱地、平原旱地和大于 20°坡度旱地）。利用旱地提取結(jié)果對(duì)所提取冬小麥面積進(jìn)行掩膜處理，以排除其他土地利用類型對(duì)冬小麥提取結(jié)果的影響。

1.2.3 農(nóng)業(yè)氣象站點(diǎn)數(shù)據(jù)

研究區(qū)內(nèi)農(nóng)業(yè)氣象站點(diǎn)數(shù)據(jù)來(lái)源于中國(guó)氣象數(shù)據(jù)網(wǎng)（http://data.cma.cn/），使用了農(nóng)氣資料中的AGME_AB2_CHN_TEN數(shù)據(jù)集（中國(guó)農(nóng)作物生長(zhǎng)發(fā)育和農(nóng)田土壤濕度旬值數(shù)據(jù)集）。本文使用了140個(gè)冬小麥站點(diǎn)，其空間分布見(jiàn)圖1。所用農(nóng)業(yè)氣象站點(diǎn)數(shù)據(jù)包括各站點(diǎn)的經(jīng)緯度、作物類型、生育期名稱和生育期日期數(shù)據(jù)。所使用的生育期僅包括冬小麥播種期和成熟期，并將其日期轉(zhuǎn)換為日序 DOY（day of year）。

圖1 研究區(qū)土地利用及所用農(nóng)業(yè)氣象站點(diǎn)分布Fig.1 Land cover of study area and distribution of agro-meteorological sites in study area

1.2.4 氣象再分析資料數(shù)據(jù)

本文采用空間分辨率為0.125°的ERA-Interim日值數(shù)據(jù)集獲取氣象要素，包括氣溫、降水和輻射。ERA-Interim是歐洲中期數(shù)值預(yù)報(bào)中心提供的全球最新的大氣再分析產(chǎn)品，較之前的產(chǎn)品精度更高[19-20]，已被廣泛應(yīng)用于氣候變化等研究[21-23]。ERA-Interim數(shù)據(jù)采用nc格式存儲(chǔ)，本研究利用IDL代碼將其讀取為img格式，并根據(jù)農(nóng)業(yè)氣象站點(diǎn)的位置提取其所在柵格的值作為站點(diǎn)的氣象要素值。

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

冬小麥播種面積統(tǒng)計(jì)數(shù)據(jù)來(lái)源于各省市統(tǒng)計(jì)年鑒[24-30]。本文收集的市級(jí)統(tǒng)計(jì)數(shù)據(jù)包括北京市、天津市、河北省9個(gè)市（承德市和張家口市除外）、河南省 18個(gè)市、山東省17個(gè)市、山西省10個(gè)市（大同市除外）、湖北省16個(gè)市（襄陽(yáng)市除外）和安徽省17個(gè)市共89個(gè)地級(jí)市的2010年和2015年冬小麥種植面積數(shù)據(jù)；縣級(jí)統(tǒng)計(jì)數(shù)據(jù)包括河北省113個(gè)縣、河南省10個(gè)縣、山西省77個(gè)縣和湖北省55個(gè)縣共255個(gè)縣的2010年冬小麥種植面積數(shù)據(jù)。

1.3 研究方法

提取冬小麥的流程如圖 2所示。首先，利用隨機(jī)分布的98個(gè)農(nóng)業(yè)氣象站的觀測(cè)播種期與成熟期數(shù)據(jù)和基本氣候要素（氣溫、降水和輻射），采用逐步回歸法，構(gòu)建多元線性回歸模型模擬冬小麥播種期和成熟期，并利用其他42個(gè)站點(diǎn)驗(yàn)證模擬結(jié)果；其次，對(duì)時(shí)序EVI數(shù)據(jù)進(jìn)行 Savitzky-Golay（S-G）濾波；然后，根據(jù)模擬得到的播種期和成熟期，截取每個(gè)像元處于播種期至成熟期之間的EVI序列，計(jì)算截取后得到的EVI序列的峰值頻數(shù)；同時(shí)，逐像元獲取成熟期前后的EVI，計(jì)算每個(gè)像元成熟期前后EVI發(fā)生突變的情況（Slope）；最后，提取峰值頻數(shù)為2[15]和Slope小于-0.02[17-18]的像元，并與黃淮海平原旱地?cái)?shù)據(jù)取交集，得到冬小麥的提取結(jié)果，同時(shí)利用統(tǒng)計(jì)數(shù)據(jù)對(duì)提取結(jié)果進(jìn)行驗(yàn)證。

圖2 冬小麥面積提取流程圖Fig.2 Flowchart of winter wheat area extraction

1.3.1 多元線性逐步回歸法

區(qū)域間不同的水熱組合條件是造成地區(qū)間物候差異的根本原因。本文采用多元線性逐步回歸法[31]擬合冬小麥播種期和成熟期日序DOY與氣溫、降水和輻射之間的關(guān)系。研究區(qū)內(nèi) 140個(gè)冬小麥農(nóng)業(yè)氣象站點(diǎn)，隨機(jī)選擇其中98個(gè)站點(diǎn)用于構(gòu)建播種期和成熟期DOY與氣候要素的關(guān)系，并利用剩余42個(gè)站點(diǎn)（保證驗(yàn)證站點(diǎn)數(shù)占總站點(diǎn)數(shù)的30%[32]）的生育期數(shù)據(jù)對(duì)擬合關(guān)系進(jìn)行驗(yàn)證。

1.3.2 S-G濾波

S-G濾波擬合算法是一種基于平滑時(shí)間序列數(shù)據(jù)和最小二乘原理的卷積算法[33-35]

式中為擬合EVI值；Yj+i為時(shí)間序列中的第(j+i)個(gè)EVI值； Ci為第i個(gè)EVI值濾波時(shí)的卷積系數(shù)；msg為濾波窗口大小的一半；Nsg為濾波器長(zhǎng)度，等于滑動(dòng)數(shù)組的寬度(2msg+1)。本文利用IDL（interactive data language）實(shí)現(xiàn)時(shí)間序列EVI影像的S-G濾波，計(jì)算卷積因子的函數(shù)為SAVGOL，采用試錯(cuò)法最終設(shè)置參數(shù)為(2,2,0,2)，分別表示濾波核中心點(diǎn)左側(cè)的數(shù)據(jù)點(diǎn)個(gè)數(shù)、濾波核中心點(diǎn)右側(cè)的數(shù)據(jù)點(diǎn)個(gè)數(shù)、導(dǎo)數(shù)階數(shù)和平滑多項(xiàng)式的次數(shù)。

1.3.3 差分法

本文采用差分法[36-37]計(jì)算生育期內(nèi) EVI的峰值頻數(shù)（peak number，PN）。該方法是一種離散點(diǎn)求取極大值的方法，曾被應(yīng)用于復(fù)種指數(shù)的提取[15,36]。差分法對(duì)小峰極為敏感，因此，計(jì)算峰值頻數(shù)之前需要對(duì) EVI數(shù)據(jù)進(jìn)行濾波以避免數(shù)據(jù)噪聲干擾。

1.3.4 光譜突變法

處于不同生長(zhǎng)階段的冬小麥會(huì)呈現(xiàn)不同的光譜特征，也會(huì)有不同的光譜植被指數(shù)形態(tài)。由文獻(xiàn)[17]中的圖1和文獻(xiàn)[18]中的圖 5-1可以發(fā)現(xiàn)，冬小麥成熟收獲后NDVI比收獲前NDVI呈現(xiàn)明顯下降的特征。而這段時(shí)間內(nèi)其他植被仍處于生長(zhǎng)期內(nèi)，NDVI不會(huì)明顯下降。故該特征使冬小麥明顯區(qū)別于其他作物[17]。冬小麥的EVI時(shí)間序列曲線也比較明顯地反映了其成熟期前后植被指數(shù)突降的特征[16,38-40]。本文采用EVI探測(cè)光譜突變的像元。利用冬小麥成熟期的前一期和后一期影像計(jì)算光譜突變斜率Slope（式（2）），判定Slope小于負(fù)0.02的像元為冬小麥[17-18]。

式中EVIm-16和EVIm+16分別為冬小麥成熟期前、后16 d的EVI值。

2 結(jié)果分析

2.1 生育期估算及其驗(yàn)證

采用98個(gè)農(nóng)業(yè)氣象站點(diǎn)的氣象要素，利用多元線性逐步回歸方法，模擬各像元的播種期和成熟期DOY。播種期和成熟期DOY模擬結(jié)果分別為

式中Ys和Ym分別為擬合得到的播種期和成熟期DOY，X1為生長(zhǎng)期內(nèi)降水（mm），X2為生長(zhǎng)期內(nèi)氣溫之和（℃），X3為生長(zhǎng)期內(nèi)總輻射（0.001 W/m2），N為樣本數(shù)（有效站點(diǎn)數(shù)）。

在研究區(qū)內(nèi)模擬冬小麥播種期和成熟期的結(jié)果見(jiàn)圖3。冬小麥表現(xiàn)出由北向南逐漸進(jìn)入播種期（圖 3a）；而在第二年冬小麥則由南向北逐漸進(jìn)入成熟期（圖3b）。同時(shí)可以看出，冬小麥由南向北或由北向南逐漸進(jìn)入生育期的規(guī)律表現(xiàn)出一定的緯度地帶性，但并不嚴(yán)格與緯度走向一致，也不是與緯度之間有嚴(yán)格的角度關(guān)系。本文得到的冬小麥播種期和成熟期表現(xiàn)出來(lái)的緯向變化規(guī)律，與前人使用反距離權(quán)重法 IDW（inverse distance weighting）對(duì)站點(diǎn)生育期數(shù)據(jù)插值得到的冬小麥播種期空間分布[41]和多元逐步回歸方法獲取的冬小麥播種期和成熟期空間分布[31]基本是一致的。這也證明了本文使用的播種期和成熟期擬合方法及所得到的結(jié)果是可靠的。

分別提取剩余30%的站點(diǎn)（42個(gè)站）的播種期和成熟期模擬值，利用觀測(cè)值對(duì)其進(jìn)行驗(yàn)證（圖3c和3d）。驗(yàn)證結(jié)果顯示，播種期和成熟期的R2分別為0.69和0.67；RMSE分別為6.12和4.88 d，均在一旬之內(nèi)；且播種期和成熟期的驗(yàn)證結(jié)果都通過(guò)了顯著性檢驗(yàn)（P<0.001）。由此可知，用氣象要素模擬得到的生育期是可靠的。

圖3 冬小麥2009年播種期和2010年成熟期模擬值及其驗(yàn)證Fig.3 Simulated sowing date in 2009 and maturity date in 2010 for winter wheat and their validations

2.2 S-G濾波及峰值頻數(shù)和Slope計(jì)算結(jié)果

圖4 a為河北省藁城農(nóng)氣站所在像元EVI經(jīng)S-G濾波前后的曲線。可以看出，濾波前的數(shù)據(jù)有一些沒(méi)有明顯規(guī)律的小峰波動(dòng)，而經(jīng) S-G濾波后的數(shù)據(jù)中處于越冬期的小波峰被平滑掉了，比較接近真實(shí)情況。S-G濾波后的EVI數(shù)據(jù)雖然小幅度地加強(qiáng)了波谷值，但是并不影響峰值頻數(shù)的統(tǒng)計(jì)。

圖4b和圖4c分別為EVI峰值頻數(shù)等于2和Slope小于負(fù)0.02的像元空間分布。從圖4b中可以看出，PN等于 2的像元分布很多，這是由于差分法對(duì)小峰極為敏感，即便計(jì)算PN前對(duì)EVI曲線進(jìn)行S-G濾波，但一些小峰仍然能夠被探測(cè)到。圖4c顯示Slope小于負(fù)0.02的像元雖然也比較多，其主要分布在河北省中南部、河南省、山東省中西部、安徽省北部及江蘇省中北部，基本上為冬小麥的主要分布區(qū)?？梢跃C合利用這兩個(gè)特征提取冬小麥種植面積。

圖4 河北省藁城農(nóng)氣站濾波前后EVI曲線、峰值頻數(shù)為2和Slope小于-0.02的像元的空間分布Fig.4 EVI curves before and after S-G filter at Gaocheng station in Hebei province, spatial distribution of pixels with PN of 2 and that of pixels with slope value less than -0.02

2.3 冬小麥提取結(jié)果及驗(yàn)證

2.3.1 冬小麥種植面積空間分布

將圖4b、圖4c與研究區(qū)旱地分布取交集，即可得到冬小麥種植面積，提取結(jié)果如圖 5所示?？梢钥闯觯←溨饕植荚诤颖笔≈心喜?、河南省、山東省中西部、安徽省北部和江蘇省北部；而在山東省東部、河北省北部和西部、北京、天津、山西省、湖北省以及安徽和江蘇南部種植面積較小。與以往學(xué)者的提取結(jié)果[2-3,12,15-16]在空間分布上也較為一致。

圖5 2010年冬小麥面積提取結(jié)果Fig.5 Extracted winter wheat area in 2010

2.3.2 驗(yàn)證結(jié)果及參數(shù)影響分析

將不同濾波參數(shù)下冬小麥面積提取結(jié)果以省為統(tǒng)計(jì)單元進(jìn)行匯總，并與統(tǒng)計(jì)年鑒數(shù)據(jù)[24-30]對(duì)比（表 1）。河南省和山東省冬小麥種植面積最大，其次是河北省和安徽省。這與統(tǒng)計(jì)年鑒數(shù)據(jù)也是相符的。利用統(tǒng)計(jì)值和遙感提取值，計(jì)算后者相對(duì)于前者的高估和低估情況[3]：（后者-前者）/前者，精度[3]計(jì)算方法為：100-|高估或低估|。結(jié)果顯示，除天津、山西、湖北和江蘇之外，其余省份提取結(jié)果精度均達(dá)到80%以上。與姜亞珍等[3]結(jié)果對(duì)比，共有省份中黃淮海平原冬小麥主產(chǎn)區(qū)如河北、河南、山東三省，本文的提取精度均有明顯提高。

在研究區(qū)內(nèi)，將冬小麥提取面積分別以市和縣為統(tǒng)計(jì)單位進(jìn)行匯總，與統(tǒng)計(jì)年鑒數(shù)據(jù)[24-30]對(duì)比結(jié)果見(jiàn)圖6。在市級(jí)和縣級(jí)水平，冬小麥面積提取結(jié)果與統(tǒng)計(jì)數(shù)據(jù)的R2（RMSE）分別為 0.91（60.08×103hm2）和 0.80（8.97×103hm2）（P<0.001）。

在各省范圍內(nèi)，以市為統(tǒng)計(jì)單位匯總冬小麥提取面積，與統(tǒng)計(jì)年鑒數(shù)據(jù)[24-30]對(duì)比驗(yàn)證，驗(yàn)證結(jié)果的統(tǒng)計(jì)值見(jiàn)表 2?？梢钥闯?，盡管山西各市的提取驗(yàn)證結(jié)果 R2比較高，但遙感提取結(jié)果明顯低于統(tǒng)計(jì)值（a<1）；湖北各市的提取驗(yàn)證結(jié)果R2低于0.5，且遙感提取結(jié)果也明顯低于統(tǒng)計(jì)值；其余省份中，河南各市的遙感提取驗(yàn)證結(jié)果稍高于統(tǒng)計(jì)值（a>1），京津冀、山東和安徽各市的遙感提取值（a≈1），且R2均高于0.90。

表1 不同S-G濾波參數(shù)下各省冬小麥提取面積與統(tǒng)計(jì)數(shù)據(jù)對(duì)比Table 1 Comparison of extracted area and statistical data of winter wheat for each province with different filter parameters of S-G

圖6 2010年冬小麥面積遙感提取結(jié)果在市級(jí)水平和縣級(jí)水平與統(tǒng)計(jì)數(shù)據(jù)對(duì)比Fig.6 Comparison between extracted winter wheat areas in 2010 based on remote sensing and statistical winter wheat areas at city level and county level respectively

2.4 提取方法普適性分析

2.4.1 不同S-G濾波參數(shù)的提取結(jié)果

使用S-G濾波參數(shù)為(3,3,0,2)得到的冬小麥提取結(jié)果以省為單位匯總后與統(tǒng)計(jì)年鑒數(shù)據(jù)[24-30]對(duì)比（表 1）。相比于使用S-G濾波參數(shù)為(2,2,0,2)得到的結(jié)果，在估算精度較低的天津、山西、湖北和江蘇省中，雖然山西、湖北和江蘇的估算精度有不同程度的下降，但天津估算精度較半波窗口設(shè)置為2時(shí)有所提高，由63.18%提高到了90.90%。將半窗口設(shè)置為 3時(shí)得到的冬小麥提取結(jié)果與市級(jí)統(tǒng)計(jì)數(shù)據(jù)對(duì)比，R2為0.87，RMSE為77.07×103hm2。與 S-G 濾波參數(shù)為(2,2,0,2)得到的結(jié)果相比，R2下降，RMSE 上升了 16.99×103hm2。

表2 各省以市為單位驗(yàn)證結(jié)果統(tǒng)計(jì)值Table 2 Statistical values of validation at city level for each province

可見(jiàn) S-G濾波的半波窗口設(shè)置對(duì)提取結(jié)果有一定的影響，但影響不是很大。因此，本文使用的提取方法具有一定的普適性。

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

為了進(jìn)一步證明本文所構(gòu)建提取方法的普適性，使用該方法提取了2015年黃淮海平原3省2市（北京、天津、河北、河南、山東）冬小麥種植面積。因本文不側(cè)重冬小麥種植面積的時(shí)間變化，故不在此展示2015年提取結(jié)果的空間分布，僅對(duì)提取及驗(yàn)證結(jié)果做文字說(shuō)明。2015年冬小麥的空間分布特征與2010年基本保持一致。以市為統(tǒng)計(jì)單位（共計(jì)44個(gè)市，不包含承德市、張家口市、安陽(yáng)市和新鄉(xiāng)市）匯總各市冬小麥面積與統(tǒng)計(jì)年鑒數(shù)據(jù)對(duì)比驗(yàn)證，R2為0.84，RMSE為14.90×103hm2。可以看出，在不同年份間，該提取方法仍然適用，且可以達(dá)到較好的精度。由此，本文構(gòu)建的提取方法具有一定的普適性。

3 討 論

提取結(jié)果中顯示，南方的省份如湖北和江蘇，還有山區(qū)分布的山西省，提取精度較低（表 1）。湖北和江蘇北部，由于地處黃淮海大平原的邊緣，相比于平原中部地區(qū)耕地較為破碎。從圖 1中也可以看出，湖北和山西是因?yàn)榱值睾秃档鼗祀s分布，江蘇則是旱地和水田摻雜分布；而河北中南部、山東和河南的大部分地區(qū)，旱地多為大面積連續(xù)分布，耕地結(jié)構(gòu)較為簡(jiǎn)單，多為冬小麥-夏玉米輪作制度[15]。當(dāng)利用250 m分辨率遙感數(shù)據(jù)提取冬小麥面積時(shí)，使用本文的提取模型，大面積連續(xù)分布的冬小麥可以比較準(zhǔn)確地提取出來(lái)，而在耕地結(jié)構(gòu)復(fù)雜、地塊破碎的地區(qū)則提取精度有限。這是因?yàn)椋谙裨旌铣潭容^高的地區(qū)，冬小麥面積所占比例較小時(shí)極易被識(shí)別為其他地類，從而使冬小麥面積提取結(jié)果低于統(tǒng)計(jì)值，影響提取精度。

采用250 m分辨率遙感數(shù)據(jù)進(jìn)行大區(qū)域冬小麥種植面積提取是可行的，但在種植結(jié)構(gòu)復(fù)雜，耕地破碎的地區(qū)應(yīng)該優(yōu)先采用更高分辨率的遙感數(shù)據(jù)或使用其他方法。

4 結(jié) 論

本文利用氣象再分析資料（氣溫、降水和輻射）、農(nóng)業(yè)氣象站生育期觀測(cè)數(shù)據(jù)和MODIS EVI數(shù)據(jù)，構(gòu)建了一種既考慮區(qū)域間物候差異，又避免閾值設(shè)置的大范圍冬小麥快速提取方法。主要結(jié)論如下：

1）采用逐步進(jìn)入法構(gòu)建生育期-氣象要素的多元線性回歸模型，在像元尺度模擬冬小麥播種期和成熟期，并用站點(diǎn)數(shù)據(jù)進(jìn)行了驗(yàn)證。驗(yàn)證結(jié)果表明，播種期和成熟期的R2分別為0.69和0.67；RMSE分別為6.12和4.88 d，均在一旬之內(nèi)；且播種期和成熟期的驗(yàn)證結(jié)果都通過(guò)了顯著性檢驗(yàn)（P<0.001）；播種期和成熟期的模擬結(jié)果是可靠的。

2）利用統(tǒng)計(jì)年鑒數(shù)據(jù)對(duì)提取結(jié)果的驗(yàn)證表明，在市級(jí)和縣級(jí)尺度R2（RMSE）分別為0.91（60.08×103hm2）和 0.80（8.97×103hm2）。本文所使用的提取方法可以達(dá)到監(jiān)測(cè)區(qū)域冬小麥面積的業(yè)務(wù)需求精度。

3）山西、湖北和江蘇省的提取精度較低，使用250 m空間分辨率的影像導(dǎo)致面積小的冬小麥種植區(qū)不容易被提取出來(lái)。應(yīng)該發(fā)展算法或選擇其他遙感數(shù)據(jù)源，以提高在南方和山區(qū)作物面積的提取精度。

4）本文構(gòu)建的方法既考慮了大范圍內(nèi)冬小麥的物候差異，也避免了經(jīng)驗(yàn)閾值的設(shè)置，更具有普適性，可用于在大范圍內(nèi)快速提取冬小麥面積。

[1] 吳炳方. 中國(guó)農(nóng)情遙感監(jiān)測(cè)研究[J]. 中國(guó)科學(xué)院院刊，2004，19(3)：202－205，241.Wu Bingfang. China crop watch system with remote sensing[J]. Bulletin of Chinese Academy of Sciences, 2004,19(3): 202－205, 241. (in Chinese with English abstract)

[2] 張佳華，胡小夏，劉學(xué)鋒，等. 基于 MODIS數(shù)據(jù)提取華北典型區(qū)冬小麥種植面積[J]. 中國(guó)科學(xué)院大學(xué)學(xué)報(bào), 2013,30(5)：637－643.Zhang Jiahua, Hu Xiaoxia, Liu Xuefeng, et al. Extraction of winter wheat planting areas based on the time-series MODIS data in typical region of North China[J]. Journal of University of Chinese Academy of Sciences, 2013, 30(5):637－643. (in Chinese with English abstract)

[3] 姜亞珍，張瑜潔，孫琛，等. 基于 MODIS-EVI黃淮海平原冬小麥種植面積分帶提取[J]. 資源科學(xué)，2015，37(2)：417－424.Jiang Yazhen, Zhang Yujie, Sun Chen, 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 with English abstract)

[4] Xu X M, Conrad C, Doktor D. Optimising phonological metrics extraction for different crop types in Germany using the moderate resolution imaging spectrometer (MODIS)[J].Remote Sensing, 2017, 9(3): 254.

[5] Atzberger C, Rembold F. Mapping the spatial distribution of winter crops at sub-pixel level using AVHRR NDVI time series and neural nets[J]. Remote Sensing, 2013, 5(3): 1335－1354.

[6] 吳炳方. 全國(guó)農(nóng)情監(jiān)測(cè)與估產(chǎn)的運(yùn)行化遙感方法[J]. 地理學(xué)報(bào)，2000，55(1)：25－35.Wu Bingfang. Operational remote sensing methods for agriculture statistics[J]. Acta Gegraphica Sinica, 2000, 55(1):25－35. (in Chinese with English abstract)

[7] 葛廣秀，李衛(wèi)國(guó)，景元書. 基于NDVI密度分割的冬小麥種植面積提取[J]. 麥類作物學(xué)報(bào)，2014，34(7)：997－1002.Ge Guangxiu, Li Weiguo, Jing Yuanshu. Area of winter wheat extracted on NDVI density slicing[J]. Journal of Triticeae Crops, 2014, 34(7): 997－1002. (in Chinese with English abstract)

[8] 黃青，吳文斌，鄧輝，等. 2009年江蘇省冬小麥和水稻種植面積信息遙感提取及長(zhǎng)勢(shì)監(jiān)測(cè)[J]. 江蘇農(nóng)業(yè)科學(xué)，2010(6)：508－511.

[9] 楊小喚， 張香平，江東. 基于MODIS時(shí)序NDVI特征值提取多作物播種面積的方法[J]. 資源科學(xué)，2004，24(6)：17－22.Yang Xiaohuan, Zhang Xiangping, Jiang Dong. Extraction of multi-crop planting areas from MODIS data[J]. Resources Science, 2004, 26(6): 17－22. (in Chinese with English abstract)

[10] Qiao Hongbo, Zhang Hui, Cheng Dengfa. Application of EOS/MODIS-NDVI of different time sequences on monitoring winter wheat acreage in Henan Province[J]. Agricultural Science & Technology, 2008, 9(3): 124－126, 132.

[11] 陳健，劉云慧，宇振榮. 基于時(shí)序 MODIS-EVI數(shù)據(jù)的冬小麥種植信息提取[J]. 中國(guó)農(nóng)學(xué)通報(bào)，2011，27(1)：446－450.Chen Jian, Liu Yunhui, Yu Zhenrong. Planting information extraction of winter wheat based on the time-series MODIS-EVI[J]. Chinese Agricultural Science Bulletin, 2011,27(1): 446－450. (in Chinese with English abstract)

[12] 姜亞珍. MODIS監(jiān)測(cè)黃淮海平原冬小麥長(zhǎng)勢(shì)與土壤濕度[D]. 北京：中國(guó)農(nóng)業(yè)科學(xué)院，2015.Jiang Yazhen. MODIS Based Monitoring of Winter Wheat Growth and Soil Moisture on HHH Plain–supporting for Effect Assessment of Dry-hot Wind[D]. Beijing: Chinese Academy of Agricultural Sciences, 2015. (in Chinese with English abstract)

[13] 張明偉，周清波，陳仲新，等. 基于物候模型的作物種植面積變化監(jiān)測(cè)方法[J]. 農(nóng)業(yè)工程學(xué)報(bào)，2006，22(10)：139－144.Zhang Mingwei, Zhou Qingbo, Chen Zhongxin, et al. Crop acreage change detection based on phenology model[J].Transaction of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2006, 22(10): 139－144. (in Chinese with English abstract)

[14] 黃青，李丹丹，陳仲新，等. 基于 MODIS數(shù)據(jù)的冬小麥種植面積快速提取與長(zhǎng)勢(shì)監(jiān)測(cè)[J]. 農(nóng)業(yè)機(jī)械學(xué)報(bào)，2012，43(7)：163－167.Huang Qing, Li Dandan, Chen Zhongxin, et al. Monitoring of planting area and growth condition of winter wheat in China based on MODIS data[J]. Transactions of the Chinese Society for Agricultural Machinery, 2012, 43(7): 163－167.(in Chinese with English abstract)

[15] 潘學(xué)鵬，李改欣，劉峰貴，等. 華北平原冬小麥面積遙感提取及時(shí)空變化研究[J]. 中國(guó)生態(tài)農(nóng)業(yè)學(xué)報(bào)，2015，23(4)：497－505.Pan Xuepeng, Li Gaixin, Liu Fenggui, et al. Using remote sensing to determine spatio-temporal variations in winter wheat growing area in the North China Plain[J]. Chinese Journal of Eco-Agriculture, 2015, 23(4): 497－505. (in Chinese with English abstract)

[16] 王學(xué)，李秀彬，談明洪，等. 華北平原2001-2011年冬小麥播種面積變化遙感監(jiān)測(cè)[J]. 農(nóng)業(yè)工程學(xué)報(bào)，2015，31(8)：190－199.Wang Xue, Li Xiubin, Tan Minghong, et al. Remote sensing monitoring of changes in winter wheat area in North China Plain from 2001 to 2011[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2015, 31(8): 190－199. (in Chinese with English abstract)

[17] 王云峰，沈潤(rùn)平，田國(guó)珍. 基于 MODIS數(shù)據(jù)光譜突變法提取冬小麥種植面積研究[J]. 內(nèi)蒙古氣象，2009(6)：18－21.Wang Yunfeng, Shen Runping, Tian Guozhen. Study on extracting winter wheat planting area based on MODIS data by spectral mutation[J]. Meteorology Journal of Inner Mongolia, 2009(6): 18－21. (in Chinese with English abstract)

[18] 孫振蓉. 京津冀地區(qū)冬小麥面積估算和旱情遙感監(jiān)測(cè)研究[D]. 北京：北京林業(yè)大學(xué), 2015.Sun Zhenrong. Estimation of Winter Wheat Area and Drought Monitoring Research Based on Remote Sensing Data in Beijing-Tianjin-Hebei Region[D]. Beijing: Beijing Forestry University, 2015. (in Chinese with English abstract)

[19] 趙靜旸，宋淑麗，朱文耀. ERA-Interim應(yīng)用于中國(guó)地區(qū)地基GPS/PWV計(jì)算的精度評(píng)估[J]. 武漢大學(xué)學(xué)報(bào)(信息科學(xué)版)，2014，39(8)：935－939，1008.Zhao Jingyang, Song Shuli, Zhu Wenyao. Accuracy assessment of applying ERA-Interim reanalysis data to calculate ground-based GPS/PWV over China[J]. Geomatics and Information Science of Wuhan University, 2014, 39(8):935－939, 1008. (in Chinese with English abstract)

[20] Dee D P, Uppala S M, Simmons A J, et al. The ERA-Interim reanalysis: Configuration and performance of the data assimilation system[J]. Quarterly Journal of the Royal Meteorological Society, 2011, 137(656): 553－597.

[21] 李川，張廷軍，陳靜. 近40年青藏高原地區(qū)的氣候變化——NCEP和ECMWF地面氣溫及降水再分析和實(shí)測(cè)資料對(duì)比分析[J]. 高原氣象，2004，23(增刊1)：97－103.Li Chuan, Zhang Tingjun, Chen Jing. Climatic change of Qinghai-Xizang Plateau Region in recent 40 year reanalysis and surface observation data: contrast of observational data and NCEP, ECMWF surface air temperature and precipitation[J]. Plateau Meteorology, 2004, 23(Supp.1): 97－103. (in Chinese with English abstract)

[22] 韋芬芬，湯劍平，王淑瑜. 中國(guó)區(qū)域夏季再分析資料高空變量可信度的檢驗(yàn)[J]. 地球物理學(xué)報(bào)，2015，58(2)：383－397.Wei Fenfen, Tang Jianping, Wang Shuyu. A reliability assessment of upper-level reanalysis datasets over China[J].Chinese Journal of Geophysics, 2015, 58(2): 383－397. (in Chinese with English abstract)

[23] 張婷，高慶九，孫玉婷. 不同時(shí)間尺度下中國(guó)東部夏季最高溫度觀測(cè)與再分析資料的比較[J]. 高原氣象，2017，36(2)：468－479.Zhang Ting, Gao Qingjiu, Sun Yunting. Comparison of reanalysis data and observation about summer maximum temperature on different time scales in Eastern China[J].Plateau Meteorology, 2017, 36(2): 468－479. (in Chinese with English abstract)

[24] 中華人民共和國(guó)國(guó)家統(tǒng)計(jì)局. 中國(guó)統(tǒng)計(jì)年鑒[M]. 北京：中國(guó)統(tǒng)計(jì)出版社，2011.

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

[26] 河南省統(tǒng)計(jì)局，國(guó)家統(tǒng)計(jì)局河南調(diào)查總隊(duì). 河南統(tǒng)計(jì)年鑒[M]. 北京：中國(guó)統(tǒng)計(jì)出版社，2011.

[27] 山東省統(tǒng)計(jì)局. 山東統(tǒng)計(jì)年鑒[M]. 北京：中國(guó)統(tǒng)計(jì)出版社，2011.

[28] 山西省統(tǒng)計(jì)局，國(guó)家統(tǒng)計(jì)局山西調(diào)查總隊(duì). 山西統(tǒng)計(jì)年鑒[M]. 北京：中國(guó)統(tǒng)計(jì)出版社，2011.

[29] 湖北省農(nóng)村統(tǒng)計(jì)年鑒編輯委員會(huì). 湖北農(nóng)村統(tǒng)計(jì)年鑒[M].北京：中國(guó)統(tǒng)計(jì)出版社，2011.

[30] 安徽省統(tǒng)計(jì)局. 安徽統(tǒng)計(jì)年鑒[M]. 北京：中國(guó)統(tǒng)計(jì)出版社，2011.

[31] 楊建瑩，梅旭榮，劉勤，等. 氣候變化背景下華北地區(qū)冬小麥生育期的變化特征[J]. 植物生態(tài)學(xué)報(bào)，2011，35(6)：623－631.Yang Jianying, Mei Xurong, Liu Qin, et al. Variations of winter wheat growth stages under climate changes in northern China[J]. Chinese Journal of Plant Ecology, 2011,35(6): 623－631. (in Chinese with English abstract)

[32] 劉勁松，陳輝，楊彬云，等. 河北省年均降水量插值方法比較[J]. 生態(tài)學(xué)報(bào)，2009，29(7)：3493－3500.Liu Jinsong, Chen Hui, Yang Binyun, et al. Comparison of interpolation method on annual mean precipitation in Hebei Province[J]. Acta Ecologica Sinica, 2009, 29(7): 3493－3500.(in Chinese with English abstract)

[33] 黃耀歡，王建華，江東，等. 利用S-G濾波進(jìn)行MODIS-EVI時(shí)間序列數(shù)據(jù)重構(gòu)[J]. 武漢大學(xué)學(xué)報(bào)(信息科學(xué)版)，2009，34(12)：1440－1443，1513.Huang Yaohuan, Wang Jianhua, Jiang Dong, et al.Reconstruction of MODIS-EVI time-series data with S-G filter[J]. Geomatics and Information Science of Wuhan University, 2009, 34(12): 1440－1443, 1513. (in Chinese with English abstract)

[34] 邊金虎，李愛(ài)農(nóng)，宋孟強(qiáng)，等. MODIS植被指數(shù)時(shí)間序列Savitzky-Golay濾波算法重構(gòu)[J]. 遙感學(xué)報(bào)，2010，14(4)：725－741.Bian Jinhu, Li Ainong, Song Mengqiang, et al.Reconstruction of NDVI time-series datasets of MODIS based on Savitaky-Golay filter[J]. Journal of Remote Sensing,2010, 14(4): 725－741. (in Chinese with English abstract)

[35] 王連喜，徐勝男，李琪，等. 基于決策樹(shù)和混合像元分解的江蘇省冬小麥種植面積提取[J]. 農(nóng)業(yè)工程學(xué)報(bào)，2016，32(5)：182－187.Wang Lianxi, Xu Shengnan, Li Qi, et al. Extraction of winter wheat planted area in Jiangsu province using decision tree and mixed-pixel methods[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CASE), 2016, 32(5): 182－187. (in Chinese with English abstract)

[36] 范錦龍，吳炳方. 復(fù)種指數(shù)遙感監(jiān)測(cè)方法[J]. 遙感學(xué)報(bào)，2004，8(6)：628－636.Fan Jinlong, Wu Bingfang. A methodology for retrieving cropping index from NDVI profile[J]. Journal of Remote Sensing, 2004, 8(6): 628－636. (in Chinese with English abstract)

[37] 范錦龍. 復(fù)種指數(shù)遙感監(jiān)測(cè)方法研究[D]. 北京：中國(guó)科學(xué)院研究生院(遙感應(yīng)用研究所)，2003.Fan Jinlong. Study on Remote Sensing Methods for Monitoring Multiple Cropping Index[D]. Beijing: Institute of Remote Sensing Applications, Chinese Academy of Sciences,2003. (in Chinese with English abstract)

[38] 蘆杰，婁玉欽，姚麗. 基于MODIS EVI時(shí)序數(shù)據(jù)的冬小麥種植面積提取[J]. 湖北農(nóng)業(yè)科學(xué)，2014，53(19)：4727－4730.Lu Jie, Lou Yuqin, Yao Li. Extraction of winter wheat planting area based on data of MODIS EVI time-series[J].Hubei Agricultural Sciences, 2014, 53(19): 4727－4730. (in Chinese with English abstract)

[39] 潘耀忠，李樂(lè)，張錦水，等. 基于典型物候特征的MODIS-EVI時(shí)間序列數(shù)據(jù)農(nóng)作物種植面積提取方法——小區(qū)域冬小麥實(shí)驗(yàn)研究[J]. 遙感學(xué)報(bào)，2011，15(3)：578－594.Pan Yaozhong, Li Le, Zhang Jinshui, et al. Crop area estimation based on MODIS-EVI time series according to distinct characteristics of key phenology phases: A case study of winter wheat area estimating in small-scale area[J]. Journal of Remote Sensing, 2011, 15(3): 578－594. (in Chinese with English abstract)

[40] 張霞，帥通，楊杭，等. 基于MODIS EVI圖像時(shí)間序列的冬小麥面積提取[J]. 農(nóng)業(yè)工程學(xué)報(bào)，2010，26(增刊1)：220－224.Zhang Xia, Shuai Tong, Yang Hang, et al. Winter wheat planting area extraction based on MODIS EVI image time series[J]. Transaction of the Chinese Society of Agricultural Engineering (Transactions of the CASE), 2010, 26(Supp.1):220－224. (in Chinese with English abstract)

[41] 王培娟，梁宏，謝東輝，等. 氣候變暖對(duì)華北冬小麥返青前熱量條件的影響[J]. 麥類作物學(xué)報(bào)，2014，34(1)：54－63.Wang Peijuan, Liang Hong, Xie Donghui, et al. Influence of climate warming on heat conditions before green-up stages for winter wheat in North China Plain[J]. Journal of Triticeae Crops, 2014, 34(1): 54－63. (in Chinese with English abstract)