Divergent flash drought risks indicated by evaporative stress and soil moisture projections under warming scenarios

Offline drought indices have been utilizable for monitoring drought conditions, but their reliability in projecting future drought risks is uncertain due to non-meteorological influences on atmospheric water demand (Ep ). This study investigated the impact of surface resistance sensitivity (rs ) to elevated CO2 (eCO2) on projections of future flash droughts (FD). We evaluated FD risks across an East Asian region during a historical period (1981–2020) and a future period (2021–2060) using two evaporative stress index (ESI) series. One series employs the conventional Penman-Monteith (PM) equation for Ep , while the other incorporates a generic rs sensitivity to eCO2 derived from advanced Earth System Models (ESMs). We compared the FD risks identified by the two ESI series with assessments based on soil moisture data from atmospheric reanalysis and multiple ESM projections under two emission scenarios linked with the Shared Socioeconomic Pathways. Results showed that the response of rs to eCO2 has had minimal influences on temporal variations of ESI for the past decades, likely due to its low sensitivity and the masking effects of other environmental factors. However, for the future decades, the ESI projected by the conventional PM equation significantly diverged from soil moisture projections, overestimating future FD risks even under a low emission scenario. We found that incorporating the generic rs sensitivity into the PM equation did not simply resolve the disparity in FD frequencies between ESI and soil moisture projections. Many associated factors contributing to stomatal responses to eCO2 complicate the understanding of future flash drought risks. This study suggests that overreliance on the conventional Ep formula, which neglects non-meteorological effects, could decrease the ability of ESI to detect future FD events under eCO2.


