CMADS在玉龙喀什河径流模拟中的适用性研究
2022-09-24骆成彦陈伏龙何朝飞龙爱华乔长录
骆成彦, 陈伏龙, 何朝飞, 龙爱华,2, 乔长录
(1.石河子大学水利建筑工程学院,新疆 石河子 832000;2.中国水利水电科学研究院流域水循环模拟与调控国家重点实验室,北京 100038)
玉龙喀什河位于和田流域,属于西北内陆干旱区,水资源严重制约着该区域的经济发展,在近50 a来和田流域的气候特征有了明显的变化,洪水和干旱等极端气候频繁发生,严重影响了人类活动[1],亟需精确的水文模型来模拟预测该流域的地表径流,为水资源和洪水风险管理提供科学指导[2-3]。通常,水文模型采用地面气象站收集的气象数据来进行水文过程和径流预测的研究[4]。但由于环境特殊、技术条件落后等因素,在西北内陆干旱区,尤其是高寒山区地面气象测站稀疏,甚至无地面气象站,即使有测站的地方也存在不同程度的缺测记录[5]。因此,地面气象数据的不完备严重制约着水文模型的模拟精度,迫切需要更高分辨率的气象数据驱动水文模型。
目前,有2种高分辨率气象数据:一种是全球气候模式或区域气候模式的降尺度数据,另一种是大气再分析数据集[6]。相比于前一种,大气再分析数据由于其覆盖范围广和更高分辨率的优点成为驱动水文模型的首选,被国内外学者广泛使用于各类水文模型模拟研究[7]。大气再分析数据是由地面实测数据和气象卫星数据再结合的气候模式,在真实的环流场强迫下,通过数据同化等技术制作得到历史气象数据集[8]。国际上使用较为广泛的主要有美国NCEP 气候再分析数据[Climate Forecast System Reanalysis(CFSR)、National Center for Atmospheric Research(NCAR)]、美国国家航空航天局再分析产品(Modern-Era Retrospective Analysis for Research and Applications,MERRA)、欧洲中期天气预报中心(European Centre for Medium-Range Weather Forecast,ECMWF)的EAR-Inter-Im、EAR-15、EAR-40和中国大气同化数据集CMADS 等[9-11]。许多学者对CMADS在中国不同地区的适用性进行了验证,在中国西部地区,孟现勇等[12]以黑河流域为例,使用CMADS、CFSR 和TWS 再分析数据驱动SWAT 水文模型,系统地分析了3 类再分析数据对黑河流域水文过程的影响,结果表明使用CMADS 再分析数据集的径流模拟效果均优于其他2 个;张利敏等[13]使用CMADS 驱动SWAT 模型研究了其在辽宁浑河流域的适用性,认为CMADS 能够很好地表现下垫面地 形 特 征;Jun 等[6]和 刘 兆 晨 等[14]对 比 分 析 了CMADS 和CFSR 的气候要素与气象站的观测值,并结合SWAT 模型模拟水文过程,研究表明CMADS+SWAT 模式的模拟效果最优。通过以上多地区、多方面的对比应用研究,基于卫星和再分析降水数据的广泛应用,极大地促进了地面实测气象站点稀少地区的水文气象研究,且CMADS 再分析数据集比其他再分析数据集表现出更高模拟潜力。
玉龙喀什河流域地形复杂,属于典型的高寒山区,流域内无气象站点分布,不能满足分布式水文模型模拟的气象驱动要求,必须利用气候模式或再分析资料作为流域内气象替代数据。CMADS 作为适应于分布式水文模型的气象驱动数据,还尚未在中国西北寒旱区玉龙喀什河流域得到充分的研究,需要评估其在研究区模拟气温和降水的能力,为玉龙喀什河流域水文模型的构建提供备选数据集,对研究该流域水资源的形成和分布规律,理解高寒山区水文循环过程具有重要科学意义,对促进玉龙喀什河流域水资源合理利用和防范夏洪春旱等自然灾害提供一定科学指导。
1 数据与方法
1.1 研究区概况
玉龙喀什河(以下简称玉河)为和田河源流(东支)之一,位于塔里木盆地南缘,发源于昆仑山北坡,自南向北流经和田绿洲,经出山口同古孜洛克水文站(79.91609°E,36.81055°N)与喀拉喀什河(西支)并行,在下游于阔什拉什汇流入和田河,最终注入塔里木河,全长504 km[15-16]。研究区位于玉龙喀什河上游河源产流区(77.25°~81.75°E,34.75°~36.25°N,图1),地势南高北低,落差极大,集水面积为14887 km2;研究区属于内陆干旱区,气候干燥、降水稀少,蒸发强烈[17],玉河出山口同古孜洛克水文站日均气温12.7 ℃,年均降水94.2 mm[18]。玉河地表径流补给主要依靠冰川积雪融水及部分高山降水,平均年径流量为22.6×108m3,年际、年内变化大,6—9 月为玉河洪水期,其地表径流量占到全年的80%,最大年径流超出平均径流64.2%,最小年径流低于平均径流34.9%。
图1 研究区地形及站点分布Fig.1 Topography and station distribution of the study area
1.2 数据来源与处理
1.2.1 气象数据 气象数据作为分布式水文模型的主要驱动数据,其精确性与完整性至关重要[19]。由于研究区位于高寒山区,地形复杂多变,研究区内没有相应的气象观测站,只有同古孜洛克水文站,但该站位于流域边缘平原处,无法替代整个流域的气温降水情况,故本文采用中国大气同化驱动数据集CMADS V1.1 版本作为SWAT 模型的气象驱动数据[20]。CMADS 数据集数据源为国家级自动气象站观测数据,包括测站的经纬度和海拔高度以及降水、气温、气压、比湿、辐射和风速等气象要素数据,
融合卫星数据产品,利用多重网格三维变分、循环嵌套、重采样等方法,并以NCEP/GFS再分析数据为背景场做地形调整并插值到0.25°×0.25°格点上[20]。CMADS 数据集已按照SWAT 模型气象数据输入格式进行了调整,不需要建立每个站点的天气发生器数据库,只需建立研究区的气象站索引表。研究区内CMADS站点分布见图1,数据来源见表1。
1.2.2 下垫面数据 SWAT 模型所需下垫面数据包括DEM 数字高程、土地利用数据、土壤数据(表1)。DEM数字高程是获取研究区地形地貌、坡度的重要依据,用以生成河网水系、划分子流域及确定流域边界[21],经过SWAT 模型处理生成的玉河水系以及流域边界(图1);土地利用数据是反映研究区某个时期的土地表面要素的空间分布状态、地表特征以及人类对土地的开发利用的数据资料[22],对于土地利用数据,需进行土地利用重分类预处理,之后与SWAT模型内部数据库的索引表进行链接。根据中
表1 数据类型及来源Tab.1 Data types and sources
国科学院资源环境科学数据中心分类指标和研究区实际情况,将玉河土地利用类型重分类为6类(图2),分别是耕地、林地、草地、水域、永久性冰川雪地、裸地,其中裸地的占比最大,为49.23%;土壤数据用以反映研究区不同土壤的空间分布以及各种土壤物理属性;通过对世界土壤数据库的研究区范围裁剪以及重分类后,得到研究区的土壤类型空间分布(图3),其中永冻薄层土占比最大,为36.04%;对于土壤物理属性,部分土壤基本数据可以从HWSD土壤数据库获得,如土壤各粒径成分的含量、碳酸盐含量、有机碳含量等,大部分数据无法直接测得,需通过SPAW 软件的Soil Water Characteristic模块计算参数,如土壤湿密度、土壤有效含水量、饱和水力传导系数等。
图2 土地利用类型Fig.2 Land use types
图3 土壤类型Fig.3 Soil types
1.3 研究方法
1.3.1 SWAT 模型设置 SWAT 模型是由美国农业部农业研究中心开发的基于流域尺度的一个长时段分布式流域水文模型,可以直接模拟水流、泥沙等物理过程[23],可以定量化的输入参数(如气候、植被等变化)以研究对径流等的影响[24-25]。根据研究区实际情况,设置流域面积划分阀值为25000 hm2,将流域划分为53 个子流域,子流域平均面积为280.88 km2;将流域划分为1863 个水文响应单元HRU(Hydrologic Response Unit)。在加载2008—2016 年CMADS 气象驱动数据之后,建立SWAT 月径流模型,设置2008—2009 年为预热期,2010—2012年为模型参数率定期,2013—2016年为模型验证期。
模型参数的率定采用SWAT-CUP(SWAT Calibration and Uncertainty Program)软件,此软件专门用于SWAT 模型的参数敏感性分析、模型的率定以及不确定性分析[26-27]。通过比较SWAT-CUP所提供的5种算法,SUFI-2算法结合了不确定性分析与优化,并使用拉丁超立方抽样方法处理参数,具有运行简单、计算效率高以及可灵活选择目标函数等优点[28-30]。故本文选用SWAT-CUP提供的SUFI-2算法,并选择纳什效率系数(Nash-Sutcliffe Efficiency Coefficient,NSE)为目标函数来率定模型参数,利用t-Stat 值和P-Value值进行模型参数敏感性检验。
1.3.2 CMADS 数据集评价方法 为评价CMADS 数据集在研究区的模拟精度,选取距离同古孜洛克水文站最近的CMADS 的148-81(80.03125°E,36.78125°N)站点进行分析。因同古孜洛克水文站仅有逐日降水以及日平均气温,故只评价CMADS数据集与观测站点的降水与平均气温。对于平均气温主要采用相关系数、均方根误差和相对误差3项指标;而对于降水为了全面表征CMADS 数据集在该地区模拟的效果,增加了探测率、误报率和临界成功指数3项指标,具体公式如表2。
1.3.3 径流模拟效果评价方法 对SWAT模型的径流模拟精度采用纳什效率系数、模型决定系数和相对误差3类指标进行评价,具体公式如表2。
表2 精度评价公式Tab.2 Accuracy evaluation formula
2 结果与分析
2.1 CMADS精度及时空分布特征评价
2.1.1 CMADS 精度及评价 降水数据的精度很大程度上决定了模型的径流模拟效果,为直观分析CMADS 降水数据与研究区同古孜洛克水文站实测降水的拟合程度,作日尺度降水过程线(图4)。由图4可以看出,CMADS降水所表现出的峰值以及峰值出现的时间与观测降水在大部分年份都有较好的对应,但在个别年份存在峰值模拟效果较差,如在2010—2011 年出现明显的错峰现象,在2013—2014 年CMADS 的降水峰值明显低于实测降水,说明CMADS 对连续的极端降水和变化剧烈的降水事件的模拟性能较差。
图4 同古孜洛克水文站与CMADS 148-81格点日降水过程对比Fig.4 Comparison of precipitation process of Tongguziluoke hydrological station between and grid 148-81 of CMADS
为定量评估CMADS 降水数据的拟合效果,作散点图进行线性拟合(图5)。CMADS 降水与实测降水的相关系数CC为0.65,两者之间具有较强的线性相关性;相对误差RE为24.67%大于0,CMADS降水与实测降水数据年均降水量分别为97.77 mm 与91.03 m,CMADS 对实测降水有所高估,但相差仅6.74 mm;均方根误差RMSE 为1.14,表明两者之间误差较小;探测率POD、误报率FAR 和临界成功指数CSI分别为0.87、0.25和0.68,说明CMADS对于实际降水事件发生的概率拟合准确性较高。总体而言,CMADS降水数据不论在降水量还是降水事件再现的模拟上,都有着较好的表现,CMADS 降水数据在研究区有较高的精度和可靠性,加上其覆盖全面,分辨率高的特点,可在研究区代替实测降水数据驱动SWAT水文模型进行径流的模拟。
图5 同古孜洛克水文站与CMADS的148-81格点降水散点图Fig.5 Precipitation scatter diagram of Tongguziluoke hydrological station between and grid 148-81 of CMADS
气温作为SWAT模型必不可少的气象数据输入项,决定着SWAT 模型降水形态的转变以及融雪速度和融雪水量。研究区属于内陆高寒山区,融雪是径流的重要组成部分,故对融雪模拟的效果也将对最终径流的模拟效果产生很大的影响。相对于降水来说,CMADS数据集对气温的拟合效果要优于降水。从气温拟合散点图(图6和图7)可知,对于最高气温和最低气温的拟合,CMADS气温与实测最高气温与最低气温之间的相关系数CC分别达到了0.998与0.995,CMADS 气温与实测气温之间具有相当高的线性相关性;相对误差RE 分别为11.22%与16.52%均大于零,CMADS 对实际气温是低估的;均方根误差RMSE 分别为1.787 和1.75,CMADS 对实际气温的模拟误差较小。通过以上指标分析,CMADS 对于气温的模拟具有相当高的精度和可靠性,可用于研究区作为SWAT 模型的气象驱动数据来模拟径流。
图6 同古孜洛克水文站与CMADS的148-81格点最高气温散点图Fig.6 Maximum temperature scatter diagram of Tongguziluoke hydrological station between and grid 148-81 of CMADS
图7 同古孜洛克水文站与CMADS的148-81格点最低气温散点图Fig.7 Minimum temperature scatter diagram of Tongguziluoke hydrological station between and grid 148-81 of CMADS
2.1.2 CMADS降水气温时空分布特征评价 SWAT模型作为分布式水文模型,不仅对输入数据精度有较高的要求,而且对气温与降水的空间分布特征也有较高的要求。从地形图(图1)可以看出,研究区地形复杂,海拔变化大,36°15′N为流域的高程划分线,其下至流域出水口,是山前冲积扇平原区,地势平坦;其上至流域边界,海拔均在3000 m以上,是冰川积雪的主要覆盖地,沟壑交错。为充分考虑地形对降水和气温的影响,采用ArcGIS中地统计处理模块的协同克里金法对CMADS 的气温与降水进行空间插值,协同克里金插值是基于高程利用区域化变量的原始数据和变异函数的结构特点,对未采样点的区域化变量的取值进行线性无偏、最优估计[31]。
对于降水,从时间分布上来看(图8),年内降水主要集中于夏季,夏季流域降水量占全年总降水量的38.03%,春季降水最少,占全年总降水量的18.52%;从空间分布来看,降水并未表现出明显的随海拔升高而降水量增加的趋势,四季的结果均表现出强烈的地形特征的降水集中区,图8 中所表现的降水集中区刚好位于流域内一处开口向东的马蹄形盆地,气流从东面开口进入,使气流辐合上升,形成地形雨。此外,四季均存在一条20~30 mm 的降水分割带,该降水带主要沿昆仑山脉走势分布。以北由于昆仑山阻挡了南来的暖空气,造成高山区降水较少;同时,该降水带的移动也表征了四季的变换。总体而言,CMADS降水数据的空间分布能够较好地反映研究区的地形特征。
对于气温,从时间分布来看(图8),日最高气温均值出现在秋季,达到28 ℃,日最低气温均值出现在春季和冬季,达到-28 ℃;以日最高气温均值为例,自春季开始,4~8 ℃气温带逐渐向南转移,秋季到达流域中部,之后再次转移到流域最北端,此为年内的气温循环。从空间分布来看,整体上呈现出随海拔升高气温降低的趋势,气温直减率约为-1.56 ℃·km-1,其符合高寒山区的气温直减率[32]。另外,四季气温空间分布均存在一个气温低值区,位置恰好位于降水集中区,并在此区域形成环状气温带,易形成气流辐合,从另一方面说明了该区域为降水集中区的原因。通过以上对降水和气温的时空分布特征分析,说明CMADS 很好地捕捉到了地形的变化。
图8 CMADS降水气温时空分布Fig.8 Temporal and spatial distribution of CMADS precipitation and temperature
2.2 模型参数率定结果
根据研究区实际情况,拟定17 个参数进行率定,因研究区位于高寒山区,径流补给主要为高山降水、地下水和冰川融雪,故参数主要根据径流补给来源选择(表3),t-Stat值表示敏感性的程度,绝对值越大越敏感;P-Value 值决定敏感性的显著性,值越接近0,越显著。SUFI-2 为SWAT 模型最常用的参数识别算法,该算法识别得到的参数最优值可能落在参数先验分布范围之外,当模型率定出现不合理参数取值时,应对这些参数取值进行分析和调整[33]。
表3 参数率定结果Tab.3 Parameter calibration results
在敏感性排序前5 的参数中,与融雪模块和地下水模块相关的参数各占2 个,其中最大融雪因子(SMFMX,北半球为6 月21 号的融雪因子)最为敏感,最优值为5.1345 mm·℃-1·d-1;融雪基温SMTMP为0.5 ℃,当子流域内的积雪温度超过此阀值将发生融雪,以上参数值与同属于寒旱区的开都河月尺度径流模拟参数取值相近[34]。对于地下水模块,基流α因子ALPHA_BF 表现最为敏感,它是地下水径流对补给量变化响应的直接指标,率定最优值为0.1100;其次为地下水延迟时间GW_DELAY,它表征了地下水流出土壤剖面补给地表水的时间延迟,根据尤扬等[35]在不同气候情景下和田河上游径流变化中,建立的SWAT月尺度模型确定参数,最终率定为329 d;TLAPS 最终率定结果为1.532 ℃·km-1,与气温的垂直分布特征所表现的结果一致。通过以上分析,SWAT-CUP 率定的各项参数和敏感性排序,与刘晓笛[15]基于SWAT 模型的和田河上游气候和土地利用变化的水文效应模拟研究中,所构建的月尺度模型的参数结果较为相近,均处于合理的物理范围内。
2.3 径流模拟结果及适用性评价
由表4 可知,模型在率定期和验证期的纳什系数分别达到0.845 和0.836,模型决定系数分别达到0.8562 和0.8153,验证期此2 项指标均有所减小。率定期和验证期实测月径流量均值较模拟均值分别高出7.464 m3·s-1和10.975 m3·s-1,相对误差较率定期增大为-12.47%,模拟的均值偏低和相对误差为负,均说明径流的模拟值低估了实测值。但总体上来看,率定期和验证期的各项指标较为可观,表明SWAT模型在该流域的径流模拟中具有较好的适应性。
表4 径流模拟评价结果Tab.4 Runoff simulation evaluation results
模拟期内模拟径流与实测径流曲线趋势基本一致(图9),在每年5 月径流逐渐增大,于8 月达到峰值,在9月骤降,从10月开始到次年4月期间均为基流,基流天数约210 d。局部来看,存在部分误差较大的时间点,其中最大负误差出现在2012 年8月,模拟值低于实测值246.8 m3·s-1,最大正误差出现在2013 年9 月,模拟值高于实测值244.3 m3·s-1,最小误差为2010年9月,误差仅为0.3 m3·s-1。大部分模拟值处于实测值的下方(图9),与评价指标中相对误差所得到的结果一致,且低估现象主要出现在基流期和洪水期。基流期一般从冬季到次年春季末,雨季已过,且气温已经降至0 ℃以下,冰川和积雪融水骤减,径流主要由地下水补给。但研究区处于高寒山区,对于地下水的监测和研究较少,缺少资料,无法对参数进行准确的率定,导致模型中地下水模块的模拟结果不够理想,造成模拟期基流偏低。
图9 月径流模拟结果Fig.9 Runoff simulation results
模拟结果表明,校准期和验证期内,玉河径流量模拟值和实测值之间的各项评价指标的结果与刘晓笛[15]、罗开盛等[36]和宋玉鑫等[37]使用SWAT 模型在西北寒旱区建立的月尺度径流模拟结果相近,均在允许范围之内。因此,SWAT 模型总体上适用于该区域高寒山区的水文过程模拟。
3 讨论
本文的研究结果表明:气象数据作为水文模型的输入数据,其不确定性对径流模拟结果有着重要的影响[25]。由于研究区的资料限制,只能使用仅有的出山口控制站同古孜洛克水文站来评价CMADS在此位置的精度,具有较大的不确定性。对于降水,其不确定性主要在于单站点的降水控制面积有限,无法代表整个流域的实际降水情况,包括降水强度变化、降水位置以及降水事件的发生[23],降水站点密度对水文建模具有一定影响[38];对于气温,研究区海拔变化幅度大,且存在冰川积雪,实际气温的变化无法得到测量,所以无法准确评估CMADS降水和气温时空变异性。在评价CMADS 的时空分布特征时,使用协同克里金插值对未采样点的区域进行取值,得到的结果虽然是一个最优估计值,但并不能反映未采样点的实际值,其空间分布特征也存在一定误差,可进一步引入海拔、坡向、坡度等因子,提高插值精度[39]。整体而言,CMADS 在精度以及时空特征分布方面均有着良好的表现,用以驱动SWAT 模型模拟研究区的径流,其结果有着较高的可靠性。
本文使用2008—2012 的实测径流率定SWAT模型,并以率定好的模型验证2013—2017的实测径流,结果表明验证期的指标较率定期的差,模拟效果较差不仅是因为产流模式与模型预设不匹配、验证期水文环境变化等,还与率定径流序列过短导致率定期模型水文过程表征不准确有关[33]。SWAT水文模型作为分布式水文模型,尽管在结构上考虑了融雪和冻土对水文循环的影响[38],拥有较多的物理参数,有着相对完善的水文循环过程,但是由于缺乏冰川模块,在夏季汛期缺少冰川融水资料,故造成模拟过程中模拟极值达不到实测极值。孙铭悦等[40]使用格点数据驱动HBV 水文模型在西北寒旱区呼图壁河流域的日尺度径流模拟中,同样存在模拟极值低于实测径流极值的现象。祖拜代·木依布拉等[41]使用单站点构建了乌鲁木齐河上游的月尺度SWAT模型,验证期决定系数仅为0.75,且相对误差高达22.20%,相比之下,本文以CMADS格点数据驱动SWAT模型在寒旱区的径流模拟效果要更好,故CMADS+SWAT 模型在寒旱区玉河流域有着较好的适用性。
4 结论
本文基于研究区DEM、2010年土地利用和土壤数据,构建了玉龙喀什河产流区的SWAT模型,验证了CMADS 再分析数据集的精度和空间分布特征,模拟了研究区的月尺度径流,得出如下主要结论:
(1)总体上CMADS 再分析数据的精度和可靠性较好,可替代无站点流域作为气象数据驱动水文模型。气温的模拟效果优于降水,降水大部分模拟值与实测值较为接近,CMADS再分析数据的时空分布特征基本符合实际,能够准确地捕捉到下垫面的地形特征变化,但CMADS 对于极端的气候变化的模拟能力有待提高。
(2)以CMADS+SWAT 模式模拟玉河流域月尺度径流,整体上模拟径流序列与实测径流变化趋势一致,在研究区具有较好的适用性,但在峰值部分出现明显的低估现象,主要原因是缺少夏汛期的冰川融水资料;另外,CMADS 数据集本身的时间跨度不够长,导致验证期的模拟效果稍差,模型的参数率定存在一定的偏差。
总体而言,CMADS再分析数据集通过研究区的验证,对于无站点流域的水文模拟效果较好,可为无降水、少降水站点或缺失降水资料的区域提供较为可靠的数据来源,为扩大水文模拟的时间和空间尺度提供了可能。但缺乏与其他再分析数据的比较,在后续的研究中将耦合其他再分析数据进行比较分析,对再分析数据进行集合平均,结合各类再分析数据的优势,减少误差;并收集研究区冰川资料,加入冰川模块,完善水文模型结构,使模拟结果精度更高。