基于Sentinel-2A NDVI时间序列数据的冬小麦识别

甄晓菊 张雪红 吴国明 傅晓艺 何泱 洪长桥

摘要：鉴于农作物类型识别中存在严重的“异物同谱”效应，基于归一化植被指数（normalized difference vegetation index，简称NDVI）时间序列数据及物候特征的农作物遥感识别已成为热点。针对现阶段NDVI时间序列数据空间分辨率普遍较低的问题，以河北省辛集市为研究区，基于Sentinel-2A数据构建了10 m高空间分辨率NDVI时间序列，并提出了积分法、斜率法和决策树法3种冬小麦识别模型，同时与传统的光谱角质图（spectral angle mapper，简称SAM）法进行了比较。结果表明，以上方法均达到了较好的识别效果，其中积分法、斜率法和决策树法的总体精度均优于97.6%，而SAM法因仅仅考虑了时间序列曲线的形态，使得稀疏林地与冬小麦之间容易误分;Sentinel-2A卫星（Sentinel-2星座重访周期为5 d）提供的高时空分辨率时间序列数据，在农作物的季相节律特征提取以及农作物的识别中具有巨大潜力。

关键词：Sentinel-2A;NDVI;冬小麦;物候期

中图分类号： TP79;S127  文献标志码： A

文章编号：1002-1302（2019）16-0239-06

收稿日期：2018-05-14

基金项目：中国博士后科学基金（编号：2017M610338）;河北省创新能力提升计划（编号：18964201H）;江苏省研究生科研与实践创新计划（编号：KYCX17_0892）;河北省气象与生态环境重点实验室开放研究基金（编号：Z201607Y）。

作者简介：甄晓菊（1981—），女，河北秦皇岛人，硕士，工程师，主要从事城市环境气象研究。

通信作者：张雪红，博士，副教授，主要從事农业遥感研究。

冬小麦是我国主要的粮食作物之一，及时、准确地获取冬小麦的空间分布信息，对于开展冬小麦长势监测和估产、区域粮食安全评估、国家/区域农业经济发展规划的制定、种植业结构调整指导等工作具有重要意义[1-3]。遥感技术具有快速、客观、探测面积大等特点，为准确、快速识别冬小麦种提供了重要的技术支撑[3-5]。

农作物的遥感识别方法主要包括两大类：单时相遥感影像识别、多时相遥感变化检测[4，6]。对于单时相遥感影像识别而言，主要基于农作物关键生长期的单期遥感影像，采用不同的分类方法，如非监督分类[7]、决策树分类[8]、支持向量机[9]等，或通过单期遥感影像提取植被指数[10]、叶面积指数[11]等，并设定阈值进行农作物信息的提取。由于“同物异谱、异物同谱”效应，给同期农作物识别造成大量的混分现象，识别精度难以保证。而农作物生长在短期内会发生强烈的变化，基于多时相遥感影像提供的农作物光谱变化及季相节律信息，可在一定程度上消除农作物相混的问题，提高作物的识别精度，保障农作物种类识别及面积的准确监测[2，4]。

