卫 炜，吴文斌，李正国，杨鹏，胡琼，周清波※

(1．中国农业科学院农业资源与农业区划研究所，北京 100081;2．农业部农业信息技术重点实验室，北京 100081)

植被指数是多光谱遥感数据经过空间变换或者不同波段之间的线性、非线性组合构成，其基本原理是利用植被的光谱特征，通过强化可见光 (主要是红光)与近红外波段反射率之间的差异来反映植被的生长状况［1］。时间序列植被指数能够较好反映植被的季节性生长发育过程，是监测不同尺度植被覆盖及其动态变化特征，提取植被生长季开始、生长最盛和结束等物候信息，描述全球气候变化及陆地生态系统响应的重要工具［2］。然而，星载传感器在数据采集和处理过程中受到各种因素 (如太阳高度角、传感器观测角、云、水汽和气溶胶等)干扰，使得植被指数时序曲线波动大并出现许多噪声［3］。因此，在进行各种趋势分析和信息提取之前，有必要先对时间序列植被指数进行去噪和平滑处理，即时序植被指数重构［4］。

尽管目前大多数植被指数产品都采用了最大值合成法来部分消除一些干扰因素的影响［5］，但是得到的时序植被指数仍然存在大量的噪声。为此，研究者们提出了许多时序植被指数重构方法，总体上可以归为滤波方法和函数拟合方法两大类。滤波方法分为时域滤波和频域滤波，时域滤波是在一个给定大小的滤波窗口内，采用某种方法对数据中的噪声进行处理，滤波窗口的大小与最终的滤波效果有关，是一种局部处理的方法;频域滤波则是通过数学变换将数据从时间域转换到频率域，噪声通常存在于高频部分，利用一个低通滤波器即可实现去噪，是一种整体处理的方法。函数拟合方法则是通过某种形式的函数对时序数据曲线进行最小二乘拟合，用拟合后的平滑曲线代替原来的时序数据曲线来实现平滑去噪。同样，函数拟合方法也存在局部处理和整体处理之分，其中经过局部拟合后的曲线还需要一个整体函数将它们连接起来以构成最终的拟合曲线。常用的滤波方法包括Savitzky-Golay滤波法、最佳指数斜率提取法、滑动中值/均值滤波法和基于傅里叶/小波变换的频域低通滤波法等;而常用的函数拟合方法包括非对称高斯函数拟合法、双Logistic函数拟合法和傅里叶函数拟合法等。这些方法有其各自的特点和区域适宜性，因此根据应用目的和研究区域对不同方法进行对比和评价具有十分重要的意义。国内外一些学者已经对部分重构方法进行了比较研究［6-9］，然而上述研究多是针对森林、草地等自然植被植被指数进行，且仅针对不同方法的重构效果进行了对比，很少考虑原始数据质量的影响和重构结果的应用情况。

耕地资源是保障国家粮食安全的基础性资源，耕地植被指数是进行农作物长势监测与估产以及农业生态系统对气候变化响应等研究的重要工具［10］。东北三省是我国重要的商品粮生产基地，在我国粮食安全保障体系和农业生产中占有重要地位［11］，同时该地区也是受气候变化影响最为显著的区域之一［12］。该研究基于耕地植被指数，在东北三省利用MODIS数据产品对非对称高斯函数拟合、双Logistic函数拟合和Savitzky-Golay滤波3种常见的时序植被指数重构方法进行比对研究。从原始数据质量对时间序列植被指数重构的影响、不同重构方法对原始数据的保真性和基于不同重构方法的耕地生长季长度空间格局等方面进行分析，以此来评价不同方法在该地区耕地时序植被指数重构中的适宜性。

1 研究区域和数据

1.1 研究区域

东北三省包括黑龙江省、吉林省和辽宁省，区域总面积79．18万km2，其中耕地面积为26．44万km2，约占全国耕地总面积的16．5%［13］(图1)。该地区属寒温带、中温带湿润 (半湿润)气候，大部分地区≥0℃积温为2 000～4 200℃，≥10℃积温为1 600～3 600℃，夏季平均气温为20～25℃，无霜期为140～170天，年降雨量500～800mm，其中60%集中在7～9月。受热量条件限制，该地区作物种植制度为1年1熟制，主要作物类型包括大豆、春玉米、春小麦和一季稻［14］。

