Main

Compound extreme weather events are hazardous and can produce higher impacts than individual or single-hazard events. Compound events may be classified into multivariate events, spatially compounding events and temporally compounding events1. Multivariate events result from multiple hazards occurring at the same time, while spatially compounding events occur when multiple locations in a region are impacted by hazards within a limited time window1. Temporally compounding events involve multiple hazardous events happening consecutively, and the later hazards can produce more damage due to the lowered resistance of infrastructure and communities. For example, tropical cyclones (TCs) and heatwaves can be temporally compounded. A TC may destroy the power system, leaving residents without air conditioning when a subsequent heatwave arrives2,3.

TCs themselves can generate multiple hazards: strong winds, heavy rainfall and storm surges, which can weaken coastal infrastructure4. These hazards can happen jointly (for example, rainfall and surge causing compound flooding5,6) and produce more damage. Previous studies have focused on the joint hazards from pairs of TC hazards5,6,7,8. A study of triple TC hazards (that is, wind, surge and rainfall) is needed to gain comprehensive understanding of coastal risk. Furthermore, previous studies have not addressed sequential hazards from TCs, which have recently impacted the United States. For example, Hurricane Ida (2021) and Hurricane Nicholas (2021) arrived in Louisiana within 15 days of each other (Supplementary Fig. 1), and the influence of Ida (for example, saturated soil conditions) led to a more severe impact from Nicholas9,10. The recent occurrence of such events is consistent with previous findings that the chance of TCs making sequential landfall is increasing11 and that TCs are becoming more hazardous6,12,13. However, it is still unclear whether the occurrence rate of sequential TC hazards in the historical period shows any trend, whether we can make robust projections of such events considering the uncertainties in TC frequency projection14 and what physical mechanisms cause the change of sequential TC hazards.

To address these questions, we investigated the change in sequential TC hazards using both historical observations and climate simulations. We chose the yearly minimal impact interval (MII; unit: days) between hazard-producing TCs as a metric to describe sequential TC hazards. We used the 95th percentiles of TC daily maximum water level, daily total rainfall and daily maximum wind as thresholds for the three hazards, respectively, and we defined storms with at least one hazard component exceeding its threshold to be hazard-producing. The MII is the yearly minimum of the intervals between the hazards of two sequential TCs (the results are not sensitive to the selected hazard thresholds). We focused on the sequential impact from storm hazards at the same location (temporally compounding events). We also present a brief analysis of the sequential impact from severe TCs on the US mainland (spatially compounding extremes). As an additional reference, we analysed marginal and joint TC wind, surge and rainfall hazards (multivariate events) to help us understand the cause of changes in sequential TC hazards.

Detecting the historical trend of sequential TC hazards directly from observations is challenging due to data limitations11. To circumvent the challenge, we expanded a probabilistic model that describes sequential TC landfalls11 and fit the model using historical observations of TC hazard events (Methods). We then used the model to perform Monte Carlo simulations for each year from 1949 to 2018 to generate a large sample of sequential TC hazard events at nine selected locations along the US East and Gulf coasts and estimated the historical trend of sequential TC hazards by calculating their annual occurrence probability from the simulations (Methods).

For the future projection (2070–2100), we applied a physics-based TC hazard analysis method to investigate sequential TC hazards impacting the US Gulf and East coasts under both the high (Shared Socioeconomic Pathway 5 8.5 (SSP5 8.5)) and moderate (SSP2 4.5) emission scenarios. For each scenario, we used synthetic TCs generated from a statistical–deterministic TC model6,15 forced by six Coupled Model Intercomparison Project Phase 6 (CMIP6) climate models. For each storm, we simulated storm tides (storm surge plus astronomical tide) with the advanced circulation (ADCIRC) hydrodynamic model16, rainfall with the physics-based TC rainfall (TCR) model17,18 and ten-minute sustained wind using the complete wind profile model19. To evaluate the impact of sea-level rise (SLR) on sequential TC hazard events, we incorporated the probabilistic SLR projection for each emission scenario from the Intergovernmental Panel on Climate Change (IPCC) six assessment20. The abovementioned probabilistic model was then fitted using the simulated TC hazards to investigate the change of the return period (that is, the reciprocal of the annual exceedance probability) of various MII levels (Methods).

Increased hazard-producing storms from 1949 to 2018

