近40年黄河流域季节性植被覆盖时空动态及驱动因子探测

任天晨，陈军锋，刘 楠，范 强

（太原理工大学水利科学与工程学院，山西 太原 030024）

0 引 言

植被作为生物圈生态可持续的重要评价因子，其生长变化不仅指示了气候变化如何影响生物圈，同时也可以通过植被蒸腾作用、地表反照率和粗糙度影响陆地和大气之间能量交换，对区域防止水土流失、调节气候等有重要影响［1-3］。然而随着气候持续变暖和人类活动增强，对植被生长产生了显著影响［2］。因此，全球已经实施了许多与植被恢复相关的项目来改善植被生长状况［2-4］。黄河流域作为黄河水系从源头到入海所影响的重要地理生态区域和我国灌溉农业发达地区，直接养活着1.07亿人口，另外还有4 亿人口依赖于黄河流域的水生活［4］。但黄河流域生态本底极为敏感脆弱，为了恢复黄河流域生态环境，自1999 年我国实施了一系列生态恢复工程并取得了阶段性成果［4-6］。但有研究指出在生态工程实施过程中没有综合考虑区域气候、人文等因子，目前高植被覆盖区域的植被（主要为林地）已经达到区域可容纳的最高阈值而呈退化趋势［7-9］，这可能会影响黄河流域生态可持续发展。因此，为了黄河流域植被恢复的长期成功，需要明确黄土高原植被覆盖变化的转折点及不同季节植被覆盖变化的主要驱动因子，为黄河流域的生态恢复工程效应的评价提供一定参考。

归一化植被指数（Normalized Difference Vegetation Index，NDVI）是目前最常用的反映大尺度地表植被状况和变化的指数［8-11］。曾有研究发现20世纪80年代以来，全球植被生长季内平均NDVI和年最大NDVI均呈增加趋势，植被覆盖增长速率、幅度均以显著增加趋势为主，这一趋势在植树造林的黄河流域最为显著［12，13］。其中，有研究发现植被覆盖增加的主要驱动力是二氧化碳施肥效应，其次为碳沉降、气候变化和人类活动，黄河流域植被变绿关键因素主要为植树造林［14-16］。然而有研究也发现随着黄土高原植被覆盖增加，植被蒸腾作用会加强水循环过程，从而影响区域土壤水分、地表径流和地表干燥化程度等，反过来抑制黄河流域植被生长［9，11］。如Fu 等［14］通过对黄土高原研究发现通过大规模的生态恢复工程减少了黄河泥沙量，生态环境也得到改善，但由于前期的植被恢复工程缺乏科学的指导，部分地区植被恢复不合理规划导致土壤水分缺失、径流减少和生态系统退化［9，14］。同时，有研究也发现新种植的植被使蒸发量增加，导致土壤水分减少从而造成土壤表层干旱化程度加剧［15-18］，这将进一步影响黄土高原植被覆盖的可持续发展。Xin 等［19］研究也发现1982-1999 年黄土高原西北部低覆盖度地区植被覆盖呈上升趋势，而部分高覆盖地区呈下降趋势，这可能是该区域植被覆盖已经达到该区域环境的承载能力从而趋向于饱和。其次，由于全球气候变暖加剧了气候系统的不稳定性，导致全球极端气候事件出现的频率和强度均呈现明显的增加趋势，这对黄河流域植被覆盖的可持续生长产生了显著影响［15，20］。如Jiang 等［21］发现在气候变暖背景下极端干旱和极端湿润交替出现对植被生产力的影响具有1 到5 年的记忆效应，弥补了过去的全球尺度研究强调极端干旱之后植被生产力下降得结论。因此，极端气候事件对黄土高原植被生长造成的影响不容忽视。虽然目前针对人类活动、平均气候和极端气候等因素对植被生长产生的影响已有很多研究，但大部分研究是基于年尺度进行，很少有研究综合分析不同季节植被生长对人类活动和气候变化的敏感性程度。

