CN109902266B - 一种基于Copula函数的河道流量演算方法 - Google Patents

Abstract

本发明公开了一种基于Copula函数的河道流量演算方法，通过收集河段上下断面的流量过程资料，在确定边缘概率分布函数的基础上，求解给定当前时刻上断面入流、下一时刻上断面入流、当前时刻下断面出流时下一时刻下断面出流的条件概率分布函数，进而进行逐时段连续演算得到下断面的出流过程。本发明河道流量演算时不需要假定河段槽蓄量与示储流量存在线性关系，同时考虑了流量的非正态特征，而且能准确地捕捉流量在前后时刻的时间相关性和上下断面的空间相关性，具有可靠的理论基础，更加符合实际情况，为河道流量演算提供一种精度高且可操作性强的新方法。

Description

一种基于Copula函数的河道流量演算方法

技术领域

本发明属于水利工程防洪减灾领域，特别涉及一种基于Copula函数的河道流量演算方法。

背景技术

河道流量演算是根据河道上断面入流过程推求下断面出流过程，在水文预报、水利工程防洪调度中应用广泛。马斯京根法是工程实际中最常采用的方法，该法将河段水流圣维南方程组中的连续方程简化为水量平衡方程，把动力方程简化为河槽蓄泄方程，对简化的方程组联解，得到流量演算方程。根据确定的演算方程，根据上断面入流过程和下断面起始出流，逐时段计算下断面的出流过程。

马斯京根法虽然计算较为简单，但是该法假定河段槽蓄量与示储流量存在线性关系，这一假定可能与实际不相符，影响河道流量演算的精度。此外，实际应用时为了满足这一线性假定通常需要对河段进行分段进行演算。为了改进这一假定，目前已提出了多种马斯京根非线性槽蓄方程或变动参数进行洪水演算。此外，考虑到人工神经网络强大的非线性模拟能力，也有学者建立了人工神经网络的河道流量演算方法。人工神经网络需要大量的参数，且网络结构的选择至今尚无统一而完整的理论指导，只能由经验选定。

基于马斯京根法的思想，河道流量演算本质上表征的是下一时刻下断面出流Q_t+1与当前时刻上断面入流I_t、下一时刻上断面入流I_t+1及当前时刻下断面出流Q_t三者之间的统计关系。采用多元随机统计理论，也可以认为是推求给定I_t、I_t+1和Q_t时，Q_t+1的条件概率分布和条件分位数。Copula函数可以构造具有任意边缘分布的多个随机变量的联合分布，进而求解条件概率的解析表达式，能很好地捕捉多个随机变量的非正态特征和非线性相关结构。目前，没有文献将Copula函数引入河道流量演算研究中。

发明内容

针对现有技术存在的不足，本发明提供了一种基于Copula函数的河道流量演算方法。

为解决上述技术问题，本发明采用如下的技术方案：一种基于Copula函数的河道流量演算方法，包括步骤：

步骤1，收集河段上下断面的流量过程资料；

步骤2，根据步骤1中的河段上下断面的流量过程资料，选取适当的边缘分布线型，估计其参数，确定最优边缘概率分布函数；

步骤3，根据步骤2估计的边缘分布函数，求解给定当前时刻上断面入流、下一时刻上断面入流、当前时刻下断面出流时下一时刻下断面出流的条件概率分布函数；

步骤4，依据步骤3所得的条件概率分布函数，进行逐时段连续演算得到下断面的出流过程。

所述步骤2中，将正态分布、对数正态分布、Gumbel分布、Gamma分布和皮尔逊III型分布作为备选边缘概率分布函数线型，并采用线性矩法估计备选边缘概率分布函数的参数。

所述步骤2中，将一维理论频率与经验频率的均方根误差最小的备选边缘概率分布函数作为最优的边缘概率分布函数。

所述步骤3中，采用非对称Gumbel-Hougaard函数构造当前时刻上断面入流、下一时刻上断面入流、当前时刻下断面出流和下一时刻下断面出流的联合概率分布函数。

所述步骤3中，采用极大似然法估计非对称Gumbel-Hougaard函数的参数。

本发明通过收集河段上下断面的流量过程资料，在确定边缘概率分布函数的基础上，求解给定当前时刻上断面入流、下一时刻上断面入流、当前时刻下断面出流时下一时刻下断面出流的条件概率分布函数，进而进行逐时段连续演算得到下断面的出流过程。

与现有技术相比，本发明的有益效果在于：

本发明河道流量演算时不需要假定河段槽蓄量与示储流量存在线性关系，同时考虑了流量的非正态特征，而且能准确地捕捉流量在前后时刻的时间相关性和上下断面的空间相关性，具有可靠的理论基础，更加符合实际情况，为河道流量演算提供一种精度高且可操作性强的新方法。

附图说明

图1是本发明方法的流程图。

图2是上断面入流、下断面出流流量过程和基于Copula函数的河道流量演算示意图。

具体实施方式

下面通过实施例，并结合附图对本发明作进一步说明。

