CN102542169A - 一种在水文频率计算过程中进行线型选择的方法 - Google Patents

Abstract

本发明公开了一种在水文频率计算过程中进行线型选择的方法，首先依据已知的地区先验信息，分别选择合理的参数先验分布类型、参数初始样本和似然函数，应用am-mcmc方法进行参数后验分布采样，得到各线型对应的参数后验分布采样结果；利用pome对各参数的后验分布采样结果的概率分布进行分析和定量描述，得到不同线型中各参数的后验分布表达式，应用参数区间近似求和的方法代替线型边缘分布积分过程，根据下式求解水文线型m_j 相对于线型m_i 的贝叶斯因子b_ji ；

在求得贝叶斯因子b_ji 的基础上，即可进行水文线型选择和综合分析。本发明合理地分析和描述参数不确定性，使得分析计算结果得到明显改善。

Description

一种在水文频率计算过程中进行线型选择的方法

 

技术领域

本发明涉及水文学及水资源领域，具体是一种在水文频率计算过程中进行线型选择的方法。

背景技术

水文频率计算过程中涉及两个关键问题：线型选择和参数估计。由于水文循环过程中许多随机和复杂因素的共同作用和影响（刘光文. 水文分析与计算. 北京：水利电力出版社, 1989；芮孝芳. 流域水文模型研究中的若干问题. 水科学进展, 1997, 8(1): 94-98；郭生练. 设计洪水研究进展与评价. 北京: 中国水利水电出版社, 2005；Sang YF, Wang D, Wu JC, et al. The relation between periods’ identification and noises in hydrologic series data. Journal of Hydrology, 2009, 368(1-4): 165-177），各地水文特性不同，水文成因各异，当掌握的地区信息较少时线型选择较为困难，而采用统一水文线型(如P-III线型)并不能很好地拟合不同地区的水文特性。此外，由于实际水文资料有限（特别当序列长度较短时），参数估计结果也存在一定误差。因此，线型选择和参数估计结果均存在不确定性，使水文设计结果存在较大不确定性和风险。

基于贝叶斯理论的水文线型选择和综合方法是解决水文频率分析不确定性的有效方法（Tung WH, Mays LW. Risk models for flood levee design. Water Resources Research, 1981, 17(4): 833-841；刘攀, 郭生练, 田向荣, 张洪刚. 基于贝叶斯理论的水文频率线型选择与综合. 武汉大学学报(工学版), 2005, 38(5): 36-40），其通过将已知地区信息（包括先验信息和专家经验等）与水文模型耦合，可实现最大程度地描述研究对象（线型、参数等）的不确定性，因此可以有效提高水文设计结果的精度和可靠性。该方法具有完善的理论依据，但求解贝叶斯因子时存在两个难点问题：参数先验分布确定和线型边缘分布数值积分。实际中常用的一些简化和近似方法均存在一定的局限和缺陷。例如BIC准则，由于极大似然法的缺陷使其具有一定局限性；且BIC准则未能分析和描述参数的不确定性。为从根本上提高贝叶斯因子求解结果及水文线型选择与综合结果的准确性，需要深入探讨关于参数先验分布确定和线型边缘分布数值积分的有效方法和手段。

贝叶斯线型选择和线型综合理论具体如下：

设备选水文线型M _i (i=1,2,.…,q)的概率密度函数为f _i (X︱θ _i )，X代表样本序列，θ _i 表示水文线型M _i 中的参数。则线型M _j 相对于线型M _i 的贝叶斯因子B _ji 求解如下：

                                                 

                            (1)

式中，p(X︱M _i )表示水文线型M _i 的边缘分布，π(θ _i )表示参数θ _i 的先验分布。在求解出贝叶斯因子B _ji 之后，一般根据表1对水文线型的优劣进行评判。

为了进一步减小线型选择结果的不确定性，可对各线型结果进行综合。设样本序列X属于线型M _i 的先验概率为P(M _i )，则根据贝叶斯方法可求得对应的后验概率为P(M _i ︱X)。对于各线型的先验概率P(M _i )，在无信息的情况下一般认为各线型的先验概率相等，即P(M _i )=1/q。

                   (2)