基于此，基于GIMMS NDVI 和MODIS NDVI 数据探究了1982-2020 年间黄河流域不同季节植被覆盖时空动态，进一步考虑到黄河流域部分地区水蚀和风蚀等影响严重并结合其他学者对黄河流域的研究（很多学者研究发现在黄河流域开展的生态恢复项目由于没有考虑到植被蒸散发等作用，导致区域蒸散发增加，这使得地表干燥化加重。因此，本文探究了气温（日最高气温最小值、日最低气温最大值、日最低气温最小值、年均气温、日最低气温最大值）、降水（每年连续五天降雨量最大值、年均累计降水）、土地利用、风速、海拔、坡度对黄河流域植被生长的影响程度，旨在回答以下问题：①自1982 年以来黄河流域植被覆盖时空动态如何？②不同季节NDVI 对极端气候、平均气候、人类活动和地形因子的响应程度如何，不同因子的相互作用对黄河流域植被生长影响如何？

1 研究区与研究方法

1.1 研究区

黄河流域是我国第二大流域，西起巴颜喀拉山，东至渤海，北临阴山，南到秦岭，流域面积达7.52×105km2。地势西高东低，从西到东横跨青藏高原、内蒙古高原、黄土高原和黄淮海平原 4个地貌单元。地势总体呈西高东低的态势，西部河源地区平均海拔超过4 000 m，中部海拔介于1 000～2 000 m［1］。地貌复杂多样，西部多山，常年积雪，中部为黄土高原地貌，水土流失严重，东部由黄河积平原组成。流域内不同地区气候差异显著，从西到南依次为半湿润气候、半干旱气候和干旱气候［3］。全年日照时间长、辐射强，气温日较差大，年较差小。降水量少而不均，湿度小，蒸发量大，多冰雹、沙暴天气。土地利用类型主要为草地和农用地，其次为未利用地、林地和水域（图1）。

图1 黄河流域海拔和坡度空间分布图Fig.1 Spatial Distribution of Altitude and Slope in the Yellow River Basin

1.2 数据来源与处理

（1）NDVI数据：本文使用了来自西部环境与生态科学数据下载中心（https：∕∕ecocast.arc.nasa.gov∕data∕）提供的1982-2015年15天8 km 的GIMMS NDVI3g数据集和NASA（https：∕∕lpdaacsvc.cr.usgs.gov∕appeears）提供的2001-2020 年的16 天250 m 的MOD13Q1 产品的NDVI数据集用于监测黄河流域植被覆盖变化。为进一步去除NDVI噪声，使用了Savitzky-Golay 滤波对GIMMS NDVI3g 和MODIS NDVI 数据进行时间序列重建，以平滑NDVI因云或其他因素导致的错误峰值［15］。并将MODIS NDVI 采样为8 km 分辨率，然后基于最大值合成法（MVC）将两种NDVI数据合成为月数据。最后用2001-2015 的MODIS NDVI数据对GIMMS NDVI重建后的数据进行建模，基于建模结果将GIMMS NDVI 数据延长至2020 年。并基于MVC 方法将12月、1月、2月的NDVI数据合成为冬季，将3月、4月、5月的NDVI数据合成为春季，将6月、7月、8月的NDVI数据夏季，将9月、10月、11月的NDVI数据合成为秋季。

（2）土地利用数据：1980、1990、2000、2005、2010、2015 和2020 年土地利用数据来源于中国科学院资源环境科学数据平台，本文主要选取一级分类，包括草地、林地、耕地、建设用地、水域和未利用地。

（3）DEM 和坡度数据：90 m 的DEM 数据源于中国科学院资源环境科学数据平台（http：∕∕www.resdc.cn）提供，为与本文NDVI数据分辨率保持一致，将90 m 的DEM 数据重采样为8 km。坡度数据是基于DEM 数据利用ARCGIS 软件提供的计算工具计算得到。

（3）气候数据：本文气候数据包括平均气候数据和极端气候数据。平均气候［气温（TEM）、降水（PER）、风速（WIN）］来源于中国气象数据网的月值数据（http：∕∕data.cma.cn∕data∕），时间跨度为1982-2020 年。极端气候数据源于HadEX3 数据集（www.climdex.org），该数据集采用综合观测资源来量化昼夜温度和降水变化，由29 个气候极端指数组成，这些指数是通过世界温度及降水站数据计算而来的，综合反映了极端温度和降水事件的频率和强度，该数据被应用于极端天气事件研究中［29，30］。根据黄河流域温度和降水的实际情况，本文选择了5个最能反映温度日变化范围［日最高气温最大值（TXX）、日最低气温最大值（TNX）、日最高气温最小值（TXN）和日最低气温最小值（TNN）］、长时降水指标［每年连续五天降雨量最大值（RX5DAY）］来分析黄河流域植被覆盖对极端气候指标的响应情况。将12 月、1 月、2 月定义为冬季，3 月、4 月、5 月定义为春季，6 月、7 月、8 月定义为夏季，9 月、10 月、11 月定义为秋季，最后利用专业气象插值软件ANUSPLINE 将季节性气象数据插值为空间分辨率为8km的栅格数据集。

