基于BMA多模型组合的疏勒河径流预测研究

doi:10.7522/j.issn.1000-0240.2022.0141

基于BMA多模型组合的疏勒河径流预测研究

周婷^,¹^,², 温小虎^,¹, 冯起¹, 尹振良¹, 杨林山¹

1.中国科学院西北生态环境资源研究院，甘肃兰州 730000

2.中国科学院大学，北京 100049

Study on runoff prediction of Shule River based on BMA multi-model combination

ZHOU Ting^,¹^,², WEN Xiaohu^,¹, FENG Qi¹, YIN Zhenliang¹, YANG Linshan¹

1.Northwest Institute of Eco-Environment and Resources, Chinese Academy of Sciences, Lanzhou 730000, Gansu, China

2.University of Chinese Academy of Sciences, Beijing 100049, China

通讯作者: 温小虎，研究员，从事水资源与水环境研究. E-mail： xhwen@lzb.ac.cn

收稿日期: 2021-07-23 修回日期: 2021-12-23

基金资助:

国家自然科学基金项目.  42001035
中国科学院前沿科学重点研究项目.  QYZDJ-SSW-DQC031
中国科学院“西部之光”项目.  Y829821001

Received: 2021-07-23 Revised: 2021-12-23

作者简介 About authors

周婷，硕士研究生，从事水文与水资源研究.E-mail：zhouting@nieer.ac.cn , E-mail：zhouting@nieer.ac.cn

摘要

准确可靠的径流预测对于水资源的科学管理与规划具有重要意义，特别是在水资源紧缺的干旱半干旱地区，径流预测对流域内水资源高效利用与水利工程经济运行具有重要现实意义。针对径流预测通常采用单一方法进行建模与预测，难以利用各预测模型优势的问题，本文利用极限学习机（ELM）模型、支持向量机（SVM）模型、多元自适应回归样条（MARS）等机器学习方法建立了疏勒河上游未来1~7日的径流预测模型。在此基础上，运用贝叶斯模型平均（BMA）方法对ELM、SVM、MARS模型的预测结果进行组合，构建了径流组合预测模型，以获取更可靠的预测结果，并采用蒙特卡洛抽样方法获取BMA的95%置信区间，对预测结果进行了不确定性分析。结果表明：ELM、SVM、MARS模型以及BMA组合模型均适用于干旱半干旱地区的中长期日径流预测；BMA的预测精度较单一模型更高，能够提供更准确的预测值；BMA的95%置信区间对实测值覆盖率高，同时能够提供较好的确定性预测和概率预测结果。表明BMA在资料有限的条件下，表现出较单一模型更高的预测性能，可以成为干旱半干旱地区中长期日径流预测的有效方法。

关键词： 日径流量预测 ; 机器学习 ; 贝叶斯模型平均 ; 不确定性分析 ; 疏勒河

Abstract

Accurate and reliable runoff prediction is of great significance for the scientific management and planning of water resources， especially in the arid and semi-arid areas where water resources are scarce. Runoff prediction has important practical significance for the efficient utilization of water resources and the economic operation of water conservation projects.In view of the problem that it is difficult to make use of the advantages of each prediction model because a single method is usually used for modeling and prediction of runoff prediction. In this paper， the Extreme Learning Machine （ELM） model， Support Vector Machine （SVM） model and Multivariate Adaptive Regression Spline （MARS） model were used to develop the runoff prediction model in the upper reaches of Shule River in 1 to 7 day. On this basis， the Bayesian Model Average （BMA） method was also used to combine the prediction results of ELM， SVM and MARS models， and a combined runoff prediction model was constructed to obtain more reliable predictions. The 95% confidence interval of BMA was obtained by Monte Carlo sampling method， and the uncertainty of the predictions was analyzed. The results show that ELM， SVM， MARS model and BMA combination model are suitable for medium and long term daily runoff prediction in arid and semi-arid areas； BMA has higher prediction accuracy than the single models and can provide more reliable and accurate predictions； The 95% confidence interval of BMA has high coverage of measured values， and can provide better deterministic and probabilistic predictions. The results suggest that BMA has better prediction performance than the single models under the condition of limited data， and can be an effective method for medium and long-term daily runoff prediction in arid and semi-arid areas.

Keywords： daily runoff ; machine learning ; Bayesian Model Average ; uncertainty analysis ; Shule River

PDF (3631KB) 元数据多维度评价相关文章导出 EndNote| Ris| Bibtex 收藏本文

本文引用格式

周婷, 温小虎, 冯起, 尹振良, 杨林山. 基于BMA多模型组合的疏勒河径流预测研究[J]. 冰川冻土, 2022, 44(5): 1606-1619 doi:10.7522/j.issn.1000-0240.2022.0141

ZHOU Ting, WEN Xiaohu, FENG Qi, YIN Zhenliang, YANG Linshan. Study on runoff prediction of Shule River based on BMA multi-model combination[J]. Journal of Glaciology and Geocryology, 2022, 44(5): 1606-1619 doi:10.7522/j.issn.1000-0240.2022.0141

0 引言

准确可靠的日径流预测是水文研究的重要内容，也是水资源管理的重要环节^［1］，可为流域内水资源的科学管理、合理利用以及水资源评估提供重要的科学依据^［2］。其中，3日以上的中长期日径流预测对于流域内的防洪、水资源调控、生态环境保护等具有十分重要的意义^［3］。特别对于我国西北干旱半干旱地区，水资源已成为限制区域社会经济发展、影响生态安全的主要因素^［4］。因此，准确的径流预测对于区域内的水资源管理与可持续利用都具有非常重要的现实意义^［5-6］。

