Copula-MCP洪水概率预报方法研究与应用

沈 婕，梁忠民，胡义明，王 军，李彬权

(河海大学水文水资源学院，南京 210098)

0 引 言

水文模型都是对自然子流域的一种概化，无论结构是否复杂、参数是否准确，都只是对水文过程的近似，无法真实地还原自然的水文过程，因此实时洪水的预报中不可避免地存在不确定性，用于描述这种不确定性的洪水概率预报逐渐得到了重视和应用[1,2]。

洪水概率预报途径可以概括为两类[3]：一是总误差分析途径，直接在确定性预报结果的基础上，量化分析最终的预报不确定性，进而实现概率预报；二是不确定性要素耦合途径，即分析预报过程中各环节的主要因子的不确定性，估计其概率分布，再将这些因子的概率分布耦合到洪水预报模型中，从而实现概率预报。从实时洪水预报角度，第一类的总误差分析途径计算相对简单，更满足时效性要求，因此在洪水作业预报中具有广泛应用。当然，不管是哪类途径，都是在确定性预报(定值预报)的基础上实现概率预报，因此，如何提高原始定值预报结果的精度，是降低不确定性、提升概率预报可靠性或合理性的关键。

在总误差分析途径框架下，由Krzysztofowicz等[1]提出的贝叶斯预报系统(BFS)是最典型的概率预报方法，其中水文不确定性处理器[2](HUP)在贝叶斯理论方法的基础上，对预报变量进行了正态化变换与线性假设，推求预报量后验分布函数。王善序[4]对贝叶斯概率水文预报理论做了简单的介绍，认为其可以不需任何附加假定与要求地定量考虑预报的不确定性，同时也指出了其中的线性-正态假设在水文过程中未必适用。邢贞相等[5]采用BP网络描述了水文过程的非线性特征，并应用AM-MCMC方法进行概率预报，在此过程中仍采用了正态性假设。Todini等[6]提出了模型条件处理器(MCP)，采用截尾正态联合分布揭示了不同流量量级时预报值与实测值的关系。梁忠民等提出了考虑误差异分布[7](Heterogeneity of Error Distributions, EHDA)方法，分析了预报误差的异分布性，进而实现洪水概率预报。

上述这些方法中，均需对原始变量进行正态化变换，并/或对正态空间中的变量或似然函数进行线性假设。然而，正态化变换过程中可能会造成信息的丢失或偏差[8]，变量或似然函数的线性假设也具有较强的主观性。为避免正态化变换和线性假设，刘章君等[9]推导了基于Copula函数的贝叶斯转移预报(BTF)方法的后验密度函数，提出了CBTF与CMHUP方法，取得了更好的概率预报效果。本文采用这种思路，将Copula函数与MCP方法结合，构建了Copula-MCP洪水概率预报模型，并以淮河王家坝断面为例进行了应用研究。

1 基于Copula-MCP的洪水概率预报模型

1.1 Copula函数

Copula函数[10,11]是将随机变量的联合分布函数与各自的边缘分布函数相连接的连接函数。设X、Y为连续随机变量，其边缘分布函数分别为Fx(x)与Fy(y)，联合分布函数为F(x,y)，若Fx(x)与Fy(y)连续，则存在唯一确定的Copula函数Cθ(u,v)使得：

F(x,y)=Cθ(u,v)

(1)

式中：θ为Copula函数待定参数，u=Fx(x)，v=Fy(y)。

Copula函数将变量的相关性结构与边缘分布分开处理，可以考虑变量间的非线性与非正态关系，没有任何限制。

本文采用Archimedean Copula函数族来构造实测流量与模拟流量的联合分布函数。常见的Archimedean Copula函数有Gumbel-Hougaard Copula函数、Clayton Copula函数和Frank Copula函数，其二维表达式如下：

Gumbel-Hougaard Copula函数：

(2)

Clayton Copula函数：

(3)

Frank Copula函数：

(4)

式中：θ为Copula函数待定参数；u与v分别为单变量的边缘分布函数。

为了对Copula函数的拟合优度进行检验，本文选用K-S检验法，其统计量D计算公式如下：

(5)

联合分布的经验累积分布与理论分布的拟合情况通过RMSE准则来评价，并根据最小RMSE值来优选Copula函数。RMSE准则计算公式如下：

(6)

式中：n为样本容量；Fc(xi)为样本集的累积分布函数；F(xi)为理论分布。

1.2 模型条件处理器(MCP)

(7)

式中：n为历史数据的数量；i为变量的排位顺序。

(8)

(9)

(10)

预测不确定性一般定义为以模型预测值为条件的预报量的条件分布，在已知变量联合分布和边缘分布时，通过计算联合分布与边缘分布的比值，就可以推得预测不确定性，即预报量的条件概率分布：

(11)

由此可知，在正态空间中，预测不确定性相当于均值和方差如下的正态分布：

