黄河上游地区土壤保持服务时空变化及归因

2022-08-16吴成永曹广超陈克龙巴丁求英刘宝康王润科石亚亚鄂崇毅

水土保持学报 2022年4期

吴成永，曹广超，陈克龙，巴丁求英，刘宝康，王润科，石亚亚，鄂崇毅

(1.天水师范学院资源与环境工程学院，甘肃天水 741001；2.青藏高原地表过程与生态保育教育部重点实验室，青海省自然地理与环境过程重点实验室，西宁 810008；3.高原科学与可持续发展研究院，西宁 810008)

土壤保持服务是生态系统提供的一项重要调节服务，包括潜在和实际供给服务。潜在供给服务是由于受到植被覆盖等生态系统结构和组分因素的影响，自然生态系统具有控制侵蚀和拦截泥沙等的潜在供给能力。当这种能力给人类带来粮食生产和清洁水质等效益时，潜在供给服务即变为生态系统为人类提供的实际土壤保持服务。气候变化，特别是温度和降雨模式的变化对生态、水文和土壤过程及土壤有机物、植物和作物等造成巨大影响，进而影响粮食、水资源等的供给能力以及区域社会经济发展和人类福祉。因此，气候变化背景下土壤保持服务相关研究受到了国际社会的广泛关注。世界气候研究计划(WCRP)将土壤保持列入重要研究内容，一些国际会议探讨了气候变化与土壤保持等相关议题，如2008年第15届国际土壤保持大会、2017年联合国粮农组织大会第40届会议、2021年土地退化、土壤保持与可持续发展国际论坛(LASOSU)。

黄河流域横跨青藏高原生态屏障、黄土高原生态屏障及北方防沙带“两屏一带”，是我国水土流失防治的重点区域。经过多年持续不懈的治理，其水土流失严重局面得到初步扭转，但监测表明，2019年水土流失面积达26.42万km，水土流失依然是构建生态安全屏障的短板。黄河上游地区，横跨干旱、半干旱、半湿润、湿润/半湿润4个气候区，作为黄河流域生态保护的关键地区，是一个发展问题和生态问题交织在一起的社会—经济—自然复合生态系统，其生态系统土壤保持服务是防止水土流失和推动高质量发展的重要保障。那么，黄河上游地区土壤保持服务基本特征及影响因素是什么？土壤保持服务在不同气候区的变化趋势如何，其主导影响因素是什么？这些科学问题的回答，将有助于深入理解土壤保持服务的变化机理，以及制定应对气候变化的策略措施。目前该区土壤保持的相关研究主要在黄河源、三江源、洮河流域、黄土高原、宁夏内蒙古河段和北方农牧交错带及其局部地区开展，全境范围内的完整性研究较少，不足以形成关于土壤保持服务的系统性认识。同时，这些研究多数使用全局性的回归模型、相关分析、趋势分析和参数控制等方法得到整个研究区的影响因子，对于宏观制定水土保持政策具有一定的指导意义，但难以提出不同地域、气候区局部性的差异化参考建议。地理加权回归(geographically weighted regression, GWR)模型，是一种探索数据空间异质性的局部技术，将空间位置信息纳入到回归方程，可以分析影响因子在空间上的差异性，其分析结果对于实施基于地理位置的差异性水土保持措施具有指导意义。

鉴于此，本文运用广泛应用的土壤保持服务评估模型RUSLE、GWR模型和趋势线等方法，明晰黄河上游地区土壤保持服务时空动态特征，分析土壤保持服务在各气候区的变化趋势，在此基础上，结合本文提出的因子影响力和影响强度等指标方法，尝试探寻土壤保持服务变化的影响因子并量化影响度，以期为该地区应对气候变化制定有针对性的分区分类宏观策略和因地制宜的具体措施，为提升土壤保持服务提供科学依据，助力生态保护和高质量发展。

1 材料与方法

1.1 黄河上游地区

黄河上游是指从河源至内蒙古自治区托克托县河口村(图1)。薛娴等以DEM数据运用河道烧录法与河道标量法划定了黄河上游流域；刘小鹏等从流域经济的角度出发，以省级行政区为基本单元，界定了黄河上游地区；杨永春等从地市行政区的角度界定了黄河上游地区。为构筑生态安全屏障和推进生态文明建设，我国以县级行政区为基本单元实施生态管理，基于此，本文将与黄河上游流域边界相交的县(市、区、旗)纳入黄河上游地区，共涉及137个县级行政区，其面积为97.94万km，占全国国土面积的10.21%。

