Atmospheric Aridity and Apparent Soil Moisture Drought in European Forest During Heat Waves

Land‐atmosphere feedbacks, in particular the response of land evaporation to vapor pressure deficit (VPD) or the dryness of the air, remain poorly understood. Here we investigate the VPD response by analysis of a large database of eddy covariance flux observations and simulations using a conceptual model of the atmospheric boundary layer. Data analysis reveals that under high VPD and corresponding high temperatures, forest in particular reduces evaporation and emits more sensible heat. In contrast, grass increases evaporation and emits less sensible heat. Simulations show that this VPD feedback can induce significant temperature increases over forest of up to 2 K during heat wave conditions. It is inferred from the simulations that the effect of the VPD feedback corresponds to an apparent soil moisture depletion of more than 50%. This suggests that previous studies may have incorrectly attributed the effects of atmospheric aridity on temperature to soil dryness.


Introduction
Over the past decades, Europe has seen a growing number of heat waves. This is in line with the model projections that heat waves will intensify and become more frequent and that European summer temperatures will increase faster than the mean global temperature (Meehl & Tebaldi, 2004;Rasmijn et al., 2018;Schär et al., 2004). Land-atmosphere feedbacks associated with soil moisture drought conditions are suggested to form the central reason of this summer temperature amplification (Kala et al., 2016;Miralles et al., 2019;Rasmijn et al., 2018;Seneviratne et al., 2010;Teuling, 2018). To better understand the processes driving heat waves and to improve the robustness of climate model projections, it is therefore essential to understand the feedbacks that regulate exchange of water and energy at the land surface. However, while the role of soil moisture in land-atmosphere interaction has been studied extensively, other processes have received less attention.
Soil moisture has been identified as a key factor driving the land surface-atmosphere interaction during heat waves. A lower availability of moisture in the root zone (induced by precipitation deficit or enhanced atmospheric demand for evaporation) reduces evaporation leading to a shift in the surface energy balance toward increasing sensible heat Miralles et al., 2019). A landmark study showed that the increase in future summer temperature variability could be attributed to the role of increasing soil dryness . Other studies have also shown the close link between observed summer temperature extremes and soil moisture conditions (Fischer et al., 2007;Hirschi et al., 2011;Quesada et al., 2012;Miralles et al., 2014;Philip et al., 2018). However, other studies have questioned the importance of soil moisture. European heat waves are often driven by atmospheric blocking circulation that steers hot and dry air northward, a process that can occur on much shorter timescales (order of days) than the development of soil moisture drought (order of weeks to months). From analysis of FLUXNET observations, Teuling et al. (2010) found that at the onset of heat waves, the increase in sensible heat flux is much larger over forest than over grass, whereas the opposite is expected based on the larger rooting depth of forest. Vegetation might thus play an important role in controlling the link between soil moisture and atmospheric conditions. Open triangles indicate negative effects, closed triangles positive effects. Each line style/color depicts a distinct feedback loop. LE is the evapotranspiration, H is the sensible heat flux, and q are the potential temperature and the specific humidity of the atmospheric boundary layer, h is the height of that layer, r s is the stomatal resistance, and SMI is the soil moisture index. On the diurnal timescale, it is assumed that SMI changes are small. trolling the opening of their leaf stomata in response to several environmental variables (including soil moisture and VPD), plants can either enhance and reduce evapotranspiration in response to an increased VPD depending on air dryness and temperature ( Figure 1). Using the simplified land-atmosphere model, van Heerwaarden and  showed that the observational results of Teuling et al. (2010) could be reproduced by assuming that forest stomata close in response to high VPD during heat waves, whereas grasses do not. Sulman et al. (2016) similarly found that observed evaporation in a mixed forest responded to both VPD and soil moisture anomalies with similar magnitudes, but the timescales of VPD response (hours) was much shorter than the soil moisture response (multiple days or weeks). It should be noted that the VPD sensitivity has been found to be variable and species dependent (Hetherington & Woodward, 2003;Gu et al., 2006;Merilo et al., 2018) and that model representations are in disagreement (Massmann et al., 2019;Zhou et al., 2019). This calls for a better understanding of the role of VPD on land surface-atmosphere interaction, in particular during episodes of warm and dry air advection typical for European heat waves.
In this study, we investigate the sensitivity of energy balance partitioning over different vegetation (grass, broadleaf forest, and needleleaf forest) to VPD, shortwave radiation, and temperature. We also quantify the impact of the VPD feedbacks on atmospheric conditions and temperature. We test the hypothesis that the short-term effect of VPD on temperatures can be (mis)interpreted as soil moisture drought due to the similar impact on evapotranspiration. First, to examine the role of land use type on surface energy balance response, we use observations of surface energy fluxes from a large network of eddy covariance sites covering different land use conditions and at sufficient temporal resolution to resolve the diurnal cycle. Second, we quantify the impact of the VPD feedbacks by modeling the diurnal evolution of the atmospheric boundary layer using a conceptual slab model of the coupled vegetation-atmosphere system (van Heerwaarden et al., 2010;Vilà-Guerau de Arellano et al., 2015), which has the advantage over more complex distributed models that initial and boundary conditions can be easily controlled. Furthermore, assumptions and results are more interpretable in a conceptual slab model.

