Large increases in extreme temperature are a robust feature of climate-change projections (Li et al 2021). Furthermore, in most of the United States, there is a greater expected rise in heat waves when considering moisture and temperature together than temperature alone (Barnston et al 2020). This translates to a growing risk of extreme heat stress (Coffel et al 2018, Li et al 2020)—the combination of temperature and other climate variables, like humidity, known to negatively affect the human body's functioning in hot weather—with major implications for human health (Hanna and Tait 2015, Kjellstrom et al 2016, Mora et al 2017). The effect of high temperatures on water and energy supply and demand, ecosystems, and economic and agricultural output can also be relatively larger or smaller according to the co-occurring humidity (Dunne et al 2013, Harpold and Brooks 2018, Haqiqi et al 2021). Even regions with comparatively mild climates can be strongly affected by heat stress due to a lack of physiological, behavioral, or infrastructural preparedness, stemming from resource limitations or perceived lack of need for investments such as air conditioning or heat-resistant materials (Ferranti et al 2016, Guirguis et al 2018). As a result of this conditioning on the historical range, heat stress's impacts (e.g. as measured by wet-bulb globe temperature) tend to increase exponentially in response to small changes near or beyond the high end of its historical values (Wu et al 2014, Cheng et al 2019).

While recent studies have considered projections of extreme heat (Broadbent et al 2020) or heat stress (Barnston et al 2020, Chen et al 2020, Li et al 2020), none have included a focus on how these increases may vary by elevation. This knowledge gap exists despite greater warming at higher elevations having been shown in observations and projections for mountain ranges including the Cascades, Sierra Nevada, Rockies, Alps, and Himalayas (Rangwala et al 2012, 2013, Pepin et al 2015, Rupp et al 2017, Minder et al 2018, Palazzi et al 2018).

Although elevation-dependent warming is relatively smaller in summer than in other seasons, it is significant: by the late 21st century, Western-US daily-maximum summer temperatures are expected to increase by approximately an additional 0.5 °C with each 1000 m of elevation gain, driven by a shift toward drier atmospheric and land-surface conditions and less cloud cover at higher elevations; little to no such effect is seen for daily minimum temperatures (Rangwala et al 2012, Palazzi et al 2018). Regional analyses have revealed the importance of additional factors, such as smaller future temperature changes over oceans affecting nearby coasts and thus creating an elevational signature (Raymond and Mankin 2019). For example, in the Northwest, projected soil drying leads to summer-mean temperatures increasing most around 1500 m, with this increase being less toward sea level (due to persistent marine cooling) and toward higher elevations (due to smaller decreases in soil moisture) (Rupp et al 2017). In the Southwest, the differential temperature increase proceeds fairly linearly from low to high elevations (Rupp et al 2017). Due to the drop-off in specific humidity with elevation, heat stress at higher elevations is more closely connected with heat than moisture, leaving open the possibility that heat stress could also have a positive elevation dependence. Conclusive analyses for elevation dependence of any kind have often been limited by the 5–30 km grid spacing of regional climate models and the sparseness of observational data in mountainous areas.

In the western US—where population centers and important transportation, water, and power infrastructure are situated at elevations ranging from sea level to above 2000 m—understanding current and future elevation profiles of extreme heat stress would aid in assessing the competing effects of projected large extreme-temperature increases (Barnston et al 2020): circulation and vegetation changes lead to strong drying and heightened wildfire risk (Mankin et al 2017, Brown et al 2021), but a weaker North American Monsoon and stronger evapotranspiration lead to increases in water vapor (Pascale et al 2017). In the eastern US, elevation variability is smaller than in the West but absolute humidities and thus baseline heat stress values are higher, creating the potential for notable divergence from regional means in mountain and coastal communities.