随着高空间分辨率及高时间分辨率遥感的快速发展，基于时间序列的农作物遥感识别已成为研究热点。张晶等基于GF-1号影像获得了归一化植被指数（NDVI）时间序列数据，构建了空间向量法、曲线积分法和坐标转换法3种冬小麦识别模型，总体精度均达到98.65%[6];王利民等基于GF-1号多时相影像，结合波段反射率、波段反射率之和、波段反射率比值、NDVI等光谱特征，采用分层决策树分类法实现了对越冬前的冬小麦面积早期提取[12];杨闫君等基于GF-1 NDVI时间序列数据，采用最大似然法、神经网络分类、支持向量机等分类方法对花生、水稻、冬小麦—夏玉米等农作物进行分类，研究发现，支持向量机分类方法总体精度最高，达到 96.33%[13];杨闫君等基于NDVI时间序列数据，通过同时考察矢量的方向和大小建立了冬小麦识别矢量分析模型[14];姜涛等利用中分辨率成像光谱仪（MODIS） NDVI时间序列数据，基于冬小麦分蘖期内的平均NDVI值、苗期-三叶期和开花期-乳熟期2个时间段的NDVI值线性拟合斜率构建了冬小麦遥感识别模型[5];张喜旺等基于MODIS NDVI时间序列数据，通过综合季相节律和特征光谱对冬小麦种植面积进行了遥感估算[2];朱爽等提出了软硬变化检测（soft and hard change detection，简称SHCD）的作物识别方法，它通过同时考虑离散变化的纯净像元区和连续变化的混合像元区，不仅能够通过地表覆盖类型状态变化来有效地识别作物，而且可识别出地表覆盖类型的状态变化程度，从而可定量表达农作物的丰度[4]。农作物及其面积的识别精度不仅与时相信息有关，同时还与空间分辨率密切相关[15]。2015年6月23日发射升空的Sentinel-2A卫星，能够以10 d的重访周期（与2017年3月7日发射的Sentinel-2B卫星组成的Sentinel-2观测星座，重访周期可缩短为5 d）提供高空间分辨率（最高达10 m）的多光谱数据。Sentinel-2A同时具有较高的空间分辨率和时间分辨率，它在农作物的季相节律特征提取以及农作物的识别中具有巨大潜力。

在以往关于时间序列的农作物识别方面，主要针对中低分辨率影像，制约了识别精度的提高。鉴于Sentinel-2A同时具备高时间和高空间分辨率的优势，本研究拟以河北省辛集市为研究区，基于10 m分辨率的Sentinel-2A NDVI时间序列数据，结合冬小麦独特的季相节律特征及其种植方式，构建冬小麦遥感识别模型，开展冬小麦的遥感识别研究。

1 研究区与试验数据

1.1 研究区概况

研究区位于河北省石家庄市辛集市，中心地理位置为115.22°E，37.93°N，占地面积约为951 km2。辛集市地势平坦，处于北半球暖温带地区，大陆性季风气候特点显著，四季分明，光照充足，适宜作物生长。年平均气温为12.5 ℃，最冷月为1月（平均气温为-3.9 ℃）;最热月为7月（平均气温为 26.5 ℃），年平均降水量为534.6 mm，年平均日照时数为 2 513.9 h，年平均无霜期为197 d。该区的主要粮食作物为冬小麦和夏玉米，一年两熟轮作制度。辛集市是河北省冬小麦主要粮食产区之一，选定为冬小麦信息提取的研究区具有一定代表性。

2.2 冬小麦遥感识别方法

2.2.1 积分法

从图2可以看出，冬小麦从返青期至孕穗期NDVI值迅速提高，在此阶段冬小麦NDVI值显著高于其他地物，冬小麦的NDVI时间序列曲线随时间的积分值相对较大。因此，本研究充分利用这一特征构建冬小麦的积分法识别模型。积分的几何含义就是求算曲线与坐标轴之间的覆盖面积，因此积分值可转化为计算曲线面积Si，公式如下：

Si=（V5，i+V6，i）（t6，i-t5，i）+（V6，i+V7，i）（t7，i-t6，i）+（V7，i+V8，i）（t8，i-t7，i）2。（2）

式中：Si为像元i在t5至t8时相的NDVI积分值（本研究中对应于冬小麦返青期至孕穗期）;V5，i、V6，i、V7，i、V8，i为像元i在t5至t8时相的NDVI值;t5，i、t6，i、t7，i、t8，i为像元i的t5至t8时相。基于每个像元的Si，构建冬小麦的积分法识别模型：

地物类型=冬小麦其他地物，Si≥a。

式中：a为阈值参数。根据公式（2）并结合1 479个冬小麦训练样本及2 455个其他地物训练样本，将a设定为1.20。

2.2.2 斜率法

研究区内的所有植被类型中，唯有冬小麦在4月中下旬生长达到鼎盛阶段，6月进入成熟期，其NDVI时间序列曲线在4月下旬至6月出现快速下降，而其他作物和植被的NDVI时间序列曲线均处于不同程度的上升阶段。因此，通过计算4月下旬至6月的NDVI时间序列曲线的变化率，对冬小麦的识别具有巨大的潜力。为此构建冬小麦的斜率法识别模型，曲线的变化率在几何上的含义即为斜率Ki，公式如下：