分布式水文模型与机器学习方法是径流模拟预测的常用工具^［7-8］。分布式水文模型可以准确地描述流域内的水文物理过程，但通常需要大量的水文、气象等数据作为输入，而西北地区由于观测站点稀疏，水文、气象资料不足，使用分布式模型会使建模过程既耗时又昂贵^［3，9］。机器学习方法仅从数学上考虑输入与输出之间的非线性关系，不需要描述水文过程的物理参数，通过挖掘数据本身的潜在规律对径流进行预测。近年来，机器学习由于算法简单、易于实现等优点，在资料有限的条件下，展现出更优的预测能力，已经广泛应用于各种条件下的径流模拟预测研究中^［10］。其中，极限学习机（Extreme Learning Machine，ELM）模型、支持向量机（Support Vector Machine，SVM）模型、多元自适应回归样条（Multivariate Adaptive Regression Spline，MARS）模型等在西北地区径流预测研究中均取得了较好的预测结果^［9-11］。虽然单一模型在特定的时空条件下都表现了较好的模拟预测性能。但是，由于径流过程的高度非线性、形成机制的复杂性以及模型结构自身的不确定性，单一模型的模拟预测结果存在较大不确定性^［12］。

针对径流预测研究通常采用单一方法进行建模与预测，难以有效利用各预测模型优势的问题，有学者提出多模型融合的方法，即将多个单一模型进行组合，通过降低单一模型的偏差，发挥各模型的优势，使预测值更加接近于实测值，进而提高径流预测结果的精确度与可靠性^［13］。Zhang等^［14］运用奇异谱分析（SSA）与自回归综合移动平均模型（ARIMA）结合对年径流量进行预测。结果表明组合模型具有更佳的预测能力。Wang等^［15］将人工神经网络模型（ANN）与集成经验模式分解（EEMD）结合用于中长期径流时间序列预测，其结果较单一ANN模型有了明显的改进。Shamseldin等^［16］运用简单平均法（SAM）、加权平均法（WAM）、神经网络方法（NNM）对5个水文模型进行组合用于日径流量预测，结果表明三组组合模型的预测精度均高于单一模型。以上方法均基于确定性的理论，集成不同模型的预测结果，能够提供更精确的综合预测结果，但无法定量评价模型结构的不确定性^［12］。基于贝叶斯理论发展起来的贝叶斯模型平均方法（Bayesian Model Averaging，BMA）通过计算单一模型的后验概率，并基于后验概率对组合模型成员进行加权平均，判断成员模型优劣。因此，BMA可以有效处理组合模型成员的不确定性，在提供精度更高的预测结果的同时，还能提供可靠的预测概率，定量评价模型结构不确定性对预测结果的影响^［17］。近年来，BMA方法也被应用于水文学领域，董磊华等^［18］运用BMA方法综合SIMHYD、SMAR（Soil Moisture Accounting and Routing），以及新安江模型的预报结果，推求出可靠的预报值的概率分布；Jiang等^［19］利用BMA方法集成多卫星降水产品，对湘江日径流进行模拟；Xu等^［20］以中国为研究区，研究了BMA在多模态水文预报后处理中的应用。以上研究均表明BMA能够提高水文预测的准确性与可靠性，还能得到更为优良的预测区间。然而，目前运用BMA方法进行中长期日径流预测的相关研究尚不多见。

疏勒河是河西走廊三大内陆河之一，其上游出山径流量制约着下游玉门、瓜州以及敦煌等地区的经济发展和生态环境变化^［21］。因此，对疏勒河上游径流进行准确地模拟预测对于加强区域内水资源管理、实现水资源高效利用以及促进人与生态可持续发展都具有重要的现实意义^［22］。因此，本文将以疏勒河上游作为研究区，利用疏勒河上游出山口昌马堡水文站2010—2017年实测日径流数据，基于SVM、ELM、MARS模型对疏勒河上游流域未来1~7日径流量进行预测，运用贝叶斯模型平均（BMA）方法对ELM、SVM、MARS模型的预测结果进行组合，构建径流组合预测模型，并提供确定性预测与概率预测结果。

1 研究方法

1.1　极限学习机（Extreme Learning Machine，ELM）

ELM通过对输入权值和偏置随机赋值，并根据“最小二乘法”的原理，运用Moor-Penrose伪逆矩阵计算输出权值，在保证学习和训练精度的前提下，克服了传统神经网络训练速度慢、容易陷入局部最优等缺点。假设n、l、m分别是网络输入层、隐含层以及输出层的节点数，对于N个样本 $(x_{i}, t_{i}), i = 1,2, \dots, N$ ，其中 $x_{i} = {[x_{i 1}, x_{i 2}, \dots, x_{i n}]}^{T} \in R^{n}$ ， $t_{i} = {[t_{i 1}, t_{i 2}, \dots, t_{i m}]}^{T} \in R^{m}$ ，ELM模型可表示为：

o_{j} = \sum_{j = 1}^{l} β_{i} g (ω_{i} \cdot x_{i} + b_{i}), j = 1,2, 3, \dots, N

（1）

其中：g（x）是第l个隐节点的阈值和非线性映射函数； $β_{i}$ 是隐含层和输出层的连接权值向量； $ω_{i}$ 是输入层和隐含层的连接权值向量； $b_{i}$ 是隐含层神经元阈值； $o_{j}$ 是网络输出。此时，上式可以改写为：

H β = T

（2）

其中： $H$ 是隐含层的输出矩阵； $β$ 是输出权值矩阵； $T$ 是样本集目标矩阵。

实现基于ELM预测径流量的过程如下。第一步：确定隐含神经元个数，本研究根据均方误差（Mean Squared Error，MSE）最小原则，设定隐藏神经元个数。随机设定输入层与隐含层的连接权值 $ω_{i}$ 和隐含神经元的阈值 $b_{i}$ 。第二步：选择一个无限可微的函数作为隐含层神经元的激活函数，进而计算隐含层输出矩阵 $H$ 。通过交叉验证方法，选取误差最小的sigmoid函数为激活函数。第三步：计算输出层的连接权值矩阵 $β$ ∶ $β$ = $T$ + $H$ 。

1.2　支持向量机（Support Vector Machine，SVM）

支持向量机是根据统计学理论提出的一种通用学习方法，通过核函数把低维空间中的非线性回归问题映射到高维特征空间中，在高维特征空间中求解凸优化的问题，能很好地解决非线性、高维数、小样本以及局部极小点等问题。主要转换过程如下：

