无定河流域土地覆被空间分异机制及相关水碳变量变化

吕锦心，梁 康，刘昌明，张仪辉，刘 璐

(1.中国科学院地理科学与资源研究所，中国科学院陆地水循环及地表过程重点实验室，北京 100101；2.中国科学院大学，北京 100049)

土地覆被作为连接土壤、大气、水文等生态要素的重要纽带，深刻影响陆气能量交换、陆地水碳循环和区域生态环境演变。作为全球环境变化的重要组成部分［1-2］，土地覆被变化不仅能反映陆地表层景观的时空动态过程，还可以在不同程度上改变全球和区域的水循环、能量循环和生物物质循环过程［3-4］。因此，国内外学者在全球、区域及流域尺度上对土地覆被动态格局、驱动机制、生态环境效应等方面做了大量研究［5-10］，对认识气候变化方向、人-地复杂系统以及地-水-碳相互作用机理等具有重要参考价值。

自20 世纪70 年代以来，中国在黄土高原地区开展了大规模生态保护和生态恢复工程，人为对当地土地覆被格局进行了干预［11-12］。黄土高原位于干旱半干旱区，是我国生态系统最脆弱的地区之一，对气候变化和人类活动十分敏感［13］。其中无定河流域作为典型的生态脆弱地区［14］，是黄土高原退耕还林还草生态恢复工程实施的重点区域，也是国家贯彻实施黄河流域生态保护和高质量发展战略的重要地区，流域内人-地相互作用剧烈，覆被变化明显，尤其植被恢复显著［15］。针对无定河流域以及整个黄土高原地区，以往与土地覆被相关的研究侧重分析土地覆被（例如NDVI）本身的格局变化［16-17］，或强调覆被变化在改善环境、水土保持等方面发挥的作用和导致河川径流下降等问题［18-20］，而考虑地-水-碳耦合机制以评估生态系统关键变量的研究相对薄弱；此外，对土地覆被空间分异机制的定性分析多，定量化研究及机理研究相对不足。因此，综合评价土地覆被的时空格局、空间分异机制、相关水碳变量变化特征具有重要意义。

以往开展土地覆被时空变化特征分析时，针对长时间连续数据序列，通常采用简单线性回归、Sen’s斜率估计、Mann-Kendall 检验等方法描述变量的时间变化［21-22］；针对不连续数据，多采用转移矩阵、图谱分析、重心迁移模型等方法分析土地覆被的转型结构、时空演变等［23-25］。在开展土地覆被变化驱动力分析时，通常将影响要素划分为自然要素和社会经济要素两大类，其中自然要素影响主要表现为累积性效应，而社会经济要素的影响效果快速而显著［26］，以往研究主要通过定性描述或使用加权回归、logistic 回归模型等方法对驱动力进行探究［27-29］，然而相关研究往往主要针对土地覆被的时间变化，对其空间分异机制的研究相对较少。近年来，地理探测器［30］为单期下垫面数据的空间驱动机制分析提供了有效工具，在量化多环境因子对土地覆被空间格局的影响方面获得广泛应用［31-33］。在水碳变量分析方面，蒸散发（Evapotranspiration，ET）连接水循环和能量循环等关键过程，是绿水资源耗散状况的重要指标［34-35］，总初级生产力（Gross Primary Productivity,GPP）是表征陆地生态系统吸收CO2的重要指标［36-38］，水分利用效率（Water Use Efficiency，WUE）是耦合生态系统碳水关系的综合指标［39-40］，这3 个水碳变量作为陆地水碳循环的关键要素和生态系统关键变量，对土地覆被变化具有敏感响应，因此常被用做土地覆被水碳效应的分析指标［41-44］。然而，由于不同源遥感ET 与GPP 数据会导致WUE 的错误估算，加上在数据分析时考虑耦合水碳过程的遥感产品有限，当前研究主要侧重对单个变量的分析，综合水碳变量的研究较少。

基于此，本文以无定河流域为研究对象，首先使用连续的土地覆被数据，利用Mann-Kendall 趋势检验、空间分析等方法，分析流域土地覆被的时空变化，在此基础上，结合气象、DEM、人口、GDP 数据，运用地理探测器探究多个自然因子和社会经济因子对土地覆被空间分异的影响，最后，依靠水碳耦合的PML_V2 陆地蒸散发与总初级生产力数据集，在同源数据基础上分析流域蒸散发、总初级生产力、水分利用效率的时空变化，以期为黄土高原地区的水土资源利用、生态文明建设等提供科学参考。

1 数据与方法

1.1 研究区概况

