植被净初级生产力（net primary productivity，NPP）是指绿色植物在单位时间、单位面积内由光合作用产生的有机物总量扣除自养呼吸后的剩余部分（Field et al.，1998）。植被NPP是研究陆地生态系统物质和能量交换的重要指标，其空间分布与区域气候、植被生长以及人类活动等因素息息相关，该指标的变化能反映植被群落的生产能力，是生态系统功能和结构变化的重要表征（朱文泉等，2007）。在“双碳”目标的大背景下，研究植被NPP时空变化特征及影响机制对于固碳增汇和生态修复治理工程等具有重要的意义。

植被NPP的早期研究主要为站点实测和统计模拟，近年来，随着新型遥感科学技术与方法的快速发展，凭借实时性强、易于获取、覆盖范围广等优点，利用遥感数据进行模型准确估算已成为区域NPP监测和时空变化趋势分析的强有力手段，如采用BIOME-BGC模型（Sun et al.，2017）、CASA模型（Liu et al.，2018）等。MOD17A3HGF是利用BIOME-BGC模型及光能利用率模型模拟得到的全球生态系统植被NPP数据，该数据集已经在全球不同地区研究中得到应用与验证（Wang et al.，2021; Venter et al.，2021; Ge et al.，2021）。国内学者从不同空间尺度探讨了我国众多区域如秦巴山区（李金珂等，2019）、西辽流域（朱丽亚等，2020）、南方农牧交错带（赵唯茜等，2018）、青藏高原（陈舒婷等，2020）、山东省（骆艳和张松林，2019）、中国东南部（崔林丽等，2016）等植被NPP时空变化及驱动因素，均发现NPP与气温、降水有一定的相关性，但不同地区气象因子对植被NPP的影响力存在明显差异; 当然，植被NPP的时空异质性除了与气象因素有关外，还受到地形、植被类型等其他自然因素（朱利欣和袁金国，2019; 刘婧等，2021）和人为因素（袁甲等，2016; Ge et al.，2021）的影响。因此，区域尺度植被NPP时空异质性的驱动机制尚未完全清晰。

近年来，海河流域由于受气候变化、人为干扰等因素影响，生态系统十分脆弱，水土流失、植被退化等生态环境问题日益突出（杨艳丽，2017）。植被NPP是判定区域生态系统变化及其可持续性的重要指标，而目前关于该区域植被NPP的研究较少，现有大多数研究多集中于京津冀地区，且更多从单一因素考虑对植被NPP的影响。因此，本研究以海河流域为研究区，利用2000—2020年MOD17A3HGF数据集，结合气温、降水、地形和土地利用等数据，采用相关分析、趋势分析以及地理探测器等方法，通过对海河流域近21年间植被NPP的时空格局动态演变特征以及从自然因素和人为因素对其变化的主要影响因素进行研究，拟探讨以下问题:（1）2000—2020年海河流域植被NPP时空变化规律;（2）海河流域植被NPP未来可能的发展趋势;（3）定量识别海河流域植被NPP时空异质性驱动因素及其影响力。本研究在弥补海河流域当前研究不足的同时为该区域植被恢复治理、生态建设工程实施提供科学参考和决策依据。

1 研究区概况与数据来源

1.1 研究区概况

海河流域位于中国华北平原北部，112°—120° E、35°—43° N之间（图1:A），包括北京和天津两个直辖市、河北的大部分地区及山东、山西、辽宁、河南和内蒙古的部分地区，地势为西北部高东南部低，气候类型为温带东亚季风气候（赵安周等，2019）。东部为平原区，地势较低; 西部为山区，主要包括太行山脉和燕山山脉，地势较高。植被类型主要有阔叶林、针叶林、草地、灌丛、农用地，其占比分别约为4.01%、0.08%、40.61%、2.76%、52.54%（图1:B）。由于海河流域受到坡度大、暴雨集中、土质疏松等自然因素和城市扩张等人为活动的强烈干扰，干旱、沙尘暴和水土流失等自然灾害严重，生态系统十分脆弱。

1.2 数据来源及处理

