Monday, November 25, 2024
No menu items!
HomeNatureMortality caused by tropical cyclones in the United States

Mortality caused by tropical cyclones in the United States

Data

Wind speed data

We measure TC incidence using the LICRICE model (version 4)28. LICRICE is a parametric wind-field model that estimates the maximum sustained winds experienced by every location throughout the lifetime of each TC recorded in the International Best Track Archive for Climate Stewardship (IBTrACS) database2,29,31,37,53,54,55,56. LICRICE uses observed maximum wind speeds to reconstruct wind fields throughout the storm based on internal storm structure, location, and each storm’s observed translational velocity. Many summaries of storm wind incidence are possible to construct using LICRICE, such as integrated power dissipation28,29; however, prior analyses have shown that the maximum wind experienced within storms is the most predictive of social and economic outcomes among numerous parsimonious metrics previously analysed. This result is consistent with most components of built structures failing catastrophically based on whether or not a threshold of stress is applied. LICRICE does not explicitly model storm surge, rainfall or flooding, however these dimensions of impact are captured in the analysis to the extent that they are correlated with integrated maximum wind speeds. For example, we find that LICRICE wind speed and NCEP rainfall from TCs are correlated within storms at the pixel level (Supplementary Fig. 3) and at the state-by-storm level (R2 = 0.31, P < 0.001; Fig. 1e) for a limited sample of storms for which granular rainfall data are available. Iterations of the LICRICE model has been used to measure various social and economic impacts of TCs, including direct deaths and damages29, changes to household income and expenditure2, infant mortality2, GDP growth28,37 and depreciation3. We also compare our measure of wind speed against total direct national economic damages (normalized by GDP) from a limited sample of TCs estimated by Nordhaus51. We find that national wind speed exposure is a meaningful predictor of total national damages (Fig. 1d and Supplementary Fig. 2), although we note that this outcome is highly uncertain and widely understood to be biased.

Here we use a new reconstruction of incidence at the sub-national-level within CONUS. We reconstruct incidence for 0.125° × 0.125° pixel of CONUS in each of 1,230 Atlantic storms between 1930 and 2015. Supplementary Fig. 1 shows all decadal averages of these output (four example maps are also shown in Fig. 1a), illustrating the TC climatology for CONUS–however these aggregates over time are not themselves used in subsequent analysis.

To match state-level mortality data, TC incidence is collapsed from pixels to states for every month. If multiple storms impact a cell within a month, the maximum incidence at the cell level is recorded, and monthly averages are computed across pixels in each state. This spatial averaging causes our measures of incidence to be substantially lower than the maximum sustained wind speed commonly reported for storms, since only a small number of pixels experience those extreme conditions within each storm. We note that TC events generate the highest average state wind speeds compared to wind speeds from other intense storm phenomena, such as tornadoes. For reference, the minimum monthly state average wind speed we compute from TCs in our sample is 3.34 × 10−4 ms−1 and the 1st percentile is 8.4 × 10−3 ms−1. By contrast, the maximum monthly state average wind speeds from tornadoes in CONUS between 1950 and 2022 is 9.6 × 10−4 ms−1 and the 99th percentile is 3.8 × 10−5 ms−1, and for non-TC and non-tornado wind/hail events the maximum is 1.1 × 10−3 ms−1 and the 99th percentile is 1.3 × 10−4 ms−1. Therefore, the maximum non-TC wind events are comparable to the minimum (non-zero) TC events; and absent a TC, states do not experience average wind speeds of a similar magnitude as those from TCs.

Prior analysis by Hsiang & Jina3 demonstrated that spatial aggregates of TC exposure can be used as independent variables in a regression framework to obtain unbiased average effects that are expressed at finer spatial resolutions (footnote 13 on pages 16–17 of ref. 3). As long as there is no systematic correlation between the average intensity of a storm and the likelihood that the most intense regions within that storm strike the most populated (or economically active or vulnerable) pixels within a state, regression coefficients will not be biased by spatial aggregations. This condition would be violated if, for example, there were systematic patterns such that the eyes of a Category 3 hurricanes tended to pass directly over dense cities, but the eyes of Category 2 hurricanes tended to miss cities. However, given that the paths of storms are primarily controlled by random steering winds at high altitude, interacting with the beta-effect induced by the Earth’s meridional vorticity gradient, we have strong reason to believe that the spatial distribution of TC incidence within each state is orthogonal to the spatial distribution of underlying populations; and further that this covariance is independent of average TC intensity. Thus far, we know of no evidence that the trajectory of stronger (or weaker) storms systematically strike more vulnerable locations on land.

Of the 1,230 TCs that we reconstruct, 501 come within 250 km of a CONUS coastline. Intersecting these storms with state boundaries generates a total of 3,317 state-by-TC events. These longitudinal data reveals rich variation in the timing and intensity of TC incidence for individual states52 (Extended Data Fig. 1). Within-state variation in incidence season-to-season and month-to-month provides substantial variation in TC impulses that enable us to identify the impulse-response of mortality empirically.

All-cause mortality data

We analyse all-cause mortality at the state-year-month level between 1930 and 2015 using data from multiple sources. Data from 1900 to 2004 were digitized and assembled by Barreca et al.31 in their report identifying the impact of temperature on mortality. According to Barreca, this is the most comprehensive data on mortality assembled in this context. The remaining data was assembled by the authors using the CDC Underlying Cause of Death database. Data prior to 1959 was digitized from the Mortality Statistics of the United States annual volumes and is not otherwise available in a machine-readable format. Therefore, data in years prior to 1959 do not include cause of death (for example, cardiovascular disease) or demographic information (for example, age 1–44, Black)31. From 1959 to 2004, the data are from the machine-readable MCOD files, which include cause of death and demographic data.