1.3 研究方法

1.3.1 NDVI突变点检验——Mann-Kendall突变点检验

本文采用Mann-Kendall 方法对1982-2020 年的春、夏、秋和冬季的NDVI进行突变点检验。首先使NDVI时间序列构造一个X秩序列记为Sk，在时间序列为随机的假设下［23］，定义统计量：

其中，E（Sk）和Var（Sk）分别是Sk的均值和方差。UFk形成UF的特征曲线，通过可信度检验来判断NDVI是否有明显的变化趋势。UFk为标准正态分布，在给定显著性水平α上，取显著性水平a=0.05，则U0.05 的临界值为1.96 的绝对值。将UF、UB的曲线和0.05 的显著水平在同一张图上描述出来，UF＞0，表示序列呈上升趋势；反之，表明呈下降趋势，大于或小于0.05 显著的临界线的，表示在95%置信水平上有显著的上升或下降趋势。NDVI突变开始的时间即UF与UB曲线出现相交且交点位于临界线之间，则交点对应的时间为NDVI突变点［23］。

1.3.2 NDVI趋势分析——Sen斜率法

采用Sen 斜率法趋势分析计算1982-2020 年春、夏、秋和冬季的NDVI变化趋势［8］。Sen 斜率法是一种非参数统计的计算方法，对测量误差和离群数据不敏感［23］。具体公式如下：

式中：β为斜率，i和j代表年份，如果β＞0，NDVI时间序列具有正趋势；如果β＜0，NDVI序列具有负趋势。

1.3.3 NDVI对不同因子的敏感性——地理探测器

地理探测器模型由王劲峰等人［24］提出，是将自变量空间分布与潜在因素分布进行比较，适合用于测量空间分层非均质性程度的空间分析方法。本文利用地理探测器度量季节性NDVI对海拔、坡度、土地利用变化、气候因子的敏感性，用q值解释这个程度。公式如下：

式中：q是影响因子对NDVI时空变化的解释力，h是不同驱动因子的分类或分区数据；L为影响因子的样本数量，Nh和N分别是h和整个区域的单元数；δ2h和δ2是h和整个区域的方差。q值越大，代表该因子对NDVI影响程度越大。并使用交互作用探测识别不同因子之间的相互作用，并评估他们的组合效应以观测任何一对因素共同作用是否会增加或降低对NDVI空间分布的解释力。

2 结果分析与讨论

2.1 GIMMS NDVI 与MODIS NDVI 的一致性分析与处理

两种原始NDVI月均数据的趋势和相关性极为显著，两者的相关系数（R2）为0.926 5，均方根误差（RMSE）为0.039 5，但GIMMS NDVI 在NDVI极高值和极低值处明显高于MODIS NDVI［图2（a）］。从MODS NDVI 对GIMMS NDVI 重建后的月均NDVI来看，两者的R2为0.980 3，RMSE为0.019 9，在NDVI极值处也较为接近［图2（b）］。从像元尺度来看，两种数据间的相关性也较好［R2= 0.917 1，RMSE=0.030 5，图2（c）］，利用MODIS NDVI 数据对GIMMS NDVI 重建后的数据的相关性显著提高［R2= 0.950 5，RMSE=0.014 9，图2（d）］。因此，可以使用MODIS NDVI数据对GIMMS NDVI数据进行延长。

图2 黄河流域 GIMMS 3g NDVI 和 MODIS NDVI 数据的对比Fig.2 Comparison of GIMMS 3g NDVI and MODIS NDVI Data in the Yellow River Basin

2.2 季节性NDVI的年际变化特征