设F(x ₀ )表示水文变量大于x ₀ 的概率，F(x ₀ ︱M _i )表示应用线型M _i 时水文变量大于x ₀ 的概率，则依据样本序列X的水文线型综合(加权平均)分析结果如下：

                          (3)

      估计方差： 

   (4)

可以看出上述贝叶斯线型选择和综合方法具有很好的理论依据。但在利用式(1)求解贝叶斯因子时，参数先验分布π(θ _i )难以确定和选取，同时线型边缘分布积分也困难。

BIC准则(式(5))是一种常用的近似处理方法（Kass R E, Raftery A E. Bayes Factors. Journal of the American Statistical Association, 1995, 90(430): 773-795），但存在一些缺陷：（1）由于其实质是通过极大似然法进行参数寻优后求解贝叶斯因子，而当P-III线型中参数C _s >2时，极大似然法无解（金光炎. 极大似然法估计τ分布参数的注记. 水资源研究, 2008, 29(2): 14-17），因此该方法存在一定局限性；（2）实际中参数同样存在不确定性，而BIC准则的实质和核心是参数寻优，并未对参数的不确定性进行分析和描述。

                              (5)

式(5)中，k _i 和

分别为线型M _i 中参数θ _i 的个数和极大似然估计值，N为样本序列的长度。

发明内容

本发明所要解决的技术问题是提供一种基于贝叶斯理论和最大熵原理(principle of maximum entropy, POME)相结合，综合确定参数后验分布表达式，最终实现在水文频率计算过程中的线型选择。

本发明所述的一种在水文频率计算过程中进行线型选择的方法，其包括以下步骤：

1）对于待分析的各水文线型，首先依据已知的地区先验信息，分别选择合理的参数先验分布类型、参数初始样本和似然函数，应用AM-MCMC方法进行参数后验分布采样，得到各线型对应的参数后验分布采样结果；

2）利用POME对各参数的后验分布采样结果的概率分布进行分析和定量描述，得到不同线型中各参数的后验分布表达式；

3）应用参数区间近似求和的方法代替线型边缘分布积分过程，根据式(14)求解水文线型M _j 相对于线型M _i 的贝叶斯因子B _ji ；

             (14)

其中，k _i 表示线型f _i (X︱θ _i )中的参数个数，A _i 表示线型f _i (X︱θ _i )对应的参数贝叶斯采样结果，θ _i 表示A _i 中的第i组参数，

表示由POME确定的第l个参数的后验分布，θ _il 表示第i组参数中的第l个参数值；

4）在求得贝叶斯因子B _ji 的基础上，即可进行水文线型选择和综合分析。

步骤1）的具体过程为：

（1）对实测水文序列进行可靠性、一致性和代表性审查，并选定合理的水文线型f(x/θ)加以描述；

（2）以线型参数随序列长度变化的分析结果作为依据，并结合实际水文特性确定各参数的初始取值范围；然后选定各参数的先验分布类型π (θ _i ), i=1,2…,n，n为参数个数；

（3）AM-MCMC采样方法初始化，i＝0；

（4）依据所研究序列进行线型参数估计，并将其结果作为采样的初始样本θ ₀ ；

（5）利用式(7)求解该参数样本的协方差矩阵，并计算生成新的参数样本θ*~N(θ _i ，C _i )。式中，i ₀ 为初始采样次数；C ₀ 为固定初始协方差矩阵；S _d 为一比例因子；ε为一较小正数，保证C _i 不为奇异矩阵； I _d 为单位矩阵；

 

                          (7)

（6）根据贝叶斯条件概率分布公式(8)，计算新样本的接受概率α；

               (8)

（7）生成一随机数u~U(0,1)，与α比较，若α>u接受θ _i+1 =θ*，否则θ _i+1 ＝θ _i ；

（8）i=i+1,重复步骤(5)－(7)，直至生成足够的样本数为止；

（9）分析各参数样本特征值随样本数增大时的变化规律，判断采样结果是否收敛，若无，则需要增加采样次数直至结果达到稳定；

