基于Copula函数的沿海城市雨潮复合灾害风险研究

许瀚卿，谭金凯，李梦雅，王 军

（1.华东师范大学地理信息科学教育部重点实验室，上海 200241；2.华东师范大学地理科学学院，上海 200241；3.华东师范大学上海城市公共安全研究中心，上海 200241；4.中山大学大气科学学院，广东珠海 519082）

引言

近年来，我国极端自然灾害频发，社会经济影响严重，已成为制约社会经济可持续发展的一个重要因素。沿海地区人口密集，经济发达，是全球城市群主要分布地带［1］。面对气候变化带来的极端海洋灾害威胁，沿海地区必须提高现有规划设计标准来应对日益加剧的灾害风险［2］。沿海城市极有可能面临暴雨、天文高潮位和风暴增水等多致灾因子相互作用造成的复合洪涝灾害威胁，其中短时强降水和风暴增水所带来的组合风险往往会造成严重损失，非常值得研究关注［3，28］。准确、客观、高效地评估多致灾因子协同作用下的灾害发生概率和风险，对有效防范极端灾害风险具有重要意义。

国内外已有学者运用统计学和数值模拟方法开展了复合洪涝灾害危险性研究，并提出复合洪涝灾害应对策略和措施［3-7］。然而，在多数灾害风险评估和城市规划中，极端和多致灾因子的复合危险性往往被忽略［8，9］。传统的单致灾因子识别和评估分析已较为成熟，但复合洪涝灾害的发生是由多致灾因子复合作用决定的，受多个并发或连续致灾因子影响，从而大大增加了其不确定性［10，31］。王璐阳等［5］提出大气-海洋-陆地相耦合的一体化风暴洪水淹没模拟方法。宋城城等综合考虑海平面上升、陆域和海域地形变化、海塘沉降等因素，构建了12种复合灾害情景，利用MIKE21模拟复合洪涝灾害对上海的影响［11，12］。多致灾因子情景设计是准确开展复合洪涝灾害危险性评估的重要基础，开展强降水和风暴增水联合概率研究，设计雨潮联合分布函数，对于提升危险性模拟精度和应对措施的有效性、减少城市复合洪涝灾害损失具有重要意义［3］。

近年来，Copula函数因其形式多样和使用灵活而被广泛应用到水文学、气候学和洪水风险分析等方面，在极端事件的联合概率分析中发挥了重要作用［13］。Copula函数不限制边缘分布类型，容易扩展到多维，能够灵活构造多维联合分布，在风雨、雨潮、雨洪和洪潮遭遇分析中具有非常大的应用潜力［14，15］。范嘉炜等［16］基于Copula函数分析了潖江河大庙峡流域洪峰流量与洪水历时的联合频率分布特征。黄锦林等［17］对降雨和潮位进行关联性分析，发现适当提高设计标准能够有效降低雨潮风险概率。陈立华等［18］以钦州市为例，利用Copula函数研究历史台风条件下风雨组合规律。陈浩等［19，20］定量评估了深圳河流域不同重现期下雨潮组合的风险率。许红师等［21］发现以单变量作为设计依据会低估具有一定严重程度的台风灾害发生频次，多致灾因子联合重现期的计算结果更加贴近实际。王一新等［22］采用Copula函数建立了太湖超汛限水位与台风发生时间的联合分布，构建了联合概率模型，并量化了太湖超汛限水位变化规律。尽管已有许多研究利用Copula函数开展了多致灾因子的联合概率研究，但现有研究大多基于“季节”或“年度”时间尺度，多采用年极值构建联合分布函数，缺乏在事件尺度量化复合洪涝灾害风险的研究。

由于单致灾因子和年极值联合分布的设计方法会对雨潮设计标准存在高估或低估的情况［17，23］，因此，本文基于二元Copula函数，以上海市1979-2014（36a）实测降水和吴淞口风暴增水资料，筛选极值累积降水量和风暴增水，通过边缘函数优选和Copula函数拟合优度检验，构建上海市年最大日降水量与风暴增水组合风险分析模型，定量评估不同联合重现期下极端降水和风暴潮的遭遇风险概率，在此基础上，计算上海市降水和风暴增水的工程设计值，这对研究区工程防范极端复合灾害具有重要实践意义。

1 数据与方法

1.1 研究区概况