Here, we focus on evaluating climatically extreme (99th percentile) heat stress over the US and its future changes under Representative Concentration Pathway 8.5 (RCP8.5) due to three variables: 2 m temperature (T), 2 m water vapor content (q), and surface downwelling solar radiation (r) (Steadman 1979, Moran et al 2001). Elevation dependence in historical values of each variable and in the processes influencing them under climate change leads to uncertainty about their combined effect on extreme heat stress. For example, elevation-dependent warming and enhanced summer drying at higher elevations would cause offsetting effects on heat stress; solar radiation represents an additional factor which affects the rate of heat transfer from the human body (Kjellstrom et al 2016) and is itself affected by vertical profiles of atmospheric thickness and optical depth (Annear and Wells 2007). Motivated by a desire to build physical understanding and lay the groundwork for assessing potential impacts of heat stress on human health and ecosystems, the aim of our study is to determine how temperature, moisture, and radiation interact to shape extreme heat stress under climate change across regions and elevations of the US.

We compare equal-length subsets of the historical (1980–2005) and RCP8.5 (2074–2099) simulations for 20 models from the Coupled Model Intercomparison Project 5 (CMIP5) that have been statistically downscaled and optimized via the Multivariate Adaptive Constructed Analogs (MACA) method to match an observational reference, the gridMET reanalysis (Abatzoglou and Brown 2012, Abatzoglou 2013). This correction is necessary because of model biases in representing historical climates, and while it results in improvements (Abatzoglou and Brown 2012), it does not address process errors, meaning that considerable model-related uncertainty remains (Lorenz et al 2018). The resultant MACA-based dataset is available at daily, ∼4 km resolution, and is a preferred choice for representing regional and subregional extreme events, including in high-elevation and coastal areas (e.g. Jiang et al 2018).

To represent extreme heat stress, we use the environmental stress index (ESI) (Moran et al 2001, 2003), which combines daily-maximum temperature (T), daily-mean specific humidity (q) (converted to relative humidity [RH]), and daily-mean shortwave radiation (r):

$$\begin{align*}{\text{ESI}} &= 0.63T - 0.03{\text{RH}} + 0.002r + 0.0054\left( {T*{\text{RH}}} \right)\nonumber\\ &\quad - 0.073{\left( {0.1 + r} \right)^{ - 1}}\end{align*}$$

for ESI in °C, with temperature T in °C, shortwave radiation r in W m⁻², and RH in %.

The timescales of the variables follows Coffel et al (2018), who argue that daily-maximum T and daily-mean q together approximate daily-maximum heat stress. As a supplementary analysis we also compute ESI using daily-minimum T and daily-mean q, reflecting the additional value this version of the metric may have for heat-stress impacts.

The ESI correlates highly (r ∼0.99) with the heat-stress-relevant but more-difficult-to-calculate wet-bulb globe temperature (Moran et al 2001), while including solar-radiation effects that are absent in wet-bulb temperature and other common heat indices and that are essential to a complete understanding of heat stress (Hanna and Tait 2015, Kjellstrom et al 2016, Gao et al 2017). It also has nearly the same values as wet-bulb globe temperature, whose impacts rise rapidly above 28 °C and especially above 32 °C (Moran et al 2001, Cheng et al 2019). The variables that constitute ESI, therefore, serve as a good approximation of potential impacts to human health, productivity, and energy demand (Kjellstrom et al 2016). However, actual heat-related impacts vary considerably according to population demographics, acclimation, socioeconomic status, physical activity, clothing, and other factors (meaning that the same impacts can occur at lower levels of heat stress depending on population vulnerability); this makes the health relevance of any particular heat threshold only approximate (Kjellstrom et al 2016, Cheng et al 2019, Vicedo-Cabrera et al 2021).

We find that, in the historical period, the MACA dataset exhibits biases relative to gridMET that are generally less than 1.5 °C for extreme T and 0.5 g kg⁻¹ for extreme q; these correspond to approximately 0.75 and 0.5 when normalized against the distribution of the 99th percentile (figures S1 and S2 (available online at stacks.iop.org/ERL/17/064046/mmedia)). Regional and global models have previously struggled to replicate the elevation dependence of observed trends in summer extreme T (Pepin et al 2015), but we find that the MACA dataset is typically at least as consistent with two observational datasets—the ERA5 reanalysis and the global weather stations of the Hadley Centre Integrated Surface Database (HadISD)—as they are with each other (figure S3) (Dunn et al 2012, Hersbach et al 2020).

