1954—2016年疏勒河上游径流变化特征及影响因素

孙美平， 张 磊， 姚晓军， 彭莉斌， 张 浩， 牛舒婷

（1.西北师范大学地理与环境科学学院，甘肃兰州730070； 2.中国科学院西北生态环境资源研究院，甘肃兰州730000）

0 引言

河川径流是重要的地表水资源，尤其是在西北干旱区，缺水与生态问题相互交织，是制约区域社会经济可持续发展的主要瓶颈［1-2］。1987 年以来西北干旱区呈现明显暖湿化趋势，致使冰川快速消融、山区降水增加，水资源时空分布和水循环过程发生了显著变化［3-5］。IPCC 第六次评估报告将气候变化对水文过程和情势的影响作为第一工作组和第二工作组的焦点领域，并重点关注气候变化背景下流域尺度上的水资源变化及其影响［6］。

疏勒河流域位于祁连山西段，其水资源是维持下游绿洲及城镇稳定和繁荣的重要保障。由于特殊的地理位置，疏勒河流域在我国西部生态安全和经济发展中具有十分重要的战略地位。在气候变暖背景下，西北干旱区冰川普遍退缩［7］，其融水补给的内流水系径流变化趋势以及对气候变化的响应规律一直备受关注。早在20 世纪90 年代初，杨针娘［8］综合应用冰川融水径流模数法、流量与气温关系法、对比观测实验法等对甘肃河西地区主要水文站点的径流变化及冰川融水补给进行了较系统分析，估算河西内流水系冰川融水径流总量约9.99×108m3，占河川径流量的3.3%。目前，关于疏勒河流域径流变化特征及规律已形成一定认识。孙栋元等［9］利用线性倾向方法分析了疏勒河干流不同时间尺度的径流量变化和突变特征，发现流域径流集中于5—9月，季节变化明显。李计生等［10］应用坎德尔秩次相关法分析了1956—2013 年疏勒河流域出山径流变化规律，并预测了其变化趋势，认为2014—2018 年疏勒河干流径流偏丰。李浩杰［11］基于多源遥感数据分析了疏勒河上游冰川面积和高程变化，发现2000—2015年期间该区域冰川呈退缩减薄趋 势，其 中 冰 川 面 积 减 少（57.6±2.68）km2（11.7%），冰川表面高程降低（2.58±0.6）m，冰储量减 少（1.271±0.307）km3。Zhang 等［12］应 用VICCAS 模型模拟了1971—2012 年期间疏勒河冰川融水对径流贡献，研究表明，流域冰川面积覆盖占比仅约4%，但冰川融水对河流径流的平均贡献约为23.6%。Wu 等［13］应用SPHY 模型模拟了疏勒河上游1971—2020年水文过程变化，得到冰川径流贡献为28%，随着冰川面积的减少和温度的升高，冰川和积雪径流对总径流的年贡献率呈下降趋势。李洪源等［14］应用寒区水文模型亦定量模拟了疏勒河流域径流过程，得出冰川融水径流占流域总径流量的30.5%，且1971—2015 年径流增加近70%。尽管学者们在不同时期分析和总结了疏勒河流域径流与冰川变化特征，但受限于模型方法、时空尺度以及可获取的资料，使得对该流域径流变化的系统认识仍有待加强。

本文基于1954—2016 年昌马堡水文站日径流资料、疏勒河流域周边气象台站和探空站观测资料以及中国两次冰川编目数据，以径流变化特征及对气候变化的响应为切入点，系统分析疏勒河上游径流长时间序列变化趋势、周期特征及关键气候要素，以期揭示气候-冰川变化背景下的冰川融水补给型河流径流变化规律，从而为西北干旱区水资源配置优化和科学管理提供参考依据。

1 研究区概况

疏勒河是河西内流水系中仅次于黑河的第二大河流，发源于祁连山脉西段托来南山与疏勒南山之间的沙果林那穆吉木岭，在青海省境内自东南向西北流经高山峡谷地区，进入甘肃省境内折向北流入昌马盆地和河西走廊冲积洪积平原（图1）。昌马堡水文站以上为疏勒河上游，集水区面积为10 946 km2，海拔介于1 990～5 826 m，平均海拔3 885 m，最高峰为团结峰（又名岗则吾结，5 826 m）。流域年平均气温为5.3°C，年平均降水量为93.1 mm，属大陆荒漠干旱型气候［15］。上游祁连山区降水较丰沛，现代冰川发育。根据中国第二次冰川编目，疏勒河流域是河西内流水系拥有冰川资源最多的二级流域，共有冰川660 条，冰川面积和冰储量分别为509.87 km2和29.66 km3，各占祁连山冰川相应总量的24.59%、31.91%和35.11%［16］。需要说明的是，昌马堡水文站以上集水区并不包括祁连山面积最大的老虎沟12号冰川（又名透明梦柯冰川，编码为5Y448D0012，面积20.42 km2）所在的老虎沟流域。除冰川外，流域上游下垫面多为寒漠和裸地。