（10）根据参数采样结果，分析各参数后验分布的特征值，包括最大值(max)、最小值(min)、

、C_υ、C_s和众数，并通过确定各置信区间的值域，定量分析各参数后验分布的不确定性；

（11）根据参数贝叶斯采样结果，求解指定频率对应水文设计值的概率分布，通过定量分析水文设计值概率分布的特征值和各置信区间，掌握水文设计结果的概率分布特性并估计风险性。

本发明基于贝叶斯理论和POME综合确定参数后验分布表达式，由于该方法的分析和计算过程能够利用已知的地区先验信息并通过贝叶斯方法与水文线型耦合，因此可实现最大程度地描述线型参数的不确定性；且由于基于POME，因此分析结果合理和客观，可避免人为选择参数先验分布时的主观性。经实例分析和Monte-Carlo统计试验，验证了该方法的有效性和实用性，并同时分析了影响水文线型选择和综合结果可靠性的主要因素。

附图说明

    图1 五强溪坝址洪峰流量频率分析结果。

具体实施方式

1、参数后验分布的确定

由于水文线型参数的先验分布未知，且不同线型参数的先验分布类型是否相同也未知，实际中当掌握的地区先验信息较少时，欲准确确定线型参数的先验分布则更加困难，且人为选取参数先验分布时亦存在较大的主观性。若首先依据已知的地区先验信息合理选择参数的先验分布后推求后验分布，然后再进行水文线型选择和综合，可提高水文频率分析结果的精度和可靠性。

具体分析思路和过程如下。

1.1 参数贝叶斯后验分布采样  

基于贝叶斯理论的参数后验分布推导公式如式(6)。在样本序列X给定的条件下，π(θ)表示参数θ的先验分布，实际中需根据已知的地区信息进行选择和确定；p(X/θ)称为似然函数，其综合了模型信息和样本信息，针对水文频率计算的特点，本发明选择似然函数为p(X/θ)=∏f(x _i /θ)；h(X/θ)则将模型信息、样本信息和参数先验分布有机联合起来；m(X)是一比例常数；π(θ/X)为求得参数的后验分布，也即条件概率分布。

         (6)

由于式(6)右端分母无法准确积分，难以获得参数后验分布的解析解，本发明采用AM-MCMC (Adaptive Metropolis -Markov Chain Monte Carlo)进行待估参数的贝叶斯后验分布采样（Haario H, Saksman E, Tamminem J. An adaptive Metropolis algorithm. Bernoulli, 2001, 7 (2): 223-242；邢贞相,芮孝芳,崔海燕,余美. 基于AM-MCMC算法的贝叶斯概率洪水预报模型. 水利学报, 2007, 38(12): 1500-1506），以获得参数的后验分布采样结果。具体采样过程中需要根据地区先验信息对参数初始取值范围确定、参数先验分布选择、参数初始样本确定、似然函数选择等关键问题进行分析，详见文献（Sang YF, Wang D. Risk estimation for design floods. Theory and Practice of Risk Analysis and Crisis Response. Paris: Atlantis Press, 2008, 363-369；桑燕芳, 王栋, 吴吉春. 基于贝叶斯理论的水文线型参数不确定性分析. 水电能源科学, 2007, 27(6), 15-19）。

AM-MCMC主要分析和计算步骤如下：

（1）对实测水文序列进行“三性”（可靠性、一致性和代表性）审查，并选定合理的水文线型f(x/θ)加以描述；

（2）以线型参数随序列长度变化的分析结果作为依据，并结合实际水文特性确定各参数的初始取值范围；然后选定各参数的先验分布类型π (θ _i ), i=1,2…,n，n为参数个数；

（3）AM-MCMC采样方法初始化，i＝0；

（4）依据所研究序列进行线型参数估计，并将其结果作为采样的初始样本θ ₀ ；

（5）利用式(7)求解该参数样本的协方差矩阵，并计算生成新的参数样本θ*~N(θ _i ，C _i )。式中，i ₀ 为初始采样次数；C ₀ 为固定初始协方差矩阵；S _d 为一比例因子；ε为一较小正数，保证C _i 不为奇异矩阵； I _d 为单位矩阵；

 

                          (7)

