(2)多变量联合分布函数Copula 三维copula重现期matlab制图代码

1 多变量联合分布函数Copula概要

目前,基于Copula函数的多变量水文频率分析主要集中于两变量,计算方法比较成熟,应用非常广泛。然而,研究中很多变量不止是两个,比如水文中的干旱包含干旱事件强度-面积-历时的三个方面的特征,而洪水过程一般包含洪峰、时段洪量和历时等多个方面特征,具有较强的复杂性。当然这些特征可以不止三个,还会有更多。这里我们以三维洪水频率分析为例。在洪水频率分析时通常将一个完整的洪水过程可划分2~3个时段,不同历时的洪量之间往往是相依的,与洪峰流量之间的相关性随着历时的增加而减弱。此时仅用二维联合分布难以描述整个洪水过程,也不能够完整地表述水文事件的特性,这需要将Copula函数向高维拓展。那么怎么计算呢?

2 多维联合分布公式解析

一个n维联合分布分布函数可表示为:
在这里插入图片描述
式中,n为变量( )的维度; 为n维变量联合分布分布函数; 为n维变量的概率密度函数。
根据Sklar定理,一个n维联合分布函数可用Copula函数表示为:
在这里插入图片描述
式中, 为第n个变量的边缘分布,表示为 ;C为Copula连接函数。当变量维度数目显著增加时,Copula函数的灵活性和准确性受到了限制[17]。对整个洪水过程而言,洪峰流量和不同历时的时段洪量之间的相关性有所不同,随着历时的增加,相关性也在逐步减弱,相关性参数的也会发生相应的变化,然而,对称性Copula函数的相关性参数是不变的,难以描述复杂的洪水过程各特征变量之间的相关性,这需构建非对称性的Copula函数结构。
Hierarchical Archimedean Copulas(HAC)在多维联合分布起着至关重要的作用,打破了变量的限制,可以通过嵌套Archimedean Copula来描述不同变量之间的相关性,其组合方式主要包括两种:完全嵌套和部分嵌套。下面以四维联合分布为例进行说明,完全嵌套模式由变量u4和u3连接成C3,u2和C3连接形成C2,u1和C2连接形成C1,其形式见图 1-a和公式(4)),部分嵌套模式由变量u4和u3连接成C3,u2和u1连接形成C2,C3和C2连接形成C1,其形式见图 1-b和公式(5)。
在这里插入图片描述
考虑到洪水过程特殊的相依性,一般而言,较长历时的洪量包含较短历时的洪量,完全嵌套方式符合洪水过程的特点。针对于此,本次研究采用该种嵌套模式来描述整个洪水的各特征变量。
在这里插入图片描述
图 1四维完全嵌套Archimedean 结构(a)和部分嵌套Archimedean 结构(b)
对于这种三维的联合重现期和同现重现期的计算公式如下:
在这里插入图片描述

3 三维代码出图效果

部分matlab代码片段
data=importdata(‘G:\copula-work\三维copula\Triaxialcopula1\copula-data.csv’);
data=data.data;

%以下是示例数据,输入经过边缘分布拟合的数据
duration=data(:,2);
peak=data(:,3);
volume=data(:,4);
%以下是示例数据,选择最优的拟合分布
[D_U1, PD_U1] = marginfitdist(duration);
[D_U2, PD_U2] = marginfitdist(peak);
[D_U3, PD_U3] = marginfitdist(volume);

% 计算拟合的边缘分布的cdf
EP1 = cdf(PD_U1{1},duration);
EP2 = cdf(PD_U2{1},peak);
EP3 = cdf(PD_U3{1},volume);
% 计算拟合边缘分布的单变量重现期
T1=1./(1- EP1);
T2=1./(1- EP2);
T3=1./(1- EP3);

EPdurpeak = [EP1 EP2];EPdurpeak(EPdurpeak == 0) = 1e-4;
EPdurvolume = [EP1 EP3];EPdurvolume(EPdurvolume == 0) = 1e-4;
EPpeakvolume = [EP2 EP3];EPpeakvolume(EPpeakvolume == 0) = 1e-4;
EPdpv= [EP1 EP2 EP3];EPdpv(EPdpv==0)= 1e-4;

[Family1, thetahat1, loglik1]=doublecopulaselect(EPdurpeak, ‘aic’);
[Family2, thetahat2, loglik2]=doublecopulaselect(EPdurvolume, ‘aic’);
[Family3, thetahat3, loglik3]=doublecopulaselect(EPpeakvolume, ‘aic’);

Pdp = copulascdf(EPdurpeak,Family1,thetahat1);
Pdv = copulascdf(EPdurvolume,Family2,thetahat2);
Ppv = copulascdf(EPpeakvolume,Family3,thetahat3);

[A,fam,theta] = mucopulaselect(EPdpv,‘aic’);
%aic 准则进行选择最优copula函数
% A 表示估计的 vinecopula的结构
%fam表示选择的最佳家族
%这里设定一个例子
threecdfvalue = vinecdf(EPdpv,A,fam,theta);
在这里插入图片描述
在这里插入图片描述
如果需要更多关于三维copula的联合重现期、同现重现期和条件重现期等的代码(MATLAB和R语言),(代码可直接运行,包含注释和示例数据,替换数据即可运行出结果)如需请关注公众号趣品科研回复 copula

4 参考文献

He, K., Chen, X., Zhou, J., Zhao, D., & Yu, X. (2024). Compound successive dry-hot and wet extremes in China with global warming and urbanization. Journal of Hydrology, 636, 131332.
He, K., Chen, X., Yu, X., Dong, C., Zhao, D., 2024. Evaluation and prediction of compound geohazards in highly urbanized regions across China’s Greater Bay Area. Journal of Cleaner Production, 449: 141641. DOI:https://doi.org/10.1016/j.jclepro.2024.141641
Jiang, T., Su, X., Zhang, G., Zhang, T., & Wu, H. (2023). Estimating propagation probability from meteorological to ecological droughts using a hybrid machine learning copula method. Hydrol. Earth Syst. Sci., 27, 559-576
Tao, Y., Wang, Y., Wang, D., Ni, L., & Wu, J. (2021). A C-vine copula framework to predict daily water temperature in the Yangtze River. Journal of Hydrology, 598, 126430.

评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值
>