图3 为黄河流域春、夏、秋和冬季NDVI年际变化及其Mann-Kendall（M-K）检验曲线。从春季NDVI年际变化图来看，1982-2020黄河流域春季NDVI值呈增加趋势，线性变化倾斜率为0.002 7 ∕a，说明黄河流域近40 年春季植被覆盖呈增加趋势。其中1982 年和2019 年分别是该地区近40 年NDVI最小和最大的年份，其NDVI值分别为0.251 和0.331。由UF曲线可以看出，春季NDVI为持续上升趋势，但在2000年之前上升趋势呈锯齿状变化，2000 年之后为稳定上升趋势。UF和UB曲线相交于2005年，说明2005年是春季NDVI突变的开始，也表明了黄河流域春季NDVI自2005 年之后NDVI上升速率显著增加。秋季NDVI的UF曲线和春季相似，线性变化倾斜率为0.001 6 ∕a。2002年为秋季NDVI突变年，从曲线来看2002年之后秋季NDVI增加速率显著大，其增加速率大于春季（UF变化值）。夏季NDVI值在近40年也呈增加趋势，线性变化倾斜率为0.001 8 ∕a。其中1982 年和2018 年分别是该地区近40 年NDVI最小和最大的年份，其NDVI值分别为0.436 和0.531。由UF曲线可以看出，黄河流域夏季NDVI为持续上升趋势。UF和UB曲线相交于2008 年，说明2008 年是夏季NDVI突变开始的年份，从UF曲线可以看出2008 年以前NDVI变化速率基本为稳定状态，但在2008 年之后夏季NDVI迅速增加。冬季NDVI值线性变化倾斜率为0.0014∕a，为四季增加速率最小的季节。其中1985 年和2017年分别是近40年冬季NDVI最小和最大的年份，其NDVI值分别为0.182和0.266。UF和UB曲线相交于2007年，说明2007年是冬季NDVI突变的开始，且在2007年之后NDVI上升速率显著增加。由UF曲线看出，冬季和夏季NDVI的UF曲线走势较为一致，但冬季NDVI的UF曲线在2007 年的变化波动性较大，其变化值小于夏季。

图3 黄河流域NDVI时序曲线和Mann-Kendall（M-K）突变点检验曲线图Fig.3 NDVI time series curve and Mann Kendall （M-K） mutation point test curve in the Yellow River basin

总体而言，夏季NDVI增加速率最快，其次为秋季，再为春季和冬季。四季的NDVI突变年份均在2000 年之后，且突变年后NDVI增加速率均为上升趋势，这也进一步说明黄河流域1999年之后实施的生态恢复工程对植被覆盖度增加作用显著。

2.3 季节性NDVI的空间分布特征

从1982-2020年黄河流域NDVI不同季节的平均值发现（图4）：干旱区NDVI值最低，其次为半干旱区，湿润区NDVI值最高，其次为半湿润区，这表明NDVI的空间分布与水热资源分布一致。春季NDVI最高值集中在半湿润区和湿润区，海拔高度基本在1 000 m 以上，NDVI 值大于0.4，植被覆盖度较高。37.33%的NDVI值小于0.2，集中分布在干旱和半干旱地区，极端最小值分布在河套平原地区（NDVI＜0.1）。但在半干旱的刘家峡地区的NDVI集中在0.3～0.4，高于周围地区。内蒙古、宁夏、山西、陕西和甘肃北部地区的NDVI值基本在0.2 以下，其植被覆盖度较低。而在其他地区NDVI值基本在0.2 以上，且植被覆盖度较高，尤其在陕西省的南部的植被覆盖有最高值，而在内蒙古的毛乌素沙漠的植被覆盖度有最低值，基本在0.2 以下。究其原因可能是因为在山西、陕西、甘肃陇东及河南地区等地区林地分布范围较大，土地荒漠化程度较黄土高原北部地区的内蒙古、宁夏、青海等地较小，植被绿化程度较高。夏季NDVI值较春季均有增加，尤其在黄河流域源头地区，NDVI基本大于0.4。半干旱区的NDVI也从春季的小于0.2 增加到了大于0.2，但在干旱区的内蒙古的NDVI仍小于0.2，河套平原部分地区仍小于0.1。秋季NDVI与春季NDVI空间分布较为一致，但在黄河源区的NDVI值大于春季，集中在0.4 以上。冬季黄河流域除半湿润和湿润区的植被覆盖度较高，57.530%的地区NDVI值在0.2以下，是四季植被覆盖度最低月。

