季烨云，刘鹏举，高 影，任 怡，余 涛，牛思圆

（1. 中国林业科学研究院 资源信息研究所，北京 100091；2. 中国林业科学研究院 国家林业和草原局林业遥感与信息技术重点实验室，北京 100091；3. 国家林业和草原局 调查规划设计院，北京 100029）

植被是生态系统的重要组成部分，对全球气候变化和生态系统的结构功能都有显著影响[1-3]。研究地表植被覆盖及其长期变化不仅可以揭示地表空间变化规律，而且有利于分析区域生态系统及环境变化情况[4]，对因地制宜开展植被的保护管理工作具有重要意义。黄土高原位于黄河中游地区，是世界上水土流失最为严重的区域之一。陕北是黄土高原的中心部分，水土流失严重、生态环境脆弱，但自三北防护林、天然林资源保护、退耕还林还草等工程实施以来，该区域生态环境得到了明显改善。因此有必要开展工程实施效果阶段性评估，便于为今后该区域的林业生态建设和保护提供理论依据。

时间序列的遥感数据可以动态监测大面积的植被生长状况[5-6]、发展变化趋势特征[7-8]、扰动情况[9-10]等。现有的大多数研究中，对归一化植被指数(normalized difference vegetation index，NDVI)时间序列的分析多针对空间和时间总体变化展开[11]，如STOW等[12]提出用一元线性回归方法来分析每个栅格像元多年的变化趋势，以此来计算植被的绿度变化率(greenness rate of change，GRC)。当植被遭遇砍伐、病虫害、火灾等自然灾害或人类活动影响时，植被绿度变化趋势可能会出现突然的转变。然而，一元线性回归方法只能表现植被变化的整体趋势，不能监测到以上趋势发生突然变化的情况，具有很大的局限性[13-14]。因此，学者们提出了专用于处理遥感时序数据的断点检测方法，如LandTrendr(Landsat based detection of Trends in Disturbance and Recovery)方法[15]、BFAST(breaks for additive seasonal and trend)方法[16]等。由于上述方法受限于卫星传感器、数据类型和数据序列长度等方面的不足，JAMALI等[17]提出了DBEST(detecting breakpoints and estimating segments in trend)算法，该算法不仅对卫星传感器类型和数据序列长度没有限制，而且支持周期性和非周期性植被指数时间序列，检测效果较好且具有普遍的适用性[18-20]。

本研究以延安市富县为例，利用2000—2019年Landsat影像构建NDVI时间序列，运用DBEST算法识别NDVI时间序列的断点，分析断点的时空分布，并研究植被变化多阶段趋势特征，为生态修复工程规划和长期生态环境保护决策提供依据。

1 研究区与数据

1.1 研究区概况

富县位于陕西省北部，延安市南部，35°44′06″～36°23′23″N，108°29′30″～109°42′54″E，总面积为4 182 km2，属于渭北黄土高原丘陵沟壑地带。该县属大陆性暖温带季风气候，四季分明，年平均气温为7～9 ℃，年平均降水量为500～600 mm。区域内分布着丰富的森林资源，乔木树种主要有栎类Quercus spp.、桦木Betula spp.、刺槐Robinia pseudoacacia等；灌木主要有胡颓子Elaeagnus pungens、马桑Coriaria nepalensis、白刺花Sophora davidii等。居民地、河流、铁路和公路分布概况如图1所示。

图1 2019 年富县基础地理信息示意图Figure 1 Basic geographic information of Fu County in 2019

1.2 数据来源及预处理