（6）根据贝叶斯条件概率分布公式(8)，计算新样本的接受概率α；

               (8)

（7）生成一随机数u~U(0,1)，与α比较。若α>u接受θ _i+1 =θ*，否则θ _i+1 ＝θ _i ；

（8）i=i+1,重复步骤(5)－(7)，直至生成足够的样本数为止；

（9）分析各参数样本特征值随样本数增大时的变化规律，判断采样结果是否收敛。若无，则需要增加采样次数直至结果达到稳定；

（10）根据参数采样结果，分析各参数后验分布的特征值(包括最大值(max)、最小值(min)、

、C_υ、C_s和众数)，并通过确定各置信区间的值域，定量分析各参数后验分布的不确定性；

（11）根据参数贝叶斯采样结果，求解指定频率对应水文设计值的概率分布(即“概率值”的概率分布)。通过定量分析水文设计值概率分布的特征值和各置信区间，掌握水文设计结果的概率分布特性并估计风险性。

1.2 参数后验分布表达式的确定  

POME以整体的观念去描述随机变量含有的不确定性，熵最大意味着因为数据不足而作的主观假定最小，因而所得的解是最客观和最似然的（Singh VP. Entropy-based parameter estimation in Hydrology. Kluwer Academic Publishers (Boston/London), 1998）。在利用AM-MCMC方法获得参数后验分布采样样本的基础上，需要进一步确定样本对应的概率密度，以便进行水文线型选择和综合。当直接对AM-MCMC采样结果进行区间划分后求解概率密度时，无法准确确定逐个样本的概率密度，且计算结果会受到区间划分大小、采样结果误差等因素的影响。为更客观地分析参数不确定性，本发明应用POME对参数采样结果的概率分布进行分析和定量描述，以得到各参数后验分布表达式。由于分析过程中仅依据采样样本本身，因此结果合理可靠。

对于参数θ的后验分布采样结果，在满足某些给定约束条件的情况下，根据最大熵原理确定其后验分布，可归为求解下述数学规划问题：

                            (9)

该规划问题为一泛函条件极值问题，通过引入Lagrange乘子λ ₀ ，λ ₁ ，…，λ _m ，应用变分法对其求解得到式(10)。对于一般工程问题，n取3可以满足精度要求（Siddall JN. 工程概率设计：原理和应用. 北京: 北京科学出版社, 1989）。

                         (10)

可以看出，只要确定各Lagrange乘子的值便可完全确定

的解析表达式。根据式(9)中的各约束条件进行数学推导，可知各Lagrange乘子应满足式(11)：

       (11)

式(11)是关于λ ₀ ，λ ₁ ，…，λ _m 的方程组。为便于数值求解，对其作变形可得：

           (12)

然后建立参数优化问题：

                              (13)

    对式(13)的参数优化问题进行求解并结合式(11)得到λ ₀ ，λ ₁ ，…，λ _m 之后，代入式(10)即可求得参数后验分布

的表达式，并用其替代式(1)中的参数先验分布。

 

2、线型边缘分布近似求和方法

在求得参数后验分布表达式之后，为解决线型边缘分布积分的困难，本发明采用参数后验分布采样结果中逐个样本近似求和的方法来代替积分过程。在近似求和过程中做如下两点近似处理：（1）假定水文线型各参数之间相互独立；（2）经分析认为参数取值范围一般取[M ₀ -5σ, M ₀ +5σ]时近似求和结果与积分结果相差无几（刁艳芳, 王本德, 刘冀. 基于最大熵原理方法的洪水预报误差分布研究. 水利学报, 2007, 38(5): 591-595），其中M ₀ 表示参数后验分布众数，σ表示参数后验分布标准差。

此外，在具体的近似求和过程中应遵循下述两点原则：（1）在近似求和之前，需要判断各线型对应前述的参数后验分布采样结果中由最大值和最小值构成的区间是否能够包含[M ₀ -5σ, M ₀ +5σ]。若能，则可利用该采样结果得到的逐个参数样本进行近似求和；若无，则需要增加采样次数，直至采样结果的取值范围能够包含[M ₀ -5σ, M ₀ +5σ]；（2）应保证在取不同线型进行分析时，各线型对应的参数后验分布采样次数（即采样样本总数）相等。