For the years 2005–2015, we analyse mortality data from the public CDC Underlying Cause of Death database (2017). The data are based on death certificates for U.S. residents, which gives a single underlying cause of death and demographic data. Cause of death prior to 2000 was indexed using the four-digit ICD-9 code and 2000 onwards the index changed to the four-digit ICD-10 code. Cause of death was indexed using a four-digit ICD-10 code. We harmonized the cause of death into five categories that matched the cause of death variables from Barreca et al. We also construct a 6th category which is the difference between all-cause mortality and the sum of the 5 cause-specific categories, called ‘other’. Notably, the change in CDC ICD code methodology resulted in a shift in the counts of deaths from specific causes, particularly infectious diseases and cardiovascular disease.

To account for differences in underlying age-specific mortality we decomposed the effect of TCs on all-cause mortality by four age groups in the data: <1, 1–44, 45–64, and 65+ years of age. We were limited to these age groups because these are the designations in our historical data. We computed mortality with respect to the underlying population by these same four age groups. We compute mortality by race with respect to the population by race in each state and year. Black and white are the only race categories available for the entire sample. Extended Data Fig. 3 shows the monthly all-cause mortality rate, and our predicted monthly mortality rate, for all the states in our sample.

Direct deaths from TCs

Direct deaths from TCs are deaths that are officially attributed to a storm by the US government. We combine official death counts from two NOAA data sources. For storms between 1950 and 1996 we use the NOAA National Hurricane Center and Central Pacific Hurricane Center’s Hurricane in History6. Storms from 1997 to 2015 are from the NOAA Storm Events Database7.

Population data

We normalize state mortality by the population (per 100,000 people) in the state each month. Similar to the all-cause mortality data, these data must be combined from multiple sources. Pre-1968 population estimates are from Haines57; estimates for 1969–2000 are from the National Cancer Institute (2008); estimates for 2000–2010 are from the US Census Bureau, Intercensal Population and Housing Unit Estimates: 2000 to 2010 (ref. 58); estimates for 2010–2017 are from the US Census Bureau, US Population Estimates59.

Temperature data

Average monthly temperature data are from Berkeley Earth Surface Temperatures (BEST) land surface air temperature. BEST provides a monthly mean of average, minimum and maximum surface air temperature over land covering 1753 to the present56,60. The temperature data are based on a large inventory of observations from over 30,000 weather stations. Using these observations gridded temperature fields are reconstructed statistically, incorporating the reliability of individual weather stations and spatial variability of temperature56,60. Gridded BEST temperature data are then spatially aggregated, weighted by population, to the state-month-level.

Analysis

The econometric approach that we apply here is a top-down strategy, commonly called a ‘reduced-form’ analysis, that describes the overall net change of an aggregate outcome y (mortality) in response to exogenous treatments z (TC incidence). Under suitable conditions, this approach can identify causal effects on the outcome y induced by exogenous changes in independent variable z without explicitly describing all underlying mechanisms that link z to y, without observing intermediary variables x (for example, retirement savings accounts or healthcare infrastructure) that might link z to y, or without explicitly tracking other determinants of y unrelated to z (such as demographic trends or health policy), denoted w (refs. 23,61,62). Let f(•) describe a complex and unobserved process that generates state-level mortality rates \({y}_{{t}_{2}}\), occurring at time t2 based on x, w and z that occur both at times t1 and t2 (t1 < t2):

$${y}_{{t}_{2}}=f({z}_{{t}_{2}},{x}_{{t}_{1}}^{1}({z}_{{t}_{1}}),\ldots ,{x}_{{t}_{1}}^{K}({z}_{{t}_{1}}),\,{x}_{{t}_{2}}^{1}({z}_{{t}_{1}},{z}_{{t}_{2}}),\ldots ,{x}_{{t}_{2}}^{K}({z}_{{t}_{1}},{z}_{{t}_{2}}),\,{w}_{{t}_{1}\,}^{1}\ldots \,{w}_{{t}_{1}\,}^{J},\,{w}_{{t}_{2}\,}^{1}\ldots \,{w}_{{t}_{2}}^{J})$$

(1)

where \({x}_{{t}_{1}}^{k}\left({z}_{{t}_{1}}\right)\) indicates that the kth factor xk, which influences mortality rates y, at time t1 is itself affected by TCs at time t1. At time t2, xk may be influenced both by TCs at t2 and those that occur in the past at t1. Here, we let there be K pathways through which y is impacted by intermediary variables (x) and J ways through which determinants unrelated to TCs (w) impact y.

In this framework, the direct mortality impact of TC incidence usually reported by government agencies are the partial derivative:

$$direct{\rm{\_}}{deaths}_{{t}_{2}}=\frac{{\rm{\partial }}{y}_{{t}_{2}}}{{\rm{\partial }}{z}_{{t}_{2}}}$$

(2)

which are the deaths that occur contemporaneously and directly as a result of the geophysical event itself, holding fixed all other factors.

In this analysis, we directly estimate the total change in mortality that results from TC incidence in both current and prior moments in time, allowing for the possibility that changes in other factors that are influenced by TCs may indirectly affect mortality. In this case, the overall total change in mortality from TC incidence in t2 is the total derivative:

