基于MODIS-EVI及物候差异免阈值提取黄淮海平原冬小麦面积

张 莎，张佳华※，白 雲 ，姚凤梅

（1. 中国科学院遥感与数字地球研究所数字地球重点实验室，北京 100094；2. 中国科学院大学地球与行星科学学院，北京 100049）

0 引 言

冬小麦是中国3大粮食作物之一[1-2]。获取准确冬小麦种植面积对于粮食估产和保障国家粮食安全有重要意义[1]。黄淮海平原是中国最大的平原区，是重要的冬小麦生产基地[3]。因此，准确估计黄淮海地区冬小麦种植面积有重要的意义。

遥感作为一种地球观测技术，在大范围面积和资源监测上有着比传统手段更为明显的优势[2]，已被广泛地应用于作物面积提取[4-5]。前人在黄淮海区域冬小麦面积的遥感提取方面已开展了大量工作。研究表明，当利用中低分辨率影像提取中国北方冬小麦面积时，对绿度指数进行图像分割就能取得较好的效果[6]。葛广秀等[7]使用归一化植被指数 NDVI（normalized difference vegetation index）密度分割方法提取了江苏省沭阳县的冬小麦面积，精度达92.37%。根据不同省市冬小麦的物候特征，使用关键生育期时相的 NDVI影像并设置阈值构建各省或各市的冬小麦提取模型，都可以达到较好的精度。黄青等[8]根据江苏省冬小麦物候特征，基于MODIS NDVI数据构建了适合于江苏省的冬小麦提取模型，总体提取精度达78%以上；杨小唤等[9]根据冬小麦的生育期特点，使用关键时相的 NDVI影像并设置阈值提取了北京市包含冬小麦在内的多种作物种植面积，总体精度达95%以上。Qiao等[10]也采用 NDVI阈值法较为准确地提取了河南省冬小麦种植面积。

可以看出，当研究区范围小于或等于省级尺度时，研究者即便使用同一套阈值参数，即忽略冬小麦在区域间的物候差异，也可以较为准确地提取冬小麦面积；而当研究范围纬度跨越较大时，如同一省份南北方向延伸较大或包含多个省份时，物候差异则不能忽略，研究者往往会选择分省或分市进行阈值设置，以区分冬小麦在不同地区的的物候差异。陈健等[11]基于时序中分辨率成像光谱仪-增强型植被指数（MODIS-EVI，moderate resolution imaging spectroradiometer-enhanced vegetation index）数据提取河北省（36°01′N—42°37′N）冬小麦种植面积时，就采用了分市设置阈值的方法。张佳华等[2]根据黄淮海各省冬小麦物候特征，采用EVI阈值模型，分别对研究区内5个省市的冬小麦生育期共10个EVI阈值进行设置，最终得到黄淮海平原典型区的冬小麦种植面积，在省级尺度验证整体精度可达 98%。当研究区包含较多省份时，这种方法无疑增加了工作量。还有学者根据纬度进行水平分带[3]和倾斜分带[12]以区分冬小麦在区域间的物候差异，根据各纬度带冬小麦生育关键期不同而使用不同时相的影像，各纬度带的提取模型使用相同的阈值。

采用植被指数阈值法，无论采用分区方案，还是根据纬度进行水平分带或倾斜分带方案，都无法避免对植被指数阈值进行确定。研究者通常依据作物在各区的物候特征根据经验值不断调整来确定阈值，分区越小提取精度相应越高[8,13-14]。再者，不同区域之间，或同一区域不同年份之间，作物长势有所差异则阈值也要相应调整。阈值的不确定性增加了冬小麦种植面积提取的不确定性。因此，找到一个避免阈值设置的提取特征十分有必要。潘学鹏等[15]利用MODIS数据，以基于复种指数的决策树分类方法提取了黄淮海平原主体地区的冬小麦面积，京津冀地区长时间序列的市级尺度验证结果 R2为0.885。王学等[16]利用MODIS数据和TM数据，根据小麦生育期内波峰和波谷的特征构建提取模型，提取了黄淮海平原主体部分冬小麦种植面积。王云峰等[17]使用MODIS NDVI数据，根据冬小麦区别于同期作物的物候特征，即在收获期植被指数突降（此阶段其他作物均为生长季），采用两景影像提取山西省运城市冬小麦面积，并取得了较为合理的结果。孙振蓉[18]也利用了类似的方法提取了京津冀地区的冬小麦种植面积。为了减小阈值设置给提取结果带来的不确定性，同时在大区域范围内既能考虑物候差异又能快速提取冬小麦面积，本文将引入波形特征。

