英语原文共 14 页,剩余内容已隐藏,支付完成后下载完整资料
一种基于Savitzky-Golay滤波重建高质量NDVI时间序列数据集的简单方法
期刊:《环境遥感》91 (2004) 332–344
作者及单位:Jin Chen(a,b,*), Per.Jonsson (c), Masayuki Tamura (b) , Zhihui Gu (a) , Bunkei Matsushita (b),Lars Eklundh (d)
a: Key Laboratory of Environmental Change and Natural Disaster, Beijing Normal University, Ministry of Education of China Beijing 100875, China
b: Social and Information System Division, Japan National Institute for Environmental Studies, 16-2 Onogawa, Tsukuba, 305-8506, Japan
c: Division of Mathematics, Natural Sciences and Language, Malmouml; University, Malmouml;, Sweden, and the Department of Physics, Lund University, Lund, Sweden
d: Department of Physical Geography and Ecosystems Analysis, Lund University, Lund, Sweden
摘要:虽然来自NOAA/AVHRR,SPOT/VEGETATION,TERRA或AQUA/MODIS的归一化差异植被指数(NDVI)时间序列数据已经成功地应用于关于全球环境变化的研究,但NDVI时间序列数据中的残留噪声即使采用严格的预处理,也会妨碍进一步的分析,并存在产生错误结果的风险。根据NDVI时间序列伴随植被生长和衰退的年周期的假设,云或恶劣的大气条件通常会降低NDVI值。本研究开发了一种基于Savitzky-Golay滤波算法的简单而强大的方法以平滑NDVI时间序列中的噪声,特别是由云污染和大气变化引起的噪声。我们的方法是为了使数据接近NDVI上限,并通过迭代过程反映NDVI模式的变化。通过将新开发的方法应用于10天MVC SPOT VGT-S产品,我们为新方法提供了优化参数,并将其与BISE算法和基于傅立叶的拟合方法进行了比较。我们的研究结果表明,新方法在获得高品质NDVI时间序列方面有更好的效果。
1.介绍
因为NDVI载有有关地表物业的有价值信息,来自NOAA / AVHRR,SPOT / VEGETATION,TERRA或AQUA / MODIS的归一化差异植被指数(NDVI)的时间序列数据已被证明适用于检测长期的土地利用/覆盖变化以及对全球、大陆和地区范围内陆地生态系统进行建模(IGBP, 1992; Justice et al., 1985; Myneni et al., 1997;Potter et al., 1993; Prince, 1991; Reed et al., 1994;Running amp; Nemani, 1988; Tucker amp; Sellers, 1986; Tucker et al., 1985及其他相关研究)。理论上,由近红外(NIR)和红外反射率的归一化变换计算的NDVI是电磁光谱红光和近红外部分植被的吸收和反射特性的指标。因此,NDVI时间序列的变化表明植被状况的变化与有效光合辐射的吸收比例(Sellers, 1985)。然而,这些时间序列几乎总是受到云污染,大气变化和双向效应的干扰。这些干扰大大影响了土地覆盖类型和陆地生态系统的监测,并显示为系统噪声(Cihlar et al., 1997; Gutman, 1991)。虽然最常用的NDVI数据集是10天最大值复合(MVC)产品(Holben,1986),例如Pathfinder土地数据集,GIMMS NDVI数据集和SPOT VGT产品,但其中仍包括很多这样的噪声。为此,在过去二十年中,许多减少噪声和构建高质量NDVI时间序列数据的方法被进一步制定和研究。这些方法可以大致分为三种一般类型:(1)基于阈值的方法,如最佳坡度指数提取算法(BISE) (Viovy et al., 1992); (2)基于傅立叶的拟合方法(Cihlar, 1996; Roerink et al., 2000; Sellers et al., 1994); (3)非对称函数拟合方法,如非对称高斯函数拟合方法(Jonsson&Eklundh,2002)和加权最小二乘线性回归方法(Swets et al., 1999)。上述方法具有各自的优势,并已成功应用于某些应用的NDVI时间序列预处理。BISE算法已被用于提取植被覆盖度季节性指标(如Reed等,1994),对植被或土地覆盖类型进行分类(Lovell&Graetz,2001; Xiao et al。,2002),并估计 初级生产力(GPP)和净初级生产力(NPP)(Ruimy等,1996)。已经采用基于傅里叶的拟合方法来获得陆生生物物理参数(如Sellers等人,1994),并评估NPP动力学(如Malmstrom等,1997)。不对称函数拟合方法主要用于提取物候学研究的季节性信息(Jonsson&Eklundh,2002)。
然而,这些方法也具有限制使用的若干缺点(Jonsson&Eklundh,2002)。例如,BISE算法需要根据通常主观的经验策略,确定滑动周期和滑动区间内NDVI再生长可接受百分比增加的阈值,结果往往取决于分析师的技能和经验。因此,与其他基于阈值的方法一样,应用BISE算法后的剩余噪声可能使提取的时间信息不可靠。当应用于不规则或非对称NDVI数据时,基于傅立叶的拟合方法可能是有问题的,因为它们主要依赖对称的正弦和余弦函数。此外,它们可能在NDVI时间序列中产生杂散振荡。与上述方法相比,非对称高斯函数拟合方法在获得高质量NDVI时间序列方面更加灵活有效。然而,可能难以确定可以拟合本地参数的合理和一致的最大值和最小值集合,特别是对于噪声数据或不存在明显季节性的区域的数据。此外,这种方法的复杂性使其更耗时。
鉴于上述缺点,本文提出了一种基于Savitzky-Golay滤的简单而强大的方法,以更有效地减少主要由于云污染和大气变化引起的NDVI时间序列中的噪声污染。我们的方法是为了使数据接近NDVI上限,并通过迭代过程描绘NDVI变化的模式。该方法用由法国、欧盟委员会、比利时、意大利和瑞典联合开发的VEGETATION计划产生的10天MVC SPOT VGT-S产品进行了测试。
2.方法
与其他减少噪声和构建高质量NDVI时间序列的策略类似,我们的方法是基于两个假设:(1)卫星传感器的NDVI数据主要与植被变化有关。 因此,NDVI时间序列遵循增长和衰退的年度周期;(2)云和恶劣的大气条件通常会使NDVI值下降,要求与植被变化过程不相符的突然下降的NDVI值被认为是噪声并予以去除。根据这两个假设,开发了基于Savitzky-Golay滤波算法的新方法,以使数据接近上NDVI包络线,并通过迭代过程得到最适合全植被季节的NDVI变化规律。该方法可应用于以不同间隔采样的NDVI数据集,包括日数据、10天或月MVC。此外,NDVI和特定传感器的缩放没有任何限制。 在下文中,首先简要介绍Savitzky-Golay滤波方法,然后根据图1所示的流程图描述实现新方法的主要步骤。
2.1. Savitzky-Golay滤波
Savitzky和Golay(1964)提出了一种简化的最小二乘拟合卷积,用于平滑和计算一组连续值(谱)的导数。卷积可以被理解为加权移动平均滤波器,其加权为一定程度给定的多项式。加权系数(以下称为系数)当应用于信号时,在滤波器窗口内执行多项式最小二乘法拟合。该多项式被设计为保留数据内的更高的值并减少由滤波器引入的偏置。当数据点沿所选横坐标处于固定且均匀的移动时,该滤波器可以应用于任何连续数据,并且通过绘制点形成的曲线必须是连续的并或多或少平滑的。NDVI时间序列清楚地满足了这些条件。NDVI时间序列平滑的简化最小二乘卷积的一般方程如下:
hellip;hellip;hellip;hellip;hellip;hellip;hellip;hellip;hellip;hellip;hellip;hellip;hellip;hellip;hellip;hellip;hellip;hellip;.(1)
其中Y是原始NDVI值,Y*是最终的NDVI值,Ci是滤波器(平滑窗口)的第i个NDVI值的系数,N是卷积整数的数量,并且等于平滑窗口大小( 2m 1)。参数j是原始纵坐标数据表的运行参数。平滑阵列(滤波器大小)由2m 1个点组成,其中m是平滑窗口的半宽度。Savitzky-Golay滤波器(Ci)的系数可以直接从Steinier等人的研究获得(1972年)作为Savitzky和Golay研究(1964年)的修正版本,或由Madden(1978)提出的方程式计算。
图1 新的滤波方法流程图
当滤波器应用于NDVI时间序列平滑时,需要根据NDVI观测值确定两个参数。第一个参数是m,即平滑窗口的半宽度。通常情况下,m的值越大,在尖锐的峰值处产生的平滑效果更为明显。第二个参数是指定平滑多项式的程度的整数(d),其通常设置在从2到4的范围内。d的值较小将产生更平滑的效果,但可能引入偏差;较高的d值会降低滤波器偏差,但可能会“过度拟合”数据,并产生出现更多噪声的结果。
2.2.新方法的实施
2.2.1.第1步:多云区NDVI值的线性插值
大多数NDVI时间序列数据集,例如Pathfinder土地数据集或SPOT VGT产品都将云标志波段作为辅助数据(Stowe等,1991);这些数据提供了时间序列中每个数据点的云状态的有价值的指标,尽管这些数据不包括NDVI数据点受云和恶劣大气条件影响的所有情况。利用这些标志数据来估计NDVI值的不确定度是非常重要的(Jonsson&Eklundh,2002)。在本研究中,云标志数据用于通过线性插值来改善NDVI时间序列。具体来说,假设有NDVI时间序列的数据点(ti,Ni,Fi),i = 1,2,3 ... n,其中ti是日期,Ni是NDVI值,Fi是云标志,如果第j个点的Fj被识别为混浊点,则使用不被识别为混浊点的相邻点将Nj替换为线性内插值。此外,20天内随机NDVI增加大于0.4的点也被剔除,并用相邻点的线性内插值代替,因为这种增加不会由于自然植被变化引起。因此,我们获得了新的NDVI时间序列数据点(ti,),i = 1,2,3 ... n,其中ti是日期,是线性插值后的新NDVI值,图2a用圆圈表示出了根据云标志的多云区NDVI值的线性插值的示例。可以看出,云标志识别的两个点显示NDVI值突然下降,而NDVI值中的其他突然下降不会被云标记指示,这意味着根据对可见光和近红外反射率(VNIR)以及热量带应用某些阈值对云标记进行分类是困难的。
图2以新开发方法(第8号测试像元的NDVI数据)的不同步骤显示NDVI时间序列的示例:(a)原始NDVI时间序列和云标志点(图中圈出);(b)由Savitzky-Golay滤波器安装的长期变化趋势曲线(粗实线);(c)Savitzky-Golay滤波器首次拟合的NDVI时间序列。将第一个拟合的NDVI时间序列绘制为粗实线,并将线性内插NDVI时间序列(ti,Ni0)绘制为实线。 (d)使用Savitzky-Golay滤镜(粗实线)的最终NDVI时间序列。
2.2.2.第2步:使用Savitzky-Golay滤波器的长期变化趋势拟合
根据上述假设,NDVI时间序列应遵循年度植被周期的变化过程,NDVI时间序列的突然下降与规律不相符,可以视为受云或具有噪声的较差的大气条件影响。因此,如果我们能够获得代表年度植被周期变化过程的长期变化趋势曲线,则有助于确定这些带有噪声的数据点,并将其作为对其他配件的重要参数(Jonsson&Eklundh,2002)。要获得令人满意的长期变化趋势曲线,需要考虑两个标准:(1)NDVI时间序列的长期变化趋势曲线应遵循年度植被周期的变化过程,并且不会有太大的时间细节损失;(2)大多数噪声点应低于长期变化趋势曲线,因为云层和恶劣的大气条件引起的噪声具有负面的偏差。基于这两个标准,利用Savitzky-Golay滤波平滑NDVI变化并获得长期变化趋势。由Savitzky-Golay滤波算法的特点可以发现,m(平滑窗口的半宽度)的设置太小可能会“过度拟合”数据点,并导致捕获长期变化趋势的困难,而m值设置太大可能忽略NDVI时间序列中的一些重要变化。因此,在4-7范围内的m的中间值可以被认为是产生长期变化趋势曲线的合适参数。 因此,当d通常设定在2到4的范围内时,存在m和d(多项式度)的4times;3组合。考虑到NDVI时间序列中时间细节对物候研究的重要性,自动选择使用所有组合的最小二乘法拟合方法,对应于最佳拟合的m和d的组合作为生成长期的最佳参数合成变化趋势曲线。根据我们的实验结果,如下所示,这种类型的参数选择可以在保留NDVI时间序列的时间细节和识别噪声点之间提供良好的权衡,图2b显示了使用Savitzky-Golay滤波获得的NDVI时间序列(ti,)的长期变化趋势曲线。很明显,长期变化趋势曲线保留了NDVI时间序列的时间细节,所有突然下降都被看作低于长期变化趋势曲线。
2.2.3.第3步:确定NDVI时间序列中每个点的权重值
在NDVI时间序列(ti,)经过平滑后,得到表示NDVI长期变化趋势的新时间序列为(ti,)。 这个新的时间序列用于通过将其与时间序列(ti,)进行比较来确定NDVI
剩余内容已隐藏,支付完成后下载完整资料
资料编号:[27349],资料为PDF文档或Word文档,PDF文档可免费转换为Word
以上是毕业论文外文翻译,课题毕业论文、任务书、文献综述、开题报告、程序设计、图纸设计等资料可联系客服协助查找。