上海地处长江入海口，且地势低平，对海平面上升极为敏感（图1）。海平面上升和地面沉降带来的相对海平面上升会放大沿海洪涝灾害风险［4］。此外，上海易受风暴潮及暴雨内涝侵袭，强降雨和外江高潮位直接影响着上海城区遭受洪涝灾害的程度［24］。2013年“菲特”台风期间发生了自1949年以来首次“风、暴、潮、洪”四碰头事件，在“菲特”台风和冷空气共同影响下，上海和周边地区普降暴雨和特大暴雨，上游洪水下泄量大，又正逢天文高潮位，在上海市造成了十分严重的洪涝灾害。

1.2 数据

全球潮汐和浪涌再分析数据集（Global Tide and Surge Re⁃analysis，GTSR）是第一个基于水动力模型对风暴增水和极端海平面进行模拟的全球再分析产品，该数据集提供了基于1979-2014年期间的极端海平面的估计值［25，26］。日累积降水数据来源于中国气象局，选取的时间为1979-2014年（36a）。降水观测的时间序列通常比风暴增水数据的时间序列更长，更完整［29］。因此，风暴增水数据的可用性决定了两套数据集重叠部分的长度。本研究中，徐家汇站1d（天）最大降水取年最大值代表上海市局地暴雨影响，吴淞口风暴增水取1d最大降水量当天和后3d中的最大值，代表风暴增水的影响。通过以上极值，构成上海市降水和风暴增水数据集，分别表征影响上海市的暴雨、风暴增水事件，从而进行雨潮复合灾害的风险分析。

1.3 研究方法

1.3.1 二元Copula函数及边缘分布函数

Sklar定理是Copula理论和应用的基础［27］。Sklar提出，将一个联合分布分解为n个边缘分布和一个Cop⁃ula函数，Copula函数描述了这n个变量间的相关结构。由此可见，Copula函数实际上是一类将变量的边缘分布函数和它们的联合分布连接在一起的函数，也称为联结函数［30］，能够根据实际情况准确计算多致灾因子遭遇组合下的风险概率，被广泛用于降水和水文事件联合概率计算中。其中，极值分析中常用的二元Copula函数有Gumbel、Clayton和Frank Copula函数，其表达形式如表1所示。

1.3.2 统计检验

Copula函数模型的检验和评价包括边缘分布模型检验和Copula函数的拟合优度评价。其中边缘分布检验主要用于评价所选模型对变量的实际拟合效果，这是构建Copula函数的关键。Kolmogorov⁃Smirnov（K⁃S）检验是Copula函数参数估计中使用最为广泛的方法之一，其主要用于估计变量的边际分布参数。为了检验各分布函数拟合的有效性，可构造赤池信息准则（Akaike Information Criteria，AIC）和贝叶斯信息准则（Bayesian Information Criteria，BIC）统计量从而实现对Copula函数拟合优度的评估。其中，各检验的K⁃S，AIC和BIC的值越小，则表明拟合效果越好。

1.3.3 雨潮组合同现和联合概率模型

风险率模型基于Copula函数，能够定量直观评估多致灾因子大于特定阈值的分布概率，或单一致灾因子大于特定阈值的分布概率［19］。本文基于上海市1979-2014（36a）年最大日降水量与吴淞口站相应风暴增水，计算同现重现期和联合重现期的风险率，用于定量评估上海市雨潮遭遇复合灾害风险。联合概率指的是包含多个条件且所有条件同时成立的概率。边缘概率是与联合概率对应，即仅与单个随机变量有关的概率。

根据Copula函数的定义，二维联合分布的公式为：

同现重现期风险率为：

通过求解同现重现期风险率的倒数即得到“同现重现期”：

联合重现期风险率为：

通过求解联合重现期风险率的倒数即得到“联合重现期”：

式（5）中，降水的边缘分布函数为u，风暴增水的边缘分布函数为v，降水和风暴增水的联合分布Copula函数为C(u,v)。

1.3.4 雨潮组合设计值

风险率无法直接作为实际工程设防标准的参考值，现探讨基于重现期的雨潮组合设计值。对于极端降水、风暴增水复合事件，在给定联合重现期条件下，设计一系列的(x d,y d)组合使得P{x＞x d,x＞y d}最大化，从而得到最优的组合设计值，其计算如下：

式（6）中α为与(x d,y d)相对应的概率临界面，即要求雨潮组合设计值(x d,y d)对应的边缘概率(u d,v d)均在概率临界面上。