The scarcity of sequential TC hazard events in the observations prevents the direct detection of any trend. However, several parameters that influence sequential TC hazards can be analysed on the basis of observations, including TC hazard frequency and duration. We analysed the trends in hazard frequency and duration for the nine coastal locations (Supplementary Fig. 2). Except for Charleston, South Carolina, and Pensacola, Florida, the hazard frequency has increased since 1949 (Supplementary Fig. 3), due to the increase in the hazard-producing capability of TCs (Supplementary Fig. 4) rather than their landfall frequency (Supplementary Fig. 5). The increased hazard-producing ability is probably a result of increased TC intensity21. Except for Charleston, South Carolina, the hazard duration has increased since 1949 (Supplementary Fig. 6), which may be a result of decreased TC translation speed22,23.

We estimated the yearly probability of experiencing sequential TC hazards using the probabilistic model with the hazard frequency and duration parameters fitted for each year (ten-year moving average; Fig. 1). The highest annual probability of experiencing sequential TC hazards with MII < 30 days is less than 10%, implying the rareness of such events. However, except for Charleston, South Carolina, and Pensacola, Florida, there is a clear increasing trend in the yearly probability of sequential TC hazards. For example, the probability doubled over the past seven decades at Sandy Hook, New Jersey.

Fig. 1: Estimated yearly probability of experiencing sequential TC hazards.
figure 1

ai, The results are shown for three thresholds of MII (10 days, 15 days and 30 days) for the nine selected US coastal locations. The dashed lines are fitted linear trends of the probabilities. An asterisk in the key indicates that the trend is significant.

Source data

Projected increase of sequential TC hazards

The climate projection shows that the return period of the MII of sequential TC hazards will substantially decrease along the US East and Gulf coasts by 2100 (Fig. 2). For example, the MII being smaller than 15 days is a 40-year event in the control climate for Texas, but under the SSP5 8.5 (SSP2 4.5) scenario, if we don’t consider SLR (or in practice if we adapt to SLR), such an event will become a 4-year (6-year) event; if we consider SLR, it will become a 2-year (3-year) event. Across the US East and Gulf coasts, under the SSP5 8.5 (SSP2 4.5) scenario, the return period of 15-day MII will decrease from 10–92 years to ~2–5 (2–11) years not considering SLR and to ~1–2 (1–3) years considering SLR.

Fig. 2: Return period of MII in the climate projection.
figure 2

ah, In each plot, the black curve indicates the return level of MII in the control period, the red solid (dashed) curve indicates the return level of MII in the SSP5 8.5 (SSP2 4.5) scenario without considering SLR, the pink solid (dashed) curve indicates the return level of MII in the SSP5 8.5 (SSP2 4.5) scenario with consideration of SLR, the blue solid (dashed) curve indicates the return level of MII in the SSP5 8.5 (SSP2 4.5) scenario without considering SLR or storm landfall frequency change, and the brown solid (dashed) curve indicates the return level of MII in the SSP5 8.5 (SSP2 4.5) scenario considering SLR but not considering storm landfall frequency change. The curves are averaged from six climate models (CanESM5, CNRM-CM6-1, UKESM1-0-LL, EC-Earth3, IPSL-CM6A-LR and MIROC6). The results were obtained under the 95th-percentile definition of hazard-producing; the results under the 90th-percentile and 99th-percentile definitions are shown in Supplementary Figs. 7 and 8, respectively. The return period for west Florida not considering SLR or storm frequency change (blue curves) is slightly lower in SSP2 4.5 than in SSP5 8.5, mainly due to higher surge-producing ability of the storms in SSP2 4.5.

Source data

The estimated change in MII between TC hazards was compared with the change in MII between TC landfalls (as in ref. 11). For example, for the 15-day MII return level under the SSP5 8.5 (SSP2 4.5) scenario, the decrease of return period defined by landfall ranges from 50.6% (48.6%) to 78.7% (74.0%), while the decrease defined by hazard ranges from 88.7% (86.6%) to 98.1% (97.1%) considering SLR and 84.7% (81.7%) to 94.4% (90.2%) not considering SLR (Supplementary Tables 1 and 2). The decrease in the return period of MII defined by landfall is caused by the projected increase in landfall frequency11, while the decrease in the return period of MII defined by hazard is due to both increased landfall frequency and increased storm hazard-producing ability.