Landsat TM/ETM+/OLI多光谱遥感影像来源于美国地质调查局提供的地表反射率(surface reflectance，SR)产品，即L2级别产品(https://earthexplorer.usgs.gov/)，空间分辨率为30 m，时间分辨率为16 d，光谱范围包括可见光、近红外、短波红外和热红外波段[21]。由于2003年之后的Landsat 7 ETM+影像出现数据条带丢失的问题，因此需先利用ENVI软件中的Landsat_gapfill插件去除此类影像中的条带[22]。然后利用Fmask云检测算法对所有影像进行云、云影识别提取[23]。最后利用红色和近红外波段计算NDVI，并基于最大值合成法[24](maximum value composite，MVC)合成了2000—2019年的NDVI时间序列数据。

本研究还收集了研究区2019年的林地类型，数据来源于国家森林资源智慧管理平台中的“森林资源一张图”数据，用于分析不同林地的植被变化情况。各类型面积如图2所示：有林地2 664.77 km2，疏林地 82.80 km2，灌木林地 554.95 km2，未成林造林地 253.85 km2，其他用地 625.63 km2，包括人造地表、水体、耕地等。

图2 2019 年富县林地类型与分布示意图Figure 2 Forest land types and distribution of Fu County in 2019

2 研究方法

2.1 断点检测方法

DBEST算法主要是由时间序列分解和趋势分割两部分组成[17]。时间序列分解主要是指植被指数时间序列数据的预处理过程，而趋势分割是实现断点检测的过程。由于本研究对生长季NDVI数据采用最大值生成年度NDVI数据，因此时间序列分解部分未进行以局部移动平均为基础的季节趋势分解(seasonal-trend decomposition procedure based on Loess，STL)，仅利用水平位移第 1 阈值、水平位移第2阈值、持续时长共3个参数检测候选断点。趋势分割(即断点检测)包括3个过程：①找到时间序列的转点(即峰/谷点、距离阈值点和水平位移点)；②通过下一个转点减去当前转点处的NDVI值来计算局部趋势变化(trend local change，TLC)；③所有转点按TLC降序重新排列，利用贝叶斯信息准则(bayesian information criterion，BIC)确定断点数量，若BIC准则得到的断点数量少于转点数量，则最终的断点数为BIC准则得到的结果。

运用DBEST算法检测断点，设置参数如下：水平位移第1阈值(为0.1，即当时序NDVI数据变化超过0.1时认为该点是候选水平位移点；水平位移第2阈值)为0.2，即当候选水平位移点前后最小步长范围内的子序列均值变化超过0.2认为候选水平位移点符合第2阈值；持续时长()为2 a，在该时段内只能有1个水平位移点，以检测NDVI短时间内的平均变化；以上3个参数可确定时间序列中的水平位移点。变化级别()、泛化百分比()、断点数目(m)共3个参数可根据需求三选一达到不同的目的。由于富县植被覆盖率较高，NDVI位于-0.2～0.2之间的变化在原始影像上的目视解译不太明显，因此本研究将设为0.2，可保留大于0.2的断点，以便于检测较为显著的变化。算法输出变化开始时间(即断点位置)、变化结束时间、变化幅度和变化持续时长。

2.2 NDVI趋势分析

2.2.1 总体趋势 本研究利用一元线性回归方法来拟合每个像元年内NDVI最大值(INDV(max))的年际变化趋势[24]：

其中：i为年序号；n为总年数；INDV(max)为年内NDVI最大值；为第i年的INDV(max)；为NDVI年内最大值的年际变化趋势。将基于一元线性回归方法得到的富县2000—2019年NDVI时间序列的趋势定义为总体趋势。依据宋怡等[25]提出的分类标准将趋势划分为3种类型：退化(≤-0.001 0)、基本不变 (- 0.001 0＜SNDVI≤0.001 0)、改善 (SNDVI＞0.0010)。

2.2.2 多阶段趋势 ① 对 NDVI时间序列进行断点检测，分析断点及 NDVI变化幅度的时空分布。② 针对未检测到断点或NDVI变化幅度小于0.2的区域，研究其总体趋势；针对断点NDVI变化幅度大于0.2的区域，对其进行多阶段趋势分析。多阶段趋势分析过程如下：以NDVI时间序列的所有断点及对应的变化结束时间对其进行分段，分别进行拟合，从而得到趋势序列。断点及对应的变化结束时间生成的拟合段定义为变化段，除了变化段以外的拟合段定义为趋势段。当前阶段趋势定义为最末一段趋势段。趋势类型按2.2.1中标准分类。对于多趋势段区域，研究其当前阶段趋势，以分析植被的现状。对于持续时长小于5 a的趋势段，由于时间太短，不进行趋势分析。③ 根据林地类型分析林地总体趋势和当前阶段趋势的变化情况。

3 结果与分析

3.1 断点检测与分析

应用DBEST算法对富县NDVI时间序列数据进行断点检测，时间序列中第1断点NDVI变化幅度的空间分布与不同幅度面积占比如图3所示。发生正向变化的区域占全县面积的79.34%；负向变化的区域占20.66%，多位于道路两侧及居民点附近(图1)。NDVI变化幅度为-0.2～0.2的区域占总面积66.47%，主要位于西部及东北地区，植被发生的变化较小；变化幅度绝对值超过0.2的区域占总面积的33.53%，位于富县的东部及河流道路周围区域，属于NDVI剧烈变化的区域。

图3 富县归一化植被指数变化幅度及其面积占比Figure 3 Change magnitude and area percent of NDVI in Fu County

图4显示：2000—2019年间富县NDVI变化幅度大于0.2的区域断点个数的空间分布及各类型的面积占比。其中，66.47%的区域断点NDVI变化幅度小于0.2，主要分布在富县西部和东北地区，说明NDVI变化幅度小，植被相对稳定，因此不作进一步分析。断点NDVI变化幅度大于0.2的区域占总面积的33.53%，其中断点数量为1或2个的区域占总面积的21.59%，位于东部地区；断点数量多于4个的区域，仅占5.88%，集中分布于道路、河流沿线，植被变化频繁，与人为活动有关。

图4 富县不同 NDVI变化幅度小于 0.2 与大于0.2的断点数目及其面积占比Figure 4 Number of breakpoints and area percent of NDVI (＜0.2 and＞0.2) in Fu County

3.2 NDVI时间序列趋势分析

3.2.1 总体趋势 富县 NDVI 总体趋势有很强的空间差异。对趋势值进行分类，得到不同类别的分布及面积(图5)。20 a间富县超过97%的区域植被呈改善趋势；退化区域占比不足2%，植被改善的面积远大于退化的面积。退化区域集中分布在富县东部区域、道路及河流沿线，在西北部散落分布，其中东部部分人造地表退化面积较大；基本不变的区域面积最小，不足1%，分布较零散。

图5 2000—2019 年富县 NDVI 时间序列总体趋势类别及其面积占比Figure 5 Overall trend level and area percent of NDVI time series in Fu County from 2000 to 2019

3.2.2 多 阶 段 趋 势 基 于 3.1 断 点 检 测 结 果 ，对2000—2019年富县NDVI时间序列进行多阶段趋势分析。由于66.47%的区域无断点，因此该区域只分析总体趋势；33.53%的区域具有1个及以上断点，对该区域时间序列进行分段，分析多阶段趋势中的当前阶段趋势。当前阶段趋势统计结果及分布(图6A)表明：27.57%的区域当前阶段趋势持续时长大于5 a，植被当前阶段改善的面积远大于退化的面积，改善面积占24.54%，分布于东部区域及道路河流沿线的周边区域；退化面积仅占2.12%，基本不变面积不足退化面积的一半，分布较零散。当前阶段趋势的开始时间(图6B)均在2014年之前，时间分布跨度大；开始于2002年的区域面积最大，占12.76%，且均为改善趋势，表明经历一次变化以后，区域植被一直呈现改善状态；其次是开始于2014年的区域，占3.68%，主要位于道路、河流沿线；其余年份所占面积均不足2%。

图6 2000—2019年富县 NDVI时间序列当前阶段趋势及其面积占比Figure 6 Current stage trend and area percent of NDVI time series in Fu County from 2000 to 2019

3.3 不同林地趋势分析

如图2所示：富县植被退化的区域大部分分布在西北、西南、东北以及中部的一些地区，如今中部地区以及一些道路沿线已经变为永久性的建设用地；其他大部分均已改善，大多变为有林地和疏林地。2019年富县林地类型中，有林地占63.72%，疏林地占1.98%，灌木林地占13.27%，未成林造林地占6.07%，4种林地占富县总面积的85.04%。如表1所示：有林地、疏林地、灌木林地、未成林造林地的改善面积占其总面积的比例分别为96.98%、94.32%、94.41%、92.97%，其他用地改善和基本不变的比例分别为86.66%和2.50%，退化的比例达10.84%，由于其他用地包括建设用地，说明大部分退化与城市和道路建设相关。整体来说，富县生态建设工程取得了良好的效果，大部分林地处于改善的趋势，特别是有林地和灌木林地。

表1 富县林地类型趋势情况统计Table 1 Statistics on trend of forest land types in Fu County

4 结论

本研究基于断点检测的多阶段趋势分析方法，研究了陕西省富县近20 a植被的时空变化特征，得到以下结论：①第1断点检测结果表明：富县发生正向变化的区域面积大，占79.34%。检测结果中NDVI变化幅度小于0.2的区域占66.47%，主要分布在富县西部和东北地区，植被相对稳定；大于0.2的区域占33.53%，其中断点数量为1或2个的区域比例为21.59%，位于东部地区；断点数量多于4个的区域，仅占5.88%，集中分布于道路、河流沿线，变化频繁，与人为活动有关。②植被的多阶段趋势表明：富县共有27.57%的区域当前阶段趋势持续时长大于5 a，当前阶段趋势为改善的面积占24.54%，分布于东部区域及道路河流沿线的周边区域；为退化的面积占2.12%，基本不变面积仅占0.91%，分布均较零散。当前阶段趋势的开始时间均在2014年之前，时间分布跨度大，空间异质性强，揭示了富县植被变化的多样性及复杂性。③近年来富县有林地、疏林地、灌木林地、未成林造林地均是改善情况明显，其中有林地和灌木林地中的改善比例最大。林地中呈改善态势的面积占全县总面积的81.83%，富县生态建设工程取得了良好的效果。同时对于检测到退化区域的林地，需要进一步分析，采取有效措施提高质量。

1999年以来，富县积极响应国家号召，实施退耕还林、荒山荒地造林、封山育林等措施。这些生态工程建设促进了植被的恢复，然而在不同地区的效果是存在差异的。本研究获得的NDVI的变化直接反映了植被恢复的效果，但DBEST算法参数设置对分析结果的影响、趋势变化与影响因素的关系尚不清楚，在后续研究中，需从多尺度应用DBEST算法，探究其应用技术，并综合考虑气温、降水、人为因素等对趋势变化的影响。