2 结果分析与讨论

2.1 雨潮联合分布分析

2.1.1 最优边缘分布选择

本文采用Burr，GEV，Gamma，Weibull和Lognormal边缘分布函数拟合上海市雨潮观测数据。为了选取降水和风暴增水的最优边缘分布，本研究进一步采用K⁃S、AIC和BIC准则分别对其进行拟合优度统计，计算结果见表2。图2为上海市1979-2014（36a）最大日降水量与吴淞口风暴增水的边缘分布函数对比。

图2 降水和风暴增水边缘分布函数Fig.2 Marginal distribution functions of precipitation and storm surge

图2结果表明，Burr、Gamma、GEV和Lognormal分布对于降水具有较好的拟合效果，而Burr和GEV分布对于风暴增水具有较好的拟合效果。根据表2的降水拟合优度统计量来看，GEV分布的K⁃S统计量最小，Lognormal分布的AIC和BIC具有较好的拟合优度。因此，本文以GEV分布作为降水的边缘分布函数。对于风暴增水而言，GEV和Burr在K⁃S统计量中是一致的，均为0.101 0，且GEV和Burr的AIC和BIC两个优度检验相差不大，综上考虑，本文选择GEV分布作为的边缘分布函数。

表2 降水和风暴增水拟合优度统计量Table 2 Goodness⁃of⁃fit statistics of precipitation and storm surge

根据上海市36a最大日降水量与相应风暴增水的最优边缘分布函数，估算了单致灾因子不同重现期下的降水量和风暴增水（图3）。从图3可以看出，降水的数据拟合情况较好，只有一次极端降水事件的拟合效果一般。对于风暴增水来说，部分风暴增水数据拟合效果一般，但总体能够反映数据的边缘分布情况。

图3 降水与风暴增水理论重现期与经验重现期Fig.3 Theoretical and empirical return periods of precipitation and stormsurge

上海市单致灾因子在5、10、20、50和100a重现期下降水和风暴增水的估算值如表3所示。当重现期为10a时，降水量为149.32 mm，相应的风暴增水为0.90 m；当重现期为50a时，降水量为236.10 mm，相应的风暴增水为1.77 m；当重现期为100a时，降水量为287.57 mm，相应风暴增水为2.39 m。

表3 上海市不同重现期降水量及相应风暴增水Table 3 Precipitation in different return periods and corresponding storm surge in Shanghai

2.1.2 最优Copula函数选择

在优选了降水量和风暴增水的最优边缘分布函数后，本文分别采用Gumbel Copula、Clayton Copula和Frank Copula函数推求降水量和风暴增水的联合分布。在推求出各个Copula函数的参数值后，将拟合后的联合分布与经验分布进行比较，从而计算出3组Copula函数的K⁃S检验、AIC准则和BIC。如表4所示，Gum⁃bel、Clayton和FrankCopula函数的τ相关系数分别为0.299 5、0.277 0、0.290 8，表明两者呈现一定正相关。通过对比3组Copula函数的3项检验统计量，发现Frank Copula函数K⁃S检验统计量小于Gumbel和Clayton Copula函数，表明Frank Copula函数更能够表达上海市雨潮灾害事件联合分布的典型特征（图4）。当优选出最优Copula函数之后，就可以得出降水量和风暴增水的联合分布函数。上海市雨潮遭遇最优联合分布为Frank Copula，参数值θ为2.814 9，K⁃S、AIC和BIC分别为3.128 6、-126.420 7和-121.253 7。

表4 Copula函数参数及检验统计量Table 4 Copula function parameters and test statistics

图4 降水与风暴增水联合概率分布图Fig.4 Joint probability distribution map of precipitation and storm surge

2.2 不同联合重现期下雨潮组合的同现及联合概率

最优边缘分布函数能够较好表征上海市最大日降水量与吴淞口风暴增水典型复合灾害事件，其中降水量和风暴增水的边缘分布均为GEV分布，其联合分布函数为Frank Copula函数。在此基础上，进一步计算上海市不同联合重现期下降水和风暴增水同现概率。如表5所示，在5、10、20、50和100a雨潮联合重现期下，雨潮遭遇的同现概率分别为0.078 0、0.023 5、0.001 1和0.000 3。