The projected increase of landfall frequency is a major driver of the shortened MII, but it has large uncertainties14. However, a sensitivity test using future storms with the control frequency shows that the return period of MII will still decrease significantly. For example, in Texas under the SSP5 8.5 (SSP2 4.5) scenario, an MII of less than 15 days would change from a 40-year event to a 22-year (28-year) event not considering SLR and a 9-year (11-year) event considering SLR. Across the US East and Gulf coasts, under the SSP5 8.5 (SSP2 4.5) scenario, the return period of the 15-day MII will decrease from 10–92 years to ~8–46 (7–70) years not considering SLR and to ~3–10 (4–12) years considering SLR. This sensitivity test implies that the enhanced hazard severity12,13, caused mainly by increased TC intensity24,25,26 and SLR, significantly increases the frequency of sequential TC hazards, which is not contingent on there being an increase in the overall frequency of TCs. The projected change in the MII return period solely due to increased hazard severity is comparable to the change caused by increased landfall frequency in our projection (Supplementary Table 1).

Mechanisms for the increase of sequential TC hazards

Sequential TC hazards were shown to significantly increase in the future due to the increased severity of TC hazards (Fig. 2). Individual and joint TC hazards may become more severe and last longer. From the control (1984–2005) to the future climate (2070–2100), the return period of joint hazards will substantially decrease for all US coastal regions (Supplementary Fig. 9). Meanwhile, the proportion of landfalling TCs that can produce joint hazards (Fig. 3) will increase from 3.6–12.9% to 7.92–26.5% (13.6–30.6% with SLR) under the SSP5 8.5 scenario and to 6.0–24.4% (10.0–28.3% with SLR) under the SSP2 4.5 scenario. The return period of TCs that can produce at least one type of hazard will also substantially decrease across the United States (Supplementary Fig. 10). The proportion of landfalling storms producing at least one hazard will increase from 25.4–43.7% to 35.9–54.0% (80.1–85.3% with SLR) under the SSP5 8.5 scenario and to 27.9–53.6% (75.8–82.9% with SLR) under the SSP2 4.5 scenario (Fig. 3).

Fig. 3: The ratio between hazard-producing and landfalling storms.
figure 3

ap, The ratios between single-hazard-producing, at-least-one-hazard-producing and joint-hazard-producing TCs (λhazard) and landfalling TCs (λlandfall) are shown. For each location, the upper plot shows the SSP5 8.5 scenario (ad,il), and the lower plot shows the SSP2 4.5 scenario (eh,mp). Each box plot shows the uncertainty among the six climate models. The azure boxes indicate the control scenario, the orange boxes indicate the future scenario without considering SLR and the red boxes indicate the future scenario with consideration of SLR. In each plot, the box extends from the first quartile to the third quartile of the data, with a line at the median. The whiskers extend from the box by 1.5× the interquartile range. Because SLR does not influence rainfall and wind in the simulation, their ratios in future scenarios with SLR are omitted. The results were obtained under the 95th-percentile definition of hazard-producing; the results under the 90th-percentile and 99th-percentile definitions are shown in Supplementary Figs. 12 and 13, respectively.

Source data

To understand which hazard component of TCs drives the change in sequential events, the ratio between single-hazard-producing TC frequency and landfall TC frequency was examined (Fig. 3). In the control simulation, surge has the highest hazard-to-landfall ratio (ranging from 16.9% to 37.6%). In the future, if SLR is not considered, the leading hazard component switches to rainfall (except in Mississippi–Alabama under the SSP2 4.5 scenario). For example, in Louisiana, the ratios for surge, rain and wind hazard are 35.0%, 29.2% and 24.2% in the control climate and 41.7% (40.5%), 48.9% (43.2%) and 36.8% (34.8%) in the SSP5 8.5 (SSP2 4.5) climate, respectively. The relative increase of the ratio of rainfall-hazard-producing TCs ranges from 67.5% to 125.4% (45.9% to 59.2%) across the coastal regions under the SSP5 8.5 (SSP2 4.5) scenario (Supplementary Tables 3 and 4). The increase of rainfall hazard is a result of the combination of TC climatology change (for example, an increase in intensity and a decrease in translational speed) and the increase of atmospheric water vapour content27,28.