Ki=V12，i-V8，it12，i-t8，i。（3）

式中：Ki为像元i在t8至t12时相的NDVI斜率;t8、t12时在本研究中对应冬小麦的孕穗期和成熟期;V8，i、V12，i为像元i在t8至t12时相的NDVI值。

基于每个像元的Ki，构建冬小麦的斜率识别模型：

地物类型=冬小麦其他地物，Ki≤a。

式中：a为阈值参数，根据公式（3）并结合1 479个冬小麦训练样本及2 455个其他地物训练样本，将a设定为-0.12。

2.2.3 决策树法

冬小麦—夏玉米属于一年两季的种植模式，其NDVI时间序列曲线具有多峰特征，并且在3至6月NDVI值出现的峰和谷是区别于其他植被和作物的显著特征之一。因此，通过建立决策树来提取此特征，从而用来识别冬小麦。结合训练样本数据，建立如下决策规则如下：

地物类型=冬小麦其他地物，V8，i≥0.6且V13，i<0.25。

2.2.4 SAM法

为了和本研究提取的以上3种方法进行对比，并鉴于冬小麦—夏玉米NDVI时间序列曲线在整个时间序列曲线与其他地物均存在明显差异的特点，也构建了常用的光谱角度制图（SAM）模型来识别冬小麦。算法原理：将每个像元的NDVI时间序列看作1个多维向量，以冬小麦训练样本的均值时间序列数据为参考向量，计算像元i与冬小麦参考向量之间的夹角（光谱角），通过光谱角θi的大小来判断待分像元i是否为冬小麦。像元i与冬小麦参考向量之间的光谱角θi定义为下式：

θi=cos∑15k=0Vk，i·Rk∑15k=0V2k，i·∑15k=0R2k。（4）

式中：Vk，i为像元i在时相K的NDVI值;RK为冬小麦在时相K的NDVI值。通过计算像元i的θi，并结合冬小麦训练样本，构建冬小麦的SAM识别模型如下：

地物类型=冬小麦其他地物，θi<0.25。

3 结果与分析

3.1 冬小麦提取结果及精度评价

基于“2.2”节中的4种冬小麦提取方法，对2016年至2017年辛集市的冬小麦进行遥感识别与提取，提取结果如图3所示。进一步采用2 790个验证样本（其中冬小麦1 253个，非冬小麦1 537个）对4种方法的提取结果进行精度评价。由表3可知，对于总体精度而言，决策树法最高，为9867%，SAM方法最低，为96.55%，而斜率法、积分法介于二者之间，分别为98.49%、97.67%;对于用户精度，積分法最高（97.83%），然后依次是决策树法（97.20%）、斜率法（96.90%）、SAM方法（93.27%）;而生产者精度中，积分法最低，为96.97%，另外3种方法均超过99%。

因此，总体而言，SAM方法因过多的非冬小麦像元被误分为冬小麦（错分误差为6.73%），从而造成总体精度及冬小麦用户精度最低;决策树法因冬小麦的生产者精度及用户精度均较高，使得其总体精度在4种方法中最高。

3.2 冬小麦识别算法对比

不同的冬小麦提取算法，其识别冬小麦的原理不同，从而导致识别的精度也存在差异。为了进一步探索导致冬小麦识别精度存在差异的原因，图4通过从研究区内截取1个典型区域，结合冬小麦识别结果，对各识别算法进行分析和比较。

对于积分法，由于其从冬小麦返青期至孕穗期（本研究从t5至t8时相）对NDVI时间序列曲线进行积分，而一些过于稀疏的冬小麦像元或冬小麦与其他地物的混合像元，其NDVI积分值偏小，容易和果树等林地相混淆，造成冬小麦像元的漏分。如图4中的红色矩形区域内，因冬小麦密度偏小，和其他3种识别算法相比，积分法（图4-a）对于稀疏的冬小麦漏分率最高，从而导致其生产者精度偏低。

