基于copula函数的山洪灾害概率雨型推求方法与流程

本发明属于山洪灾害研究技术领域，具体涉及一种基于copula函数的山洪灾害概率雨型推求方法。

背景技术：

受全球气候变化的影响，极端暴雨天气事件发生频率和强度逐渐加剧，全球或局地的降水结构发生了显著变化，导致水汽循环发生改变，由此引发一系列气象灾害。而近几十年中国极端暴雨天气的强度和频次亦呈现逐渐增加趋势，由此造成我国小流域山洪灾害肆虐严重，对人类生活、社会发展和自然生态系统产生了严重影响。

极端暴雨天气是由降雨在时空分布上的不确定性所造成的，降雨作为山洪灾害预警指标计算的关键输入数据，其不确定性造成了山洪灾害预警的不确定性。降雨总量、峰值雨强、降雨时程分布等雨型基本因素均对预警指标有较大影响，由于在现实情境中，降雨具有随机性，基于单因素发生概率推求的山洪灾害暴雨雨型已经不能满足实际需求，如何基于多因素的联合发生概率，科学合理地确定一定概率下的山洪灾害暴雨雨型对山洪灾害预警预报具有重要意义。

常用的雨型确定方法有模式雨型、芝加哥雨型、三角形雨型、p&c法雨型、huff雨型等，但这些雨型推求方法大多是基于城市洪涝灾害而提出的，对于这些方法应用到小流域山洪灾害暴雨雨型的确定有一定的局限性。如：p&c法雨型虽然考虑了一定的概率分布，但是仅仅只是将雨峰放在最可能出现的位置，其它雨峰位置没有考虑。并且在其余各时段分配比例偏向于极端化，致使降雨过程陡涨陡落，不符合实际降雨过程；huff雨型根据最大雨强发生的具体时间段,对其作出多种不同频率的无因次时间分配过程，推求不同概率下降雨时程分配，但是该方法仅仅将最大雨强位置按照降雨历时分为四种，没有考虑雨峰位置的随机性影响；三角形雨型和芝加哥雨型在应用时雨峰位置固定且不考虑随机概率分布，虽然应用简单边界，但不满足降雨过程的随机性和不确定性，易与实际降雨过程产生较大出入，且芝加哥法雨型不能用于历时大于2h的降雨时程分配，不满足小流域山洪灾害暴雨特性。

技术实现要素：

本发明为了解决现有技术中的不足之处，提供一种基于copula函数的山洪灾害概率雨型推求方法；该方法通过耦合总雨量与峰值雨强构建联合分布合数，同时结合雨峰位置系数函数，充分考虑总雨量、峰值雨强、雨峰位置的随机性和不确定性，推求不同概率下总雨量、峰值雨强、雨峰位置的联合雨型。并且该方法基于大量实测资料，充分考虑小流域降雨特性，将除雨峰位置外其余时段分配比例平均化，避免降雨过程中出现极端化的陡涨陡落现象。

本发明采用如下技术方案：基于copula函数的山洪灾害概率雨型推求方法，包括以下步骤，

s1、以小流域为单元收集典型小流域降雨资料；

s2、确定所求时段雨型的雨峰位置系数函数；

s3、采用copula函数构造时段峰值雨强和时段总雨量的联合分布函数，并估计copula函数参数；

s4、基于联合分布函数推求一定总雨量所对应的峰值雨强条件概率分布；

s5、根据条件概率分布推求相同概率下峰值雨强；

s6、确定所求时段某一概率下除雨峰位置外，其余时段的降雨分配比例；

s7、确定研究区某一概率下的雨型。

优选地，s1步骤中需根据典型小流域降雨资料依据一定原则进行雨场划分，整理出一定时段的降雨场次；雨场划分原则如下所述：

a：长时间间隔的场次降雨分割采用前后连续24小时无降雨的原则，即一场连续降雨过程，其前后连续24小时均无降雨；

b：中等时间间隔的场次降雨分割采用前后连续12小时无降雨的原则，即一场连续降雨过程，其前后连续12小时均无降雨；

c：短时间间隔的场次降雨分割采用以4mm为降雨场次的原则，即一场连续降雨过程中，以小时雨量大于4mm时相应时刻为该场次降雨的开始时刻，以小时雨量连续3小时小于4mm的开始时刻为该场次降雨的结束时刻。