图4 黄河流域NDVI空间分布和空间分布区间统计图Fig.4 Spatial distribution and spatial distribution interval statistics of NDVI in the Yellow River basin

从40 年间黄河流域不同季节NDVI变化趋势的空间分布、变化趋势等级所占面积百分比（图5）发现：近40年四季NDVI均为增加趋势，尤其在陕西北部、山西西部及甘肃陇东地区增加速率较大，但不同季节的NDVI增加速率存在显著的空间差异。春季NDVI仅有9.2%的区域为减小趋势，集中在小浪底、河套平原和毛乌素沙地区。1.6%的区域NDVI增加速率大于0.005 ∕a，集中在黄河流域下游地区。其他地区的NDVI变化速率为干旱区＜半干旱＜湿润地区，变化趋势与水热资源显著相关。夏季NDVI的在黄土高原地区的变化速率集中在0.001 ∕a 以上，但在黄河源区的NDVI变化速率集中在-0.001 ∕a～0.001 ∕a，且大部分区域NDVI为减小趋势。在小浪底、三门峡和壶口地区的植被NDVI也为减小趋势。秋季8.1%的区域NDVI为减小趋势，集中分布在三门峡地区，其他区域的变化情况与夏季较为一致。冬季24.9%的区域NDVI为减小趋势，主要分布在干旱和半干旱地区，是四季中NDVI呈减小区域最多的季节，仅有16.4%的区域NDVI增加速率高于0.001 ∕a，其他区域NDVI基本为不变趋势。

图5 黄河流域NDVI变化趋势图Fig.5 NDVI Change Trend in the Yellow River Basin

2.4 黄河流域土地覆盖变化

从1980、1990、2000、2005、2010、2015 和2020 年土地利用分布情况（图6、图7）发现：1980-1900 年间草地和水域面积的减少主要转换为耕地和建设用地，林地和未利用地基本没有变化。1990-2000 年间草地面积有一定的增加（增加了2.649%），耕地、林地和水域面积分别减少了2.798%、0.045%和0.073%，建设用地面积仍为持续增长趋势，增长了0.14%。2000-2005年间减少的耕地面积主要转换为草地、建设用地和未利用地。这主要可能与1999开始的退耕还林（还草）政策有关，将坡耕地退耕还林还草、宜林荒山荒地造林，因此在2000 年之后四季植被覆盖均出现转折点。2005-2015 年间草地面积减少了1.72%，林地和耕地面积分别增加了0.42%和1.8%。2015-2020年间虽然林地和建设用地面积分别增加了0.15%和1.24%，但耕地和未利用地面积减少了2.90%和1.24%，使草地面积增加了0.87%。据统计数据显示这与目前的农村人口减少和保护生态环境的意识增加，使植被覆盖在2015-2021 年间有所改善［16］。有研究发现这是因为特定生态环境变化和人类活动导致部分草地转为耕地、林地和其他土地［18］，这主要是由于人口快速增长和三北防护林体系工程的实施，建设用地的发展，导致许多草地遭到破坏，并且长期过度放牧加上气候变化和啮齿动物破坏的影响，长江源区天然草地退化严重，将不可避免地导致长江源区荒漠化加剧，且研究发现黄河源区植被覆盖度呈微弱的减小趋势［14］。然而，随着政府实施了一系列移民安置项目，导致黄河流域建设用地持续占用农田、草地和湿地。因此，小浪底等地区植被覆盖为退化趋势。

图6 黄河流域土地利用变化图Fig.6 Land use change map of different land use areas in the Yellow River basin

图7 黄河流域不同土地利用面积百分比统计图Fig.7 Percentage statistics of different land use areas in the Yellow River basin

2.5 季节性NDVI对地形因子的响应