SAM法是将整个NDVI时间序列作为1个向量，通过计算该向量的方向（与地物的类型密切相关）与参考向量的一致度来识别冬小麦，即考察待分像元NDVI时间序列曲线形态与典型地物的变化趋势一致性。本研究中SAM法提取的冬小麦用户精度最低（为93.27%），可能的原因是从t1至t8时相，林地NDVI时间序列曲线变化趋势同冬小麦的差异较小，再加上不同时相卫星空间采样的微小差异，使得部分稀疏林地或林地斑块边缘在t13时相附近的NDVI值偏低（如图4中的椭圆区域，误分为冬小麦的林地NDVI时间序列曲线如图5所示），从而导致其与冬小麦—夏玉米NDVI时间序列曲线的变化趋势相似，导致容易被误分为冬小麦。

对于斜率法及决策树法，前者利用冬小麦在t8至t12时相NDVI时间序列曲线陡然下降的独特特征，后者则结合冬小麦NDVI时间序列曲线在t8和t13时相分别出现峰和谷的特点提出冬小麦的识别算法，因此二者均充分利用了冬小麦与夏玉米的一年两季种植模式的特点，使得它们提取冬小麦时的各项精度指标均较高。

4 讨论

冬小麦独特的物候期以及一年两季的冬小麦—夏玉米种植模式，是基于NDVI时间序列数据进行冬小麦遥感识别的理论基础。针对冬小麦的物候期及种植模式，本研究构建了积分法、斜率法及决策树法来识别冬小麦，其中积分法、斜率法及决策树法主要是基于冬小麦2—6月的物候特征而建立的，此物候特征在NDVI时间序列上表现为显著的波峰。

积分法通过计算特定时间段内NDVI时间序列曲线所包围的面积来判定是否为冬小麦，对于冬小麦，在积分区间内其积分面积通常大于其他农作物及植被类型。由本研究可知，林地在10月及11月份（即冬小麦播种期和苗期）的NDVI值往往高于冬小麦[5-6]。因此，如果积分时段选择全生育期，则导致冬小麦与林地的积分面积差异变小，从而降低了该方法对冬小麦识别的敏感性。为此，结合Sentinel-2A的高时间分辨率的优势，本研究选择冬小麦关键生育期（返青期至孕穗期）进行积分，以提高积分法识别冬小麦的精度。

冬小麦NDVI时间序列曲线表现出从4月下旬至6月急剧下降的特点，而其他植被则刚好处于快速增长阶段。因此本研究通过计算此阶段NDVI值的变化率，提出了斜率法冬小麦识别模型。本研究表明，此方法的总体精度达到9849%，能很好地识别冬小麦。姜涛等利用MODIS NDVI时间序列数据，基于冬小麦分蘖期的平均NDVI值、苗期-三叶期和开花期-乳熟期2个时间段的NDVI值线性拟合斜率构建的冬小麦遥感识别模型也能有效对冬小麦进行识别[5]。相比较而言，本研究的斜率法仅要计算抽穗期-成熟期的斜率，算法更为简便，所需的NDVI时间序列更短。同理，决策树法仍然是基于冬小麦NDVI时间序列曲线4月至6月出现显著的“一峰一谷”特征，通过判断峰谷的位置及幅度来识别冬小麦。研究表明，该算法的总体精度最高（达98.67%）。

为了与本研究提出的冬小麦识别算法进行对比，同时采用传统的SAM法识别与提取冬小麦。由于SAM法主要考察的是向量的方向，即NDVI时间序列曲线的形状，而没有考虑NDVI的幅度。因此，部分稀疏的林地容易与冬小麦相混淆。针对SAM的局限性，张晶等通过同时考虑向量的方向及大小，以提供识别精度[6，14]。

农作物及其面积的识别与精度提取同时还与空间分辨率密切相关，王利民等对冬小麦面积识别精度与遥感数据分辨率之间的关系研究表明，当空间分辨率过小（如低于 100 m）或较高（优于10 m）时，识别精度变化较为缓慢，而在中等分辨率（15 m或30 m）区间时，用户精度变化明显[15]。本研究表明，Sentinel-2A同时具有较高的空间分辨率和时间分辨率的特点，为实现农作物的准确识别提供了可能。

5 结论

本研究以河北省辛集市为研究区，基于10 m分辨率的Sentinel-2A NDVI时间序列数据，结合研究区物候历，分别构建积分法、斜率法、決策树法以及SAM法等识别模型对冬小麦信息进行提取，得出以下结论：