$$\frac{{\rm{d}}}{{\rm{d}}{z}_{{t}_{2}}}({mortality{\rm{\_}}rate}_{{t}_{2}})=\frac{{\rm{d}}{y}_{{t}_{2}}}{{\rm{d}}{z}_{{t}_{2}}}=\mathop{\underbrace{\frac{{\rm{\partial }}{y}_{{t}_{2}}}{{\rm{\partial }}{z}_{{t}_{2}}}}}\limits_{{\rm{d}}{\rm{i}}{\rm{r}}{\rm{e}}{\rm{c}}{\rm{t}}\,{\rm{d}}{\rm{e}}{\rm{a}}{\rm{t}}{\rm{h}}{\rm{s}}}+\mathop{\underbrace{\mathop{\sum }\limits_{k=1}^{K}\frac{{\rm{\partial }}{y}_{{t}_{2}}}{{\rm{\partial }}{x}_{{t}_{2}}^{k}}\frac{{\rm{\partial }}{x}_{{t}_{2}}^{k}}{{\rm{\partial }}{z}_{{t}_{2}}}}}\limits_{{\rm{i}}{\rm{n}}{\rm{d}}{\rm{i}}{\rm{r}}{\rm{e}}{\rm{c}}{\rm{t}}\,{\rm{d}}{\rm{e}}{\rm{a}}{\rm{t}}{\rm{h}}{\rm{s}}}$$

(3)

which includes both direct deaths and deaths that result from any of the K possible pathways that depend on the intermediate variables xk at time t2. Empirically, we find that direct deaths are much smaller than indirect deaths in CONUS.

In addition, we also account for the possibility that deaths are delayed relative to TC incidence. Because direct deaths are usually tabulated immediately following storms, there are negligible direct deaths that are delayed. However, once we begin considering indirect deaths, it becomes possible for substantial delays to emerge due to the dynamics of different pathways xk. In our analysis, we also estimate the total deaths that occur at time t2 as a result of TC incidence that occurs at an earlier time t1, which is the total derivative

$$\frac{{\rm{d}}}{{\rm{d}}{z}_{{t}_{1}}}({mortality{\rm{\_}}rate}_{{t}_{2}})=\frac{{\rm{d}}{y}_{{t}_{2}}}{{\rm{d}}{z}_{{t}_{1}}}=\mathop{\underbrace{\mathop{\sum }\limits_{k=1}^{K}\frac{{\rm{\partial }}{y}_{{t}_{2}}}{{\rm{\partial }}{x}_{{t}_{2}}^{k}}\frac{{\rm{\partial }}{x}_{{t}_{2}}^{k}}{{\rm{\partial }}{z}_{{t}_{1}}}}}\limits_{{\rm{i}}{\rm{n}}{\rm{d}}{\rm{i}}{\rm{r}}{\rm{e}}{\rm{c}}{\rm{t}}\,{\rm{v}}{\rm{i}}{\rm{a}}\,{x}_{{t}_{2}}}+\mathop{\underbrace{\mathop{\sum }\limits_{k=1}^{K}\frac{{\rm{\partial }}{y}_{{t}_{2}}}{{\rm{\partial }}{x}_{{t}_{1}}^{k}}\frac{{\rm{\partial }}{x}_{{t}_{1}}^{k}}{{\rm{\partial }}{z}_{{t}_{1}}}}}\limits_{{\rm{i}}{\rm{n}}{\rm{d}}{\rm{i}}{\rm{r}}{\rm{e}}{\rm{c}}{\rm{t}}\,{\rm{v}}{\rm{i}}{\rm{a}}\,{x}_{{t}_{1}}}$$

(4)

This expression does not contain a term for direct deaths, but it contains two summations which capture the effects of past TC incidence \(({{z}_{t}}_{1})\) on current mortality \((\,{y}_{{t}_{2}})\) via past intermediate variables \(({x}_{{t}_{1}})\) and current intermediate variables \(({x}_{{t}_{2}})\). In practice, we explore the possibility of indirect effects that emerge over the course of 240 months following TC incidence, one could generalize this framing to a corresponding number of summations.

The possibility of delayed indirect deaths has two major implications regarding how indirect mortality is estimated and how those results are interpreted. First, because TC incidence at multiple points in the past, as well as the present, might affect current mortality, we must account for both the present and past influence of TC incidence simultaneously for each instance of the outcome. This is accomplished via deconvolution25,26,27,56,63,64, implemented here using a distributed lag-model solved via ordinary least squares, detailed below.

Second, each TC event affects mortality outcomes at multiple points in time, thus computing the full impact of a TC event requires summing these impacts that might emerge gradually. In the simplified two-period framework above, the total impact from TC incidence at t1 is then

$$\frac{{\rm{d}}}{{\rm{d}}{z}_{{t}_{1}}}(mortality{\rm{\_}}{rate}_{{t}_{1}}\,+\,mortality{\rm{\_}}{rate}_{{t}_{2}})\,=\,\frac{{\rm{d}}{y}_{{t}_{1}}}{{\rm{d}}{z}_{{t}_{1}}}\,+\,\frac{{\rm{d}}{y}_{{t}_{2}}}{{\rm{d}}{z}_{{t}_{1}}}$$

(5)

which can be expanded further by substituting from the equations above. These terms, if plotted separately, characterize the impulse response of y in reaction to the TC impulse \({z}_{{t}_{1}}\). In our actual analysis, we compute the average cumulative impact of a single TC that occurs at t0 over 240 subsequent months (t1 − t240). Following substitution and simplification, this can be expressed as