图1 疏勒河上游及冰川、水文站与气象站分布Fig. 1 The upper Shule River and distribution of glaciers，hydrological and meteorological stations

2 数据与方法

2.1 数据来源

本文所用数据包括疏勒河上游昌马堡水文站及流域周边气象台站、探空站长时间观测资料（表1）。其中，昌马堡水文站1954—2016年日径流数据为甘肃省水文水资源局的实测整编资料，用于分析疏勒河上游径流变化特征；流域周边气象台站包括托勒、酒泉、玉门镇、安西、敦煌、冷湖和大柴旦7 个国家标准气象站点，观测要素为1954—2016年降水量、平均气温、潜在蒸发、相对湿度、风速等日数据，用于分析疏勒河上游径流变化特征及与气候变化间的关系；疏勒河流域周边探空站极少，仅有马鬃山、敦煌和酒泉3 个探空站，考虑到2010 年后探空观测仪器更新，本文选用1971—2010 年夏季（6—8月）700 hPa、500 hPa 和400 hPa 等压面的气温、位势高度等日数据，用于计算0°C 层高度及对汛期疏勒河上游径流的影响。上述7个国家标准气象站点和3 个探空站数据从国家气象科学数据中心（http：//data.cma.cn）获得。此外，疏勒河上游冰川数据来自中国第一次冰川编目［17］和中国第二次冰川编目［18］，前者数据源为1956—1963 年的1∶5 万航测图并经过人工修订，后者数据源为2005—2010 年的Landsat TM 遥感影像。

表1 疏勒河上游水文站、周边气象站和探空站信息Table 1 Information of hydrological station in upper Shule River，meteorological and radiosonde stations

2.2 研究方法

2.2.1 线性倾向估计与滑动平均

为揭示1954—2016 年疏勒河上游径流变化趋势及年代际波动特征，采用线性倾向估计和10 a 滑动平均方法，计算径流随时间变化相关系数并检验其显著性水平，从而判断径流变化趋势是否显著，线性倾向估计和滑动平均计算方法见文献［19］。

2.2.2 Mann-Kendall趋势检验

Mann-Kendall（M-K）趋势检验是世界气象组织（World Meteorological Organization，WMO）推荐并广泛应用于气象、水文变化趋势与突变检验的一种非参数统计方法，其优点在于样本不需服从一定的分布，且不受少数异常值的干扰，定量化程度高。在M-K 趋势检验方法中，其趋势分析和突变检测通过统计量UFk和UBk进行判别［20-21］。在给定的显著性水平α下，当M-K 检验曲线中UFk线在临界线内变动，表明变化趋势和突变不明显。当UFk值大于0意味着序列呈上升趋势，反之为下降趋势；当UFk或UBk值超过临界值时，表明上升或下降趋势显著。若UFk和UBk两条曲线出现交点，且交点在临界线内，则该交点对应的时刻即为突变开始年。

2.2.3 小波分析与经验模态分解

小波分析（wavelet）亦称多分辨分析，可以将信号和图像分解成交织在一起的多尺度成分，得到时间系列相关参数（如平均值、频率、振幅）的变化特征，从而能够不断聚焦到所研究对象的微小细节。相较于小波分析，经验模态分解法（empirical mode decomposition，EMD）能更直观、直接地得到有物理意义的频率，其实质是通过特征时间尺度识别数据资料中所内含的所有振动模态［22］。对于给定的一时间序列X（t），EMD 方法通过不断剔除极大值和极小值连接上下包络的均值将原数据信号分解为：

式中：cj（t）为一个本征模函数（intrinsic mode function，IMF）分量；Rn（t）为残余分量，一般为数据序列的平均趋势。EMD方法的关键是经验模式分解，使复杂资料分解为有限个IMF，所分解出来的各IMF分量不断从高频到低频进行提取，最终得到一个频率接近为0的残余分量项。小波分析和经验模态分解两种方法常用于水文序列的时频结构分析，以探求径流阶段性和周期性的变化规律，二者具体算法见文献［23-24］。