2000—2020年植被NPP数据为美国国家宇航局NASA（https://lpdaac.usgs.gov/）提供的MOD17A3HGF数据产品，该数据来自给定年份的所有8 d净光合作用（net photosynthesis，PSN）产品MOD17A2HGF数据之和，空间分辨率为500 m，时间分辨率为1 a。结合产品说明，取数据集有效值为-3 000~32 700，剔除范围以外的值。气象数据选取中国气象数据网（https://data.cma.cn/）提供的2000—2020年海河流域及其周边36个气象站点的月均气温、月降水量、平均风速、日照时数以及平均相对湿度数据，利用反距离加权法（inverse distance weighted，IDW）将其插值为空间分辨率为500 m的栅格数据。DEM数据为地理空间数据云（http://www.gscloud.cn/）提供的SRTM 90 m数据，坡度数据为通过高程数据计算所得，根据研究区实际情况（申丽娜，2018），基于高程信息将研究区域划分为平原、丘陵、低山、中山、高山5类，坡度划分为微坡、缓坡、斜坡、陡坡、急坡5类。土地利用数据选取2000、2010、2020年3期30 m全球地表覆盖Globeland30数据，下载地址为http://www.globallandcover.com/，植被类型数据源于NASA提供的2000—2019年MCD12Q1数据产品，空间分辨率为500 m，时间分辨率为1 a，根据产品数据说明，将土地利用和植被类型数据分别重分类为6类（耕地、林地、草地、水域、建设用地、未利用土地）、5类（针叶林、阔叶林、灌丛、草地、农用地）。将高程、坡度和土地利用等重采500 m的数据，使上述数据具有相同的像元大小。

2 研究方法

2.1 变异系数

采用变异系数对植被NPP随时间变化的稳定性进行分析（刘宪锋等，2015），按如下公式计算。

式中:CV为变异系数; CV值越小，表示植被NPP波动性越小，反之则越大; NPP_i为第i年的植被NPP值; $\overline{NPP}$为海河流域2000—2020年NPP的平均值，n为21。

2.2 趋势分析

采用一元线性回归法对植被NPP的变化趋势进行分析（穆少杰等，2013），按如下公式计算。

式中:θ_Slope为变化率，θ_Slope＞0，表示植被NPP呈上升趋势，反之则下降; n为21; NPP_i为第i年的NPP值。查F分布临界值表，将植被NPP变化趋势分为5类，即极显著减少（θ_Slope＜0，P＜0.01）、显著减少（θ_Slope＜0，0.01≤P＜0.05）、变化不显著（P≥ 0.05）、显著增加（θ_Slope＞0，0.01≤P＜0.05）、极显著增加（θ_Slope＞0，P＜0.01）。

2.3 相关性分析

采用逐栅格空间分析法，探讨植被NPP与气温、降水之间的相关性（骆艳和张松林，2019），利用t检验法对结果进行显著性检验，按如下公式计算。

式中:R_xy为相关系数; x_i为第i年的植被NPP值; y_i为第i年的气温或降水值; $\stackrel{-}{x}$为21年间植被NPP平均值; $\stackrel{-}{y}$为21年间气温或降水平均值; m为控制变量个数; n为样本数。

2.4 Hurst指数

利用Hurst指数探究研究区植被NPP的未来变化趋势（任涵玉等，2021），基本原理为假定一个植被NPP时间序列NPP_t，对于任意正整数，则定义以下内容。

均值序列:

累计离差:

极差:

标准差:

若存在$R/S\propto {\tau }^H$，表明植被NPP时间序列存在Hurst现象，称H值为Hurst指数，包括持续性（0.5＜H＜1）、随机性（H=0.5）、反持续性（0＜H＜0.5）3种情况。

2.5 地理探测器

地理探测器是探测空间分异性，以及揭示其背后驱动因素的一种统计学方法（王劲峰和徐成东，2017）。本文选取分异及因子探测器对海河流域2000—2020年植被NPP的影响因素进行分析，按如下公式计算。

式中:h为层数; N_h和N分别为层h和研究区的样本数; σ_h和σ分别为层h和研究区植被NPP的方差; q的值域为 [0，1]，q值越大表示该因素对海河流域植被NPP影响力越大，反之则越小。

