基于水文连通性的流域输沙量评估——以堵河流域为例

张媛平, 刘成帅, 肖海兵, 王 剑, 史志华

(华中农业大学 资源与环境学院, 武汉 430070)

流域输沙会导致土壤养分流失、污染物扩散和河道淤积，对农业生产和水质健康具有重要影响[1-2]。针对流域输沙过程及其影响因素[3-5]，国内外学者已经开展了大量研究。例如，Ranzi等[6]评估了水库修建对越南Lo River流域输沙量的影响，指出修建水库显著削弱了流域出口的泥沙产量。Wang等[7]利用泥沙归因诊断，探究了黄河流域输沙量变化的原因，分析了降水、径流系数和含沙量等因素对黄河泥沙量减小的贡献。虽然以往研究已取得较大进展，但对于流域内部的输沙过程仍然缺乏了解。近年来，水文连通性的概念逐渐被研究者关注，它为解决流域输沙过程及其时空变异特征提供了可能，逐步成为研究热点[8]。

水文连通性是指泥沙、营养盐、生物体等各种物质以水为载体，在空间异质景观或斑块间进行迁移、传输或扩散的便利程度[9]。根据其功能定位，可分为结构连通性和功能连通性。结构连通性反映了流域景观单元的空间分布特征，可通过连通性指数定量表征。目前应用较为普遍的有汇流路径长度指数Flowlength[10]、方向性渗透指数DLI[11]、地形湿度指数TWI[12]和连通性指数IC等[13]。由于IC指数计算简单、数据要求低且精度高，应用最为广泛[14-15]。近年来，国内外学者尝试将IC指数与土壤侵蚀模型相结合，以评估流域输沙量[8]。例如，Zhao等[16]通过应用IC指数与泥沙输移比的经验关系，准确评估了黄土高原延河流域输沙量。Sean等[17]结合RUSLE模型与IC指数，绘制了流域输沙空间分布图，强调基于水文连通性评估流域输沙量的可行性。虽然，相关研究已证实基于水文连通性可有效评估流域输沙量，但不同流域自然背景差异较大，受人类影响的程度也不同。因此，模型的适用性仍需进一步研究。

堵河流域位于湖北省西北部汉江南岸，流域内山脉连绵、坡陡、谷深、地势落差大，是水土流失多发地区[18]。流域经历了耕作和再造林时期，输沙量年变化显著。因此，本文以堵河流域为研究对象，分析流域退耕前后土壤侵蚀模数和水文连通性的时空特征，并基于IC指数进一步估算流域输沙量，评估基于水文连通性估算流域输沙量的可行性。研究结果对于提高模型的适用性以及理解流域泥沙输移过程具有重要意义。

1 材料和方法

1.1 研究区概况

本研究以竹山站为界的堵河流域(31°30′—32°37′，109°11′N—110°25′E)为研究对象(图1)，流域面积为8973 km2，地势西南高东北低，以亚高山地貌为主。该地区属于北亚热带季风气候，年均气温14℃，年均降雨量810 mm，主要集中在6—10月。土壤类型主要以黄棕壤为主(占流域总面积70%以上)，有机质含量丰富，易被降水冲刷和淋溶[19]。域内植被覆盖率高，类型多样，主要有亚热带常绿阔叶林和常绿落叶阔叶混交林。年均径流量约为5.03×109m3，多年平均输沙量为3.62×106t，年平均侵蚀模数为4.04万 t/hm2。为有效防控水土流失，当地政府先后实施了一系列水土流失防治措施，如“丹江口水库水源区水土保持重点防治工程”等。结合堵河流域的历史背景，将研究时段划分为P1(1990—2000年，退耕前)和P2(2000—2010年，退耕后)两个时期。

图1 堵河流域地理位置

1.2 数据来源

流域1990年、1995年、2000年、2005年、2010年的遥感影像数据，来源于中国科学院计算机网络信息中心地理空间数据云平台(http:∥www.gscloud.cn)。数字高程模型(DEM)来源于国家基础地理信息中心，比例尺为1∶5万。土地利用主要为耕地、林地、草地、水域和居民用地，来源于中国科学院资源环境科学数据中心(http:∥www.resdc.cn)。降雨数据和水文数据来源于湖北省水利厅，包括流域内7个气象站日尺度降雨量数据及竹山水文站日径流泥沙数据。土壤可蚀性因子(K)源于国家科技基础条件平台—国家地球系统科学数据中心—地理资源分中心。本研究所用数据空间分辨率均为30 m，与国内外其他类似研究一致[4-6, 16]。

1.3 侵蚀模数估算

采用RUSLE模型[20]估算土壤侵蚀模数：

A=R·K·LS·C·P