3、贝叶斯因子求解新方法

根据上述关于确定参数后验分布和线型边缘分布近似求和这两个关键问题的讨论和分析结果，提出贝叶斯因子求解新方法，具体求解步骤如下：

（1）对于待分析的各水文线型，首先依据已知的地区先验信息，分别选择合理的参数先验分布类型、参数初始样本和似然函数，应用AM-MCMC方法进行参数后验分布采样，最终得到各线型对应的参数后验分布采样结果；

（2）利用POME对各参数的后验分布采样结果的概率分布进行分析和定量描述，得到不同线型中各参数的后验分布表达式；

（3）应用参数区间近似求和的方法代替线型边缘分布积分过程，最终根据式(14)求解水文线型M _j 相对于线型M _i 的贝叶斯因子B _ji 。

            (14)

其中，k _i 表示线型f _i (X︱θ _i )中的参数个数，A _i 表示线型f _i (X︱θ _i )对应的参数贝叶斯采样结果，θ _i 表示A _i 中的第i组参数，

表示由POME确定的第l个参数的后验分布，θ _il 表示第i组参数中的第l个参数值。在求得贝叶斯因子B _ji 的基础上，即可进行水文线型选择和综合分析。

算例分析及结果讨论

1、实例分析

以沅水五强溪坝址的50年实测年最大洪水流量序列为例进行分析。将水文频率计算中常见的皮尔逊-III型(P-III)、对数皮尔逊-III型(LP-III)、三参数对数正态分布(LN3)和广义极值分布(GEV)共四种线型作为待分析的水文线型，分别采用BIC准则和本发明求解贝叶斯因子，然后进行线型综合分析。分析过程中采用Weibell公式计算经验频率，选择线性矩法估计参数。为便于结果对比分析，选用χ ² 拟合度和均方误差(RMSE)两个定量指标对结果优劣进行判别。相关分析结果分别见表2和图1。

 

该实例分析结果显示：（1）P-III线型对应的各设计值求解结果最小，BIC准则对应的各设计值求解结果最大，而本发明对应的结果介于上述二者之间；（2）在考虑历史洪水时该序列对应0.5%的洪峰流量为43400m³/s（刘光文. 水文分析与计算. 北京：水利电力出版社, 1989），可以看出BIC准则对应的设计值44880 m³/s结果偏大，该结果可能会过高估计设计洪水值；（3）本发明对应的贝叶斯综合线型的χ ² 拟合度和均方误差均为最小，分别为0.103和713.87，表明该新方法的线型综合结果对经验频率点据的拟合度最好，因此其适线结果更为合理。

2、Monte-Carlo统计试验

2.1 试验方案设计 

考虑到实测水文序列服从的真实水文线型未知，且利用各种经验频率分布进行计算时结果均存在不同程度的偏差等问题，为进一步验证本发明所提贝叶斯因子求解新方法的实用性和有效性，同时也为便于比较不同贝叶斯因子求解方法之间的优劣，此处采用Monte-Carlo统计试验作进一步分析和验证。

同样选用P-III、LP-III、LN3和GEV共四种水文线型进行Monte-Carlo统计试验。为分析样本序列长度、参数取值大小等因素对贝叶斯因子求解和线型后验概率计算结果准确性的影响，生成了不同长度和不同参数值的随机序列分别进行统计分析。其中，序列长度和参数取值分别设定3种不同情况，共计4类各9种统计试验方案(表3-表6)。对每种试验方案下生成的服从某指定水文线型分布的随机序列，应用BIC准则和本发明分别计算各线型之间的贝叶斯因子B _ji ，并利用式(2)计算各线型的后验概率值。根据结果的收敛性确定每种方案的统计试验次数为200，然后取平均值。最后对各种方案的统计试验结果进行对比分析，以验证新方法的实用性和有效性。

2.2 统计试验结果

