CN117634325B - 资料受限河口区域极值事件识别和复合洪水灾害分析方法和系统 - Google Patents
资料受限河口区域极值事件识别和复合洪水灾害分析方法和系统 Download PDFInfo
- Publication number
- CN117634325B CN117634325B CN202410110058.5A CN202410110058A CN117634325B CN 117634325 B CN117634325 B CN 117634325B CN 202410110058 A CN202410110058 A CN 202410110058A CN 117634325 B CN117634325 B CN 117634325B
- Authority
- CN
- China
- Prior art keywords
- typhoon
- data
- flood
- estuary
- typhoons
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 116
- 239000002131 composite material Substances 0.000 title claims abstract description 87
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 claims abstract description 106
- 238000004458 analytical method Methods 0.000 claims abstract description 50
- 150000001875 compounds Chemical class 0.000 claims abstract description 47
- 238000004088 simulation Methods 0.000 claims abstract description 29
- 238000011144 upstream manufacturing Methods 0.000 claims abstract description 20
- 238000012216 screening Methods 0.000 claims abstract description 6
- 238000012502 risk assessment Methods 0.000 claims abstract description 4
- 241000039077 Copula Species 0.000 claims description 53
- 238000010801 machine learning Methods 0.000 claims description 34
- 238000004364 calculation method Methods 0.000 claims description 24
- 238000012937 correction Methods 0.000 claims description 23
- 241000345998 Calamus manan Species 0.000 claims description 16
- 235000012950 rattan cane Nutrition 0.000 claims description 16
- 230000006870 function Effects 0.000 claims description 14
- 238000007781 pre-processing Methods 0.000 claims description 14
- 238000004140 cleaning Methods 0.000 claims description 11
- 230000002159 abnormal effect Effects 0.000 claims description 10
- 238000012549 training Methods 0.000 claims description 10
- 238000007405 data analysis Methods 0.000 claims description 8
- 230000007774 longterm Effects 0.000 claims description 8
- 238000005457 optimization Methods 0.000 claims description 7
- 238000009499 grossing Methods 0.000 claims description 6
- 238000012360 testing method Methods 0.000 claims description 6
- 230000009467 reduction Effects 0.000 claims description 3
- 230000000875 corresponding effect Effects 0.000 description 17
- 238000011160 research Methods 0.000 description 13
- 230000008859 change Effects 0.000 description 9
- 238000007637 random forest analysis Methods 0.000 description 9
- 230000000694 effects Effects 0.000 description 8
- 238000005516 engineering process Methods 0.000 description 7
- 230000035945 sensitivity Effects 0.000 description 7
- 230000008569 process Effects 0.000 description 6
- 238000000354 decomposition reaction Methods 0.000 description 4
- 238000005259 measurement Methods 0.000 description 4
- 230000003416 augmentation Effects 0.000 description 3
- 230000001186 cumulative effect Effects 0.000 description 3
- 230000001419 dependent effect Effects 0.000 description 3
- 238000001914 filtration Methods 0.000 description 3
- 238000011835 investigation Methods 0.000 description 3
- 230000009022 nonlinear effect Effects 0.000 description 3
- 230000000737 periodic effect Effects 0.000 description 3
- 238000010206 sensitivity analysis Methods 0.000 description 3
- 238000007619 statistical method Methods 0.000 description 3
- 230000009466 transformation Effects 0.000 description 3
- 238000012795 verification Methods 0.000 description 3
- 238000006424 Flood reaction Methods 0.000 description 2
- 241000282414 Homo sapiens Species 0.000 description 2
- 230000008901 benefit Effects 0.000 description 2
- 238000013145 classification model Methods 0.000 description 2
- 238000007621 cluster analysis Methods 0.000 description 2
- 238000010219 correlation analysis Methods 0.000 description 2
- 238000002790 cross-validation Methods 0.000 description 2
- 230000003247 decreasing effect Effects 0.000 description 2
- 238000011161 development Methods 0.000 description 2
- 230000009021 linear effect Effects 0.000 description 2
- 230000007246 mechanism Effects 0.000 description 2
- 238000012544 monitoring process Methods 0.000 description 2
- 238000000513 principal component analysis Methods 0.000 description 2
- 238000013179 statistical model Methods 0.000 description 2
- 238000012300 Sequence Analysis Methods 0.000 description 1
- 230000009471 action Effects 0.000 description 1
- 238000013528 artificial neural network Methods 0.000 description 1
- 238000010224 classification analysis Methods 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 230000002596 correlated effect Effects 0.000 description 1
- 230000008878 coupling Effects 0.000 description 1
- 238000010168 coupling process Methods 0.000 description 1
- 238000005859 coupling reaction Methods 0.000 description 1
- 238000003066 decision tree Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 230000003993 interaction Effects 0.000 description 1
- 238000005065 mining Methods 0.000 description 1
- YHXISWVBGDMDLQ-UHFFFAOYSA-N moclobemide Chemical compound C1=CC(Cl)=CC=C1C(=O)NCCN1CCOCC1 YHXISWVBGDMDLQ-UHFFFAOYSA-N 0.000 description 1
- 230000002265 prevention Effects 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 238000000611 regression analysis Methods 0.000 description 1
- 230000002787 reinforcement Effects 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 230000007704 transition Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/27—Design optimisation, verification or simulation using machine learning, e.g. artificial intelligence, neural networks, support vector machines [SVM] or training a model
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
- G06Q50/00—Information and communication technology [ICT] specially adapted for implementation of business processes of specific business sectors, e.g. utilities or tourism
- G06Q50/10—Services
- G06Q50/26—Government or public services
- G06Q50/265—Personal security, identity or safety
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/02—Reliability analysis or reliability optimisation; Failure analysis, e.g. worst case scenario performance, failure mode and effects analysis [FMEA]
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A10/00—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE at coastal zones; at river basins
- Y02A10/40—Controlling or monitoring, e.g. of flood or hurricane; Forecasting, e.g. risk assessment or mapping
Landscapes
- Engineering & Computer Science (AREA)
- Business, Economics & Management (AREA)
- Theoretical Computer Science (AREA)
- Physics & Mathematics (AREA)
- Evolutionary Computation (AREA)
- General Physics & Mathematics (AREA)
- Tourism & Hospitality (AREA)
- Strategic Management (AREA)
- Development Economics (AREA)
- Human Resources & Organizations (AREA)
- Marketing (AREA)
- Primary Health Care (AREA)
- Economics (AREA)
- Health & Medical Sciences (AREA)
- General Business, Economics & Management (AREA)
- Educational Administration (AREA)
- General Health & Medical Sciences (AREA)
- Artificial Intelligence (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Computer Security & Cryptography (AREA)
- Medical Informatics (AREA)
- Software Systems (AREA)
- Computer Hardware Design (AREA)
- Geometry (AREA)
- General Engineering & Computer Science (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
本发明公开了一种资料受限河口区域极值事件识别和复合洪水灾害分析方法和系统,包括如下步骤:采集并修正台风最大风速,构建适用于台风模拟的驱动风场;确定河口区域的台风影响圆;采用调和分析法基于已有实测潮位数据计算风暴增水;筛选出影响河口区域的所有台风路径,形成台风路径集合;构建河口风暴潮数值模型,计算所有在河口区域产生影响的台风造成的风暴增水,并获得风暴高潮位时间;根据风暴高潮位时间,基于实测流量在台风期间变异性确定上游流量值,采用滑动相关性法确定日台风降雨量及其对应时间;分析河口区域复合洪水概率并输出分析结果。能够更全面和综合地分析复合洪水的成因和影响,提高了复合洪水的风险评估的精度和稳定性。
Description
技术领域
本发明涉及洪水灾害仿真技术,尤其是一种资料受限河口区域极值事件识别和复合洪水灾害分析方法和系统。
背景技术
河口地区是海陆交互作用的复杂区域,受到多种气象、水文、海洋等因素的影响,容易发生洪涝灾害,对人类的生命财产安全和社会经济发展造成严重威胁。在气候变化的背景下,极端天气和气候事件的频率和强度可能增加,导致河口地区的洪涝灾害风险加剧。因此,研究河口地区的洪涝灾害机理和评估方法,提高河口地区的洪涝灾害防御能力,具有重要的理论意义和实际价值。
河口地区的洪涝灾害往往并非由单一因素引起,而是由多种因素相互作用和叠加的结果。这些因素主要包括海洋因素(如天文潮、风暴潮和海浪)、河流因素(如上游河道流量)和大气因素(如降雨造成的城市内涝)。当这些因素同时或相继发生时,可能导致河口地区的水位超过防洪标准,引发洪涝灾害。这种由多种因素组合导致的洪涝灾害称为复合洪水(Compound Flooding,CF)。复合洪水的危害程度往往高于单一因素的洪水,因为复合洪水的发生概率较低,防洪措施往往没有考虑到复合洪水的可能性,而且复合洪水的持续时间和范围较大,造成的损失较难恢复。
目前,对河口地区的复合洪水研究主要有两类方法,一类是基于统计模型的方法,另一类是基于动力数值模型的方法。这两类方法各有优缺点,基于统计模型的方法可以较好地理解复合洪水的时空变化和气候因素的影响,但需要较长的观测数据。基于动力数值模型的方法可以较好地刻画复合洪水的演变过程和物理机制,但需要较高的计算能力和边界条件,而且难以考虑不确定性和误差的影响。
目前,对河口地区的复合洪水研究还存在一些问题和挑战,主要包括数据、模型和方法等多个方面。对于数据问题,复合洪水的研究需要多种类型的数据,如潮位、流量、降雨、风场、气压等,而这些数据往往分布在不同的机构和平台,难以获取和整合。而且,由于观测时间较短或部分站点数据不公开等原因,一些国家和地区的数据长度受限,难以满足复合洪水研究的要求。数据的质量和完整性也影响复合洪水研究的可靠性和有效性。
因此,河口地区的复合洪水研究是一个重要而复杂的课题,需要在数据、模型和方法等方面进行深入的探索和创新。
发明内容
发明目的,提供一种资料受限河口区域极值事件识别和复合洪水灾害分析方法和系统,以解决现有技术存在的上述问题。
技术方案,根据本申请的一个方面,提供一种资料受限河口区域极值事件识别和复合洪水灾害分析方法,包括如下步骤:
步骤S1、采集预定区域的台风历史数据并预处理和修正台风最大风速,构建适用于台风模拟的驱动风场;
步骤S2、确定河口区域的台风影响圆;采用调和分析法基于已有实测潮位数据计算风暴增水;
步骤S3、基于台风影响圆确定的范围,筛选出影响河口区域的所有台风路径,形成台风路径集合;
步骤S4、构建河口风暴潮数值模型,进行数值模拟,计算所有在河口区域产生影响的台风造成的风暴增水,并获得风暴高潮位时间;
步骤S5、根据风暴高潮位时间,基于实测流量在台风期间变异性确定上游流量值,调取中国地面气候资料日值数据集,采用滑动相关性法确定日台风降雨量及其对应时间;
步骤S6、基于步骤S3至S5的输出结果,采用R-vine Copula对河口区域复合洪水概率进行分析并输出分析结果。
根据本申请的一个方面,所述步骤S1进一步为:
步骤S11、从数据网站获取预定区域的台风历史数据,并进行清洗和预处理,剔除异常值和缺失值,统一数据格式和单位,台风历史数据包括每一台风的编号、名称、时间、位置、最大风速和最低气压;
步骤S12、从预定数据库获取再分析数据集,提取与台风历史数据对应的风场和气压场数据,并进行插值和平滑;
步骤S13、采用公式法和预训练的Light GBM机器学习模型对台风最大风速进行修正,获得修正后的台风最大风速;其中公式法包括如下修正公式:
Vm1=a(φ,4971)(1010-P0)b(φ,4971);Vm2=a(φ,7220)(1010-P0)b(φ,7220);
Vm'=Vm-(Vm1-Vm2);其中,a(φ,4971)、b(φ,4971)、a(φ,7220)、b(φ,7220)分别为1949~1971年和1972~2020年时段内采用纬度为φ的台风数据拟合得到的参数;Vm和P0分别为待订正的风速和对应气压;Vm'为修正后的风速;
步骤S14、将修正后的台风最大风速与再分析数据集中的风场和气压场数据结合,构建台风模拟的驱动风场。
根据本申请的一个方面,所述步骤S2进一步为:
步骤S21、针对河口区域的特点进行概化,获取圆心坐标位置,并以N倍的台风最大风速半径为区域半径,获得台风影响圆,获取台风影响圆圆周上若干点的坐标并存储;N为正整数;
步骤S22、获取河口区域的河口潮数据并利用调和分析法对河口潮的实测潮位数据进行分解,得到天文潮和每个潮周期内的偏峰增水;
步骤S23、修正偏峰增水并计算风暴增水以及风暴增水的极值和频率分布。
根据本申请的一个方面,所述步骤S3进一步为:
步骤S31、构造机器学习模型的输入数据特征,包括台风影响圆圆心、台风距离特征、台风中心最低气压及最大风速、台风最大风速半径和台风路径方位特征;台风距离特征是每场台风路径与台风影响圆圆心的最小距离;
步骤S32、构建机器学习模型并采用输入数据进行训练和测试,同时对机器学习模型的参数进行调优;
步骤S33、采用训练完成后的机器学习模型对台风进行识别,获得所有影响河口区域的台风路径,形成台风路径集合。
根据本申请的一个方面,所述步骤S4进一步为:
步骤S41、获取再分析风场并构建影响河口区域的台风风暴潮的驱动风场集;
步骤S42、以驱动风场集作为驱动,采用三维有限体积法FVCOM进行批量数值模拟,计算近M年所有在河口区域产生影响的台风造成的风暴增水;M为正整数;
步骤S43、提取计算结果中的增水、天文高潮位和风暴高潮位数据,包括场次、时间和高潮位。
根据本申请的一个方面,所述步骤S5进一步为:
步骤S51、获取河口上游预定地点的实测流量数据,包括逐日流量值和极值流量值;
步骤S52、获取并从中国地面气候资料日值数据集中提取河口区域附近地区的日降雨量数据,包括逐日降雨量值和极值降雨量值;
步骤S53、根据得到的风暴高潮位时间,基于实测流量在台风期间的变异性确定上游流量值,作为极端复合洪水的一个灾害因子;基于中国地面气候资料日值数据集,采用滑动相关性法,确定日台风降雨量及其对应时间,作为极端复合洪水的另一个灾害因子;基于修正后的台风最大风速和计算得到的风暴增水,作为极端复合洪水的另外两个灾害因子;利用调和分析得到的短期天文潮,作为极端复合洪水的最后一个灾害因子;
步骤S54、将得到的上游径流、降雨、天文潮、风暴潮和增水数据集按照时间对齐,构建河口区域极端复合洪水灾害因子数据集。
根据本申请的一个方面,所述步骤S6进一步为:
步骤S61、构建藤Copula函数集合并筛选,所述藤Copula函数至少包括C藤、D藤和R藤;
步骤S62、读取复合洪水的时间序列数据,根据台风编号和台风名称进行分组,得到一个按台风分组的列表,每个分组包含一个台风在河口区域各个潮位站造成的复合洪水的时间序列;
步骤S63、利用R-vine Copula模型对每个台风分组进行数据分析,建立复合洪水多灾害因子之间的联合概率分布模型;利用AIC和BIC指标,对R-vine Copula模型进行选择和验证,选择最优的模型参数和结构;
步骤S64、利用参数优化后的R-vine Copula模型,进行复合洪水的概率分析和风险评估,计算不同复合洪水等级的发生概率和超越概率,分析复合洪水的风险等级和影响范围。
根据本申请的一个方面,所述步骤S13还包括:
步骤S13a、分别采用公式法和预训练的Light GBM机器学习模型对台风最大风速进行修正,获得修正后的台风最大风速;
步骤S13b、构建风速窗口,并分别计算每个风速窗口内公式法和Light GBM机器学习模型计算得到的台风最大风速的均方根误差;
步骤S13c、根据均方根误差,确定每个风速窗口的修正模型,获得各个风速窗口的修正数据,形成复合修正后的台风最大风速。
根据本申请的一个方面,所述步骤S64还包括:
步骤S641、根据防洪标准和历史洪水资料,确定不同复合洪水等级对应的水位高度和流量大小;
步骤S642、利用参数优化后的R-vine Copula模型,对每个台风分组的复合洪水数据进行模拟,得到不同复合洪水等级的发生概率和超越概率,以及其置信区间;
步骤S643、根据复合洪水等级的发生概率和超越概率,评估不同台风分组的复合洪水风险等级,分成若干个等级,根据复合洪水的水位高度和流量大小,分析不同台风分组的复合洪水影响范围;
步骤S644、根据复合洪水的风险等级和影响范围,制定长期的防洪规划和减灾策略。
根据本申请的另一个方面,提供一种资料受限河口区域极值事件识别和复合洪水灾害分析系统,其特征在于,包括:
至少一个处理器;以及
与至少一个所述处理器通信连接的存储器;其中,
所述存储器存储有可被所述处理器执行的指令,所述指令用于被所述处理器执行以实现上述任一项技术方案所述的资料受限河口区域极值事件识别和复合洪水灾害分析方法。
有益效果,本发明针对数据长度受限的问题,提出了一种极值事件识别和复合洪水灾害研究的新方法,旨在为河口地区的复合洪水研究提供一种新的思路和技术支撑。相关技术优势,将在具体实施方式进行详细说明。
附图说明
图1是本发明的流程图。
图2是本发明步骤S1的流程图。
图3是本发明步骤S2的流程图。
图4是本发明步骤S3的流程图。
图5是本发明步骤S4的流程图。
图6是本发明步骤S5的流程图。
图7是本发明步骤S6的流程图。
图8是不同相位差下日累计降雨量与最高风暴潮位相关性曲线图。
图9是各节点密度分布结构图。
具体实施方式
如图1所示,根据本申请的一个方面,提供一种资料受限河口区域极值事件识别和复合洪水灾害分析方法,包括如下步骤:
步骤S1、采集预定区域的台风历史数据并预处理和修正台风最大风速,构建适用于台风模拟的驱动风场;
步骤S2、确定河口区域的台风影响圆;采用调和分析法基于已有实测潮位数据计算风暴增水;
步骤S3、基于台风影响圆确定的范围,筛选出影响河口区域的所有台风路径,形成台风路径集合;
步骤S4、构建河口风暴潮数值模型,进行数值模拟,计算所有在河口区域产生影响的台风造成的风暴增水,并获得风暴高潮位时间;
步骤S5、根据风暴高潮位时间,基于实测流量在台风期间变异性确定上游流量值,调取中国地面气候资料日值数据集,采用滑动相关性法确定日台风降雨量及其对应时间;
步骤S6、基于步骤S3至S5的输出结果,采用R-vine Copula对河口区域复合洪水概率进行分析并输出分析结果。
在本实施例中,为了分析河口区域极端复合洪水和灾害因子之间的依赖性,首先提出了基于机器学习Light GBM模型和经验公式的复合法对历史台风最大风速进行校正,采用随机森林法对可能在河口区域产生增水的台风进行识别;其次开发了RHSS工具箱对近60余年河口区域历史风暴潮的批量数值模拟,结合实测降雨和径流资料构建了河口区域极端复合洪水灾害因子基础数据集;最后引入R-vine Copula函数构建了河口区域多灾害因子极端复合洪水模型,对多灾害因子依赖结构和敏感因子进行了分析和识别。
在本实施例中,利用公式法和Light GBM机器学习模型对台风最大风速进行修正,消除历史台风风速数据的非一致性现象,提高风暴潮的模拟精度和实时性。利用机器学习模型对台风路径进行识别,获得所有影响河口区域的台风路径集合,提高台风路径的筛选效率和准确性。利用三维有限体积法FVCOM进行风暴潮的数值模拟,计算所有在河口区域产生影响的台风造成的风暴增水,并获得风暴高潮位时间,提高风暴潮的模拟精度和稳定性。利用滑动相关性法确定日台风降雨量及其对应时间,反映复合洪水事件的相关性和时空变化,提高复合洪水的数据输入的灵敏度和可靠性。利用R-vine Copula对河口区域复合洪水概率进行分析,描述多变量之间的复杂依赖关系,提高复合洪水的概率分析的精度和稳定性。
如图2所示,根据本申请的一个方面,所述步骤S1进一步为:
步骤S11、从数据网站获取预定区域的台风历史数据,并进行清洗和预处理,剔除异常值和缺失值,统一数据格式和单位,台风历史数据包括每一台风的编号、名称、时间、位置、最大风速和最低气压;
步骤S12、从预定数据库获取再分析数据集,提取与台风历史数据对应的风场和气压场数据,并进行插值和平滑;
步骤S13、采用公式法和预训练的Light GBM机器学习模型对台风最大风速进行修正,获得修正后的台风最大风速;其中公式法包括如下修正公式:
Vm1=a(φ,4971)(1010-P0)b(φ,4971);Vm2=a(φ,7220)(1010-P0)b(φ,7220);
Vm'=Vm-(Vm1-Vm2);其中,a(φ,4971)、b(φ,4971)、a(φ,7220)、b(φ,7220)分别为1949~1971年和1972~2020年时段内采用纬度为φ的台风数据拟合得到的参数;Vm和P0分别为待订正的风速和对应气压;Vm'为修正后的风速;
步骤S14、将修正后的台风最大风速与再分析数据集中的风场和气压场数据结合,构建台风模拟的驱动风场。
台风中心最低气压和最大风速半径处风速是台风强度的两个重要指标,且风速和气压的变式具有较强的线性关系,尤其是风速在50m/s以内时。20世纪70年代以前,主要通过航空器飞行观测和海洋测量船监测的方法获取海面台风数据。随着气象卫星的发射和数据图像处理技术的发展,风廓线雷达的使用,台风监测数据的数量和精度都得到较大的提高。观测方法的改进和数据整编方式的不同致使台风风压历史数据出现非一致性现象。为了解决这一问题,提供了新的修正公式,具体如上所示。
该方法的主要思路是分别求出1972~2020年和1941~1971年两时段内的参数,将气压代入,以上述两式结果的差值作为待修正风速和真值的差值,从而消除台风早期风速数据在不同年和不同纬度上的误差。以a0=0.6,b0=0.6为初值,采用最小二乘法对上式进行求解。经研究发现,当纬度从10~40°N增大时,系数a先变小后变大,并在25°N左右达到最小;指数b的变化规律恰恰相反,并在25°N左右达到最大。台风中心最低气压是台风最大风速的主导因素,纬度是最大风速的次要因素。
Light GBM模型采用1972~2007年台风路径数据为训练集,2008~2020年数据为测试集。在经过RSCV(Randomized Search CV)和GSCV(Grid Search CV)方法组合对参数进行优化和KFold交叉验证后,Light GBM模型最优参数确定为:最大深度max_depth=13,叶节点数量num_leaves=18,迭代次数num_iterations=700,学习率learning_rate=0.01。其余参数主要维护模型稳定性,保持默认。
在研究中发现,在一般风速情景下,Light GBM模型预测结果和公式订正结果的效果相近;在高风速情景下,Light GBM模型预测结果更好,更接近真实值。Light GBM模型仅在样本数为340~350之间出现轻微过拟合现象,预测风速偏大。对于整个测试集的验证结果,Light GBM模型的均方根误差为1.86m/s,公式法的均方根误差为1.97m/s,Light GBM模型的整体误差更小。仅从两种方法的预测风速的值来看,Light GBM模型预测风速略大于公式法,在构造风暴潮计算驱动风场时,风速偏大比风速偏小得到的数值模拟结果更保守。
在本实施例中,采用Light GBM机器学习模型和公式法相结合的方法对台风最大风速进行修正,均方根误差可下降到1.69m/s,相较于“公式法”校正结果误差下降了12.0%,相较于传统Light GBM模型结果预测风速误差减小9.2%。
根据本申请的一个方面,所述步骤S13还包括:
步骤S13a、分别采用公式法和预训练的Light GBM机器学习模型对台风最大风速进行修正,获得修正后的台风最大风速;
步骤S13b、构建风速窗口,并分别计算每个风速窗口内公式法和Light GBM机器学习模型计算得到的台风最大风速的均方根误差;
步骤S13c、根据均方根误差,确定每个风速窗口的修正模型,获得各个风速窗口的修正数据,形成复合修正后的台风最大风速。
在风速较小(<20m/s)或较大(>55m/s)时,公式法的均方根误差较大;风速在20~40m/s区间内时,两者均方根误差相当;风速在41~55m/s区间内时,公式法误差小于LightGBM模型(最大相差约0.5m/s)。因此,采用公式法与Light GBM模型相结合的方法,即当风速在41~55m/s区间内时采用公式法,其余情况采用Light GBM模型,最终均方根误差可下降到1.69m/s。
在本申请的另一实施例中,还包括利用多种数据源和方法对台风最大风速进行修正的步骤:包括利用卫星遥感数据、风速公式和机器学习模型等,综合考虑台风的位置、气压、风场和历史数据等因素,提高台风最大风速的修正效果。
从数据网站获取预定区域的台风历史数据,并进行清洗和预处理,剔除异常值和缺失值,统一数据格式和单位,台风历史数据包括每一台风的编号、名称、时间、位置、最大风速和最低气压;
从预定数据库获取再分析数据集,提取与台风历史数据对应的风场和气压场数据,并进行插值和平滑;
利用卫星遥感数据,如风云卫星、海洋卫星和气象卫星等,获取台风的实时位置、风速和气压等信息,并与台风历史数据进行对比和校验;
采用公式法和预训练的Light GBM机器学习模型对台风最大风速进行修正,获得修正后的台风最大风速;
将修正后的台风最大风速与再分析数据集中的风场和气压场数据结合,构建台风模拟的驱动风场。
在本申请的另一实施例中,还包括利用误差分析和敏感性分析对台风最大风速的修正结果进行评估,计算修正后的台风最大风速与实测风速的均方根误差,分析台风最大风速对台风模拟的驱动风场的影响,选择最优的修正方法和参数。具体步骤如下:
从数据网站获取预定区域的台风实测风速数据,并进行清洗和预处理,剔除异常值和缺失值,统一数据格式和单位,台风实测风速数据包括每一台风的编号、名称、时间、位置和风速;
分别采用公式法和Light GBM机器学习模型对台风最大风速进行修正,获得修正后的台风最大风速;
构建风速窗口,并分别计算每个风速窗口内公式法和Light GBM机器学习模型计算得到的台风最大风速与实测风速的均方根误差(RMSE),作为误差分析的指标;
利用河口风暴潮数值模型FVCOM,分别以公式法和Light GBM机器学习模型计算得到的台风最大风速为驱动风场,进行台风模拟,计算河口区域的风暴增水,并与实测增水进行对比,计算相关系数(R),作为敏感性分析的指标;
根据误差分析和敏感性分析的结果,选择RMSE最小且R最大的修正方法和参数,作为最优的修正方案。
如图3所示,根据本申请的一个方面,所述步骤S2进一步为:
步骤S21、针对河口区域的特点进行概化,获取圆心坐标位置,并以N倍的台风最大风速半径为区域半径,获得台风影响圆,获取台风影响圆圆周上若干点的坐标并存储;N为正整数;
步骤S22、获取河口区域的河口潮数据并利用调和分析法对河口潮的实测潮位数据进行分解,得到天文潮和每个潮周期内的偏峰增水;
步骤S23、修正偏峰增水并计算风暴增水以及风暴增水的极值和频率分布。
台风在平面结构上是一个近似圆形的、影响半径为几百到上千公里的旋转风场。台风气旋的强弱和平面尺度大小,台风路径与研究地点距离的远近等因素都会对增水高度产生影响。台风与研究区域之间的距离是影响风暴增水大小的重要指标。在一定范围内,台风越接近研究区域引起强增水的可能性越大。但在没有潮位资料的前提下,仅根据台风路径等气象信息无法准确判定台风在研究区是否产生增水。
因此,为从数量上缩小台风研究范围,提出台风影响圆的假设。将河口区域被台风可能影响的范围概化为一个以某岛东岸为圆心,10倍台风最大风速半径为半径的圆形区域。当台风路径与长江台风影响圆相交或相切时,则认为该台风可能在河口区域引起风暴增水。
接着,采用调和分析的方法从河口潮中提取天文潮,计算每个潮周期内的偏峰增水。并根据台风发生时间,对偏峰增水不合理的结果进行修正。
在气候变化和海面升高的背景下,从长期角度看调和常数并非常数,而是随着时间的推移也在发生着改变。在河口地区,由于上游径流的非线性作用,一年内甚至一天之间的“平均水位”都可能发生变化。调和分析计算的本质是通过数学方法(如最小二乘法)求解未知参数,“平均水位”的变化导致了振幅和迟角的非平稳性。
在本实施例中,引入了非平稳潮调和分析方法NS_TIDE。NS_TIDE认为每个时刻的平均海面,振幅和迟角都是变动的,换句话说,该方法考虑了上游径流对潮汐的非线性作用,具体如下:
其中,η(tj)为第tj时刻的潮位,η0(tj)为模型中假定的第tj时刻的平均海平面高度,代表NS_TIDE中的亚潮模型的水位波动;∑n k=1(·)为径潮模型,代表了日潮及以上的高频潮族水位波动;QR(tj)为低通滤波后的上游径流值,低通滤波是为了平滑数据,过滤掉下游潮波或上游大坝运行对径流的影响;Rt为下游参考站一个太阴日内的最大潮差。NS_TIDE模型中要求下游参考站受径流影响较小。为了确定潮差极值并消除潮差的不连续性,将小时间隔的潮位数据在通过高通滤波后采用样条曲线插值至6分钟间隔,随后采用步长1h、跨度27h的滑动窗口计算每个小时的最大潮差。(ps,qs,rs)和(pf,qf,rf)为每个站点和频率区间的待定指数。
NS_TIDE认为分潮的频率是已知的并假设相近频率分潮的待定指数相同,因而将这些分潮划分至同一频率区间。tQ和tR分别为上流径流和下游潮波传播到研究站点的时间差,两者都通过相应的互相关系数最大时的时间差确定。Si,ci,j和sk,i(i=0,1,2)为模型中待定系数;Hk(tj)和φk(tj)为tj时刻的振幅和迟角。
本实施例的主要目的是将天文潮从实测潮位中分离出来并计算风暴增水,更关注台风发生时的调和分析效果,而对全年调和分析的整体精度并没有过高要求。影响河口区域的强台风一般发生在夏、秋两季,而每年夏季东南季风带来的太平洋暖湿气流恰好会在长江中下游地区造成长时间、大范围的降雨和洪水。NS_TIDE在洪水期可以得到更准确的水位回报结果,在某站点采用NS_TIDE回报的天文潮计算风暴增水精度更高,结果更为可信。
在本实施例中,基于实测潮位的调和分析结果,采用随机森林法对影响河口区域的台风路径进行了识别,极大程度上减少了后续数值模拟计算组次和工作量,相较于仅基于距离的“最大影响圆法”可节约30%以上的数值模拟计算资源。
通过NS_TIDE非平稳潮调和分析可分离出1988~2020年共33年的天文潮时间序列。计算每个潮周期内最高天文潮位和实测最高潮位的差,再通过历史台风路径信息进行比对和筛选,即可得到台风偏峰增水时间序列。由于亚潮潮族和其他高频潮族在某些情景下并不能做到完全分离,导致调和分析回报的天文潮相位存在微小偏差。由研究结果可知,河口潮汐相位差对风暴增水(非潮残差)的计算结果有较大影响。偏峰增水只计入风暴潮和天文潮在一个潮周期内高潮位差。相比之下,采用偏峰增水作为河口风暴增水的指标,可以很好的避免相位偏移引入的误差的影响。
在本申请的另一实施例中,确定圆心和半径的过程如下:
步骤一、利用历史数据和统计分析确定台风影响圆的范围,根据河口区域历史发生的复合洪水事件,分析与之相关的台风的位置和特征,确定一个合理的台风影响圆的半径和圆心,使其能够覆盖河口区域的主要风暴潮来源。在获取复合洪水和台风数据后,进行相关性分析。
获取复合洪水事件数据。从数据网站获取河口区域历史发生的复合洪水事件的数据,并进行清洗和预处理,剔除异常值和缺失值,统一数据格式和单位,复合洪水事件的数据包括每一次事件的时间、地点、水位、流量、降雨等;
获取相关台风数据。从数据网站获取与复合洪水事件相关的台风的数据,并进行清洗和预处理,剔除异常值和缺失值,统一数据格式和单位,台风的数据包括每一场台风的编号、名称、时间、位置、最大风速和最低气压等;
分析复合洪水事件与台风的关系。利用统计分析技术,包括相关性分析、回归分析和聚类分析方法,分析复合洪水事件与台风的关系,确定影响河口区域的主要台风的位置和特征,如距离、方位、强度等;
确定台风影响圆的范围。利用历史数据和统计分析的结果,确定至少一个合理的台风影响圆的半径和圆心,使其能够覆盖河口区域的主要风暴潮来源,如根据河口区域的地理位置和台风的分布特征,选择一个代表性的河口潮位站作为圆心,选择一个能够包含所有影响河口区域的台风的距离作为半径。
步骤二、利用动态调整的台风影响圆的范围,根据台风的实时位置和特征,动态调整台风影响圆的半径和圆心,使其能够适应台风的变化,减少不必要的计算和误差。
获取台风的实时位置和特征的数据。从数据网站获取台风的实时位置和特征的数据,并进行清洗和预处理,剔除异常值和缺失值,统一数据格式和单位,台风的实时位置和特征的数据包括每一场台风的编号、名称、时间、位置、最大风速和最低气压等;
分析台风的实时位置和特征的变化。利用数据分析技术,包括时间序列分析、趋势分析和预测分析,分析台风的实时位置和特征的变化,预测台风的未来走向和强度,可以利用ARIMA模型、指数平滑法或神经网络等方法,根据台风的历史和当前的数据,建立台风的预测模型,输出台风的未来位置和特征的预测值;
动态调整台风影响圆的范围。利用数据分析和预测的结果,动态调整台风影响圆的半径和圆心,使其能够适应台风的变化,减少不必要的计算和误差,如根据台风的预测位置和特征,更新台风影响圆的半径和圆心,使其能够包含所有可能影响河口区域的台风,同时排除一些不相关或不重要的台风。
步骤三、利用多层次的台风影响圆的范围,根据台风的不同等级和强度,设置不同的台风影响圆的半径和圆心,使其能够区分台风的影响程度,提高计算的效率和精度。
从数据网站获取台风的不同等级和强度的标准和定义,并进行清洗和预处理,统一数据格式和单位,台风的不同等级和强度的标准和定义包括每一种等级和强度的台风的名称、最大风速和最低气压的范围等;
利用数据分析技术,包括分类分析、决策树分析和规则挖掘等,分析台风的不同等级和强度的特征,如台风的位置、方位、形状、移速、持续时间等;
利用数据分析的结果,根据台风的不同等级和强度,设置不同的台风影响圆的半径和圆心,使其能够区分台风的影响程度,提高计算的效率和精度,如根据台风的等级和强度,按照预定的比例,缩小或扩大台风影响圆的半径,或者根据台风的方位和形状,调整台风影响圆的圆心,使其更贴近台风的实际位置。
如图4所示,根据本申请的一个方面,所述步骤S3进一步为:
步骤S31、构造机器学习模型的输入数据特征,包括台风影响圆圆心、台风距离特征、台风中心最低气压及最大风速、台风最大风速半径和台风路径方位特征;台风距离特征是每场台风路径与台风影响圆圆心的最小距离;
步骤S32、构建机器学习模型并采用输入数据进行训练和测试,同时对机器学习模型的参数进行调优;
步骤S33、采用训练完成后的机器学习模型对台风进行识别,获得所有影响河口区域的台风路径,形成台风路径集合。
在一些场景中,西太平洋每年约30场台风生成,向北移动并影响某河口区域的台风每年仅2~3场。研究长江历史风暴潮首先需要确定哪些台风在河口区域造成影响,但目前尚无可以快速准确地识别出这些历史台风的方法,尤其是年代久远的台风。通过水动力数值模拟的方法对每一场台风进行计算,虽可以准确识别出台风是否在河口区域产生增水,但这种方法计算量太大,耗时长。通过实测水位和调和分析回报天文潮差值的方法可以较快的计算台风增水,但也存在两个缺陷。其一在于长周期实测潮位资料较难获取;其二在于我国各区域潮位观测时间序列长度一般短于气象资料。比如在河口区域河口段的潮位站多数设立于上世纪80年代初期,而公开的台风资料一般可追述到1949年。在上述实施例中,采用了ISYE-R500(impact scope of the Yangtze Estuary,台风影响圆)的方法通过距离对台风是否影响河口区域进行判定。但这种方法偏于保守,会将一些距离河口区域较近(但风速较弱)、没有在河口区域产生增水的台风计入在内。
针对上述问题,基于随机森林机器学习模型的河口区域影响台风识别方法。基于上述实施例的调和分析研究成果,1989~2020年河口区域增水是可以明确计算出来的。随机森林模型以此为标签数据,对1949~1988年西太平洋发生的台风进行分类(台风路径中风速为修正后的风速)。
在进行模型训练之前,首先需要基于物理过程和实际意义对数据特征进行构造。一般台风距离目标地点越近,台风强度越高,造成增水的可能性越大。因此本实施例将河口区域的位置概化为一个点,计算每场台风路径与该点的距离并以最小值作为距离特征,此时台风采样点的中心最低气压和最大风速作为气压和风速特征,台风最大风速半径为台风半径特征。台风右侧风速一般大于左侧,因而从河口区域南侧登陆的台风引起的增水往往更大。因此增加台风路径方位特征,在河口区域北侧登陆的台风路径赋值为1,其他类台风路径赋值为0。构造时间特征R400和R500,分别代表了台风进入河口区域台风最大影响圆400km和500km范围内的时长。若台风路径距离较远,在其生命过程中未进入上述两区域,则该特征的值为0。
从研究结果可知,台风中心与河口区域的距离对分类结果影响最大,最小距离河口区域台风影响圆R400和R500三特征重要度合计74.4%。最大风速半径特征的重要程度为7.4%,台风是否在河口区域北侧登陆特征的重要程度为6.7%。重要程度最低的反而是气压和风速,占比合计约11.5%。事实上,在其他特征相同的条件下,台风最大风速半径特征越大,表明台风离目标区域越近。台风最大风速半径特征也应计入距离特征内,“距离特征”的重要程度合计超过了80%。因此,在没有其他依据的情况下,根据台风路径到目标地点的距离大致判断其是否会影响到目标区域具有一定合理性。
随机森林机器学习分类模型采用了1989~2020年的台风为训练数据,测试集占比15%。训练数据中包含了没有编号的台风,这样可以帮助模型强化学习无编号台风的特征并对含有该特征的所有台风予以排除。
采用GridSearchCV工具对主要参数进行调优得:弱学习器的最大迭代次数n_estimators=16,最大深度max_depth=5,最大特征数max_features=4,叶子节点最少样本数min_samples_leaf=3,内部节点再划分所需最小样本数min_samples_split=4,判断准则criterion='gini',模型分类准确率为96.6%(通过参数为5的交叉验证)。本实施例研究的主要目的在于从西太平洋台风中尽可能多的识别出影响河口区域的台风,以减少后续的风暴增水数值模拟计算量。因为风暴潮数值模拟可以对增水二次进行验证,所以本工作对模型错误率有一定容忍度。但要求错误最好为第二类错误,减少第一次错误。即可以接受“假台风”(台风没有影响河口区域,但是模型认为影响了),但不能拒绝“真台风”(台风影响河口区域,但模型未识别出)。因此,当概率超过50%,则认为该台风可能在河口区域产生增水。模型验证结果中仅有3场未影响河口区域台风被误识别在内,满足要求。根据研究结果可以知道,1949~1988年西太平洋共发生1283场台风,符合要求,需参与后续数值模拟计算的台风共计90场,后续的数值模拟计算工作量减少了93%。随机森林法识别的台风路径都被ISYE-R500(impact scope of the Yangtze Estuary,台风影响圆)捕获,但相较于仅采用最大影响圆法(共112场台风),随机森林机器学习法仍可使数值模拟最终计算量减少约30%。
综上,采用随机森林机器学习分类模型对西太平洋的台风进行识别,可大幅提高河口区域风暴潮数值模拟的计算效率。
在本申请的另一实施例中,利用多尺度分析和聚类分析对台风路径进行优化,包括利用小波变换、经验模态分解或主成分分析等方法,将台风路径的时空数据分解为不同的尺度和频率,分别分析台风路径的长期趋势、周期变化和短期波动等特征。然后利用聚类分析,包括K均值、层次聚类或模糊聚类,根据台风路径的多尺度特征,将台风路径分为不同的类别,优化台风路径的表示和计算。具体步骤如下:
从数据网站获取预定区域的台风路径的数据,并进行清洗和预处理,剔除异常值和缺失值,统一数据格式和单位,台风路径的数据包括每一场台风的编号、名称、时间、位置、最大风速和最低气压等;
利用多尺度分析技术,包括小波变换、经验模态分解或主成分分析等方法,将台风路径的时空数据分解为不同的尺度和频率,分别分析台风路径的长期趋势、周期变化和短期波动等特征,如利用小波变换,将台风路径的数据分解为近似分量和细节分量,分别反映台风路径的平滑趋势和局部波动;
利用聚类分析技术,包括K均值、层次聚类或模糊聚类等方法,根据台风路径的多尺度特征,将台风路径分为不同的类别,优化台风路径的表示和计算,如利用K均值,根据台风路径的长期趋势和周期变化,将台风路径分为直线型、弯曲型、回旋型等类别,分别采用不同的参数和方法进行台风路径的表示和计算。
如图5所示,根据本申请的一个方面,所述步骤S4进一步为:
步骤S41、获取再分析风场并构建影响河口区域的台风风暴潮的驱动风场集;
步骤S42、以驱动风场集作为驱动,采用三维有限体积法FVCOM进行批量数值模拟,计算近M年所有在河口区域产生影响的台风造成的风暴增水;M为正整数;
步骤S43、提取计算结果中的增水、天文高潮位和风暴高潮位数据,包括场次、时间和高潮位。
基于沿海气象站实测风压资料的台风驱动风场校正方法,提高了风暴潮模型驱动风场的计算效率和近岸风准确度;开发了RHSS工具包,提高了区域历史风暴潮数值模拟的计算效率;结合实测降雨和径流数据重构了近60余年河口区域极端复合洪水灾害因子数据集。
结合1959~2020年era5再分析风场对影响河口区域的台风风暴潮驱动风场集进行了构建,并在此基础上采用FVCOM数值模型计算了河口区域近60余年风暴潮。从研究结果可以看出,天文高潮位一般在3~6m区间内,其10%分位数为5.0m;偏峰增水值一般在1.7m内,其10%分位数为0.83m,超过1.0m的偏峰增水仅占总量的6.8%。风暴高潮位一般在4~7m区间内,其10%分位数为5.4m,前10%大值对应的天文高潮位和偏峰增水均值分别为5.06m和0.92m。可见,极端风暴高潮位一般由极端天文大潮和中等或强增水构成。需指出的是,在自然演变和人类活动的共同作用下,长江下游河道经历了复杂的历史变迁。本实施例仅以2016~2020年的实测地形数据为基础对河口区域水动力进行建模计算。在河口区域现有岸线地形条件和防护水平的基础上,对历史台风再次发生时可能产生的风暴灾害后果和风险进行研究。
在本申请的另一实施例中,还包括:根据河口区域的不同特点,设置不同的模型分辨率和参数,如在河口入海口和河口内部,采用高分辨率和高精度的模型,如三维有限体积法FVCOM;在河口外部和远海,采用低分辨率和低精度的模型,如二维有限差分法ADCIRC。利用多层次的模型,可以兼顾计算的精度和速度。
获取河口区域的地形数据。从数据网站获取河口区域的地形数据,并进行清洗和预处理,剔除异常值和缺失值,统一数据格式和单位,河口区域的地形数据包括河口的水深、宽度、长度、坡度等;
设置不同的模型分辨率和参数。根据河口区域的不同特点,设置不同的模型分辨率和参数,如在河口入海口和河口内部,采用高分辨率和高精度的模型,如三维有限体积法FVCOM,设置较小的网格尺寸和较高的时间步长,考虑河口的三维水动力和海陆交互等因素;在河口外部和远海,采用低分辨率和低精度的模型,如二维有限差分法ADCIRC,设置较大的网格尺寸和较低的时间步长,忽略河口的垂向变化和非线性效应等因素;
构建多层次的河口风暴潮数值模型。利用不同的模型分辨率和参数,构建多层次的河口风暴潮数值模型,如利用FVCOM和ADCIRC的耦合模式,将河口入海口和河口内部的高分辨率模型与河口外部和远海的低分辨率模型进行连接,实现模型的多层次结构,提高模型的计算效率和精度。
如图6所示,根据本申请的一个方面,所述步骤S5进一步为:
步骤S51、获取河口上游预定地点的实测流量数据,包括逐日流量值和极值流量值;
步骤S52、获取并从中国地面气候资料日值数据集中提取河口区域附近地区的日降雨量数据,包括逐日降雨量值和极值降雨量值;
步骤S53、根据得到的风暴高潮位时间,基于实测流量在台风期间的变异性确定上游流量值,作为极端复合洪水的一个灾害因子;基于中国地面气候资料日值数据集,采用滑动相关性法,确定日台风降雨量及其对应时间,作为极端复合洪水的另一个灾害因子;基于修正后的台风最大风速和计算得到的风暴增水,作为极端复合洪水的另外两个灾害因子;利用调和分析得到的短期天文潮,作为极端复合洪水的最后一个灾害因子;
步骤S54、将得到的上游径流、降雨、天文潮、风暴潮和增水数据集按照时间对齐,构建河口区域极端复合洪水灾害因子数据集。
台风在登陆期间往往伴随着大范围的暴雨,导致径流暴涨,加剧了城市洪水淹没风险。本实施例根据风暴潮高潮位的出现时间对河口区域暴雨和洪水时间进行筛选。日累积降雨量数据来自国家气象科学数据中心的中国地面气候资料日值数据集,采用某市附近三个雨量站点日累积降雨作为台风降雨参考站点。
如图7、8和9所示,根据本申请的一个方面,所述步骤S6进一步为:
步骤S61、构建藤Copula函数集合并筛选,所述藤Copula函数至少包括C藤、D藤和R藤;
步骤S62、读取复合洪水的时间序列数据,根据台风编号和台风名称进行分组,得到一个按台风分组的列表,每个分组包含一个台风在河口区域各个潮位站造成的复合洪水的时间序列;
步骤S63、利用R-vine Copula模型对每个台风分组进行数据分析,建立复合洪水多灾害因子之间的联合概率分布模型;利用AIC和BIC指标,对R-vine Copula模型进行选择和验证,选择最优的模型参数和结构;
步骤S64、利用参数优化后的R-vine Copula模型,进行复合洪水的概率分析和风险评估,计算不同复合洪水等级的发生概率和超越概率,分析复合洪水的风险等级和影响范围。
当联合分布的变量较多时(>3),传统Copula函数求解复杂。配对(Pair)Copula模型采用一系列双变量Copula组成的层次结构来描述高维联合分布,其原理是基于条件概率公式将d维的多元函数概率密度函数分解为d(d-1)/2个二维Copula密度函数与边缘概率密度的乘积。其中前n-1个无条件Copula,其余为条件Copula。藤(Vine)Copula是基于正则树结构的多变量Copula建模函数。VineCopula模型中可以衡量直接和间接的灾害因子传导渠道,衡量不同灾害因子之间的依赖性。由于分解过程较为复杂,一般借助图论进行分解。
根据本申请的一个方面,所述步骤S64还包括:
步骤S641、根据防洪标准和历史洪水资料,确定不同复合洪水等级对应的水位高度和流量大小;
步骤S642、利用参数优化后的R-vine Copula模型,对每个台风分组的复合洪水数据进行模拟,得到不同复合洪水等级的发生概率和超越概率,以及其置信区间;
步骤S643、根据复合洪水等级的发生概率和超越概率,评估不同台风分组的复合洪水风险等级,分成若干个等级,根据复合洪水的水位高度和流量大小,分析不同台风分组的复合洪水影响范围;
步骤S644、根据复合洪水的风险等级和影响范围,制定长期的防洪规划和减灾策略。
具体的,包括如下步骤:
构建模型,站点径流(简写F,编号1,下同),日累积降雨量(R,2),偏峰增水(S,3),风暴潮(ST,4)和天文潮(T,5)五个河口极端复合洪水灾害因子,根据AIC准则从GaussianCopula、Studentt Copula、Frank Copula和Clayton Copula、Gumbel Copula、Joe Copula及其旋转90°、180°、270°共15个Copula函数中选择最佳Copula构建基于R-vine Copula架构的河口区域多因子极端复合洪水风险模型。
模型验证,借鉴自助法的思想,将估计的Copula参数代入R-vine Copula五维河口区域多因子极端复合洪水风险模型中,生成1000组符合分布结构的拟合数据(维度为1000×5),计算其F,R,S,ST,T两两之间的Kendall相关系数。将以上工作重复100次,并把计算结果绘制箱线图与原始数据Kendall系数进行比较。可以看出,除站点流量(F)相关性较小外,多数变量间的原始Kendall系数均落在箱型图上、下四分位数之间,并且大部分在中位数附近,说明R-vine Copula模型对变量的相关结构刻画较为可靠。
结果分析,首先对R-vine Copula结构非条件Copula组分进行分析。从空间关系上看,站点流量的增减主要受长江中上游降雨汇流和三峡出流量控制的影响,与下游的台风风暴潮动力系统相关性不强;从时间关系上看,站点水文站水流达到河口区域需4~6天时间,此时的流量并不受台风系统影响。
其次对R-vine Copula结构中单条件Copula组分进行分析,当增水为先验条件时,日累计降雨量和风暴高潮位有0.05的kendall秩相关性,说明了日累计降雨量和风暴高潮位间接依赖结构为正,但相关性较弱。当风暴高潮位为先验条件时,站点流量和天文高潮位有-0.04的kendall秩相关性,说明了站点流量和天文高潮位间接依赖结构为负,但相关性同样不明显。当风暴高潮位为先验条件时,增水和天文高潮位有-0.74的kendall秩相关性,说明了增水和天文高潮位间接依赖结构为负,相关性较强,一个增长则另一个相应减弱。
再次对R-vine Copula结构中两条件Copula组分进行分析。Copula结构的先验条件越多,其灾害依赖结构更为间接,理论上其相关性会逐渐递减。当风暴高潮位和增水为先验条件时,天文高潮位和日累计降雨量相关性接近0,相关性极弱。当风暴高潮位和天文高潮位为先验条件时,增水和站点流量相关性为负,且绝对值小于0.05,说明了增水和天文高潮位间接依赖结构为负,一个增长对应另一个减弱,但相关性微弱。
最后对R-vine Copula结构中三条件Copula组分进行分析。当增水,风暴高潮位和天文高潮位为先验条件时,站点径流和降雨的相关性为0,表明间接依赖结构在三条件Copula组分中不明显。
综上,河口区域风暴潮灾害因子间存在明显的依赖结构,具体表现为:天文高潮位和风暴高潮位强相关(相关性0.73),历史上极端风暴高潮位基本都是在台风遭遇天文大潮的基础上产生的;增水和风暴高潮位上尾相关,对极值风暴潮位具有一定贡献(相关性0.34);日累计降雨量和增水之间表现为弱相关(相关性0.26),且不具备尾部相关性。站点径流与其余四灾害因子的相关性极弱。当风暴高潮位为先验条件时,增水和天文高潮位间接依赖结构为负(相关性-0.74)。其余条件依赖关系均较弱(相关性绝对值低于0.05)。
在本申请的另一实施例中,当站点径流,日累计降雨量,增水,风暴高潮位和天文高潮位的基准重现期为五年一遇时,各灾害因子对应的联合分布概率密度值相对其他工况均较大,但各因子之间的差别较小,说明各因子对于联合风险均较为敏感。当单一灾害因子重现期增大时,各因子对于联合风险敏感性均发生不同程度的降低。极端复合洪水风险对站点流量和增水敏感性大幅下降,站点流量敏感性减小幅度最大,从最敏感的因子之一变为敏感性最弱的因子,即站点流量的增加几乎不会对五年基准重现期下河口区域极端复合洪水风险产生影响,可将其从关键风险因子中排除;日累计降雨量在其重现期为十年一遇时的概率密度远大于其他变量,说明五年基准重现期与单变量十年重现期组合情节下,日累计降雨量为最敏感风险因子;随着重现期的增加,风暴高潮位和天文高潮位下降较为缓慢,说明两者为极端复合洪水主要灾害因子。
当基准重现期为十年一遇时,单变量重现期为五年一遇的站点流量和增水的概率密度最大,说明该情景组合下站点流量和增水对联合风险率的影响最大,为敏感风险因子。但是当单变量重现期超过五十年,站点流量和增水联合风险率的影响大幅降低;随着单变量重现期的增加,日累计降雨量的概率密度经历了先增大后减小的过程,并在二十年一遇时达到最大(且大于其他灾害因子),说明在十年一遇基准重现期、单变量二十年重现期情景下,日累计降雨量为最敏感风险因子;随着重现期的增加,风暴高潮位和天文高潮位随重现期的增加下降较为缓慢,说明两者为极端复合洪水主要灾害因子。
当基准重现期为二十年一遇时,站点流量(十年一遇)为最敏感因子,增水(五年、十年一遇)次之。当单一灾害因子重现期增大至五十年一遇时,极端复合洪水风险对站点流量和增水敏感性几乎降至0,说明单一灾害因子重现期在五年~五十年一遇范围内时,站点流量和增水对联合风险起主导作用;而超过五十年一遇后,日累计降雨量,天文高潮位和风暴高潮位起主导作用。随着重现期的增加,风暴高潮位和天文高潮位随重现期的增加下降较为缓慢,说明两者为极端复合洪水主要灾害因子。
当基准重现期为五十年一遇时,站点流量和增水(二十年一遇)为最敏感因子,天文高潮位次之。随着重现期的增加,站点流量和增水的概率密度都经历了“先增大,后减小”的过程,并在二十年一遇重现期时达到最大。表明站点流量和日累计降雨量在“小重现期”(重现期年)时对极端复合洪水风险的影响是不断增强的,在“大重现期”(重现期为50年和100年)时对极端复合洪水风险的影响不断减弱的。随着重现期的增加,风暴高潮位和天文高潮位随重现期的增加下降较为缓慢,说明两者为极端复合洪水主要灾害因子。
当基准重现期为百年一遇时,天文高潮位(五年一遇)为最敏感因子,增水(二十年一遇)和站点流量(五十年一遇)次之。从整体上看,站点径流和降雨对联合灾害风险率的影响十分有限,增水,天文高潮位和风暴高潮位起主导作用。
综上,对于河口区域多因子极端复合洪水风险模型,多数情景下风暴高潮位和天文高潮位都是复合洪水主要灾害因子,但不同基准重现期下复合洪水灾害的敏感风险因子不同:五年一遇基准重现期下,各因子均为敏感因子;十年、二十年和五十年一遇基准重现期下,站点流量和增水为敏感因子;百年一遇基准重现期下,天文高潮位为敏感因子。
为分析河口区域极端复合洪水灾害中站点径流,日累计降雨量,增水,风暴高潮位和天文高潮位的敏感性,识别关键风险因子,在完成对R-vine Copula函数建模的基础上,采用五年、十年、二十年、五十年和百年基准重现期,改变其中某个灾害因子的重现期,检验其对风险率的贡献值,识别出关键风险因子。
采用五维风暴潮灾害因子联合分布的概率密度值分析极端复合洪水灾害风险对各因子的敏感性。在各因子处于同一重现期水平时,其所对应的联合分布概率密度值越大,说明该因子对联合分布的变化率影响越大,即联合风险对其越敏感。
根据本申请的另一个方面,提供一种资料受限河口区域极值事件识别和复合洪水灾害分析系统,其特征在于,包括:
至少一个处理器;以及
与至少一个所述处理器通信连接的存储器;其中,
所述存储器存储有可被所述处理器执行的指令,所述指令用于被所述处理器执行以实现上述任一项技术方案所述的资料受限河口区域极值事件识别和复合洪水灾害分析方法。
以上详细描述了本发明的优选实施方式,但是,本发明并不限于上述实施方式中的具体细节,在本发明的技术构思范围内,可以对本发明的技术方案进行多种等同变换,这些等同变换均属于本发明的保护范围。
Claims (5)
1.资料受限河口区域极值事件识别和复合洪水灾害分析方法,其特征在于,包括如下步骤:
步骤S1、采集预定区域的台风历史数据并预处理和修正台风最大风速,构建适用于台风模拟的驱动风场;
步骤S2、确定河口区域的台风影响圆;采用调和分析法基于已有实测潮位数据计算风暴增水;
步骤S3、基于台风影响圆确定的范围,筛选出影响河口区域的所有台风路径,形成台风路径集合;
步骤S4、构建河口风暴潮数值模型,进行数值模拟,计算所有在河口区域产生影响的台风造成的风暴增水,并获得风暴高潮位时间;
步骤S5、根据风暴高潮位时间,基于实测流量在台风期间变异性确定上游流量值,调取中国地面气候资料日值数据集,采用滑动相关性法确定日台风降雨量及其对应时间;
步骤S6、基于步骤S3至S5的输出结果,采用R-vine Copula对河口区域复合洪水概率进行分析并输出分析结果;
所述步骤S1进一步为:
步骤S11、从数据网站获取预定区域的台风历史数据,并进行清洗和预处理,剔除异常值和缺失值,统一数据格式和单位,台风历史数据包括每一台风的编号、名称、时间、位置、最大风速和最低气压;
步骤S12、从预定数据库获取再分析数据集,提取与台风历史数据对应的风场和气压场数据,并进行插值和平滑;
步骤S13、采用公式法和预训练的Light GBM机器学习模型对台风最大风速进行修正,获得修正后的台风最大风速;其中公式法包括如下修正公式:
Vm1=a(φ,4971)(1010-P0)b(φ,4971);Vm2=a(φ,7220)(1010-P0)b(φ,7220);
Vm'=Vm-(Vm1-Vm2);其中,a(φ,4971)、b(φ,4971)、a(φ,7220)、b(φ,7220)分别为1949~1971年和1972~2020年时段内采用纬度为φ的台风数据拟合得到的参数;Vm和P0分别为待订正的风速和对应气压;Vm'为修正后的风速;
步骤S14、将修正后的台风最大风速与再分析数据集中的风场和气压场数据结合,构建台风模拟的驱动风场;
所述步骤S2进一步为:
步骤S21、针对河口区域的特点进行概化,获取圆心坐标位置,并以N倍的台风最大风速半径为区域半径,获得台风影响圆,获取台风影响圆圆周上若干点的坐标并存储;N为正整数;
步骤S22、获取河口区域的河口潮数据并利用调和分析法对河口潮的实测潮位数据进行分解,得到天文潮和每个潮周期内的偏峰增水;
步骤S23、修正偏峰增水并计算风暴增水以及风暴增水的极值和频率分布;
所述步骤S3进一步为:
步骤S31、构造机器学习模型的输入数据特征,包括台风影响圆圆心、台风距离特征、台风中心最低气压及最大风速、台风最大风速半径和台风路径方位特征;台风距离特征是每场台风路径与台风影响圆圆心的最小距离;
步骤S32、构建机器学习模型并采用输入数据进行训练和测试,同时对机器学习模型的参数进行调优;
步骤S33、采用训练完成后的机器学习模型对台风进行识别,获得所有影响河口区域的台风路径,形成台风路径集合;
所述步骤S4进一步为:
步骤S41、获取再分析风场并构建影响河口区域的台风风暴潮的驱动风场集;
步骤S42、以驱动风场集作为驱动,采用三维有限体积法FVCOM进行批量数值模拟,计算近M年所有在河口区域产生影响的台风造成的风暴增水;M为正整数;
步骤S43、提取计算结果中的增水、天文高潮位和风暴高潮位数据,包括场次、时间和高潮位;
所述步骤S5进一步为:
步骤S51、获取河口上游预定地点的实测流量数据,包括逐日流量值和极值流量值;
步骤S52、获取并从中国地面气候资料日值数据集中提取河口区域附近地区的日降雨量数据,包括逐日降雨量值和极值降雨量值;
步骤S53、根据得到的风暴高潮位时间,基于实测流量在台风期间的变异性确定上游流量值,作为极端复合洪水的一个灾害因子;基于中国地面气候资料日值数据集,采用滑动相关性法,确定日台风降雨量及其对应时间,作为极端复合洪水的另一个灾害因子;基于修正后的台风最大风速和计算得到的风暴增水,作为极端复合洪水的另外两个灾害因子;利用调和分析得到的短期天文潮,作为极端复合洪水的最后一个灾害因子;
步骤S54、将得到的上游径流、降雨、天文潮、风暴潮和增水数据集按照时间对齐,构建河口区域极端复合洪水灾害因子数据集。
2.如权利要求1所述的资料受限河口区域极值事件识别和复合洪水灾害分析方法,其特征在于,所述步骤S6进一步为:
步骤S61、构建藤Copula函数集合并筛选,所述藤Copula函数至少包括C藤、D藤和R藤;
步骤S62、读取复合洪水的时间序列数据,根据台风编号和台风名称进行分组,得到一个按台风分组的列表,每个分组包含一个台风在河口区域各个潮位站造成的复合洪水的时间序列;
步骤S63、利用R-vine Copula模型对每个台风分组进行数据分析,建立复合洪水多灾害因子之间的联合概率分布模型;利用AIC和BIC指标,对R-vine Copula模型进行选择和验证,选择最优的模型参数和结构;
步骤S64、利用参数优化后的R-vine Copula模型,进行复合洪水的概率分析和风险评估,计算不同复合洪水等级的发生概率和超越概率,分析复合洪水的风险等级和影响范围。
3.如权利要求2所述的资料受限河口区域极值事件识别和复合洪水灾害分析方法,其特征在于,所述步骤S13还包括:
步骤S13a、分别采用公式法和预训练的Light GBM机器学习模型对台风最大风速进行修正,获得修正后的台风最大风速;
步骤S13b、构建风速窗口,并分别计算每个风速窗口内公式法和Light GBM机器学习模型计算得到的台风最大风速的均方根误差;
步骤S13c、根据均方根误差,确定每个风速窗口的修正模型,获得各个风速窗口的修正数据,形成复合修正后的台风最大风速。
4.如权利要求3所述的资料受限河口区域极值事件识别和复合洪水灾害分析方法,其特征在于,所述步骤S64还包括:
步骤S641、根据防洪标准和历史洪水资料,确定不同复合洪水等级对应的水位高度和流量大小;
步骤S642、利用参数优化后的R-vine Copula模型,对每个台风分组的复合洪水数据进行模拟,得到不同复合洪水等级的发生概率和超越概率,以及其置信区间;
步骤S643、根据复合洪水等级的发生概率和超越概率,评估不同台风分组的复合洪水风险等级,分成若干个等级,根据复合洪水的水位高度和流量大小,分析不同台风分组的复合洪水影响范围;
步骤S644、根据复合洪水的风险等级和影响范围,制定长期的防洪规划和减灾策略。
5. 一种资料受限河口区域极值事件识别和复合洪水灾害分析系统,其特征在于,包括:
至少一个处理器;以及
与至少一个所述处理器通信连接的存储器;其中,
所述存储器存储有可被所述处理器执行的指令,所述指令用于被所述处理器执行以实现权利要求1至4任一项所述的资料受限河口区域极值事件识别和复合洪水灾害分析方法。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202410110058.5A CN117634325B (zh) | 2024-01-26 | 2024-01-26 | 资料受限河口区域极值事件识别和复合洪水灾害分析方法和系统 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202410110058.5A CN117634325B (zh) | 2024-01-26 | 2024-01-26 | 资料受限河口区域极值事件识别和复合洪水灾害分析方法和系统 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN117634325A CN117634325A (zh) | 2024-03-01 |
CN117634325B true CN117634325B (zh) | 2024-04-02 |
Family
ID=90036038
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202410110058.5A Active CN117634325B (zh) | 2024-01-26 | 2024-01-26 | 资料受限河口区域极值事件识别和复合洪水灾害分析方法和系统 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN117634325B (zh) |
Families Citing this family (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN118396031A (zh) * | 2024-04-22 | 2024-07-26 | 浙江省海洋监测预报中心 | 一种风暴潮位实时校正方法 |
Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN108320055A (zh) * | 2018-01-26 | 2018-07-24 | 中国海洋大学 | 多个河口风暴潮增水联合重现期的确定方法 |
CN112036691A (zh) * | 2020-07-24 | 2020-12-04 | 东南大学 | 一种基于jpm-os-q模型的河口地区极端风暴潮水位计算方法 |
CN114756817A (zh) * | 2022-02-22 | 2022-07-15 | 南方科技大学 | 一种基于Copula函数的复合洪涝灾害联合概率分析方法 |
CN116522764A (zh) * | 2023-04-17 | 2023-08-01 | 华中科技大学 | 一种考虑气候变化影响下的热浪-洪水复合灾害评估方法 |
Family Cites Families (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN112966926B (zh) * | 2021-03-02 | 2022-04-22 | 河海大学 | 一种基于集成学习的洪水敏感性风险评估方法 |
-
2024
- 2024-01-26 CN CN202410110058.5A patent/CN117634325B/zh active Active
Patent Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN108320055A (zh) * | 2018-01-26 | 2018-07-24 | 中国海洋大学 | 多个河口风暴潮增水联合重现期的确定方法 |
CN112036691A (zh) * | 2020-07-24 | 2020-12-04 | 东南大学 | 一种基于jpm-os-q模型的河口地区极端风暴潮水位计算方法 |
CN114756817A (zh) * | 2022-02-22 | 2022-07-15 | 南方科技大学 | 一种基于Copula函数的复合洪涝灾害联合概率分析方法 |
CN116522764A (zh) * | 2023-04-17 | 2023-08-01 | 华中科技大学 | 一种考虑气候变化影响下的热浪-洪水复合灾害评估方法 |
Non-Patent Citations (3)
Title |
---|
基于云模型的汉江上游安康市洪水灾害风险评价;石晓静 等;水利水电科技进展;20170510(第03期);全文 * |
沿海地区复合洪水危险性研究进展;方佳毅 等;气候变化研究进展;20210531;第17卷(第3期);全文 * |
洪水对长江口风暴潮的影响;罗俐雅 等;水运工程;20230331(第3期);全文 * |
Also Published As
Publication number | Publication date |
---|---|
CN117634325A (zh) | 2024-03-01 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
US9262723B2 (en) | Predicting climate data using climate attractors derived from a global climate model | |
CN117634325B (zh) | 资料受限河口区域极值事件识别和复合洪水灾害分析方法和系统 | |
CN114049545B (zh) | 一种基于点云体素的台风定强方法、系统、设备及介质 | |
CN115423163A (zh) | 一种流域短期洪水事件预测方法、装置及终端设备 | |
CN114935759A (zh) | 基于高频地波雷达观测的波浪场缺失值填补方法及系统 | |
CN116777079A (zh) | 一种基于贝叶斯层间结构模型的沙漠化侵扰灾害预测方法 | |
CN117113854A (zh) | 一种基于ConvLSTM和三维数值模拟的咸潮预报方法 | |
CN117556197A (zh) | 一种基于人工智能的台风涡旋初始化方法 | |
Jeon et al. | Characterization of extreme precipitation within atmospheric river events over California | |
Lubkov et al. | Method for reconstructing the monthly mean water transparencies for the northwestern part of the Black Sea as an example | |
CN114676621A (zh) | 基于深度学习权重负载提高陆地水储量异常准确性方法 | |
Baordo et al. | Intercomparison and assessement of wave models at global scale. | |
Hema et al. | Reconstructing missing hourly real-time precipitation data using a novel intermittent sliding window period technique for automatic weather station data | |
Rajabi-Kiasari et al. | Forecasting of absolute dynamic topography using deep learning algorithm with application to the Baltic Sea | |
CN116976227A (zh) | 一种基于lstm机器学习的风暴增水预报方法及系统 | |
CN116882731A (zh) | 一种基于斜坡单元的地质灾害危险性评估方法及系统 | |
Majdzadeh Moghadam | Neural network-based approach for identification of meteorological factors affecting regional sea-level anomalies | |
CN117688853B (zh) | 基于短期潮位和长期气象资料的区域风暴潮破坏性评估方法及系统 | |
CN110751398A (zh) | 一种区域生态质量评价方法及装置 | |
Shi et al. | A Grey Model for Short-Term Prediction of the Ionospheric TEC | |
CN117709242B (zh) | 一种河口咸潮上溯强度评估方法 | |
CN116842766B (zh) | 一种基于时空组合分解电离层异常的全球tec经验模型 | |
CN118194235B (zh) | 一种多模式融合的降水干湿时序变化重构方法及系统 | |
CN116204836B (zh) | 一种河口滩涂红树林宜林区筛选方法 | |
Olsen | Long Term Trends in Lake Michigan Wave Climate |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |