1. 引言
水文模型通过确定模型参数来模拟降雨径流关系,通常将实测资料年份分为率定期和检验期,通过参数优化方法优化出一组最优参数用于径流模拟。如王中根[1]等利用SWAT模型,模拟了莺落峡流域的日、月径流,模型的效率系数达到了0.83。江燕[2]等采用新安江模型,通过对粒子群优化算法进行改进,对天福庙水库径流进行模拟,率定期和检验期效率系数均超过0.82。赵彦增[3]等利用HBV模型,对淮河官寨流域径流进行了模拟。由于水文机理的复杂性、流域的多样性和模型的局限性等原因,降雨径流关系并非唯一不变,不同时期的降雨径流关系会随着下垫面等情况的改变而发生变化。本研究通过构建汉江上游三水源新安江模型,采用实测资料和遗传算法及SCE-UA算法分时期率定模型,研究水文模型参数在不同年代间的变化规律;建立基于贝叶斯理论的新安江参数不确定性分析模型,开展水文模型的后验分布和不确定性分析;通过水文模型参数的后验分布,分析水文模型参数在不同年代的时变规律,为变化情境下径流模拟提供参考。
2. 研究概况和研究方法
2.1. 研究概况
汉江发源于秦岭南麓,襄阳以上河流总体流向东,襄阳以下转向东南,于武汉市注入长江,干流全长1577 km,流域面积约15.9万km2。汉江干流丹江口以上为上游,河长925 km,占干流总长的59%,控制流域面积9.52万km2。丹江口至钟祥为中游,河长270 km,占干流总河长的17%,控制流域面积4.68万km2。钟祥以下为下游,长382 km,占干流总长的24%,控制流域面积1.7万km2。汉江流域水系和站网示意图见图1。
流域多年平均气温14~16℃,极端最高气温42℃,极端最低气温−13℃,多年平均相对湿度74%,最大风速21 m/s,多年平均蒸发量为848 mm,最大月蒸发量出现在6月或7月,最小出现在12月或1月。多年平均无霜期212~254 d。
Figure 1. Diagram of the river network and hydrological stations in Hanjiang River Basin
图1. 汉江流域水系站网图
2.2. 研究模型
本次研究采用国内外广泛使用的三水源新安江模型[4]。最初的模型为两水源——地表径流和地下径流,20世纪80年代,相继提出了三水源和四水源的新安江模型。新安江模型在国内外湿润半湿润地区得到了广泛的应用。
Figure 2. Structure diagram of Xin’anjiang Model
图2. 新安江模型结构图
新安江模型的主要核心是通过流域蓄水容量曲线来反映流域下垫面土壤的不均匀性对产流变化的影响,模型计算主要包含蒸散发计算、产流计算、分水源计算以及汇流计算四个部分,模型结构如图2所示。
2.3. 研究方法
研究分别采用遗传算法[5]和SCE-UA [6]算法对水文模型参数进行优选。遗传算法具有多点寻优、并行处理等一系列特点。遗传算法的搜索过程是从初始解群开始的,以模型对应的适应函数作为寻优判据,适者生存、劣者淘汰,从而可以直接对解群进行操作,而与模型的具体表达方式无关。SCE-UA算法在综合了随机搜索和确定性搜索,聚类分析法以及生物群体竞争演化等多种较优的全局优化方法的基础上进行改进,可以快速高效地搜索到全局最优解,广泛应用于国内外水文模型参数优选。
本次参数后验分布研究采用DREAM算法[7]。DREAM算法融合了自适应Metropolis算法和传统差分进化算法的优势,与其它马尔科夫链蒙特卡罗方法不同,DREAM算法从多个不同的搜索起点同时产生多条平行的马尔科夫链,使得马尔科夫链充分遍历参数空间,从而搜索到全局最优解。同时DREAM算法在搜索的过程中能够根据计算结果自动调整寻优的方向和步长,提高了算法的运算效率。凭借算法超强的全局收敛能力和鲁棒性,DREAM算法被广泛应用于复杂系统的估计以及后验分布的推求等各个领域。
2.4. 评价指标
选用Nash-Sutcliffe效率系数(NSE)和水量相对误差(RE)来评估模型模拟效果,其中Nash-Sutcliffe效率系数越接近1,同时水量相对误差越接近0,说明模拟效果越好。效率系数和水量误差计算公式如下所示。
(1) 水量相对误差(RE):
(1)
(2) Nash-Sutcliffe效率系数(NSE):
(2)
式中:
为第i天的实测值;
为实测系列的均值;
为第i天的模拟值;
为系列的总长度。
3. 水文模型模拟结果
3.1. 模型构建
以汉江上游1961~1990年降雨数据、流量数据和蒸发数据作为输入,建立汉江上游三水源新安江模型。其中,将汉江上游多个雨量站点的降雨数据由泰森多边形加权平均得到面均降雨数据;蒸发数据采用彭曼公式根据流域临近站点的气象数据计算得到;输入流量为白河站径流过程,模型采用日时段。
3.2. 模拟结果分析
分别采用SCE-UA算法和遗传算法来优选出汉江上游新安江模型的最优参数,表1列出了汉江上游新安江模型的模型效率系数(NSE)与水量相对误差(RE)。
Table 1. Comparison of simulation results for Xin’anjiang Model in the upper reaches of Hanjiang River
表1. 汉江上游新安江模型模拟效果对比表
时段 |
NSE |
RE |
遗传算法 |
SCE-UA |
遗传算法 |
SCE-UA |
1961~1970 |
0.817 |
0.841 |
−4.26% |
−1.17% |
1971~1980 |
0.821 |
0.853 |
−4.16% |
−1.12% |
1981~1990 |
0.812 |
0.833 |
−4.94% |
−1.26% |
模型效率系数表示模拟径流与实测径流吻合程度。模拟和实测径流拟合程度高,模型效果较好。由表1中的统计数据看出,运用新安江模型在汉江上游模拟日径流,两种参数优化方法得到的不同时段模型确定性系数均超过0.80,水量相对误差最大值为4.94%,新安江模型在汉江上游的模拟效果较好。
两种参数优化方法对比发现,SCE-UA的效率系数要高于遗传算法,同时水量相对误差小于遗传算法,SCE-UA优化方法结果优于遗传算法。1981~1990实测和模拟径流过程见图3。因此,选择SCE-UA算法优化汉江上游新安江模型,进而开展参数的不确定性分析。
Figure 3. Comparison of simulated and observed flow using Xin’anjiang Model and SCE-UA in the Upper Reaches of Hanjiang River (1981~1990)
图3. 汉江上游新安江模型SCE-UA优化方法实测径流与模拟径流过程(1981~1990)
4. 参数年际变化规律研究
4.1. 参数设置
本次采用SCE-UA全局优化算法在参数可行域范围内对不同时期内的新安江模型进行参数优化率定,为了比较不同时期的优化结果,本文采用相同的参数可行域范围对不同时期的水文模型参数进行优选,各参数的范围见表2。
Table 2. Upper and lower limits of parameters in Xin’anjiang Model
表2. 新安江模型参数取值上下限表
参数 |
WM |
X |
Y |
KE |
B |
SM |
EX |
KI |
KG |
IMP |
C |
CI |
CG |
N |
NK |
上限 |
300 |
0.4 |
0.6 |
1.3 |
0.8 |
80 |
1.5 |
0.45 |
0.45 |
0.1 |
0.2 |
1 |
1 |
12 |
25 |
下限 |
100 |
0.01 |
0.2 |
0.6 |
0.1 |
10 |
1 |
0.01 |
0.01 |
0.01 |
0.15 |
0.7 |
0.98 |
0.5 |
0.8 |
4.2. 年际变化规律
将不同时期下的降雨、蒸发、径流系列作为新安江模型的输入,根据贝叶斯理论,采用DREAM算法推求参数的后验分布。考虑到模型参数的相关性和不确定性以及贝叶斯方法的计算效率,本次仅仅对模型敏感参数(KE、SM、KI、KG、CI、CG、N和NK)进行后验分布的推求,各敏感参数的可行域范围见表2。对剩下的不敏感参数(WM、X、Y、B、EX、IMP、C)则采用各时期下SCE-UA算法优化得到的模型参数值的平均值进行固定。待8条马尔科夫链充分混合并收敛后,取每条链收敛后的1000组参数,共计8000组参数统计模型各参数的分布,得到不同时期下各敏感参数的后验分布如图4所示。
(a) KE (b) SM
(c) KI (d) KG
(e) CI (f) CG
(g) N (h) NK
Figure 4. Posteriori distribution boxplots of parameters in different periods of the Xin’anjiang Model in the upper reaches of the Hanjiang River
图4. 汉江上游新安江模型不同时期各参数的后验分布箱线图
从模型参数的中值来看,KI、KG、CI、CG的参数值随不同时期逐渐增加,而SM、N的参数值随不同时期逐渐减小,KE和NK没有出现趋势性变化。分析参数的物理意义可以发现,KE影响着水量平衡,SM和KI、KG作为三个分水源参数,共同决定了地面径流、壤中流和地下径流的比例,其相关性较强,对于汇流参数CI、CG、N和NK,影响着各水源的汇流过程,最终决定了水文过程线的形状。
从参数的分布来看,对于参数KI和KG,随时段的改变,其95%置信区间的宽度逐渐增大,对于参数SM、CI、CG和N,随时段的改变,其95%置信区间的宽度逐渐减小,而参数KE和NK随着时段的改变,参数未发生趋势性变化。但参数NK具有非常强的敏感性,参数值微小的变化将引起模型模拟精度的强烈变化。
将8000组收敛于后验分布的水文模型参数放入新安江模型中计算,得到8000组模拟径流序列,并根据模拟径流序列与实测径流序列计算其效率系数NSE和水量相对误差RE,得到后验分布参数的模拟结果如图5所示。
(a) NSE
(b) RE
Figure 5. Simulation results of posterior distribution parameters at different periods
图5. 不同时间尺度下后验分布参数模拟结果图
从后验分布参数的模拟结果可以看出,虽然后验分布参数模拟径流序列的精度要比与SCE-UA算法优化得到的径流序列的精度略低,但是后验分布参数模拟的径流的效率系数(NSE)值均在0.74~0.84之间,水量相对误差(RE)在−4%~−1%之间,说明后验分布参数能够较好地模拟实测径流。从两个评价指标的分布来看,1981~1990时段模型模拟的不确定性区间要大于1961~1970和1971~1980时段模拟的不确定性区间,主要原因是随着人类活动作用的增强,降雨径流关系发生变化,虽然考虑了径流的还原,但难以准确还原到天然过程。
计算上述8000组模拟径流序列对应的每一时刻的95%置信区间,得到模型参数不确定性导致的水文模拟径流的不确定性区间,然后在该8000组模拟径流序列的基础上,计算出每一时刻对应的95%置信区间,结果如图6所示。
从图中可以看出,大部分实测径流都包含在模拟径流的95%置信区间内,后验分布参数能够较好地模拟实测径流,通过开展参数的不确定性分析,可以得到不同置信区间的范围值,进而可以更好地为决策提供依据。
Figure 6. 95% confidence interval obtained from simulated flow by posterior distribution of parameters (1981~1990)
图6. 采用参数后验分布中抽取参数模拟径流过程得到的95%置信区间(1981~1990)
5. 结论
本次研究选取汉江上游白河站作为控制点,建立了汉江上游三水源新安江模型,分别采用遗传算法和SCE-UA算法对新安江模型的参数进行优选,用1961~1990年逐日实测径流资料分不同时期对三水源新安江模型进行率定,应用结果表明:两种参数优化方法得到的不同时段模型确定性系数均超过0.80,水量相对误差最大值为4.94%,新安江模型在汉江上游的模拟效果较好;SCE-UA的效率系数要高于遗传算法,同时水量相对误差小于遗传算法,SCE-UA优化方法结果优于遗传算法;但两种方法均对于高水部分有些低估。
根据贝叶斯理论,采用DREAM算法抽样,通过研究不同时期水文模型参数的后验分布,进而分析参数的时变规律和不确定性。结果表明:后验分布参数模拟的径流的NSE均在0.74~0.84之间,RE在−4%~−1%之间,说明后验分布参数能够较好地模拟实测径流。1981~1990时段模型模拟的不确定性区间要大于1961~1970和1971~1980时段模拟的不确定性区间,说明随着人类活动作用的增强,水文模型参数的不确定性增加。
基金项目
长江科学院开放研究基金资助项目(编号:CKWV2021866/KY);水利部重大科技项目(SKS-2022038);国家重点研发计划项目(2021YFC3200304)。