本发明属于水文相依结构建模方法，具体涉及一种基于互信息和vinecopula的水文相依结构建模方法。
背景技术：
水文事件都是由多个变量所共同表征的，例如洪水由洪峰、历时、洪水总量三个水文变量表征，干旱由干旱烈度、干旱历时、干旱峰值三个变量表征，因此水文风险分析和管理需要对多变量的相关性进行分析，并对水文变量的相依结构进行建模。vine-copulas已经越来越多地应用于多变量模型中，其中vine-copula的结构选择对多变量相依结构建模起着关键作用。水文系统中存在着大量的不确定性，在对水文相依结构的建模中，应该分析考量变量的不确定性，并实现对其量化。
互信息(mutualinformation,mi)来自信息理论，它不仅是一种度量变量之间相关关系的方法，还能度量知道其中一个水文变量的信息对另一个水文变量的不确定性减少的程度。同时，互信息是负的copula熵，因此互信息也可以看做是变量间相关关系的不确定性的度量。条件互信息(conditionalmutualinformation,cmi)则是以给定变量为条件的互信息，同样可度量变量的条件相关关系以及给定变量的不确定性的减少程度。通常vine-copula的构造需要三个主要步骤：(1)选择vine的结构，即选择哪两个变量成为一对；(2)为步骤(1)中确定的每一对变量确定双变量的copula类型；(3)估计每个copula的相应参数。由于步骤(2)和步骤(3)是取决于vine-copula的结构的，因此选择适当的树结构对多变量的建模，至关重要。
随着变量的数量的增加，可能的vine的构造数量呈指数增长，因此使用穷举法从所有可能的vine-copula来挑选最优结构并不现实。最常用的方法是按照一棵树紧接着一棵树的顺序依次完成树的结构选择，首先先确定第一棵树，依次确定第二棵，第三棵……现在的vine-copula选择方法有基于经验tau,aic准则，自由度，拟合度(goodness-of-fit,gof)，以及p值。基于p值，自由度和aic的选择需要事先对copula进行估计，过程较为繁琐。
目前还未有将互信息和条件互信息应用基于互信息和vinecopula的水文相依结构建模方法，且传统的方法在水文相依结构建模时不能同时考量变量的相关性以及不确定性，而利用基于互信息和copula的水文相依结构建模方法则实现对相关关系和不确定性的同时考量。
技术实现要素：
发明目的：本发明的目的在于解决现有技术中存在的不足，提供一种基于互信息和vinecopula的水文相依结构建模方法。
技术方案：本发明一种基于互信息和vinecopula的水文相依结构建模方法，依次包括如下步骤：
(1)从vinecopula的第一棵树开始，根据公式(1)分别计算两两水文变量间的互信息，结合相关性最强以及不确定性最小的原则，选择使得互信息之和最大的配对方式作为第一棵树的边；
其中，x和y是任意两两水文变量，p(x,y)的x和y联合概率密度，p(x)和p(y)分别是x和y的边缘概率密度；且若总计有n个水文变量，则计算n(n-1)/2个两变量对的互信息；
(2)结合确定的前一棵树，根据公式(2)计算可能配对的变量间的条件互信息，通过比对条件互信息的值，结合水文条件相关性最强以及水文变量信息的不确定性最小的原则，选择使得条件互信息之和最大的配对方式作为vinecopula第二棵树的边；
其中，xs是变量的集合；xi和xj是可能配对的水文变量；
(3)重复步骤(2)，直到所有的vinecopula树的边都确定；
(4)将p-iii型分布，gamma分布，极值型分布，对数正态分布等分布作为边缘分布的候选分布线型，利用单变量水文数据，应用极大似然法等参数估计方法分别对候选分布线型进行参数估计，利用aic准则选择最优分布，并利用kolmogorov-smirnov检验对选择的单一变量分布进行拟合优度检验；
(5)根据步骤(1)-(3)确定的vinecopula树的构造，从第一棵树开始，将常用的椭圆copula(例如高斯和tcopula)、archimedeancopula(例如clayton,gumbel,frank,和joe)以及他们旋转90°，180°，270°后的类型作为每个边的候选copula类型，利用极大似然估计等参数估计方法估计相应的参数，利用aic准则选择每个边的copula的类型，并对选择的最优copula做cm拟合优度检验；
(6)利用h函数(即公式(3))，计算变量的条件边缘分布，估计候选copula参数，利用aic准则选择copula类型，并，最后做cm检验；
其中，xi,xj,是变量，xs是条件变量集合，xs-j表示不包含xj的条件变量集合h(fi|s-j(xi|xs-j),fj|s-j(xj|xs-j))是h函数，ci,j；s-j(fi|s-j(xi|xs-j),fj|s-j(xj|xs-j))是copula函数，fi|s-j(xi|xs-j)，fj|s-j(xj|xs-j)是边缘分布函数；表示对c做f的偏导数；
(7)重复步骤(6)，直到所有边的copula都确定；
(8)根据vinecopula定义，连接所有树以及树的边即可完成对水文相依结构的建模。
有益效果：本发明同时实现了变量间的相关关系的度量以及水文变量的不确定性的考量，与现有技术相比，本发明具有以下优点：
(1)本发明利用互信息以及条件互信息作为变量间相关关系的度量，并选择最大化互信息和条件互信息的变量配对方式作为每层树的边。这个选择方法的每一个步骤都是基于原始数据的，在选择过程中不需要对数据进行转变，是独立于模型的，因此在选择过程中不需要事先对copula进行选择和估计，使得建模过程更加简便快捷。
(2)本发明利用互信息和条件互信息，能度量水文变量间的不确定性，并且根据互信息和条件互信息的特性，能够实现水文变量不确定性的量化。
综上所述，本发明利用互信息和条件互信息，既能度量多变量间的相关关系，又能同时考量变量间以及变量间相关关系的不确定性，具有合理性和有效性。
附图说明
图1为本发明的流程示意图；
图2为实施例中vinecopula的结构示意图；
图3为实施例中单一变量拟合分布的概率示意图；
图4为实施例中相应的copula函数以及估计的参数示意图；
图5为实施例中多站点径流的相依结构图；
图6为实施例中联合分布示意图；
图7为实施例中的散点对比图。
其中，图3(a)是a站点的拟合的理论概率分布与经验概率分布的点图(probability-probabilityplot,pp图)；图3(b)是b站点的pp图；图3(c)是c站点的pp图；图3(d)是d站点的pp图；图7(a)是站点a和b的实测数据和模拟数据的散点图；图7(b)是站点a和c的实测数据和模拟数据的散点图；图7(c)是站点a和d的实测数据和模拟数据的散点图；图7(d)是站点b和c的实测数据和模拟数据的散点图；图7(e)是站点b和d的实测数据和模拟数据的散点图；图7(f)是站点c和d的实测数据和模拟数据的散点图。
具体实施方式
下面对本发明技术方案进行详细说明，但是本发明的保护范围不局限于所述实施例。
便于理解本发明，做以下说明：
r-vine:
定义1：树是每两个节点被不同边连接起来的图，根据树的基础定义，则v＝(t1,t2,…,tn)在满足以下条件时就是一个包含了n个变量的r-vine:
(1)t1是具有n1＝(1,2,…,n)个节点和边为e1的一棵树；
(2)对于i＝2,…n-1,ti是具有节点为ni＝ei-1和边为ei的一棵树；
0)对于i＝2,…n-1,{a,b}∈ei，其中a＝{a1,a2},b＝{b1,b2},并且a和b满足#(a∩b)＝1，#表示集合的势。
copula:
定义2：h是具有边缘分布为fi的联合分布函数，若xi为随机变量，那么存在copula函数c，使得：h(x1,x2,...,xn)＝c(f1(x1),f2(x2),...,fn(xn))。
vine-copula：
定义3：一个具有相同边缘分布且联合分布服从[0,1]n多元分布的copula称之为vine-copula。在(f,v,b)中，如果f＝(f1,f2,…,fn)是一系列连续且可逆的分布函数，v是一个n维的r-vine结构，b＝{be|i＝1,…n-1；e∈ei}是一个copula函数集(其中be为二元copula函数)，则(f,v,b)为r-vinecopula。
如图1所示，本发明一种基于互信息和vinecopula的水文相依结构建模方法，依次包括如下步骤：
(1)从vinecopula的第一棵树开始，根据公式(1)分别计算两两水文变量间的互信息，结合相关性最强以及不确定性最小的原则，选择使得互信息之和最大的配对方式作为第一棵树的边；
其中，x和y是两两水文变量，p(x,y)的x和y联合概率密度，p(x)和p(y)分别是x和y的边缘概率密度；且若总计有n个水文变量，则计算n(n-1)/2个两变量对的互信息；
(2)结合确定的前一棵树，计算可能配对的变量间的条件互信息，通过比对条件互信息的值，结合水文条件相关性最强以及水文变量信息的不确定性最小的原则，选择使得条件互信息之和最大的配对方式作为vinecopula第二棵树的边；
其中，xs是变量的集合；xi、xj是可能配对的水文变量；
(3)重复步骤(2)，直到所有的vinecopula树的边都确定；
(4)将p-iii型分布，gamma分布，极值型分布，对数正态分布等分布作为边缘分布的候选分布线性，利用单变量水文数据，应用极大似然法等参数估计方法分别对候选分布线性进行参数估计，利用aic准则选择最优分布，并利用kolmogorov-smirnov检验对选择的单变量分布进行拟合优度检验；
(5)根据步骤(1)-(3)确定的vinecopula树的构造，从第一棵树开始，将常用的椭圆copula(高斯和tcopula),archimedeancopula(clayton,gumbel,frank,和joe)以及他们旋转90°，180°，270°后的类型作为每个边的候选copula类型，利用极大似然估计等参数估计方法估计相应的参数，利用aic准则选择每个边的copula的类型，并对选择的最优copula做cm拟合优度检验；
(6)利用h函数(公式(3))，计算变量的条件边缘分布，估计候选copula参数，利用aic准则选择copula类型，并，最后做cm检验；
其中，xi,xj,是变量，xs是条件变量集合，xs-j表示不包含xj的条件变量集合h(fi|s-j(xi|xs-j),fj|s-j(xj|xs-j))是h函数，ci,j；s-j(fi|s-j(xi|xs-j),fj|s-j(xj|xs-j))是copula函数，fi|s-j(xi|xs-j)，fj|s-j(xj|xs-j)是边缘分布函数；表示对c做f的偏导数；
(7)重复步骤(6)，直到所有边的copula都确定；
(8)根据vinecopula定义，连接所有树以及树的边即可完成对水文相依结构的建模。
实施例：
本实施例以colorado流域的多站点径流相关结构建模作为实际应用
以colorada流域四个站点的径流数据为例，应用基于互信息和vinecopula的水文相依结构建模方法对四个站点的径流相关结构进行建模。
(1)实测数据
数据来源于colorada流域从1906年-2015年四个站点的年径流数据，四个站点分别为：maybell附近的yampariver(usgs编码：9251000)，lily附近的littlesnakeriver(09260000)，randlett附近的duchesneriver(09302000)和greenriveratgreenriver(09315000)，分别表示为站点a，b，c，d。
关于数据库的更多细节可以从以下网站获取：
(http://www.usbr.gov/lc/region/g4000/naturalflow/previous.html)。
(2)第一棵树的互信息
根据上表，可知变量对(a,d),(b,d),(c,d)之间的相关性最强，同时，若已知d变量的信息，则a,b,c变量的不确定性由原来的2.30分别可减少0.98,1.06,1.08，不确定性减少程度最高，因此，结合相关性最强和不确定最小的原则，选择(a,d),(b,d),(c,d)变量对作为vinecopula第一棵树的边。
(3)第二棵树的条件互信息
根据上表，可知在给定条件变量为d的情况下变量对(a,b),(b,c)的条件相关性最强，同时，若给定d变量条件下，已知b变量的信息，则a,b变量的不确定性由原来的1.32,1.22分别可减少0.64,0.62，不确定性减少程度最高，因此，结合相关性最强和不确定最小的原则，选择(a,b；d),(b,c；d),变量对作为vinecopula第二棵树的边。
(4)选择的r-vine结构，最终选择的vinecopula的结构如图2所示。
(5)a,b,c,d的单一变量拟合分布的概率，如图3所示；图3(a)是a站点的拟合的理论概率分布与经验概率分布的点图(pp图)，拟合分布为参数为(12.15,1.14,-1.46)的gamma分布；图3(b)是b站点的pp图，拟合分布为参数(7.90,0.62,-0.22)的gamma分布；图3(c)是c站点的pp图，拟合分布为参数(8.02,1.04,-0.52)的gamma分布；图3(d)是d站点的pp图，拟合分布为参数(15.41,4.17,-10.74)的gamma分布；
(6)相应的copula函数以及估计的参数如图4所示，其中，*表示该copula函数了旋转180°，#表示旋转了90°，※表示旋转了270°；
(7)连接所有vinecopula树和边，得到多站点径流的相依结构如图5所示；
(8)径流的联合分布的pp图，如图6所示；
(9)实测数据(observed)和模拟数据(simulated)的散点比较如图7所示，图7(a)是站点a和b的实测数据和模拟数据的散点图；图7(b)是站点a和c的实测数据和模拟数据的散点图；图7(c)是站点a和d的实测数据和模拟数据的散点图；图7(d)是站点b和c的实测数据和模拟数据的散点图；图7(e)是站点b和d的实测数据和模拟数据的散点图；图7(f)是站点c和d的实测数据和模拟数据的散点图；
通过上述实施例看出，本发明利用互信息和条件互信息来对多站点径流相依结构进行建模，能够很好地度量变量间的相关关系，同时考量了系统中不确定性，更适合于水文相依结构的建模，且更为简单和直观。