Urban vegetation cooling potential during heatwaves depends on background climate

The capacity of vegetation to mitigate excessive urban heat has been well documented. However, the cooling potential provided by urban vegetation during heatwaves is less known even though heatwaves have been projected to be more severe with climate change. Across 24 global metropolises, we combine 30 m resolution satellite observations with a theoretical leaf energy balance model to quantify the change of the leaf-to-air temperature difference and stomatal conductance during heatwaves from 2000 to 2020. We found the responses of urban vegetation to heatwaves differ significantly across cities and they are mediated by climate forcing and human management. During heatwaves, vegetation in Mediterranean and midlatitude-humid cities shows a significant decrease in cooling potential in most cases due to large stomatal closures, while vegetation in arid cities shows a cooling enhancement with an unmodified stomatal opening likely in response to intense irrigation. In comparison, the cooling potential of vegetation in high-latitude humid cities does not show significant changes. These responses have implications for future urban vegetation management strategies and urban planning.


Introduction
Heatwaves have been observed to increase both in frequency and severity due to climate change (Skinner et al 2018, Perkins-Kirkpatrick andLewis 2020). Nearly one-third of the world's population is currently exposed to deadly heat for at least 20 d per year and an increasing number of people are likely to experience such climatic conditions in the future under all emission scenarios (Mora et al 2017). For people living in urban areas, the threat of this impending heat stress is greater due to localized synergies of heatwaves and urban heat islands (UHIs) in different climates and urban development patterns (Li and Bou-Zeid 2013, Wouters et al 2017, Khan et al 2020, He et al 2021. The case in metropolises could be even worse because of the positive relationship between the level of heat stress and population density (Luo andLau 2018, Zander et al 2018).
An often-proposed mitigation strategy to counteract extreme urban heat is the increase of vegetation in cities (Gaffin et al 2012, Wong et al 2021. It has been well-documented that vegetation reduces the local temperature in multiple ways such as by providing shade and transpiration, which replaces sensible heat with latent heat (Gunawardena et al 2017, Lai et al 2019, Meili et al 2021. However, how urban vegetation behaves during heatwaves when its cooling is most needed is uncertain but crucial for vegetation survival and its ecological benefits in cities. During heatwaves, raised temperature and vapor pressure deficit (VPD) typically increase transpirative demand but induce stomatal closure to avoid water loss, which in turn restrains transpirative cooling (Grossiord et al Figure 1. Example simulation of (a) the relationships between Tv and Ta for one single leaf receiving different levels of radiation (Qa) and (b) how Tv − Ta changes with VPD under different levels of water stress, where β = 0% indicates the leaf is well-watered and β = 50% and 80% represent the level of stomatal closure when compared to a well-watered leaf. When not specified Qa = 1000 W m −2 , Ta = 25 • C, relative humidity = 65%, wind speed = 1.5 m s −1 , atmospheric CO2 concentration = 400 ppm and atmospheric pressure = 1013 hPa. The dashed black line in (a) is the 1:1 standard line. 2020, Kimm et al 2020. This is even more the case when the heatwave coincides with a drought that largely decreases soil moisture (He et al 2022). However, for urban vegetation, irrigation can potentially maintain evapotranspiration elevated and therefore enhance the cooling potential Santamouris 2019, Gao et al 2020). An opposite behavior was also observed in some species showing that plants open or retain their stomatal opening if they receive irrigation or can access water with deep roots (Drake et al 2018, Aparecido et al 2020. A glasshouse experiment in Australasia indicates that this is the case for plants with inherently low stomatal conductance (g s ) and which typically experience droughts because they are likely to be at a higher risk of heat-related leaf damage and therefore, need transpiration to reduce their leaf temperature (Marchin et al 2022). In the US, satellitebased evidence also shows a prevalent enhanced cooling capacity of urban trees during heatwaves (Wang et al 2019) but the magnitude of such enhancements varies across cities. Thus, whether plants provide more or less cooling during heatwaves and the magnitude of cooling they provide are uncertain and they will depend on the environmental conditions (soil moisture, VPD, wind speed, etc) and the extent to which plants adjust their stomatal opening. All these factors vary geographically and therefore require a more extensive investigation across cities in a wide range of climates.
Vegetation regulates its leaf canopy temperature (T v ) through transpiration to maintain functional biochemical and physiological processes (TESKEY et al 2015, Muller et al 2021. Depending on leaf traits and the leaf atmospheric coupling, T v can be higher or lower than air temperature (T a ) (Leuzinger et al 2010, Feng and Zou 2019). The leaf-to-air temperature difference (T v − T a ) is a key variable describing vegetation cooling potential that affects the sensible heat flux between plant surfaces and the air (Moran et al 1994, Muller et al 2021 and changes with meteorological conditions. T v − T a can be computed by concurrently solving leaf photosynthesis, g s , and the energy balance for a given leaf or a canopy. In the following, we show an example of how climatic and environmental conditions such as available radiation, VPD and water availability can influence the relationship between T v and T a using a mechanistic modeling approach (Bonan 2019) (figure 1). In the given example, a well-watered leaf receiving a moderate amount of radiation (Q a = 1000 W m −2 ), where Q a is the sum of net shortwave radiation and incoming longwave radiation, has a theoretically higher T v than T a when T a is lower than 25 • C, above which T v becomes lower than T a and the leaf starts cooling the ambient air ( figure 1(a)). However, when Q a increases to 1200 W m −2 , the transpirative cooling cannot offset the radiative warming and higher T v than T a is observed across all analyzed T a values. The opposite occurs at Q a = 800 W m −2 . All else being equal, T v − T a significantly decreases (more cooling) as VPD increases but the rate of the decrease can be notably suppressed by water stress which limits g s ( figure 1(b)). Although theoretically these mechanisms are well described and can be used to evaluate vegetation's cooling potential (i.e. T v − T a in this study), measuring them in real conditions is challenging, and little is known about how heatwaves modify T v − T a of urban plants in different cities and climates, especially considering contrasting stomatal responses.
To quantitatively derive the cooling potential and stomatal behavior of urban vegetation, we used 30 m satellite-retrieved land surface temperatures (LSTs) and meteorological variables from the ERA5-Land reanalysis product to compute T v − T a during heatwave periods since 2000 across 24 global cities located in climates ranging from humid to semi-arid to arid (figure S1). The value of T v − T a is used to infer g s by inverting a theoretical formulation of the canopy energy budget and to answer the questions: do urban plants close or open their stomata during heatwaves in comparison to normal climatic conditions? Is there a correlation between stomatal behavior and city background climate? Addressing these questions can shed light on the potential of vegetation to provide cooling during heatwaves, when it is likely most needed, and concurrently guide urban greening strategies.

Data and methods
Heatwave days in this study are defined as five or more consecutive with a maximum T a above its 90th percentile in the city during summer (from June to early September) for the period 2000-2020 (figure S2). Normal summer days are defined as the remaining summer days during this period. Days with daily precipitation of more than 2 mm were excluded to remove effects associated with recent rainfall and interception. To analyze the response of urban vegetation to heatwaves, we selected several typical urban green spaces (UGSs) that are present from 2000 to 2020 and are fully covered by vegetation in each city. Most UGSs are fully covered by trees while UGSs in a few arid cities are covered by grasses mixed with trees. The boundaries of the UGSs are determined using multiple datasets including the 10 m European Space Agency (ESA) WorldCover, the 30 m National Land Cover Database (Yang et al 2018), the long-term global land-cover product with fine classification system at 30 m (GLC_FCS30) (Zhang et al 2021), and high-resolution google historical images. The 30 m Landsat thermal images (from all of Landsat 5, Landsat 7 and Landsat 8) were used to retrieve the LST (T v ) of these UGSs and the 9 km ERA5-Land hourly reanalysis product (globally available at hourly scale since 1981) (https://cds.climate.copernicus.eu/#!/home) was used to extract the 2 m air temperature (T a ) that is spatially and temporally matched to each of the T v observations to calculate T v − T a . Other meteorological variables used in this study are also extracted from the ERA5-Land hourly product. More details on the selection criteria of UGSs and description of the T v and T a data are provided in supporting information S1.
Note, while there is a mismatch in the spatial resolution between the two data sets, the spatial heterogeneity of T a in urban areas is generally much smaller than that of T v (Eliasson 1990). For example, data from urban weather stations in Shenzhen show that within-city daytime LSTs ranging from 25 • C to 40 • C only result in a T a spatial variability of 30 • C-31.5 • C (Cao et al 2021). A comparison of the ERA5-Land Hourly T a values with a high-resolution (100 m) simulation product of T a generated for European cities with the urban climate model Urb-Clim (Hooyberghs et al 2019) shows that the maximum difference of T a at 11:00 a.m. local time, which is close to the Landsat overpass times (∼10:30 a.m. local time), on one typical summer day in Paris is less than 1 • C (figure S3). Furthermore, throughout July 2015, there is a high correlation between T a at different resolutions in vegetated areas (R 2 = 0.93, figure  S4) for a humid city (Paris) and semi-arid city (Madrid), selected as examples.
A single leaf absorbs incoming solar radiation R ↓ sw and long-wave radiation L ↓ from the atmosphere and surrounding surfaces and emits long-wave radiation L ↑ as a function of its leaf temperature. The net absorbed radiation is then partitioned into sensible heat H and latent heat λE according to the energy balance as follows where ρ a is the air density, C p is the specific heat capacity of the air, r b is the leaf boundary layer resistance and r a is the where γ ≈ 67 is the psychrometric constant, r s is the stomatal resistance, e sat is the saturation vapor pressure and e a is the actual vapor pressure. Equation (1) can be rewritten as where Q a is the total available energy for the leaf which is the sum of absorbed solar radiation and incoming long-wave radiation, g b = 1/r b is the boundary layer conductance, g a = 1/r a is the aerodynamic conductance, and g s = 1/r s is the stomatal conductance. Here g b and g a follow the parameterization as used in the Urban Tethys-Chloris (UT&C) model (Meili et al 2020). A more detailed description of the calculation of g b and g a is provided in supporting information S2.
Since remotely sensed T v measures the vegetation canopy temperature rather than the temperature of a single leaf, it is necessary to upscale g b and g s from the leaf to the canopy scale. This is done by simply multiplying the conductance g b and g s by the leaf area index (LAI), which was calculated by an artificial neural network trained on a radiative transfer model (PROSAIL) inversion that can predict LAI from the Landsat surface reflectances (Martínez-Ferrer et al 2022). Then, equation (2) becomes Equation (4) is non-linear in T v . To obtain an analytical solution, we use a Taylor's expansion to approximate T v around T a , which gives which gives an analytical expression for T v − T a : For each observation of T v − T a , we numerically solve equation (6) to obtain the value of g s that satisfies equation (6). All other terms in equation (6) are either measured (T v , T a , Q a , VPD, etc) or calculated (g b , g a , etc), based on aerodynamic considerations. In this way, we approximate g s during heatwaves and normal summer days over UGSs in the analyzed cities. All variables/parameters and their units used in this study are listed in table S1.

Results
We first compared the changes in meteorological variables during heatwaves and normal climatic conditions in all 24 selected cities (figure 2). By definition, T a is statistically larger during heatwaves ( figure 2(b)). We also found that VPD significantly (p < 0.05, t-test) increased on heatwave days compared to normal summer days in most cities ( figure 2(d)). Additionally, half of the cities also showed significantly stronger solar radiation during heatwaves ( figure 2(a)). T a can exceed 40 • C for some arid cities such as Baghdad, Dubai and Phoenix with VPD exceeding 6 kPa. Seven cities showed a significantly decreased wind speed during heatwaves ( figure 2(c)). Since heatwave durations in most cities are relatively short (less than 10 d, figure  S5), UGSs only showed a slightly decreased greenness (indicated by normalized difference vegetation index, NDVI) and the decreases are only significant in Houston and Denver (figure 2(e)). A few cities showed a significant increase of NDVI (e.g. Wuhan and Shanghai), which might be related to modified radiative conditions. Despite the relatively consistent changes in meteorological variables, we observed distinct patterns of the leaf-to-air temperature difference T v − T a (figure 2(f)), which in most conditions is a positive number, i.e. vegetation is warmer than the surrounding air. For several high-latitude humid cities including Baltimore, Moscow, London, Berlin, Oslo and Chicago, T v − T a showed almost no change between heatwave and normal summer days, while for Mediterranean cities such as Madrid, Rome, Barcelona and Lisbon, T v − T a increased (although not always statistically significant) on average by 0.8 • C-2.6 • C, when compared to T v − T a on normal summer days. The increasing T v − T a indicates that the cooling potential of plants in these cities decreased during heatwaves as sensible heat increases driven by the positive T v − T a . Some midlatitude humid cities such as Wuhan, Paris and Shanghai also showed an increased T v − T a . In contrast, the average T v − T a decreased by 1.4 • C-4.8 • C in six arid cities (Baghdad, Dubai, Phoenix, Las Vegas, Abu Dhabi and Denver).
Within a given city, the response of urban vegetation to heatwaves is spatially consistent across UGSs (figure 3). For example, the T v − T a of UGSs in Paris universally increased during heatwaves however with different magnitudes in various Figure 2. Comparison of meteorological variables, NDVI and Tv − Ta between normal summer days (green box) and heatwave days (yellow box). The asterisk indicates that a significant (p < 0.05, t-test) difference exists between the mean value of the two groups. According to the Köppen-Geiger climate classification (figure S1), we group these cities into high-latitude humid cities (in purple), predominantly midlatitude humid cities (in dark green), Mediterranean cities (in light green), arid cities (in tan), and cities in other climates (in black).

UGSs (figures 3(a)-(c)
). Small UGSs showed a higher increase in T v − T a compared to larger UGSs and the T v of small UGSs can be up to 10 • C higher than T a on heatwave days in Paris (figure 3(c)). The effect of heatwaves was more pronounced in Madrid, where T v − T a increased up to 20 • C for most grasslands, which are likely wilted or highly water-stressed and almost completely lose their cooling function (figures 3(d)-(f)). However, tree-covered UGSs in downtown Madrid were less affected by the heatwaves compared to the grasslands. This discrepancy is likely related to different rooting depths and also irrigation practices. In contrast to Paris and Madrid, all UGSs in Phoenix, which are mostly grasslands and irrigated, showed a consistently decreased T v − T a during heatwaves (figures 3(g)-(i)).
To explain the different responses of T v − T a to heatwaves across cities, we numerically computed the g s for each observation of T v − T a (figure 4). For the same value of g s , T v − T a is commonly lower during heatwaves than during normal summer days, indicating a higher cooling potential of vegetation in hotter and higher VPD conditions if plants are not water stressed, which is the case where plants are intensively irrigated in arid cities. In this case, results suggest that vegetation keeps its stomata relatively open during heatwaves (Baghdad, Abu Dhabi, Dubai, Las Vegas, Phoenix and Denver).
However, for Mediterranean cities experiencing a hot and dry summer such as Lisbon, Rome, Barcelona and Madrid, T v on normal summer days was prevalently higher than T a , resulting in a positive T v − T a larger than 3 • C. This indicates that plants in these cities are already under water stress on normal summer days ( figure 1(b)). Their g s is generally between 0.2 and 0.3 mol H 2 O m 2 s −1 ( figure 5(a)). During heatwaves, the plants further close their stomata and therefore show a decreased g s of around 0.15 mol H 2 O m 2 s −1 ( figure 5(b)) and raised T v − T a , which approaches 6 • C in some cities ( figure 4).
For humid cities, plants also commonly closed their stomata and show a decreased g s (although not always statistically significant) during heatwave days ( figure 4). However, T v − T a of high-latitude cities (London, Oslo, Berlin, Moscow, Chicago and Baltimore) showed no significant difference compared to normal summer days even when g s significantly decreased (figures 2(f) and 4). Most midlatitude humid cities including Paris, Wuhan, Shanghai and Houston showed increases of T v − T a ranging from 1.1 • C to 2.7 • C. We also found that T v − T a is positively related to heatwave intensity (defined as the degree-hours that are the sum of hours with T a above the 90th percentile of daily maximum T a in the city during summer, weighted by the departure of hourly T a from the 90th percentile) in a few cities (e.g. Moscow, Chicago, Wuhan, Houston, Atlanta, Baghdad, Las Vegas, Phoenix, figure S6) and thus an increase in T v − T a is observed with increasing heatwave intensity. Meanwhile, g s mostly follows the change of T v − T a and shows a negative relationship with heatwave intensity in some cities leading to larger stomatal closure with higher intensity (figure S7). . Dots indicate the observation of Tv − Ta and the corresponding computed gs (green ones for normal summer days and purple ones for heatwave days). Tv − Ta and gs of the two groups are statistically compared in the form of the box-and-whisker plot and the asterisk indicates that a significant (p < 0.05, t-test) difference exists between the mean gs of the two groups. To avoid overlap over individual dots, these boxes were moved down by nine units ( • C) and readers can check the secondary y-axis for specific values. According to the Köppen-Geiger climate classification (figure S1), we group these cities into high-latitude humid cities (in purple), predominantly midlatitude humid cities (in dark green), Mediterranean cities (in light green), arid cities (in tan), and cities in other climates (in black). The unit of gs is converted from m s −1 to mol H2O m −2 s −1 to be comparable with plant physiological literature (supporting information S3).

Discussion and conclusion
We combined 30 m resolution satellite observations with a theoretical leaf energy balance model to analyze how urban vegetation responds to heatwave conditions across 24 global metropolises. Unavoidable uncertainties in LST from remote sensing observations (Cook 2014, Laraby 2017 and T a from reanalysis (Hersbach et al 2020, Araújo et al 2022) are filtered by using a large number of time steps for each city and by fitting a theoretical T v − T a vs. g s model in each time step. Although T a in cities experiencing strong UHIs is likely underestimated in the used dataset, the canopy T a UHI is on average much smaller than the surface temperature UHI (Venter et al 2021), especially during daytime and it is unlikely to generate any major effect on the results. To quantify the potential impact of variations or uncertainty in T a , we also tested the sensitivity of the g s results to the variability in T a by adding an error term to T a for each observation. Even with a large error term of ε ∼ N (3 • C, 1 • C), that assumes a large uncertainty in T a , results show that the T v − T a and g s response to heatwaves (i.e. how T v − T a and g s changes during heatwaves) in each city remain the same (figures 4 and S8). Hence, the changes of T v − T a and g s are robust at the city scale even without accounting for the within-city spatial heterogeneity of T a .
Our results highlight that the response of urban vegetation has a geographical divergence which is largely related to background climate forcing and to vegetation management, i.e. irrigation. In highlatitude humid cities, there are no significant changes in the leaf-to-air temperature difference (T v − T a ) (figures 2(f) and 4), which indicates that the vegetation cooling potential of these cities was not affected by the heatwave conditions. This is because of (a) similar g s during normal summer days and heatwave days in some cities such as London (figure 4) and (b) relatively higher g s (figure 5(a)) and overall less sensitivity of T v − T a to higher g s than to lower g s (see the curvature change of the theoretical relationship between T v − T a and g s in figure 4). For instance, a significantly decreased g s does not change T v − T a in Moscow as lower g s is likely compensated by changes in meteorological conditions (figures 2 and 4). However, in midlatitude humid cities such as Wuhan, Shanghai and Houston, we found significant stomatal closure which led to increased T v − T a during heatwaves. Due to generally abundant summer rainfall in these cities, urban vegetation is mainly rainfed and expected to still provide substantial cooling even during heatwaves or extended dry periods when cooling is likely most needed. However, our results do not support this assumption and suggest that urban vegetation in such cities, with the current management strategies, provides less rather than more cooling under extreme heat conditions. This effect is even exacerbated in vegetation with lower canopy heights, which showed higher T v − T a than that of vegetation with higher canopy heights for similar g s (see an example in Houston, figure S9). This suggests that increased irrigation is likely needed to fulfill the water demand of urban vegetation and to maintain its cooling potential during heatwaves.
Vegetation in the Mediterranean is widely known for experiencing and being adapted to low water availability (Galmés et al 2007, Rana et al 2020. We found that heatwaves exacerbated plant water stress as stomatal closure is more pronounced than during normal conditions and such closure limits transpiration, which is originally already at low levels ( figure 5). During heatwaves, as plants close stomata further, we observe a much higher T v than T a indicating that Mediterranean urban plants are actually unlikely to considerably cool urban air. However, plant T v could be still lower than those of impervious surfaces and their exact cooling potential during heatwaves requires further investigation.
UGSs in arid cities, which are mostly covered by grasses, were found to retain stomatal opening during heatwaves at a similar level than during normal days. Given the higher atmospheric water demand, T v − T a for a similar g s shows lower values during heatwaves than normal summer days in arid climates (figure 4). This contrasts with the reduced green roof cooling potential during heatwave/drought conditions found in other studies (Speak et al 2013, Zhang et al 2020. However, by keeping stomata relatively open, these plants have a considerably enhanced cooling potential, which given the extremely high VPD, is counterintuitive, and not captured by most stomatal conductance models (Leuning 1995, Damour et al 2010, Medlyn et al 2011, Meili et al 2021 that will predict stomatal closure at such VPD levels independently from water availability (Yang et al 2019, Meili et al 2021. This result sheds light on the importance of watering vegetation during extreme weather conditions (Zhang et al 2020).
Regardless of heatwave conditions, T v − T a generally follows the gradient of vegetation NDVI (figures 2(e) and (f)), i.e. dense and healthy UGSs lead to higher cooling potential. However, all our results combined show the critical role of stomatal behavior and background climate in the responses of the cooling potential to heatwaves (figure 4), especially given the remarkable difference of such responses, for example, in Phoenix and Barcelona with both having a decreased NDVI (figures 2(e) and (f)). A large decrease of g s and stomatal closure during heatwaves can directly suppress plant transpiration and increase leaf temperature (T v ), which might prevent any transpirative cooling and puts plants at risk of lethal overheating if they fail to keep T v below the leaf critical temperature, generally 46 • C-49 • C (Hüve et al 2011, O'sullivan et al 2017). Nevertheless, even in the dry Mediterranean cities experiencing severe water stress, plants close their stomata as these critical temperatures are likely not reached even during heatwaves in these cities ( figure 2(b)), or herbaceous vegetation is already wilted. Conversely, when plants have enough water available (in humid and irrigated cities) they keep stomata open even when atmospheric water demand is significant.
Our analysis has unavoidable limitations in terms of data and methodological choices. The rare occurrence of heatwaves and cloud contamination in some cities results in a limited number of observations during heatwave days and induces uncertainty in cooling potential and g s estimation. The use of the firstorder Taylor's expansion can also cause some bias in estimating g s , especially when leaf temperature strongly deviates from air temperature. Beyond this, even though most UGSs chosen in this study are fully covered by trees, in a few arid cities (i.e. Abu Dhabi, Dubai, Las Vegas and Phoenix) UGSs are mainly covered by grasses and contain only a few trees, in which case our modeled g s could have larger uncertainties. Irrigation on urban grasslands can strongly change soil moisture and can lead to evaporation from water intercepted on the grass leaves or from the soil underneath, which can make T v − T a not representative of leaf temperature only. Compared to the low-frequency Landsat observations used here, future studies could focus on daily observations to quantify the change of stomatal behavior throughout the whole heatwave and potentially associated drought period.
In summary, we explore for one of the first times how T v − T a and g s , which are good proxies for the cooling potential of urban vegetation, respond to current heatwaves in different climates and cities. Results highlight the crucial role of human intervention through irrigation in all those climates where vegetation might undergo partial or substantial water stress. Mediterranean and midlatitude humid cities are shown to experience a significantly suppressed cooling potential as plants are mostly rainfed and they largely decrease g s when cooling (transpiration) will be most needed (e.g. during heatwaves). However, cooling enhancement is observed in arid cities because irrigation is largely shifting the behavior of g s leading to equal or higher g s than on normal summer days. As a result, while urban greening is a desirable strategy to achieve multiple ecosystem services (Haase et al 2014, Richards et al 2022), its capability to reduce temperature during heatwaves might not match expectations derived from normal summer days and cannot be thought of disjointly from strategies to supply necessary water requirements. This is becoming even more relevant in a rapidly changing climate.

Data availability statements
All data that support the findings of this study are included within the article (and any supplementary files).