$$\frac{{\rm{d}}}{{\rm{d}}{z}_{{t}_{0}}}\left(\mathop{\sum }\limits_{\ell =0}^{240}{mortality{\rm{\_}}rate}_{{t}_{\ell }}\right)=\mathop{\underbrace{\frac{{\rm{\partial }}{y}_{{t}_{0}}}{{\rm{\partial }}{z}_{{t}_{0}}}}}\limits_{{\rm{d}}{\rm{i}}{\rm{r}}{\rm{e}}{\rm{c}}{\rm{t}}\,{\rm{d}}{\rm{e}}{\rm{a}}{\rm{t}}{\rm{h}}{\rm{s}}}+\,\mathop{\underbrace{\mathop{\sum }\limits_{\ell =0}^{240}\mathop{\sum }\limits_{k=1}^{K}\left(\mathop{\sum }\limits_{j=\ell }^{240}\frac{{\rm{\partial }}{y}_{{t}_{j}}}{{\rm{\partial }}{x}_{{t}_{\ell }}^{k}}\right)\frac{{\rm{\partial }}{x}_{{t}_{\ell }}^{k}}{{\rm{\partial }}{z}_{{t}_{0}}}}}\limits_{{\rm{i}}{\rm{n}}{\rm{d}}{\rm{i}}{\rm{r}}{\rm{e}}{\rm{c}}{\rm{t}}\,{\rm{d}}{\rm{e}}{\rm{a}}{\rm{t}}{\rm{h}}{\rm{s}}}$$

(6)

which describes the overall total impact of a storm through all pathways across all possible delays ℓ = [0, 240]. Note that neither K nor xk need ever be specified explicitly in our estimation below. This expansion reveals that, accounting for numerous possible pathways operating over different delays, a single TC event can potentially generate a total mortality impact much larger than the direct deaths traditionally reported.

To compute the overall mortality burden imposed by the TC climate of CONUS, we compute the full mortality response across all age groups in each state, accounting for the incidence of each storm on the state:

$$\begin{array}{l}mortality{\rm{\_}}burden\,=\,\mathop{\sum }\limits_{{\ell }=0}^{240}\sum _{t\in {\rm{m}}{\rm{o}}{\rm{n}}{\rm{t}}{\rm{h}}{\rm{s}}}\sum _{s\in {\rm{s}}{\rm{t}}{\rm{o}}{\rm{r}}{\rm{m}}{\rm{s}}}\sum _{i\in {\rm{s}}{\rm{t}}{\rm{a}}{\rm{t}}{\rm{e}}{\rm{s}}}{population}_{i,t+{\ell }}\\ \,\,\,\,\,\cdot {z}_{sit}\cdot \frac{{\rm{d}}}{{\rm{d}}z}{(mortality{\rm{\_}}rate)}_{i,t+{\ell }}\end{array}$$

(7)

where zsit is the TC incidence of storm s on state i in month t, populationi,t+ℓ is the population in state i in month t + ℓ, and \(\frac{{\rm{d}}}{{\rm{d}}z}{(mortality{\rm{\_}}rate)}_{i,t+{\ell }}\) is our estimate for the total impact of TC incidence on the mortality rate in state i in month t + ℓ. In practice, this effect is nonlinear, but it is expressed linearly here for simplicity. Information on state i affects the impulse responses used in these calculations because TC risk is computed by state and affects the structure of the impulse-response function.

Econometric implementation

Identification

Our econometric analysis exploits the quasi-random variation in the location and intensity of TC incidence to estimate the impact of TCs on mortality separately from other known and unknown factors that affect mortality across locations and over time. As described above, this reduced-form approach captures the effect of all possible channels of influence that may increase mortality after a TC33,65. Because the location, timing and intensity of TC incidence is determined by oceanic and atmospheric conditions that are beyond the control of individual states, we assume mortality TC incidence is as good as randomly assigned23,61. For reference, Extended Data Fig. 1 shows the sequence of monthly TC incidence by state for all the states in our sample.

We note that some early analyses of natural disaster impacts utilized social outcomes (for example, direct economic damage66 or direct mortality) as a proxy measure of physical hazard severity. However, it is now understood that use of these metrics as independent variables may confound estimated treatment effects, since they are endogenously determined by many of the same underlying covariates (for example, healthcare, infrastructure, inequality and institutions) that mediate other outcomes from disasters22,29,33,34,35,67. Thus, use of these proxy measures for hazard severity exposes analyses to selection biases, since population characteristics may cause observational units to ‘select’ into more or less severe treatment23. We therefore focus this analysis strictly on independent variables that are physical measures of TC incidence (wind speed), because they are exogenous and cannot be influenced by the populations that are impacted22.

Deconvolution

In considering the long-run impact of TCs on mortality, we hypothesize that there may be a delay between the geophysical event and components of the mortality response. Because TCs are regular events that occur frequently in CONUS, the possibility of this delay means that the time series of mortality outcomes we observe in data may be the result of overlapping responses from multiple storms. Extended Data Fig. 2 displays a cartoon of this data-generating process. In such a context, the empirical challenge is isolating the impact from individual storms which might be partially confounded by the overlapping TC signals from earlier or later storms. We use the well-established signal-processing approach of deconvolution25,26,27,33,63,64 to recover the characteristic impulse-response function for a TC impulse. Conceptually, this approach searches for an impulse-response function that, if applied to all TCs in the data simultaneously, best fits the observed outcome data. Stated another way, this approach estimates the effect of each TC accounting for the potential overlapping impact of all other TCs, subject to the constraint that TCs share a characteristic impulse-response function.

This method assumes that the overlapping responses influencing mortality at a moment in time are additively separable, an assumption that we think is reasonable given the overall small impact that any individual storm event has mortality rates at a moment in time in a particular region (0.019% on average, 0.04% for storms at the 95th percentile). We solve for the structure of the impulse response, characterized by a set of coefficients β, using ordinary least squares (OLS). This is a standard procedure that is commonly applied in a wide range of disciplines26. In some fields, such as econometrics, deconvolution is frequently described as estimation of distributed lags27.

Baseline specification

Our main results are based on a linear model of TC incidence on mortality rate. Indexing states by i and month of sample by t, we solve the model

