基于贝叶斯网络的流域内水文事件丰枯遭遇研究

吕振豫　穆建新　王富强刘姗姗

摘要：基于贝叶斯网络理论结合Copula函数建立了东江流域上、中、下游三个站点降雨、径流丰枯遭遇的风险管理模型，直观地描述了各个站点间降雨、径流的相互关系。利用Copula函数建立联合分布模型计算了站点间不同丰枯组合状态调水不利情况的风险概率。通过贝叶斯网络结构的反向推理功能，以后验知识作为输入，进一步对流域内未来调水可能面临的丰枯遭遇情况进行了仿真计算。结果表明，系统中一个节点的变化将会对其他节点的概率值产生巨大影响，以上、中游降雨为枯水情况作为后验输入，站点间丰枯遭遇调水不利风险概率增幅至55%以上；仿真结果可以为流域内调水方案的制定提供理论支撑。

关键词：东江流域；Copula函数；贝叶斯网络；丰枯遭遇；风险分析

中图分类号：P333.9；X820.4 文献标志码：A 文章编号：

16721683（2016）05001808

水文事件丰枯变化的差异性及不确定性，直接影响流域内调、受水的可控水量，对工程调水风险、流域水资源配置产生重大影响。近年来，不同区域水文事件的丰枯遭遇问题引起了广泛关注[14]，分析方法层出不穷。郑红星[5]和韩宇平[6]等通过建立联合分布模型采用统计方法计算了南水北调调、受水区水文变量间的丰枯遭遇概率，这种统计方法虽然直观易懂，但不适用于多维随机变量遭遇研究，没有考虑各变量间的相关关系。Copula函数作为一种新兴的方法弥补了统计方法的不足，是一种将随机变量联合分布与各自边缘分布相结合的理论联合分布。目前基于二维及多维Copula函数建立丰枯遭遇联合分布模型已得到大量实践应用[710]。此外，传统风险分析方法只考虑采用降雨、径流事件的先验概率计算风险，并没有考虑后验知识；贝叶斯网络作为一种基于贝叶斯条件概率的风险分析方法，既考虑了先验风险概率，又可以利用后验信息进行仿真模拟[11]，在水文事件丰枯遭遇风险分析中有其显著的优势。康玲等[12]运用贝叶斯网络理论建立了南水北调中线水源区与受水区降水丰枯遭遇风险管理模型，充分利用贝叶斯网络的情景仿真和后验推理功能，模拟了不同丰枯组合对调水的影响；Daniel[13]利用并验证了非参数蒙特卡洛贝叶斯理论在洪水频率分析中的独特优势；H.van de Vyver[14]利用贝叶斯网络的后验推理功能对降雨极值强度持续时间频率相关关系进行了分析，并与传统方法进行对比验证了贝叶斯方法的优越性。本文以东江流域上游龙川站、中游河源站及下雨博罗站实测降雨、径流数据为研究对象，运用Copula函数结合贝叶斯网络建立流域内丰枯遭遇风险管理模型，为流域内水资源优化配置决策提供理论支持，具有重要的理论意义和实践意义。

1 数据和方法

1.1 数据来源

本文所用数据均由广东省气象局提供，包括东江流域上游龙川站、中游河源站以及下游博罗站三个代表性水文站点1956年-2005年年尺度实测降雨、径流资料。流域概况及水文站点分布情况见图1。

1.2 研究方法

1.2.1 二维Copula函数理论

Copula函数是定义在[0，1]区间上均匀分布的多维联合分布函数，其主要构造形式如下：

1.2.3 多维Copula函数

三维及多维Copula函数的构造形式与二维情况类似，篇幅限制这里不再赘述，其构造形式详见文献[17]。对于多维Copula函数的参数估计（以三维为例），不能直接利用参数与变量的Kendall秩相关[CM（22]系数的关系直接计算，一般采用两阶段极大似然进

1.2.4 贝叶斯网络概述

贝叶斯网络[18]简称BN，是基于贝叶斯定理和条件概率建立的有向无环图，它由代表变量的节点和连接各节点的有向边构成，以图形化的形式直观地表达系统内各元素之间的相互影响关系。一个完整的贝叶斯网络模型由网络结构和网络参数两部分组成。图2所示为一个6节点贝叶斯网络模型，可用N=〈〈V，E〉，P〉表示[19]，其中：