Methods
The sensitivity of energy balance partitioning over different vegetation to VPD is investigated using the FLUXNET2015 database (FLUXNET, 2018). These data, measured at eddy covariance flux towers, were quality checked and gap filled by the FLUXNET community. We selected only sites located in Europe (68 • N to 33 • E to 41 • N to 0 • E) and with grass, broadleaf forest, or needleleaf forest land cover. This resulted in a selection of 30 flux tower sites, with data covering a total of 309 summers in the period 1996-2014. Only summer (June-August) observations that were labeled as "good quality gap filled" were used. Subsequently, we averaged the half hourly measurements between 9:00-13:00 local standard time following Teuling et al. (2010), assuming that the heating of the land surface is maximum during this time period and precedes the daily maximum temperature. Finally, days with precipitation were removed to isolate the plant evaporation response. The resulting data were analyzed using a combination of bootstrapping and locally weighted polynomial regression (LOESS). Bootstrapping was conducted by randomly sampling the averaged half hourly measurements. LOESS was applied to the latent and sensible heat fluxes with incoming shortwave radiation, temperature, and VPD as independent variables. It should be noted that soil moisture is not routinely observed at FLUXNET sites. However, taking into account that annual precipitation of the grass sites is similar to the forest sites (supporting information) and trees generally have deeper root systems than grass, it seems reasonable to assume that forest sites were not more water stressed than the grassland sites. So, the responses of the different vegetation types to atmospheric conditions could be compared.
To disentangle the impact of the VPD response from other feedbacks, we conducted experiments with a conceptual model of the coupled land-vegetation-atmosphere system (van Heerwaarden et al., 2010; CLASS: https://classmodel.github.io/). This model simulates how vegetation and the diurnal boundary layer interact by making use of idealized atmospheric profiles. We maintained early morning relative humidity across a range of early morning temperatures, such that the specific humidity profile changed with temperature (van Heerwaarden . A latitude of 50 • N was used as representative for European midlatitude heat waves. CLASS calculates LE (evapotranspiration), which originates from vegetation, according to Penman-Monteith. LE is largely determined by the surface resistance (r s ), in which the VPD response is included. The surface resistance was parameterized following the original ideas by Jarvis (1976) and the implementation by van Heerwaarden et al. (2010): where r s,min (s m −1 ) is the minimal stomatal resistance, LAI the leaf area index, and n are dimensionless stress functions. These functions account for incoming short wave radiation (SW in ), soil moisture index (SMI), VPD, and temperature (T) (Vilà-Guerau de Arellano et al., 2015). The soil moisture index is a linear scale and ranges 0 (wilting point) to 1 (field capacity). The soil moisture and VPD functions are of particular relevance for this study, since they represent the biotic response of vegetation to soil moisture supply and atmospheric moisture demand. Both functions are taken from ECMWF IFS. The response to soil moisture is where SMI wilt is soil moisture index at wilting point and SMI fc is soil moisture index at field capacity. The VPD function is where g D (hPa −1 ) is the empirical constant that determines the strength of the response to VPD (Balsamo et al., 2009).
To compare the impact of increased VPD with the impact of soil moisture depletion, we created a new variable: the apparent soil moisture index (SMI * ). If soil moisture would enhance stomatal resistance as severely as both correction functions of soil moisture and VPD do together, apparent extractable soil moisture equals where 2 and 3 represent equation (2) and (3). SMI * accounts for both the reducing impact of soil moisture and VPD on stomatal resistance and thereby on LE. So, SMI * represents soil moisture in an alternate reality where VPD has no influence over LE, to give the same LE as one would model with VPD's impact included.
We simulated only the response of forest, because the data analysis suggests that grass responds much less to VPD than both forest types ( Figure 2). Also, the ECMWF (2018) parameters lead to no stomatal response of grass to VPD. We took the values of the model parameters for forest from van Heerwaarden and Teuling (2014) and ECMWF (2018), with g D = 0.03 hPa −1 . We ran model experiments with and without the VPD response for a range of soil moisture and initial atmospheric conditions, and analyzed midday values (9:00-13:00) similar to the FLUXNET analysis.

Analysis of Flux Observations
Using locally weighted polynomial regression (LOESS), we determined the statistical dependency of latent heat to the main environmental variables: shortwave radiation, temperature, and VPD. When these environmental variables are relatively low, grass, broadleaf, and needleleaf forest sites show similar latent heat fluxes. However, the fluxes start to diverge at certain levels. From an incoming short wave radiation of 650 W m −2 (Figure 2a), a temperature of 20 • C (2b), and a VPD of 8 hPa (2c), LE of broadleaf and needleleaf forest becomes substantially smaller than LE fluxes of grass. Although potential evapotranspiration should increase with increasing incoming short wave radiation, the observed LE of broadleaf and needleleaf forest does not show this continuous increase in contrast to grass LE, which shows proportionality to shortwave radiation over the whole range. This suggest a reduction of LE not related to soil moisture, since forests generally have deeper root systems. Above a VPD of 15 hPa, LE of broadleaf forest shows a clear negative rather than positive trend with increasing VPD, pointing at VPD as having an important reducing rather than stimulating impact on LE.
The response to LE to environmental conditions has direct implications for the relations between sensible heat (H) and shortwave radiation (Figure 2d), temperature (Figure 2e), and VPD (Figure 2f). In contrast to the LE fluxes, the H fluxes of broad-and needleleaf forest are substantially greater than the H flux of grass over the whole measured range of environmental variables. Especially with increasing VPD, the difference between both forest types and grass becomes prominent. At 20 hPa the difference is as large as 120 W m −2 . At elevated VPD conditions, forests, and especially needleleaf trees, emit generally more than twice the amount of H flux emitted by grass, thus acting as important sources for atmospheric heat.
The opposite responses of forests and grass are also reflected in Bowen ratio (H/LE). With increasing temperature, forests show increasing Bowen ratio's, whereas grass shows decreasing ratio's ( Figure 2h). Increasing VPD makes the trends of forest and grass even more distinctive (Figure 2i).

Model Simulations
Next, we simulate the interaction between a land surface with dual sensitivity to VPD (broadleaf forest) and the diurnal atmospheric boundary layer. When the VPD response is turned on (g D > 0 in Equation 3), forest increases stomatal resistance and therefore reduces LE for higher values of VPD. Figure 3a shows the reference situation of forest without VPD response. When the soil is close to wilting point (SMI close to 0) and heat wave temperatures prevail, temperatures exceed 32 • C. However, when the stomatal resistance of forest does respond to VPD, temperature increases further: by up to 2.2 • C (Figure 3b). This maximum increase occurs around 0.5 SMI and heat wave temperatures. The direct effect of the VPD response is a decrease of LE (Figure 3c, depicted in green) and indirectly the enhancement of H (Figure 3c, white lines), which explains the temperature increase of Figure 3b. It should be noted that these differences are comparable in magnitude to the observed differences in Figure 2. When the LE reduction of Figure 3c is (erroneously) attributed to soil moisture, it becomes clear that the atmospheric temperature matches a soil that appears much drier when the VPD feedback is not accounted for. In terms of absolute apparent soil moisture, this effect is strongest for wet soils and high morning temperatures (>0.4). When considered relative to soil moisture under conditions of VPD feedback, the effect is strongest for dry soils and heat wave temperatures

Discussion and Conclusion
In this study, we find from analysis of flux observations that under high VPD, forest in particular reduces evaporation and emits more H (all other factors being equal). In contrast, grassland increases LE and emits less H in response to high VPD. Simulations with a slab model of the atmospheric boundary layer show that this VPD feedback can induce significant temperature increases over forest of up to 2 K during heat wave conditions. It is inferred from the simulations that the effect of the VPD feedback corresponds to an apparent soil moisture depletion of more than 50%. This suggests that previous studies that have attributed heat wave temperatures to soil moisture (e.g., Fischer et al., 2007;Vogel et al., 2017) may have overlooked the role of atmospheric aridity on temperature.
Atmospheric aridity can limit LE as long as vegetation can take up soil moisture. During a persistent heat wave, the soil moisture reservoir will be depleted. Since we did not encounter clear declining evaporation rates with increasing radiation, soil moisture drought has likely not been dominant. This is consistent with positive monthly LE anomalies observed during European summer droughts (Teuling et al., 2013). Therefore, the strong VPD response found in this study is characteristic for warm conditions without strong soil moisture depletion, so during the onset of a heat wave or during short warm extremes due to advected heat (Schumacher et al., 2019). During persistent dry summers, reduced LE due to soil moisture depletion will start to interact with the VPD response. The timescales associated with this drying are much longer than the duration of a typical heat wave with values ranging from 2 weeks to over 1 month depending on land cover (Boese et al., 2019;Teuling et al., 2006).
During heat waves, LE can be reduced by concurrent soil moisture drought as well as atmospheric aridity (high VPD) (e.g., Sulman et al., 2016). To exclude that broadleaf and needleleaf forest sites experienced more soil moisture drought than grass sites and therefore showed smaller LE, we compare available moisture. Although soil moisture conditions are not only dependent on precipitation, the latter can be used as an approximation. The annual precipitation data (supporting information) gives no reason to assume a systematic difference in soil moisture conditions between grassland (on average 994 mm year −1 ) and broadleaf or needleleaf forest sites (892 and 915 mm year −1 , respectively). Also, trees generally have deeper rooting systems than grass and can therefore access a larger soil moisture pool. Both annual precipitation and rooting depth make it unlikely that broadleaf and needleleaf forest experienced more soil moisture drought than grass sites. This points at a distinctive response of the different vegetation types to VPD rather than soil moisture.
Vegetation adapted to drought conditions is able to regulate stomatal resistance and close stomata in case of water stress or high VPD (isohydric). Less adapted species continue high transpiration rates even when soil moisture is limited or VPD is high (anisohydric) Tardieu and Simonneau (1998). Grassland is found to be more anisohydric than woody vegetation such as forests Walther et al. (2019). This confirms our results with markedly larger LE over grassland than over forests with similar VPD conditions (Figure 2c). The highest VPD values were measured above broadleaf forests sites. This is partly due to site selection and data availability and also to the natural occurrence of tree species with needleleaf species often occupying colder and northern regions. The difference between the two forest types shows mostly in the relation between H and shortwave radiation, where the mostly darker colored needleleaf trees emit more H.
The stomatal response in CLASS is empirically parametrized following Jarvis (1976), which is used in the majority of numerical weather prediction models. The parametrization strongly simplifies the mechanisms within plants that drive the VPD response. However, these mechanisms are still poorly understood and contradicting hypotheses exist about how plants sense and respond to VPD (Streck, 2003). The Jarvis parametrization is based on four correction functions, which are assumed to be independent, although one of the correction functions is based on VPD and another one on temperature, which strongly covary in reality. Under high-temperature conditions, VPD is mainly governed by temperature (van Heerwaarden et al., 2010) and thus questions the independency. Whereas this inconsistency can lead to serious errors on the vegetation scale van de Boer et al. (2014), it has been shown to perform well on the landscape scale under various conditions (e.g., Noilhan & Mahfouf, 1996;van Heerwaarden et al., 2010) and therefore seems to be suitable for our study. Furthermore, CLASS has been demonstrated to reproduce observations of the boundary layer well during conditions typical for heat waves Miralles et al. (2014). The rather generic model simulates single idealized days. The VPD response reduced LE mostly at high temperatures. This confirms our FLUXNET data analysis (Figure 2b). Moreover, CLASS showed that LE can be equally reduced by high VPD as by soil moisture drought, a result also found in a 13 year observational study in the United States (Sulman et al., 2016). This emphasizes the importance of atmospheric aridity, which can be confused with soil moisture drought, and may play an even larger role in a future warmer climate. This work used eddy covariance data (https://fluxnet.fluxdata.org/data/ fluxnet2015-dataset/) acquired by the FLUXNET community and in particular by the CarboEuropeIP network. The FLUXNET eddy covariance data processing and harmonization was carried out by the European Fluxes Database Cluster, AmeriFlux Management Project, and Fluxdata project of FLUXNET, with the support of CDIAC and ICOS Ecosystem Thematic Center, and the OzFlux, ChinaFlux, and AsiaFlux offices.