If SLR is considered, the coastal extreme water level will be the leading hazard, and the ratio between surge-hazard-producing storms and landfalling storms will reach 79.3–84.7% (74.5–82.5%) under the SSP5 8.5 (SSP2 4.5) scenario (Fig. 3). The relative increase in the rate of surge-hazard-producing TCs ranges from 126.1% to 379.1% (113.6–352.9%), and the relative increase in the rate of hazard-producing TCs ranges from 93.6% to 223.0% (83.9% to 215.9%) across the coastal regions (Supplementary Table 2). These different increases imply that weak TCs will be capable of producing high peak water levels under SLR, though these TCs may not be intense enough to produce extreme rainfall and wind. In addition to the increased rates of hazard-producing TCs in the future climate, the individual hazard levels will shift to higher values (Supplementary Fig. 11).

The differences between the SSP2 4.5 and SSP5 8.5 scenarios in the change in hazard-producing storm ratios are relatively small (Fig. 3 and Supplementary Tables 3 and 4), especially for the Gulf Coast and for surge-producing storm ratios (not considering SLR). In Texas and West Florida, the surge-producing storm ratio is even larger under the SSP2 4.5 scenario (median 24.5% and 17.7%) than under the SSP5 8.5 scenario (18.6% and 14.3%). A possible reason is that TCs with higher intensity (which are more frequent under SSP5 8.5) are associated with a smaller radius of maximum wind (given the climate-invariant outer size as assumed in this study) and therefore may impact smaller areas. The hazard-producing storm ratio considering SLR is around 80% along the Gulf and East coasts in both scenarios, with SSP5 8.5 being slightly higher (Fig. 3).

Besides making weak TCs hazard-producing, SLR also increases hazard duration (Fig. 4). Compared with the projection without considering SLR, the projection considering SLR shows a drastic increase in the average hazard duration. The increase in hazard duration caused by SLR indicates that, due to SLR, TCs can be hazardous even when they are relatively weak or distant, especially under the SSP5 8.5 scenario. Longer hazard duration under the SSP5 8.5 scenario is another reason that the occurrence rate of sequential TC hazards is higher in the SSP5 8.5 scenario.

Fig. 4: Exceedance probability of TC hazard duration.
figure 4

ah, The hazard duration is defined as the number of days when at least one hazard component occurs during a storm. The exceedance probability is defined as one minus the cumulative probability function. The mean of the hazard duration is presented in the key. An asterisk in the key indicates that the change in the mean is significant compared with the control simulation. The solid black curve represents the control simulation, the solid (dashed) blue curve represents the SSP5 8.5 (SSP2 4.5) simulation and the solid (dashed) red curve represents the SSP5 8.5 + SLR (SSP2 4.5 + SLR) simulation. The results were obtained under the 95th-percentile definition of hazard-producing; the results under the 90th-percentile and 99th-percentile definitions are shown in Supplementary Figs. 14 and 15, respectively.

Source data

Previous research29,30 argued that the increase in peak storm surge caused by storm climatology change may be comparable to or dominant over the effect of SLR at some coastal locations, which seems contradictory to the findings of this study. However, previous research focused on extreme events, such as 100-year events. Here, we define lower hazard thresholds that are not extreme but can still cause flooding and impacts. The differences between the findings of this study and those of previous research imply that SLR is more capable of changing weak TCs from non-surge-producing to surge-producing than of increasing the chances of the most extreme surge events.

Grey swan sequential hazard-producing storms

In previous sections, we examined the probability of two hazard events sequentially impacting the same location. Here, we consider the possibility of two extreme “grey swan” events31 sequentially impacting the United States, in which case dispatching limited rescue resources to the affected areas may be difficult (e.g., Hurricanes Harvey, Irma and Maria in 2017). As an example, we investigated the chance of a Katrina-like storm (causing a water level > 8 m in at least one coastal location) and a Harvey-like storm (causing total rainfall > 1,000 mm in at least one coastal location) impacting the US mainland sequentially (for example, with an impact interval < 15 days). Such sequential events cannot be found in the 1,375-year control simulation. However, such events may occur in the future. The return period of ‘Katrina’ and ‘Harvey’ impacting the US coastline within 15 days is around 250 (1,300) years without considering SLR in the future climate and 85 (650) years considering SLR under SSP5 8.5 (SSP2 4.5) (Fig. 5). The chance of occurrence of extreme sequential events is much higher in the SSP5 8.5 scenario than in the SSP2 4.5 scenario, though the difference is relatively small for general sequential events (Fig. 2).