式中：A为土壤侵蚀模数t/(hm2·a)；R为降雨侵蚀力因子[MJ·mm/(hm2·h·a)]，采用章文波等[21]提出的基于日降雨量的半月降雨侵蚀力模型计算；K为土壤可蚀性因子[t·h/(MJ·mm)]，依据国家地球系统科学数据中心绘制的全国30 m分辨率的土壤类型分布图进行计算；LS是坡长坡度因子(无量纲)，采用Wischmeier等[22]和刘宝元等[23]的方法计算；C为植被覆盖与管理因子(无量纲)，采用Durigon等[24]提出的基于归一化植被指数(NDVI)的方法计算；P为水土保持措施因子(无量纲)，采用Fu等[25]基于坡度的方法计算。计算公式如下：

P=0.2+0.03S

式中：S为坡度(%)。根据水利部颁布的《土壤侵蚀分级标准》(SL190-2007)，将流域土壤侵蚀划分为6个等级：微度、轻度、中度、强烈、极强烈和剧烈。

1.4 水文连通性评估

采用IC指数评估流域水文连通性：

式中：W为上坡贡献区的平均权重因子(无量纲)，用RUSLE模型中的C因子定量表征；S是上坡贡献区的平均坡度(m/m)；A是上坡贡献区的面积(m2)；di是沿着水流路径第i个单元的长度(m)；Wi是第i个单元的权重(无量纲)；Si是第i个单元的坡度(m/m)。IC范围为[-∞, +∞]，其值越大，表明水文连通性越大。为了更好地解释结果，采用Crema和Cavalli等[26]的方法将IC值分为4个等级：低、中低、中高和高连通性。

1.5 泥沙输移比(SDR)和输沙量估算

采用SDR分布式模型计算泥沙输移比，公式如下：

式中：SDRmax为流域内泥沙输移比最大值，取值范围为[0, 1]；IC0,i及KIC,i是校正参数，分别取0.5，2[27]，ICi为第i个单元的IC值。

流域输沙量计算如下：

式中：Q为输沙量(t/a)；Ei为第i个单元的侵蚀量(t/a)；SDRi为第i个单元的泥沙输移比(无量纲)；n是流域像元总数。

1.6 统计分析

利用SPSS 22.0软件进行输沙量估算值和实测值间的相关性分析；利用MATLAB计算均方根误差(RSME)和平均绝对误差(MAE)，评估实测输沙量和估算值间的相对偏差，确定模拟的准确性。

2 结果与分析

2.1 土壤侵蚀模数时空分布特征

流域土壤侵蚀模数呈现先增大后减小的趋势。2000年流域土壤侵蚀模数最大，为17.17 t/(hm2·a)；2010年最小，为3.48 t/(hm2·a)。P1时期，流域土壤侵蚀模数增加了215.05%，其中轻度侵蚀面积增加了41.86%，中度侵蚀和强烈以上侵蚀面积呈指数增加。P2时期，土壤侵蚀等级主要由轻度向微度侵蚀转变(表1，图2)。流域东北部、北部和河道附近土壤侵蚀模数大，西南、南部和中北部平原地区侵蚀模数小。总的来说，土壤侵蚀空间分布呈现边缘大，中间小的特点。微度侵蚀集中于流域中北部平原丘陵区，在坡度为15°～35°的区域零散分布；轻度侵蚀与中度侵蚀主要分布在流域北部和东北部坡度>25°的区域(图3)；强烈侵蚀及以上侵蚀主要以点和线状分布在流域东、北及河道景观附近。

表1 流域土壤侵蚀等级面积占比及其变化

图2 堵河流域土壤侵蚀等级分布

2.2 水文连通性的时空变化规律

由表2看出，2000年流域IC均值最高，为0.08；2010年最低，为-1.11。P1时期，IC均值增加了114.81%。其中，中高和高连通性面积分别增加了24.99%和16.37%。P2时期流域IC均值减小了105.56%，中高连通性及以上的面积减小了66.38%，而低和中低连通性面积分别增加了54.9%和12.0%(图4—5)。总的来说，1990—2010年流域IC值呈现先增大后减小的变化趋势。此外，IC值呈现“靠近河道大、远离河道小”的空间分布特征。流域低连通性主要分布在居民用地和林地，中低和中高连通性集中分布在坡度较陡的坡耕地与沟谷附近，高连通性则主要分布在流域坡度>25°的地区和靠近河道的区域。

表2 堵河流域IC值统计

2.3 流域输沙量估算

流域输沙量估算值呈现先增大后减小的年际变化趋势。2000年输沙量估算值最大，为9.42×106t/a，2010年最小，为1.27×106t/a。输沙量实测值与估算值呈现明显的线性关系(R2=0.92)(图6)，均方根误差(RSME)与平均绝对误差(MAE)表明估算输沙量与实测输沙量间差异较小(RSME=2.81；MAE=2.01)，模拟结果可接受。但是，流域输沙量估算值均大于实测值，其中2000年流域输沙量估算甚至达到实测值的2.5倍。研究结果表明基于IC指数与SDR的关系模型估算流域输沙量可能会在一定程度上造成数值偏高。