（1）研究区内不同植被类型的物候期存在显著差异，不同农作物的种植模式不同，其中只有冬小麦—夏玉米为一年两季种植，使得其NDVI时间序列曲线在年周期内呈现“三峰三谷”的独特特征，尤其是2—6月的波峰是识别冬小麦的关键特征。

（2）构建的积分法、斜率法、决策树法以及SAM法等冬小麦识别模型，均取得了较高的分类精度，总体精度都高于96.5%，表明覆盖作物完整生长期的Sentinel-2A NDVI时间序列可以实现冬小麦的高精度分类与识别。其中斜率法、决策树法精度较高，总体精度均超过98%，这表明二者更适用于基于NDVI时间序列的作物分类。

（3）Sentinel-2A卫星能以10 d的重访周期提供高空间分辨率影像数据，其像元纯度高，有效地缓解了作物识别中混合像元的问题，识别精度较高，在农作物的季相节律特征提取以及农作物的识别中具有巨大潜力;此外，该数据的免费共享可为我国农情遥感监测业务提供高时间分辨率、高空间分辨率的遥感数据。

参考文献：

[1]Lobell D B，Asner G P. Cropland distributions from temporal unmixing of MODIS data[J]. Remote Sensing of Environment，2004，93（3）：412-422.

[2]张喜旺，秦耀辰，秦 奋. 综合季相节律和特征光谱的冬小麦种植面积遥感估算[J]. 农业工程学报，2013，29（8）：154-163，295.

[3]权文婷，王 钊. 冬小麦种植面积遥感提取方法研究[J]. 国土资源遥感，2013，25（4）：8-15.

[4]朱 爽，张锦水，帅冠元，等. 通过软硬变化检测识别冬小麦[J]. 遥感学报，2014，18（2）：476-496.

[5]姜 涛，朱文泉，詹 培，等. 一种抗时序数据噪声的冬小麦识别方法研究[J]. 遥感技术与应用，2017，32（4）：698-708.

[6]张 晶，占玉林，李如仁. 高分一号归一化植被指数时间序列用于冬小麦识别[J]. 遥感信息，2017，32（1）：50-56.

[7]邹金秋， 陈佑启，Uchida S，等. 利用Terra/MODIS数据提取冬小麦面积及精度分析[J]. 农业工程学报，2007，23（11）：195-200.

[8]Pal M，Mather P M. Support vector machines for classification in remote sensing[J]. International Journal of Remote Sensing，2005，26（5）：1007-1011.

[9]Xu M，Watanachaturaporn P，Varshney P K，et al. Decision tree regression for soft classification of remote sensing data[J]. Remote Sensing of Environment，2005，97（3）：322-336.

[10]王長耀，林文鹏. 基于MODISEVI的冬小麦产量遥感预测研究[J]. 农业工程学报，2005，21（10）：90-94.

[11]李开丽，蒋建军，茅荣正，等. 植被叶面积指数遥感检测模型[J]. 生态学报，2005，25（6）：1491-1496.

[12]王利民，刘 佳，杨福刚，等. 基于GF-1卫星遥感的冬小麦面积早期识别[J]. 农业工程学报，2015，31（11）：194-201.

[13]杨闫君，占玉林，田庆久，等. 基于GF-1/WFVNDVI时间序列数据的作物分类[J]. 农业工程学报，2015，31（24）：155-161.

[14]杨闫君，占玉林，田庆久，等. 利用时序数据构建冬小麦识别矢量分析模型[J]. 遥感信息，2016，31（5）：53-59.

[15]王利民，刘 佳，高建孟，等. 冬小麦面积遥感识别精度与空间分辨率的关系[J]. 农业工程学报，2016，32（23）：152-160.

[16]Liang S L. Quantitative remote sensing of land surfaces[M]. New Jersey：John Wiley & Sons，2003.

[17]Main-Knorn M，Pflug B，Debaocker V，et al. Calibration and validation plan for the l2A processor and products of the Sentinel-2 mission[C]//International Archives of the Photogrammetry Remote Sensing & S，2015，XL-7/W3（7）：1249-1255.