$$\begin{array}{l}mortality{\rm{\_}}{rate}_{it}=\mathop{\sum }\limits_{{\ell }=-72}^{240}({\beta }_{{\ell }}\cdot wind{\rm{\_}}{speed}_{i,t-{\ell }})\\ \,\,\,\,\,\,\,\,+{\delta }_{1,i}\cdot {temp}_{it}\cdot {s}_{i}+\,{\delta }_{2,{\rm{i}}}\cdot {temp}_{it}^{2}\cdot {s}_{i}+{\mu }_{1}\cdot {m}_{it}\\ \,\,\,\,\,\,\,\,+\,\mathop{\sum }\limits_{n=1}^{8}({\eta }_{n}\cdot {s}_{i}\cdot {t}^{n})\,+\,{\mu }_{2}\cdot {m}_{it}\cdot t+{\mu }_{3}\cdot {h}_{t}+{{\epsilon }}_{it}\end{array}$$

(8)

via OLS. Here wind_speedi,t−ℓ is TC maximum wind speed ℓ months prior to month t, si is a vector of state-specific dummies, mit are state-by-month-of-the-year dummies (for example, an indicator variable for whether state = Florida and month = January), ht are month-of-sample dummies (for example, an indicator variable for whether the month = January, 1974), \({temp}_{it}\cdot {s}_{i}\) is month-of-sample temperature interacted with state dummies, and \({temp}_{it}^{2}\cdot {s}_{i}\) is squared month-of-sample temperature (also interacted with state). Each coefficient βℓ measures the marginal effect of an additional ms−1 of wind speed incidence on mortality ℓ months after a TC conditional on the effect of any prior TC. We include 72 lead terms in equation (8) as a falsification test, also known as negative exposure controls, since idiosyncratic future TC incidence should not alter current health outcomes.

This model accounts for state-specific quadratic effects of temperature on mortality based on prior literature, which has shown that very hot and very cold temperatures cause higher levels of mortality relative to more moderate temperatures31,33,34,35. Each state is allowed to express a different mortality response to temperature extremes, implemented via interaction with the state dummy variable si. Supplementary Fig. 6 shows the state-specific shape of the quadratic functions we estimate for the temperature-mortality response. Consistent with prior findings studying patterns of adaptation31,33,34,35, we observe that some states have a flatter response at temperatures that are more common for that state (for example, cold in Minnesota) while other states have steeper curves at those same temperatures if they are less common (for example, cold in Florida). In an effort to balance parsimony with model richness, we omit extended lags of temperature based on prior literature demonstrating that impacts on mortality dissipate within a month33,68.

This model also non-parametrically accounts for:

  • State-by-month-specific constants (fixed effects) that capture average differences between states, as well as unique seasonable patterns within states (\({\mu }_{1}\cdot {m}_{it}\)). These terms will account for differences in mortality driven by unobserved factors at the state level, such as health policies, as well as factors that cause seasons within a state to exhibit higher mortality, such as holidays.

  • State-specific nonlinear trends in mortality, captured by eighth-order polynomials in month-of-sample interacted with state fixed effects (\({\sum }_{n=1}^{8}({\eta }_{n}\cdot {s}_{i}\cdot {t}^{n})\)). These trends account for unobserved factors that have caused mortality within states to change over time, such as changing health policies or demographic trends.

  • Trends in state-specific seasonal patterns of mortality, captured by a linear trend in month-of-sample interacted with state and month fixed effects (\({\mu }_{2}\cdot {m}_{t}\cdot t\)). These trends are additive to the state-specific polynomial and allow for the model to express gradual convergence or divergence in the seasonality of mortality within a year, and allows for these changes to differ by state. These trends account for unobserved factors that drive gradual changes over time that may cause mortality in certain times of year (for example, January) to change relative to other times of year (for example, June). For example, if adoption of safety standards has reduced wintertime mortality from motor vehicle accidents or improvements in medical care have reduces summertime deaths from infectious diseases. Extended Data Fig. 4a illustrates the combined effect of these state-specific seasonal trends, state-specific polynomials, and state-by-month- specific constants on model predictions for Florida and New Jersey. For example, the seasonality of mortality in Florida has lessened over time, in conjunction with other nonlinear trends.

  • National month-of-sample fixed effects that capture nonlinear and/or discontinuous changes in mortality rates nation-wide (\({\mu }_{3}\cdot {h}_{t}\)). These terms are particularly important for capturing idiosyncratic spikes in mortality that result from nation-wide conditions, such as influenza outbreaks, as well as any systematic changes in the accounting methodology of mortality by the CDC. Comparisons of Extended Data Fig. 4b,c illustrates how the inclusion of these terms in the model alters the ability of the model to capture unusual spikes in mortality that are not captured by other model elements, including the trends listed above, TCs, and temperature.

Overall, the fit for this model is high (in-sample adjusted R2 = 0.93 with 25,062 degrees of freedom). Extended Data Fig. 3 overlays predictions with observations for all states (same as in Fig. 1c). We find that all of the non-parametric controls listed above are important for passing standard specification checks. For example, failure to account for trends flexibly enough causes estimated leads to deviate from zero or randomization-based placebo tests (described below) to recover non-zero central estimates. These results are unchanged if we use a Poisson regression specification (Extended Data Fig. 6a).

In a robustness test, we interact the month-of-sample fixed effects with 3 region indicator variations. We continue to obtain our main findings after introducing these additional 2,062 parameters to the model, although the estimates become much noisier and attenuate slightly. Both of these effects are well understood results of including a large number of highly flexible variables that absorb a meaningful fraction of the true variation in the independent variable69.