（1）〈V，E〉表示网络结构有向无环图。图中节点V={V1，V2，…，V6}表示变量，节点间的有向边E代表变量间的相关关系。对于有向边（Vi，Vj），Vi称为Vj的父节点，Vj则为Vi的子节点，没有父节点的称为根节点（V1），没有子节点的称为叶节点（V6）。规定Vi的父节点集合及非后代节点集合分别用fa（Vi）和A（Vi）表示，则贝叶斯网络系统包含如下条件假设：

P（Vi|fa（Vi），A（Vi））=P（Vi|fa（Vi））[JY]（6）

即，在给定父节点情况下，子节点与其非父节点条件独立。

（2）P表示系统根节点概率和非根节点条件概率。由条件独立性假设可知，非根节点的条件概率分布可用P（Vi|pa（Vi））表示，表达了节点与其父节点的相关关系。给定根节点先验概率及非根节点条件概率分布，可据此计算包含所有节点的联合概率分布，图2所示包含全部节点的联合概率分布函数为：

2 结果分析

2.1 多情境丰枯遭遇组合概率

流域内降雨、径流的丰枯变化受地区气候特性

及下垫面等的影响，不同区间的遭遇概率往往是随机的。上游出现某一量级枯水年情况时，中、下游可能出现不同量级的丰、枯情况。为定量分析流域内不同区间降雨、径流的丰枯遭遇情况，计算其遭遇概率及条件概率，文中采用Copula函数，通过建立丰枯遭遇理论联合分布模型进行研究。频率分析中我国一般采用PⅢ曲线作为降雨、径流的边缘分布，其构造形式如下：

利用线性矩法计算得到流域上、中、下游三个站点降雨、径流的参数估计结果见表2。分析可知，流域内年降雨量上游龙川站最小为1 61018 mm，中游河源站最大达1 85065 mm；径流量从上游到下游程递增趋势，变异系数分别为04、062和022，说明径流变化在空间上属于中等变异。

采用单参数ArchimedeanCopula函数建立降雨、径流丰枯遭遇的联合分布函数，不同Copula函数的拟合优度情况与随机变量的相关性有关[20]。考虑到流域内降雨、径流存在较强正相关性，结合几种Copula函数的构造形式及其参数估计方法，文中选取Clayton Copula函数作为两站点降雨、径流丰枯遭遇联合分布模型的构造函数；三站点丰枯遭遇组合则选取GumbelHougaard Copula函数构造；各站点组合copula联合分布函数参数估计结果见表3。

2.2 丰枯遭遇条件概率

基于贝叶斯网络建立丰枯遭遇的风险管理模型，需要确定站点间降雨、径流丰枯遭遇的条件概率作为网络的初始参数。以已知上游来水情况，中、下游遭遇条件概率计算为例，公式如下：

采用GumbelHougaard Copula函数构造三站点降雨、径流丰枯遭遇理论联合分布，其与经验频率的拟合情况见图5、图7，观察可知，理论分布与经验频率拟合情况良好。图6、图8所示为上游龙川站降雨、径流为枯水情况下，中游河源站与下游博罗站不同量级降雨、径流遭遇的条件概率等值线，可从图中直接定量出某一固定值丰枯遭遇条件概率大小。如，上游径流为枯水情况，中、下游同为枯水时的概率P（Y≤y37.5%，Z≤z37.5%|X≤x37.5%）=P（Y≤128.03，Z≤210.66|X≤55.36）=0.685。此外，对比分析图4和图8，已知上游径流为枯水情况，中、下游河源、博罗站不同量级径流遭遇条件概率等值线有明显的前移趋势，等值线密集程度增加，中、下游径流同丰概率由0260降低到0003，同枯风险概率由原来的0293增加到现在的0685，增长了近25倍，调水不利风险急剧增大。

2.3 流域内丰枯遭遇风险管理模型建立

利用贝叶斯网络构造流域内丰枯遭遇的风险管理模型一般包括三个步骤：（1）确定网络结构；（2）确定初始网络参数；（3）根据后验知识进行仿真推理。

2.3.1 初始网络结构确定