优选地，s2步骤具体包括如下步骤：

ⅰ：统计雨峰位置系数；

ⅱ：选取适当的概率分布函数线型,估计概率分布函数的参数，确定雨峰位置系数函数，并推求不同概率下雨峰位置。

优选地，确定雨峰位置系数具体包含以下步骤，

（1）：确定场次降雨中雨峰出现时间；

（2）：确定场次降雨总历时；

（3）：确定雨峰出现时间与降雨总历时的比值，即为雨峰位置系数。

优选地，s3步骤具体包括如下步骤：

a：选取适当的边缘概率分布函数线型,确定峰值雨强边缘分布函数；

b：选取适当的边缘概率分布函数线型，确定时段总雨量边缘分布函数；

c：采用copula函数构造峰值雨强和时段总雨量的联合分布函数，并估计copula函数相关参数。

优选地，步骤a、b中边缘概率分布函数线型均采用p-ⅲ型分布，并采用矩法估计边缘分布函数的相关参数。

优选地，步骤c中copula函数采用二维g-hcopula函数，并采用kendall秩相关性系数法估计二维g-hcopula函数的参数。

优选地，s6步骤具体包括如下步骤：

①：依据雨峰位置系数函数，确定某一概率下的雨峰位置。

②：选取样本资料中与①中雨峰位置相同的样本。

③：确定各样本中除雨峰位置外，其余各时段的降雨分配比例。

④：将相同时段降雨量比例平均化，即为某一概率下除雨峰位置外，其余各时段的降雨分配比例，且各时段的降雨分配比例和为1。

优选地，确定某一概率下的雨峰位置；确定某一概率下的峰值雨强和和总雨量；确定除雨峰位置外其余各时段的降雨分配比例；将总雨量按比例分配到各个时段，即为某一概率下的暴雨雨型。

采用上述技术方案，与现有研究相比，本发明具有以下优点：

本发明所述的一种基于copula函数的山洪灾害概率雨型确定方法，其优点在于，

通过收集以小流域为单元的典型小流域的降雨资料，统计出所需时段降雨的雨峰位置系数、峰值雨强、总雨量，在此基础上构建雨峰位置系数概率密度函数和基于copula函数的峰值雨强和总雨量的联合分布函数，充分考虑了雨型确定过程中雨峰位置、峰值雨强、总雨量三个不确定性因素，并将三者以概率密度函数有机结合起来，得以推求出某一概率下的山洪灾害暴雨雨型，为山丘区小流域概率雨型确定提供参考和依据。

附图说明

图1是本发明所述基于copula函数的山洪灾害概率雨型确定方法的流程示意图；

图2是本发明所述的前后24小时以上无降雨的一场连续降雨分割示意图。

图3是本发明所述的前后12小时以上无降雨的一场连续降雨分割示意图。

图4是本发明所述的以降雨4mm为界的一场连续降雨分割示意图。

图5是本发明所述的实测6小时总雨量与相应场次降雨峰值雨强散点图。

图6是本发明所述的经验联合分布函数与理论联合分布函数拟合图（θ=1.36）。

图7是本发明所述的当总雨量为50%时，峰值雨强的条件概率分布图。

具体实施方式

如图1所示，本发明的基于copula函数的山洪灾害概率雨型推求方法，包括以下步骤：

s1、以小流域为单元收集典型小流域降雨资料；

s2、确定所求时段雨型的雨峰位置系数函数；

s3、采用copula函数构造时段峰值雨强和时段总雨量的联合分布函数，并估计copula函数参数；

s4、基于联合分布函数推求一定总雨量所对应的峰值雨强条件概率分布；

s5、根据条件概率分布推求相同概率下峰值雨强；

s6、确定所求时段某一概率下除雨峰位置外，其余时段的降雨分配比例；

s7、确定研究区某一概率下的雨型。

优选地，s1步骤中需根据典型小流域降雨资料依据一定原则进行雨场划分，整理出一定时段的降雨场次；雨场划分原则如下所述：

a：长时间间隔的场次降雨分割采用前后连续24小时无降雨的原则，即一场连续降雨过程，其前后连续24小时均无降雨；如附图2所示；

b：中等时间间隔的场次降雨分割采用前后连续12小时无降雨的原则，即一场连续降雨过程，其前后连续12小时均无降雨；如附图3所示；