2.2.4 分层回归

分层回归（hierarchical multiple regression）又称层次回归，其本质是对所建立的两个及两个以上回归模型进行比较，根据两个模型所解释的变异量差异判定自变量的重要性，且差异量可用统计显著性来估计和检验。通常情况下，采用容忍度和方差膨胀因子进行自变量之间的多重共线性诊断，诊断标准为容忍度小于0.1 且方差膨胀因子大于10 时，则认为存在严重的多重共线性。在分层回归方法中，决定性系数（R2）是一个重要指标，可反映自变量解释因变量变异的程度，从其改变量可以检验增加的自变量是否具有统计学意义［25］。

3 疏勒河上游径流长时间序列观测事实

3.1 径流变化趋势

图2为1954—2016 年昌马堡水文站日径流序列和年内径流分布及变化趋势。从日尺度径流序列［图2（a）］来看，疏勒河上游径流量随着时间推移呈增大趋势，尤其是2000 年后径流量增加明显，较2000 年前增加50.5%。1954—2016 期间，昌马堡水文站年径流最大值与最小值之比为4.18，年径流变差系数仅0.29，表明疏勒河上游径流量年际变化小，即径流较为稳定，这与流域内冰雪融水补给调蓄有关。然而，疏勒河上游径流年内分配极不均衡［图2（b）］，汛期（5—9月）来水量占全年径流总量的57%～79%，月径流量的极值分别达到0.09×108m3和5.95×108m3。从年内径流变化趋势［图2（c）］看，各月径流自有观测以来均呈增加趋势，且主要出现在冰雪消融期，尤其是7—8 月增加幅度最大，达0.24×108m3·（10a）-1；而1—4 月和10—12 月径流年内变化微弱，其中1—3 月径流增加幅度仅0.03×108m3·（10a）-1，仅是7—8月径流增幅的1/8。

图2 昌马堡水文站日径流序列及年内径流分布与变化Fig. 2 Daily series，intra-annual distribution and variation of runoff at Changmabao hydrological station

3.2 径流变化特征检测

1954—2016年昌马堡水文站实测年径流量呈明显的上升趋势［图3（a）］，变化速率为1.00×108m3·（10a）-1，且通过0.01 显著性检验。从10 a 滑动平均曲线来看，1970s—1990s中期疏勒河上游来水偏少，年平均径流量为9.84×108m3；自21 世纪初以来进入丰水期，并维持在高位波动，年平均径流量为12.60×108m3，较2000 年 前 增 加 近28%。Mann-Kendall检验统计［图3（b）］同样表明，疏勒河上游年径流量自1960年以来呈上升趋势，且该趋势在1998年超过显著性水平0.05临界线，1999年径流发生突变，这与径流变化曲线一致。无论是汛期还是非汛期，1954—2016 年昌马堡水文站观测到的径流均呈显著上升趋势［图3（c）、3（e）］，变化率分别为0.73×108m3·（10a）-1和0.31×108m3·（10a）-1，其中汛期多年平均径流量为7.39×108m3，远高于非汛期多年平均径流量（2.44×108m3）。与年径流突变年份相同，汛期径流突变亦出现在1999年［图3（d）］，而非汛期径流量不存在明显的突变年份［图3（f）］。

图3 昌马堡水文站径流变化及Mann-Kendall趋势检验Fig.3 Runoff variation and Mann-Kendall trend test at Changmabao hydrological station

3.3 径流周期分析

由小波实部系数［图4（a）］和小波方差［图4（b）］可知，1954—2016 年期间疏勒河上游年径流量呈现3 个时间尺度的周期变化规律：10～16 a 尺度上表现为枯、丰交替的准3 次振荡，且具有全域性；6～10 a尺度上存在准4次丰、枯振荡；2～6 a尺度上存在一定交替振荡，但较为紊乱。从年径流小波方差判断，12～15 a尺度上振荡能量最强，7～9 a尺度上振荡能量较弱，即15 a为第一主周期，7 a为第二主周期。疏勒河上游汛期径流量的小波实部系数表现出与年径流量相一致的周期情况，分别为准16 a、8 a 和5 a；非汛期小波方差在小波影响锥有效情况下存在2个明显的振荡，分别为12～13 a和6～7 a。