对于网络结构的确定有两种方法[22]，一种是利用大量的实测数据通过不同的优化算法进行结构学习，确定最优网络结构模型，这种方法是在大量数据的基础上实现的；第二种方法是根据专家知识结合数据间的相关关系，直接勾画出网络结构，这种方法对数据的多少要求不高，计算简便。本文使用第二种方法构造流域内降雨、径流丰枯遭遇的贝叶斯网络结构模型，结合流域上、中、下游站点降雨、径流的丰枯遭遇关系，确定贝叶斯网络结构模型见图9。网络结构图充分说明了上、中、下游降雨、径流之间的内在联系，上游龙川站降雨为根节点A1，其子节点包括河源站降雨B1、博罗站降雨C1以及龙川站径流A2；叶子节点D1、D2、D3、D4为各站点径流丰枯遭遇对调水不利的风险概率。

2.3.2 确定网络参数

贝叶斯网络结构的参数包括各节点的先验概率及非根节点的条件概率。根据实测资料计算得各站点降雨、径流发生丰、平、枯情况的概率，结合Copula函数计算得到各种组合情况的丰枯遭遇概率，输入到网络结构中得到贝叶斯网络初始模型详见图10。鉴于东江流域内已建成枫树坝水库、新丰江水库等具有年调节功能的大型水库，规定，两个或两个以上站点同为枯水年是对调水构成风险的情况，其中，三站同枯为调水最不利情况。分析图10可知，河源站降雨枯水情况发生概率为40%，易于丰水情况发生；博罗站降雨枯水情况发生概率36%，难于丰水情况发生，说明流域内降雨在空间上存在差异性；流域内径流三站点丰、平、枯水情况发生概率大致相同；龙河径流遭遇调水不利风险概率27%，龙博、河博径流遭遇调水不利风险概率均为29%，调水最不利情况（三站同枯）发生概率26%，表明流域内基本可以实现正常调水。

2.3.3 仿真模拟

利用贝叶斯网络构造流域内降雨、径流的丰枯遭遇风险管理模型，不仅可以通过实测数据计算各站点丰枯遭遇调水不利情况的发生概率（先验概率），还可以通过预测某一节点或者某几个节点的发生情况，作为后验信息输入到网络模型中，利用贝叶斯网络的反向推理功能，推测出这一节点变化对其他节点条件概率的影响，为决策者制定应急方案提供数据支撑。为充分体现贝叶斯网络仿真模拟在流域内降雨、径流丰枯遭遇分析中的优越性，本文以龙川站降雨为枯水和龙川、河源降雨均为枯水两种情境作为后验信息输入到贝叶斯网络中进行仿真研究。

（1）第一次仿真模拟：将龙川站降雨为枯水情况这一后验信息输入到网络结构中得到仿真计算结果见图11。对比图10分析可知，当输入龙川站降雨为枯水这一后验信息后，其他节点各种情况的发生概率均发生较大变化。以三站点降雨及调水风险发生概率为例，只输入先验知识情况下，河源站降雨丰、枯水情况发生概率分别为34%、40%；博罗站降雨丰、枯水发生概率分别为42%、36%；几种调水不利情况风险概率均低于30%。输入龙川站降雨为枯水这一后验知识后，河源站丰、枯水情况发生概率分别为8%、75%；博罗站发生概率为17%、63%；调水不利风险概率基本在50%以上；各站点枯水情况发生概率有大幅度增加。

（2）第二次仿真模拟：贝叶斯网络构造风险管理模型对输入后验知识的节点个数没有限制，可能出现同时输入多节点后验知识的情况。以输入两个节点的后验知识为例，将龙川、河源站降雨均为枯水这一后验知识输入到初始网络结构中，得到第二次仿真模拟结果见图12。对比初始网络模型（图10）分析可知，输入龙川、河源站为枯水情况这一后验知识后，博罗站降雨及三个站点径流枯水情况的发生概率均有大幅度增加，达到65%以上；博罗站径流枯水情况发生概率增幅最大，由原来的28%，增加到现在的82%，受上、中游降雨变化影响最为严重。就调水风险来看，各站点丰枯遭遇调水不利风险概率均超过55%，其中河源博罗站调水不利情况最易发生，风险概率达67%；三站点同枯风险发生概率63%。针对这一情况，当气象部门侦测到流域上、中游降雨同为枯水情况时，相关部门需要制定出具有针对性的应急预案措施。

3 结论