Consistent with our focus on extreme events, we compute the 99th percentile of daily ESI for each model over the extended summer season (May–September, MJJAS), averaged over 1980–2005 and separately 2074–2099. The selection of the 99th percentile follows previous studies (Chen et al 2020, Vicedo-Cabrera et al 2021), but we do not intend to comprehensively characterize drivers of heat-stress risk; a more epidemiologically oriented study might consider multiple percentiles, or both daytime and nighttime values. For both ESI and its constituent variables, we take the mean of models within the central two quartiles of the distribution, as a balance between representing the model spread and ensuring spatial consistency. We term this quantity the 'extreme ESI' for a time period. Model versions produced by the same modeling center, or whose outputs are highly correlated, have weights reduced such that their sum equals that of a single more-independent model (table S1; Knutti et al 2013). We consider extreme ESI for individual gridpoints and as a mean value over the seven National Climate Assessment regions of the contiguous US (USGCRP 2017): Northwest (NW), Southwest (SW), Northern Great Plains (NGP), Southern Great Plains (SGP), Midwest (MW), Southeast (SE), and Northeast (NE) (figure 1(a)). Results are computed for elevation bins in 50 m increments, with regional means computed in all cases where there are at least 25 gridpoints in a bin. While we use only data from the RCP8.5 scenario, to ease inter-scenario comparisons we include supplemental figures showing changes per degree of change in global-mean surface air temperature (figures S4 and S5; Collins et al 2013).

Zoom In Zoom Out Reset image size

Historical extreme ESI values are highest in the SW deserts, SGP, MW, and SE, consistent with previous work connecting heat stress most closely to extreme q (figure 1(b)) (Raymond et al 2017). Lower values are apparent in the mountains of both the western and eastern US, as is a general latitudinal gradient. Projected extreme heat stress increases are larger at higher latitudes and in interior areas of the western US (figure 1(c)), broadly matching patterns for extreme T (Wuebbles et al 2015, Vose et al 2017, Raymond and Mankin 2019). Coastal areas including California, the Gulf Coast, and Florida experience smaller increases. The influence of elevation on projected changes is moderate and can be seen most clearly in the Great Basin area of the SW.

The projected increases in extreme ESI are more dramatic when the future climate is compared against the historical 99th percentile (figure 1(d)): this threshold is expected to be exceeded on 15% to 60% of future summer days, with a spatial pattern largely controlled by historical variability but also bearing some elevational signature, such as in the central Sierra Nevada and southern Appalachian Mountains. This result adds new detail to previous global-scale studies (Dosio et al 2018), and in the regional means is consistent with findings from dynamically downscaled simulations (Zobel et al 2017).

We decompose the overall pattern of ESI changes shown in figure 1 by considering the contributions from changes in each variable (figure 2). T increases account for 60%–90% of the ESI increase, a proportion which is largest in the NW; the proportion attributable to q increases ranges from less than 10% in parts of the West to more than 40% in the SE; and the r proportion is less than 5% everywhere (figures 2(c), (f) and (i)). We find the largest changes in extreme T, of up to 8 °C, in the MW, NGP, and NW, while coastal areas of California, Texas, and Florida see changes of only around 4 °C, in agreement with prior results (figure 2(a)) (Wuebbles et al 2015). However, the latter regions experience the largest relative increases, due to their smaller historical variability. Few, if any, topographic effects are seen. In contrast, extreme q (figure 2(d)) increases sharply in the eastern US but also in elevation-specific portions of the West that are comparatively drier than the East, such as the California Central Valley and Idaho's Snake River Plain. Extreme radiation has positive changes in the greater Appalachians, northern Rockies, and interior NW, and negative changes in Florida, the Gulf Coast, and southern California (figure 2(f)); in the former regions, model output may reflect connected biases in temperature, cloudiness, and land-atmosphere coupling (Cheruy et al 2014).

Zoom In Zoom Out Reset image size

