ABSTRACT
Objectives A nationwide study of the impact of high temperature on respiratory disease hospitalizations among older adults (65+) living in large urban centers.
Methods Daily rates of short-stay, inpatient respiratory hospitalizations were examined with respect to variations in ZIP-code-level daily mean temperature in the 120 largest US cities between 2000-2017. For each city, we estimated cumulative associations (lag-days 0-6) between warm-season temperatures (June-September) and cause-specific respiratory hospitalizations using time-stratified conditional quasi-Poisson regression with distributed lag non-linear models. We estimated nationwide associations using meta-regression and updated city-specific associations via best linear unbiased prediction. With stratified models, we explored effect modification by age, sex, and race (Black/white). Results were reported as percent change in hospitalizations at high temperatures (95th percentile) compared to median temperatures for each outcome, demographic-group, and metropolitan area. Excess hospitalization rates were estimated for days above median temperatures.
Results At high temperatures, we observed increases in the percent of all-cause respiratory hospitalizations [1.2 (0.4, 2.0)], primarily driven by an increase in respiratory tract infections [1.8 (0.6, 3.0)], and chronic respiratory diseases/respiratory failure [1.2 (0.0, 2.4)]. East North Central, New England, Mid-Atlantic, and Pacific cities accounted for 98.5% of the excess burden. By demographic group, we observed disproportionate burdens of heat-related respiratory hospitalizations among the oldest beneficiaries (85+ years), and among Black beneficiaries living in South Atlantic cities.
Conclusion This study found robust impacts of high temperature on respiratory failure and chronic inflammatory and fibrotic diseases among older adults. The geographic variation suggests that contextual factors account for disproportionate burdens.
INTRODUCTION
Exposure to extreme heat and humidity has a direct impact on human health and well-being (1–3), and is projected to result in a 370% increase in heat-related mortality among those 65 years and over by mid-century (2). Robust associations between heat and excess mortality among older adults are consistently reported in the literature. However, associations between heat and health care utilization, such as emergency department visits and hospitalizations, are often weaker with inconsistent findings across study locations (4–9). Previous studies of heat and respiratory morbidity have primarily focused on all-cause respiratory events or major categories of respiratory disease (chronic obstructive pulmonary disease, asthma, and respiratory tract infections), with little attention given to less common respiratory diseases. The impact of extreme heat on cause-specific morbidity is an important knowledge gap that must be addressed to optimize patient outcomes and increase the resilience of health care systems for future extreme heat events.
Older adults are among the most physiologically sensitive to prolonged heat exposure because of the effects of aging on thermoregulation and evaporative heat loss. In addition, older adults may have impaired mobility, pre-existing diseases, or take medications that impact body fluid balance or interact with hemodynamic regulation (3, 10). When prolonged heat exposure results in hyperthermia, older adults may experience myocardial strain due to increased cardiac output, interstitial fluid pooling, as well as systemic inflammation and reactive oxygen species production (3, 11–13). Hyperthermia can also lead to increases in ventilation rate, tidal volume, and respiratory rate (14), which may damage the lung parenchyma and exacerbate underlying respiratory illnesses. Epidemiologic and mechanistic studies provide evidence that some heat-related pulmonary injury may also be mediated by inhaling hot air (11, 15–18). Proposed pathways include activation of bronchopulmonary vagal afferent C fibers, increased cholinergic responses, and stimulation of specific heat shock proteins leading to both epithelial barrier dysfunction and airway inflammation (11, 16, 17).
In addition to advanced age, numerous studies have demonstrated that disproportionate exposure and vulnerability to hazardous heat is influenced by complex interactions between municipal intervention efforts, socioeconomic inequality (19), structural racism (20, 21), access to health care, adaptive capacity, risk perception, and physiological acclimatization and susceptibility (6, 19–23). Urban populations have also been shown to be at an increased risk of heat-related health outcomes due to the urban heat island effect (24) in which the heat retaining and generating properties of the urban form (low albedo surfaces, building configurations that trap heat, diverse sources of heat generation, decreased vegetation, and less evapotranspiration) result in higher air temperatures in urban environments compared to rural areas (25).
In this large, nationwide study of respiratory hospitalizations among Medicare beneficiaries (aged 65 years and older) living in large urban centers, we assessed the role of high ambient temperatures on cause-specific respiratory morbidity between 2000 and 2017. We further quantified the heat-related respiratory burden for specific disease outcomes and assessed whether the burden of all-cause respiratory hospitalizations varied by demographic group and/or geographical setting. A novel contribution of this work defines the relationship between heat and previously unexplored respiratory diagnoses, as well as evaluating the variation in heat-attributable burden and risk across 120 large urban centers.
METHODS
Outcome Ascertainment
We obtained short-stay, inpatient hospitalization data for all Medicare beneficiaries aged 65-114 years during 2000-2017 from the Centers for Medicaid and Medicare Services. We excluded hospitalizations for Medicare beneficiaries who were older than 114 years, as well as “long-stay” or “skilled nursing facility” hospitalizations. Medicare records included date of admission, International Classification of Diseases, Ninth and Tenth Revisions diagnosis codes (ICD-9; ICD-10), age, race, sex, indicator of short stay/long stay/skilled nursing facility, and ZIP code of patient residence. Health outcomes of interest were identified using ICD-9 and ICD-10 diagnosis codes in the billing claim (Supplemental Table E1). Five cause-specific respiratory groupings were considered: all-cause, asthma, chronic obstructive pulmonary disease (COPD), respiratory tract infection (RTI), and all other respiratory diseases that did not fit into the previous categories. We refer to this last group as chronic respiratory disease/respiratory failure (CRD/RF) as nearly 90% of CRD/RF hospitalizations were for diagnoses of respiratory failure and progressive inflammatory and fibrotic pulmonary diseases (Supplemental Table E2). Health endpoints were also distinguished based on whether respiratory events were reported in the principal diagnostic code position (principal diagnoses) or the first three diagnostic code positions (first-three diagnoses). Respiratory hospitalizations that occurred within four days of discharge from a prior hospitalization were assumed to be related to the prior event and excluded from the analysis.
Following previous work (24), we included hospitalizations of beneficiaries living in Metropolitan Statistical Areas (MSAs) within the contiguous US with a total 2010 population of 500,000 or more. MSAs with a population greater than 2.5 million were subdivided into metropolitan divisions, as defined by the US Office of Management and Budget, and analyzed as 120 separate study areas, hereafter referred to as large urban centers (LUCs) (Figure E1).
Ambient meteorological data
Weather station observations of hourly ambient temperature and relative humidity, obtained from the National Oceanic and Atmospheric Administration, were interpolated to census tract centroids using thin-plate spline regressions and aggregated to the daily ZIP code-level (24).
Meteorological variables used in analyses included daily mean (24-hour average), minimum, and maximum temperature, daily mean relative humidity, and daily mean and maximum heat index calculated using the weathermetrics package in R. The analysis was restricted to the warm season, defined as June through September, the four warmest months in all large urban centers (LUCs) during the study period (26).
Statistical Analyses
We applied a two-stage modeling approach to estimate associations between daily ZIP-code level ambient temperature and respiratory hospitalizations within each LUC and nationwide. In the first stage, associations were estimated for each LUC using time-stratified conditional quasi-Poisson regression with distributed lag nonlinear models (27, 28), conditioning on ZIP code of residence, year, month, and day of the week of hospitalization (29). In primary analyses, we defined the exposure-response function with a natural cubic spline and two internal knots placed at the 33rd and 66th percentiles (REFS) of the LUC warm-season temperature distribution based on ZIP code level temperatures within an LUC across all years. The lag-response association was defined using a natural cubic spline with three internal knots equally spaced on the log-scale. A maximum lag window of seven days was chosen a priori to account for the effects of the exposure distributed over time (30–35). All models controlled for Federal holidays, same-day (lag day 0) average humidity using a natural cubic spline with four degrees of freedom, and residual temporal trends using natural cubic splines on day of the warm season (four degrees of freedom). Separate models were fitted for each outcome, and for all-cause respiratory outcomes across the following population groups: sex (male, female); age group (65–74, 75–84, 85–114); and race (Black, white).
In the second stage, we pooled the risk estimates using multivariate meta-regression to obtain nationwide effect estimates (36, 37). The meta-regression included LUC average ambient temperature and temperature range (difference between minimum and maximum temperature) as fixed effects to explain variation across LUCs, as done previously (24, 38, 39). We used the best linear unbiased prediction (BLUP) to update the estimates of the LUC-specific associations (24, 36). Sensitivity analyses assessed whether results were robust to model specification (additional details are available in the supplemental material).
For all models, the cumulative effect is the sum of the relative risk (RR) estimates across lag days 0-6 associated with a specific temperature value relative to the LUC-specific warm-season median temperature (reference exposure). We report lag-response associations, cumulative percent change in hospitalizations (%Δ= (RR-1)*100), and 95% confidence intervals (CIs) comparing high warm-season temperatures (95th ambient temperature percentile) to median warm-season temperatures (40–43).
Attributable Burden Analyses
Attributable burden was quantified at the national, LUC, and subpopulation level using the forward perspective as previously described (44) for the following attributable measures: (1) attributable number (AN), defined as the number of respiratory hospitalizations attributable to temperatures above the 50th percentile; and (2) annual attributable rate (AR), defined as the average annual AN per 100,000 Medicare beneficiaries within subgroup. For the attributable measures, interval uncertainty was obtained empirically through Monte Carlo simulations (44) and are reported as 95% empirical confidence intervals (eCIs). Additional details on the attributable burden calculation are reported in the supplemental material.
RESULTS
Descriptive results
In an open cohort of 24.4 million Medicare beneficiaries aged 65-114 living in the 120 largest US LUCs during 2000-2017, we identified 3,275,033 hospitalizations with primary respiratory diagnoses and 8,374,625 hospitalizations with a respiratory diagnosis in the first three diagnostic code positions (Table 1). Pneumonia accounted for over 94% of all respiratory tract infections (Supplemental Table E1). Black beneficiaries comprised the second largest racial/ethnic demographic in the study cohort. Although other race/ethnicities were reported, daily hospitalizations among Asians, Hispanics, and North American Natives were not sufficiently numerous for separate analyses by these race/ethnicities. Cause-specific respiratory hospitalizations had similar proportions by age, sex, and race, except for asthma, which had proportionally more women and Black beneficiaries compared to other outcomes.
Across the study period and all LUCs, warm-season daily mean temperatures ranged from 0.1°C to 42.6°C (Supplemental Table E3). LUCs in the Pacific and Mountain US Divisions experienced the widest warm-season temperature ranges, while temperatures were higher, on average, in LUCs of the South Atlantic and West South Central Divisions (Supplemental Table E3). Across all locations, daily mean warm-season temperatures were highly correlated with other temperature metrics (Supplemental Table E4). Spearman correlations between warm-season mean temperature and warm-season relative humidity within each LUC were more variable and ranged from −0.71 to 0.10 (mean [SD] ρ= −0.28 [0.25], Supplemental Table E4). LUC average warm-season temperature and temperature range helped explain the heterogeneity across locations, as meta-regression I2 statistics were <15% for all outcomes and all population subgroups (Supplemental Table E5). LUC-specific daily mean temperature summary statistics are reported in Supplemental Table E6.
Nationwide Risk and Burden
Nationwide, 7-day cumulative associations at high warm-season temperatures resulted in a 1.2% (95% CI: 0.4, 2.0) increase in hospitalizations for primary diagnoses of all-cause respiratory diseases, primarily driven by increases in RTI [1.8% (95% CI: 0.6, 3.0)] and CRD/RF hospitalizations [1.2% (95% CI: 0.0, 2.4)] (Table 2, Figure 1). For these outcomes, associations monotonically increased with increasing temperatures above the 50th percentile of the warm-season temperature distribution (Figures 1-2; Supplemental Figure E2). We did not observe associations between high warm-season temperatures and asthma or COPD hospitalizations (Figure 1; Supplemental Figure E2). For all outcomes, principal diagnoses had a stronger association with high warm-season temperatures compared to the first-three diagnoses (Figure 1-2; Supplemental Figure E2). Given this finding, only principal diagnosis outcomes were considered in subgroup and LUC specific analyses.
For all-cause respiratory diseases, RTIs, COPD, and CRD/RF, the pooled lag-response associations demonstrated a robust increase in risk on lag day 0, followed by a decline in risk on lag day 1 possibly due to morbidity displacement, and weak or slightly negative risk across later lag periods (Figure 2; Supplemental Figure E2). For all-cause respiratory hospitalizations, we observed a 2.1% (1.5, 2.8) increase in hospitalizations on lag day 0, with attenuation on the following 6 days (Figure 2). Across all study locations and years, we estimated 11,710 (95% eCI: 8,290-14,670) excess all-cause respiratory hospitalizations due to warm-season temperatures above the 50th percentile (Table 2). Results from Sensitivity Analyses are reported in the supplemental material.
Subgroup-Specific Risk and Burden
We did not observe effect measure modification in relative risk (RR) by age group, sex, or race/ethnicity at the national level (Table 2; Supplemental Figure E5). Due to very few observed hospitalizations among Black beneficiaries, we were unable to estimate associations for three LUCs (Provo-Orem, UT; Boise City, ID; and Portland-South Portland, ME). The burden of heat-related respiratory hospitalizations was largest among beneficiaries aged 85 and over, as reflected in the ARs (Table 2). Attributable burden at the national-level was similar among males and females (Table 2).
LUC-Specific Risk and Burden
We observed considerable geographic variation in the magnitude and direction of heat-related relative risk (RR) and attributable risk (AR) across the study areas. For heat-related all cause respiratory diseases, percent change in hospitalizations ranged from a 3.3% decrease to a 4.7% increase, and ARs ranged from −7.6 to 11.9 annual excess hospitalizations per 100,000 beneficiaries (Figure 3). Across the nine US Census divisions, the highest attributable burden rates were observed in LUCs located in the Pacific, East North Central (a subregional division of the Midwest), New England, and Mid-Atlantic (Figure 3; Supplemental Figure E6; Supplemental Table E8). Although the 58 LUCs in these divisions only accounted for 56% (1,826,777) of the total number of respiratory hospitalizations, 98.5% (11,539) of all heat-related excess respiratory hospitalizations between 2000-2017 were reported in these LUCs. The lowest RRs and ARs were observed in LUCs in the South Atlantic region, especially Florida (Figure 3; Supplemental Figure E6; Supplemental Table E8). RTI and CRD/RF accounted for most of the heat-related respiratory hospitalization burden across all LUCs (Supplemental Figure E7).
At the metropolitan-level we observed a disproportionate burden of heat-related all-cause respiratory hospitalization across age groups, sex, and race. In contrast to the overall population, the highest burden among Black beneficiaries was observed in LUCs of the South Atlantic and East South-Central divisions (Figure 4; Supplemental Figure E8). In the mid-Atlantic and Pacific divisions, the burden was more evenly distributed among white and Black beneficiaries. In the Mountain, West North Central, and New England divisions - which have much fewer Black beneficiaries – the white population accounted for most of the burden (Supplemental Figure E8). Across age groups, the 85 and over population accounted for most of the burden in almost every metropolitan area (Supplemental Figure E9). On average, burden rates were similar between males and females across most US divisions. However, we consistently observed higher ARs among males living in LUCs of the Mountain and West North Central divisions (Supplemental Figure E10).
DISCUSSION
In this nationwide study, we characterized the relationship between short-term exposure to high warm-season temperature and cause-specific respiratory morbidity among older populations living in the 120 largest US urban centers. Notably, we observed that elevated ambient temperature led to excess hospitalizations for respiratory failure and chronic inflammatory and fibrotic diseases such as pneumonitis, pleural effusion, and interstitial pulmonary diseases.
Across demographic groups, most of the burden was attributed to the oldest age group (85+), and most heat-related excess hospitalizations were due to exacerbation of RTI and CRD/RF. Finally, demographic and LUC results suggest heat-related health disparities were likely driven by complex social and contextual forces acting on individuals to intensify their exposure and respiratory health risk.
To our knowledge, this study is the first large population-based study in the US to define a relationship between ambient high temperature and hospitalizations for respiratory failure and chronic inflammatory and fibrotic respiratory diseases (pneumonitis, pleural effusion, and interstitial pulmonary diseases). A nascent literature base on the biological effects of high ambient temperature on respiratory injury points to a few key mechanisms involving inflammatory responses and epithelial barrier permeability that may act independently or in concert to initiate and sustain pulmonary injury (11, 16, 17, 45, 46). Possible mechanisms being investigated include activation of specific heat-shock proteins, which have been associated with respiratory disease progression of both asthma and fibrosis and may function through dysregulation of tight junction permeability. However, so few mechanistic studies have been performed that the collective biological effects of hot, humid air on the respiratory tract are still to be discovered.
The results of this study underscore the importance of considering geographic and socio-demographic factors when assessing heat-health risks. At the national level we did not observe effect modification by sex, age, or race; however, LUC and regional analyses suggest a disproportionate burden of heat-related respiratory hospitalizations among Black beneficiaries compared to their white counterparts across LUCs in the South Atlantic and East South Central Divisions. In these locations, it is possible that structural racism manifests differently than in other parts of the country and intersects with physiological and contextual factors (e.g., comorbidities and cumulative impacts of environmental and social stressors) to disproportionately burden Black beneficiaries. The increasingly warming climate is likely to deepen existing disparities resulting from past and present structural racism. We also observed higher heat-related excess hospitalization rates among beneficiaries of advanced age (85 years and older) compared to other age groups in most LUCs. While age-related susceptibility to extreme heat is well established, this study is among the first to demonstrate the disproportionate burden of heat on the respiratory health of older adults across nearly every major US metropolitan area.
In the overall Medicare population, risks were generally lowest in the Southeast and Southwest, possibly due to higher prevalences of air conditioning (47, 48) and greater physiological adaptation to heat. Findings of lower heat-related mortality and morbidity risk in communities with hotter summers compared to milder summers are commonly reported in multi-city, multi-national studies (4, 15, 43). However, the southwestern and southeastern US geographies are expected to experience more severe extreme heat over the coming decades and human heat acclimatization capacity has a physiological ceiling that may not be able to overcome projected temperature increases without behavioral and technological adaptations (12, 13). This study also identified LUCs in the Pacific, the East North Central, New England, and the Mid-Atlantic that are particularly vulnerable to the effects of extreme heat exposure, with high risk and high burden rates among the Medicare population (Figure 3). On average, LUCs in these US divisions have historically lower prevalence of central air conditioning (47, 48), mild summers, and a wide range of warm-season temperature (Supplemental Table E4, Supplemental Table E6). In addition, a recent study suggests that high temperature warning thresholds may be too high in colder US locations and may not take into consideration the sensitivity of the local population (49).
Overall, the nationwide pooled risk estimates for all-cause respiratory disease were in the range of risk estimates reported in other heat-hospitalization studies among Medicare beneficiaries. In this study, we observed a 1.22% increase in all-cause respiratory hospitalizations across lag days 0-6. Other Medicare based studies, with similar cumulative lag considerations and outcome definitions, report heat-related increases in all-cause respiratory hospitalizations ranging from 0.0-4.0% for a change in the previous week’s ambient heat (6, 15). In contrast to previous studies, we did not detect a cumulative association between heat and COPD. However, we observed elevated COPD exacerbation risk on lag day 0 [1.70 (1.04, 2.36)] with a robust negative association possibly due to morbidity displacement on lag day 1 [−1.47 (−2.12, −0.80)], and no association across later lag periods. Differences in study locations, study-period, heat exposure metrics, the functional form of temperature, and how we defined a health relevant change in temperature (95th percentile to median temperature) could explain differences in reported risk estimates.
The ability to evaluate heat-related respiratory hospitalizations among older adults across 120 LUCs was facilitated by an 18-year study period with rich patient level data and fine-scale meteorological data. However, additional considerations should be acknowledged when interpreting results. First, this study is not able to account for air conditioning prevalence and other potentially important housing characteristics due to limited data availability. Populations that do not have access to air conditioning or cannot afford to adequately cool their homes likely experience higher levels of heat exposure and could be at greater risk. As other studies have indicated, heat-related mortality and morbidity risk may be in competition in some locations.
Thus, high rates of heat-related mortality could artificially create a protective effect when examining heat-related hospitalizations (6, 8, 50). Given the use of spatially interpolated weather station data aggregated to the ZIP code level, there could be exposure misclassification error that would likely result in attenuation toward the null. Finally, reported risks for the entire Medicare population are largely driven by the strength of the association among white beneficiaries living in large LUCs (Supplemental Figure E10) and may not be generalizable to other races/ethnicities or to rural populations.
CONCLUSION
This study considerably extends the current understanding of the relationship between heat and respiratory hospitalizations among older adults by characterizing novel respiratory endpoints that have not been previously examined in a large, multi-city epidemiologic study. The findings raise new questions on the role that heat may play in relation to respiratory failure and airway diseases like interstitial lung disease, pneumonitis, and pleural effusion. Additionally, we observed disparities in the impact of heat on the respiratory health of Black beneficiaries that were geographically dependent. These results suggest that the drivers of heat-related respiratory morbidity are likely due to the joint effects of physiological susceptibility and contextual forces acting on individuals to intensify their exposure and risk. Additional studies to further explore these disparities, as well as the effects of heat on understudied airway diseases, are needed to corroborate our findings. Mechanistic studies, especially controlled exposures, will be essential to determine whether the associations identified in this study are causative in nature. Findings from this study and similar ones can be used at the local level to motivate intervention studies/efforts and repeated studies can aid in monitoring progress towards reducing heat-related health disparities.
Data Availability
Medicare hospitalization data are restricted by IRB protocols and data use agreements with Centers for Medicaid and Medicare (CMS) but other researchers may obtain the same Medicare data directly from CMS. Daily meteorological data are publicly available from the National Oceanic and Atmospheric Administrations National Climatic Data Centers Global Surface Summary of the Day database. American Community Survey (ACS) data are publicly available through the Census Bureau.
Supplemental Material
ADDITIONAL DETAIL OF THE METHODS
Study Areas
To investigate the health effects of extreme heat on urban populations, this study included the 120 largest US metropolitan areas (the names of each metropolitan area are reported in Table E6). Metropolitan areas were identified using US Office of Management and Budget definitions for Metropolitan Statistical Areas (MSA). Similar to Cleland et al., 2023(E1), we used the 2010 Decennial Census survey to restrict our analysis to MSAs within the contiguous US with a total 2010 population of 500,000 or more. Metropolitan Statistical Areas with a population greater than 2.5 million were subdivided into metropolitan divisions, according to OMB definitions, and analyzed as separate study areas. For example, the metropolitan divisions, Anaheim-Santa Ana-Irvine, CA and Los Angeles-Long Beach-Glendale, CA are considered two separate study areas in our analyses; however, these metropolitan divisions make up the Los Angeles-Long Beach-Anaheim, CA MSA. We refer to our 120 study areas as Large Urban Centers (LUCs). In our analyses, all residential ZIP codes that fully or partially overlapped with Metropolitan Statistical Area boundaries were included. ZIP codes were excluded if they were classified as Post Office Boxes, businesses, or universities. The resulting study included 11,924 ZIP codes. For regional analyses, large urban centers (LUCs) were grouped according to the nine US Census divisions.
Additional methodological details
Extracting hospitalization records
We extracted patient records and summary files from Medicare billing claims data for short stay, in-patient hospitalizations for any respiratory event reported in the first 3 diagnostic code positions from 2000-2017. Respiratory events were identified using International Classification of Diseases, Ninth and Tenth Revisions codes, specifically, 460-519 (ICD-9) and J00-J99 (ICD-10). We extracted data only for the 11,924 ZIP codes included in our study locations.
Outcomes
“Principal” respiratory hospitalizations were explored based on ICD codes reported in the principal diagnosis code position. “First-three” respiratory hospitalizations were explored based on ICD codes reported in the first three diagnosis code positions. We created several cause-specific outcomes using ICD codes (Table E1).
Warm-Season
We defined the warm-season as June-September as these were the four warmest months on average from 2000-2017 in each study location. To account for N/A values when creating the 7-day lags, we utilized hospitalization data from May 25 to September 30th.
Reference exposure
The LUC-specific median temperature was used as the reference exposure. We defined a health relevant change in temperature as a change from the 50th to 95th percentile of the warm-season temperature distribution. In this study, we used median temperature as the reference exposure instead of the minimum hospitalization percentile to facilitate risk communication and to allow for easier comparisons across outcomes, demographic groups, and alternative exposure metrics as the minimum hospitalization percentile varied across these factors (Table E3).
Statistical Analyses
Stage 1 – Large Urban Area specific models
In the first stage of our statistical approach we implemented time-stratified conditional quasi-Poisson regression for each location matching on day, month, year, and ZIP code of Medicare beneficiary residence
Conditional Poisson modeling was performed for each study location for each outcome and for each subgroup
Non-linear temperature exposure-response associations were explored by placing 2 internal knots at 33rd and 66th percentiles of the warm-season temperature distribution
Cumulative effects were explored across lag days 0-6 using DLNM with 3 internal knots equally spaced on the log scale
We specified quasi-Poisson (family=”quasipoisson”)
We excluded empty strata (subset=keep)
We excluded missing data (na.action = “na.exclude”)
The general form of the Stage 1 model for each study location is as follows: where, Yzt indicates the daily number of hospitalizations in ZIP code z on day t, and t indexes the event (hospitalization) and control days. ωk denotes the indicator variables that distinguish the case–control sets for each ZIP code, x is the total number of case–control sets (11,924), and δk denotes parameters specific to the case–control sets (which are not estimated in conditional Poisson regression). cb(tempz,t,l) is the crossbasis function for daily ambient temperature in ZIP code z at time t and lag l (in this case 7 days). Other model covariates included a natural cubic spline with 4 degrees of freedom for same day relative humidity (rh), a natural cubic spline with 4 degrees of freedom for day of the warm-season (DOWS), and a binary indicator variable for Federal holidays.
Stage 2 - Multivariate meta regression pooling 120 LUC-specific estimates
In the second stage of our statistical approach, estimates from LUC-specific stage 1 models were reduced to the cumulative risk during the lag period for each location.
We pooled 120 LUC-specific estimates using multivariate meta-regression with fixed effects for LUC-specific average warm-season temperature and LUC-specific temperature range to account for possible effect modification on the pooled relationship.
For all models, the cumulative effect is the sum of the relative risk (RR) estimates across lag days 0-6 associated with a specific temperature value relative to the LUC-specific warm-season median temperature (reference exposure). We report lag-response associations, cumulative percent change in hospitalizations (%Δ= (RR-1)*100), and 95% confidence intervals (CIs) comparing high warm-season temperatures (95th ambient temperature percentile) to median warm-season temperatures.
We report the I2 heterogeneity statistic for all meta-regressions (Table E5).
Best Linear Unbiased Predictions (BLUP) were extracted for each city and used to report LUC-specific RRs and 95% CIs that were “shrunk” based on the underlying distribution across all 120 study locations.
Attributable burden and excess hospitalizations
BLUPs were used to estimate the number of excess hospitalizations (attributable number) within each ZIP code over the next 7 days for each day of the series above the ZIP code-specific 50th percentile of the warm-season distribution. Attributable numbers were summed across all days and ZIP codes to derive LUC-level burdens for all outcomes and subgroups, as well as a nationwide burden estimate, for all outcomes and subgroups.
Specifically, the following equation was used:
Where n is the total number of respiratory hospitalizations in ZIP code z on day t. βxt−1,l represents the LUC-specific log relative risk (based on BLUP estimates) between ambient temperature (xt−1) and respiratory hospitalizations at lag l. The attributable burden calculation was restricted to days when temperatures were above the 50th percentile of the warm-season temperature distribution for each ZIP code. Empirical 95% CIs for the burden estimates were obtained through Monte Carlo simulations, using 1,000 iterations and assuming a multivariate normal distribution of the location-specific BLUP of reduced coefficients (E2). To estimate the excess heat-related respiratory hospitalizations, for each metropolitan area we summed ANzt and the Monte Carlo estimates across all days and ZIP codes within each metropolitan area. To derive the nationwide estimate, we summed ANzt and the Monte Carlo estimates across all days and all ZIP codes.
The total attributable number for each LUC and all 120 LUCs was divided by 18 and by the annual average number of beneficiaries in each LUC and across all locations to estimate an annual attributable rate that is adjusted for the Medicare beneficiary population. Annual attributable rates are expressed as the average annual number of excess hospitalizations per 100,000 beneficiaries overall and within each subgroup.
We calculated empirical confidence intervals (eCI) using Monte Carlo simulations assuming a multivariate normal distribution of the BLUP of reduced coefficients.
Sensitivity Analyses - Approach
Sensitivity analyses were performed for three purposes
To determine appropriate temporal and confounder control;
To evaluate whether results were robust to model specification; and
To evaluate whether results and interpretations were consistent across different exposure metrics.
To determine appropriate temporal and time-varying confounder control, we evaluated whether inclusion of time splines during the warm season improved the model fit. We also explored the impact of inclusion of an interaction between year and day of the warm season (used to capture between-year differences in previous studies) (E3, E4). Finally we evaluated model performance using humidity with lagged effects verses same-day humidity.
To evaluate whether results were robust to model specification, we examined knot placement in the exposure-response function comparing models with two evenly spaced internal knots to models with one internal knot placed at either the 50th or 75th percentile of the warm-season temperature distribution, or two internal knots placed at the 50th and 90th percentile. Sensitivity analyses did not explore three internally placed knots. We also estimated associations using unconstrained lag specifications, and constrained lagged effects with two and three internally placed knots.
Finally, we estimated associations across four additional exposure metrics (minimum temperature, maximum temperature, mean heat index, maximum heat index) and lag days 0-3 instead of 0-6 to further explore the sensitivity of our model. When estimating associations between respiratory disease and heat index metrics, we did not include additional control for relative humidity since it was used in the calculation of the heat index. Results of the sensitivity analyses are reported in Supplemental Table E7.
Sensitivity Analyses - Results
Sensitivity analyses demonstrated that results and interpretations from main analyses were robust to changes in knot placement and lag structure (Table E7). Risk estimates were attenuated when defining heat exposure using minimum and maximum ambient temperature and maximum heat index (Figures E3-4; Table E7). Associations were stronger when using mean heat index as the exposure (%Δ: 1.49% [95%CIs: 0.76, 2.23]) and when estimating associations across lag days 0-3 ([%Δ: 1.54% (95% CIs: 0.86%, 2.23%) compared to lag days 0-6 (Figures E3-4; Table E7).
Acknowledgements
The authors would like to acknowledge the important contributions of members of Clinical Research Branch (Wei-Lun Tsai, Corinna Keeler, Kathryn Burns, Riley Short, William Steinhardt, Cavin Ward-Caviness) of the US EPA who have facilitated our ability to examine climate influenced exposures on health. The successes of this manuscript and its future impacts are largely due to the hard work of behind-the-scenes researchers whose combined contribution to this study has been invaluable.
Footnotes
Contributions: CRO, SEC, MWT, LMN, IJ, MER, AGR designed the study and directed its implementation. LMN, AGR, KLH, SEC provided exposure data, analytical design and modelling assistance. CRO analyzed the data. CRO, SEC, MWT, LMN, EMM, MER, IJ, AGR, AJG, and MER interpreted the results. SEC and EMM provided analytical and visualization support. SEC developed many of the R programs used in analysis. CRO, SEC, MWT, LMN, IJ, AGR, AJG, MER wrote the manuscript. All co-authors reviewed and edited the manuscript.
Funding: This research was supported by cooperative agreement, Grant CR84033801, between the U.S. EPA and University of North Carolina-Chapel Hill. The content of this manuscript is solely the responsibility of the authors and does not necessarily reflect the views and policies of the U.S. EPA. Further, U.S. EPA does not endorse the purchase of any commercial products or services mentioned in this manuscript.
Data sharing statement: Medicare hospitalization data are restricted by IRB protocols and data use agreements with Centers for Medicaid and Medicare (CMS) but other researchers may obtain the same Medicare data directly from CMS. Daily meteorological data are publicly available from the National Oceanic and Atmospheric Administration’s National Climatic Data Center’s Global Surface Summary of the Day database (National Oceanic and Atmospheric Administration, 2022). American Community Survey (ACS) data are publicly available through the Census Bureau.
Ethics approval: Ethics approval was obtained by the EPA Institutional Review Board (IRB number: IRB00006222).
Declaration of competing financial interests (CFI): The authors declare they have no actual or potential competing financial interests.