如图1-图2所示，一种基于Copula函数的河道流量演算方法，收集河段上下断面的流量过程资料，在确定边缘概率分布函数的基础上，求解给定当前时刻上断面入流、下一时刻上断面入流、当前时刻下断面出流时下一时刻下断面出流的条件概率分布函数，进而进行逐时段连续演算得到下断面的出流过程。图1是本实施例的计算流程图，按照以下步骤进行：

1.收集河段上下断面的流量过程资料。

本具体实施中流量的时间尺度为6小时。上下断面的实测流量过程取自同一场洪水过程，从水文年鉴中摘录，上下断面的实测流量样本数量均为N＝36。当前时刻上断面入流、下一时刻上断面入流、当前时刻下断面出流和下一时刻下断面出流分别记为I_t、I_t+1、Q_t和Q_t+1。如图2所示，给出了上断面入流(实线)和下断面出流(虚线)的流量过程示意图。

2.确定上断面入流和下断面出流的边缘概率分布函数。

根据步骤1中的河段上下断面的流量过程资料，选取适当的边缘分布线型，并估计其参数，最后确定最优边缘概率分布函数，本步骤包括三个子步骤：

2.1备选边缘概率分布函数线型

考虑到上断面入流I_t、I_t+1的边缘概率分布函数是相同的，下断面出流Q_t和Q_t+1的边缘概率分布函数也是相同的。本具体实施中，只需要估计上断面入流I和下断面出流Q的边缘概率分布函数。I_t、I_t+1、Q_t和Q_t+1的边缘概率分布函数分别记为u₁＝F_I(i_t)、u₂＝F_I(i_t+1)、u₃＝F_Q(q_t)和u₄＝F_Q(q_t+1)。

由于河段上下断面的流量的总体分布频率线型是未知的，通常选用能较好拟合多数水文样本资料系列的线型。本具体实施中采用将正态分布、对数正态分布、Gumbel分布、Gamma分布和皮尔逊III型分布作为备选边缘概率分布函数线型。

2.2估计边缘分布线型的参数

当频率分布线型选定后，接下来需要进行估计频率分布的参数。目前常用的方法主要有矩法、极大似然法、适线法、概率权重矩法、权函数法和线性矩法等。其中，线性矩法的特点是对序列的极大值和极小值没有常规矩那么敏感，求得的参数估计值比较稳健，是目前国内外公认的高精度参数估计方法。

本具体实施中采用L-矩法估计备选边缘概率分布函数的参数。

2.3最优边缘概率分布函数确定

采用均方根误差(Root Mean Square Error，RMSE)准则评价边缘分布的一维理论频率与经验频率拟合情况，RMSE值越小，说明拟合效果越好。

式中：F(x_i)为观测值x_i的理论频率；m(i)为实测系列中满足x≤x_i的观测值个数，n为样本长度。

本具体实施中，采用RMSE值最小的备选边缘概率分布函数作为最优的边缘概率分布函数。

3.求解上断面入流和下断面出流的条件概率分布函数。

根据步骤1中得到的当前时刻上断面入流I_t、下一时刻上断面入流I_t+1、当前时刻下断面出流Q_t和下一时刻下断面出流Q_t+1以及步骤2中估计的边缘概率分布函数，选取适当的Copula函数作为连接函数构造I_t、I_t+1、Q_t和Q_t+1的联合概率分布函数，并估计其参数，求解给定I_t、I_t+1和Q_t时，Q_t+1的条件概率分布函数，本步骤包括三个子步骤：

3.1选择Copula函数

令I_t、I_t+1、Q_t和Q_t+1的概率密度函数分别为f_I(i_t)、f_I(i_t+1)、f_Q(q_t)和f_Q(q_t+1)。由Copula函数理论可知，I_t、I_t+1和Q_t的联合概率分布函数可以用一个三维Copula函数C表示：

F_3d(i_t,i_t+1,q_t)＝C_3d(F_I(i_t),F_I(i_t+1),F_Q(q_t))＝C_3d(u₁,u₂,u₃) (5)

一般而言，I_t、I_t+1和Q_t之间存在正相关关系且相关性不对称。本具体实施中，采用三维非对称Gumbel-Hougaard Copula函数构造I_t、I_t+1和Q_t的联合概率分布函数，表达式如下：

其中，参数θ_3d＝{θ₂,θ₁}为三维Copula函数的参数，且满足θ₂≥θ₁≥1。

同理，借助Copula函数，可以将I_t、I_t+1、Q_t和Q_t+1的联合概率分布函数写为下式：

F_4d(i_t,i_t+1,q_t,q_t+1)＝C_4d(F_I(i_t),F_I(i_t+1),F_Q(q_t),F_Q(q_t+1))＝C_4d(u₁,u₂,u₃,u₄) (7)

同样地，I_t、I_t+1、Q_t和Q_t+1之间存在正相关关系且相关性不对称。本具体实施中，采用四维非对称Gumbel-Hougaard Copula函数构造I_t、I_t+1、Q_t和Q_t+1的联合概率分布函数，其表达式如下：