3 结果与分析

3.1 植被NPP时间变化特征

海河流域2000—2020年植被NPP总量波动范围为72.02~117.41 Tg C（1 Tg C=10¹² g C，图2:A），多年均值为98.73 Tg C。研究期内植被NPP总量呈上升趋势，其中2020年NPP总量达到最大值，为117.41 Tg C; 2001年NPP总量达到最小值，为72.02 Tg C。就不同年份NPP总量偏离其均值的程度而言，2001、2007、2020年偏离均值的程度较大，NPP变化显著，而2005—2006年、2009—2010年偏离均值的程度较小，NPP变化不显著。

将海河流域不同年份的NPP均值划分为6个区段，对各区段所占面积的百分比进行统计（图2:B）。由图2:B可知，2001年NPP均值最低，为238.44 g C·m^-2·a^-1，2020年NPP均值最高，为388.54 g C·m^-2·a^-1，2000—2020年NPP均值为326.75 g C·m^-2·a^-1。NPP主要集中于200~500 g C·m^-2·a^-1，面积占比为69.49%~94.74%，为海河流域NPP均值集中区。

利用2000—2019年植被类型数据对2000—2019年NPP数据进行分区统计，对于缺少的2020年植被类型数据，考虑到相邻年份差别不大，因此，采用2019年植被类型数据对2020年NPP数据进行统计，得到不同植被类型NPP均值年际变化趋势（图3）。由图3可知，不同植被类型的年均NPP存在一定的差异性，不同年份NPP变化趋势与总趋势（图2:A）基本一致。不同植被类型多年NPP均值排序为阔叶林（385.28 g C·m^-2·a^-1）＞灌丛（353.03 g C·m^-2·a^-1）＞草地（320.12 g C·m^-2·a^-1）＞针叶林（313.59 g C·m^-2·a^-1）＞农用地（295.22 g C·m^-2·a^-1）。阔叶林和灌丛NPP多年均值明显高于全区多年均值（326.75 g C·m^-2·a^-1），因此，阔叶林和灌丛为该区固碳、生产力最强的植被。

3.2 植被NPP空间变化特征

3.2.1 NPP年均值及变异系数

海河流域2000—2020年不同地区年平均植被NPP波动范围为29.32~746.84 g C·m^-2·a^-1（图4: A），空间差异性显著，整体呈现西高东低、北高南低的空间格局，高值区主要分布于长治市、晋中市及承德市等地，而低值区主要分布于邯郸市、邢台市及衡水市等地。结合海河流域自然条件、地形地貌及人为活动等因素分析可知，西北部地区主要为太行山脉、燕山山脉，植被覆盖度高，植被类型主要为阔叶林和灌丛，受人为影响小，植被长势较好，植被NPP值高，东南部地区主要为平原，大多为耕地和城镇区，植被较为稀疏，受人为影响大，植被NPP值低。

NPP变异系数为0.05~1.80（图4:B），均值为0.17，水平相对较低，表明研究期内海河流域植被NPP空间上整体处于相对稳定状态，植被生态系统处于良性循环。高值区域主要集中分布在北京市中西部以及承德市西南部，该地区近些年来大规模进行退耕还林还草、封山育林等生态保护建设工程，使得NPP波动性相对较大。

3.2.2 趋势分析

海河流域2000—2020年植被NPP的θ_Slope为-23.27~30.10（图5:A），呈增长、减少的占比分别为97.67%、2.33%，表明研究区NPP空间上整体呈现增长趋势，与上述时间上变化趋势结论相同。由图5:B可知，该地区植被NPP整体变化显著，显著变化面积占比为71.04%，其中显著增加比例为70.78%，显著减少比例为0.26%。就空间分布而言，NPP极显著减少和显著减少地区位于东部沿海; 极显著增加和显著增加地区位于西北部山区; 变化不显著地区位于东南部平原。总体而言，海河流域NPP呈显著增加趋势，结合研究区自然条件和一些生态保护政策导向可知，西北部海拔高，生态环境脆弱，建议应持续加强该地区生态管控和巩固现有成果，而东南部海拔较低，多为平原，适合开展农业活动和城镇化建设，但考虑到整个海河流域的生态环境情况，应坚持保护优先的发展战略。

