多维联合分布计算方法及其在水文中的应用                                                                                 葛维亚       

       我在水文水资源职场工作几十年，对于多维概率分布在水文中的应用，和其他同事一起作过一些具体研究工作，取得一些进展。

        1多维联合分布计算方法

目前，国内外多维联合分布的计算方法研究主要集中于分布为正态分布，对数正态分布,Gumbe1分布,指数分布等，但这些方法大多局限于二维分布，对三维以上联合分布的求解较为困难.本文介绍3种具有一般性的多维联合分布计算方法，即将原始数据变换为正态的Moran方法，并将多维分布转化为一维分布的FEI方法，以及基于经验频率分析的EFM方法。

1.1       正态变换方法(Moran方法)

       由于原始变量的多维联合分布的解析求解比较困难，因此通常将原始变量变换成正态分布变量，继而求解多维正态分布概率。多维正态分布的密度函数如下

,(…,=ex一号(一(—tt)】                               (1)

式中:.,,…,为正态分布变量;为协方差矩阵;p=(.,,…,),=(.,,…,).

采用Moran方法计算多维联合分布的一个关键问题是采用何种正态化变换工具。水利工程领域常用的正态化变换工具主要有Box.Cox变换和多项式正态变换(PNT).研究表明,Box.Cox变换后样本的正态性优于PNT法，因此,本研究中正态化变换工具采用Box—Cox变换。

       对一组原始样本,Box—Cox变换可通过式(2)将其变换成近似正态分布样本，(2)式如下：

         Ya.=,(≠o);     Ya.=In(:),(=0)                         (2)

式中:.是原始样本;Ya是变换后的样本;是变换参数；Ya服从正态分布，N(/1,);/1为均值标准差。

       本研究采用极大似然法结合遗传算法求解变换参数，当原始变量换成正态分布变量后，利用数值积分方法计算多维正态联合分布概率。

1.2将多维转化为一维的FEI方法

利用事件积先后推导出二维及多维随机变量转化为一维随机变量的计算方法，并用来求解丰枯遭遇频率。多维联合分布的计算公式如下：

           P(XI≥I,X2≥2,…,X≥)=P(Z≥.)                     (3)

       此时，随机变量z的系列为:=min(…:+a:,…,+a)。另一种形式的计算公式

          P(Xl≤l,X2≤2,…,Xn≤n)=P(Z≤:)                    (4)

       此时,随机变量z的系列为:=max(2+a2,…,+a)，其中,a=l，z=lp为随机变量频率为P的设计值。

       在实际工作中，随机变量X1,X2…,Xn的分布一般都假设为已知的某种分布，而随机变量z的分布一般与x的分布类型不相同，属于总体分布未知的水文变量，因此,P(Z≥)的计算可采用经验频率或者非参数统计方法，研究表明，多项式正态变换方法(PNT)不需要假设总体分布，是一种统计性能优良的频率分析方法，因此本研究中P(z≥)的计算采用PNT法。一个偏态分布的变量，可以通过三阶的多项式做正态变换，即

           X=a0+a1Z+a2Z+aZ                                 (5)

式中:X是偏态分布的变量；Z是标准正态分布的变量; a.,a.,a,a 是多项式系数，多项式系数a.,a.,a:,a可采用概率权重矩(PWM)估计.

          按照PWM的定义rI卢=l(u)udu                  (6)

J0将式(5)代人式(6)可得：

             卢()=aoC,o+alC,l+a2C,2+a3C.3               (7)

式中:C…:lgrt(z)dz,，(),，庐(z)分别为正态分布的分布函数与概率密度函数。

        当系数a0,a1,a2,a3已知时，就可利用两种分布变换后的分位点对应相等的原则，给定一个正态分布的分位点(分位数)Z，代人式(5)就可推出偏态分布的对应设计值。

1.3基于经验频率分析的EFM方法

当随机变量X.,X:,…,X的系列较长时，且变量个数较少时，可采用经验频率方法(EFM方法)进行计算。以下就二维变量的联合分布概率计算进行说明。首先构造一个X,Y的二维表，X和y的观测值分别按照升序排列，则有样本观测值(,Y)的概率为

              P(f,Y)=P(X=.,Y=Y)=n/(Ⅳ+1)                   (8)

       2统计试验研究

             V(x,),)=P(≤Xi,Y≤),)=∑∑几/(N+1)m=ln=1 (9)