无定河流域位于37°09′～39°03′N、107°45′～110°34′E 之间，跨内蒙古、陕西两省，面积约为3.6×104km2，流域高程在834～2001 m 之间，西北高，东南低，其地理位置及高程如图1所示，西北部位于毛乌素沙地南缘，东南部位于黄土高原丘陵沟壑区，生态脆弱。气候属温带大陆性季风气候，流域多年平均气温约为8.9 ℃，多年平均降水量约为369.9 mm，年内分配不均，主要集中在7—9 月，降水空间分布为东南多，西北少［45］。无定河是黄河中游河口至龙门区段最大的支流，发源于陕西省榆林市定边县白于山北麓，流经定边、靖边、米脂、绥德和清涧县，由西北向东南注入黄河，干流全长491 km，主要支流有北岸的纳林河、海流兔河、榆溪河和南岸的芦河、大理河、淮宁河等。

图1 无定河流域地理位置及高程Fig.1 Location and elevation of the Wuding River Basin

1.2 数据来源

土地覆被数据（1990—2019年，共30期）来自中国土地覆盖产品［46］，空间分辨率为30 m，数据基于5463个目视解译样本，总体精度达80%。本文主要关注草地、耕地、林地、建设用地、荒地等5 大类型，故对源数据进行重分类，其中耕地为源数据中的农田，林地为源数据中的森林和灌木，建设用地为源数据中的不透水区，草地与荒地和源数据类型相同。

年平均气温、年降水量、人口密度、GDP、高程等数据来自中国科学院资源环境科学与数据中心（https://www.resdc.cn/），空间分辨率均为1 km。年平均气温、年降水量、人口密度、GDP等数据的年份为1995 年、2000 年、2005 年、2010 年、2015 年、2019年（共6 期）［47-48］，坡度、坡向数据基于高程数据，通过ArcGIS软件计算得来。