注：黄河上游流域边界，来源于国家青藏高原科学数据中心；生态地理分区边界，来源于中国科学院资源环境数据中心，其中IIC3为内蒙古高平原东部，半干旱区；IID1为内蒙古高平原西部及河套，干旱区；IID2为阿拉善及河西走廊，干旱区；HIID1为柴达木盆地，干旱区；IIIC1为晋中陕北甘东高原丘陵，半干旱区；HIIC1为青东祁连山地，半干旱区；HIC1为青南高原宽谷，半干旱区；HIB1为果洛那曲丘状高原，半湿润区；HIIA/B1为川西藏东高山深谷，湿润/半湿润区。

1.2 数据来源与处理

本文使用的主要数据见表1，其中归一化植被指数(NDVI)、降水量和土地覆被产品数据时段为2001—2015年。所有数据在ArcGIS 10软件平台下统一处理为Albers等积投影。

表1 主要数据与来源

2 研究方法

2.1 土壤保持服务评估模型

土壤保持服务通常采用土壤保持量进行评估，基于RUSLE方程的土壤保持服务评估模型为：

=-=×××(1-×)

(1)

式中：为土壤保持量(t/(hm·a))，由潜在侵蚀()与实际侵蚀()之差决定；为降雨侵蚀力因子((MJ·mm)/(hm·h·a))；为土壤可蚀性因子((t·hm·h)/(hm·MJ·mm))；为地形因子；为植被覆盖因子；为水土保持措施因子。本文中，使用修正的降雨侵蚀力公式；与采用中国环境保护部和发展改革委员会联合发布的《生态保护红线划定指南》(环办生态[2017]48号)推荐的计算办法；采用王丽霞等的计算办法；采用地类赋值法。

2.2 趋势线法

使用趋势线法，即以土壤保持量的变化斜率定量分析土壤保持服务的年际变化。

(2)

式中：为年序号；为年数；为土壤保持量变化斜率，代表土壤保持服务变化趋势。当>0，表明研究期间土壤保持量上升，土壤保持服务提升；当<0，表明土壤保持量降低，土壤保持服务下降；当=0，则土壤保持服务保持不变。

2.3 地理加权回归模型

GWR模型考虑数据的空间位置，因逐点进行参数估计，其分析结果是局部性的，而非全局性的，其模型为：

(3)

式中：、、(，)和(，)分别为空间位置处的因变量、解释变量、地理坐标和回归系数；(，)为截距；为解释变量的个数；为误差修正项，即残差。本文以土壤保持量为因变量，降水、海拔、坡度和植被覆盖度为解释变量在ArcGIS 10软件中进行GWR分析。

2.4 因子影响度指标及其计算

借鉴已有研究，本文尝试提出影响强度、影响率和影响力3个指标以度量影响因子对土壤保持服务的影响。

2.4.1 影响强度 GWR分析的回归系数反映了解释变量与因变量之间的关系强度，如果关系较强，则系数相对较大，如果关系较弱，则系数接近于零。因此，回归系数的绝对值可以定义为因子的影响强度(influence intensity,)：

=|(,)|

(4)

式中：和(，)分别为因子在空间位置处的影响强度和GWR分析的回归系数。

2.4.2 影响率残差可用于确定模型的拟合程度，较小的残差表明模型拟合效果较好。ArcGIS 10平台的GWR分析结果中，标准化残差绝对值的范围为0～2.5，本文将其<1.5的模型视为较优模型，并将对应的解释变量个数()作为因子对土壤保持服务的影响数量，于是，因子的影响率(influence ratio,)可定义为：

=()×100

(5)

式中：为因子对土壤保持服务的影响率；为解释变量的总数，即参与GWR分析的因子总个数；为因子的影响数量。

2.4.3 影响力影响因子对土壤保持服务的影响力大小，一是取决于影响强度，二是取决于影响率。于是，因子影响力(influence level, IL)可定义为：

=×

(6)

式中：为因子对土壤保持服务的影响力；为因子对土壤保持服务的影响率；为(因子影响数量)个因子的影响强度均值。

3 结果与分析

3.1 土壤保持服务基本特征

3.1.1 土壤保持服务分区将2001—2015年间(称“研究期间”)土壤保持量变化斜率，在ArcGIS 10软件平台下，使用自然断裂法分为5类，据此将土壤保持服务变化划分为明显下降区、略微下降区、基本不变区、略微上升区和明显上升区(图2)。研究期间，基本不变、下降和上升区面积分别为64.08万，28.06万，5.54万km(表2)，面积比约为66∶29∶5。