对各种方案的统计试验结果（表3-表6）进行综合分析和整理，可得到以下结论：

（1）样本序列长度对各线型后验概率求解结果的精度有很大影响。当序列长度较短时，两种方法对各线型后验概率计算结果的偏差均相对较大。如表3和表4中序列长度为10时，BIC准则结果中真实线型的后验概率均不超过0.5；随着序列长度的增大，两种方法均能更好地识别出真实线型。例如表3三参数为(100，0.5，1.0)的统计方案中，当序列长度为10、50和100时，BIC结果中P-III线型的后验概率分别为0.2147、0.7094和0.8358。

（2）参数取值大小对BIC准则的分析计算结果有很大影响。以P-III线型的统计试验结果为例，当序列的偏态系数C _s 取较小值时（例如C _s =0.2），BIC准则难以识别出真实的线型。同样当偏态系数较大值时（例如C _s =3），BIC准则也不能识别出真实的线型。分析原因认为是由参数估计结果存在较大误差使得贝叶斯因子求解结果不准确造成。

（3）各种线型在被识别过程中均受到样本序列长度、参数取值大小等因素的影响，但由于不同水文线型之间存在本质差异，因此其受这些因素影响程度的大小不同。本统计试验结果表明，在序列较短时更容易识别LN3和GEV线型，而欲准确识别P-III和LP-III线型则需要更长的序列。

（4）对于各统计试验方案，本发明提出的贝叶斯因子求解新方法的分析结果较BIC准则分析结果均有明显改善。同一试验方案中，本发明所求得真实线型的后验概率值均明显高于对应BIC准则的结果。特别在序列长度偏短、参数取较大值或较小值等不利情况下，本发明分析结果的准确性更加明显高于BIC准则。例如表3和表4中当序列长度为10时，本发明结果中真实线型的后验概率均超过0.6。由于新方法合理地分析和描述了参数不确定性，因此可有效提高各线型后验概率计算结果的准确性和可靠性。

表3   P-III型模拟序列统计试验对应真实线型后验概率值结果

注：估计P-III线型中参数值时采用线型矩法

表4   LP-III型模拟序列统计试验对应真实线型后验概率值结果

注：表3和表4中括弧里三个参数分别指P-III线型和LP-III线型的均值、变差系数（C_v）和偏态系数（C_s）。

表5   LN3型模拟序列统计试验对应真实线型后验概率值结果

注：括弧中的三个参数分别指LN3线型的均值、变差系数（C_v）和位置参数。

表6   GEV型模拟序列统计试验对应真实线型后验概率值结果

注：括弧中的两个参数分别指GEV线型的位置参数和刻度参数，生成每个模拟序列时，GEV线型中的参数K随机产生。

3、两种方法对比分析

结合文献（Sang YF, Wang D. Risk estimation for design floods. Theory and Practice of Risk Analysis and Crisis Response. Paris: Atlantis Press, 2008, 363-369；桑燕芳,王栋,吴吉春. 基于贝叶斯理论的水文线型参数不确定性分析. 水电能源科学, 2007, 27(6), 15-19）中关于水文线型参数不确定性的研究结果进行综合分析，可以认识到本发明与BIC准则之间的本质不同在于：

（1）根据式(5)可知， BIC准则的实质是通过寻求一组最优参数值对贝叶斯因子进行求解，未考虑参数的不确定性。因此当参数值求解存在较大误差时，相应的贝叶斯因子求解结果也同样会存在较大误差。特别是当序列较短、参数取较大值或较小值等情况下，参数估计结果往往会存在较大误差，相应BIC准则的线型后验概率计算结果也会存在较大偏差。

（2）本发明贝叶斯因子求解新方法较BIC准则的分析结果有明显改善，这是因为该法通过贝叶斯方法和POME合理地确定参数的后验分布表达式，实现了对参数不确定性的定量描述，可避免BIC准则单一参数最优值造成的不确定性和计算误差，因此可以明显改善线型后验概率计算结果的精度。对各模拟试验方案，本发明的分析结果均明显优于BIC准则；