海拔和坡度会影响区域水热条件、植被类型等，因此基于参考文献［22］将黄河流域海拔划分为极低海拔区（＜1 000 m）、低海拔区（1 000～2 000 m）、较低海拔区（2 000～3 000 m）、中低海拔区（3 000～4 000 m）和高海拔区（＞4 000 m），将坡度分为平地（0°～2°）、平坡（2°～6°）、缓坡（6°～15°）、斜坡（15°～25°）和陡坡（＞25°）。然后统计了黄河流域四季NDVI变化速率在不同分区的分布状况，发现NDVI变化速率在不同高程梯度中存在显著差异，总体表现为随着海拔高度增加，NDVI增加速率减小。其中，春季NDVI变化速率在极低海拔区和低海拔区分布相对分散，均值分布在0.001 8 ∕a。海拔高于2 000 m 之后，春季NDVI变化速率集中在0.001 6～0.001 8 ∕a，分布相对集中。夏季NDVI变化速率在小于3 000 m 的地区均分布较春季分散，且随着海拔高度增加其增加速率逐渐降低。秋季NDVI变化速率在海拔低于3 000 m 的时候下降较慢，在海拔为3 000～4 000 m 时，NDVI变化速率突然下降，之后随着海拔高度增加NDVI变化速率呈稳定态势。主要原因是因为海拔高于3 000 m 的区域人类活动较少，植被覆盖度相对较高，但在海拔超多4 000 m 时，已远越过林线，植被覆盖度较低，且该区域的生态恢复工程实施困难，因此NDVI变化相对稳定，且增加速率最低。海拔小于2 000 m的区域NDVI增加速率较大，植被覆盖改善效果明显，主要是因为自1999 年之后在该区域大面积实施退耕还林还草等政策，使植被改善效果最为显著。冬季NDVI增加速率在海拔小于1 000 m 时较高，均值为0.000 3 ∕a。在海拔大于1 000 m后，NDVI变化速率基本为0，且变化趋势较为稳定［图8（a）］。主要是因为黄河流域大部分地区冬季较为干旱，植被覆盖度相对较低，因此在海拔高于1 000 m 的地区NDVI的增加速率基本0，植被覆盖度改善不明显［25］。

图8 黄河流域NDVI变化速率随海拔和坡度变化图Fig.8 Variation of NDVI change rate with altitude and gradient in the Yellow River basin

从不同坡度区域来看，坡度对春季和冬季NDVI变化速率影响不大，但夏季和秋季NDVI变化速率随着坡度增加而显著减小，尤其在坡度大于6°之后夏季NDVI变化速率减小最明显。主要原因是因为秋冬季人类活动相对较夏秋季少。其次，在坡度较缓区域的人类活动较为密集，坡度较大区域人类活动受到地形限制逐渐减少，植被受到的干扰也就越小，且坡度大于15°的区域多为林地分布区域，土地利用不易发生转换，这就使随着波段升高，NDVI的变化速率较为稳定且较小，而在坡度较小的区域的NDVI变化速率较大且分散［图8（b）］。

2.6 季节性NDVI对所有因子的敏感性探测

通常认为人类活动（如土地利用变化）、地形（如坡度、海拔）和气候因子（平均气候和极端气候）的共同或单独作用导致植被覆盖发生变化。虽然目前针对上述的单一因子对植被覆盖的影响的研究很多，但定量描述不同因子对植被覆盖度的贡献大小和因子之间的相互作用对植被覆盖度的影响的研究较少。因此，本文基于地理探测器的因子探测器和交互探测器分析了年均和季节性气候、土地利用和地形对植被覆盖的驱动机制。从因子探测器的结果可知（图9），近40年不同因子对NDVI的决定力q值大小为：日最高气温最小值（TXN）的q决定力＞日最低气温最大值（TNX）＞日最低气温最小值（TNN）＞年均气温（TEM）＞土地利用（LUCC）＞日最低气温最大值（TX）＞每年连续五天降雨量最大值（RX5DAY）＞风速（WIN）＞海拔（DEM）＞坡度（SLOPE）＞年均累计降水（PER）。由此可见，温度（平均温度和极端温度）和土地利用变化对NDVI的影响大于其他因子，年均降水量对植被NDVI的影响较小，其次为地形。

图9 不同因子对年均NDVI变化的q决定力值分布图Fig.9 Distribution of q driving force value of different factors on annual average NDVI change