3.2.3 NPP未来变化趋势

通过Matlab编程计算得到Hurst指数，结果见图6:A。由图6:A可知，海河流域2000—2020年植被NPP的Hurst指数位于0.10~0.96，均值为0.44。将Hurst指数结果划分为强反持续（0＜H≤0.25）、弱反持续（0.25＜H＜0.5）、不确定（H=0.5）、弱持续（0.5＜H≤0.75）、强持续（0.75＜H≤0.96）5个等级，并统计各等级的占比（图6:B）。由图6:B可知，海河流域未来植被NPP反持续性大于不确定和持续性，其中强反持续、弱反持续占比分别为0.51%、71.39%，说明该区植被NPP未来变化趋势整体将呈现反持续性变化，即与过去变化趋势相反。就空间分布来看，强反持续性分布较为零散，弱反持续性主要集中分布于西部山区和承德等地，主要在于西部地区主要为燕山和太行山脉，地势较高，生态环境脆弱，而承德位于坝上草原地区，受人为影响较大，因此未来治理应重点关注该地区，提高其生态系统稳定性。

3.3 植被NPP驱动因素分析

3.3.1 NPP与气象因素的关系

海河流域2000—2020年年均气温的变化范围为10.31~11.72℃，多年均值为11.15℃，呈上升趋势; 年降水量的变化范围为367.14~617.35 mm，多年均值为498.81 mm，呈上升趋势（图7）。植被NPP与气温的相关系数为-0.85~0.80（图8:A），均值为-0.04，正相关的比例为42.97%，负相关的比例为57.03%，表明植被NPP与气温整体表现为负相关关系，其中，极显著负相关和不显著负相关面积占比分别为1.11%和53.46%（图9:A）。由图9:A可知，植被NPP与气温呈不显著负相关关系; NPP与降水的相关系数为-0.64~0.86（图8:B），均值为0.47，正相关的比例为99.40%，负相关的比例为0.60%，表明植被NPP与降水整体表现为正相关关系，其中，极显著正相关和不显著正相关的面积占比分别为37.76%和40.16%（图9:B）。由图8:B可知，植被NPP与降水呈不显著正相关关系。就海河流域植被NPP与气象因素的相关系数均值和显著性面积大小而言，降水均高于气温，表明该地区植被NPP更易受到降水影响，降水是影响植被NPP变化的主要气象因素（图10）。

3.3.2 NPP与地形因素的关系

对高程、坡度不同分级的植被NPP均值及总量进行了统计，结果见表1。由表1可知，海河流域以平原、中山为主，高程在0~1 500 m占比约为93.55%，其中，平原所占比例比最大，高山所占比例较小。就不同高程等级的植被NPP均值而言，存在一定差异，总体变化趋势为随着海拔的升高先增大后减小，具体表现如下:平原地区容易受到城镇化建设、农业耕作等人类活动干扰，植被生产力处于较低水平，但因其面积占比高，所以植被NPP总量较高; 植被类型为以林地为主的丘陵和低山地带，植被生产力处于较高水平，但因其面积占比低，所以植被NPP总量较低; 植被生长环境恶劣且较为稀少的高山地带，植被生产力水平较低，同时由于面积占比小，因此植被NPP总量较小。坡度方面，海河流域主要以微坡、斜坡为主，坡度在0°~15°的面积占比约为97.94%，其中，微坡所占比例最大，急坡最小。就不同坡度等级的植被NPP均值而言，变化趋势与高程一致，即随坡度的增加先增大后减小，具体表现为:微坡植被NPP总量较高，而当坡度逐渐增加，坡面侵蚀强度逐渐增大，植被生长受到一定的限制，植被覆盖度降低，NPP总量下降。

3.3.3 NPP与土地利用变化的关系