（3）研究结果显示，参数不确定性分析结果会受到参数灵敏度、样本序列长度和样本序列统计特性等多方面因素的影响。应用本发明求解贝叶斯因子时同样会受到这些因素的影响。特别当序列偏短，序列的离散性和偏态性明显、参数的灵敏度较小时，相应参数后验分布的不确定性会较大，因此贝叶斯因子计算结果的误差也会有所增大，但整体上该新方法的分析结果仍能够很好地满足实际应用时所需的精度要求。

水文频率计算对许多实际工作具有重要的意义。通过水文线型选择和综合方法可以有效地提高水文设计结果的精度和可靠性。本发明在研究和解决参数后验分布确定和线型边缘分布近似求和这两个难点问题的基础上，提出了贝叶斯因子求解新方法。通过实例分析和Monte-Cralo统计试验，并与BIC准则进行对比，验证了新方法的实用性和有效性。但在实际应用中，如何更好地利用已知的地区先验信息，以提高参数后验分布分析结果的准确性和合理性需要进一步研究。此外，对于短序列情况下如何进一步提高水文线型选择和综合结果的精度也仍需作深入研究。

Claims (2)

1.一种在水文频率计算过程中进行线型选择的方法，其特征在于包括以下步骤：

1）对于待分析的各水文线型，首先依据已知的地区先验信息，分别选择合理的参数先验分布类型、参数初始样本和似然函数，应用AM-MCMC方法进行参数后验分布采样，得到各线型对应的参数后验分布采样结果；

2）利用POME对各参数的后验分布采样结果的概率分布进行分析和定量描述，得到不同线型中各参数的后验分布表达式；

3）应用参数区间近似求和的方法代替线型边缘分布积分过程，根据下式求解水文线型M _j 相对于线型M _i 的贝叶斯因子B _ji ；

其中，k _i 表示线型f _i (X︱θ _i )中的参数个数，A _i 表示线型f _i (X︱θ _i )对应的参数贝叶斯采样结果，θ _i 表示A _i 中的第i组参数，

表示由POME确定的第l个参数的后验分布，θ _il 表示第i组参数中的第l个参数值；

4）在求得贝叶斯因子B _ji 的基础上，即可进行水文线型选择和综合分析。

2.根据权利要求1所述的在水文频率计算过程中进行线型选择的方法，其特征在于步骤1）的具体过程为：

（1）对实测水文序列进行可靠性、一致性和代表性审查，并选定合理的水文线型f(x/θ)加以描述；

（2）以线型参数随序列长度变化的分析结果作为依据，并结合实际水文特性确定各参数的初始取值范围；然后选定各参数的先验分布类型π (θ _i ), i=1,2…,n，n为参数个数；

（3）AM-MCMC采样方法初始化，i＝0；

（4）依据所研究序列进行线型参数估计，并将其结果作为采样的初始样本θ ₀ ；

（5）利用式(7)求解该参数样本的协方差矩阵，并计算生成新的参数样本θ*~N(θ _i ，C _i )；

 式中，i ₀ 为初始采样次数；C ₀ 为固定初始协方差矩阵；S _d 为一比例因子；ε为一较小正数，保证C _i 不为奇异矩阵； I _d 为单位矩阵；

 

                          (7)

（6）根据贝叶斯条件概率分布公式(8)，计算新样本的接受概率α；

               (8)

（7）生成一随机数u~U(0,1)，与α比较，若α>u接受θ _i+1 =θ*，否则θ _i+1 ＝θ _i ；

（8）i=i+1,重复步骤(5)－(7)，直至生成足够的样本数为止；

（9）分析各参数样本特征值随样本数增大时的变化规律，判断采样结果是否收敛，若无，则需要增加采样次数直至结果达到稳定；

（10）根据参数采样结果，分析各参数后验分布的特征值，包括最大值(max)、最小值(min)、

、C_υ、C_s和众数，并通过确定各置信区间的值域，定量分析各参数后验分布的不确定性；

（11）根据参数贝叶斯采样结果，求解指定频率对应水文设计值的概率分布，通过定量分析水文设计值概率分布的特征值和各置信区间，掌握水文设计结果的概率分布特性并估计风险性。

Images