Introduction
Drought, traditionally considered a slow-developing phenomenon with subtle onsets and termination (Wilhite 2000, Mishra andSingh 2010), is now being influenced by anthropogenic climate change, becoming quicker and more intense (Trenberth et al 2014, Yuan et al 2023. Recent decades have witnessed the occurrence of impactful drought events characterized by rapid development, giving rise to a new concept known as flash drought (FD; Otkin et al 2018a).
Despite its short duration, an FD event can cause large agricultural and ecological losses due to its suddenness and unpredictable nature (Ford andLabosier 2017, Jin et al 2019). The scientific communities, thus, have increasingly paid attention to assessing the risks associated with FDs worldwide to formulate proactive measures (e.g. Anderson et al 2016, Apurv and Cai 2020, Liu et al 2020. FDs have been characterized by drought metrics sensitive to changes in atmospheric evaporative demand (E p ) and terrestrial evaporation (E).
Abruptly rising surface radiation and vapor pressure deficit (VPD) can accelerate the moisture flux to the atmosphere and deplete soil moisture (Pendergrass et al 2020). Water stress over land surfaces has been objectively measured by comparing E and E p (Vicente-Serrano et al 2018, Kim et al 2019, allowing the identification of FDs across the globe. Examples include assessments of the rapid soil moisture depletion during the 2019 Australian drought (Nguyen et al 2019), the global frequency and trends of FDs (Christian et al 2021), and the catchmentscale spatio-temporal dynamics of FDs (Li et al 2020) among others. Even solely tracking E p could detect FD onsets in a continental area and their relationship to the atmospheric teleconnections (Parker et al 2021).
While the E p -based indices have been highly utilizable in diagnosing FDs, the definition of E p is still unclear and many different methods are available to calculate it (Vicente-Serrano et al 2020). For instance, E p is traditionally defined as the E rate that should occur with unlimited water supply (e.g. Thornthwaite 1948); however, a regional-scale E p can be substantially adjusted by land-atmosphere interactions even under the same radiation and wind speed (Szilagyi et al 2022). Furthermore, one can employ a semiempirical method oversensitive to changes in air temperature or surface radiation (e.g. Priestley andTaylor 1972, Hargreaves andSamani 1985), potentially leading to incongruous E p estimates and drought trends (e.g. Vicente-Serrano et al 2014, Bai et al 2016).
Though how to define and calculate E p may lead to significantly different FD identifications, the standard FAO-56 Penman-Monteith (PM) equation (Allen et al 1998) has been commonly used in previous E p -based FD identifications (e.g. Basara et al 2019, Christian et al 2021. The physics-based equation can capture changing surface radiation, air temperature, VPD, and wind speed together. While it assumes the idealized grass surface, the effect of heterogeneous land properties could be nullified by normalizing the resulting reference E p , providing a reasonable description of changes in atmospheric water demand across large areas (e.g. Hobbins et al , McEvoy et al 2016, Otkin et al 2018b. Its broad applicability can be partly attributed to the widespread negative correlations between VPD and soil moisture (Zscheischler et al 2018, Zhou et al 2019. As a limited latent heat flux amplifies the sensible heat flux and VPD, an abnormally high reference E p can become an indicator of rapidly developing drought events (Hobbins et al 2016).
However, it is important to consider nonmeteorological factors that can influence land-surface water conditions, but cannot be sensed by changes in the reference E p . Several studies have recently raised a concern that elevated atmospheric CO 2 concentrations (eCO 2 ) could alter plants' stomatal responses and thereby disrupt the global hydrologic cycle (Milly and Dunne 2016, Swann et al 2016, Yang et al 2019, Berg and McColl 2021. The stomatal conductance under eCO 2 has shown to decrease in many experimental studies and meta-analyses (e.g. Long et al 2004, Ainsworth and Rogers 2007, Gimno et al 2016. Although plants' physiological responses to eCO 2 are highly uncertain (McDermid et al 2021, those model and field experiments suggest that non-meteorological effects may significantly alter future FD risks if the reference E p is continuously used for the E p -based indices. The conventional FAO PM equation, which prescribes non-meteorological factors (e.g. a fixed surface resistance; r s ), is unlikely to capture effects of significant environmental changes on E p .
Hence, if FDs are identified by drought indices relying on the reference E p , the associated risks may largely deviate from actual soil moisture conditions in the future. Kim et al (2021) found that future drought risks assessed by the E-E p differences could be significantly offset by the increasing r s with eCO 2 . While such non-meteorological factors may become increasingly important in operational FD monitoring, previous studies have mostly focused on characterizing FD episodes in past decades using the standard FAO PM equation (e.g. Ngyuen et al , Christian et al 2021. In this work, we hence investigated how much the FD risks identified by the reference E p can deviate from indications from soil moisture deficits in the coming decades. We assessed FD risks in the past and the future decades using anomalies of the reference E p and soil moisture datasets from an advanced reanalysis system and Earth System Models (ESMs). Additionally, we employed an alternative PM equation to examine whether considering a sensitivity of r s to eCO 2 can close the gap in the future FD risks assessed by E p -and soil moisturebased indicators.

Reanalysis datasets and CMIP6 projections
We assessed FD risks within an East Asian domain spanning from 15 • to 60 • N and 75 • to 150 • E, for two 40 year periods of 1981-2020 and 2021-2060. The study area is characterized by the monsoonal climate that creates a strong southeast-northwestern vegetation gradient (Peng et al 2018). The southern and the southeastern parts, which are under warm and (semi-)humid climates, are densely vegetated. Conversely, the northern and the northwestern parts have low vegetation due to the cold and semi-arid climates. While monsoonal P is the major control of surface water availability, severe drought event during the past decade have often been triggered by increased heating from anomalous atmospheric circulations (Yuan et al 2015, Ha et al 2020, suggesting that FD events in the study area could increase under global warming. For an observational analysis for 1981-2020, we collected the atmospheric forcing and soil moisture datasets from the state-of-the-art ERA5 reanalysis archive (Hersbach et al 2020). The 2 m air temperature, 2 m dew point temperature, and 10 m wind speed (u 10 ) at 0.25 • × 0.25 • were downloaded from the Copernicus Climate Data Store (https://cds.climate. copernicus.eu/, last access on 11 May 2022). We also collected the near-surface latent heat and sensible heat flux datasets and the top-soil-layer water content data together. All of the hourly datasets were aggregated to the pentad scale to reduce computational costs while being within the subseasonal timescale of FDs (Pendergrass et al 2020).
To project the reference E p from 1981 to 2060, the same atmospheric and soil moisture data were collected from the Coupled Model Intercomparison Project Phase 6 (CMIP6; Eyring et al 2016) archive (https:// esgf-node.llnl.gov/projects/cmip6/, last access on 2 February 2023). We considered the seven ESMs (table S1) that produced daily simulations under the two emission scenarios combined with the Shared Socioeconomic Pathways (SSP; O'Neill et al 2016), namely, SSP1-2.6 and SSP 5-8.5. The SSP1-2.6 is a sustainable scenario that assumes a relatively low greenhouse gas radiative forcing of +2.6 W m −2 by 2100. The SSP5-8.5, on the other hand, is a businessas-usual scenario that assumes continued reliance on fossil fuels and a much higher greenhouse gas radiative forcing increase of +8.5 W m −2 by 2100. The historical simulations and the projections of the six variables were aggregated to the same pentad scale and their grid resolutions were bilinearly unified to 1 • × 1 • .
We first input the ERA5 and the CMIP6 forcing data into the conventional FAO PM equation (Allen et al 1998): where, E p is the atmospheric evaporative demand on the idealized grass surface (mm d −1 ), ∆ is the slope of the saturation vapor pressure curve (kPa • C −1 ) at the near-surface air temperature T ( • C), R n is the surface net radiation minus the soil heat flux (MJ m −2 d −1 ), which can be equated with the sum of latent and sensible heat fluxes, γ is the psychrometric constant (kPa • C −1 ), u 2 is the 2 m wind speed transformed by u 2 = u 10 (2/10) 1/7 , where u 10 is the 10 m wind speed (m s −1 ), and VPD is the difference between the saturation and the actual vapor pressure (kPa).
Since it is formulated with a constant r s at 70 s m −1 , equation (1) could overestimate E p if plants reduce their stomatal conductance under rising eCO 2 (Milly and Dunne 2016). Yang et al (2019) remedied this problem by quantifying a generic sensitivity of r s to eCO 2 . They inverted the PM equation at nonwater-limited locations with CMIP5 ESM projections and related the consequent r s to eCO 2 , suggesting that a reasonable sensitivity of r s to eCO 2 is around 0.09% ppm −1 . By incorporating this value into equation (1), the reference E p can be alternatively estimated by: where, [CO2] signifies the atmospheric CO2 concentration (ppm).
After the temporal and spatial gap filling, we compared the global mean of the [CO2] estimates with observations from the National Oceanic and Atmospheric Administration (https://gml.noaa.gov/ ccgg/trends/, last access on 3 June 2022). Figure  S1 illustrates that the estimated global mean [CO2] agrees with the observations, increasing to 647 ppm by 2060 under the SSP5-8.5 scenario. The spatial variation of the grid [CO2] estimates was also found to be plausible, because it indicates much higher CO2 emissions in major cities in the northern hemisphere. In East Asia, high [CO2] levels were observed in densely populated areas such as eastern China, Korea, and Japan.

Flash drought identification and frequency assessment
In the early 'bucket model' theory (Manabe 1969), actual E is estimated by multiplying Ep to a limiting factor, which is often linearly related with soil moisture availability (Seneviratne et al 2010). Conversely, the ratio of E to Ep (r E ) serves as an indicator of water stress, allowing users to identify anomalously low soil moisture conditions. This concept has proven to be effective in assessing drought risks at local, regional, and global scales (e.g. Anderson et al 2013, Christian et al 2021). However, for r E to purely reflect water stress, it is essential that the other factors influencing E and Ep remain the same. Whereas actual E observed or synthesized by ESMs could be significantly offset by non-meteorological factors (e.g. changes in rs and leaf areas), the reference Ep from the 'offline' FAO PM equation cannot account for them, potentially introducing biases to r E .
To quantify how much the generic sensitivity of rs to eCO2 can change the FD risks indicated by r E , we generated the evaporative stress index (ESI; Otkin et al 2014) using the two Ep estimates from equations (1) and (2) with the collected E data. Simply, ESI is the z-scored r E anomalies as: where, Φ −1 [·] is the inverse function of the standard normal distribution, and F [·] is the cumulative probability density function at each pentad of a year. To preserve the empirical distributions of r E , the cumulative probability was calculated by the nonparametric kernel density function as proposed by Kim et al (2003). We used the normal kernel density function, and the bandwidths were fitted to the empirical distributions of Ep during the historical reference period of 1981-2020.
We also standardized the collected soil moisture data in the same manner, to compare them with the two ESI series as: where, SSI is the standardized soil moisture index that describes anomalous changes in the volumetric soil water content θ.
To capture the variation of SSI with ESI, we determined the optimal timescale for resampling E and Ep. In most locations below 50 • N, the Pearson correlation coefficients (r) between ESI and SSI were highest when the pentad-scale E and Ep were directly used. However, they needed to be resampled on a longer timescale to better capture underlying soil moisture variations at higher latitudes and elevations, meaning that earlier E and Ep may contribute significantly to the current moisture status (figure S2).
We generated ESI series using the spatially changing timescales, and hereafter the two ESI series from the original FAO PM equation (1) and the modified equation (2) will be referred to as ESIo and ESIm, respectively. Using the ESIo, ESIm, and SSI series from the reanalysis ERA5 datasets, we investigated correlations between evaporative stress and soil moisture anomalies and their long-term trends over 1981-2020, as well as frequency of FD events. We then used the CMIP6 projections to assess the FD frequency and areal extent during 2021-2060 relative to 1981-2020.
FD events were identified based on the U.S Drought Monitor operational drought classification (https://droughtmonitor.unl.edu/About/ AbouttheData/DroughtClassification.aspx). As proposed by Pendergrass et al (2020), an FD event was defined as a severe drought (ESI or SSI <−1.3) that worsened by two or more categories within three pentads (about two weeks) and lasted at least another three pentads. To focus on FD risks during cropping seasons, we only counted FD events whose onsets were between April and September. The FD events were assumed to terminate when ESI and SSI values recovered above −1.3.
Since ESI is a monitoring tool of low soil moisture conditions, we evaluated the true positive rate (TPR) of ESI by treating the indications by SSI as true water-stress conditions. The TPR was simply calculated by dividing the number of time steps where both ESI and SSI indicated FDs by the total number of time steps where ESI identified FDs. The TPR of the two ESI series was evaluated at the locations where SSI and ESI were fairly correlated in the historical period. Figures 1(a) and (b) illustrate the Pearson r of ESIo and ESIm to SSI from the ERA5 datasets over the study domain. The Pearson r values exhibited a range of 0.60 ± 0.37 (median ± interquartile range), indicating nearly identical distributions for both ESIo and ESIm. The anomaly correlations were strong below 30 • N and in eastern China, South Korea, and southern Japan, and were also relatively high in some northern parts of China. However, they tended to decline in higher latitudes and elevations, where the infiltration and evaporation processes are influenced by the snowmelt process and low energy availability.

Observational analysis
The distribution of Pearson r was consistent with the studies by Li et al (2017) and Lesk et al (2021), which showed the strong negative correlations between soil moisture and air temperature in the Indian peninsula, Southeast Asia, and eastern China, and the weak correlations around the Tibetan plateau and in high latitudes. The mechanisms underlying the strong correlations differ between dry and wet regions. In (semi-)arid western areas, limited soil moisture can enhance the sensible heat flux, amplifying VPD and Ep. On the contrary, in the wetter eastern regions, increased cloud cover reduces surface radiation and subsequently lower Ep, while soil moisture is replenished by ample P (Li et al 2017). Lesk et al (2021) referred to the two coupling mechanisms as the 'land-atmosphere interaction coupling' and the 'atmospheric circulation coupling' , respectively. The  ESI and SSI series averaged over the locations with Pearson r > 0.5 had similar temporal variations with declining trends (figure 1(c)). The identical temporal variations of ESIo and ESIm also suggest that the generic sensitivity of rs to changing CO2 was not great enough to produce significant differences in Ep with the 40 year rise of eCO2 (∼75 ppm).
Using the nonparametric MannKendall tests, we found that the trends of ESIo, ESIm, and SSI were consistent at the locations with the coupling mechanisms. In figure 2, more than 70% of the chosen locations showed negative trends in the three indices during 1981-2020 (figure 2). Specifically, the ESI and SSI series have decreased in Mongolia, northern and eastern China, Korea, and southern Japan, which suggests increasing FD risks. On the other hand, FD risks have decreased around South Asia, partly due to enhanced moisture transport during monsoon seasons (Roxy et al 2017). The long-term trends of ESI and SSI were found to be −0.056 ± 0.16 and −0.062 ± 0.17 per decade, respectively. The scatter plots between the point-scale trend estimates (figure S3) demonstrate that tracking r E with the reference Ep could provide consistent indications of FD risks, aligning with observations from soil moisture anomalies. The identical trend distributions between ESIo and ESIm also suggest that long-term changes in evaporative stress over the past decades were unlikely to be influenced by the rs sensitivity to eCO2. Figure 3 illustrates the FD frequencies identified by ESIo, ESIm, and SSI for the historical period. The SSI series shows relatively high FD frequencies in the southeastern region, indicating that aberrations in soil moisture levels are primarily influenced by precipitation (P) variability. This finding aligns with Koster et al (2019) and Parker et al (2021). In areas of southeastern China with humid climates, upper-layer soil moisture tends to fluctuate in response to P, often exceeding the field capacity. Consequently, belowaverage soil moisture levels occur frequently due to the high variability of P (Mo and Lettenmaier 2016). In contrast, the FD frequency tended to be low in the arid northwestern part. The normal soil moisture level under the arid climate is already very low, and thus, instances of extreme dryness would be rarely driven by exceptionally long dry spells and/or heatwaves (Mo and Lettenmaier 2015). The decrease in FD risk towards the arid northwestern part of China has been consistently observed in several previous studies (Hoffmann et al 2021, Wang and Yuan 2022).
The two ESI series also indicate that the FD frequencies are high along the southeastern coastline, gradually decreasing towards the northwestern arid areas (figures 3(a) and (b)). In comparison to the SSI (1.3 ± 1.0 events per decade), the ESI series exhibited slightly higher FD frequencies, with 1.5 ± 1.3 events per decade. This difference can be attributed, in part, to two factors: abrupt changes in Ep and slow soil moisture responses. Ep can increase independently of soil moisture changes, driven by factors such as wind speed, thereby triggering FD events on its own (Hoffmann et al 2021). Low soil moisture level can persist in certain regions in subsequent years, leading to a reduction in the number of FD events. Nevertheless, both the ESI and SSI series consistently demonstrate that the areas affected by severe FDs have increased over the past four decades. The mean trend of FD areas identified by the SSI was +0.16% per decade, while the ESIo and ESIm showed trends of +0.17% and +0.15% per decade, respectively. The slight discrepancy in FD-area trends between ESIo and ESIm suggests that the sensitivity of rs to eCO2 has had minimal influence to the areal extent of past FDs.
As expected, the TPRs of the two ESI series in detecting low SSI (0.50 ± 0.43) were insignificantly different. The two ESI series can better detect soil moisture anomalies when their correlations to SSI become stronger (figure S4). They were generally high in eastern China, India, and the Indochinese peninsula where the two coupling mechanisms exist. Despite the relatively high correlations between ESI and SSI, the detectability of ESI could be limited when the frequency of FD events was very low.

Projections of future FD risk under the SSP scenarios
By applying the same methods to the CMIP6 ESMs projections, we generated the ESIo, ESIm, and SSI series at 1 • × 1 • from 1981 to 2060, considering the two warming scenarios. The strongest correlations of ESIo to SSI were found primarily at 1 or 2 pentad timescales during the reference period of 1981-2020, and the timescales of ESI tended to become longer in higher latitudes and elevations (figure S5). The optimal timescales of ESIo to detect low soil moisture seemed to be slightly longer for the seven ESM projections compared to the observational analysis. The spatial patterns of the highest Pearson r between ESIo and SSI, derived from the seven ESMs for 1981-2020, were consistent with the observational analysis. The range of the Pearson r values was 0.57 ± 0.33, there were no notable differences between the two scenarios during the early phase of the projection period. The ESM-generated ESIo and SSI exhibited reasonably strong correlations below 45 • N, with notable emphasis on regions such as India and the southeastern area (figure S5).
At locations where the Pearson r exceeded 0.5, the temporal variations of the median ESIo and ESIm were closely linked to the median SSI (figure 4). Up until 2020, the difference between the median ESIo and ESIm series appeared to be negligible. However, they began to diverge after 2020 under both warming scenarios. The Pearson r of the median ESIo to SSI was significantly decreased in 2021-2060 under the SSP5-8.5 scenario, while the median ESIm was projected to maintain a stable correlation ( figure 4(b)). The divergence between the median ESIo and ESIm could potentially be attributed to enhanced water use efficiency (WUE) in plants' photosynthesis, because the only difference between the two was the rs sensitivity to eCO2 incorporated in ESIm. However, the grid-scale correlations between ESIm and SSI were also projected to decline in many locations, although not as severely as between ESIo and SSI (figure S6). Hence, considering the rs sensitivity to eCO2 alone may not fully close the trend gap between ESIo and SSI.
Using the multi-model median of the seven ESM projections, we found that the FD frequency identified by ESIo and ESIm was 1.5 ± 1.0 events per decade in locations with the coupling mechanisms during 1981-2020, while the SSI series estimated it to be 1.3 ± 0.8 events per decade ( figure  S7). The pattern correlation between the two frequency estimates was 0.70, indicating a reasonably high level of agreement. It is worth noting that the FD identifications during 2015-2020 did not exhibit significant differences between the SSP1-2.6 and SSP5-8.5 scenarios, because the greenhouse gas emissions under the two scenarios had not yet diverged.
However, looking ahead to the coming decades (2021-2060), atmospheric CO2 was projected to increase by 135 ppm from 1981 to 2060 even under the 'sustainability' SSP1-2.6 scenario. This increase leads to a 17% smaller increase in Ep when using the modified PM equation compared to the conventional one. As a result, the ESIm projections yielded a 12% smaller FD frequency (2.5 ± 1.8 events per decade) compared to the ESIo series (2.8 ± 1.8 events per decade). Under the SSP5-8.5 scenario, the differences in FD frequency became more pronounced between ESIo (3.3 ± 2.7 events per decade) and ESIm (2.8 ± 2.3 events per decade), while the ESM-generated SSI projected a less frequent occurrence of FD events in the future (2.5 ± 2.0 events per decade) (figure 5).
While all the three indices indicated that FD risks would increase in over 75% of the chosen locations during the coming decades, there were notable differences in the frequency estimates between ESIo and ESIm due to the constant rs assumed in the FAO PM equation. However, it is important to note that considering the sensitivity of rs to eCO2 did not fully close the frequency gaps between ESIo and SSI. Under the SSP5-8.5 scenario, the largest changes in FD frequency were projected to be 7.8 and 4.5 events per decade for the ESIm and the SSI series, respectively. This suggests that other non-meteorological factors could play a significant role in the FD risks overestimated by ESIo projections.  The annual maximum FD areas identified by the three indices consistently show divergent FD risks in the coming decades ( figure 6). Under the SSP5-8.5 scenario, the ESIo series projected the FD areas expanding to more than 13% in the last decade, while the assessment by SSI mostly indicated values below 10%. The annual FD areas estimated by ESIm fell between the two. To further understand the contribution of climate controls to the interannual variability of FD areas, we conducted a multiple regression analysis. The annual P and Rn data from the seven ESMs were averaged over the coupled locations and then regressed with the medians of maximum annual FD area estimates. The regression models revealed that the two climatic controls only explain 39% of the FD area variation by SSI under SSP5-8.5. On the other hand, they explained 63% and 48% of the FD area variations indicated by ESIo and ESIm, respectively. This implies that the FD risks assessed by evaporative stress are more sensitive to changing climates under global warming compared to the indication by soil moisture projections.
The histograms in figure 7 compare the FD detectability of ESIo and ESIm for the periods 1981-2020 and 2021-2060, when assuming that the identifications by SSI were assumed to be true events. The distribution of the median TPRs from the seven ESMs indicates that the detectability of ESIo has shifted further to the left compared to ESIm (p-value < 0.05). The overestimated Ep can undermine the ability of ESIo to identify actual soil moisture deficiency, even in the locations where the land-atmosphere and circulation couplings exist.

Discussion and conclusions
Our comparative assessments imply that evaporative stress projections may deviate from actual water stress levels, when there are inconsistencies between E and Ep estimates. The CMIP6 ESMs do not use the concept of Ep. Instead, they employ the Monin-Obukhov similarity theory to simulate terrestrial E by directly adjusting rs to changes in various factors such as water stress, VPD levels, and stomatal physiology (Vicente-Serrano et al 2022). To purely indicate surface water stress by ESI, the ESM-generated E should be compared with Ep estimated under conditions where other controlling factors remain constant. However, the FAO PM equation neglects the non-meteorological factors that may alter the atmospheric water demand, and thus potentially introducing biases to future E-Ep ratios (r E ) under eCO2.
Although it has not posed practical problems in detecting land-surface water stress coupled with r E (e.g. Nguyen et al 2019), the historical insensitivity of rs to eCO2 can be attributed to the relatively low eCO2 and the masking effects of other factors (e.g. nutrient deficiency and cold temperatures), despite the potential enhancement of carbon assimilation (Higgins et al 2023). Through a Free Air Carbon dioxide Enrichment (FACE) system combined with a lysimeter, Vremec et al (2023) found that the stomatal conductance decreased over a grass surface with insignificant changes in leaf areas. The ESM-based rs sensitivity to eCO2 of Yang et al (2019) fell within the narrow range from the recent FACE experiment. Hence, it is arguable that the assumption of rs = 70 s m −1 in the reference Ep is outdated. A drought index based on offline Ep estimates could overpredict atmospheric aridity and land-surface drying (e.g Milly and Dunne 2016, Swann et al 2016, Bonfis et al 2017. This issue could become particularly significant in the case of ESI, because the non-meteorological effects would be more pronounced when comparing the reference Ep to E, which is less variable than P. Thus, the rapid soil moisture depletion may not occur as frequently as indicated by ratios between the reference Ep and E under warming scenarios (e.g. Christian et al 2023).
However, it is also crucial to acknowledge that the response of rs to eCO2 alone could not fully resolve the discrepancy between ESIo and SSI projections. This finding aligns with model experiments conducted by Scheff et al (2021), which revealed divergent changes in soil moisture and an offline index even after switching off the stomatal effects of eCO2 in ESM simulations. The projected ESIm, which only considers the rs sensitivity to eCO2, also indicated higher FD frequencies than did the SSI projections. This implies that additional mechanisms should be at play to further reduce E under non-water-limited conditions. These mechanisms include the direct stomatal closure in response to extreme air temperatures and high VPD (Novick et  In addition, plants' physiological responses to eCO2 (e.g. root behaviors, stomatal density and forest mortality among many others) have not been sufficiently understood (De Kauwe et al 2021, Vicente-Serrano et al 2022). While the stomatal reductions in ESM simulations were found to be oversensitive to eCO2 (Forzieri et al 2020), observational studies have suggested that the enhanced WUE under eCO2 was driven by improved photosynthetic efficiency rather than stomatal reduction (e.g. Gardner et al 2023). Since the direct stomatal responses to eCO2 are different across plant functional types and ambient air conditions (Vicente-Serrano et al 2022), the generic sensitivity of rs to eCO2 derived from fractional wet areas is unlikely representative of the heterogeneous land surfaces under diverse climates. Zhou et al (2022) attributed the substantial discrepancies in hydrological projections of ESMs to the limited understanding in fast physiological responses to eCO2.
Thus, it remains challenging to conclude that future FD risks would be less severe than what is indicated by ESIo projections. If the complex physiological responses to eCO2 counteract each other, the radiative effects of eCO2 on Ep could become the primary determinant of future evaporative stress. It is worth noting that the areal extent of projected FDs identified by ESIo was strongly correlated with annual P and Ep. Recent studies have also underscored the growing impact of Ep on the incidence of FDs (Parker et al 2021, Noguera et al 2022. To better detect onsets of future FDs under rising eCO2, there is a compelling need to refine the conventional Ep formula to appropriately account for significant non-meteorological factors.
In summary, based on observational analysis and future projections of flash drought, we can draw the following conclusions: • Inconsistencies between actual evaporation and atmospheric water demand (Ep) can lead to significant deviation of evaporative stress from soil moisture projections, making flash drought events less detectable. • The stomatal response to elevated atmospheric CO2 (eCO2) has had negligible influence on evaporative stress in East Asia, likely due to its low sensitivity and the masking of other environmental factors. • When incorporating a generic stomatal sensitivity to eCO2 in the reference Ep, future flash drought frequency in East Asia was projected to decrease by 12% for 2021-2060 under a low emission scenario, compared to the case that this factor was not considered. • The disparity between conventional evaporative stress and soil moisture projections was not simply reconciled by taking into account the generic stomatal sensitivity. Deep uncertainty remains in plant physiological responses to eCO2, which complicates the understanding of the associated flash drought risks.

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