对1954—2016 年昌马堡水文站年径流量时间序列进一步分解，可得到7 个具有不同波动周期的IMF 分量和一个残余趋势项Residual［图4（c）］，表明疏勒河上游年径流量具有多时间尺度的波动特性，且总体表现为随着阶数增加，IMF分量振幅减小而周期增长。其中，IMF1、IMF2 为高频分量，IMF1波动频率最高，保持着3 a左右的波动周期，而IMF2分量对应准6 a 周期，二者在1980s 和1990s 波动幅度小且稳定，1980 年前和2000 年后波动幅度较大；IMF3 和IMF4 为中频分量，对应对波动周期分别为9～10 a、15～20 a；IMF5～IMF7 为低频分量，三个分量对应周期分别为准40 a、50 a 和60 a。残差项（Residual）表明1954—2016 年疏勒河上游年径流量亦呈增加趋势。

图4 昌马堡水文站年径流Morlet小波分析及经验模态分解［（a）为小波实部图，小波影响椎（图中弧线）以内的小波方差对应周期为有效周期；（b）为小波方差图，红色实线为方差曲线，灰色虚线为傅立叶谱，黑色虚线为95%置信水平；（c）为年径流量经验模态分解］Fig.4 Morlet wavelet analysis and empirical modal decomposition of annual runoff at Changmabao hydrological station（subfigure（a）is the real part of Morlet wavelet and the curve corresponding to wavelet variance in centrum represents the effective period；subfigure（b）is wavelet variance，and the red solid line，gray dashed line and black dashed line denote variance，Fourier spectrum and confidence level with 95%，respectively；subfigure（c）is empirical mode decomposition of annual runoff）

4 疏勒河上游径流变化影响因素

4.1 气候变化对径流变化的影响

昌马堡水文站以上人类活动较少，其径流变化可真实反映疏勒河上游水文自然变化过程。图5反映了疏勒河上游及周边7 个气象站1956—2016 年年均气温、正积温、降水与潜在蒸发量变化。显然，1956—2016 年期间疏勒河上游流域年平均气温呈显著升高趋势，增温速率为0.29°C·（10a）-1；与之类似，年正积温也总体呈现显著增大趋势，尤其是1990s 末以来，年正积温超过3 200 °C。同期，年降水量表现为在波动中增加趋势，变化速率为4.9 mm·（10a）-1，且通过0.001 显著性检验。与上述三个要素变化趋势不同的是，疏勒河上游潜在蒸发量在研究时段内整体表现为下降趋势，这可能与这一时期的降水量增加和相对湿度增大（0.22%·（10a）-1，P＞0.05）导致饱和水汽压差减小抑制地表蒸发有关，也与同时期的风速显著减小（-0.09 m·s-1·（10a）-1，P＜0.001）有直接关系。总体上，疏勒河上游流域在过去60 a间经历了“暖干→冷湿→暖湿”3个明显的气候变化过程，分别对应1956—1974年、1975—1999年和2000—2016年。

图5 疏勒河上游年平均气温、正积温、年降水与潜在蒸发量变化Fig. 5 Changes in annual mean temperature，cumulative positive temperature，annual precipitation and potential evapotranspiration in upper Shule River

为辨析影响疏勒河上游径流的关键气候要素，采用分层回归模型检验上述气象要素对昌马堡水文站径流的影响程度及重要性（表2）。若仅考虑年降水影响（模型1），决定性系数R2为0.442，表明降水单一要素仅能解释径流贡献的44%左右；当增加正积温要素（模型2）后，二者对年径流均有显著影响，且通过0.001显著性水平检验，决定性系数R2升至0.807，表明疏勒河上游80%以上径流可由降水和正积温来解释。根据标准化系数，正积温对径流影响强于降水，这与疏勒河上游流域发育有面积较大的冰川、积雪和冻土密切相关，因正积温这一热量指标能较好地反映冰冻圈要素的消融强度，即冰川、积雪融水和冻土冻融补给大于大气降水影响。当增加年平均气温（模型3）和年蒸发（模型4）两个要素时，二者均未通过显著性检验，且两个模型中的决定性系数较模型2 变化微小，这意味着年降水和年正积温是决定疏勒河上游径流变化的关键气候因子，尤其是1990s 末以来二者的同步上升导致昌马堡水文站径流快速增大。

表2 疏勒河上游径流与各气象要素的分层回归模型Table 2 Hierarchical multiple regression of runoff with each meteorological element in upper Shule River

4.2. ℃层高度对汛期径流影响