针对流域内不同区间降雨、径流的丰枯遭遇组合状况，以东江流域上游龙川站、中游河源站及下游博罗站年尺度降雨、径流为研究对象，运用Copula函数方法建立了流域内不同站点降雨、径流丰枯遭遇组合的联合分布模型，结合贝叶斯网络理论以实测资料及Copula函数计算概率作为先验知识输入，构建了流域内三个站点降雨、径流丰枯遭遇的风险管理模型。通过模型计算，不考虑后验知识的情况下，各种组合调水不利风险发生概率在25%～30%之间。利用贝叶斯网络的反向推理功能，分别以上游龙川站降雨为枯水和上游龙川站、中游河源站降雨同为枯水两种情境作为后验知识输入网络结构，对流域内可能发生的丰枯组合状态进行仿真模拟计算。计算结果显示，①输入上游龙川站为枯水情况时，各站点丰枯遭遇组合调水不利风险概率达到50%左右，其中三站同枯调水最不利风险发生概率达54%；②输入上游龙川站、中游河源站降雨同枯情况时，站点间降雨、径流丰枯遭遇调水不利风险概率增加到60%左右，其中河源博罗站调水不利风险概率最大，达67%；文中不同情境下降雨、径流丰枯遭遇概率和条件概率的计算结果以及基于贝叶斯模型的仿真模拟结果，为流域内水资源的合理调配及可持续利用提供理论支撑，具有理论和实践意义。

参考文献（References）：

[1] 冉啟香，张翔.多变量水文联合分布方法及Copula函数的应用研究[J].水电能源科学，2010，28（9）：811.（RAN Qixiang，ZHANG Xiang.Review on methods of multivariate hydrological joint distribution and Copula function[J].Water Resources and Power，2010，28（9）：811.（in Chinese））

[2] 费永法.多元随机变量的条件概率计算方法及其在水文中的应用[J].水利学报，1995（8）：6066.（FEI Yongfa.A method for estimating the conditional probability of multirandom variables and its application in hydrology[J].Journal of Hydraulic Engineering，1995（8）：6066.（in Chinese））

[3] Grimaldi S，Serinaldi F.Asymmetric copula in multivariate flood frequency analysis[J].Advances in Water Resources，2006，29（8）：11551167.

[4] 戴昌军，梁忠民.多维联合分布计算方法及其在水文中的应用[J].水利学报，2006，37（2）：160165.（DAI Changjun，LIANG Zhongmin.Computation methods of multivariate joint probability distribution and their applications in hydrology[J].Journal of Hydraulic Engineering，2006，37（2）：160165.（in Chinese））

[5] 郑红星，刘昌明.南水北调东中两线不同水文区降水丰枯遭遇性分析[J].地理学报，2000，55（5）：523532.（ZHENG Hongxing，LIU Changming.Analysis on asynchronismsynchronism of regional precipitation in planned SouthtoNorth Water Transfer Areas[J].Acta Geographica Sinica，2000，55（5）：523532.（in Chinese））

[6] 韩宇平，蒋任飞，阮本清.南水北调中线水源区与受水区丰枯遭遇分析[J].华北水利水电学院学报，2007，28（1）：811.（HAN Yuping，JIANG Renfei，RUAN Benqing.Analysis on wetnessdryness encountering of runoff flow between water source region and receiving water region in the Middle Route of the SouthtoNorth Water Transfer Project[J].Journal of North China Institute of Water Conservancy and Hydroelectric Power，2007，28（1）：811.（in Chinese））

[7] 闫宝伟，郭生练，肖义.南水北调中线水源区与受水区降水丰枯遭遇研究[J].水利学报，2007，38（10）：11781185.（YAN Baowei，GUO Shenglian，XIAO Yi.Synchronousasynchronous encounter probability of richpoor precipitation between water source area and water receiving areas in the Middle Route of SouthtoNorth Water Transfer Project[J].Journal of Hydraulic Engineering，2007，38（10）：11781185.（in Chinese））

[8] 冯平，李新.基于Copula函数的非一致性洪水峰量联合分析[J].水利学报，2013，44（10）：11371147.（FENG Ping，LI Xin.Bivariate frequency analysis of nonstationary flood time series based on Copula methods[J].Journal of Hydraulic Engineering，2013，44（10）：11371147.（in Chinese））