蒸散发（Evapotranspiration，ET）、总初级生产力（Gross Primary Productivity, GPP）数据来自“中国区域PML-V2 陆地蒸散发与总初级生产力数据集（2001—2019 年，共19 期）［49］”(http://data.tpdc.ac.cn/zh-hans/)，数据空间分辨率为500 m，基于包含中国9种植物功能类型的26 个涡动协方差通量塔的观测值校准显示，ET 数据的Nash-Sutcliffe效率（NSE）为0.75，均方根误差（RMSE）为0.69 mm·d-1；GPP 数据的NSE 为0.82，RMSE 为1.71 g C·m-2·d-1。根据定义，水分利用效率（Water Use Efficiency，WUE）为GPP与ET之比。

本文从以上数据集中提取了无定河流域范围内的数据，用于分析研究。

1.3 研究方法

1.3.1 趋势、突变检验方法 本文利用线性倾向估计、Mann-Kendall 趋势检验、Pettitt 突变检验来分析土地覆被及水碳变量序列的趋势性与突变点。其中线性倾向估计通过求变量关于时间的线性回归来表征变量的变化趋势，其中回归系数可以反映趋势方向和倾向程度［50］。Mann-Kendall趋势检验不需要样本数据遵从特定的分布，能够有效的剔除异常值，适用于非正态分布数据的趋势检验［51-52］。检验统计量Z表征数据序列的趋势性，|Z|值越大，趋势越明显；Z＞0，说明有增大趋势，Z＜0，说明有减小趋势。当给定显著性水平α后，若|Z|大于临界值，则序列的趋势显著，本文统一选用α=0.05 的显著性水平，所对应区间为±1.96。Pettitt 突变检验为非参数检验方法，最初由Pettitt［53］提出，用于检验突变点。

1.3.2 地理探测器 地理探测器［30］是探测空间分异性，以及揭示其背后驱动力的一组统计学方法，包括因子探测、交互作用探测、风险区探测和生态探测4部分。其中因子探测可用于探测某因子X对属性Y的空间分异的解释力度，该力度用q值（0 ≤q≤1）度量，q值越大表示自变量X对属性Y的解释力越强，计算公式如下：

式中：h=1,···，L；L为变量Y或因子X的分层，即分类或分区；Nh和N分别为层h和全区的单元数；σh2和σ2分别是层h和全区的Y值的方差。

本文选取年平均气温、年降水量、坡度、坡向、高程等自然因子，人口密度、GDP等社会经济因子，运用地理探测器探究不同自然、社会经济因子对土地覆被空间分异的影响。首先，利用最优离散化方法［54］，将年降水量、年平均气温、高程、坡度、坡向、人口密度、GDP 等数值型数据进行离散化，转化为类型变量，其次，利用最邻近内插法将土地覆被数据的空间分辨率重采样至1 km，最后利用地理探测器对土地覆被空间分异驱动进行计算。

2 结果与分析

2.1 土地覆被时空变化

草地和耕地是无定河流域最主要的土地覆被类型，1990—2019 年期间，流域整体草地、林地、建设用地显著增加，耕地、荒地显著减少，林草面积增加区域主要集中在流域下游及无定河沿岸地区，建设用地在榆林市周边明显增加，耕地主要减少在流域下游地区，而荒地主要减少在流域中上游地区。整体看来，1990—2019 年，无定河流域退耕还林还草生态恢复工程成效显著，生态环境得到改善。

具体而言，如图2、图3所示，草地在全流域内广泛分布，多年平均面积占流域总面积的68.62%，草地面积在1990—2019 年间呈显著增加趋势（Z=6.49，P＜0.05），并在2002 年存在突变点，增加区域主要分布在流域下游地区以及沿无定河水系附近；耕地是流域内面积第二的土地覆被类型，主要分布在流域中下游，多年平均面积占比为19.59%，虽然流域中上游局部地区出现新增耕地，但下游的耕地大面积减少，导致流域整体耕地面积在1990—2019年面积显著减少（Z=-3.39，P＜0.05）；荒地面积在流域内占到11.21%，主要分布在流域东北部，在1990—2019 年期间，荒地呈显著减少趋势（Z=-6.42，P＜0.05），减少地区主要位于无定河上游支流河岸附近，与新增草地区域有较好的一致性；建设用地面积占流域面积的0.33%，近年来面积呈增加趋势（Z=7.74，P＜0.05）,尤其是在2004年之后迅速增加，新增地区主要集中在榆林市；林地在流域内占比较少，为0.03%，但呈显著增加趋势（Z=3.96，P＜0.05），尤其在2006 年后迅速增加，新增林地主要集中在流域的最南部。

图2 1990—2019年无定河流域各土地覆被类型面积变化Fig.2 Area changes of land cover in the Wuding River Basin,1990-2019

图3 无定河流域1990—2019年土地覆被类型变化情况Fig.3 Changes of land cover in the Wuding River Basin,1990-2019

2.2 土地覆被空间分异驱动力探测

综上所述，无定河流域各土地覆被类型分布存在明显的空间分异性，并且不同区域土地覆被类型还存在不同方向的转变。因此，进一步从自然、社会经济等多角度出发，利用地理探测器对无定河流域6 期土地覆被的空间分布驱动进行了探测，结果如图4所示。

图4 土地覆被空间分异因子探测结果Fig.4 Factor detection results of land cover spatial differentiation

从多年平均解释力看，在各自然因子中，降水是影响无定河流域土地覆被空间分布机制的主控因素，在社会经济因子中，人口密度比GDP 的影响力大。具体而言，降水、气温、高程、坡度、坡向对无定河流域土地覆被空间分异的多年平均解释力q值分别为0.09、0.08、0.07、0.06、0.02，人口和GDP 的多年平均解释力q值分别为0.12、0.07，以上因子对研究区土地覆被空间分布的影响均通过显著性检验。

此外，本文将自然、社会经济方面解释力最大的前两项进行加和，得到分别代表自然和社会经济因子相对解释力大小结果的时间变化曲线。由图4可知，除2000年、2005年外，社会经济因素对土地覆被的空间分异解释力大于自然因素，但在1995—2019年间，社会经济因素的解释力存在微弱下降趋势，而自然因素的解释力存在微弱上升趋势，相比1995年，二者间的差距在减小。

2.3 相关水碳变量的时空变化

在自然与社会经济等多方面因素的影响下，无定河流域的土地覆被分布情况在2000 年后发生了明显的改变，这在一定程度上势必会影响到无定河流域水碳变量的时空变化。因此，本文选取GPP、ET、WUE这3个重要水碳变量展开分析，结果如图5所示。

图5 2001—2019年无定河流域GPP、ET、WUE的变化Fig.5 Changes of GPP,ET and WUE in the Wuding River Basin,2001-2019

从空间分布上看，如图5a、图5b、图5c 所示，无定河流域多年平均GPP、ET、WUE 的空间分布与土地覆被空间分异有较好的对应性，大体存在由东南向西北递变的规律。东南部以耕地、林地、草地为主要覆被的地区其GPP、ET、WUE相对偏高，西北部以草地、荒地为主要覆被的地区其GPP、ET、WUE相对偏低。此外，水系附近的GPP、ET、WUE 相对偏高。

从时间变化上看，如图6、图5e、图5f、图5g 所示，2001—2019 年间，流域整体GPP、ET、WUE 均呈增加趋势，GPP、WUE 在流域绝大部分区域均显著增加，而ET主要在流域中下游地区显著增加。进一步分析不同土地覆被类型的水碳变量及其变化可以发现，林地、耕地、草地、荒地等土地覆被类型的GPP、ET、WUE 也呈增长趋势，林地的GPP、ET、WUE 在各土地覆被类型中占据第一，耕地、草地次之，荒地最小。此外，不同土地覆被类型的GPP、WUE差别较大，并且不同土地覆被类型间差距还在扩大，而不同土地覆被类型间的ET差别较小。

图6 2001—2019年无定河流域整体及不同土地覆被类型的GPP、ET、WUE变化Fig.6 Changes in GPP,ET and WUE of the whole and different land cover in the Wuding River Basin,2001-2019

3 讨论

无定河流域整体GPP、ET、WUE均呈增加趋势，主要受到人为改变土地覆被与气候变化的影响。植被可通过光合作用吸收降低空气中CO2浓度，与其他土地利用类型相比，植被覆盖较大的林草地的固碳能力也相对较强。1999 年，陕西作为全国试点，首先开展了退耕还林还草生态恢复工程，故此大大促进了当地GPP的增长，这类碳固存增加主要体现林草地恢复明显的流域南部地区。而全流域大范围GPP增长可能由气候变化驱动，依据主要有两点，一是空间驱动力分析结果显示降水和气温对流域土地覆被格局的影响在逐渐扩大，毛乌素沙地的暖湿化［55］在一定程度上驱动了流域西北部低GPP的荒地向高GPP的草地的转变；二是不同土地覆被类型的GPP均在增加，在气候变暖推动下所有土地覆被类型中绿色植物的生存环境得到改善，碳固存能力均得到提高。与GPP相比，流域内ET显著增加的区域较少，主要分布在中下游地区，而且不同土地覆被类型间的ET差别小，这与流域的自然气候条件密切相关，无定河流域位于干旱半干旱区，当地的实际蒸散发量受到水分供应的限制［56］，此外，水分还是制约植被恢复的主要因素。因此，水分条件相对优渥的流域中下游地区植被恢复明显，同时也促进了该区域蒸散发的增加。WUE 耦合了GPP 与ET 的变化特征，在GPP 与ET 均以增加为主的情形下，流域绝大多数区域的WUE 仍呈显著增加趋势，说明该区域GPP的显著增加主导了WUE的增加，在植被恢复耗散水分增加的情况下，流域整体的水分利用效率是提高的。

在无定河流域所在的水资源不足的黄土高原区以及广大西北干旱半干旱区开展生态保护工程实践时，需要以水资源植被承载力为依据，权衡植被恢复的碳效益和水效应。从生态环境角度看，植树造林会增加区域固碳能力，改善生态环境；而从水资源耗散角度看，森林覆被增加会导致更强的蒸散发，水资源耗散显著增加，根据水量守恒定律，那么人类赖以利用的径流量就会减少。因此，耦合水碳机制来对区域植树造林的可行性和退耕退林还草效果的可持续性进行科学评估具有重要意义。在追求林草植被恢复的固碳功能的同时，有必要充分考量当地水分条件、水分情势变化和未来可持续性等，真正做到以水定绿、以水定碳，因地因水开展适宜性生态建设工作，切不可为了短期内提升碳汇能力而忽视水分条件的限制。

本文在无定河流域开展土地覆被的变化特征分析、空间分异因子探测、相关水碳变量变化分析的“三段式”综合研究模式可为黄土高原区、黄河流域及其他区域的地-水-碳耦合研究提供科学参考，对支撑区域水土资源保护与规划、实施水土保持工程、辅助实现“双碳”目标等工作也具有重要的实践参考意义。但本研究存在土地覆被空间分异机制分析不连续、相关水碳变量分析数量不全面等不足，今后将进一步获取更高精度、逐年连续的社会经济数据、拓展更多生态环境变量等来展开更深入的研究。

4 结论

本文通过分析1990s以来无定河流域土地覆被的基本特征、空间分异驱动力以及相关水碳变量（GPP、ET、WUE）的变化特征，得出以下结论：

（1）近30 a 来流域林草覆被得到良好恢复，退耕还林还草生态恢复工程的实施以及气候变化均是其背后重要的推力。作为流域内最主要的土地覆被类型，草地于2002 年左右出现明显的增长，在流域中下游及无定河沿岸地区大面积增加，其中下游地区草地主要由耕地转化而来，中上游河岸附近的草地主要由荒地恢复而来；而林地明显增加出现在2006年左右，增加地区也仅局限在流域最南部。

（2）流域土地覆被在空间上主要呈现由东南到西北逐渐变化的分布特征，人口密度和降水对该格局形成的影响力最大，GDP、气温、高程、坡度、坡向等因素也具有重要作用。此外，降水、气温等变率较大的自然因素对流域土地覆被空间分异的影响呈增强趋势。

（3）流域GPP、ET、WUE 的空间分布与土地覆被空间分异有较好的对应性，并且在林草恢复为主的覆被变化背景下，流域整体的GPP、ET、WUE明显增加。在流域绝大部分区域碳固存增加、水分利用效率提高的同时，伴随着中下游地区水分耗散加剧的问题，但总体而言，流域生态环境偏向转好。