冰川融水是疏勒河上游径流的重要补给源。伴随气温上升，冰川消融区范围和强度不断扩大和加剧，其融水对河流径流补给持续加强继而进入汛期。冰川消融与0 ℃层高度变化存在密切关系，后者又称冰冻层高度、消融层高度，强烈影响着高海拔地区的冻融过程，对冰雪消融有明确的物理意义［26］。根据马鬃山、酒泉和敦煌三个探空站700 hPa、500 hPa 和400 hPa 等压面的气温（T）和位势高度（H）数据，二者分别存在以下函数关系：T=-0.00702H+32.63（马鬃山），T=-0.00721H+35.15（酒泉），T=-0.00684H+33.06（敦煌）。由图6（a）可知，1971—2010 年三个探空站夏季0 ℃层高度均呈显著升高趋势，变化速率分别为77.3 m·（10a）-1、52.7 m·（10a）-1和48.6 m·（10a）-1，这与该流域冰川中值面积海拔呈上升趋势一致［16］。夏季0 ℃层高度的上升意味着暴露在0 ℃层之下的冰川消融区范围增大，冰川物质亏损使得更多的冰川冰转换为液态融水。由昌马堡水文站夏季径流与0 ℃层高度散点拟合可知，疏勒河上游夏季径流量与三个探空站的0 ℃层高度均存在显著的正相关关系，相关系数在0.7 以上，且均通过0.01 显著性检验［图7（b）］，这与陈忠升等［27］对西北干旱区夏季径流量对大气0 ℃层高度的响应研究结果相吻合，表明西北干旱区夏季0 ℃层高度升降已成为影响冰川融水补给型河流径流量变化的关键指示因子。

图6 疏勒河流域夏季0°C层高度变化及与径流量变化关系Fig.6 The change of 0 ℃layer height and relation between runoff in summer and 0 ℃layer height in upper Shule River basin

4.3 冰川变化对疏勒河上游径流影响

疏勒河上游流域冰川融水主要源自祁连山脉的疏勒南山，1956—2010 年祁连山冰川自西向东呈加速退缩趋势［16］。由昌马堡水文站以上疏勒河流域两次冰川编目数据（表3）可知，1966—2006 年期间流域冰川整体呈退缩态势，冰川面积减少79.17 km2（-19.06%）。随着冰川萎缩，间接导致较大规模冰川转为小规模冰川，如面积＜0.1 km2的冰川数量增加63 条，面积＜1.0 km2的冰川数量占比由78.02%上升为81.44%，使得该流域冰川数量以小规模冰川为主这一特征更加凸显，而冰川面积减少则以≥1.0 km2的冰川为主，其面积减少（-65.12 km2）远大于＜1.0 km2的冰川（-14.05 km2）。相比河西内流水系其他流域，疏勒河流域冰川平均规模最大（0.78 km2）［16］，其中面积介于5～10 km2、10～20 km2的冰川分别有13 条和4 条，这也是该流域冰川融水对径流贡献相较其他流域大的原因，已有研究得出疏勒河流域冰川融水约占出山径流比重的23.6%～30.5%［12-14］，表明冰川变化对于疏勒河上游径流变化具有举足轻重的作用。

表3 昌马堡水文站以上疏勒河流域不同面积等级冰川变化Table 3 Changes in glaciers of different area sizes in the Shule River basin above Changmabao hydrological station

5 结论

本文基于昌马堡水文站1954—2016 年径流观测资料及疏勒河流域周边气象观测资料，结合两次冰川编目数据，较系统分析了疏勒河上游长序列径流变化事实及其影响因素，得出以下主要结论：

（1）1954—2016 年疏勒河上游年径流量总体呈显著上升趋势，变化速率为1.00×108m3·（10a）-1，其中1970s—1990s 初径流偏枯（平均年径流量为8.86×108m3），2000 年之后径流处于丰水期（平均年径流量为13.19×108m3），径流突变发生在1999 年；汛期和非汛期径流亦均呈显著增加趋势，汛期径流突变时间与年径流一致，而非汛期径流不存在明显的突变年份。疏勒河上游径流周期振荡特征显著，年径流、汛期和非汛期径流大体都存在15a 和7a 左右的显著主周期。

（2）正积温和降水是决定疏勒河上游径流变化的主导气候因子，二者可解释80%以上的径流变化，且正积温对径流的影响强于降水，这与疏勒河上游流域冰川、积雪和冻土普遍发育有关。疏勒河上游夏季径流与0 °C 层高度存在显著正相关关系（r≥0.7），表明夏季0 ℃层高度升降是影响疏勒河上游径流量变化的关键指示因子，可间接表征疏勒河上游汛期径流量变化。1966—2006 年疏勒河上游流域冰川冰储量减少约5.77 km3，冰川变化对于疏勒河上游径流变化具有举足轻重的作用。