本研究采用蒙特卡罗模拟分析技术，模拟多维(二维)对数正态分布的随机数，比较了Moran方法，FEI方法，EFM方法计算多维对数正态联合分布的效果。

    2.1试验方案设计

  2.1.1总体参数选择对二维联合分布来说，总体参数的组合情况较多,考虑到本文在方法比较的基础上，需进行南水北调东线工程不同水文区年径流丰枯遭遇分析，因此,总体参数选择为以下4种类型:

长江黄河组合型:CI=0.13,Csl=0+26,C=0.24,Cs2=0.96; 长江淮河组合型:Cl=0.13,Csl=0.26,C=0.63,Cs2=1.26; 黄河淮河组合型:C,=0.13,C.=0.96,C=0.63,C:=1.26;淮河沂沭泗组合型:Cyl=0.63,Cs1.26,C=0.87,Cs2=1.74.

       由于均值与结果无关，故统一选择=1000，此外由于各水文区的相关系数较为接近，统计试验中也都采用了p=0.3,模拟组数取=1000.

  2.1.2计算频率选择计算频率分为以下9种情况:

丰丰型:P1=P(≥,y≥y); 丰平型:p2=P(≥,<y<y); 丰枯型:p3=P(≥,y≤y); 平丰型:p4=P(<<,y≥y); 平平型:p5=

P(<<,y<y<y); 平枯型:p6=P(<<,y≤yp); 枯丰型:p7=P(≤,y≥); 枯平型:p8=P(≤,Yl,<y<y);枯枯型:p9=P(≤y≤y)。其中∥,为丰枯划分频率,本文选择∥=37.5%,pk=62.5%.

  2.1.3统计量比较准则统计量(记为0)的比较准则为

,k,厂—————————～

bO=(百1∑/0.一1)×100%:SO=√[∑(一Oo)]IKIO.×100%(10)～i=Ii=I

其中II越小,表示估计的无偏性越好，越小表示有效性越好, 0为二维对数正态分布的理论值.本文比较了9种计算频率的无偏性与有效性.

    2.2统计试验结果分析

  2.2.1无偏性比较无偏性统计试验结果如表1所示.分析表明,Moran方法与EFM方法的无偏性比较接近，偏差相对较小，FEI方法的无偏性，对长江～黄河的径流组合型情况，与另两个方法相当，对其它径流组合情况,相对较差，其原因可能在于组合中某变量的C较大时，采用非参数的PNT法估计频率的效果随之降低, 从而影响了FEI方法的性能。

  2.2.2有效性比较有效性统计试验结果如表2所示.分析表明,就有效性而言,Moran方法,EFM方法, FEI方法的有效性能递减,4种组合的平均有效性指标约为25%,45%和70%。

    2.3统计试验结论由各种计算频率的无偏性及有效性的比较可知：

(1)FEI方法的无偏性能与有效性能均最差，其原因是FEI方法原始的公式只能求解丰丰型或者枯枯型的频率，对其他组合,如:P(≥,y≤y),需要先计算P(≥)及P(≤,y<y),再用P(≥,y≤y肚)=P(≥)一P(≤,y<y)进行求解。由于FEI方法转换后系列服从的分布未知，采用非参数估计的PNT法对其求解也就存在误差,即P(≥)及P(≤,y<y)的计算

存在误差，再通过代数和计算误差累计，从而造成FEI方法在丰枯遭遇计算时成果不稳定，当变量增加时，这种误差将进一步加大。