鉴于此，本研究以黄淮海平原为研究区，拟提出一种既能考虑区域间物候差异又具有普适性的提取方法。首先，采用基本气候要素（气温、降水和辐射）模拟冬小麦的播种期和成熟期；然后，使用滤波后的MODIS EVI数据，利用差分法逐像元提取播种期至成熟期内的 EVI峰值频数；逐像元计算成熟期前后 EVI的下降速率；最后，利用统一的提取规则提取研究区内冬小麦种植面积并用统计数据验证。旨在为大区域范围内快速、准确地提取冬小麦种植面积提出一套既能满足作物种植面积监测精度需求又具有普适性、能在大范围区域使用的提取模型，服务于国家冬小麦面积监测的业务化运行。

1 材料与方法

1.1 研究区概况

本文研究区为黄淮海平原，位于中国东部地区（32°~42°N，113°~120°E），覆盖北京市、天津市、河北省、河南省、山东省、山西省、湖北省、安徽省和江苏省共计 9个省市。黄淮海平原是中国最大的平原区，也是中国冬小麦的主产区。该区范围大，在区域范围内辐射和气温等气候要素有明显的区域变化，在植被上变现出明显的物候差异。

1.2 数据来源与处理

1.2.1 MODIS EVI数据

所采用的MODIS EVI 数据为MOD13Q1的C6数据（https://search.earth data.nasa.gov/），空间分辨率为250 m。本文获取从2009年第241天至2010年第209天每16 d合成一期的影像，共22期。这些影像涵盖了研究区（h26v04，h26v05，h27v04，h27v05，h27v06，h28v05，h28v06）冬小麦的全部生育期。MODIS EVI数据有效值为-3 000～10 000，比例因子为0.000 1。为去除非植被像元干扰，本文提取EVI大于0的像元以进行下一步研究。

1.2.2 土地利用数据和旱地数据