c：短时间间隔的场次降雨分割采用以4mm为降雨场次的原则，即一场连续降雨过程中，以小时雨量大于4mm时相应时刻为该场次降雨的开始时刻，以小时雨量连续3小时小于4mm的开始时刻为该场次降雨的结束时刻，如附图4所示。

优选地，s2步骤具体包括如下步骤：

ⅰ：统计雨峰位置系数；

ⅱ：选取适当的概率分布函数线型,估计概率分布函数的参数，确定雨峰位置系数函数，并推求不同概率下雨峰位置。

优选地，确定雨峰位置系数具体包含以下步骤，

（1）：确定场次降雨中雨峰出现时间；

（2）：确定场次降雨总历时；

（3）：确定雨峰出现时间与降雨总历时的比值，即为雨峰位置系数。

优选地，s3步骤具体包括如下步骤：

a：选取适当的边缘概率分布函数线型,确定峰值雨强边缘分布函数；

b：选取适当的边缘概率分布函数线型，确定时段总雨量边缘分布函数；

c：采用copula函数构造峰值雨强和时段总雨量的联合分布函数，并估计copula函数相关参数。

优选地，步骤a、b中边缘概率分布函数线型均采用p-ⅲ型分布，并采用矩法估计边缘分布函数的相关参数。

优选地，步骤c中copula函数采用二维g-hcopula函数，并采用kendall秩相关性系数法估计二维g-hcopula函数的参数。

优选地，s6步骤具体包括如下步骤：

①：依据雨峰位置系数函数，确定某一概率下的雨峰位置。

②：选取样本资料中与①中雨峰位置相同的样本。

③：确定各样本中除雨峰位置外，其余各时段的降雨分配比例。

④：将相同时段降雨量比例平均化，即为某一概率下除雨峰位置外，其余各时段的降雨分配比例，且各时段的降雨分配比例和为1。

优选地，s7步骤包括确定某一概率下的雨峰位置；确定某一概率下的峰值雨强和和总雨量；确定除雨峰位置外其余各时段的降雨分配比例；将总雨量按比例分配到各个时段，即为某一概率下的暴雨雨型。

下面结合实例对本发明做出进一步的详细说明。

以下基于某研究区小流域某雨量站1982-2013年实测降雨资料，推求某典型小流域50%概率下的6小时山洪灾害概率雨型为例，阐述基于copula函数的山洪灾害概率雨型确定方法。具体包括如下步骤：

1.以研究区典型小流域为单元进行降雨资料的收集，本发明以a省b县的典型小流域c为研究对象，选取某雨量站1982-2013年实测降雨资料，按照短时间间隔的场次降雨分割原则划分场次暴雨，如附图4所示。

1.1．根据各场次暴雨的降雨历时选择样本，本例中,选择降雨历时在6小时左右的场次暴雨作为样本,共选择了20场次降雨样本作为基础资料。

2.确定该20场6小时降雨样本中雨峰位置系数函数。

2.1按照6小时时段雨量中最大降雨量所处的位置统计20场暴雨样本中雨峰位置。

2.2.计算样本中雨峰位置系数，即计算各场次暴雨中最大降雨量所处时段与场次降雨历时的比例。

2.3.选择概率分布函数线型。

20世纪60年代以来，我国经过多年分析比较规定p-ⅲ型分布作为我国水文分析计算的线型。本发明采用p-ⅲ型分布作为雨峰位置系数的分布线型。

2.4.估计概率分布函数的参数。

估计概率分布函数的方法有矩法、适线法、线性矩法、优化适线法、概率权重矩法等，基于各方法的优缺点及适用性，本发明采用矩法估计边缘分布函数的参数。并在此基础上确定雨峰位置频率曲线。

2.5.确定雨峰位置系数函数，并推求50%概率下雨峰位置

在参数确定的基础上，将雨峰位置系数按照降序排列后，通过水文频率计算绘制雨峰位置频率曲线。并确定50%概率下的雨峰位置系数为0.58。

3.确定6小时峰值雨强和时段总雨量的联合分布函数。

3.1.统计样本中6小时时段雨量峰值雨强（最大1h降雨量）。

3.2.统计样本中6小时时段总雨量。

3.3.选择边缘概率分布函数线型。