3 讨 论

研究结果表明退耕前(P1时期)堵河流域主要土壤侵蚀等级由微度逐渐转变为轻度，退耕后(P2时期)则主要由轻度转变为微度。以往研究指出土壤侵蚀年际变化主要受土地利用、降雨与地形等因素的影响[28]。自20世纪70—80年代“家庭联产承包责任制”在全国范围的推广，堵河流域大面积林地和草地被开垦为耕地，这一举措加剧了水土流失与土地退化，导致P1时期土壤侵蚀模数显著增加。90年代末国家“退耕还林”政策的施行，以及当地政府先后启动了“丹江口水库水源区水土保持重点防治工程”、“生态县”和“农发基金”等建设项目。流域林地和草地面积得以恢复，土壤侵蚀模数和侵蚀等级不断下降(2000年>2005年>2010年)。其次，降雨量是流域土壤侵蚀变化的关键影响因子[29]。2000年流域降雨量最大，其年降雨侵蚀力最高，因而其对表层土壤的剥蚀和运输能力最大。与之相反，2010年流域降雨量大于1995年和1990年，但其土壤侵蚀模数却最小，这可能是由于土地利用变化对土壤侵蚀的影响远大于降雨所造成的。此外，坡度[23]对土壤侵蚀也具有重要影响，主要表现为坡度越陡，侵蚀发生的潜力越大。以上结果与黄萱等[19]在该流域的研究结果一致。

图3 不同坡度土壤侵蚀等级面积占比

水文连通性与土壤侵蚀的年际变化趋势一致，主要受流域土地利用、植被和地形的影响。其中土地利用变化可能是影响水文连通性的关键因素。土地利用通过改变地表粗糙度、土壤有机质含量、土壤结构和入渗速率[30]，进而影响水文连通性。相比耕地，林地和草地具有较高的保水能力和泥沙拦截系数，水文连通性较低。因而，退耕前流域连通程度增大与耕地面积增加密切相关。其次，植被主要通过植被密度和群落结构影响水文连通性。通过冠层截留，枯落物蓄水、改善土壤结构延缓和阻碍径流，降低水文连通性[31]。袁亚男等[32]分析了昕水河流域的土地利用/覆被与水文连通性的关系，指出IC指数受土地利用/覆被的影响较为明显，不同土地利用/覆被下的水文连通性差异显著，而同一土地利用/植被覆盖条件下的水文连通性则变化较小。此外，流域地形(坡度等)会影响地表径流产生和发展。通常坡度越陡，地表径流产生的概率越大，水文连通性越强[33]。堵河流域东南部坡度普遍大于西北部，因而该区域水文连通性较大。

图4 堵河流域IC等级分布

注：横坐标分别表示1990年、2000年的连通性等级，纵坐标分别表示P1和P2时期连通性等级的面积变化百分比(分别基于1990年、2000年)。

图6 输沙量估算值和实测值

本研究基于IC指数进一步估算了流域输沙量，结果显示流域输沙量估算值与实测值具有显著的线性关系(R2=0.92)，但估算值要高于实测值。这表明该模型会高估流域实测输沙量。本文的研究结果与Zhao等[16]在黄土高原延河流域类似研究的结果相似，他们认为模型总体上可以较好的估算出低值年输沙量，而对高输沙量会产生错误估计。主要原因是水库、大坝的运行，会导致下游实测泥沙量出现异常变化，从而导致模型对泥沙量的高估或低估。本研究认为产生该结果的原因如下：(1) IC指数是基于流域土地利用和地形参数建立起来的，并未考虑降雨，入渗等的影响，因而导致估算值与实测值存在偏差；(2) 流域内水库和湖泊的建造，会影响河流流速，降低水体含沙量，进而拦截相当数量的泥沙，降低流域输沙量。本研究结果证明了基于IC指数评估流域输沙量的可行性，为该模型的广泛应用提供了依据。后续研究需充分考虑沟渠、水库、湖泊等人为工程设施对流域输沙的影响，改进模型以更准确评估流域输沙量。

4 结 论

流域土壤侵蚀模数和IC值均呈现先增加后减小的年际变化趋势。2000年，流域土壤侵蚀模数最大，为17.17 t/(hm2·a)，侵蚀强度以轻度侵蚀为主；平均IC值为0.18，水文连通性最强。2010年，流域土壤侵蚀模数最小，为3.48 t/(hm2·a)，侵蚀强度以微度侵蚀为主；IC均值为-1.11，连通程度最弱。此外，流域水文连通性呈现“靠近河道大，远离河道小”的空间分布特征。流域输沙量估算值与实测值均呈现先增大后减小的年际变化规律。RSME(2.81)和MAE(2.01)表明两者差异较小，模拟结果可接受。但由于模型尚未考虑水库等对流域输沙的影响，导致输沙量估算值大于实测值。因此，后续研究中应进一步改进模型以准确评估流域输沙量。