中国2010年土地利用数据来源于中国科学院资源环境科学数据中心(http://www.resdc.cn)，空间分辨率为1 km。将该数据重采样为250 m，并用矢量图裁剪，得到研究区土地利用（图1）。除耕地按二级分类体系显示外，其余土地利用类型均将编码合并处理按一级分类体系显示。本文提取了土地利用类型中的旱地（编码为121，122，123和124，分别为山地旱地、丘陵旱地、平原旱地和大于 20°坡度旱地）。利用旱地提取结果对所提取冬小麦面积进行掩膜处理，以排除其他土地利用类型对冬小麦提取结果的影响。

1.2.3 农业气象站点数据

研究区内农业气象站点数据来源于中国气象数据网（http://data.cma.cn/），使用了农气资料中的AGME_AB2_CHN_TEN数据集（中国农作物生长发育和农田土壤湿度旬值数据集）。本文使用了140个冬小麦站点，其空间分布见图1。所用农业气象站点数据包括各站点的经纬度、作物类型、生育期名称和生育期日期数据。所使用的生育期仅包括冬小麦播种期和成熟期，并将其日期转换为日序 DOY（day of year）。

图1 研究区土地利用及所用农业气象站点分布Fig.1 Land cover of study area and distribution of agro-meteorological sites in study area

1.2.4 气象再分析资料数据

本文采用空间分辨率为0.125°的ERA-Interim日值数据集获取气象要素，包括气温、降水和辐射。ERA-Interim是欧洲中期数值预报中心提供的全球最新的大气再分析产品，较之前的产品精度更高[19-20]，已被广泛应用于气候变化等研究[21-23]。ERA-Interim数据采用nc格式存储，本研究利用IDL代码将其读取为img格式，并根据农业气象站点的位置提取其所在栅格的值作为站点的气象要素值。

1.2.5 统计数据

冬小麦播种面积统计数据来源于各省市统计年鉴[24-30]。本文收集的市级统计数据包括北京市、天津市、河北省9个市（承德市和张家口市除外）、河南省 18个市、山东省17个市、山西省10个市（大同市除外）、湖北省16个市（襄阳市除外）和安徽省17个市共89个地级市的2010年和2015年冬小麦种植面积数据；县级统计数据包括河北省113个县、河南省10个县、山西省77个县和湖北省55个县共255个县的2010年冬小麦种植面积数据。

1.3 研究方法

提取冬小麦的流程如图 2所示。首先，利用随机分布的98个农业气象站的观测播种期与成熟期数据和基本气候要素（气温、降水和辐射），采用逐步回归法，构建多元线性回归模型模拟冬小麦播种期和成熟期，并利用其他42个站点验证模拟结果；其次，对时序EVI数据进行 Savitzky-Golay（S-G）滤波；然后，根据模拟得到的播种期和成熟期，截取每个像元处于播种期至成熟期之间的EVI序列，计算截取后得到的EVI序列的峰值频数；同时，逐像元获取成熟期前后的EVI，计算每个像元成熟期前后EVI发生突变的情况（Slope）；最后，提取峰值频数为2[15]和Slope小于-0.02[17-18]的像元，并与黄淮海平原旱地数据取交集，得到冬小麦的提取结果，同时利用统计数据对提取结果进行验证。

图2 冬小麦面积提取流程图Fig.2 Flowchart of winter wheat area extraction

1.3.1 多元线性逐步回归法

区域间不同的水热组合条件是造成地区间物候差异的根本原因。本文采用多元线性逐步回归法[31]拟合冬小麦播种期和成熟期日序DOY与气温、降水和辐射之间的关系。研究区内 140个冬小麦农业气象站点，随机选择其中98个站点用于构建播种期和成熟期DOY与气候要素的关系，并利用剩余42个站点（保证验证站点数占总站点数的30%[32]）的生育期数据对拟合关系进行验证。

1.3.2 S-G滤波

S-G滤波拟合算法是一种基于平滑时间序列数据和最小二乘原理的卷积算法[33-35]

式中为拟合EVI值；Yj+i为时间序列中的第(j+i)个EVI值； Ci为第i个EVI值滤波时的卷积系数；msg为滤波窗口大小的一半；Nsg为滤波器长度，等于滑动数组的宽度(2msg+1)。本文利用IDL（interactive data language）实现时间序列EVI影像的S-G滤波，计算卷积因子的函数为SAVGOL，采用试错法最终设置参数为(2,2,0,2)，分别表示滤波核中心点左侧的数据点个数、滤波核中心点右侧的数据点个数、导数阶数和平滑多项式的次数。

1.3.3 差分法

本文采用差分法[36-37]计算生育期内 EVI的峰值频数（peak number，PN）。该方法是一种离散点求取极大值的方法，曾被应用于复种指数的提取[15,36]。差分法对小峰极为敏感，因此，计算峰值频数之前需要对 EVI数据进行滤波以避免数据噪声干扰。

1.3.4 光谱突变法

处于不同生长阶段的冬小麦会呈现不同的光谱特征，也会有不同的光谱植被指数形态。由文献[17]中的图1和文献[18]中的图 5-1可以发现，冬小麦成熟收获后NDVI比收获前NDVI呈现明显下降的特征。而这段时间内其他植被仍处于生长期内，NDVI不会明显下降。故该特征使冬小麦明显区别于其他作物[17]。冬小麦的EVI时间序列曲线也比较明显地反映了其成熟期前后植被指数突降的特征[16,38-40]。本文采用EVI探测光谱突变的像元。利用冬小麦成熟期的前一期和后一期影像计算光谱突变斜率Slope（式（2）），判定Slope小于负0.02的像元为冬小麦[17-18]。

式中EVIm-16和EVIm+16分别为冬小麦成熟期前、后16 d的EVI值。

2 结果分析

2.1 生育期估算及其验证

采用98个农业气象站点的气象要素，利用多元线性逐步回归方法，模拟各像元的播种期和成熟期DOY。播种期和成熟期DOY模拟结果分别为

式中Ys和Ym分别为拟合得到的播种期和成熟期DOY，X1为生长期内降水（mm），X2为生长期内气温之和（℃），X3为生长期内总辐射（0.001 W/m2），N为样本数（有效站点数）。

在研究区内模拟冬小麦播种期和成熟期的结果见图3。冬小麦表现出由北向南逐渐进入播种期（图 3a）；而在第二年冬小麦则由南向北逐渐进入成熟期（图3b）。同时可以看出，冬小麦由南向北或由北向南逐渐进入生育期的规律表现出一定的纬度地带性，但并不严格与纬度走向一致，也不是与纬度之间有严格的角度关系。本文得到的冬小麦播种期和成熟期表现出来的纬向变化规律，与前人使用反距离权重法 IDW（inverse distance weighting）对站点生育期数据插值得到的冬小麦播种期空间分布[41]和多元逐步回归方法获取的冬小麦播种期和成熟期空间分布[31]基本是一致的。这也证明了本文使用的播种期和成熟期拟合方法及所得到的结果是可靠的。

分别提取剩余30%的站点（42个站）的播种期和成熟期模拟值，利用观测值对其进行验证（图3c和3d）。验证结果显示，播种期和成熟期的R2分别为0.69和0.67；RMSE分别为6.12和4.88 d，均在一旬之内；且播种期和成熟期的验证结果都通过了显著性检验（P<0.001）。由此可知，用气象要素模拟得到的生育期是可靠的。

图3 冬小麦2009年播种期和2010年成熟期模拟值及其验证Fig.3 Simulated sowing date in 2009 and maturity date in 2010 for winter wheat and their validations

2.2 S-G滤波及峰值频数和Slope计算结果

图4 a为河北省藁城农气站所在像元EVI经S-G滤波前后的曲线。可以看出，滤波前的数据有一些没有明显规律的小峰波动，而经 S-G滤波后的数据中处于越冬期的小波峰被平滑掉了，比较接近真实情况。S-G滤波后的EVI数据虽然小幅度地加强了波谷值，但是并不影响峰值频数的统计。

图4b和图4c分别为EVI峰值频数等于2和Slope小于负0.02的像元空间分布。从图4b中可以看出，PN等于 2的像元分布很多，这是由于差分法对小峰极为敏感，即便计算PN前对EVI曲线进行S-G滤波，但一些小峰仍然能够被探测到。图4c显示Slope小于负0.02的像元虽然也比较多，其主要分布在河北省中南部、河南省、山东省中西部、安徽省北部及江苏省中北部，基本上为冬小麦的主要分布区。可以综合利用这两个特征提取冬小麦种植面积。

图4 河北省藁城农气站滤波前后EVI曲线、峰值频数为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 冬小麦提取结果及验证

2.3.1 冬小麦种植面积空间分布

将图4b、图4c与研究区旱地分布取交集，即可得到冬小麦种植面积，提取结果如图 5所示。可以看出，冬小麦主要分布在河北省中南部、河南省、山东省中西部、安徽省北部和江苏省北部；而在山东省东部、河北省北部和西部、北京、天津、山西省、湖北省以及安徽和江苏南部种植面积较小。与以往学者的提取结果[2-3,12,15-16]在空间分布上也较为一致。

图5 2010年冬小麦面积提取结果Fig.5 Extracted winter wheat area in 2010

2.3.2 验证结果及参数影响分析

将不同滤波参数下冬小麦面积提取结果以省为统计单元进行汇总，并与统计年鉴数据[24-30]对比（表 1）。河南省和山东省冬小麦种植面积最大，其次是河北省和安徽省。这与统计年鉴数据也是相符的。利用统计值和遥感提取值，计算后者相对于前者的高估和低估情况[3]：（后者-前者）/前者，精度[3]计算方法为：100-|高估或低估|。结果显示，除天津、山西、湖北和江苏之外，其余省份提取结果精度均达到80%以上。与姜亚珍等[3]结果对比，共有省份中黄淮海平原冬小麦主产区如河北、河南、山东三省，本文的提取精度均有明显提高。

在研究区内，将冬小麦提取面积分别以市和县为统计单位进行汇总，与统计年鉴数据[24-30]对比结果见图6。在市级和县级水平，冬小麦面积提取结果与统计数据的R2（RMSE）分别为 0.91（60.08×103hm2）和 0.80（8.97×103hm2）（P<0.001）。

在各省范围内，以市为统计单位汇总冬小麦提取面积，与统计年鉴数据[24-30]对比验证，验证结果的统计值见表 2。可以看出，尽管山西各市的提取验证结果 R2比较高，但遥感提取结果明显低于统计值（a<1）；湖北各市的提取验证结果R2低于0.5，且遥感提取结果也明显低于统计值；其余省份中，河南各市的遥感提取验证结果稍高于统计值（a>1），京津冀、山东和安徽各市的遥感提取值（a≈1），且R2均高于0.90。

表1 不同S-G滤波参数下各省冬小麦提取面积与统计数据对比Table 1 Comparison of extracted area and statistical data of winter wheat for each province with different filter parameters of S-G

图6 2010年冬小麦面积遥感提取结果在市级水平和县级水平与统计数据对比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滤波参数的提取结果

使用S-G滤波参数为(3,3,0,2)得到的冬小麦提取结果以省为单位汇总后与统计年鉴数据[24-30]对比（表 1）。相比于使用S-G滤波参数为(2,2,0,2)得到的结果，在估算精度较低的天津、山西、湖北和江苏省中，虽然山西、湖北和江苏的估算精度有不同程度的下降，但天津估算精度较半波窗口设置为2时有所提高，由63.18%提高到了90.90%。将半窗口设置为 3时得到的冬小麦提取结果与市级统计数据对比，R2为0.87，RMSE为77.07×103hm2。与 S-G 滤波参数为(2,2,0,2)得到的结果相比，R2下降，RMSE 上升了 16.99×103hm2。

表2 各省以市为单位验证结果统计值Table 2 Statistical values of validation at city level for each province

可见 S-G滤波的半波窗口设置对提取结果有一定的影响，但影响不是很大。因此，本文使用的提取方法具有一定的普适性。

2.4.2 不同年份的提取结果

为了进一步证明本文所构建提取方法的普适性，使用该方法提取了2015年黄淮海平原3省2市（北京、天津、河北、河南、山东）冬小麦种植面积。因本文不侧重冬小麦种植面积的时间变化，故不在此展示2015年提取结果的空间分布，仅对提取及验证结果做文字说明。2015年冬小麦的空间分布特征与2010年基本保持一致。以市为统计单位（共计44个市，不包含承德市、张家口市、安阳市和新乡市）汇总各市冬小麦面积与统计年鉴数据对比验证，R2为0.84，RMSE为14.90×103hm2。可以看出，在不同年份间，该提取方法仍然适用，且可以达到较好的精度。由此，本文构建的提取方法具有一定的普适性。

3 讨 论

提取结果中显示，南方的省份如湖北和江苏，还有山区分布的山西省，提取精度较低（表 1）。湖北和江苏北部，由于地处黄淮海大平原的边缘，相比于平原中部地区耕地较为破碎。从图 1中也可以看出，湖北和山西是因为林地和旱地混杂分布，江苏则是旱地和水田掺杂分布；而河北中南部、山东和河南的大部分地区，旱地多为大面积连续分布，耕地结构较为简单，多为冬小麦-夏玉米轮作制度[15]。当利用250 m分辨率遥感数据提取冬小麦面积时，使用本文的提取模型，大面积连续分布的冬小麦可以比较准确地提取出来，而在耕地结构复杂、地块破碎的地区则提取精度有限。这是因为，在像元混合程度较高的地区，冬小麦面积所占比例较小时极易被识别为其他地类，从而使冬小麦面积提取结果低于统计值，影响提取精度。

采用250 m分辨率遥感数据进行大区域冬小麦种植面积提取是可行的，但在种植结构复杂，耕地破碎的地区应该优先采用更高分辨率的遥感数据或使用其他方法。

4 结 论

本文利用气象再分析资料（气温、降水和辐射）、农业气象站生育期观测数据和MODIS EVI数据，构建了一种既考虑区域间物候差异，又避免阈值设置的大范围冬小麦快速提取方法。主要结论如下：

1）采用逐步进入法构建生育期-气象要素的多元线性回归模型，在像元尺度模拟冬小麦播种期和成熟期，并用站点数据进行了验证。验证结果表明，播种期和成熟期的R2分别为0.69和0.67；RMSE分别为6.12和4.88 d，均在一旬之内；且播种期和成熟期的验证结果都通过了显著性检验（P<0.001）；播种期和成熟期的模拟结果是可靠的。

2）利用统计年鉴数据对提取结果的验证表明，在市级和县级尺度R2（RMSE）分别为0.91（60.08×103hm2）和 0.80（8.97×103hm2）。本文所使用的提取方法可以达到监测区域冬小麦面积的业务需求精度。

3）山西、湖北和江苏省的提取精度较低，使用250 m空间分辨率的影像导致面积小的冬小麦种植区不容易被提取出来。应该发展算法或选择其他遥感数据源，以提高在南方和山区作物面积的提取精度。

4）本文构建的方法既考虑了大范围内冬小麦的物候差异，也避免了经验阈值的设置，更具有普适性，可用于在大范围内快速提取冬小麦面积。

[1] 吴炳方. 中国农情遥感监测研究[J]. 中国科学院院刊，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] 张佳华，胡小夏，刘学锋，等. 基于 MODIS数据提取华北典型区冬小麦种植面积[J]. 中国科学院大学学报, 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]. 资源科学，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] 吴炳方. 全国农情监测与估产的运行化遥感方法[J]. 地理学报，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] 葛广秀，李卫国，景元书. 基于NDVI密度分割的冬小麦种植面积提取[J]. 麦类作物学报，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年江苏省冬小麦和水稻种植面积信息遥感提取及长势监测[J]. 江苏农业科学，2010(6)：508－511.

