塞罕坝植被覆盖时空分异与驱动机制研究*

王晓红,辛守英,阳丽虹,马明浩,焦琳琳

(1.华北理工大学 矿业工程学院,河北 唐山 063210;2.河北省矿区生态修复产业技术研究院,河北 唐山 063210;3.河北省矿山绿色智能开采技术创新中心,河北 唐山 063210;4.国家知识产权局专利局专利审查协作广东中心,广东 广州 510700)

植被是陆地生态系统的重要组成部分,在气候调节、水土保持和生态系统稳定方面发挥着关键作用[1]。近些年,气候变化和人类活动对植被的影响不断加剧,而植被对气候变化和人类活动的响应也变得更加敏感[2]。因此,监测植被覆盖的时空分布变化以及了解气候因素和人类活动的驱动机制,对改善生态环境和提升生态服务功能具有重要的意义。

遥感技术以其覆盖范围广、时空连续性强的优势,成为常用于监测植被覆盖分布变化的重要工具[3]。归一化差异植被指数(normalized difference vegetation index,NDVI)作为一种高分辨率且具有良好时效性的遥感数据,被广泛用于区域植被覆盖变化的动态监测[4]。如贾云飞等[5]利用2000—2019年的NDVI数据对黄土高原延河流域及其沟壑区植被覆盖变化进行分析。吴晶晶等[6]利用长时序NDVI数据,对祁连山地区的生态修复前后植被覆盖时空变化进行分析。同时,针对NDVI对植被覆盖变化的驱动机制的研究,如于璐等[7]运用回归分析和残差分析法,分离出了不同时期京津风沙源气候变化和人类活动对植被变化的影响。石淞等[8]运用趋势分析、Hurst指数、偏相关分析、残差分析及相对作用分析法,深入解析黑龙江省植被时空变化特征及其对气候变化和人类活动的响应机制。由此,植被覆盖时空变化主要受气候和人类活动两大因素的影响[9]。其中,气温和降水是能够量化影响植被覆盖变化的主要气候因素[10]。而人类活动因素对植被覆盖变化的影响较难直接定量,主要通过获取人类活动的相对贡献率来评估,目前定量区分人类活动相对贡献率应用最广泛的是残差趋势法[11],其能够从气候因素中有效地分离出人类活动对植被覆盖变化的影响。此外,对植被覆盖的NDVI进行植被覆盖等级划分,这有助于揭示在不同干扰情况下区域植被的动态变化特征,并更好地探索植被覆盖的动态变化的驱动机制[12]。

塞罕坝机械林场作为京津冀乃至华北地区防风固沙和涵养水源的重要生态屏障,是集森林、草原、草甸、水体、湿地等为一体的复杂生态系统,其植被覆盖变化多受自然气候环境和人类活动等的综合影响[13]。然而,以往对塞罕坝的研究主要关注植被景观格局时空尺度效应,对气候变化和人类活动与不同等级植被覆盖区相互作用的研究相对较少[14]。基于此,本文采用气温、降水和MODIS-NDVI数据,通过趋势分析、偏相关性分析、残差分析和相对作用分析等方法,逐像元地分析2000—2020年塞罕坝林场的生长季植被NDVI的空间分布,并进行植被覆盖等级划分,在不同植被覆盖等级的基础上分析气候因素和人类活动的相对贡献率,揭示气候变化和人类活动与植被覆盖变化的相互作用机制,阐释塞罕坝植被覆盖现状特征和趋势变化驱动机制,为塞罕坝植被保护与森林经营管理提供科学的理论指导。

1 材料与方法

1.1 研究区概况

塞罕坝机械林场(42°02′～42°36′ N、116°51′～117°39′ E)位于河北省承德市围场满族蒙古族自治县的最北部(以下称为塞罕坝林场),如图1所示。该林场地处内蒙古高原和冀北山地的接壤地带,是省属的大型国有人工林林场[15],总面积926.347 km2,森林覆盖率为75.5 %[16];总场下有三道河口、千层板、北曼甸、阴河、第三乡和大唤起6个分场。林场属于寒温带半干旱半湿润大陆性季风气候[17],年均气温-1.4℃,年降水量490 mm;冬季受内蒙古高压控制,盛行西北风且风力大,气候寒冷;夏季受太平洋气压影响,多东南风,风力小,降雨频繁。塞罕坝林场地形复杂,分为坝上和坝下两部分。其中,北部地区属于坝上内蒙古高原的一部分,土壤类型主要为风沙土,伴随有草甸土和沼泽土;中部和南部地区属于坝下山地,土壤类型以棕壤和灰色森林土为主[14]。

图1 研究区位置与地形