[9] 杨志勇，袁喆，方宏阳，等.基于Copula函数的滦河流域旱涝组合事件概率特征分析[J].水利学报，2013，44（5）：556562.（YANG Zhiyong，YUAN Zhe，FANG Hongyang，et al.Study on the characteristic of multiply events of drought and flood probability in Luanhe River Basin based on Copula[J].Journal of Hydraulic Engineering，2013，44（5）：556562.（in Chinese））

[10] 冯平，牛军宜，张永，等.南水北调西线工程水源区河流与黄河的丰枯遭遇分析[J].水利学报，2010，41（8）：900907.（FENG Ping，ZHU Yijun，ZHANG Yong，et al.Analysis of wetnessdryness encountering probability among water source rivers and the Yellow River in the Western Route of SouthtoNorth Water Transfer Project[J].Journal of Hydraulic Engineering，2010，41（8）：900907.（in Chinese））

[11] 康玲，何小聪.南水北调中线降水丰枯遭遇风险分析[J].水科学进展，2011，22（1）：4450.（ZHANG Ling，HE Xiaocong.Risk analysis of synchronousasynchronous encounter probability of richpoor precipitation in the Middle Route of Southto North Water[J].Advances in Water Science，2011，22（1）：4450.（in Chinese））

[12] 康玲，何小聪，熊其玲.基于贝叶斯网络理论的南水北调中线工程水源区与受水区降水丰枯遭遇风险分析[J].水利学报，2010，41（8）：908913.（KANG Ling，HE Xiaocong，XIONG Qiling.Risk analysis for precipitation richpoor encounter between source area and receiving area of the Middle Route of SouthtoNorth Water Transfer Project Based on Bayesnet theory[J].Journal of Hydraulic Engineering，2010，41（8）：908913.（in Chinese））

[13] O′Connell D R H.Nonparametric Bayesian flood frequency estimation[J].Journal of Hydrology，2005，313（s 12）：7996.

[14] H.Van de Vyver.Bayesian estimation of rainfall intensitydurationfrequency relationships[J].Journal of Hydrology，2015：14511463.

[15] 张雨，宋松柏.基于Archimedean Copula的三维干旱特征变量联合分布研究[J].中国农村水利水电，2011（1）：6568.（ZHANG Yu，SONG Songbai.Research on threedimensional joint distribution of drought characteristics based on Archimedean Copulas[J].China Rural Water and Hydropower，2011（1）：6568.（in Chinese））

[16] 张冬冬，鲁帆，严登华，等.基于Archimedean Copula函数的洪水多要素联合概率分布研究[J].中国农村水利水电，2015（1）：6874.（ZHANG Dongdong，LU Fan，YAN Denghua，et al.Research on multidimensional joint distribution of flood characteristics based on Archimedean Copula[J].China Rural Water and Hydropower，2015（1）：6874.（in Chinese））

[17] NELSON R B.An introduction to copulas[M].New York：Springer，2006.

[18] Murphy K P.A brief introduction to graphical models and bayesian networks[J].Borgelt Net，1998.

[19] 周忠宝，马超群，周经伦，等.基于贝叶斯网络的多态故障树分析方法[J].数学的实践与认识，2008，38（19）：8995.（ZHOU Baozhong，MA Chaoqun，ZHOU Jinglun，et al.Multistate fault tree analysis method based on Bayesian Networks[J].Mathematics in Practice and Theory，2008，38（19）：8995.（in Chinese））

[20] 谢华，罗强，黄介生.基于三维copula函数的多水文区丰枯遭遇分析[J].水科学进展，2012，23（2）：186193.（XIE Hua，LUO Qiang，HUANG Jiesheng.Synchronous asynchronous encounter analysis of multiple hydrologic regions based on 3 D copula function[J].Advances in Water Science，2012，23（2）：186193.（in Chinese））

[21] Singh V P，Zhang L.Bivariate Flood Frequency Analysis Using the Copula Method[J].Journal of Hydrologic Engineering，2006，11（2）：150164.

[22] 张兵利，裴亚辉.贝叶斯网络模型概述[J].电脑与信息技术，2008，16（5）：4142.（ZHANG Bingli，PEI Yahui.Summary of Bayesian Network[J].Computer and Information Technology，2008，16（5）：4142.