We evaluate the distribution of the unmodelled variation represented by the error term ϵit and find that it essentially follows a Normal distribution except with slightly positive kurtosis (Supplementary Fig. 5a). The distribution of these residuals appears stationary throughout the sample period and independent over time (Supplementary Fig. 5b). The consistency of the distribution of these errors is attributed to the high degree of flexibility in the non-parametric terms of our econometric specification, which are able to capture those components of the data-generating process that would otherwise appear as auto-correlated errors. On the basis of this evaluation, we construct OLS standard-error estimates as the underlying assumptions for these estimates appear to be reasonably satisfied. In addition, we find strong support for this modelling choice when we conduct a variety of permutation tests for statistical significance (Extended Data Fig. 5), all of which indicate that our asymptotic estimates for confidence intervals are correctly (and possibly conservatively) sized and our tests for statistical significance correctly powered. Notably, these permutation tests do not rely on the assumptions used to estimate these confidence intervals, thus they can be considered independent corroboration for the validity of this approach.

Cumulative effects

To compute the total effect after the TC makes landfall, we estimate the cumulative sum of βℓ for each ℓ ∈ [−72, 240]. We compute \({\Omega }_{{\ell }}={\sum }_{k=o}^{{\ell }}{\beta }_{k}\), which denotes the cumulative impact of an additional 1 ms−1 wind speed incidence on mortality ℓ months after a TC event. We account for the estimated covariance of β when estimating uncertainty in Ω. We normalize the sums relative to the impact one month prior to the TC, ℓ = −1 such that Ω−1 = 0.

Randomization-based placebo tests

Given the complexity of our model, the long delays we study, and the absence of prior analyses of long-run total mortality from TCs, it is not possible to subjectively evaluate our econometric analysis against any prior benchmark. In such a context, there is risk of unknowingly recovering a spurious estimate generated as an artefact of our model specification. A strong test designed to avoid such artefacts is to ensure that model estimates of TC impacts on mortality are unbiased in a variety of situations where the structure of the association has been manipulated. In four tests, we shuffle the true TC data in different ways. In each case, this shuffling should break any correlation between TC incidence and mortality such that an unbiased estimate of the effect of shuffled TCs on mortality is zero. However, in each case, some of the structure in the original TC data are allowed to remain in the shuffled TC data. For example, randomization within a state over time retains the average cross-sectional patterns of TC incidence, but destroys any time- series structure. Thus, these tests allow us to examine whether, in each case, the remaining structure generates artefacts in the model that would produce a spurious result, also known as a negative exposure control70. Any non-zero correlation, on average, would indicate a biased model where the bias is driven by the non-randomized components of the original TC data.

Within each type of randomization we scramble TC assignment 1,000 times and run the linear version of the model (equation (8)) on each re-sampled version of the data. Our four randomizations are illustrated graphically in Supplementary Fig. 7 and described below:

  • Total randomization shuffles TC events across all state-by-month observations. This tests whether the unconditional marginal distribution of TC events, which has a long right tail, could generate bias. Results are shown by light blue boxes in Fig. 1g.

  • Within-state randomization shuffles the sequencing of TCs that a state experiences over time. TCs are always assigned to the correct state, but the month and year assigned to each storm is random. The cross-sectional average pattern of storm incidence is preserved in the data. Thus, this tests whether time-invariant cross-sectional patterns across states generate spurious correlations. Results are shown by dark blue boxes in Fig. 1g.

  • Within-month randomization shuffles the TC incidence across states within each month-of-sample. TCs are always assigned to the correct month and year, but the state assigned to each storm is random. The average time-series structure of TC incidence nation-wide is preserved. Thus, this tests whether national or seasonal trends, which are nonlinear, could bias estimates produced by this model. Results are shown by maroon boxes in Fig. 1g.

  • Across-state shuffles complete TC times-series across states, keeping the timing and sequence of storms correct as blocks. TCs are always assigned to the correct month and year, and the sequence of storms experienced by a state is always a continuous sequence that is observed in the data. However, the state that is assigned that sequence is randomly chosen. This tests whether trends within a state and within the sequence of storms that a state experiences could generate bias. This test differs from the within month randomization because state-level trends often differ across states (see Extended Data Fig. 1 and Extended Data Fig. 3) and there are complex seasonal patterns that could potentially affect estimates. Results are shown by red boxes in Fig. 1g.

The estimated impact of TCs in each of these placebo tests is zero on average, to within a high degree of precision. Extended Data Fig. 5 illustrates distributions of estimates for all lags. These results demonstrate that non-exchangeability across states within a month, across months within a state, or across states (conditional on month of sample) does not confound our analysis; indicating that the rich set of fixed effects and trends successfully adjust for many patterns of TC incidence and/or mortality such that the remaining conditional variation is as good as random.

Permutation tests for statistical significance