1.2 数据来源与预处理

1.2.1 MODIS-NDVI数据

本文采用的是MOD13Q1 NDVI(https://ladsweb.modaps.eosdis.nasa.gov)遥感数据,时间跨度为2000—2020年;时间分辨率是16 d,空间分辨率为250 m,行列号为h26v04。原始数据通过MRT(MODIS Reprojection Tools)工具批量进行格式转换、投影转换、裁剪等操作,即将HDF格式转换为GeoTIFF,坐标投影为WGS84坐标系,裁剪为塞罕坝林场区域;并使用ArcGIS软件剔除NDVI数据中的无效值,留下NDVI真值。为了消除云层、大气、卫星轨道偏移与太阳高度角等的影响,利用最大值合成(MVC)方法合成月尺度数据。根据以往研究,将每年的4—10月作为生长季[18]。为了更好地反映研究区植被覆盖生长状况,计算生长季的NDVI均值进行研究[19]。

1.2.2 气象数据

本文使用的气温和降水气象数据来源于国家地球系统科学数据中心——黄土高原分中心(http://loess.geodata.cn),在ArcGIS软件中根据塞罕坝林场矢量边界裁剪获得2000—2020年分辨率为1 km逐月平均气温和降水数据,并计算每年的生长季平均气温和降水数据;利用ArcGIS软件对气气温和降水数据重采样为与MODIS-NDVI数据像元大小和位置相同的栅格数据。

1.3 研究方法

1.3.1 趋势分析

本研究使用SPSS软件对NDVI栅格数据逐像元进行一元线性回归分析2000—2020年塞罕坝林场的时空变化趋势。其具体计算公式如下:

①

式中:Slope为NDVI变化趋势斜率即NDVI年际变化速率;i为研究时段内第i个年份;n为研究时段的年数(n=21);NDVIi为第i年的生长季NDVI均值。根据回归特性,当Slope<0,表示NDVI随时间减少,即区域植被覆盖呈减少趋势;当Slope>0,表示NDVI随时间增加,即区域植被覆盖呈增加趋势;总之Slope的绝对值越大,NDVI随时间变化越快,即区域植被覆盖变化趋势越快。此外,通过F显著性检验方法检验NDVI的变化趋势的可置信程度,检验的具体计算公式如下:

②

1.3.2 偏相关性分析

本研究利用偏相关分析是逐像元计算2000—2020年塞罕坝林场年均NDVI与气温和降水的气候因子的偏相关系数。其具体计算公式如下:

③

式中:rxy,z为在控制变量z的情况下x和y之间的偏相关系数,即xy相关分析中剔除z的影响;rxy、rxz和ryz为不同变量间的相关系数;x、y、z为NDVI与气温与和降水的气候因子。其不同变量间的相关系数具体计算公式如下:

④

⑤

式中:r为NDVI与气温和降水的气候因子的偏相关系数;n为研究时段的年数;k为自变量的个数。

1.3.3 残差分析

残差分析是区分气候变化和人类活动对NDVI的影响[21]。该方法通过对NDVI与气候因素的气温和降水逐像元进行多元线性回归分析,拟合得到NDVI的预测值,将其视为气候因素对NDVI的影响,进而利用NDVI的预测值与从遥感影像获取的NDVI观测值做差,差值即为残差,是表示人类活动对NDVI的影响。其具体计算公式如下:

NDVICC=a×T+b×P+c

⑥

NDVIHA=NDVIobs-NDVICC

⑦

式中:NDVICC为多元线性回归模型中NDVI预测值;NDVIobs为从遥感影像获取的NDVI观测值;NDVIHA为残差;T和P分别为平均气温和降水;a和b分别为回归系数,c为回归方程的常数项。

1.3.4 相对作用分析

Sun等应用相对作用分析方法确定了气候因素和人类活动在华北地区植被NDVI变化过程中的相对贡献率[22]。该方法是基于一元线性回归趋势分析和残差分析得到NDVICC和NDVIHA的变化趋势,在不同情景下,气候因素和人类活动对生长季NDVI相对贡献率的计算方式不同[23],如表1所示。情景1和4表示气候因素和人类活动的共同作用导致的NDVI变化,情景2和5表示仅气候因素导致的NDVI变化,情景3和6表示仅人类活动导致的NDVI变化。

表1 不同情景下植被NDVI变化的驱动因素相对贡献率计算

2 结果与分析

2.1 2000—2020植被覆盖分级及其空间分布

根据塞罕坝林场植被覆盖特征,本文将2000—2020年生长季NDVI均值划分为3个等级。具体划分标准为:低于0.55的区域被归为低等植被覆盖区,高于0.65的区域归为高等植被覆盖区,而0.55～0.65的区域则划分为中等植被覆盖区;如图2展示了这些不同等级植被覆盖区在塞罕坝林场内的空间分布。

图2 不同植被覆盖等级空间分布

根据图2可知,低等植被覆盖区主要分布在三道河口林场的边缘地带以及少量分布在千层板、北曼甸、阴河和大唤起分林场;中等植被覆盖区主要分布在西部的三道河口和千层板、东部的北曼甸和阴河林场、南部的第三乡和大唤起分林场;高等植被覆盖区主要分布在南部的第三乡和大唤起分林场以及少量分布在其他分林场。因此,塞罕坝林场植被NDVI覆盖具有明显的东高西低、南高北低的空间分布特征。这种空间分布受多方面原因的综合影响:一方面受到土壤类型的影响,研究区东部主要为黑土,为冷湿草原草甸植被下形成的土壤,土壤较湿润,有机质含量多,土层深厚,适宜植被生长;西部土壤类型主要是风沙土,土壤养分低,植被生长条件严苛[24]。另一方面水热条件也在植被覆盖变化中发挥重要作用,东部降水明显多于西部,东南部气温高于东北部[25],导致了研究区植被覆盖东高西低、南高北低的分布特点。此外,根据植被覆盖分级的像元空间分布占比可知,低等植被覆盖区面积远远小于中高等植被覆盖区,中高等植被覆盖率高达94.62%,而低等植被覆盖率仅为5.38%,如表2所示。

表2 不同植被覆盖等级像元占比

2.2 2000—2020年植被NDVI的时空变化分析

本文对2000—2020年的塞罕坝林场植被NDVI覆盖等级化数据逐像元进行趋势分析,并通过F检验的分析结果,如图3所示。

图3 低(a)、中(b)、高(c)植被覆盖区NDVI变化趋势

根据图3可知,从不同植被覆盖等级的NDVI变化速率来看,低等植被覆盖区NDVI增长速率较快(图3a),平均增长速率为0.005 6/a,中等植被覆盖区NDVI增长速率存在着空间差异(图3b),西部NDVI增长速率较快,而中部和东部增长速率较缓,平均增长速率为0.004 3/a。这可能是由于西部早期植被覆盖度不高,地势平坦,造林工作主要致力于西部沙地[26]。高等植被覆盖区NDVI增长速率较缓(图3c),平均增长速率为0.003 1/a。总之,NDVI变化速率与植被覆盖情况大致相反,高植被覆盖区NDVI增长速率慢,低植被覆盖区NDVI增长较快。此外,通过2000—2020年逐年NDVI植被覆盖等级化的数据逐像元进行趋势分析,并通过F显著性检验的结果如图4所示。

图4 低(a)、中(b)、高(c)植被覆盖区归一化差异植被指数年平均变化

根据图4可知,从塞罕坝林场逐年平均植被NDVI的变化情况来看,生长季年均NDVI呈显著波动上升趋势(Slope>0,P<0.05)。近20年以来,塞罕坝林场植被NDVI增长趋势可分为2个阶段,2000—2007年为显著上升阶段(Slope>0,0.010,P<0.01),即2007年之后,塞罕坝林场全域植被NDVI平均增长速率明显高于2000—2007的植被NDVI平均增长速率。从不同植被覆盖等级的植被的变化速率来看,低植被覆盖区NDVI增长速率快(图4a),平均增长速率为0.005 1/a,中植被覆盖区NDVI平均增长速率为0.003 3/a(图4b),高植被覆盖区NDVI增长较缓(图4c),平均增长速率为0.002 7/a。

2.3 2000—2020植被NDVI变化的驱动机制分析

2.3.1 气候因素对植被NDVI的影响

通过逐像元计算2000—2020年的生长季平均NDVI与平均气温和降水的气候因素的偏相关系数并得到通过F显著性检验的结果,如图5所示。

图5 植被NDVI与气温(a)和降水(b)的偏高相关系数空间分布

根据图5可知,植被NDVI与气温和降水的相关性主要呈现正相关,且降水与NDVI的相关性更为显著,其偏相关系数大都大于0.5。此外,在不同植被覆盖区统计了NDVI与气温和降水的平均偏相关系数,如表3所示。

表3 不同植被覆盖区NDVI与气温和降水的平均偏相关系数

根据表3可知,随着植被覆盖程度的改善,NDVI与气温的相关性越高,与降水的相关性越低。此外,从不同植被覆盖等级方面,评估气候因素对塞罕坝林场植被NDVI覆盖变化的相对贡献率,如图6所示。

图6 低(a)、中(b)、高(c)植被覆盖区气候因素的相对贡献率

根据图6可知,低植被覆盖区气候因素对NDVI的相对贡献率多低于50%(图6a),其平均相对贡献率约为39.22%。中高等植被覆盖区东北部气候因素的相对贡献率相对较高(图6b和6c),其相对贡献率均在50%以上。这主要是因为塞罕坝林场东北部属坝上地区,海拔较高,人类活动较少。

2.3.2 人类活动对植被NDVI的影响

从不同植被覆盖等级方面,评估人类活动对塞罕坝林场植被NDVI覆盖变化的相对贡献率,如图7所示。

图7 低(a)、中(b)、高(c)植被覆盖区人类活动的相对贡献率

由图7可知,低植被覆盖区的人类活动对NDVI相对贡献率多在50%～75%之间(图7a),其平均相对贡献率约为60.78%。中高等植被覆盖区的西部和东南部人类活动的相对贡献率相对较高(图7b和7c),其相对贡献率均在50%以上。

3 讨论与结论

3.1 讨论

自2000年以来,塞罕坝林场的不同植被覆盖区的NDVI呈上升趋势,其中低等植被覆盖区的NDVI逐年平均增长速率为0.005 1/a,中等植被覆盖区的NDVI逐年平均增长速率为0.003 3/a以及高等植被覆盖区的NDVI逐年平均增长速率为0.002 7/a。与贺军亮等对承德市全域的2000—2018年植被的NDVI整体呈上升趋势研究结果[27]和王倩等对张家口和承德区域的2001—2020年植被的NDVI呈增加趋势研究结果[28]一致。此外,2000—2020年,塞罕坝林场植被NDVI增长趋势由2007年转折分成2个阶段,2000—2007年为显著上升阶段;2008—2020年为极显著快速上升阶段,这很大程度上是与2007年塞罕坝林场晋升成了国家级自然保护区、进入了塞罕坝林场国家级自然保护区一期总体规划项目实施有关[29]。

根据偏相关性分析可知,塞罕坝林场植被的NDVI时空变化是受气候因素和人类活动的共同影响。在气候因素中,降水和温度是与植被覆盖增长直接相关的气候因素,低等植被覆盖区的气候因素对NDVI的平均相对贡献率约为39.22%,中高等植被覆盖区东北部气候因素的相对贡献率相对较高,这主要因为塞罕坝林场东北部属坝上地区,海拔较高,人类活动相对较少。赵明伟等[30]发现,降水是驱动华北平原北部和内蒙古大部分区域植被覆盖时空动态变化的重要因素。邵雅琪等[31]研究结果表明,NDVI与降水呈显著正相关,与气温呈显著负相关。而本研究中塞罕坝林场植被的NDVI与降水和气温的相关性均呈正相关,且降水与NDVI的相关性更为显著,即降水对植被的NDVI的影响大于气温,这可能是因为塞罕坝林场气温空间分布差异小、降水空间分布差异大所致。本研究在气候因子中仅考虑了降水和气温,并未考虑相对湿度和极端气候等其他气候因素的影响,因此,在今后的研究中应考虑该类因素影响使研究内容更加准确充实。研究期间,人类活动在低植被覆盖区的NDVI平均相对贡献率约为60.78%。中高等植被覆盖区的西部和东南部人类活动相对贡献率相对较高,这是因为西部地势平坦,是塞罕坝林场造林工程的核心区域,东南部海拔较低,属于坝下山地区,易受到人类活动的影响。

3.2 结论

本文基于2000—2020年长时序NDVI数据划分植被覆盖等级进行塞罕坝林场植被覆盖变化监测,同时以气温和降水数据为基础分析不同植被覆盖等级的植被覆盖变化受气候和人类活动的主要驱动机制,得如下结论:

根据塞罕坝林场植被覆盖特征,将林场内NDVI低于0.55的区域归类为低植被覆盖区,NDVI高于0.65的区域归类为高植被覆盖区,NDVI为0.55～0.65的区域为中等植被覆盖区。在不同等级植被覆盖等级区,植被NDVI的增长速率有所差异,低等植被覆盖区NDVI增长速率较快,其次是中等植被覆盖区和高等植被覆盖区。尤其是2007年之后,塞罕坝林场全域植被NDVI平均增长速率明显高于2000—2007的植被NDVI平均增长速率。此外,在不同等级植被覆盖等级区,NDVI的主要驱动因素也有所不同,低等植被覆盖区NDVI主要受到人类活动的影响,中高等植被覆盖区主要驱动因子存在着空间异质性,西部和东南部低海拔区的人类活动是主要的驱动因素,东北部高海拔区植被覆盖变化则主要受到气候因素的影响。