Fig. 5: Sequential occurrence of Katrina-like and Harvey-like storms.
figure 5

a, Return period of the MII of ‘Katrina’ and ‘Harvey’ impacting the US coastline sequentially. The black curve (the return period of the control simulation) does not exist because no such sequential event can be identified in the control simulation. b,c, Selected storm track of ‘Katrina’ (Storm 4562 in IPSL under SSP5 8.5) (b) and ‘Harvey’ (Storm 4472 in IPSL under SSP5 8.5) (c) from the future synthetic dataset; the two representative storms were chosen because they have landfall locations similar to those of historical Katrina (2005) and Harvey (2017). d, The water level caused by the selected ‘Katrina’. e, Total rainfall caused by the selected ‘Harvey’.

Source data

Discussion

This study demonstrates that the chances of sequential TC hazards along the US East and Gulf coasts have increased and may continue to increase. The projected increase is largely driven by the increased hazard-producing ability of TCs, especially for rainfall and surge hazard (also considering SLR). Previous studies have projected average TC rain rates to increase by 10–32%6,27,28,32,33 by 2100. The increase of TC rain rates can drive both coastal rainfall extremes6,12,34 and the proportion of storms that become hazard-producing to increase, which increases the chance of sequential TC hazards. Although TC climatology change dominates the change in extreme surge and compound flood hazard6, SLR has a stronger impact on the change in moderate events. SLR is the leading factor that causes the increase of hazard-producing TCs to be around 80% of landfalling TCs (under both the SSP5 8.5 and SSP2 4.5 scenarios). As SLR projections have uncertainties35,36, their role in causing weak TCs to become hazard-producing may be further investigated.

This study compared the increase of sequential TC hazards under high and moderate climate change scenarios. Previous research on TC hazards mainly focused on the high emission scenario5,11,13,29,30,31,34, but given current policies, this scenario is considered unlikely to happen, whereas the moderate emission scenario is considered likely37. Our analysis shows that although there are some differences in terms of single or joint hazards (Supplementary Figs. 911), the return period of sequential TC hazards is similar under the two scenarios (for example, an 89–98% versus 87–97% decrease from the control climate for MII < 15 days). These similarities are not very sensitive to the uncertainties in SLR projection and storm frequency projection, as the change in storm hazard-producing ability is substantial under SSP2 4.5 (although it is more severe under SSP5 8.5). However, SSP5 8.5 will see more sequential extreme events, as demonstrated by the grey swan sequential hazard analysis. Overall, the results indicate that the main benefit of adopting the SSP2 4.5 pathway over the SSP5 8.5 pathway is to reduce the likelihood of extreme events.

The findings in this study are inevitably subject to the uncertainties of the observational data and climate downscaling methods. The data quality of TC hazard observations in earlier decades is always a concern38. However, a clear trend of increased probability of sequential hazard-producing TCs from 1979 (the satellite era) could be found for most locations (Fig. 1), indicating that the possible missing TC observations in early decades38 would not substantially alter the conclusions of this study. Also, the general increase of TC intensity21 and the decrease of TC translation speed22,23 in the past few decades supports the conclusion of increased hazard-producing capability of landfalling TCs. For the climate projection, the major concern was that large uncertainties exist in TC frequency projection. Our synthetic TC model projects increased TC frequency25, which contradicts several previous studies projecting the opposite14. However, the influence of increased TC hazard-producing capability is dominant over (when considering SLR) or comparable to (when not considering SLR) the influence of the projected increase in landfall frequency on the increase of sequential TC hazards.

This research provides the following takeaway messages. First, this research urges the consideration of back-to-back TC impacts in the development of resilience strategies. Second, TC rainfall and SLR have substantial influence on sequential TC hazard events; thus, the resilience of coastal infrastructure should be upgraded, targeting future extreme rainfall, prolonged surge and flooding. Finally, compounding between grey swan extreme events (even when they impact different locations) will stretch the emergency response systems in unprecedented ways, so preparation is needed.

Methods

Historical observations