(12)

1.3 Copula-MCP

(13)

(14)

条件概率分布函数对应的概率密度函数为：

(15)

2 应用实例

如图1所示。王家坝站是淮河上游的总控制站，其洪水由上游干流的息县、左右岸控制站班台和潢川、息-潢-班至王家坝区间共四部分组成。息潢班-王家坝区间集水面积为7 110 km2，干流长度为360 km。本文进行王家坝洪水预报时，先采用经验降雨径流模型API，对息潢班-王家坝区间洪水进行预报，再与经马斯京根法演算至王家坝的息县、潢川、班台实测流量过程叠加，得到确定性的洪水预报结果。在此基础上，采用Copula-MCP模型推求以确定性预报结果为条件的预报变量的概率分布，实现王家坝站的洪水概率预报。

图1 淮河息潢班～王家坝区间流域示意图

2.1 确定性预报

2.1.1 区间洪水预报

采用加权平均法计算王家坝、光山等15个雨量站的区间平均雨量，并采用流量起涨前30天雨量计算流域前期影响雨量Pa，径流深R采用降雨径流相关图法进行计算，降雨径流相关图如图2所示。

图2 淮河息潢班～王家坝区间降雨径流相关图

以降雨中心位置为依据，采用不同的单位线进行区间汇流计算。单位线参数如图3所示。其中，a线适用于暴雨中心在淮北；b线适用于降雨均匀；c线适用于降雨中心在淮南。

图3 淮河息潢班～王家坝区间单位线图

2.1.2 河道洪水预报

采用马斯京根法，选用先演后合的方式对息县、潢川、班台(大洪河)上游来水进行河道洪水演算。相关参数见表1。

表1 河道洪水演算参数表

2.2 API模型的确定性预报结果

对息潢班～王家坝区间1990-2010年间的25场洪水进行模拟，其中20场作为率定，5场作为验证，计算时段Δt=6 h，精度统计见表2。可以看出，确定性预报的结果良好，洪峰的合格率(相对误差绝对值在20%以内)为68%，确定性系数的优秀率(确定性系数大于0.85)为84%。

表2 API模型精度统计表

2.3 Copula-MCP方法应用

2.3.1 边缘分布的确定

对实测及模型计算的洪水序列(或两个随机变量)，一般可假设其边际分布服从Log-Weibull分布[13]。本次采用线性矩法估计两个变量边际分布的统计参数，为检验参数估计结果是否合理，采用K-S检验法分别对两个随机变量经验累积分布与理论分布进行拟合优度检验，边际分布参数估计值与K-S检验统计量D值结果见表3。由表3可知，在5%的显著性水平下，两个变量的K-S检验统计量D均小于相应的临界值(0.039 5)，均通过拟合优度检验。实测流量 与预报流量 样本点据与其边际分布拟合情况如图4所示。由图4可知，理论边际分布与样本点据拟合较好，则使用该分布作为边际分布是合理的。

表3 y和边际分布参数估计值和K-S检验

图4 理论Weibull与样本点据拟合图

2.3.2 联合分布的建立

(16)

2.3.3 洪水概率预报结果及分析

表4 联合分布参数估计、检验及优选结果

采用Copula-MCP方法进行洪水概率预报，其结果见表5，其中提供了90%置信区间覆盖率、洪水概率分布期望值过程的确定性系数(DC)、洪峰期望值相对误差等结果；为对比分析，表中亦列出了API及MCP模型的结果。图5为其中一场洪水(19980509)的API、MCP、Copula-MCP预报结果与实测流量的对比图。由表5可知，Copula-MCP模型的洪峰相对误差合格率为92%，高于MCP模型的84%与API模型的68%；Copula-MCP模型的平均确定性系数为0.93，高于MCP模型的0.92与确定性模型的0.91，Copula-MCP模型确定性系数的优秀率为92%，亦高于MCP模型的88%与API模型的84%；从90%置信区间的覆盖率来看，Copula-MCP模型也优于MCP模型。

表5 API、MCP与Copula-MCP三个模型结果对比表

图5 19980509场次洪水Copula-MCP概率预报结果

3 结 论

(1)与目前常用的洪水概率预报面向相较，采用Copula函数直接构建实测流量与预报流量的联合分布函数，不需要进行变量的正态化变换及线性假设，避免变量间相关信息的丢失和洪水过程非线性特征的改变，由此构建的洪水概率预报模型，理论上应更具优势。

(2)洪水过程的不同发展阶段或不同流量量级的洪水，其预报误差规律可能不同，将描述这种异误差分布的MCP模型与Copula函数结合，构建的Copula-MCP洪水概率预报模型具有更好的适用性。本文研究实例表明，Copula-MCP模型整体上优于MCP模型，可进一步提高洪水预报精度、增加概率预报的可靠性。

□