Historical elevation profiles of extreme ESI illustrate the differences in vertical gradients between regions and demonstrate, such as for the NW, SW, and NE, notable departures from a steady linear drop-off with elevation (figure 3(a)). The moderating effect of water's high specific heat, accentuated by coastal upwelling, is particularly apparent in the lower extreme ESI below 500 m in the NW, and to a lesser degree the SW. Across the US, extreme ESI drops off at a mean rate of around 4 °C per 1000 m, ranging from about a 4 °C difference over only 400 m in the MW to a 6 °C difference from 250 to 2250 m in the NW. This greater elevation dependence in the eastern regions is correlated with their greater summertime RH (figure 3(a), scatterplot). Our hypothesis is that this follows from higher RH increasing the chance of precipitation; that precipitation removes moisture (the most important element for heat stress as per Raymond et al 2017); and this process is more and more likely to have occurred as one moves to higher elevations.

Zoom In Zoom Out Reset image size

Looking more closely at patterns of extreme-ESI warming, we find approximately a 25% variation in the increase across regions at a given elevation, and about a 15% variation across elevations in a given region (figure 3(b)). A latitudinal gradient is apparent, with the northern US having larger increases especially below 750 m, driven by T (figure 2(a); Wuebbles et al 2015). This gradient may be a consequence of confounding between latitude and RH during extreme heat events, i.e. changes in extreme ESI are smaller where extreme ESI is more closely correlated with q—a relationship that broadly strengthens with latitude (Raymond et al 2017).

Echoing previous studies of both means and extremes, we find across regions sharply smaller increases at sea level than at higher elevations (figure 3(b)) (Rupp et al 2017, Zobel et al 2017). The eastern regions (MW, NE, SE) share a relatively small increase at the lowest elevations relative to 250 m (a coastal effect of about 0.25 °C), comparable to that found in Raymond and Mankin (2019) using datasets produced with different downscaling techniques (Pierce et al 2014, Zobel et al 2017). The SW and SGP are expected to see the least extreme-ESI increase, which is especially small at sea level, then rises about 0.4 °C from 0 to 1000 m, and in the SW a further 0.4 °C from 1000 to 2000 m (figure 3(b), red and gold lines) (Rangwala et al 2012). The SW's T and r increases are also the smallest of any region, possibly due to model projections of continuing or even strengthening marine-layer intrusions along the coast that offset any reductions in coastal cloudiness (Lebassi-Habtezion et al 2011). Above about 800 m, where the California marine-layer influence fades away, ESI changes in the SW are more similar to those in the NW and NGP (figure 3(b)). Models tend to underestimate the temperature difference between the interior and the coast on hot days by a factor of about 2 in the eastern US (Raymond and Mankin 2019), suggesting (if this model bias applies also in the western US) that the actual coastal increases may be even smaller than indicated here.

Spatially varying increases of extreme T and q contribute to shaping the elevational profiles of extreme ESI change (figure 4). Changes in extreme shortwave radiation exhibit almost no elevation dependence, so we exclude them from figure 4 for simplicity. For all regions and elevations, increases in T are the primary factor underlying increases in ESI; however, in coastal areas below 250 m (i.e. for all regions except the NGP and MW), q increases have an ESI effect of similar magnitude (figures 2(c) and (f)). Extreme-T increases have relatively greater importance at higher elevations, most notably in the West, and are nearly constant with elevation otherwise, except for a sharply smaller increase near sea level (figure 4). This apparent coastal effect is less noticeable in the SW than the NW; the latter has greater topographically induced localization of cooling to the coast, even on extreme heat stress days (Rupp et al 2017). This phenomenon bears further examination in targeted studies. The greater q contribution at the lowest NGP elevations may be related to these being in the region's east and therefore substantially more humid than the central and western parts (figures 1(a) and S6(b)). In general, q contributions to ESI change tend to be largest near sea level and gradually smaller above 1000 m. However, because our method does not consider the dependence of q on T, it tends to underestimate (overestimate) the contribution of T (q), particularly at higher elevations where there is a greater ratio of warming to moistening.

Zoom In Zoom Out Reset image size