TCs that made landfall in the United States from 1949 to 2018 were used to analyse the historical trend of the annual occurrence probability of sequential hazard events (the probability of at least one pair of sequential events happening in a given year). The information on historical TCs was obtained from the International Best Track Archive for Climate Stewardship39, which provides six-hourly TC locations and intensities. To quantify and examine the change in impacts from historical landfalling TCs, we obtained observations of hourly water levels, daily rainfall and daily maximum wind at nine locations across the US coastal areas (the locations are shown in Supplementary Fig. 1) from the Center for Operational Oceanographic Products and Services and the National Centers for Environmental Information. If the tidal gauge sites had no information for rainfall and wind, we found the closest weather stations (within 100 km) to the tidal gauge location and used the rainfall and wind observations from these stations to represent the wind and rainfall hazards at the tidal gauge station.

We found that there were more missing data for wind observations from 1949–1979 than in later decades, while the missing data for surge and rainfall observations had no substantial decadal variations. To eliminate the possibility that the trends were caused by the more frequent missing wind observations in the earlier decades, we conducted a similar analysis for only surge and rainfall hazards, and we reached the same conclusion. The fact that using only surge and rainfall hazards can obtain a similar conclusion as using the triple hazards does not imply that the wind hazard is not important; rather, it implies that the three hazards are physically correlated.

Synthetic TCs

Synthetic TCs generated with a statistical–deterministic TC model15,25 were used to study the change in sequential hazard-producing TCs in the United States. The storms were generated under the environment of six CMIP6 climate models (CanESM5, CNRM-CM6-1, UKESM1-0-LL, EC-Earth3, IP-SL-CM6A-LR and MIROC6) for the control (1984–2005), SSP5 8.5 (2070–2100) and SSP2 4.5 scenarios (2070–2100); 4,400, 6,200 and 6,200 US landfalling storms were generated for each model under each scenario, respectively. The synthetic storms for the control and SSP5 8.5 scenarios were the same as in ref. 6. Neither the difference in the simulation period nor the difference in the number of synthetic storms in the datasets influences our analysis, as the storm frequency for each of the climate states (1984–2005 and two scenarios for 2070–2100) is separately determined in the TC model and used in the analysis. In total, 8,250 years of simulation were performed for the 1984–2005 period, and 5,821 (5,388) years of simulation were performed for the 2070–2100 period under the SSP2 4.5 (SSP5 8.5) scenario. The thousand-year simulations enable us to obtain confident estimation of extremes, and the larger dataset of the control climate provides an even more reliable baseline.

The synthetic storm model assumes that the storms in the same year are conditionally independent (given the same environmental forcing), as the feedback of storms to the environment is not captured in the one-way coupled TC modelling system. Previous studies11 have examined the independence assumption of TC genesis and landfall and find that the assumption is valid for hazard analysis.

Wind modelling

The complete wind profile model that merges the TC wind profiles of the inner core and outer radii (C15; ref. 19) was used to simulate the wind hazard associated with the synthetic TCs. To obtain the surface wind affected by environmental winds, we added a correction to the simulated wind profile following ref. 40. The C15 model was also used to prepare the required wind input for the TCR simulation.

Hydrodynamical modelling

The ADCIRC model was used to simulate the storm tides produced by the synthetic TCs. The unstructured computational mesh developed in ref. 41, which has a spatial coverage of the entire North Atlantic basin, was used in this study. Eight tidal constituents42 were applied as boundary conditions at the ocean boundary of the mesh. The wind and pressure fields associated with the synthetic TCs were needed for the storm tide simulations and were obtained using physics-based parametric models43. This wind model was used in ref. 6 to drive surge simulation, as it is simpler and has similar performance to the C15 wind model. Further details regarding the ADCIRC model and simulation setups can be found in ref. 41. In the main text, we call the extreme water level (storm tide + possible SLR) related to TCs ‘surge hazard’, though the word ‘surge’ itself means the increase in water level due to atmospheric forcing.

Rainfall modelling

The physics-based TCR model12 was used to simulate rainfall associated with the synthetic TCs. Detailed formulation on the TCR can be found in ref. 18. The simulation setup of the model followed ref. 44, and the environmental field needed for the TCR simulation was obtained following ref. 12. The C15 model was used to drive the TCR simulation following ref. 6.

SLR projection