图1 研究区域耕地范围及分布

1.2 数据来源

该研究采用的遥感数据为搭载在EOS/Terra卫星上的MODIS传感器MOD13A2数据产品，版本为v005。涉及到的数据内容包括基于最大值合成 (MVC)方法合成的16天1km空间分辨率的NDVI及其质量控制数据，该数据产品已经完成了几何校正和大气校正等处理。数据从美国NASA LP DAAC(Land Processes Distributed Active Archive Center)工作组官方网站下载 (https://lpdaac．usgs．gov/products/modis_products_table)，时间范围为2010年1月～2012年12月，3年时间共69期影像。数据产品中NDVI的DN值范围为-3 000～10 000，其中-3 000为填充值，将DN值转换成标准植被指数值的公式为:VI=DN*0．000 1。该研究使用的质量控制数据为MOD13A2中的Data Pixel Reliability数据集，其在时间和空间分辨率上与NDVI数据集相匹配。

该研究使用的土地利用数据来自中国1:25万土地覆盖遥感调查与监测数据库的2005年土地覆盖数据，该数据是在全国1:10万土地利用数据和遥感分类数据的基础上，并结合辅助资料、野外考察记录、遥感图像分析等通过全数字制图完成的。数据内容包括森林、草地、农田、聚落、湿地与水体、荒漠等6个一级类型，空间分辨率为100m。该研究将东北三省的农田数据提取出来作为耕地图层，并将空间分辨率重采样到1km使之与MOD13A2的NDVI数据相匹配，然后基于该图层提取东北三省耕地的植被指数进行分析。

2 研究方法

2.1 时序植被指数重构方法

星载传感器在获取地表信息时会受到诸多因素的影响，即使在经过了一定处理的遥感数据产品中仍然残留了很多噪声，导致时间序列植被指数曲线的周期性变化趋势不明显［15］。因此，对时序植被指数进行去噪和平滑等预处理，重构其时间序列曲线至关重要［16］。该研究主要对3种常见的时序植被指数重构方法:非对称高斯函数拟合法、双Logistic函数拟合法和Savitzky-Golay滤波法进行对比和分析。

2．1．1 非对称高斯函数拟合法

非对称高斯函数 (Asymmetric Gaussian，AG)拟合法是一种由局部拟合到整体拟合的方法，使用高斯函数来分段模拟植被的生长过程，最后通过平滑连接各段高斯拟合曲线实现时间序列的重构［17］。其主要过程大致可以分为区间提取、局部拟合和整体连接3个步骤。首先提取原始时序数据曲线中的谷值和峰值，采用高斯函数分别拟合曲线的左右部分。针对曲线突出部分拟合效果欠佳的问题，将曲线划分成左边谷值区、中部峰值区与右边谷值区，分别用不同的局部拟合函数进行描述，最后，再利用各局部拟合函数构建整体拟合函数。局部拟合函数为:

其中g(t;a1，…，a5)为高斯函数，式中c1和c2控制曲线的基准和幅度，a1决定峰值和谷值的位置;a4、a5和a2、a3分别控制曲线左、右部分的宽度和陡峭度。整体拟合函数为:

其中 ［tL，tR］是时序数据中待拟合部分的变化区间，fL(t)、fC(t)和fR(t)分别代表 ［tL，tR］区间内左边谷值、中间峰值及右边谷值对应的局部拟合函数，α(t)和β(t)为介于0和1之间的剪切系数。通过整体拟合函数将局部拟合函数连接起来是该方法的关键之一，这种从局部到整体的拟合策略避免了整体数据对局部拟合的干扰，拟合后的曲线更加接近真实情况。

2．1．2 双Logistic函数拟合法

双Logistic函数 (Double Logistic，DL)拟合法同样是基于从局部到整体的拟合思想，处理过程与AG拟合方法相类似，区别之处在于局部拟合函数为双Logistic形式且比AG函数少一个参数［18］。DL拟合的局部拟合函数为:

其中，g(t;a1，…，a4)为双Logistic函数，式中c1和c2控制曲线的基准和幅度;a1、a2和a3、a4分别控制曲线左、右部分的拐点位置及拐点处的变化速率，整体拟合函数与AG拟合相同。

2．1．3 Savitzky-Golay 滤波法

Savitzky-Golay(S-G)滤波法是Savizky和Golay于1964年提出的，它是一种通过局部多项式回归模型来平滑时序数据的时域低通滤波方法［19］。S-G滤波的基本思想是基于多项式，在滤波窗口内利用最小二乘法对数据进行最佳拟合，可以简单的理解为一种加权平均算法，权重由在一个滤波窗口范围内做最小二乘拟合的多项式次数来决定［20］。S-G滤波公式为:

其中，Y为原始数据，Y'为拟合值，Ci为第i个点的权重，N=2n+1为滤波窗口的大小。S-G滤波法对滤波窗口的大小非常敏感，滤波窗口的宽度设置偏小容易产生大量冗余数据;反之则可能遗漏一些细节信息［21］。另外，拟合多项式的次数也会对平滑效果产生影响，次数较低时结果较为平滑，当次数较高时则会导致过度拟合。

2.2 不同重构方法对比策略

该研究主要从3个方面对AG拟合、DL拟合和SG滤波3种时间序列植被指数重构方法进行比较。首先对NDVI及其质量控制数据进行分析，探讨原始数据质量对时序植被指数重构的影响;然后选择若干随机抽样点，从定性和定量两个角度来评价不同方法对时序植被指数的重构效果;最后借助研究区的耕地物候特征，间接对不同的时序植被指数重构方法进行对比说明。

2．2．1 原始数据质量的影响分析

根据Data Pixel Reliability数据集的DN值可将所有像元的数据质量划分为5个等级 (表1)［22］。由表1可知，MOD13A2数据产品中植被指数的数据质量主要受到冰雪覆盖和云掩的影响。在时序植被指数进行重构过程中，这些数据质量受到影响的像元点可能会使真实的植被指数曲线发生扭曲或偏移。本研究以Data Pixel Reliability数据集的DN值为依据对各个数据点赋予不同的权重，使得重构后的植被指数时序曲线最大限度的接近真实数据。具体的赋权规则如下:DN值为0的点权重设为1，DN值为1的点权重设为0．8，而DN值为2和3的点权重设为0．2［22］。对某一像元的NDVI数据分别利用AG拟合法、DL拟合法和SG滤波法进行时序植被指数重构，对比未引入质量控制数据和引入质量控制数据的拟合曲线，以此来评价原始数据质量对不同方法重构时序植被指数的影响。

表1 MOD13A2 v005产品Data Pixel Reliability数据集DN值说明［22］

2．2．2 抽样点数据重构效果对比

在耕地范围内随机选择100个抽样点并获取其原始的时间序列植被指数数据，分别用AG拟合、DL拟合和SG滤波3种方法对植被指数时序曲线进行重构，将不同方法的重构结果以及原始数据进行定性对比分析。在由时序植被指数提取物候信息和应用全年或生长季植被指数平均状况的研究中，要求重构算法在去除噪声的同时最大限度的保留原始高质量数据的真实性。因此，该研究从重构后的植被指数时序曲线保持原始植被指数曲线特征的能力及其与原始数据之间的平均差异程度两个方面进行分析。其中高质量数据点是进行时序植被指数重构的重要依据，它们代表了真实的植被指数时序曲线的波形特征。为了定量评价不同方法的重构效果，计算了重构前后的植被指数时序曲线之间的回归系数 (R2)和均方根误差(RMSE)，其中R2反映了二者之间的相关程度，而RMSE则反映了二者之间的平均差异程度。具体做法是:分别计算每个样本像元利用不同方法重构的植被指数时序曲线与原始植被指数曲线之间的R2和RMSE，对所有像元统计这两个指标的平均值，当R2值越高或RMSE值越低时表明重构的植被指数时序曲线对原始数据的保真性越好。考虑到原始数据的质量问题，在定量评价时分全年数据和高质量数据进行比较，其中高质量数据是指在Data Pixel Reliability数据集中DN值为0和1的像元比重较高时段的数据。

2．2．3 研究区耕地物候特征对比

一般而言，耕地的植被指数在以年为周期的变化中会经历一个从升高到峰值再到降低的动态过程［23］。这种总体的变化特征可以描述耕地作物的生长过程。因此，时序植被指数能够用来提取耕地物候信息。对时序植被指数进行去噪和平滑等预处理之后，即得到重构的植被指数时序曲线，在此基础上就可以开展物候信息提取的工作［24］。该研究基于AG拟合、DL拟合和SG滤波3种方法重构的植被指数时序曲线，分别提取研究区的耕地物候参数及其空间格局，通过耕地物候特征间接对比不同的时序植被指数重构方法，作为前面抽样点对比的补充。提取的耕地物候参数包括耕地的生长季开始期和结束期并以此来计算生长季长度，物候参数提取方法采用比例阈值法［17］。比例阈值法能够消除不同土壤背景和植被类型的影响，但是在比例参数的选择上具有一定的经验性［25］。针对农作物等植被生长季NDVI曲线不对称的问题，可以将生长季开始和结束的比例参数设置为不同值。对于东北三省而言，理论上耕地生长季开始期可以对应作物的出苗期，而生长季结束期则对应作物的成熟期［26］。该研究结合地面观测数据将耕地生长季开始和结束期的阈值分别设置为10%和50%。

3 结果与分析

3.1 原始数据质量对植被指数时序曲线重构的影响

图2是对东北三省2011年全年23期Data Pixel Reliability数据的DN值进行统计的结果，其中DN值等于0和1的点被认为是高质量数据，而DN值等于2和3的点被认为是低质量数据。从中可以发现，前5期和后3期影像中除第21期低质量数据比重略小于高质量数据之外，其它各期影像低质量数据的比重都大于高质量数据，特别是前4期影像中低质量数据比重远远大于高质量数据。即这些时段的影像整体上属于低质量数据，可能会对时间序列植被指数重构产生影响。而从第6期影像起高质量数据的比重开始远远大于低质量数据，从第7期一直到第20期影像高质量数据的比重始终保持在一个较高的水平，同时低质量数据比重维持在一个非常低的水平。因此，这些时段的影像整体上属于高质量数据，是进行时间序列植被指数重构的重要依据，也是该研究中对3种方法的重构效果进行比较时需要重点分析的数据。

图2 东北三省MOD13A2数据产品2011年各期影像数据质量统计

图3为分别利用AG拟合、DL拟合和SG滤波3种方法在无质量控制和引入质量控制的情况下得到的单个像元在2011年的时序NDVI重构结果以及原始数据的对比。从图中可以看出，前5期和后3期的原始植被指数明显偏低，与前面对Data Pixel Reliability数据DN值的统计结果分析相对应。3种方法引入质量控制前后的重构曲线之间差异明显，主要表现为:AG拟合和DL拟合的重构结果在低质量数据时段更加偏离原始曲线，而在高质量数据时段更加符合原始曲线，且二者的重构曲线除后3期外在引入质量控制后几乎趋于一致;SG滤波的重构结果在低质量数据时段的情况与前两种方法类似，而在高质量数据时段引入质量控制前后变化不大。原始数据质量的影响使得无质量控制的重构曲线在上升和下降阶段的最小值明显低于引入质量控制的重构曲线，这样就可能会给下一步的趋势分析和信息提取带来误差。考虑到原始数据质量影响的问题，该研究在后续利用3种方法对NDVI时序数据进行重构的过程中统一引入质量控制数据。

图3 单个像元引入质量控制前后3种方法重构的时序植被指数以及原始数据对比

3.2 3种重构方法对原始数据的保真性对比

图4中4条曲线分别为100个抽样点2011年的原始植被指数和经过3种方法重构后的植被指数取平均值得到的时间序列曲线。从整体上看，AG拟合和DL拟合得到的结果非常接近，而SG滤波得到的结果则与上述两种方法存在差异。在低质量数据时段，即影像的前5期和后3期时段，3种方法得到的重构曲线相对于原始的植被指数曲线都有所提升，但SG滤波结果的提升幅度小于前两种方法。在高质量数据时段，AG拟合和DL重构曲线的上升和下降阶段与原始曲线吻合度较高，而在曲线峰值附近拟合曲线较原始曲线有所提升;SG滤波得到的结果在重构曲线的上升和下降阶段相对原始曲线分别出现左偏和右移，而在曲线峰值附近较原始曲线有所降低。其中SG滤波重构曲线的左偏和右移将会导致利用比例阈值法提取的生长季开始期提前和结束期推迟，从而使得计算的生长季长度延长，另外拟合峰值的降低则进一步扩大了这种误差的范围。

图4 3种方法重构得到的植被指数均值时序曲线与原始植被指数时序曲线对比

图5为分别用AG拟合、DL拟合和SG滤波3种方法对100个抽样点2011年全年23期以及代表高质量数据的第6～20期的植被指数进行重构前后的时序曲线之间的回归系数和均方根误差对比。从全年数据来看，SG滤波的R2值最高，AG拟合次之，DL拟合最低;RMSE值SG滤波最低，AG拟合居中，DL拟合最高。说明SG滤波法对原始数据的保真性最好，AG拟合次之，DL拟合最差，但值得注意的是这3种方法的R2和RMSE值之间差距并不大。对高质量数据而言，整体上AG拟合和DL拟合的R2和RMSE值非常接近，而SG滤波与上述两种方法之间差距明显。其中AG拟合和DL拟合的R2值明显高于SG滤波，且AG拟合比DL拟合略高;而AG拟合和DL拟合的RMSE值则明显低于SG滤波，且AG拟合比DL拟合略低。由此说明，利用SG滤波重构的植被指数时序曲线对原始高质量数据的保真性最差，AG拟合与DL拟合重构的植被指数时序曲线对原始高质量数据的保真性基本上相接近且AG拟合略占优势。

图5 3种方法重构前后植被指数时序曲线高质量数据的回归系数和均方根误差对比

3.3 基于不同重构方法的耕地生长季长度空间格局

该研究分别利用AG拟合、DL拟合和SG滤波3种方法对整个东北三省2010～2012年共69期的NDVI时序数据进行重构，然后由比例阈值法提取了2011年该地区耕地的生长季开始期和结束期并以此来计算生长季长度。

3．3．1 耕地生长季长度的空间分布特征比较

图6为分别基于3种时序植被指数重构方法得到的东北三省2011年耕地生长季长度的空间分布特征。可以看出，基于AG拟合和DL拟合提取的耕地生长季长度空间分布基本一致，而基于SG滤波的提取结果与前两种方法差异明显。在基于AG拟合和DL拟合的提取结果中，吉林省西部、黑龙江省松嫩平原西北部以及三江平原的部分地区耕地生长季长度一般为90～120天，三大平原其它大部分地区的耕地生长季长度为120～135天，另外黑龙江省东部和南部的零星区域、吉林省中南部以及辽宁省东部和西部的部分地区耕地生长季长度为135～150天，该结果与李正国等［26］的研究结果基本一致;而基于SG滤波的提取结果中耕地生长季长度较前两种方法则明显延长，大部分区域的耕地生长季长度在150天以上，只有少部分区域为120～150天。这与之前定性分析3种方法对原始高质量数据的保真性时的发现相对应，即SG滤波重构曲线在上升和下降阶段相对于原始曲线分别出现左偏和右移以及在峰值附近有所降低，由此导致利用比例阈值法提取的生长季开始期提前和结束期推迟，进而使得计算的生长季长度延长。

图6 东北三省2011年耕地生长季长度空间分布特征(从左到右依次为基于AG拟合、DL拟合和SG滤波得到的结果)

3．3．2 不同生长季长度区域的重构曲线比较

根据东北三省2011年耕地生长季长度的空间分布情况，从中分别选择生长季短、中和长3种类型的随机抽样点，将原始NDVI曲线和由3种方法得到的重构曲线进行比较分析，如图7。可以看出，在耕地生长季长度较短的区域3种方法得到的重构曲线在低质量数据时段都较原始曲线有所抬升，其中AG拟合和DL拟合在前5期略有差异但后3期几乎一致，SG滤波与上述两种方法则差异明显;在高质量数据时段AG拟合和DL拟合结果与原始曲线比较吻合，但SG滤波重构曲线在上升和下降阶段较原始曲线分别出现左偏和右移。在耕地生长季长度居中的区域总体上与短生长季的区域相类似，不同的是低质量数据时段AG拟合和DL拟合在前5期几乎一致而在后3期略有差异，同时SG滤波与上述两种方法的差异增大;高质量数据时段AG拟合和DL拟合的结果更加符合原始曲线，SG滤波重构曲线在上升阶段的左偏程度减小，拟合峰值也有所降低。在耕地生长季长度较长的区域原始数据质量影响不大，3种方法得到的重构曲线非常接近且与原始曲线的吻合程度较高。

图7 不同生长季长度的随机抽样点3种方法重构的时序植被指数及原始数据对比

4 结论与讨论

4.1 结论

该研究以东北三省为研究区域，利用MODIS植被指数产品对比了AG拟合、DL拟合和SG滤波3种常见的时序植被指数重构方法。主要从原始数据质量对时间序列植被指数重构的影响、不同重构方法对原始数据的保真性以及基于不同重构方法的耕地生长季长度空间格局3个方面进行分析。研究结果表明:

(1)原始数据质量对利用3种方法进行时序植被指数重构都存在影响，并且对AG拟合和DL拟合重构结果的影响大于SG滤波。主要表现为冰雪覆盖和云掩造成的植被指数值降低，导致无质量控制的重构曲线在上升和下降阶段的最小值明显低于引入质量控制的重构曲线，这样就可能给下一步的趋势分析和信息提取带来误差。

(2)从整体上看，AG拟合和DL拟合得到的结果非常接近，SG滤波得到的结果则与上述两种方法存在差异。从全部原始数据来看，3种方法中SG滤波重构的时序植被指数的保真性最好，AG拟合次之，DL拟合最差;但是对高质量数据而言SG滤波的保真性则最差，AG拟合与DL拟合的保真性相接近且AG拟合略占优势。

(3)基于AG拟合和DL拟合的耕地生长季长度在空间分布特征基本上一致，并且与前人在该地区利用VGT植被指数得到的研究结果相符合。基于SG滤波的耕地生长季长度空间格局与上述两种方法差异显著，大部分区域的耕地生长季长度明显延长。在生长季较短和居中的区域SG滤波重构曲线在上升和下降阶段相对原始曲线分别出现左偏和右移，而在生长季长度较长的区域与原始曲线的吻合程度较高。

总体上看，AG拟合和DL拟合重构结果在各个方面都表现出非常相似的特点，且基于这两种方法得到的耕地物候特征更加符合实际情况，而SG滤波则与前两者明显不同。通过上述研究发现，在东北三省利用AG拟合、DL拟合和SG滤波进行耕地时序植被指数重构时需要引入质量控制数据，AG拟合和DL拟合在该地区耕地时序植被指数重构时表现相似且比SG滤波更加适宜。

4.2 分析和讨论

造成AG拟合、DL拟合和SG滤波3种方法重构结果差异性的原因主要与这些方法的自身特点密切相关。AG拟合和DL拟合属于函数拟合方法，是从局部到整体的处理方法，受到拟合函数形式的影响较大且得到的重构曲线更加平滑;SG滤波属于时域滤波方法，是一种局部处理的方法，重构曲线的平滑程度与滤波窗口的大小有关并且受原始数据的影响较大。AG拟合和DL拟合在原理上基本相同，只是局部拟合函数形式有所差异，因此这两种方法的重构曲线比较接近，由此得到的耕地物候空间格局也非常相似。SG滤波的原理与前两种方法截然不同，因此导致了重构曲线以及耕地物候空间格局的明显差异，其中SG滤波重构曲线在上升和下降阶段的左偏和右移是所提取的耕地生长季长度明显延长的直接原因。除此之外数据情况也与3种方法重构结果的差异性有关，该研究使用的植被指数来自Terra卫星上的MODIS传感器，虽然经过了MVC方法合成，但残余的云掩和冰雪覆盖影响仍然存在，另外1km的空间分辨率和16天的时间分辨率也造成了部分细节信息的损失。

AG拟合、DL拟合和SG滤波3种时序植被指数重构方法各具优势，但是在不同地区和不同研究目的时仍然需要进行适宜性分析。理论上以强调局部细节为特点的SG滤波法对原始数据的保真性应该最高，然而，在东北三省由于受到冰雪覆盖及云掩导致的原始数据质量下降影响，使其对时序植被指数重构的效果大打折扣，进而影响到物候信息的提取，以整体拟合为特点的AG拟合和DL拟合则能够在一定程度上避免这种影响。另外在用比例阈值法提取耕地物候信息时比例参数的选择非常关键，20%是目前应用比较普遍的比例参数，然而大部分研究都是以自然植被为对象的。尽管该研究通过结合地面观测数据将耕地生长季开始期和结束期的比例参数分别设置为10%和50%，但在缺乏地面数据的情况下如何为受人为因素影响较大的耕地选择一个能够准确获取耕地物候信息的比例参数仍然值得商榷。

［1］王正兴，刘闯，Alfredo Huete．植被指数研究进展:从AVHRR-NDVI到MODIS-EVI．生态学报，2003，23(5):979～987

［2］Zhang X，Friedl M A，Schaaf C B，et al．Monitoring vegetation phenology using MODIS．Remote Sensing of Environment，2003，84(3):471～475

［3］李儒，张霞，刘波，等．遥感时间序列数据滤波重建算法发展综述．遥感学报，2009，13(2):335～341

［4］吴文斌，杨鹏，唐华俊，等．两种NDVI时间序列数据拟合方法比较．农业工程学报，2009，25(11):183～188

［5］Taddei R．Maximum Value Interpolated(MVI):A Maximum Value Composite method improvement in vegetation index profiles analysis．International Journal of Remote Sensing，1997，18(11):2365～2370

［6］Beck P S A，Atzberger C，H?gda K A，et al．Improved monitoring of vegetation dynamics at very high latitudes:A new method using MODIS NDVI．Remote Sensing of Environment，2006，100(3):321 ～334

［7］Hird J N，McDermid G J．Noise reduction of NDVI time series:An empirical comparison of selected techniques．Remote Sensing of Environment，2009，113(1):248～258

［8］曹云锋，王正兴，邓芳萍．3种滤波算法对NDVI高质量数据保真性研究．遥感技术与应用，2010，25(1):118～125

［9］宋春桥，游松财，柯灵红，等．藏北地区三种时序NDVI重建方法与应用分析．地球信息科学学报，2011，13(1):133～143

［10］李正国，唐华俊，杨鹏，等．植被物候特征的遥感提取与农业应用综述．中国农业资源与区划，2012，33(5):20～28

［11］任国玉，初子莹，周雅清，等．中国气温变化研究最新进展．气候与环境研究，2005，(4):701～716

［12］刘志娟，杨晓光，王文峰，等．气候变化背景下我国东北三省农业气候资源变化特征．应用生态学报，2009(09):2199～2206

［13］李正国，杨鹏，唐华俊，等．气候变化背景下东北三省主要作物典型物候期变化趋势分析．中国农业科学，2011，44(20):4180～4189

［14］李正国，唐华俊，杨鹏，等．东北三省耕地物候期对热量资源变化的响应．地理学报，2011，66(7):928～939

［15］李杭燕，颉耀文，马明国．时序NDVI数据集重建方法评价与实例研究．遥感技术与应用，2009，24(5):596～602

［16］孙华生，徐爱功，林卉，等．基于不同算法的时间序列植被指数去噪效果分析．江苏农业科学，2012，40(5):375～379

［17］Jönsson P，Eklundh L．Seasonality extraction by function fitting to time-series of satellite sensor data．IEEE Transactions on Geoscience and Remote Sensing，2002，40(8):1824～1832

［18］Jönsson P，Eklundh L．TIMESAT—a program for analyzing time-series of satellite sensor data．Computers ＆ Geosciences，2004，30(8):833～845

［19］Chen J，Jönsson P，Tamura M，et al．A simple method for reconstructing a high-quality NDVI time-series data set based on the Savitzky-Golay filter．Remote Sensing of Environment，2004，91:332 ～344

［20］蔡天净，唐瀚．Savitzky-Golay平滑滤波器的最小二乘拟合原理综述．数字通信，2011，(1):63～68

［21］顾娟，李新，黄春林．NDVI时间序列数据集重建方法述评．遥感技术与应用，2006，21(4):391～395

［22］宋春桥，柯灵红，游松财，等．基于TIMESAT的3种时序NDVI拟合方法比较研究——以藏北草地为例．遥感技术与应用，2011，26(2):147～155

［23］吴文斌，杨鹏，唐华俊，等．过去20年中国耕地生长季起始期的时空变化．生态学报，2009，29(4):1777～1786

［24］胡砚霞，黄进良，王立辉．基于MODIS NDVI时序数据的丹江口库区物候格局研究．中国农业资源与区划，2012，33(4):25～30

［25］吴文斌，杨鹏，唐华俊，等．基于NDVI数据的华北地区耕地物候空间格局．中国农业科学，2009，42(2):552～560

［26］李正国，唐华俊，杨鹏，等．基于时序植被指数的东北地区耕地生长季特征识别与应用研究．北京大学学报 (自然科学版)，2011，47(5):882～892