The high-resolution downscaled model output we employ retains biases in complex terrain and near large water bodies, such that higher-resolution regional modeling experiments might yield significant upgrades in accuracy. Land-surface representations are well-known as potential sources of error, especially for extreme events (Cheruy et al 2014). Humidity and radiation are also broadly less reliable than temperature in both observations and simulations (Abatzoglou 2013). To identify patterns of extreme-heat-stress change, we have used only the RCP8.5 scenario due to its large warming signal; however, this limits our study's applicability to less severe scenarios, as there might be physical non-linearities that prevent our findings from being extrapolated. The choice of data resolution also shapes the results; considering ESI computed with daily-minimum temperature yields slightly different values and spatial patterns, with generally a greater increase in low-latitude coastal areas (figures S7 and S8). Framing our study using National Climate Assessment regions also means our conclusions may partly reflect uncorrected-for correlations across gridpoints within a region—for example, between elevation and latitude. We omit separate analysis of any urban effects the models may contain, because of their likely biases and because urbanized areas comprise only 1%–6% of the regional area fractions. But large urban areas considerably affect their local climate (as seen in both observational and model evidence) and tend to be located at low elevations within each region, so this omission may contribute to the non-linearities observed in figure 3(a), particularly in the more heavily urbanized Northeast.

Our results highlight several key attributes of projected extreme heat stress change in the US. The dominant feature is an increase of 3 °C–5 °C, leading to the late 20th-century top 1% summer heat stress occurring on about 25%–50% of summer days by the late 21st century. The majority of this relative increase is due to rising temperatures, although moisture increases are responsible for nearly half in the southern and eastern regions. Elevational variations in temperature and moisture changes add further detail to this picture. Compared to moderate-elevation areas, ocean coastlines experience smaller temperature increases not fully compensated by larger moisture increases. For example, along the Atlantic coast, nearly all models agree that extreme-ESI increases will be greater at 250 m than at sea level and at 500 m than at 250 m (figures 3(b) and S9). In the mountains of the West around 2000 m, extreme heat stress increases are greater than at 1000 m, but this differential increase is either small (approximately 0.3 °C) or statistically insignificant. Although heat stress intensity drops off with elevation, estimates of warming rate and the diagnosis of factors that might enhance it are crucial for accurately understanding risk—which is often connected more closely with relative extremes and rates of warming than with actual values—and for assessing the likelihood of unprecedented extreme events (Ward et al 2020, Fischer et al 2021).

Analysis of high-resolution statistically downscaled model projections indicates large increases in extreme heat stress across the US from the late 20th to the late 21st centuries under the RCP8.5 scenario, and that this increase is generally greater at higher latitudes and elevations. Our results aid in quantifying this broad picture and also newly reveal how heterogeneous spatial patterns of temperature and specific-humidity changes interact to shape it, including identifying hotspots of heat stress warming in the Sierra Nevada, central Rockies, and southern Appalachian Mountains. We find that in many regions large temperature increases cause extreme heat stress to rise at a somewhat faster rate at higher elevations, despite moisture increasing more at low elevations. As heat stress impacts are a function of both the rate of warming and the actual value, locations where one (or both) are especially large present a particular challenge for human health and for natural and managed systems. We find that the elevation dependence of observed heat stress exhibits regional differences that are correlated with coastal proximity and regional humidity, which bear an imprint on current and future values. Further observational and high-resolution modeling work would aid in determining how the regional patterns that we describe vary with local microclimates (especially urban-inflected ones) and specific events, while investigation of the human and landscape context in which intensified extreme heat stress events will occur is a necessary prerequisite for anticipating and mitigating their potential harm.

A portion of this work was carried out at the Jet Propulsion Laboratory, California Institute of Technology, under a contract with the National Aeronautics and Space Administration (80NM0018D0004).

The data that support the findings of this study are openly available. The MACAv2-METDATA dataset can be accessed at https://developers.google.com/earth-engine/datasets/catalog/IDAHO_EPSCOR_MACAv2_METDATA. ERA5 reanalysis is available at www.ecmwf.int/en/forecasts/datasets/reanalysis-datasets/era5, and HadISD is available at www.metoffice.gov.uk/hadobs/hadisd/. Code necessary to reproduce the figures can be found at www.github.com/cr2630git/raymondetal_heatstresselev.

The authors declare no conflicts of interest.