Localized probabilistic SLR projections in ref. 20 under the SSP5 8.5 and SSP2 4.5 emission scenarios were incorporated in this analysis. The projection of SLR was developed for tidal gauge locations; for each point along the coastline, we applied the projected SLR at the nearest tide gauge location. We followed ref. 45 to sample SLR time series over 2070–2100 from the SLR projection, and we sampled TC hazard events over 2070–2100 from the generated TC hazard datasets (assuming Poisson arrivals of storms) and combined the TC hazard events and SLR for each year.

Joint hazard analysis

Statistical analyses of the triple hazard including the modelled maximum storm maximum water level (L), maximum daily rainfall (R) and maximum wind speed (W) across each coastal segment of Texas, Louisiana, Mississippi–Alabama, West Florida, East Florida, Georgia, South Carolina and North Carolina11 were performed. The marginal distributions of rainfall, surge and wind were fitted by generalized Pareto distributions to characterize the long tails that corresponded to the extreme events6,29. The generalized Pareto distributions were fitted at each coastal segment, and the threshold was set by minimizing the mean squared error between empirical quantiles and the theoretical quantiles5. There was no parametric probabilistic distribution that described the trivariate generalized Pareto distributions, so nested Gumbel copulas46 were used to represent the dependent structure of the three individual hazards. Gumbel copulas were used here because previous studies6,7,47 showed that each pair of hazards (surge–rainfall, surge–wind and wind–rainfall) are correlated especially at tails, and Gumbel copulas are often used to quantify the tail-dependent structure48.

After fitting the marginal and joint distributions of the triple TC hazard, we calculated the probability of the three hazards jointly exceeding their respective thresholds (joint exceedance probability (JEP)) and the probability of at least one of the three hazards exceeding its respective threshold (at least one exceedance probability (OEP)). The thresholds LT, RT and WT are the marginal return levels of surge, rainfall and wind hazard with a return period of T years in the control simulation. The JEP and OEP are thus mathematically defined as:

$${\mathrm{JEP}}(T) = {\Bbb P}\left( {L > L_T \cap R > R_T \cap W > W_T} \right)$$
(1)

and:

$${\mathrm{OEP}}(T) = 1 - {\Bbb P}(L \le L_T \cap R \le R_T \cap W \le W_T)$$
(2)

To explicitly discuss the effect of SLR, we calculated the JEP and OEP with and without the consideration of SLR. To explore the causation of the change in JEP and OEP, we also investigated the changes in the exceedance probabilities of the single hazards (the marginal cumulative density function). The return period (TJE) of the three hazards jointly exceeding their respective return level thresholds (under return level T) can be calculated as:

$$T_{{\mathrm{JE}}}(T) = \frac{1}{{\lambda \times {\mathrm{JEP}}(T)}}$$
(3)

where λ is the arrival rate of the storms.

Similarly, the return period (TOE) of at least one of the three hazards exceeding its respective threshold can be calculated as:

$$T_{{\mathrm{OE}}}(T) = \frac{1}{{\lambda \times {\mathrm{OEP}}(T)}}$$
(4)

Probabilistic model for sequential hazard-producing TCs

The Poisson–Gaussian model for sequential landfalling TCs developed in ref. 11 was extended to capture the storms’ hazard-producing abilities and hazard durations. The arrival of hazard-producing storms was modelled as a non-stationary Poisson process, with arrival rate ν:

$$\nu \left( {t,s} \right) = \lambda _{{\mathrm{hazard}}}\left( t \right)S_{{\mathrm{hazard}}}\left( s \right)$$
(5)
$$S_{{\mathrm{hazard}}}\left( s \right) = \frac{1}{{{{{\mathrm{\sigma }}}}\sqrt {2{{{\mathrm{\pi }}}}} }}{{{\mathrm{e}}}}^{ - \frac{1}{2}\left( {\frac{{{{{\mathrm{s}}}} - \mu }}{{{{\mathrm{\sigma }}}}}} \right)^2}$$
(6)