(2)传统的基于频率分析的EFM方法无偏性较好，但有效性较差, 这是因为EFM方法对资料系列长度较为敏感, 对某一组样本, 计算的频率相差较大, 但模拟的组数增加时, 就趋向于真值,因此总体的无偏性能较好。但由于各组样本的计算结果变化较大, 从而有效性较差。

(3)Moran方法的无偏性及有效性均较优, 因此在进行多维联合分布计算时,建议采用Moran方法进行频率分析。表1为无偏性统计试验结果(单位:%)

       3实例应用

       根据不同水文区内具有代表性的水文站的1956—1997年的年径流资料, 采用Moran方法分析南水北调东线工程不同水文区年径流丰枯遭遇频率, 并与EFM方法进行比较。研究的水区包括:长江中下

游流域、淮河流域、黄河下游流域。其中, 长江中下游流域采用大通站径流资料; 淮河流域含有淮河及沂沭泗两个独立水系(本文以下凡涉及淮河流域均为淮河及沂沭泗两个水系总称), 淮河水系采用蚌埠站

径流资料,沂沭泗水系采用骆马湖入湖径流资料; 黄河下游流域采用花园口站径流资料。

    3.1长江黄河丰枯遭遇分析对研究站点的年径流系列，分别采用皮尔逊Ⅲ型分布(P.Il1)和对数正态分布进行了拟合，结果表明，两种统计分布对年径流系列的拟合效果十分接近。由于对多维P.Il1分布

的计算问题目前尚缺乏解析表达式，无法得到组合概率的理论值，所以本文假定长江、黄河的径流系列服从对数正态分布，采用Moran方法与EFM方法计算遭遇频率大小，并与二维对数正态分布的理论值

进行比较，统计9种丰枯组合的离差绝对值之和与离差平方和，计算结果如表3所示。分析表明,Moran方法更接近于理论值，平均离差绝对值约为1%，优于传统的EFM方法，后者平均离差绝对值约为3%  。

       分析可知,Moran方法更优，这与前面统计试验的结论是一致的。

    3.2有利于调水至黄河的频率分析，本文根据文献[1]对有利于调水频率的定义，对有利于调水至黄

河的频率做如下界定:(1)当长江处于丰水或平水状态时,黄河处于枯水或平水状态，对黄河调水有利；(2)当长江为枯水时,如果淮河水系与沂沭泗同丰或丰平遭遇，同时黄河为非丰，则也是有利于调水至黄河。

                 表3长江黄河丰枯遭遇频率(单位:%)  此表省略。

       研究结果表明，对某一组合而言，变量较少时(如两变量)Moran方法与EFM方法计算结果较为接近，变量较多时(如四变量)，计算结果相差显着。EFM方法由于资料系列较短，组合概率会出现为0的

情况，从理论来说，这是不合理的。因此，对多维联合分布汁算，Moran方法的计算结果更加合理可靠。

     对总的有利于调水的频率而言，两种方法结果接近，Moran方法与EFM方法计算的总的有利于调水至黄河的频率分别约为46%和50%.

       4主要结论

       本研究采用蒙特卡罗模拟分析技术，比较研究了Moran方法，EFM方法和FEI方法3种多维联合分布计算方法的统计性能，并将之应用于南水北调丰枯遭遇分析，得到如下主要结论: 

   (1)统计试验结果表明，在对比的3个方法中，基于正态化变换的Moran方法均具有最优的无偏性和有效性；

  (2)通过南水北调东线工程不同水文区年径流丰枯遭遇组合分析表明，与其他方法比较，Moran方法的结果更合理.因此，对多变量水文事件的概率计算问题，建议采用Moran方法。计算结果初步表明，南水北调东线工程有利于调水至黄河的概率约为46%；

  (3)本文的结论是针对水文变量服从多维对数正态分布的情况，对其他类型的多维分布是否有类似的结论，还有待进一步研究。