表2 土壤保持服务分区及其特征

图2 土壤保持服务分区

下降区，主要位于黄河上游东北部的干旱、半干旱气候区。其中，明显下降区集中分布在陇东黄土高原丘陵沟壑区的半干旱气候区，县域上，环县、彭阳县、原州区东部土壤保持服务下降明显；略微下降区，集中分布在黄土高原丘陵区西部、内蒙古高平原西部及河套地区、西北内流河片区的阿拉善及河西走廊东南部交汇地带的干旱气候区，大致是积石峡以东15 km处至河口村以西50 km处，黄河两岸宽约150 km的带状区域。县域上，乌拉特中旗、杭锦旗、杭锦后旗、鄂托克旗、鄂托克前旗、盐池县、阿拉善左旗东南部、靖远县、白银区、皋兰县、永登县东部、景泰县、古浪县、永昌县东部、凉州区、白银区、中宁县、永宁县、和政县、民勤县南部、兴庆区、金凤区和平川区等县(市、区、旗)土壤保持服务呈略微下降。

上升区，主要位于黄河上游西部和南部的半干旱、湿润/半湿润和半湿润气候区。其中，明显上升区，主要分布在半干旱的青东祁连山地和湿润/半湿润的川西藏东高山深谷，特别是祁连山、疏勒南山、岷山北段等山地，这与刘海等的研究结果类似，即土壤保持量呈现山区高、平原低的分布格局。县域上，肃南裕固族自治县西北中部、祁连县西北中部、宕昌县中南部、松潘县土壤保持服务上升明显；略微上升区，主要分布在半湿润的果洛那曲丘状高原西部、湿润/半湿润的松潘高原北部，县域上主要是称多县南部、久治县、班玛县和迭部县。

3.1.2 各分区土壤保持服务年际变化特征明显下降区，土壤保持量从2001年的4 116.92 t/(hm·a)，波动下降至2015年的1 349.34 t/(hm·a)，年均下降184.51 t/(hm·a)(图3a)；略微下降区，土壤保持量从2001年的599.81 t/(hm·a)，波动下降至2015年的345.93 t/(hm·a)，年均下降16.92 t/(hm·a)(图3b)；明显上升区，土壤保持量从2001年的11 015.00 t/(hm·a)，波动增加至2015年的11 926.65 t/(hm·a)，年均增加60.78 t/(hm·a)(图3c)；略微上升区，土壤保持量从2001年的4 998.19 t/(hm·a)，波动增加至2015年的5 332.55 t/(hm·a)，年均增加22.29 t/(hm·a)(图3d)。

图3 各分区土壤保持量年际变化

3.2 土壤保持服务变化归因分析

结合前人选用的影响因子，考虑数据的可得性，本文选取降水、海拔、坡度、植被覆盖作为土壤保持服务的影响因子。将研究期间4个因子的均值归一化至0～1以消除不同量纲间的差异，依据影响因子及其指标特征(表3)，探寻各分区的影响因子，进而分析土壤保持服务变化的原因。

表3 各分区土壤保持服务影响因子及其指标特征

3.2.1 下降区土壤保持服务影响因子明显和略微下降区，各因子对土壤保持服务的影响强度和影响力排序一致，从大到小依次为降水>海拔>植被覆盖度>坡度(表3)。降水是影响下降区土壤保持服务的主导因子，其像元尺度的影响强度见图4a。

Tang等研究发现，黄土高原1965—2014年间降水量减少。研究期间，土壤保持服务下降连片区(积石峡至河口村的黄河沿岸)，其年均降水量呈略微下降态势(图5a)。地处陇东黄土高原的甘肃省环县，研究期间的降水量亦呈现略微下降特征，其土壤保持量下降最为明显(图5b)。此外，信忠保等研究发现，黄土高原地区在降水减少的过程中，侵蚀性降水减少更为剧烈。降水量的减少，特别是侵蚀性降水的减少，降低降雨侵蚀力，致使土壤保持服务下降。

3.2.2 上升区土壤保持服务影响因子依据各因子影响强度和影响力大小(表3)，明显、略微上升区土壤保持服务的主导因子分别是植被覆盖和降水，其像元尺度的影响强度见图4b和图4c。

图4 各分区主导因子影响强度空间分布