[9] 杨小唤， 张香平，江东. 基于MODIS时序NDVI特征值提取多作物播种面积的方法[J]. 资源科学，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] 陈健，刘云慧，宇振荣. 基于时序 MODIS-EVI数据的冬小麦种植信息提取[J]. 中国农学通报，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监测黄淮海平原冬小麦长势与土壤湿度[D]. 北京：中国农业科学院，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] 张明伟，周清波，陈仲新，等. 基于物候模型的作物种植面积变化监测方法[J]. 农业工程学报，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数据的冬小麦种植面积快速提取与长势监测[J]. 农业机械学报，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] 潘学鹏，李改欣，刘峰贵，等. 华北平原冬小麦面积遥感提取及时空变化研究[J]. 中国生态农业学报，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] 王学，李秀彬，谈明洪，等. 华北平原2001-2011年冬小麦播种面积变化遥感监测[J]. 农业工程学报，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] 王云峰，沈润平，田国珍. 基于 MODIS数据光谱突变法提取冬小麦种植面积研究[J]. 内蒙古气象，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] 孙振蓉. 京津冀地区冬小麦面积估算和旱情遥感监测研究[D]. 北京：北京林业大学, 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应用于中国地区地基GPS/PWV计算的精度评估[J]. 武汉大学学报(信息科学版)，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年青藏高原地区的气候变化——NCEP和ECMWF地面气温及降水再分析和实测资料对比分析[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] 韦芬芬，汤剑平，王淑瑜. 中国区域夏季再分析资料高空变量可信度的检验[J]. 地球物理学报，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] 张婷，高庆九，孙玉婷. 不同时间尺度下中国东部夏季最高温度观测与再分析资料的比较[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] 中华人民共和国国家统计局. 中国统计年鉴[M]. 北京：中国统计出版社，2011.