\{\begin{array}{l} m i n (\frac{1}{2} {‖ω‖}^{2} + C \sum_{i = 1}^{l} (ξ_{i} + ξ_{i}^{*})) \\ s . t . \{\begin{array}{l} y_{i} - ω φ (x_{i}) - b \leq ε + ξ_{i} \\ - y_{i} + ω φ (x_{i}) + b \leq ε + ξ_{i}^{*}, i = 1,2, \dots, l \\ ξ_{i} \geq 0, ξ_{i}^{*} \geq 0 \end{array} \end{array}

（3）

最终回归函数为：

f (x) = \sum_{i = 1}^{l} (a_{i} - a_{i}^{*}) K (x_{i}, x) + b

（4）

式中： $\emptyset (x)$ 为非线性映射函数； $ω$ 为超平面法向量；C为惩罚参数；b为超平面偏移量； $ε$ 为线性不灵敏损失函数； $K (x_{i}, y_{i}) = \emptyset (x_{i}) \emptyset (y_{i})$ 是满足Mercer条件的核函数； $a_{i}$ 和 $a_{i}^{*}$ 为二次规划中Lagrange乘子。

SVM核函数的选择直接影响模型径流量预测结果好坏，高斯核函数是一种局部性极强的核函数，面对小样本数据可以使用较少的参数取得较好的分类性能，因此本研究选取高斯径向基核函数（radial）作为核函数。在此基础上，通过交叉验证网格搜索算法确定最优惩罚系数C与核函数参数γ。

1.3　多元自适应回归样条（Multivariate Adaptive Regression Spline，MARS）

多元自适应回归样条模型是一种非参数和非线性回归方法，其原理是进行局部回归建模，通过样条函数估算非线性模型的变化趋势，将数据整体划分为几个子区域，在每个特定的子区域，响应变量用线性回归进行拟合。MARS的优点是能够估算基函数的贡献，并通过模拟变量之间的交互影响和加性来确定对应变量。MARS把基函数作为基本单元，可表示为：

h_{i} (x) = m a x (0, x - c)

（5）

h_{i} (x) = m a x (0, c - x)

（6）

其中： $h_{i} (x)$ 为第i个基函数， $x$ 为输入变量； $c$ 为输入参数的阈值向量。MARS的结构为：

f (x) = β_{0} + \sum_{i = 1}^{m} β_{i} h_{i} (x)

（7）

其中： $f (x)$ 为输出结果； $β_{0}$ 为初始常数； $β_{i}$ 为第 $i$ 个基函数的系数；m为基函数的个数； $h_{i} (x)$ 为第i个基函数。

在建模过程中，模型使用常数值对目标变量进行估计，得到目标变量的平均值，根据最大基函数的个数m，所有可能的基函数都会加入到模型中，从而导致结果出现过拟合现象。为避免这一局限，通过遵循广义交互验证（Generalized Cross Validation，GCV）原则，识别重要性较低的基函数，并去除这些基函数，获得最佳MARS模型。

G C V = \frac{\frac{1}{N} \sum_{i = 1}^{N} {[y_{i} - f (x_{i})]}^{2}}{{[1 - \frac{(H + 1) + d H}{N}]}^{2}}

（8）

式中：N为数据总数；H为基函数的个数；d为惩罚因子。

实现基于MARS预测径流量的过程关键主要为惩罚因子d个数的选择，惩罚因子数量决定了基函数的个数，d越小，生成的MARS模型中基函数越多。d越大，被放置的节点数越少，函数就会更加平滑，但同时被排除在外的基函数也就越多。本研究通过反复训练模型，根据MSE最小原则确定惩罚因子个数。

1.4　贝叶斯模型平均（Bayesian Model Averaging，BMA）

BMA方法对单个模型进行加权平均，生成一个整体的概率分布函数（PDF），通过PDF集成不同模型的预测结果得到的更可靠的综合预测值^［23］。该方法的基本原理为：

设y为径流模拟值， $D = [d_{1}, d_{2}, \dots, d_{r}]$ 为径流实测值， $f = [f_{1}, f_{2}, \dots, f_{k}]$ 为所取K个水文模型组成的模型空间。根据总概率法则，BMA模拟变量y的概率密度函数可以表示为：

p (y | D) = \sum_{k = 1}^{K} p (f_{k} | D) \cdot p_{k} (y | f_{k}, D)

（9）

其中： $p (f_{k} | D)$ 为模型所模拟序列f_k 的后验概率，即在给定的数据条件下模型为最优模型的概率，反映了f_k 与实际径流量的吻合程度。实际上， $p (f_{k} | D)$ 是BMA的权重 $ω_{k}$ ，预测精度越高的模型分配到的权重越大，并且 $\sum_{k = 1}^{K} ω_{k} = 1$ 。本质上，BMA均值预测序列是对权重为 $ω_{i}$ 的不同模型的最优预测序列进行加权平均。 $p_{k (y | f_{k}, D)}$ 是给定模型预测f_k 与数据 $D$ 条件下的预测值Q的后验分布。BMA平均预测值与方差公式如下：

E (y | D) = \sum_{k = 1}^{K} p (f_{k} | D) \cdot E [p_{k} (y | f_{k}, D)] = \sum_{k = 1}^{K} ω_{k} f_{k}

（10）

V a r (y | D) = \sum_{k = 1}^{K} ω_{k} {(f_{k} - \sum_{k = 1}^{K} ω_{k} f_{k})}^{2} + \sum_{k = 1}^{K} ω_{k} σ_{k}^{2}

（11）

其中： $σ_{k}^{2}$ 为给定观测数据D和模型f_k 的条件下模拟变量的方差。BMA预测变量方差包括模型间误差和模型自身误差，式（11）中 $\sum_{k = 1}^{K} {(f_{k} - \sum_{k = 1}^{K} ω_{k} f_{k})}^{2}$ 为不同模型间误差， $\sum_{k = 1}^{K} ω_{k} σ_{k}^{2}$ 为模型自身误差。同单一模型相比，BMA模拟变量方差可以更好地描述预测变量的不确定性。

本研究运用蒙特卡洛组合抽样的方法生成BMA在任意时刻t的预测值的不确定性区间^［24］。步骤如下：

（1）通过单一模型权重［ $ω_{1}, ω_{2}, \dots, ω_{k}$ ］，在［1，2，…，K］中随机生成整数k以抽选模型。具体步骤为： $①$ 设累计概率 $ω_{0}^{'} = 0$ ，计算 $ω_{k - 1}^{'} + ω_{k}$ （k=1，2，…，K）； $②$ 在0~1之间随机生成小数u； $③$ 若满足 $ω_{k - 1}^{'} \leq u \leq ω_{k}^{'}$ ，则选择第k个模型。

（2）通过第k个模型在t时刻的概率分布 $g (Q_{t} | f_{k}^{t} ∙ σ_{k}^{2})$ 随机生成流量值Q_t， $g (Q_{t} | f_{k}^{t} ∙ σ_{k}^{2})$ 是均值为 $f_{k}^{t}$ ，方差为 $δ_{k}^{2}$ 的正态分布。

（3）对步骤（1）（2）重复M次，M是任意t时刻的样本容量，本文M=10 000。

BMA在任意t时刻得到M个样本后，将样本从小到大排序，BMA的95%预测区间即为2.5%与97.5%分位数之间的部分。

BMA均值预测过程主要依赖于Zellner的g先验下的线性模型，为减少先验信息的主观性，本文将参数先验概率分布设定为单位信息先验分布（UIP）。参数先验概率分布设定完成，还需要设置模型先验概率分布，本文选取随机分布（Random Distribution）为先验概率分布。为减少对所有模型平均计算的繁杂过程，本文选取蒙特卡洛抽样方法抽取后验概率较高的模型。

2 研究区概况与模型构建

2.1　研究区概况

疏勒河发源于托勒南山与疏勒南山之间的沙果林那木吉木岭。疏勒河上游（图1）主要包含了出山口昌马堡以上的山区（96.6°~99.0° E，38.3°~39.9° N），流域面积1.14×10⁴ km²，涵盖了甘肃省酒泉市肃北蒙古族自治州县、青海省西蒙古族自治州天骏县等地区。地形上由托勒南山、疏勒南山以及疏勒河谷组成。河谷地区相对平坦，山区地势高大险峻、地形陡峭，海拔约在2 100~5 750 m，4 500 m以上广泛分布着现代冰川，冰川资源丰富^［25］。

图1

新窗口打开| 下载原图ZIP| 生成PPT

图1 研究区位置

Fig. 1 Location of the study area

疏勒河出山径流的年内分配差异较大，受气候影响，径流多集中在5月至9月，占全年径流量的80%左右，最大日径流可达415 m³·s^-1，10月到次年4月径流量则明显偏少。其中，冬季降水常以固态形式存在，每年12月到次年3月日径流量较少，甚至无径流。4、5月气温逐渐回升，积雪和冰川融化使径流量增加，随之进入一年中降水最为集中的季节，加上高山冰雪融水的增加，7、8月份出山径流量达到最大值，10月之后随着气温的降低和暖湿气流影响的削弱，径流量逐渐减少^［26］。

2.2　数据获取与处理

本文选用昌马堡水文站2010—2017年日径流量作为疏勒河出山径流分析数据。其中，2010—2014年数据用于模型训练，2015—2017年数据用于模型测试。疏勒河上游出山口径流量数据基本统计资料见表1。由表1可知，径流量数据呈现明显的高偏态分布（训练数据集为2.39，测试数据集为2.81），并且不同数据集日径流数据的变异系数均大于1，这可能会对预测精度产生一定的影响。为了保证对模型的有效训练以及有效特征的提取，通常对输入数据进行归一化处理。本文采用最小-最大归一化方法对径流数据时间序列进行归一化处理。公式为：

Q^{*} = \frac{Q - Q_{m i n}}{Q_{m a x} - Q_{m i n}}

（12）

式中：Q^* 为经标准化处理后的径流序列；Q为原始径流序列；Q_min为原始最小径流量；Q_max为原始最大径流量。

表1 昌马堡日径流量基本统计值

Table 1 Basic daily runoff statistics of Changmapu

样本	平均值/（m³·s^-1）	最大值/（m³·s^-1）	标准差/（m³·s^-1）	变异系数	偏度
总体	45.45	415.00	51.44	1.13	2.59
训练期	42.39	348.00	47.98	1.13	2.39
测试期	50.77	415.00	56.47	1.11	2.81

新窗口打开| 下载CSV

2.3　模型构建

径流时间序列能够揭示其随时间变化的规律，并将这种规律延伸到未来。因此，径流时间序列可以作为径流预测的基本输入因素。本文通过计算2010—2017年实测日径流的自相关函数（ACF）和偏自相关函数（PACF）来确定模型的输入。如图2所示，径流时间序列的ACF与PACF均表现出拖尾特征，ACF逐渐衰减并在95%的置信区间内表现出明显的自相关性，特别是在时滞数小于5时，ACF大于0.8，表现出显著的自相关性。同时，PACF显著不为0的时滞数为1，2，3。因此，综合ACF与PACF分析结果，确定Q_t_-3、Q_t_-2、Q_t_-1为模型输入量。

图2

新窗口打开| 下载原图ZIP| 生成PPT

图2 实测径流时间序列的ACF和PACF

Fig. 2 ACF and PACF of measured runoff time series

构建一个包含7列数据的表格，表格形式为（Q_t_-3，Q_t_-2，Q_t_-1，Q_t₊₁，Q_t₊₃，Q_t₊₅，Q_t₊₇），其中，Q_t_-3、Q_t_-2、Q_t_-1为输入变量，分别表示t-3、t-2、t-1日的径流量组合；Q_t₊₁、Q_t₊₃、Q_t₊₅、Q_t₊₇表示模型在t+1、t+3、t+5、t+7日的预测径流量，即模型的输出。

对于单一模型的构建主要分为两部分：第一部分为确定模型的输入向量；第二部分为利用确定的输入向量设置模型运行的最优基本参数，建立中长期日径流预测模型。对于BMA组合模型的构建，主要是利用各模型后验概率将单一模型的预测结果组合起来得到集合预测结果，建立中长期日径流组合预测模型。此外，本文采用马尔科夫蒙特卡洛抽样对集合预测结果进行10 000次抽样，获取集合预测结果的95%置信区间，对预测结果进行不确定性分析。为验证BMA集合预测的有效性，在相同预见时间的条件下，与单一ELM、SVM、MARS模型进行对比研究。本文主要研究技术路线如图3所示。

图3

新窗口打开| 下载原图ZIP| 生成PPT

图3 研究技术路线

Fig. 3 Research technical route

2.4　模型评价指标

模拟效率可体现模型在研究区的适用情况，本文选取相关系数R、纳什效率系数NSE（Nash-Sutcliffe efficient）以及均方根误差RMSE（Root-Mean-Square Error）对预测结果进行评价，表达式如下：

R = \frac{\frac{1}{n} \sum_{i = 1}^{n} (Q_{o b s} - {\bar{Q}}_{o b s}) (Q_{f} - {\bar{Q}}_{f})}{\sqrt[]{\sum_{i = 1}^{n} {(Q_{o b s} - Q_{o b s})}^{2}} \sqrt[]{\sum_{i = 1}^{n} {(Q_{f} - Q_{f})}^{2}}}

（13）

N S = 1 - \frac{\sum_{i = 1}^{n} {(Q_{o b s} - Q_{f})}^{2}}{\sum_{i = 1}^{n} {(Q_{o b s} - {\bar{Q}}_{o b s})}^{2}}

（14）

R M S E = \sqrt[]{\frac{1}{n} \sum_{i = 1}^{n} {(Q_{o b s} - Q_{f})}^{2}}

（15）

式中：Q_obs为实测值（m³·s^-1）；Q_f为预测值（m³·s^-1）； ${\bar{Q}}_{o b s}$ 为实测均值（m³·s^-1）； ${\bar{Q}}_{f}$ 为预测值均值（m³·s^-1）；n为实测数据数量。

R是判断实测值与预测值的线性相关程度。R取值范围是-1至1，R的绝对值越接近于1，线性相关程度越高；NSE代表预测值与实测值在1∶1水平上的相似度，NSE在 $- \infty$ 至1.0之间取值，NSE≥0.5时预测结果可接受，NSE越接近1，拟合效果越好；RMSE反映预测值与实测值之间的偏差，该值越接近于0，模型的拟合度越高。当R=1、NSE=1、RMSE=0时，认为该模型为最佳模型^［27］。

评价预测结果的不确定性区间的优良性可以体现模型的可靠性。本文采用覆盖率（CR）、平均带宽（B）以及平均偏移幅度（D）三个主要指标来分析比较BMA的不确定性区间^［28］。

C R = \frac{\sum_{t = 1}^{N T} I [Q_{o b s} (t)]}{N T}

（16）

B = \frac{1}{N T} \sum_{T = 1}^{N T} (q_{u}^{t} - q_{l}^{t})

（17）

D = \frac{1}{N T} \sum_{t = 1}^{N T} |\frac{1}{2} (q_{u}^{t} + q_{l}^{t}) - Q_{o b s}^{t}|

（18）

其中： $I [Q_{o b s} (t)] = \{\begin{array}{l} 1, & q_{l}^{t} \leq Q_{o b s} \leq q_{u}^{t} \\ 0, & 其他 \end{array}$ ；N为实测数据数量；T为测量的数据所对应的时间； $q_{u}^{t}$ ， $q_{l}^{t}$ 分别为t时刻的预测区间的上界和下界； $Q_{o b s}^{t}$ 为t时刻的实测径流。

覆盖率是评价模拟区间最常用的指标，指预测区间所覆盖的实测流量数据的比率。CR值越大，表明模拟区间的覆盖率越高，预测结果包含的真实信息越多；平均带宽B是在指定的置信水平，以及较高的覆盖率保证的前提下，预测区间的平均带宽越窄预测区间的不确定程度越小，反之，预测的不确定性程度越大；平均偏移幅度D是评价预测区间的中心线相较实际径流过程线偏离程度的指标，理论上，平均偏移幅度越小，表明预测区间的对称性越好。

3 结果分析

3.1　单一模型预测精度分析

为了验证单一模型预测的有效性，表2~4给出了最优参数下三个单一模型在不同时间的预测结果评价。由表2~4可知，随着时间的增加，三个单一模型的预测精度逐渐降低。但总体上三个单一模型的预测值与实测值的R值均在0.75以上，表明预测值与实测值之间存在着很高的线性相关关系。NSE值均在0.55以上，说明预测值与实测值之间的拟合程度是可以接受的。因此，说明所构建的单一模型的预测结果合理、有效，可以作为BMA集合预测的有效成员。

表2 ELM模型参数及预测结果评价

Table 2 Parameter and evaluation of prediction results of ELM model

模型	参数	训练集			测试集
模型	神经元个数	R	NSE	RMSE/（m³·s^-1）	R	NSE	RMSE/（m³·s^-1）
ELM_（_t_+1）	15	0.959	0.919	13.642	0.964	0.928	15.135
ELM_（_t_+3）	50	0.865	0.748	24.061	0.890	0.784	26.225
ELM_（_t_+5）	50	0.809	0.654	28.209	0.823	0.669	32.463
ELM_（_t_+7）	50	0.781	0.610	29.959	0.781	0.603	35.584

新窗口打开| 下载CSV

表3 SVM模型参数及预测结果评价

Table 3 Parameter and evaluation of prediction results of SVM model

模型	参数		训练集			测试集
模型	核函数参数（γ）	惩罚因子（C）	R	NSE	RMSE/（m³·s^-1）	R	NSE	RMSE/（m³·s^-1）
SVM_（_t_+1）	0.250	16.0	0.971	0.919	13.223	0.967	0.935	14.435
SVM_（_t_+3）	0.250	0.5	0.868	0.746	24.180	0.886	0.758	27.768
SVM_（_t_+5）	0.004	128.0	0.811	0.640	28.792	0.837	0.655	33.157
SVM_（_t_+7）	0.004	256.0	0.778	0.595	30.543	0.783	0.582	36.494

新窗口打开| 下载CSV

表4 MARS模型参数及预测结果评价

Table 4 Parameter and evaluation of prediction results of MARS model

模型	参数	训练集			测试集
模型	惩罚因子（C）	R	NSE	RMSE/（m³·s^-1）	R	NSE	RMSE/（m³·s^-1）
MARS_（_t_+1）	3	0.972	0.944	11.336	0.980	0.961	11.131
MARS_（_t_+3）	4	0.874	0.764	23.309	0.842	0.708	30.491
MARS_（_t_+5）	4	0.811	0.658	28.074	0.802	0.639	33.927
MARS_（_t_+7）	4	0.786	0.612	29.657	0.781	0.604	35.503

新窗口打开| 下载CSV

3.2　BMA均值预测与单个模型的预测值对比

本研究得到的BMA方法在测试期的（t+1） d、（t+3） d、（t+5） d与（t+7） d的预测结果如图4与表5所示。通过水文过程曲线可以看出（图4），BMA的整体预测效果比较好，预测值与实测值变化趋势一致，洪峰出现时间同实际观测情况一致。BMA对径流低值的预测结果较好，但对径流高值的预测有一定的误差，并且随着预测时间的增加，峰值预测误差逐渐增大。这被认为是目前径流研究工作中的一个重要局限性，在其他径流模拟预测研究中也出现了同样的局限性^［9-10］。尽管如此，从散点图可以看出，BMA在（t+1）~（t+7） d上的中长期日径流的R均大于0.8，表明预测精度是可以接受的^［29］。另外，由表5可知，BMA的预测精度随着预测时间的增加而逐渐减小，但BMA在（t+1）~（t+7） d的R值均大于0.78，表明BMA的预测值与实测值之间存在较高的线性相关关系；NSE值均大于0.6，表明预测值与实测值之间的拟合效果较好。综上，说明BMA可以用于（t+1）~（t+7） d的中长期日径流预测，并具有较高的预测精度。然而，BMA在（t+1） d后训练期的R值不如测试期的R值，而训练期的RMSE值却优于训练期。由表1可知，训练期的变差系数高于测试期，训练期数据离散程度大，导致训练期R值低于测试期；而测试期数据的偏度大于训练期，分布不均匀，产生的较大误差。

图4

新窗口打开| 下载原图ZIP| 生成PPT

图4 BMA在测试期（t+1）~（t+7）日径流量实测值与预测值对比与95%置信区间

Fig. 4 The measured and predicted runoff of （t+1）~（t+7） days were compared with 95% confidence interval by BMA

表5 BMA在（t+1）~（t+7） d日径流量预测结果评价

Table 5 Evaluation of daily runoff predicted results of BMA in （t+1）~（t+7） days

时间/d	训练期			测试期
时间/d	R	NSE	RMSE/（m³·s^-1）	R	NSE	RMSE/（m³·s^-1）
t+1	0.979	0.959	9.775	0.982	0.964	10.708
t+3	0.874	0.764	23.294	0.894	0.799	25.306
t+5	0.814	0.661	27.898	0.841	0.708	30.570
t+7	0.786	0.618	29.654	0.800	0.640	33.878

新窗口打开| 下载CSV

由表2~5可知，BMA与组成BMA的3个单一模型均可适用于（t+1）~（t+7） d的中长期日径流预测，但进一步对比分析可知，BMA的预测结果较单一模型具有更高的预测精度。以训练期和测试期的（t+1） d与（t+7） d预测结果为例，训练期（t+1） d时BMA的R值均为0.979，高于SVM（0.971）ELM（0.959）和MARS（0.972）的R值。BMA的NSE值为0.982，高于单一ELM（0.919）、SVM（0.919）、MARS（0.944）的NSE值。BMA的RMSE值较单一ELM、SVM、MARS的RMSE值分别提高了28.35%、26.08%、13.77%。在训练期（t+7） d时，BMA与MARS的R值（0.786）相同，高于单一ELM（0.781）和SVM模型（0.778）。BMA的NSE值为0.618，分别高于三个单一模型的NSE值（ELM为0.610、SVM为0.595、MARS为0.612）BMA的RMSE值较单一ELM、SVM、MARS模型分别提高了4.79%、7.17%、4.58%。

在测试期（t+1） d时，BMA与SVM的R值相同（0.967），较ELM、MARS的R值分别提高了1.83%、0.6%；NSE值较单一ELM、SVM、MARS模型分别提高了3.7%、3.1%、0.3%；RMSE值较单一ELM、SVM、MARS模型分别降低了33.7%、25.8%、3.8%。在测试期（t+7） d时，ELM、SVM、MARS的R值分别为0.781、0.783、0.781，BMA的R值为0.800，高于单一模型。BMA的NSE值为0.64，较单一模型ELM、SVM、MARS的NSE值分别提高了5.8%、9.1%、5.6%。另外，ELM、SVM、MARS的RMSE值分别为35.548、36.494、35.503，均大于BMA的RMSE值33.878。因此，综合R、NSE、RMSE的评价结果可以看出，BMA较单一模型无论是在（t+1）~（t+3） d的短期日径流，还是在（t+5）~（t+7） d的中长期日径流预测均具有较高的预测精度，能够提供更准确的预测结果。

对不同等级的径流量进行预测对及时发布极端水文事件预警信息以及制定防治预案具有重要意义。本文参考董磊华等^［18］研究成果，根据疏勒河的径流量特征，将径流量分为3个等级：流量从大到小排序，将前10%的大流量定为高水，中间50%的流量定为中水，后40%的小流量定为低水。对比分析BMA与单一模型ELM、SVM、MARS在不同等级径流量的预测结果（表6）。从表6可以看出，BMA与单一模型对中水预测效果最好，BMA与单一模型均可预测出未来5~7 d中水径流量；低水次之，除单一SVM模型可以预测未来3日的中水径流量，其余模型均只能预测未来1日低水径流量；高水最差，仅BMA与SVM模型能预测未来1日的高水径流量。但相比之下，BMA较单一模型的预测效果更好，这是因为BMA方法综合发挥了单一模型各自的优势，在一定程度上克服了单一模型无法避免的局部极值问题^［28］。

表6 BMA和组成它的3个模型的预测值在3个不同流量等级的统计结果

Table 6 Statistical results of BMA and the forecast values of the three models that comprised it at three different flow levels

时间	指标	模型	训练期			测试期
时间	指标	模型	高水	中水	低水	高水	中水	低水
t+1	R	ELM	0.741	0.930	0.779	0.781	0.952	0.842
		SVM	0.830	0.935	0.764	0.781	0.946	0.838
		MARS	0.750	0.924	0.766	0.689	0.944	0.817
		BMA	0.835	0.940	0.765	0.949	0.958	0.962
	NSE	ELM	0.425	0.857	0.561	0.506	0.894	0.671
		SVM	0.625	0.874	0.584	0.645	0.874	0.694
		MARS	0.412	0.851	0.411	0.033	0.887	0.667
		BMA	0.835	0.941	0.769	0.523	0.891	0.832
	RMSE/（m³·s^-1）	ELM	36.657	9.659	1.971	40.675	9.319	2.126
		SVM	30.624	8.796	2.228	1.964	18.664	21.803
		MARS	34.611	9.897	2.021	46.892	9.368	2.210
		BMA	28.506	8.312	1.987	38.659	9.021	2.012
t+3	R	ELM	0.359	0.856	0.632	0.501	0.879	0.701
		SVM	0.944	0.843	0.800	0.811	0.870	0.899
		MARS	0.465	0.856	0.626	0.236	0.887	0.701
		BMA	0.459	0.891	0.834	0.478	0.871	0.676
	NSE	ELM	0.109	0.637	0.122	0.234	0.729	0.357
		SVM	0.718	0.615	0.335	0.656	0.446	0.647
		MARS	0.171	0.644	0.345	0.036	0.693	0.443
		BMA	0.172	0.670	0.405	0.228	0.729	0.375
	RMSE/（m³·s^-1）	ELM	64.431	18.944	3.381	65.412	17.287	3.966
		SVM	20.648	25.508	23.309	3.955	31.584	24.800
		MARS	59.235	18.651	3.231	82.461	18.015	4.012
		BMA	58.375	17.998	3.114	62.357	16.986	3.554
t+5	R	ELM	0.321	0.798	0.556	0.367	0.848	0.642
		SVM	0.285	0.786	0.513	0.321	0.855	0.577
		MARS	0.315	0.801	0.549	0.235	0.850	0.642
		BMA	0.326	0.812	0.561	0.392	0.862	0.621
	NSE	ELM	0.101	0.534	0.288	0.152	0.656	0.161
		SVM	0.081	0.554	0.263	0.103	0.652	0.163
		MARS	0.111	0.532	0.279	0.059	0.641	0.375
		BMA	0.126	0.611	0.306	0.158	0.734	0.191
	RMSE/（m³·s^-1）	ELM	70.603	22.687	4.302	76.542	21.015	4.736
		SVM	80.281	22.225	1.146	98.263	21.808	5.121
		MARS	69.375	22.325	4.261	92.162	19.663	5.012
		BMA	70.501	17.865	3.996	75.826	15.568	4.346
t+7	R	ELM	0.271	0.752	0.365	0.278	0.823	0.601
		SVM	0.310	0.733	0.528	0.237	0.812	0.599
		MARS	0.315	0.755	0.532	0.316	0.757	0.534
		BMA	0.331	0.864	0.641	0.290	0.827	0.610
	NSE	ELM	-2.312	0.494	0.012	-1.013	0.819	0.606
		SVM	-3.894	0.470	-1.436	-2.060	0.655	-0.084
		MARS	-2.892	0.452	0.286	-0.401	0.581	0.305
		BMA	-2.007	0.534	0.108	-1.397	0.581	0.320
	RMSE/（m³·s^-1）	ELM	74.316	22.732	4.301	85.671	22.125	5.046
		SVM	82.903	24.524	4.071	106.930	24.187	4.855
		MARS	73.015	23.752	4.406	95.587	20.989	5.405
		BMA	72.309	20.067	4.112	84.032	18.032	5.062

注：黑体为最优值。

新窗口打开| 下载CSV

3.3　BMA径流预测不确定性分析

对预测结果的不确定性分析可以判断模型结构的可靠性，对提高水资源管理的科学性、极端水文事件风险防范能力具有重要意义^［29］。本文基于覆盖度（CR）、平均区间宽度（B）以及平均偏移幅度（D）3个指标对BMA的95%置信区间的径流预测结果的不确定性进行分析。从图4可看出，大多数实测值都在不确定性区间之内，测试期内BMA在（t+1）~（t+7） d的95%置信区间的覆盖度分别为达到了94.7%、94.1%、93.8%、92.9%，几乎覆盖了整个实测径流序列，表明95%置信区间预测效果好，不确定性较小。

由于BMA对高水和低水的预测出现了较大偏差，存在较大的不确定性。因此，为进一步研究模型结构的不确定性，对每个等级径流预测结果进行评价。表7综合了BMA在3个不同流量等级的预测区间优良性，可以看出BMA预测区间覆盖了整个低水部分，中水部分除在（t+7） d时刻上覆盖率较低，其余时间均在95%以上。但是测试期与训练期在高水部分的覆盖率均不足65%，其中测试期（t+7） d的高水部分最低覆盖率仅43.52%，拉低了整体覆盖率。同时，BMA 95%置信区间的平均带宽在低水部分最小，在高水部分最大，表明BMA在高水部分的预测不确定性更大。就平均偏移幅度来说，BMA在高水部分的平均偏移幅度最大，中水部分最小，说明BMA在高水部分的预测值偏离实测值的程度最大，中水部分的预测值偏离实测值的程度最小。综上，BMA的预测不确定性主要来源于高水部分，主要原因在于高水一般出现在汛期，径流量受降水、融水影响较大，预测难度较大，因此高水预测的不确定性较大。另外，各流量等级的区间不确定性均随着时间增加而增大。但总体而言，BMA的95%置信区间效果良好，对实测径流量覆盖率高。

表7 BMA 95%在3个不同流量等级的统计结果

Table 7 Statistical results of BMA 95% confidence interval at three different traffic levels

时间/d	径流等级	训练期			测试期
时间/d	径流等级	CR/%	B/（m³·s^-1）	D/（m³·s^-1）	CR/%	B/（m³·s^-1）	D/（m³·s^-1）
t+1	高水	60.01	51.21	23.89	45.51	50.98	32.87
	中水	97.42	49.63	5.06	98.62	51.32	5.32
	低水	100	40.72	5.64	100	41.03	5.98
t+3	高水	65.03	91.05	43.01	57.61	90.66	58.84
	中水	95.63	78.45	14.32	97.66	83.52	11.64
	低水	100	61.62	15.61	100	61.34	15.99
t+5	高水	57.14	98.68	55.43	47.53	110.52	71.56
	中水	95.03	92.34	19.72	97.42	97.32	15.78
	低水	100	72.67	20.43	100	72.59	21.23
t+7	高水	61.86	113.87	57.69	43.52	114.32	78.49
	中水	74.32	94.32	21.55	97.65	104.89	18.56
	低水	100	73.01	20.81	100	74.68	21.83

注：黑体为最优值。

新窗口打开| 下载CSV

4 讨论

李洪源等^［30］在运用分布式水文模型SPHY对疏勒河上游径流进行了逐日与逐月径流预测，并取得了较好的预测精度。但SPHY在率定期的日径流NSE值为0.62，低于本文BMA在训练期（t+1） d的NSE值0.943。SPHY在测试期中的日径流的NSE值为0.79，低于本文BMA在测试期（t+1） d的NSE值（0.935），表明BMA方法具有较高的精度，这可能是因为BMA方法通过对各单一模型的后验概率进行加权平均，有效解决了单一模型的不确定性，把多个单一模型的优点组合起来，从而获得了更准确的预测值。

本文通过计算径流的ACF与PACF确定模型的输入，但如何利用径流本身的自相关关系确定模型输入暂无统一标准。对于自相关关系差的小流域来说，可以针对流域自身特点来确定模型的输入。如于海姣等^［31］在对小流域排露沟日降水-径流模拟研究中，针对降水-径流的非线性关系，选择相对百分比误差方法来确定模型的输入。Gavin等^［32］根据输入与输出之间的显著线性关系，运用互相关法确定模型输入，对澳大利亚Murray河流的含盐量进行模拟。

BMA方法通过集成各个模型的自身优势，弥补单个模型的不足获得较单一模型更高的预测精度。在疏勒河上游流域（t+1）~（t+7） d的中长期日径流预测中，BMA方法较单一模型的预测效果更好，精度更高，为缺乏资料的干旱半干旱地区的径流预测提供了有效的预测工具。但是，BMA在中长期日径流预测中仍然存在着一定的不足。BMA在高径流量和低径流量预测上表现出较单一模型更高的精度，但是预测值仍然小于实测值，存在着较大的误差。误差原因可能是：疏勒河上游流域冰川积雪、冻土分布广，由此形成了独特的寒区水文过程，流域内的温度、蒸发、冻融作用都是影响高径流量时期径流变化的重要因素，所以，只将径流数据作为模型输入过于单一，导致预测序列出现偏差^{［31，33］}。低径流时期的径流量数值较小且存在0值，经归一化处理后的径流数据产生较多的0值，导致模型在运行过程中出现一定的误差^［31］。因此，在日后的工作中可以尝试从增加模型输入、优化模型参数计算等方面提高模型预测精度，更真实地反映径流变化规律。

5 结论

本文利用BMA方法对疏勒河上游2010—2017年的（t+1）~（t+7） d中长期日径流进行了预测，对BMA与组成BMA的ELM、SVM、MARS模型的预测结果进行分析对比，并对BMA预测结果的不确定性进行定量分析。得到以下结论：

（1）BMA在（t+1）~（t+7） d的径流预测结果比单一ELM、SVM、MARS结果准确度更高，不同时期径流量的预测值与实测值的相关系数均在0.7以上，纳什系数均在0.55以上，能够再现径流量的变化趋势，说明BMA方法具有较好的实用性，在中长期日径流预测中可以提供更准确的径流量预测值。BMA对高水、低水的预测值与实测值的相关系数、纳什系数在（t+1） d后均低于0.5，但仍较单一模型保持较高精度。

（2）贝叶斯模型加权平均（BMA）方法是一种通过综合多个模型预测结果的后验分布来推断预测结果的可靠概率分布分析的工具。BMA对（t+1）~（t+7） d的95%置信区间实现了对低水径流量的100%全覆盖，对中水径流量的覆盖率在70%以上，对高水径流量的覆盖率较低，但整体覆盖率均在90%以上。表明BMA不仅可以提供更准确的综合预测结果，还能提供一个综合预测区间来进行模型结果不确定性分析，可成为资料有限条件下的径流模拟预测的有效工具。

参考文献

原文顺序

文献年度倒序

文中引用次数倒序

被引期刊影响因子

[1]

Wagena

M B

， Goering

， Collick

A S

， et al.

Comparison of short-term streamflow forecasting using stochastic time series， neural networks， process-based， and Bayesian models

［J］. Environmental Modelling & Software， 2020， 126： 104669.