In addition to establishing the unbiasedness of our main point estimates, the four randomizations above can be utilized to serve a second function: estimating statistical significance of our estimates. These randomizations enable approximate permutation tests71, allowing for different types of autocorrelation to remain in the TC data. We use these randomizations to examine the likelihood of randomly obtaining an entire impulse-response function similar our actual estimate, if in reality no such relationship exists in the data. To do this, we jointly test the significance of all true cumulative estimates Ωℓ against the null hypothesis that a similarly extreme sequence of estimates is generated randomly. Extended Data Fig. 5 overlays the true estimates for Ωℓ on distributions of similar estimates from each randomization. P values for individual lag terms \(\left({p}_{{\ell }}=\Pr \left(\left|{\Omega }_{{\ell }}^{{\rm{randomized}}}\right| > \left|{\Omega }_{{\ell }}\right|\right)\right)\) are plotted in the right subpanels and are all individually statistically significant (P < 0.05) for ℓ < 150 months in each randomization. However, the significance of the complete sequence of coefficients that together compose the entire impulse-response function is far greater. We compute a joint P value for the full impulse response between 0 and 172 months \(({p}_{{\ell }}=\Pr ({\cap }_{{\ell }=0}^{172}(| {\Omega }_{{\ell }}^{{\rm{randomized}}}| > | {\Omega }_{{\ell }}| ))\) that ranges from P = 0.012 to P = 0.0012 across the randomization approaches (Extended Data Fig. 5). We conclude that it would be extremely unlikely to obtain an impulse-response function as extreme as our main result due to chance.

Subsamples by age, race and cause of death

In addition to the all-cause mortality rate for the entire population, we also present the results stratified by age, race, and cause of death. For the six cause-specific mortality rates we compute mortality per 100,000 of the total population in that state in the time period (for example, total number of deaths from cardiovascular disease divided by total population times 100,000). ‘Other’ mortality is the difference between the total deaths and the sum of all the other causes. For age groups and race we report mortality risk as the outcome of interest. For example, for the Black population, we construct the mortality rate for Black people as the number of deaths of Black people in the state divided by the Black population. We do the same procedure by age group. We also report the mortality by these strata as a proportion of the total deaths traceable to TC incidence, see Extended Data Fig. 7.

Subsamples by average TC risk

To evaluate whether there is heterogeneity in the mortality response of states that are frequently exposed to TCs compared to those infrequently exposed, we stratify the sample by the average TC incidence they experience. We allow the mortality impulse response to differ based on quartiles of states, sorted by their average TC incidence. We implement this by including and interaction with an indicator variable for the quartile of their average wind speed incidence, following the general approach for modelling adaptation developed in refs. 29,33. Average wind incidence is a measure of the expected TC risk a population bears, which informs preventive risk reduction investments, behaviours, or other adaptive actions they take to reduce the expected harm from TCs. We approximate this measure by computing as the mean wind speed in each state i across the period t in our sample

$$wind{\rm{\_}}{speed}_{i}=\frac{{\sum }_{t=1}^{1,032}wind{\rm{\_}}{speed}_{it}}{\mathrm{1,032}}$$

and assigning the quartiles of these means to each state. We estimate a model that allows each quartile to express a different impulse response to TCs and observe little difference in the impact of TCs on mortality between the second through fourth quartile (Extended Data Fig. 6c). The effect for the second through fourth quartile are not statistically significantly different than the effect for the second through fourth quartile, combined (P = 0.38). Thus, to improve the efficiency of our model and limit unnecessary noise in our estimates, we pool quartiles 2–4 to create the ‘high average incidence’ group in the main results (shown in Fig. 2f) (‘high incidence’ in equation (9)). The ‘low average incidence’ group is the first quartile of average wind speed alone (‘low incidence’ in equation (9)). We additionally evaluate whether the spatial distribution of populations relative to the coast, within states, alters the mortality response to TCs. Stratifying states on the basis of the average fraction of the population that lives in coastal counties, we fail to find evidence that states with high concentrations of coastal populations are systematically different from states with little or no coastal population (Extended Data Fig. 6d). Lastly, we evaluate whether the overall spatial correlation between populations and average wind speed incidence, within each state, alters the mortality response to TCs. Stratifying states on the basis of within-state spatial correlation across 0.125° × 0.125° pixels, we fail to find evidence that states with higher spatial correlations are systematically different from states with little or negative correlations (Extended Data Fig. 6e).

Nonlinear effects of TCs

We evaluate whether the mortality impact of a TC is nonlinear in the physical intensity of the event. This could occur, for example, if more extreme TC events generate exponentially more physical damage51,72 or if they elicit different government responses18,73. Empirically, we find that excess mortality 180 months after a TC is well approximated by a linear function of max wind incidence, particularly for TCs with area-average max wind speeds between 0 and 20 ms−1, which is the majority of events (93%) in our sample (Extended Data Fig. 8). However, for the most extreme events (>30 ms−1, 1.4% of events) we find that excess mortality is generally lower than a linear function would predict, although these nonlinear effects are not themselves statistically significant. We lack the data to fully evaluate the underlying causes of this nonlinearity, but believe it is an important topic for future study. For example, it is possible that societal responses to the most extreme events (for example, disaster relief) are more effective at alleviating mortality impacts of TCs because these events attract a disproportionate quantity of attention, compared to less extreme events that are also harmful but less salient74. Regardless of their cause, we account for these non-linearities in calculations below as they contain information on how populations in CONUS have adapted to their TC climates5,73.

To estimate nonlinear effects of TCs, we estimate a model that is identical to the benchmark linear model in equation (1), but it allows the magnitude of the TC mortality impulse-response function to be cubic in TC incidence. The motivation for this approach is the possibility that the relationships between wind speed and long-run mortality does not increase linearly. For example, very high wind speeds may cause extreme damages and/or elicit greater governmental and humanitarian responses, which would mean that a unit of increase from 40 to 41 ms−1 wind speed may have a larger or small mortality impact compared to an increase from 5 to 6 ms−1. Since the nonlinear impact of TCs may be influenced by the historical TC experience and baseline TC risk of populations, this model also allows for the nonlinear response to differ based on the risk categorization for each state:

$$\begin{array}{l}mortality{\rm{\_}}{rate}_{it}=\\ \,\mathop{\sum }\limits_{{\ell }=-72}^{240}\mathop{\sum }\limits_{r=1}^{3}({\theta }_{r,{\ell }}^{low{\rm{\_}}incidence}\cdot wind{\rm{\_}}{speed}_{i,t-{\ell }}^{r}\cdot {Q}_{i}^{low{\rm{\_}}incidence}\\ \,+\,{\theta }_{r,{\ell }}^{high{\rm{\_}}incidence}\cdot wind{\rm{\_}}{speed}_{i,t-{\ell }}^{r}\cdot {Q}_{i}^{high{\rm{\_}}incidence})\\ \,+\,controls+{{\epsilon }}_{it}\end{array}$$

(9)

where r indicates an exponent and the controls are identical to those in equation (8). \({Q}_{i}^{low{\rm{\_}}incidence}\) (\({Q}_{i}^{high{\rm{\_}}incidence}\)) is an indicator variable that is set to one if state i is in the low incidence (high incidence) group. The coefficients θ1ℓ, θ2ℓ and θ3ℓ separately capture the cubic relationship between wind speed incidence and mortality for low and high-risk states in each lag period. Extended Data Fig. 8 displays the cumulative impact estimated using both the linear and nonlinear models after 180 months. These results are unchanged if we alternatively use a cubic spline regression specification (Extended Data Fig. 8). The low-risk response is slightly convex relative to the linear estimate, while the high-risk response is slightly concave. In both cases, impacts are relatively well approximated by the linear version of the model and only diverge (insignificantly) at very high levels of incidence that are rare in sample. Distributions for in sample frequency are shown in lower panels of Extended Data Fig. 8.

Computing mortality burdens

The total impacts of all TCs on mortality are estimated using each version of the model, presented in Supplementary Table 1. We compute the excess mortality from TCs by state, month, and TC. These estimates are presented in Figs. 3 and 4 and Supplementary Tables 1 and 2.

Figure 4 displays the estimated excess mortality from TCs by state, age and race, computed using equation (8) applied to equation (7). Figure 4a presents our estimated full TC mortality burden, similar to equation (7) but by state (mortality_burdenit) as an average proportion of total deaths (mortalityit) in each state between 1950 and 2015:

$${proportion}_{i}=\frac{{\sum }_{t\in {\rm{m}}{\rm{o}}{\rm{n}}{\rm{t}}{\rm{h}}}mortality{\rm{\_}}{burden}_{it}}{{\sum }_{t\in {\rm{m}}{\rm{o}}{\rm{n}}{\rm{t}}{\rm{h}}}{mortality}_{it}}$$

Similarly, we estimate the proportion by state and age group (proportioni,a), shown in Fig. 4d,f. Proportion and total excess mortality for the Black population is based on mortality burden estimated with the mortality risk for the Black population, therefore \({proportion}_{Black}=\frac{{\sum }_{t\in {\rm{m}}{\rm{o}}{\rm{n}}{\rm{t}}{\rm{h}}}{\sum }_{i\in {\rm{s}}{\rm{t}}{\rm{a}}{\rm{t}}{\rm{e}}}mortality{\rm{\_}}{burden}_{it,{\rm{B}}{\rm{l}}{\rm{a}}{\rm{c}}{\rm{k}}}}{{\sum }_{t\in {\rm{m}}{\rm{o}}{\rm{n}}{\rm{t}}{\rm{h}}}{\sum }_{i\in {\rm{s}}{\rm{t}}{\rm{a}}{\rm{t}}{\rm{e}}}{mortality}_{it,{\rm{B}}{\rm{l}}{\rm{a}}{\rm{c}}{\rm{k}}}}\).

Figure 4g illustrates the impact of the TC climate on the geographic differences in average annual all-cause mortality rate between states that do not experience TCs (‘non-TC states’) and states that do (‘TC states (actual)’), in this context. We also subtract the average annual TC mortality burden from the actual average annual mortality for each TC-impacted state (‘TC states without TCs’).

Decomposing trends in mortality burden

We examine the differences in TC events and population distribution before 2001 and after 2001 in order to understand why the mortality burden after 2001 is increasing more rapidly than it did before 2001 (Fig. 4h–j). Figure 4h shows the distribution of the number of TCs that made landfall each year before 2001 and after 2001. Similarly, Fig. 4i shows the maximum annual wind speed per year and Fig. 4j plots the average wind speed per year as experienced by a proportion of the CONUS population. The changes in the distribution of TC events affecting CONUS after 2001 were themselves probably caused by a combination of factors, including warmer sea surface temperatures11,75 and reductions of anthropogenic aerosol emissions76,77 (which create an environment more amenable to TC intensification); and shifts in steering winds78 (which direct a larger fraction of TCs to landfall in CONUS after formation). We note that identifying factors driving the TC climate remains an active area of research56,79.

To understand the 1950 to 2015 trend in the national aggregate mortality burden, we re-estimate mortality_burdent for each month of the sample with different populations to decompose the long-term trend based on various population patterns. We first replace populationi,t+ℓ from equation (7) with fixed 1950 or 2015 populations (Fig. 4k, red and maroon lines). To generate the yellow line in Fig. 4k, we replace populationi,t+ℓ with an estimate of the 2015 population ‘deflated’ to 1950 levels. Specifically, we first compute a national population deflation fraction,

$$\Delta =\frac{{\sum }_{i\in {\rm{s}}{\rm{t}}{\rm{a}}{\rm{t}}{\rm{e}}}{population}_{i,2015}-{\sum }_{i\in {\rm{s}}{\rm{t}}{\rm{a}}{\rm{t}}{\rm{e}}}{population}_{i,1950}}{{\sum }_{i\in {\rm{s}}{\rm{t}}{\rm{a}}{\rm{t}}{\rm{e}}}{population}_{i,1950}}$$

where populationi,2015 is the state-specific population in 2015 and populationi,1950 is the state population in 1950. We then calculate

$$deflated\_\,{population}_{i,2015}=\frac{{\sum }_{i\in {\rm{s}}{\rm{t}}{\rm{a}}{\rm{t}}{\rm{e}}}{population}_{i,2015}}{1+\Delta }$$

and apply deflated_populationi,2015 to equation (7). This value is an adjusted state-level population that allows the total national population to match 1950 level but have a relative spatial distribution that reflects 2015.

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.

RELATED ARTICLES

Most Popular

Recent Comments