where λhazard(t) is the annual frequency of hazard-producing TCs in year t. Shazard(s), representing the seasonal variation, is the likelihood of a hazard starting to occur on day s and μ is the mean and σ is the standard deviation of the impact time of storms that make landfall at the selected location. The difference between the Poisson–Gaussian model in this study and that in ref. 11 is that this study uses the annual frequency and seasonality of hazard-producing storms instead of those of landfalling storms, and the ratio between the annual frequency of hazard-producing storms and that of landfalling storms can be viewed as a metric of storm hazard-producing ability. The Poisson–Gaussian model was shown to be capable of capturing the relationship between TC landfall climatology and the minimal landfall intervals11, and a larger frequency (λhazard) and a smaller seasonal span (σ) favour sequential hazard-producing events.

We define the MII mathematically as:

$${{{\mathrm{MII}}}}({{t}}) = \mathop {{\min }}\limits_{i > j} \left( {{{{\mathrm{max}}}}\left( {B_i\left( {{t}} \right) - E_j\left( {{t}} \right),0} \right)} \right)$$
(7)

where Bi(t) is the hazard beginning time of the ith TC in year t, and Ej(t) is the hazard end time of the jth TC in year t; if the impacts of the two TCs overlap, we set the MII equal to 0. The TCs in a single year were ordered by the starting times of hazards associated with the TCs. In the Monte Carlo simulations, for each individual TC, the beginning times of hazards were randomly drawn from the probabilistic model of equations (5) and (6), and the beginning and end times of hazards were connected as:

$$E_i\left( {{t}} \right) = B_i\left( {{t}} \right) + D_i({{t}})$$
(8)

where Di(t) is the total hazard duration of the ith TC in year t. As the impact duration of TC-related hazards may change in the future, a probability distribution that describes the duration of TC hazard impacts at each selected location is needed. The probability distributions of duration in different coastal locations under different climate scenarios do not share the same parametric probability distribution. Thus, in our analysis and simulation, we applied kernel density estimation to fit the non-parametric probability distribution to the duration. For the climate simulation analysis, the physical modeling of each hazard component was performed before the probabilistic model was fitted with the simulated synthetic events.

Definition of hazard-producing TCs and hazard duration

There is no universal definition of ‘hazard-producing’ as different infrastructure systems and communities may respond differently to sequential TC events11. A statistically reasonable, physically meaningful and engineeringly applicable definition of ‘hazard-producing’ should be able to both categorize the landfalling storms into ‘hazard-producing’ and ‘non-hazard-producing’ categories and separate the days of impact from a single ‘hazard-producing’ TC into ‘hazard days’ and ‘non-hazard days’ (such as the days when TCs are too weak or too far away from the point of interest to produce hazards).

To do so, we used specific percentiles of daily maximum water level, total rainfall and maximum wind speed when a TC is within 250 km of the study point as the thresholds to define ‘hazard-producing’ TCs. In the historical analysis, for each individual point of interest, the daily maximum water level, total rainfall and maximum wind from every TC that ever approached 250 km from the location were collected, and the percentiles of each hazard component were calculated on the basis of the observations of all TCs impacting this location. For the climate simulation using CMIP6 models and the synthetic storm model, the same method was applied to calculate percentiles for the control simulation (historical period) of each climate model (so the percentiles differ between climate models). We used the thresholds in the control period instead of using the thresholds in the future period to better investigate the changes in TC hazard in the future relative to the historical hazard levels.

The specific percentile chosen as the threshold should be both high enough to eliminate nuisance TC events and low enough to include some non-extreme events that are still hazardous. In this main text of this study, the 95th percentile of each hazard was chosen as the threshold, and the days when at least one hazard component exceeded the threshold were defined as ‘hazard days’. The threshold was spatially varied by using the 95th percentiles specific to each coastal location to account for spatial variation in the preparedness for or awareness of hazards. The storms that caused at least one hazard day for a point of interest were defined as ‘hazard-producing’ storms for this location. The selection of the 95th percentile for the threshold fulfils the requirements mentioned in the beginning of this paragraph. For example, the mean of the 95th percentiles of six climate model simulations of daily accumulated TC rainfall (maximum tide level) in the control climate is 98.8 mm (1.4 m) near New Orleans, which is approximately ¼ (1/6) of the total rainfall (maximum tide level) that Hurricane Katrina produced in this location in 2005. Admittedly, the selection of the 95th percentile was ad hoc, so the results of climate simulations of sequential TC-related hazards under other thresholds (the 90th and 99th percentiles) are shown in Supplementary Figs. 7, 8 and 1215 as a sensitivity test of the results of this study.