将土地利用数据进行叠加分析，得到土地利用变化转移矩阵（表2）。由表2可知，2000—2020年期间海河流域有20.48%的土地利用类型（65 076 km²）发生了变化，具体表现如下:耕地净转出16 403.25 km²，主要转向草地（5 501.5 km²）和建设用地（18 583.5 km²）; 草地净转出3 055.25 km²，主要转向耕地（4 826.25 km²）和林地（11 823.75 km²）; 林地净转入2 573 km²，主要来自耕地（1 906.25 km²）和草地（11 823.75 km²）; 水域净转入627.75 km²，主要来自耕地（2 340.25 km²）和草地（819 km²）; 建设用地净转入15 868.5 km²，主要来自耕地（18 583.5 km²）和草地（1 240 km²）。总体而言，2000—2020年期间，海河流域土地利用变化表现为耕地、草地的转出以及林地、水域、建设用地和未利用土地的转入。

将不同土地利用类型转移量乘以单位面积植被NPP变化量，得到植被NPP变化转移矩阵（表3）。由表3可知，2000—2020年期间，海河流域耕地转为草地造成的植被NPP增量和转为建设用地造成的植被NPP减量分别为732.22×10^-3、294.65×10^-3 Tg C; 草地转为耕地和林地造成的植被NPP增量分别为62.94×10^-3、47.18×10^-3 Tg C; 耕地和草地转为林地造成的植被NPP增量分别为317.69×10^-3、47.18×10^-3 Tg C; 耕地和草地转为水域造成的植被NPP减量分别为292.59×10^-3、75.08×10^-3 Tg C; 耕地和草地转为建设用地造成的减量分别为294.65×10^-3、92.95×10^-3 Tg C。耕地转为草地导致的植被NPP变化量（732.22×10^-3 Tg C）高于其他转换方式带来的NPP变化量。综上表明，耕地转为草地是海河流域土地利用类型变化引起NPP变化的主要地类转换，该转换方式在一定程度上促进了研究区植被NPP的增加。

3.3.4 驱动因素影响力分析

植被NPP的时空分异性是多种因素综合作用的结果。基于上述研究结果分析，从自然因素和人为因素两个方面，选取高程、坡度、降水、气温、湿度、风速、日照时数和土地利用共8个因素，利用生成的31 023个3 km × 3 km等间距网格点（采样点）分别提取植被NPP数据和不同驱动因素数据的信息，采用分异及因子探测器对海河流域2000、2010、2020年植被NPP进行单因素影响力探测，结果见图10。就得到的驱动因素q值而言，值越大表示该因素对植被NPP的影响力越大，则认为其为该年NPP空间分异的主导因素。由图10可知，海河流域植被NPP不同年份的q值大小排序依次为降水（0.252）、坡度（0.228）、高程（0.235），说明该研究区没有绝对的主导因素。不同年份不同因素的q值存在明显差异性，其中波动性最大为降水，该现象在一定程度上体现了研究区区域降水的不均匀性。不同因素3个年份q值均值大小排序依次为降水（0.215）＞高程（0.212）＞坡度（0.209）＞湿度（0.198）＞气温（0.194）＞日照时数（0.192）＞土地利用（0.190）＞风速（0.185），就q值每年均值与多年均值（0.20）关系而言，可将海河流域植被NPP变化影响因素划分为主导因素和重要因素两种类型，即降水、高程和坡度是主导因素，而湿度、气温、风速、日照时数和土地利用则是重要因素。

注: —表示该方式下不转换。

Note: — indicates no conversions in this mode.

4 讨论与结论

近年来，由于受气候变化及人类活动的影响，海河流域面临着水土流失、水资源短缺等生态压力。本研究表明，时间上，2000—2020年海河流域植被NPP整体呈增加趋势，主要得益于退耕还林还草、封山育林政策的实施以及一些生态工程的建设，如京津冀防风固沙工程、平原绿化工程等，在一定程度上促进了海河流域植被恢复，生态环境得到了一定程度改善（杨艳丽，2017）; 空间上，植被NPP均值呈现由西北向东南递减的趋势，主要在于西北部山区植被以林地为主，而东南部平原地区则多为耕地。另外，不同植被类型年均NPP存在一定的差异，阔叶林、灌丛为研究区固碳、生产力最强的植被，这与李肖等（2018）研究结果一致，同时也应注意到由于研究时段、数据来源等不同，与潘竟虎和文岩（2015）的研究结果存在一定的差异。