明显上升区的植被覆盖度从2001年的65.12%波动增加至2015年的72.84%(图5c)。地处波密—川西东北部的四川省松潘县，其土壤保持服务上升最为明显，植被覆盖度从2001年的78.73%波动增加至2015年的86.35%(图5d)。杨达等研究也发现，2001—2018年间波密—川西东北部的NDVI(与植被覆盖度正相关)显著增加。植被覆盖度上升能够有效遏制土壤侵蚀的发生，进而提升土壤保持服务。

略微上升区降水量从2001年的532.88 mm缓慢波动增加至2015年的533.72 mm(图5e)。温煜华等研究发现，1961—2017年间祁连山平均降水和极端降水强度增加。研究期间，地处祁连山地区的肃南裕固族自治县，其西北部土壤保持量呈略微上升态势，年均降水也呈略微上升态势，从2001年的272.28 mm波动增加至2015年的325.68 mm(图5f)。降水量的增加，致使降雨侵蚀力增强，土壤保持服务提升。

图5 各分区及其典型县域主导因子年际变化

4 讨论

本文通过GWR模型及其因子影响度指标得到影响黄河上游地区土壤保持服务变化的主导因子为降水和植被覆盖。分析研究区局部地区土壤保持服务/土壤侵蚀相关研究案例(表4)可知，虽然其影响因子探寻方法有所不同，但影响因子与本文探寻的基本一致，表明本文影响因子探寻方法合理，探寻结果可靠。

表4 研究区局部地区土壤保持服务/土壤侵蚀研究案例

黄河上游东北部的干旱、半干旱气候区，特别是积石峡至河口村的黄河沿岸地区，影响土壤保持服务变化的主导因子是降水。降水量的变化致使水文、土壤和生态过程及植被覆盖等生态系统结构功能和组份发生相应的变化，进而影响生态系统土壤保持服务。气候模式预测认为，未来研究区降水变化存在较大的不确定性，魏洁等基于BCC-CSM1.1情景分析认为，2011—2050年间黄河上游降水增多4.31%～5.74%，而Yang等研究发现，在RCP4.5情景下，2021—2050年间黄河流域降水减少，减少率为2.99 mm/10年。青东祁连山地半干旱气候区和川西藏东高山深谷湿润/半湿润气候区，特别是祁连山、疏勒南山和岷山北段等山地，影响土壤保持服务变化的主导因子是植被覆盖。较高的植被覆盖能有效遏制土壤侵蚀的发生，从而保障较高的土壤保持供给服务。研究期间，上述植被覆盖为主导因子的地区，其植被覆盖度呈动态上升趋势，相关研究也表明，近年来黄河上游一些地区的植被覆盖得到了显著改善。考虑到降水是影响黄河流域植被覆盖的决定性因素，本文认为面对未来降水的不确定性，提高植被覆盖度是提升黄河上游地区土壤保持服务的基本策略。

结合GWR模型，本文尝试定义了因子“影响强度”指标。该指标反映了土壤保持服务与影响因子之间的关系强度，若某地理位置上某因子的影响强度值较大，则表明该因子对该地理位置处的土壤保持服务影响也较强。因此，可根据主导因子影响强度空间分布图，因地制宜地制定具体措施提升土壤保持服务，如依据图4，植被覆盖因子影响强度值较大的地理位置处，应加强植被保护等措施以进一步提高植被覆盖度，而降水因子影响强度值较大的位置处，应实施坡沟治理、改变耕作方式等保水措施以增加该位置降水的利用率。

本文探寻得到的主导因子能够解释研究区绝大多数地域的土壤保持服务变化。但位于黄河源的曲麻莱县西部，研究期间其土壤保持服务呈现下降趋势，而主导因子降水则是略微增加，这难以阐释降水与土壤保持服务的关系。曲麻莱县作为三江源国家级自然保护区和可可西里国家级自然保护区的核心区域，生态保护修复项目的实施，以及当地居民的生态环保理念等均影响土壤保持服务。事实上，土壤保持服务受到自然环境因子以及人类活动、社会经济和政治文化等人文因子直接或间接的影响，且各因子之间彼此联系、互相促进与制约，通过复杂的物理、生物化学过程影响土壤保持服务。由于数据源的限制，本文GWR分析中仅考虑降水、海拔、坡度和植被覆盖自然环境因子，土壤性状(类型、机械组成等)、雨强等也是决定区域土壤保持服务/土壤侵蚀自然环境的主要因子。未来收集或采集相关自然环境因子数据，纳入人文因子，并选取合适的方法量化因子间的相互影响是今后的研究方向。