通过不同因子交互探测结果我们发现（表1），任何两种因子的交互作用都大于或等于单个因子对NDVI的影响。其中，RX5DAY∩DEM（0.72）、RX5DAY∩SLOPE（0.70）、TXN∩DEM（0.74）、TXN∩TEM（0.71）、TNN∩DEM（0.27）的交互作用对NDVI的影响最大，影响因子大于0.7。DEM∩SLOPE（0.27）的交互作用对NDVI的影响最小。极端降雨和地形、极端气温和地形的共同作用对NDVI的影响大于其他因子的交互作用对NDVI的影响。这主要是因为黄河流域的大面积的土壤为黄土层，黄土土质松散，垂直节理发育，遇水则容易溶解，因此在坡度较大的地区极容易造成泥石流、滑坡等自然灾害，这对海拔较高的地区的林草地破坏极为严重（表1）［21］。其次，IPCC 第五次报告指出1998年之后全球地表温度由快速增温进入停滞期，这对地处干旱、半干旱和半湿润的黄河流域植被生长来说是有利的，地表温度增温的停滞有利于减小土壤表层水分的蒸散发，从而使土壤含水量提高［8，11，22］。因此，土地利用变化和气候等的综合影响使NDVI发生了突变。

表1 黄河流域年间NDVI系数Tab.1 NDVI Interaction in the Yellow River Basin

从不同季节的驱动因子随时间的变化特征来看（图10），发现TEM和PER在1980-2015 年间对植被NDVI的q决定力呈下降趋势，不同季节的极端气候指标对NDVI的影响情况在2000年达到最低，2010年次之，其他年份的影响均较高。而LUCC和SLOPE的影响力呈增加趋势，且LUCC的q值增加幅度较SLOPE大。以上原因主要是因为2000年之前黄河流域NDVI的变化很大程度上决定于气候因子的变化，但从1999年之后开展了一系列工程（如退耕还林、营造林、草地治理、耕地治理和移民搬迁等）使LUCC发生了显著变化，且其中的很多工程属于生态恢复工程，这些工程的实施极大的改变了该地区植被覆盖状况。因此在2000 年之后的植被NDVI变化对LUCC的变化较为敏感。其次，夏季的SLOPE的q值大于春、秋和冬季，有研究发现主要是因为黄河流域的黄土高原大部分地区夏季容易形成及时雨、大暴雨等天气，而由于黄土层的土壤性质，导致在SLOPE较大的地区发生滑坡、泥石流等自然灾害，从而对植被覆盖产生显著的影响。以上结论与李晓丽等［25］研究结论一致。

图10 不同因子对季节性NDVI变化的q决定力值分布图Fig.10 Distribution of q driving force values of different factors on seasonal NDVI changes

3 结 论

基于GIMMS NDVI 数据分析了1982-2020 年间黄河流域植被覆盖的季节性变化和土地利用、地形、极端气候和平均气候对季节性植被NDVI的贡献度和不同因子的交互作用对NDVI的影响。主要结论如下：

（1）1982-2020 年黄河流域干旱区NDVI值最低，其次为半干旱区，湿润区NDVI值最高，其次为半湿润区，即NDVI的空间分布与水热资源分布一致。

（2）1982-2020 年间黄河流域整体以增加趋势为主。春季NDVI增加速率为0.002 7 ∕a，其中1982 年和2019 年分别是该地区近40年NDVI最小和最大的两个年份，2005年NDVI变化的突变年。夏季NDVI增加速率0.001 8 ∕a，2008 年是NDVI突变年，且2008 年后夏季NDVI迅速增加。秋季NDVI增加速率为0.001 6 ∕a，2002 年为秋季NDVI突变年。冬季NDVI增加速率为0.001 4 ∕a，2007年是NDVI突变年。

（3）黄河流域NDVI随着海拔高度增加，NDVI增加速率减小，且夏季的这种趋势最明显。坡度对春、冬季NDVI变化速率影响不大，但夏、秋季NDVI变化速率随着坡度增加而显著减小，尤其在坡度大于6°后夏季NDVI变化速率减小最明显。

（4）近40 年不同因子对NDVI的q影响大小为：日最高气温的最小值（TXN）的q决定力＞日最低气温最大值（TNX）＞日最低气温最小值（TNN）＞年均气温（TEM）＞土地利用（LUCC）＞日最低气温最大值（TXX）＞每年连续五天降雨量最大值（RX5DAY）＞风速（WIN）＞海拔（DEM）＞坡度（SLOPE）＞年均累计降水（PER）。其中，任何两种因子的交互作用都大于或等于单个因子对NDVI的影响。且TEM和PER在1980-2020 年对植被NDVI的决定力呈下降趋势，而LUCC 和SLOPE的影响力呈增加趋势，且LUCC的q值增加幅度较SLOPE大。极端温度和降水对其影响没有明显的时间趋势。