Global Projections of Storm Surges Using High‐Resolution CMIP6 Climate Models

In the coming decades, coastal flooding will become more frequent due to sea‐level rise and potential changes in storms. To produce global storm surge projections from 1950 to 2050, we force the Global Tide and Surge Model with a ∼25‐km resolution climate model ensemble from the Coupled Model Intercomparison Project Phase 6 High Resolution Model Intercomparison Project (HighResMIP). This is the first time that such a high‐resolution ensemble is used to assess changes in future storm surges across the globe. We validate the present epoch (1985–2014) against the ERA5 climate reanalysis, which shows a good overall agreement. However, there is a clear spatial bias with generally a positive bias in coastal areas along semi‐enclosed seas and negative bias in equatorial regions. Comparing the future epoch (2021–2050) against the historical epoch (1951–1980), we project ensemble‐median changes up to 0.1 (or 20%) in the 1 in 10‐year storm surge levels. These changes are not uniform across the globe with decreases along the coast of Mediterranean and northern Africa and southern Australia and increases along the south coast of Australia and Alaska. There are also increases along (parts) of the coasts of northern Caribbean, eastern Africa, China and the Korean peninsula, but with less agreement among the HighResMIP ensemble. Information resulting from this study can be used to inform broad‐scale assessment of coastal impacts under future climate change.


Plain Language Summary
In the next few decades, coastal flooding is expected to become more frequent due to rising sea levels and changes in storms.To understand and prepare for these changes, we use a global hydrodynamic and multiple climate models to investigate how storm surges (a temporarily rise in water level during a storm) will respond to a warmer climate.These are the first global projections of storm surges based on highly detailed climate models (25-50 km).We compare the model output for recent years  with real-world weather data, known as reanalysis data, and find a good overall agreement.However, the comparison also shows relatively large differences with larger or smaller values in certain areas.For the future , results show storm surge changes up to 20%.These changes vary around the world.Some areas (Mediterranean, northern Africa, and southern Australia) might experience lower storm surges, whereas other areas (South Australia, Alaska, the northern Caribbean, eastern Africa, China, and the Korean Peninsula) might experience higher storm surges.By using advanced computer models as done here, we are able to better understand how climate change could impact coastal area and thus make informed decisions for the future. et al., 2018).Even if global warming is limited to 1.5°C above pre-industrial levels, as outlined in the Paris Agreement, approximately half of the world's coastline will experience the present-day 100-year ESL at least once a year by 2080 (Tebaldi et al., 2021).Besides higher mean sea levels, the warming of the climate and change in the large-scale atmospheric circulation can modulate the occurrence of tropical cyclones and extra-tropical storms (Seneviratne et al., 2021).The (thermo)dynamic response of storms to a warmer climate is a complex process (Shaw et al., 2016).While the total number of tropical cyclones may decrease or remain unchanged (Sobel et al., 2021), tropical cyclone intensities are projected to increase in response to warmer sea surface temperatures (Knutson et al., 2010(Knutson et al., , 2019(Knutson et al., , 2020;;Walsh et al., 2016).Most projections show an increase in the proportion of intense tropical cyclones (category 4-5) and higher maximum wind speeds (Knutson et al., 2020).In a warmer climate, tropical cyclone tracks will shift poleward with implications for the mid-latitude regions that may experience more post-tropical storms (Haarsma, 2021;Haarsma et al., 2013).Future projections of mid-latitude storms generally show small changes, although due to a poleward shift in the extra-tropical storm tracks, some regions may see substantial changes (Morim et al., 2019;Priestley & Catto, 2022a;Seneviratne et al., 2021).Changes in storm surge frequencies depend on the combined effect of changes in storm tracks, frequency and intensity.Since the response to a warming climate is influenced by many opposing factors both increases and decreases may occur.
Continental to global-scale ESL models have evolved rapidly in recent years, and now include more physical processes at increasing resolution (Bouwer, 2018;Wahl, 2017).This has led to an enhanced understanding of ESLs and associated coastal flooding (Muis et al., 2016), and how flood hazard and risks may be impacted by climate change and adaptation (Kirezci et al., 2020;Mori et al., 2019;Tiggeloven et al., 2020;Vitousek et al., 2017;Vousdoukas et al., 2018).Mean sea level will be the largest driver of ESL changes (Vousdoukas et al., 2018), and most global studies have focused on future ESL changes resulting from mean sea level changes (Frederikse et al., 2020;Kirezci et al., 2020;Rasmussen et al., 2018;Tebaldi et al., 2021;Vitousek et al., 2017).This approach neglects the high potential but not yet fully understood contribution from changes in tropical cyclones and extratropical storms (Morim et al., 2019;Seneviratne et al., 2021).There are some studies have explored how storm surge frequencies may change under climate change, often at regional to national scale (e.g., Colberg et al., 2019;De Dominicis et al., 2020;Garner et al., 2017;Howard et al., 2019;Lin et al., 2019;Little et al., 2015;Marsooli et al., 2019a), but more recently also at the continental to global scale (Colberg et al., 2019;Mori et al., 2019;Muis et al., 2020;Vousdoukas et al., 2016Vousdoukas et al., , 2018)).However, studies that have provided global to continental-scale projections of storm surges are generally based on Global Climate Models (GCMs) from the Coupled Model Intercomparison Project -Phase 5 (CMIP5) experiments (Colberg et al., 2019;Lin et al., 2019;Muis et al., 2020;Vousdoukas et al., 2018).A major limitation of the CMIP5 experiments is that the GCM's spatial resolution (∼150 km) is insufficient to resolve localized climate extremes, such as tropical cyclones (Camargo, 2013;Hodges et al., 2017;Schenkel & Hart, 2012).With further GCM development and increases in resolution (Bauer et al., 2015), there is a growing number of GCMs that can generate a credible climatology of both TC numbers and intensities in the current climate (Murakami et al., 2012;Roberts et al., 2020;Walsh et al., 2016).The High Resolution Model Intercomparison Project (HighResMIP), which is part of the CMIP6 framework, provides climate projections of models with a spatial resolution of ∼25-50 km (Haarsma et al., 2016).Roberts et al. (2020) has shown that TC frequencies are comparable to observations, and that increased horizontal resolution improves the representation of TC intensity and tracks, although the strongest TCs may still be underestimated.For extra-tropical cyclones increased resolution also leads to a better representation of the structure and intensity (Priestley & Catto, 2022b).Hence, the HighResMIP experiments provide the opportunity to improve upon existing global storm surge projections and enhance the understanding of the response the climate change.
Here we present novel global storm surge projections that are derived from a global hydrodynamic model forced with the HighResMIP climate model ensemble.While the full data set includes total water levels, tides and storm surges, here we focus on storm surges.First, we evaluate the model performance of storm-surge levels derived from the HighResMIP ensemble for the present epoch  against those derived from the ERA5 climate reanalysis.Next, we assess projected storm surge changes for the future epoch  and assess the intermodel agreement across ensemble.