[25] 河北省统计局，河北省人民政府办公厅. 河北农村统计年鉴[M]. 北京：中国统计出版社，2011.

[26] 河南省统计局，国家统计局河南调查总队. 河南统计年鉴[M]. 北京：中国统计出版社，2011.

[27] 山东省统计局. 山东统计年鉴[M]. 北京：中国统计出版社，2011.

[28] 山西省统计局，国家统计局山西调查总队. 山西统计年鉴[M]. 北京：中国统计出版社，2011.

[29] 湖北省农村统计年鉴编辑委员会. 湖北农村统计年鉴[M].北京：中国统计出版社，2011.

[30] 安徽省统计局. 安徽统计年鉴[M]. 北京：中国统计出版社，2011.

[31] 杨建莹，梅旭荣，刘勤，等. 气候变化背景下华北地区冬小麦生育期的变化特征[J]. 植物生态学报，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]. 生态学报，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滤波进行MODIS-EVI时间序列数据重构[J]. 武汉大学学报(信息科学版)，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] 边金虎，李爱农，宋孟强，等. MODIS植被指数时间序列Savitzky-Golay滤波算法重构[J]. 遥感学报，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] 王连喜，徐胜男，李琪，等. 基于决策树和混合像元分解的江苏省冬小麦种植面积提取[J]. 农业工程学报，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] 范锦龙，吴炳方. 复种指数遥感监测方法[J]. 遥感学报，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] 范锦龙. 复种指数遥感监测方法研究[D]. 北京：中国科学院研究生院(遥感应用研究所)，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时序数据的冬小麦种植面积提取[J]. 湖北农业科学，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] 潘耀忠，李乐，张锦水，等. 基于典型物候特征的MODIS-EVI时间序列数据农作物种植面积提取方法——小区域冬小麦实验研究[J]. 遥感学报，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图像时间序列的冬小麦面积提取[J]. 农业工程学报，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] 王培娟，梁宏，谢东辉，等. 气候变暖对华北冬小麦返青前热量条件的影响[J]. 麦类作物学报，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)