表5 不同重现期下降水和风暴增水同现概率Table 5 Co⁃occurrence probability of precipitation and storm surge in different return periods

如表6所示，在5、10、20、50和100a雨潮联合重现期下，雨潮遭遇的联合概率分别为：0.321 9、0.176 4、0.093 4、0.038 8和0.019 7。

表6 不同重现期下降水和风暴增水联合概率Table 6 Joint probability of precipitation and storm surgein different return periods

雨潮主要表现为上海在遭遇台风影响时，发生城市区域强降水，大风导致风暴增水。从表5和6可以看出，不同联合重现期下的雨潮同现概率远小于联合概率。通过对比可以发现，在5、10、20、50和100a雨潮联合重现期下，联合概率是同现概率的4.12、7.51、14.21、34.27、67.72倍。通过对比可以发现，不通联合重现期下的雨潮同现概率远小于联合概率。虽然同现概率小于联合概率，但是沿海城市遭遇暴雨和高潮位时所造成的城市洪涝灾害更为严重，经济损失和人员伤亡更为巨大。因此，雨潮碰头事件非常值得上海市防汛关注。

2.3 工程雨潮组合设计值

在实际的设防标准中，可根据雨潮重现期制定相应的设防标准参考。表7为在不同联合重现期下强降水与风暴增水不同组合的设计值。可以看出，Gumbel函数在5a至100a的重现期下其对应的降水与风暴增水均小于Clayton函数与Frank函数，尤其是在50a与100a重现期下更明显；总体上，Clayton函数与Frank函数对应的雨潮组合设计值相差不大，其中在50a重现期下，降水与风暴增水约为233 mm和2.7 m；在100a联合重现期下，降水与风暴增水约为276 mm和3.5 m，说明了上海市沿海地区要预防百年一遇的强降水或风暴增水，须在考虑天文大潮的基础上，设计构筑至少3.5 m的防汛墙。

表7 不同重现期下的雨潮组合工程设计值Table 7 Design values of combined rainfall and surge under different return periods

3 结论

在全球气候变暖与人类活动的双重影响下，沿海城市面临短时强降水和下游风暴增水顶托的共同作用，极易引发局部的城市内涝或流域性的复合洪涝灾害事件，进而造成巨大的社会经济损失。本研究利用上海市36a最大日降水量与吴淞口风暴增水数据，优选边缘分布函数，采用二元Copula函数，计算了上海地区降水和风暴增水的不同联合重现期下的同现概率和联合概率，并从工程角度，给出了不同联合重现期下的设计降水量和风暴增水。此研究表明二元Copula函数能够较为准确地计算出不同联合重现期下的设计降水和风暴增水，为优化防汛工程的设计标准，设计科学复合洪涝灾害情景提供了可能。根据分析结果，得到以下结论：

（1）基于K⁃S检验、AIC准则和BIC准则计算其检验统计量，分别用Gamma、Lognormal、GEV、Weibull和Burr分布对降水和风暴增水进行拟合，发现上海地区36a最大日降水量和吴淞口相应风暴增水均适合运用GEV分布来进行拟合。通过使用二元Copula函数构造了上海地区降水量和风暴增水的联合分布，发现Frank Copula函数能够较好地描述上海地区典型雨潮遭遇事件的联合分布。

（2）基于最优边缘分布函数GEV和Frank Copula联合分布函数，定量评估了上海地区最大日降水量和风暴增水在不同联合重现期下的同现概率和联合概率，结果发现：上海地区在5、10、20、50和100a雨潮联合重现期下，联合概率是同现概率的4.12、7.51、14.21、34.27、67.72倍。

（3）传统的单变量边缘分布，会对降水和风暴增水存在高估或低估的情况。通过Copula函数选择联合重现期，能够更加准确计算得到联合重现期下的工程设计值。在100a联合重现期下，降水与风暴增水约为276 mm和3.5 m，说明了上海市要预防百年一遇的雨潮复合洪涝灾害，须在考虑天文大潮的基础上，额外设计构筑至少3.5 m的防汛墙。

本文在重现期分析结果中，利用多种精度评价方法进行拟合，降低了边缘分布和Copula函数拟合的误差。但由于风暴增水数据来源于全球风暴潮模型模拟结果（即利用气候模式驱动模型生成的逐日风暴增水数据），模式模拟数据和实际观测还是存在一定误差，可利用区域水动力模型进一步优化。