Methods and Data
Figure 1 summarizes the three-step methodology that was used to generate the global storm surge projections.First, we simulate 100-years of continuous time-series of total water levels (defined as tides and storm surges) by forcing the Global Tide and Surge Model (GTSMv3.0) with sea-level pressure and wind speed and direction from an ensemble of HighResMIP models.We also do a tide-only simulations and compute storm surge levels as the difference between total water level and tidal level.Second, we split the time-series into three 30-year periods (historical epoch: 1951-1980; present epoch: 1985-2014; and future epoch: 2021-2050) and perform extreme value analysis (EVA).Last, we validate the present epoch and project how the magnitude and frequency of storm surges may be impacted by future climate change by comparing the future and present epoch against the historical epoch.Below each of the steps is explained in more detail.

Global Modeling of Storm Surges
We use the Global Tide and Surge Model (GTSMv3.0) to simulate total water levels resulting from tidal and meteorological forcing (Muis et al., 2020;Wang et al., 2021).GTSM is a depth-averaged hydrodynamic model based on the Delft3D Flexible Mesh modeling suite (Kernkamp et al., 2011).The grid resolution ranges from 2.5 km along the coast (1.25 km in Europe) to 25 km in the ocean.We obtain time-series of storm surge levels by computing the difference between a tide-only simulation and a total water level simulation.This way, non-linear interactions between tides, storm surges and mean sea levels are dynamically included.For the total water level simulation, meteorological forcing (i.e., 10 m elevation wind speed and sea level pressure from the climate models) and tidal forcing (i.e., the gravitational attraction and centrifugal forces of the Earth, Moon, and Sun) are applied.The tide-only simulation excludes meteorological forcing.Both simulations include a spatially varying sea-level rise field based on CMIP5 models for the entire simulation period (see Muis et al. (2020) for methodology; and Text S1 and Figure S1 in Supporting Information S1 for the mean-sea level projections).However, when computing the storm surge levels, mean sea level changes are removed.Time series at 10-min resolution are produced for a set of 43,119 output points (Muis et al., 2020).GTSM has been extensively validated in previous studies (Bloemendaal et al., 2018;Dullaart et al., 2020;Muis et al., 2016Muis et al., , 2019Muis et al., , 2020)), which have shown good agreement between modeled and observed water levels.Using a global set of ∼300 tide gauge stations, Muis et al. (2020) showed that the mean bias of modeled water levels (i.e., tides and storm surge) with a 10-year return period for GTSMv3.0 with the ERA5 climate reanalysis are underestimated by −0.10 m with a standard deviation of 0.32.The tidal performance is comparable to other non-assimilating models (Muis et al., 2022a(Muis et al., , 2022b;;Wang et al., 2021;Wang, Verlaan, Apecechea, & Lin, 2022;Wang, Verlaan, Veenstra, & Lin, 2022), and considered to be sufficient for this type of application.

HighResMIP Climate Projections
To analyze the storm surge changes under climate change, we force GTSM with meteorological fields of the multi-model HighResMIP ensemble (Haarsma et al., 2016), consisting of five ∼25 km resolution GCMs.The HighResMIP experiments span the period 1950-2050.For the historic period , the forcing fields of the HighResMIP simulations are similar to those used in CMIP6 (Eyring et al., 2016), although at higher resolution.For the future period (2015-2050), the forcing fields of the HighResMIP simulations are based on the SSP5-8.5 greenhouse gas concentration scenario (O'Neill et al., 2016).We use a mix of coupled atmosphere-ocean and atmosphere-only simulations; specifically, three coupled simulations (HadGEM3-GC31-HM (Roberts, 2017), Earth3P-HR (EC-Earth Consortium, 2018), CMCC-CM2-VHR4 (Scoccimarro et al., 2017)) and two atm-only simulations (HadGEM3-GC31-HM-SST (Roberts, 2017), GFDL-CM4C192-SST (Zhao et al., 2018)).The atmosphere-only simulations are forced by observed estimates of SST and sea ice concentration in the historic period, and a combination of observed historical variability and modeled future trends in the future period (Haarsma et al., 2016), with the same radiative forcing as the coupled models.TCs derive their energy from the upper ocean, and by including atmosphere-ocean interactions the coupled simulations may have more realistic storm intensities and duration and climate variability (Li & Sriver, 2018).At the same time, coupled models are characterized by large SST biases and generally show too cold in the northern hemisphere and too warm in the southern hemisphere (He & Soden, 2016;C. Wang et al., 2014), which contributes to the projection uncertainty.Ideally, we would have used a larger model ensemble based on coupled simulations only.However, at the time of carrying out the GTSM simulations these were the only simulations available which produced the required output fields with a minimum temporal resolution of at least 6 hr.

Extreme Value Statistics
We split the storm surge time-series into three epochs : 1951-1980 (reference), 1985-2014 (present), and 2021-2050 (future).Subsequently, we compute various percentiles (1, 5, 10, 25, 50, 75, 90, 95, 99, 99.5, 99.9) at each output location for each of these three epochs.Next, we apply EVA to estimate the exceedance probabilities of storm surge levels.Following Wahl et al. (2017), we use the Peak Over Threshold method (POT) and we fit the Generalized Pareto Distribution (GPD) on the peaks that exceed the 99th percentile storm surge level.Next, we derive estimates for various return periods (e.g., 10-, 25-, 50-, and 100-year).Following previous large-scale studies (Cid et al., 2015;Marcos & Woodworth, 2017;Pineau-Guillou et al., 2023;Vousdoukas et al., 2016Vousdoukas et al., , 2018;;Wahl et al., 2017), we ensure independence between storm events by using a 72 hr period for the de-clustering of the peaks (Vousdoukas et al., 2018).We consider this an appropriate threshold given that the duration of a storm duration is typically less than 3 days.Because our data have a 10-min temporal resolution, we define a minimum storm duration (defined as the time above the threshold) of 60 min.We use the Maximum Likelihood Estimation to fit the GPD parameters, however, we set a starting estimate of zero for the shape parameter to prevent overshooting due to outliers.Other commonly used extreme value distributions other than GPD were also fitted (Arns et al., 2013).This includes the Exponential distribution following POT method, and Gumbel and Generalized Extreme Value distribution following annual maxima method.We found little variance in the 1 in 10-year return levels among the various methods, which shows that the return period estimates provided by GPD are robust.

Data Analysis
As a next step, we evaluate the performance of the storm surges derived from forcing GTSM with the HighRe-sMIP ensemble against storm surges derived from forcing GTSM with reanalysis data (Section 3.1).We use the ERA5 climate reanalysis (Copernicus Climate Change Service (C3S), 2017; Hans Hersbach et al., 2020), which has hourly fields with a spatial resolution of 0.25° × 0.25° (∼31 km at the equator).The computed performance metrics include Pearson's correlation coefficient, Root-Mean-Squared Error (m), mean bias (m), and mean relative bias (%).Next, we assess the influence of climate change by analyzing the changes in the 1 in 10-year storm surge levels by comparing the future and present epoch against the historical epoch (Section 3.2).Changes are computed for each individual HighResMIP model, as well as the median over the multi-model ensemble.With such a small ensemble, mean values are easily distorted by outliers; for this reason, we focus on the median values.Moreover, we mainly focus on the 10-year return period (and not the 100-year return period) in order to minimize the uncertainty of the extreme value fit.We assess the inter-model spread by computing the standard deviation of the projected changes across the different climate models.We also assess the inter-model agreement on the sign of change by computing the number of ensemble members indicating an increase or decrease in the return period values.
For both the validation and analysis of changes, we use the reference regions of the IPCC AR6 report (Iturbide et al., 2020) to evaluate regional differences.The reference regions consist of 33 coastal regions that represent consistent regional climate features.

Validation of Storm Surges Based on the Historical Simulations
Figure 2a displays a global map of the storm surge level corresponding to 10-year return period for the present epoch , showing the median of the 5-member model ensemble.Storm surge levels exceeding 1.5 m are seen in coastal areas located in regions with a shallow continental shelf (depth < 100 m), this includes the North Sea (Northwest Europe), Hudson Bay (Canada), the Arctic Sea (Alaska), Gulf of Carpentaria (North Australia), east Patagonia (Argentina), and the Yellow Sea (China).Surge heights below 0.2 m are found in coastal areas where the shelf is narrow and steep, such as Chile or the U.S. east coast, and oceanic islands, such as Hawaiian or Canary Islands.Variations in surge height also reflects differences in climatology; where mid-latitude regions generally experience more variable and more stormy weather than equatorial regions (Merrifield et al., 2013;Rueda et al., 2017).In tropical regions the storm surges are caused by tropical cyclones which have a much smaller spatial extent, shorter duration, higher amplitude and lower probability than those produced by mid-latitude storms (Dullaart et al., 2021;von Storch & Woth, 2008).Figure 2a shows that the storm surge level corresponding to the 10-year return period is below 0.5 m in regions such as the southern Gulf Coast (U.S. and Mexico) and near Mozambique and Madagascar.These regions have seen very high storm surges in the past, but those storms have a low probability and are not well-represented by decadal simulations.The spread of the model ensemble, as indicated by the standard deviation in Figures 2b, is correlated to the value of the ensemble median.This means that the spread is large in regions with high storm surge levels and small in regions with low storm surge levels.
To assess the performance of the HighResMIP ensemble, we validate the storm surge levels corresponding to the 1, 10, and 100-year return period against the ERA5 climate reanalysis (Table 1).The multi-model median shows a good global agreement with ERA5 as indicated by a high correlation coefficient (r > 0.9) and a small bias (<0.1 m).The performance in coastal areas with storm surge levels exceeding 1 m is lower than the global average, but reasonably accurate with a relative bias of −5.9% the multi-model median.These aggregated values however Note.
The HighResMIP ensemble is compared against storm surge levels derived with the ERA5 reanalysis.For the bias, values are averaged across all output locations.We show the values for all output locations (S columns), as well as for output locations where storm surge levels higher than 1 m (S > 1 columns).This selection is based on the 1 in 100-year storm surge levels derived with ERA5 and includes approximately 25% of the locations that are typically locations located along the coast in shallow regions.

Table 1
Model Performance of the HighResMIP Ensemble for Storm Surge Levels Corresponding to a Return Period of 1, 10, and 100 Years mask considerable regional variations in model performance.Figures 2c and 2d show that the ensemble median has a clear spatial bias, with mostly negative values at the tropics and positive values at mid-to high-latitude areas.The overestimation compared to ERA5 is most profound in coastal regions along semi-enclosed seas that frequently experience high storm surges, which suggests localized weather effects or overestimated storm intensities could play a role.Figure S2 and Table S1 in Supporting Information S1 show the model performance for selected regions that regularly experience storm surge events, including North Sea (Northwest Europe), Hudson Bay (Canada), the Arctic Sea (Alaska), Gulf of Carpentaria (North Australia), east Patagonia (Argentina), and the Yellow Sea (China), we refer to Figure S2 and Table S1 in Supporting Information S1.Also, for these regions the average relative bias is reasonably small (within 10%), however, there are coastal locations with biases that exceed 30%.We further analyze regional difference in model performance in Figure 3, which shows the relative bias for the 1 in 10-year storm surge level aggregated to the IPCC reference regions.This clearly shows the spatial coherency of the bias with regions showing a tendency to over-or underestimate compared to ERA5.The spatial bias pattern is consistent for the 1, 10 and 100-year return periods, but also present for the bias for the 75th, 90th, and 95th percentiles of the surge level time-series (Figure S3 in Supporting Information S1).The fact that this bias is already present under not very extreme conditions suggest that systematic errors in the large-scale ocean and atmosphere fields could possibly play a role.
The ensemble median shows a better overall agreement with ERA5 than the individual models (Table 1), which is in favor of the use of an ensemble approach.Comparing the performance of the individual models, there are notable differences for the different return periods.The EC-Earth3P-HR and HadGEM3-GC31-HM have a comparable and relatively good performance, which is consistent for the different return periods (Table 1).The performance of HadGEM3-GC31-HM is worse for the atmosphere-only simulation than for the coupled simulation, with the mean relative bias doubling from 5% to 10% for the storm surge levels with a 10-year return period.While the atmosphere-only simulation may be less affected by a spatial SST bias, this could be counteracted by the fact that atmosphere-ocean coupling may result in a more realistic representation of extreme winds.The mean bias of CMCC-CM2-VHR4 is higher than for the other models, especially for the 100-year return period with a 0.14 m bias.GFDL-CM4C192-SST has a negative mean bias, which is most profound for the low return periods.
Although the magnitude differs, most of the individual HighResMIP models also show a clear spatial bias (Figure S4 in Supporting Information S1).Similarly, as for the ensemble median, this bias is not just caused by deficiencies in the extremes.For example, in GFDL-CM4C192-SST the normal climate conditions already introduce a negative bias, while extremes seem to be overestimated.Moreover, the distribution of the changes of the future epoch is asymmetrical with 26 of the 32 regions having a distribution that is highly skewed toward storm surge increase.While further research is needed, this suggest that the climate warming could possibly drive larger changes than currently observed.Generally, across larger spatial scales storm surge changes tend to cancel each other out as most regions see both increases and decreases.

Changes in
For almost all regions, the median change is around zero and the interquartile (25th-75th percentile) range of the changes is within 0.05 m.This means that for 50% of the output locations the projected changes are smaller than 0.05 m.An exception is for North and Central Australia where the interquartile range reaches up to an increase in surge level of 0.1 m.
The HighResMIP ensemble has considerable intermodel variability.Whilst the ensemble only contains five models, we quantify the consistency of the ensemble by showing the intermodel agreement on the sign of change.
Figures 6a and 6b display the number of models that project an increase or decrease in the 1 in 10-year storm surge levels for the present and future epoch, relative to the historical epoch.We find that there are a few regions where the majority of models (more than 3) agree on the sign of change.However, the changes relative to the historical epoch are more consistent across the ensemble for the future epoch than for the current period.   .The projected changes computed based on the individual HighResMIP models (and not the ensemble median).
based on SSP5-8.5.Regions where the ensemble consistently shows an increase in the 1 in 10-year storm surge levels include parts of the North Sea (Northwest Europe), Yellow Sea (Eastern China), Patagonia, and northern Australia.For the southern part of the Mediterranean the models agree on a decrease.For many other regions, such as the Persian Gulf (Saudi Arabi, Iran) and the Gulf Coast (southern United States and Mexico) there is no clear majority.Figures 6c and 6d shows the intermodel spread expressed as the standard deviation of the changes across the ensemble.As expected from the larger magnitude of changes, the spread is generally large in areas with high storm surges, like the North Sea and Gulf of Thailand.While having higher consistency, for the future epoch the intermodel spread is larger than for the present epoch.Figure S5 in Supporting Information S1 shows the projected changes of the individual members of the HighResMIP ensemble.
Because of the large uncertainties associated with EVA, we focus on the 10-year return period.Projected changes for the 25-, 50-, and 100-year return periods are presented in Figure S6 in Supporting Information S1.It shows that, although the magnitude of changes is slightly different and there are some outliers, the spatial pattern of increases and decreases is largely consistent for different return periods.

Understanding the HighResMIP Ensemble Biases
Our results show that the GCMs of the HighResMIP ensemble exhibit large biases in the simulation of storm surge in current climate conditions, contributing to the projection uncertainty for future changes.Understanding and addressing the biases is a key topic for future research.There are three important directions that could be explored.
First, moving forward, it is important to attribute the storm surge bias to underlying drivers in the GCMs (e.g., large-scale ocean and atmospheric processes, too many or too strong storms).The increase in model resolution from 150 km (CMIP5) to 25-50 km (HighResMIP) results in a more realistic representation of important physical processes.High-resolution GCMs have a better performance than low-resolution GCMs (Haarsma et al., 2016; M. J. Roberts et al., 2018), especially with regards to the representation of extremes that drive storm surges.For extra-tropical storms, CMIP5 models have been shown to underestimate the intensity (Zhang & Colle, 2018).The model resolution of HighResMIP compared to CMIP5, shows great improvements in storm structure, intensity and frequency (Catto et al., 2010;Priestley & Catto, 2022b;Zappa et al., 2013).For tropical cyclones, similar improvements are reported although at the resolution of HighResMIP the intensity of the strongest cyclones are still underestimated (Murakami et al., 2015;Roberts et al., 2020).Moreover, there are several long-standing climate biases that remain even at high resolution and that will require improvements in model physics (Moreno-Chamarro et al., 2022).Many climate biases are both depending on the specific climate model and region and it is not straightforward to evaluate how these reported biases will affects storm surge modeling at the global scale.While not extensively documented in literature, global storm surge projection based on CMIP5 models (Muis et al., 2020;Vousdoukas et al., 2018) are expected to exhibit similar or larger biases as the HighResMIP ensemble.
Second, a better understanding of the GCM biases could also help in developing appropriate bias correction methods.While statistical methods cannot address fundamental problems of poorly represented physical processes in GCMs, bias correction is commonly applied in climate impact studies (Lemos et al., 2020;Navarro-Racines et al., 2020).Here we minimized the effect of the bias by comparing the historical and future period for the same GCM (assuming the biases are time-invariant), also referred to as the delta-change method (Eisner et al., 2012;Hemer et al., 2013;Morim et al., 2019).The application of more advanced methods, such as quantile mapping (Lemos et al., 2020;Parker et al., 2023;Tebaldi & Knutti, 2007), could be explored in future research.The impact of the bias could be further reduced by computing a weighted ensemble mean where models that have a large bias in a specific region are given less weight (Marsooli et al., 2019a).
Third, while computationally demanding, additional GTSM simulations could give valuable insights.Currently we use a heterogenous ensemble that consists of a mix of 5 coupled and atmosphere-only simulations.The statistical uncertainty could be reduced by enlarging our ensemble, which can be done in two ways: (a) by adding more HighResMIP models to better assess the intermodel agreement; and (b) by adding additional members of the existing HighResMIP models to better assess the climate variability.The first would allow to investigate model disagreement as the result of atmosphere-ocean coupling, while the second could help in disentangling the signal of anthropogenic-forced climate change from the climate variability.An alternative direction is to extend the ERA5 simulation back to 1950 (Bell et al., 2021).This would make it possible to validate the storm surge trends from 1950 to now and evaluate which HighResMIP simulations are consistent with the reanalysis.This could be supplementation with a comparison against observed trends from GESLA tide gauge stations with long enough records (Haigh et al., 2022).

Comparison of Projected Changes With Previous Studies
Our results provide useful information for broad-scale studies focusing on future flood risks and climate adaptation, which generally do not account for changes in climate extremes.While sea-level rise will be the main driver of changes in extreme water levels (Kirezci et al., 2020;Pringle et al., 2021;Vousdoukas et al., 2018), changes in storm surge could contribute as well.The ensemble-median of our projections suggest larger changes than previously reported by global studies based on CMIP5 (Muis et al., 2020;Vousdoukas et al., 2018), even though we focus on the near-future and the end of century.In comparison with projected storm surge changes from regional studies, our projections are largely consistent.The storm surge decrease for the Mediterranean has also been reported by other studies (Makris et al., 2023;Nissen et al., 2014), as a result of a decrease in storm frequency.The increase at the mid-at high-latitudes in coastal regions in semi-enclosed seas could possibly be linked with a poleward shift of storm tracks and increased risk of post-tropical cyclones (Haarsma, 2021).In agreement with our results, Vousdoukas et al. (2016) also showed an increase in the 1 in 10-year storm surge level for the North Sea and Baltic Sea, and a decrease in the Mediterranean Sea.For Australia the largest increases along the Gulf of Carpentaria and southern Tasmania and decreases for the south coast are largely consistent with the findings from Colberg et al. (2019), although their projections show on average more decreases.For East Asia, we show storm surge increases may occur in the East China Sea and Yellow Sea.This is not fully in agreement with Yasuda et al. (2014) who projects decreases in the Yellow Sea, however, the reported decrease is based on a single climate model.For the Pacific coast of Canada four regional climate models project a slight but widespread increase in the 99th percentile of daily maximum storm surges between 2006 and 2025 to 2031-2050 for RCP8.5 (Cousineau & Murphy, 2022).The ensemble mean value of +5% is comparable to our results.The small increase to no change in 1 in 10-year storm surge level along the Atlantic coast of northeastern United States coastline is comparable to regional studies for extra-tropical storms (Lin et al., 2019;Pringle et al., 2021).For tropical cyclones in the same area, a comparison is difficult since studies tend to focus on higher return levels using thousands of stochastic events (Marsooli et al., 2019).

Limitations and Recommendations of the Used Modeling Approach
Our hydrodynamic modeling approach is based on GTSM, a depth-averaged barotropic model that is designed to globally model storm surges with high accuracy in coastal areas.The resolution of 2.5 km has shown to be able to accurately represent the storm surges of extra-tropical cyclones as well as tropical cyclones (Bloemendaal et al., 2018;Dullaart et al., 2020;Muis et al., 2019), although compared to local observations an underestimation of the highest extremes is expected.One aspect of our approach that could be improved is the ocean-atmosphere coupling of the hydrodynamic model.Currently, we use a constant value for the wind drag based on Charnock (1955), although the efficiency of the air-sea momentum transfer will depend on the sea state and the wind speed (Powell et al., 2003;Ridder et al., 2018).The implementation of a spatially and temporally varying wind drag parameterization, that should be consistent with the formulation for sea-air momentum transfer in the climate models, is a focus of future research.Another aspect of our modeling approach that could be improved is the fact that wind-waves are not included.Water level interactions can be significant (Arns et al., 2017;Wolf, 2008) and for coastal flooding it has been shown that the contribution of waves can be important (Bertin et al., 2014;Staneva et al., 2016).Moreover, waves are correlated with storm surges (Marcos et al., 2019) and future projections of waves show significant changes due to climate change (Meucci et al., 2020;Morim et al., 2019).For future research it would be valuable to extent our projections with consistent wind-wave projections also based on HighResMIP ensemble and analyze the combined effects of waves and storm surges on the frequency of coastal flooding (Kirezci et al., 2020;Parker et al., 2023;Vousdoukas et al., 2018).Finally, by calculating storm surges as the residual water levels (i.e., difference between tide-only and total water-level simulation), non-linear interactions are considered to be part of the storm surge level.This experimental setup was chosen because it has been shown that tide-surge interaction is important to consider at global scale (Arns et al., 2020).This is the case especially for shallow areas with a large tidal range, which are generally areas that also experience large storm surges.However, this setup also means that errors in tide simulations may propagate into the tide-surge interactions and consequently the storm surge levels.Given the reasonable performance of GTSM with regards to tidal amplitudes in coastal areas (Wang et al., 2021;Wang, Verlaan, Apecechea, & Lin, 2022;Wang, Verlaan, Veenstra, & Lin, 2022), this is expected to be a minor influence.However future research could aim to quantify the effect of this on the storm surge levels.
Our statistical modeling approach is based on EVA based on 30-year time periods.This is generally associated with large uncertainties (Wahl et al., 2017), which make it challenging to detect trends in extremes.Although we focus on the 10-year return period, climate variability may have a significant influence.Calafat et al. (2022) have provided observational evidence that, in Europe, climate variability can considerably affect the storm surge probabilities over periods as long as 60 years.The small sample size is especially problematic for regions prone to tropical cyclones, which will be partially masked by the EVA.Tropical cyclones generally have low probabilities that can only be robustly assessed with stochastic modeling of events representing thousands of years (Bloemendaal et al., 2020;Dullaart et al., 2021;Haigh et al., 2013;Lin & Emanuel, 2015;Marsooli et al., 2019;Orton et al., 2016).In future research the confidence of the projected changes could be improved by enlarging the sample size by adding more HighResMIP models, but also by exploring alternative statistical modeling approaches.Promising directions include non-stationary extreme value statistics that account for interannual climate variability (Rohmer et al., 2021), Bayesian modeling that enables sharing of information across space (Calafat & Marcos, 2020) and the pooling of extremes from an intermodel ensemble of climate models (Meucci et al., 2020).

Conclusions
Using the CMIP6 HighResMIP climate model simulations, we developed global multi-model projections of extreme sea levels (ESLs) from 1950 to 2050.The ∼25 km resolution of the CMIP6 HighResMIP ensemble represents a step-change compared to the previous coarse resolution CMIP5-based simulations that poorly represent climate extremes such as extra-tropical storms and tropical cyclones.Comparison of the 1 in 10-year storm surge levels derived from HighResMIP ensemble against those derived from ERA5 shows a good overall performance but also a clear spatial bias.For the same return period, the projected changes for the future epoch (2021-2050) compared to the historical epoch  are up to −0.1 m (or up to 20%).These changes are not uniform across the ensemble with large intermodel differences.The HighResMIP ensemble agrees upon a decrease along the coast of Mediterranean and northern Africa and southern Australia and an increase along the south coast of Australia and Alaska.For other regions, such as northern Caribbean, eastern Africa, China and the Korean peninsula, there is large intermodel disagreement.Overall, the projected storm surge changes are small compared to climate model bias, decadal variability, and statistical uncertainties, and future research is needed to further enhance the understanding of how ESL will change in response to anthropogenic forcing.

Figure 1 .
Figure 1.Flowchart of the modeling approach that was used to generate the global storm surge projections.

Figure 2 .
Figure 2. Global map of storm surge levels corresponding to a 10-year return period.Panel a shows the ensemble median of the HighResMIP models for present epoch (1985-2014), while panel b shows the standard deviations.Panel c and d show, respectively, the bias (m) and the relative bias (%) of multi-model median in comparison with the ERA5 reanalysis for the same epoch.

Figure 3 .
Figure 3. Relative bias (%) of storm surge levels corresponding to a 10-year return period between the HighResMIP ensemble median and ERA5, aggregated into the Intergovernmental Panel on Climate Change reference regions.Panel a shows the regions with the colors corresponding to the colors of the boxplots in panel (b) The extent of each box shows the interquartile range across locations within the region (25th-75th percentile), the vertical black line indicates the median value (50th percentile) and the extent of the whiskers indicates the range of the 5-95th percentiles.
Storm Surges for the HighResMIP EnsembleNext, we analyze the storm surge changes for the periods 1985-2014 (present epoch) and2021-2050 (future  epoch) relative to 1951-1980 (historical epoch).Figures4a and 4bshows that for most of the world's coastline the storm surge changes for the present epoch are relatively small.For approximately 10% of the output locations, the relative changes are larger than 5% and absolute changes are larger than 0.1 m.The largest increases in the 1 in 10-year storm surge levels for the present epoch occur in the southern part of the North Sea, the Gulf of Carpentaria (Australia), the Bering Sea (Russia and Alaska), and the South China Sea (Vietnam and China), while the largest decreases occur near the Gulf of St. Lawrence (Quebec, Canada) and the Yellow Sea (Northeast China and North-and South Korea).For the future epoch, the projected changes in surge levels under SSP5-8.5 are larger.For approximately 25% of the coastal output locations, the relative change is larger than 5% and the absolute changes are larger than 0.1 m.Compared to the present epoch, the change signal for the Bering Sea and Gulf of Carpentaria becomes more coherent, while for the North Sea and South China Sea the signal becomes more scattered.The largest and most coherent storm surge decrease occurs in the Mediterranean Sea.For some regions the sign of change flips when comparing the projected changes between the present and future epoch.This indicates that storm surge trends from 1950 to 2050 do not always have a consistent up-or downward trend but are possibly influenced by decadal climate variability.Next, we analyze the projected change in the 1 in 10-year storm surge levels aggregated to the IPCC-AR6 reference regions (Figure5).The changes for the future epoch are larger than for the present epoch in all regions.

Figure 4 .
Figure 4. Ensemble-median absolute (panels a and b) and relative (panels c and d) changes in storm surge levels for the 10-year return period.We show changes between the present epoch (1985-2014) and future epoch (2021-2050) compared to the historical epoch (1951-1980).

Figure 5 .
Figure 5. Statistics of projected changes in the storm surge levels corresponding to a 10-year return period between the future epoch (2021-2050) and the historical epoch, aggregated into the IPCC-AR6 regions.The colored boxes show the interquartile range (25th-75th percentile) across locations within the region, the vertical black line indicates the median value (50th percentile) and the extent of the whiskers indicates the range of the 5-95th percentiles.The color of the bars is corresponding with the map in Figure3a.The gray boxes show the interquartile range of changes between the present epoch and the historical epoch.The projected changes computed based on the individual HighResMIP models (and not the ensemble median).

Figure 6 .
Figure 6.Robustness of the projected changes across the model ensemble for the storm surge levels corresponding to a 10-year return period, showing the intermodel agreement (panel a and b) as expressed by the number of models agreeing upon the sign of change, and the model spread (panel c and d) as expressed by the standard deviation.We show the changes for the present epoch (1985-2014) and future epoch (2021-2050) compared to the historical epoch (1951-1980).