其中，参数θ_4d＝{θ₃,θ₂,θ₁}为四维Copula函数的参数，且满足θ₃≥θ₂≥θ₁≥1。

3.2估计Copula函数的参数

目前估计Copula函数的参数的常用方法有Kendall相关系数法、极大似然法、边际推断法和核密度估计法等。其中，极大似然法的思想是将似然函数关于参数最大化，从而得到参数向量的估计值，目前广泛应用于三维及以上Copula函数的参数估计中。本具体实施中，采用极大似然法估计三维非对称和四维非对称Gumbel-Hougaard Copula函数的参数。

3.3给定I_t、I_t+1和Q_t，求解Q_t+1的条件概率分布函数。

给定I_t、I_t+1和Q_t取值时，对应Q_t+1的取值存在无数种可能性，只是出现不同取值的概率有所不同，存在着一个条件概率分布函数

F(q_t+1|i_t,i_t+1,q_t)＝P_c(Q_t+1≤q_t+1|I_t＝i_t,I_t+1＝i_t+1,Q_t＝q_t) (9)

其中，P_c代表事件发生的概率值。

借助Copula函数，条件概率分布函数F(q_t+1|i_t,i_t+1,q_t)可重新写为：

其中，

为三维非对称Gumbel-Hougaard Copula函数的密度函数。

4.进行逐时段连续演算得到下断面的出流过程。

得到Q_t+1的条件概率分布函数F(q_t+1|i_t,i_t+1,q_t)后，根据数理统计原理，可以计算得到中位数q_t+1,m作为Q_t+1的点估计值，通过下式求解：

F(q_t+1,m|i_t,i_t+1,q_t)＝0.5 (11)

本具体实施中采用二分法试算求解式(11)得到数值解。

进行逐时段连续演算时，共有(N-1)个时段，即t＝1,…,N-1。本具体实施中，从t＝1开始，令第1时刻下断面出流等于上断面入流，即q₁＝i₁。已知第1时刻上断面入流i₁、第2时刻上断面入流i₂和第1时刻下断面出流q₁，就可以根据公式(11)求解第2时刻下断面出流Q₂的条件概率分布函数F(q₂|i₁,i₂,q₁)，进而推求Q₂的中位数q_2,m作为第2时刻下断面出流流量值。以此类推，最后t＝N-1时，推求得到第N时刻下断面出流Q_N的中位数q_N,m作为第N时刻下断面出流流量值。这样就可以实现由河道上断面入流过程I_t推求下断面出流过程Q_t。如图2所示，给出了基于Copula函数的河道流量演算示意图，其中，点划线表示基于Copula函数演算法推求的下断面出流过程。

综上，本发明通过收集河段上下断面的流量过程资料，在确定边缘概率分布函数的基础上，求解给定当前时刻上断面入流、下一时刻上断面入流、当前时刻下断面出流时下一时刻下断面出流的条件概率分布函数，进而进行逐时段连续演算得到下断面的出流过程。本发明河道流量演算时不需要假定河段槽蓄量与示储流量存在线性关系，同时考虑了流量的非正态特征，而且能准确地捕捉流量在前后时刻的时间相关性和上下断面的空间相关性，具有可靠的理论基础，更加符合实际情况，为河道流量演算提供一种精度高且可操作性强的新方法。

Claims (5)

1.一种基于Copula函数的河道流量演算方法，其特征在于包括以下步骤：

步骤1，收集河段上下断面的流量过程资料；

步骤2，根据步骤1中的河段上下断面的流量过程资料，选取适当的边缘分布线型，估计其参数，确定最优边缘概率分布函数；

步骤3，根据步骤2估计的边缘分布函数，求解给定当前时刻上断面入流、下一时刻上断面入流、当前时刻下断面出流时下一时刻下断面出流的条件概率分布函数；

步骤4，依据步骤3所得的条件概率分布函数，进行逐时段连续演算得到下断面的出流过程。

2.如权利要求1所述的一种基于Copula函数的河道流量演算方法，其特征在于：所述步骤2中，将正态分布、对数正态分布、Gumbel分布、Gamma分布和皮尔逊III型分布作为备选边缘概率分布函数线型，并采用线性矩法估计备选边缘概率分布函数的参数。

3.如权利要求1所述的一种基于Copula函数的河道流量演算方法，其特征在于：所述步骤2中，将一维理论频率与经验频率的均方根误差最小的备选边缘概率分布函数作为最优边缘概率分布函数。

4.如权利要求1所述的一种基于Copula函数的河道流量演算方法，其特征在于：所述步骤3中，采用非对称Gumbel-Hougaard函数构造当前时刻上断面入流、下一时刻上断面入流、当前时刻下断面出流和下一时刻下断面出流的联合概率分布函数。

5.如权利要求1所述的一种基于Copula函数的河道流量演算方法，其特征在于：所述步骤3中，采用极大似然法估计非对称Gumbel-Hougaard函数的参数。

Images