本研究结果表明，气候变化、地形条件以及人类活动共同影响着海河流域植被NPP的变化。气温、降水的增加在一定程度上促进植被生长，但激烈的气候变化也会对植被生长产生抑制作用，进而影响植被NPP的积累，如气温过高会导致植物气孔关闭，从而影响植被的光合作用。研究期内，海河流域气温和降水呈逐步上升趋势，植被生长环境比较适宜，进而促进了研究区植被NPP增加。由于研究区地势高差大，不同高程范围内植被NPP存在一定差异，整体呈现先增后减的趋势，500~1 000 m处植被NPP值最高，主要在于该范围内受人类活动影响较小，水热条件优越，植被多为碳吸收能力较强的林地。除上述因素外，人类活动也会对区域植被NPP的积累产生影响，最显著的为土地利用的改变，2000—2020年期间，研究区建设用地面积增加了15 868.5 km²，造成的植被NPP损失量最大，为387.6×10^-3 Tg C; 与此同时，退耕还草使得植被NPP增加了732.22×10^-3 Tg C; 本世纪初期，政府在全面促进社会经济发展的同时，海河流域内植被受到较大影响，造成植被退化，而近年来在国家大力推进生态文明建设的背景下，研究区持续开展生态保护、修复和建设工程，植被得到了很大程度恢复，促进了植被NPP的增加。

驱动因素影响力上，降水q值均值最高，表示降水对海河流域植被NPP的时空分异性解释力最强，主要在于海河流域处于半湿润半干旱气候过渡区，降水的增加会促进植被的光合作用，从而使植被NPP增加。地形因素中，高程解释力大于坡度，主要是高程可直接影响植被类型的空间分布，而坡度则是通过坡面侵蚀强度间接影响植被的生长状况，由此可见，地形往往通过控制水热和土壤条件，影响其他环境变量，进而对区域植被格局产生影响（Mao et al.，2014; 潘洪义等，2019），与当前研究成果（朱利欣和袁金国，2019）具有较好的一致性。虽然过去21年间海河流域土地利用方式变化较为剧烈，但其类型并未发生较大改变，仍主要为耕地、林地和草地，相比于气象和地形因素，对植被NPP产生的影响不大。

未来海河流域生态保护与建设工作应重点考虑降水、地形以及人类活动等方面的影响。例如，多关注低海拔地区植被修复，合理分配水资源，持续推进退耕还林还草、封山育林等政策和生态修复治理工程，同时需从多角度、多维度持续加大生态保护力度。另外，需要指出的是，在人为因素对区域植被NPP影响方面，目前仅考虑到土地利用类型整体变化，后续研究应定量和细化区分城市扩张、生态建设等人类活动的影响，以降低植被NPP影响因素的不确定性; 土壤类型可制约植物能否有效吸收水分，进而影响植被的生长，海河流域土壤类型丰富，不同土壤类型对该区域植被NPP时空分异性的影响有待进一步的探讨。

综上所述，得出如下结论。（1）2000—2020年期间，海河流域植被NPP整体呈显著上升趋势，NPP总量均值为98.73 Tg C，年均值为326.75 g C·m^-2·a^-1，阔叶林、灌丛为固碳、生产力最强的植被。（2）空间上，植被NPP整体呈现西北部高、东南部低的空间格局。NPP空间波动性不大，变化趋势呈增长和减少的区域占比分别为97.67%和2.33%。未来变化趋势呈反持续性，植被NPP可能下降。（3）植被NPP与气温呈不显著负相关关系，与降水呈不显著正相关关系，降水为影响植被NPP变化的主要气象因素。植被NPP随着高程和坡度增加均呈“增加—减少”的变化趋势。耕地转为草地一定程度上促进了研究区植被NPP增加。（4）驱动因素影响力上，不同驱动因素对植被NPP影响力存在差异，主导因素为降水、高程和坡度，重要因素为湿度、气温、风速、日照时数和土地利用，人为因素影响力低于自然因素。

参考文献