边缘函数的选择一般要求适用性较强、拟合程度较高、应用简单便捷、有充分的理论支撑。基于此，选择我国经过多年不同线型的分析比较后，适合我国大部分地区、认可程度较高、且拟合程度较高的p-ⅲ型分布作为峰值雨强和总雨量的边缘分布线型。

3.4．估计边缘分布函数的参数及确定函数。

采用步骤2中雨峰位置概率分布函数相同的参数估计方法，即采用矩法估计边缘分布函数的参数。在此基础上，计算1982-2013年共20组实测数据6小时降雨总雨量和峰值雨强的边缘分布函数，并由边缘分布函数确定50%概率下的总雨量为79.44mm。

3.5.利用copula函数构建峰值雨强和总雨量的概率分布函数并计算。

3.5.1.选择copula函数

假设p、s分别表示总雨量和峰值雨强，p、s分别表示相应的实现值，fp(p)、fs(s)分别为其相应的边缘概率分布函数，相应的概率密度函数为fp(p)、fs(s)。由sklar定理可知，p、s的联合分布函数可以用一个二维copula函数来表示：

fp,s(p,s)=cθ(fp(p),fs(s))=cθ(u,v)（1）

其中θ为copula函数的参数，u=fp(p)，v=fs(s)为边缘概率分布函数。

在本发明实例过程中，采用二维g-hcopula函数建立峰值雨强和总雨量的联合分布函数，以及采用kendall秩相关系数法估计参数，其公式如下：

其中样本的kendall秩相关系数τ通过下式计算：

其中：n为样本长度，。

3.5.2.计算理论频率和经验频率。

经验分布函数：将1982-2013年共20组实测数据6小时降雨总雨量p和其对应降雨过程中的峰值雨强s绘制成散点图，则其经验分布函数根据散点图计算，如f(63,33)=图上线框包围的点数/总点数。总雨量和其对应的峰值雨强散点图如附图5所示。

理论频率采用g-hcopula函数求解，参数计算取值为：τ=0.57，θ=2.33，所求理论联合频率分布函数值和经验联合概率分布函数值的拟合程度良好，如附图6所示。

3.6.采用copula函数，计算当一定概率总雨量下的峰值雨强条件概率分布函数。

给定总雨量p取值p时，所对应的峰值雨强s的取值并非唯一，只是出现的概率有所不同，因此存在着两者的条件概率分布函数，以下给出当p取值一定时的条件分布函数。

将式（2）代入得：

首先由边缘分布函数fs(s)确定50%概率下的总雨量为79.44mm。在此基础上，采用条件概率分布函数及copula联合分布函数计算相同概率下峰值雨强，当总雨量概率为50%时，峰值雨强的条件概率分布图如附图7所示，由此得到相同概率下峰值雨强为28.39mm。

4．确定50%概率下6小时雨型中除雨峰位置外，其余时段的降雨分配比例；

第2步中已经确定了50%概率下雨峰所在位置为0.58，即在6小时降雨中所占时段为第3时段。然后统计实测资料中与该概率下雨峰位置相同的样本，计算各样本中除雨峰位置外其余各时段的降雨比例并求平均，即50%概率下除雨峰位置外，其余各时段的分布比例，且满足比例之和为1。据实测资料统计有2场降雨雨峰位置在第3时段，求平均后得出其余各时段比例如下表所示。

表16小时暴雨雨型雨峰位置外其余各时段比例

5.确定研究区小流域c的50%概率下的暴雨雨型雨型。

将50%概率下的峰值雨强28.39mm放到50%概率的雨峰位置上。将50%概率下的总雨量减去相同概率下的峰值雨强得到其余5个时段的总雨量，再乘以相应比例，即可得到除雨峰位置外其余时段的降雨量。综合雨峰位置、峰值雨强、总雨量即为该小流域c的50%概率下6小时的暴雨雨型，如表2所示。

表26小时暴雨雨型

综上，本发明通过构建小流域雨峰位置概率分布函数，以及基于copula函数的总雨量和峰值雨强联合概率分布函数，推求出某一概率下的山洪灾害暴雨雨型。本发明综合考虑了雨峰位置、峰值雨强、总雨量三者之间的不确定性影响及内在关联性，能够客观的反应实际降雨过程中降雨随机性较大这一特征，为小流域山洪灾害暴雨雨型的确定提供参考和依据。

本实施例并非对本发明的形状、材料、结构等作任何形式上的限制，凡是依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化与修饰，均属于本发明技术方案的保护范围。