Global observational diagnosis of soil moisture control on the land surface energy balance

An understanding of where and how strongly the surface energy budget is constrained by soil moisture is hindered by a lack of large‐scale observations, and this contributes to uncertainty in climate models. Here we present a new approach combining satellite observations of land surface temperature and rainfall. We derive a Relative Warming Rate (RWR) diagnostic, which is a measure of how rapidly the land warms relative to the overlying atmosphere during 10 day dry spells. In our dry spell composites, 73% of the land surface between 60°S and 60°N warms faster than the atmosphere, indicating water‐stressed conditions, and increases in sensible heat. Higher RWRs are found for shorter vegetation and bare soil than for tall, deep‐rooted vegetation, due to differences in aerodynamic and hydrological properties. We show how the variation of RWR with antecedent rainfall helps to identify different evaporative regimes in the major nonpolar climate zones.


Introduction
Soil moisture (SM) plays a central role in the partition of available energy at the land surface. As the soil dries out, during extended periods without rain, land evapotranspiration (ET) changes from an energy-limited to a water-limited regime [Seneviratne et al., 2010]. The characteristics of this transition toward greater sensible heat (H) as ET diminishes determine the influence of SM on air temperature and precipitation. Precipitation deficits have been linked to hot extremes [Hirschi et al., 2011;Mueller and Seneviratne, 2012;Vautard et al., 2007], with established SM deficits interacting with the large-scale circulation to amplify summertime temperature variability [Haarsma et al., 2009;Miralles et al., 2014;Quesada et al., 2012]. This atmospheric warming can in turn affect SM through changes in humidity deficit, cloud cover, and precipitation [Fischer et al., 2007;Taylor et al., 2012]. Such soil moisture-atmosphere interactions are projected to strengthen under future climate change [Dirmeyer et al., 2013].
Systematic land hydrological and climate biases have been identified in state-of-the-art global climate model (GCM) simulations, particularly in water-limited ET regimes [Mueller and Seneviratne, 2014], which have been related to model SM dynamics [Cheruy et al., 2013]. Models with biases in historical climate simulations often have stronger warming projections [Cheruy et al., 2014;Christensen and Boberg, 2012]. International modeling efforts like the Global Soil Wetness Project [Dirmeyer, 2011] and the Global Land-Atmosphere Coupling Experiment  have made important contributions to our understanding of land surface hydrology and its relation to surface fluxes. However, evaluating how SM controls the energy budget within GCMs globally is still a major challenge, owing to the scarcity and uncertainty of observational data sets of land surface fluxes and SM at the appropriate scale [Badgley et al., 2015;Dirmeyer et al., 2006;Mueller et al., 2011].
In this paper, we present a unique global observationally based diagnostic to evaluate SM effects on fluxes at scales relevant for GCM evaluation. We follow the approach of Teuling et al. [2006] in assessing land surface behavior during rain-free periods or dry downs. However, rather than focus on a limited number of sites for which flux measurements are available, we use satellite data to quantify dry spell behavior at large spatial scales. We consider the evolution of remotely sensed land surface temperature (LST) with respect to nearsurface air temperature (T a ) as an indirect measure of surface energy partition. Expanding on a recent study that applied this approach in Europe [Folwell et al., 2015] (hereafter F15), we assess the utility of our method globally. We examine the dry spell temperature dynamics for the major nonpolar climate zones. In particular, we analyze the sensitivity of the surface warming to land cover and antecedent rainfall and compare the observations with a simple land surface model.

Methods and Data
We use satellite and reanalysis data to quantify the dry spell evolution of LST relative to T a globally. To characterize the evaporative behavior of the surface, we compute the mean rate of surface warming relative to the atmosphere during dry spells, hereafter referred to as RWR (relative warming rate, K d À1 ).
We have analyzed clear-sky daytime LST from the Moderate Resolution Imaging Spectroradiometer (MODIS) on board Terra and Aqua satellites available since 2000 and 2002, respectively. The daytime equatorial crossing time is approximately 1030 LT for Terra and 1330 LT for Aqua. Levels 1 and 2 swath data were acquired from the NASA Land Processes Distributed Active Archive Center [LPDAAC, 2001]. Daytime acquisitions were regridded to a 0.01°equal-angle LST product; where multiple observations exists on any given day and grid box, only the observation with the lowest satellite viewing angle was used. Only cloud-and aerosol-free observations with view angles of less than ±40°were retained [Trigo et al., 2008].
From the 0.01°data set, we computed daily mean LST anomalies at the 0.5°scale. Using anomalies reduces spatial sampling errors when regions are partially cloud-free and allowed us to combine data from Aqua and Terra. The anomalies were calculated relative to LST climatologies on the 0.01°grid, which were constructed separately for Aqua and Terra. The climatological value for each day of the year was computed from linear interpolation between adjacent monthly climatological means. Mean LST anomalies were calculated for four general cover classes provided there were at least 500 0.01°observations per 0.5°grid box. A land cover class was assigned to each 0.01°grid box by grouping the following classes of European Space Agency GlobCover map [Arino et al., 2008]; all classes, bare: 200, grassland: 14, 20, 30, 120, and 140, and forest: 40-100.
Global subdaily T a is needed to derive RWR. We used 3-hourly T a from European Centre for Medium-Range Weather Forecasts (ECMWF) ERA-Interim reanalysis [Dee et al., 2011]. The data were regridded (bilinear interpolation) to a regular 0.5°grid and linearly interpolated to the average satellite scan time of each 0.5°grid box. In the same way as for LST, a clear-sky (i.e., days with 0.5°LST observations) climatology was constructed for the period 2000-2014 from which daily T a anomalies were computed. Recognizing there are uncertainties in this data set, we repeated the analysis with three other reanalyses: Japanese 55-year Reanalysis (JRA-55) [Kobayashi et al., 2015], Modern Era Retrospective-Analysis for Research and Applications (MERRA) [Rienecker et al., 2011], and National Centers for Environmental Prediction Climate Forecast System Reanalysis (CFSR) [Saha et al., 2010].
To identify dry spells, we used three satellite precipitation data sets: Tropical Rainfall Measuring Mission (TRMM) [Huffman et al., 2007], Climate prediction center MORPHing product (CMORPH) [Joyce et al., 2004], and Precipitation Estimation from Remotely Sensed Information using Artificial Neural Networks (PERSIANN) [Hsu and Sorooshian, 2008], available at 0.25°spatial and 3-hourly temporal resolution. We define a dry day when the mean daily rainfall of the data sets is below 0.5 mm at 0.5°resolution and we define a dry spell as a period of at least 10 consecutive dry days. The choice of 10 days reflects the trade-off between having a long enough period to observe LST changes and maintaining a sufficient quantity of dry spells to enable the analysis globally. The combination of three data sets provides a more conservative identification of dry days, reflected by a less noisy evolution of LST during dry spells and higher RWRs. Finally, in the absence of a globally consistent SM data set, we used antecedent rainfall accumulations as a proxy for SM. Our results are based on 30 day accumulations, though periods from 10 to 90 days produced qualitatively similar results.
To derive RWR, we first calculated the composite mean difference between LST and T a anomalies for each dry spell day (LST anom À T a,amon , equation (1)).
LST anom À T a; anom ¼ where j is the dry spell day, n is the number of dry spell events, i a single event with a valid LST observation, LST c and T c a are the climatological LST and T a , and w is the number of 0.01°observations used to construct the 0.5°LST. Values of T a are only selected on days with a valid 0.5°LST value, so LST c and T c a represent clear-sky climatologies. We removed dry spell events with days in which T a was below 10°C at satellite overpass time to avoid artifacts from soil freezing. Then RWR was computed by linear regression of LST anom À T a,amon on days 2 to 10. The first day of each dry spell was discarded because of the possible low numbers of LST observations and precipitation errors [F15].

10.1002/2016GL068178
Each 0.5°grid box was assigned to one of six climatic zones based on the Köppen-Geiger classification [Kottek et al., 2006]: equatorial (As, Aw, and Cw), arid cold (Bk) and hot (Bh), temperate humid (Cf) and dry (Cs), and continental hot/warm summer (Dfa and Dfb; hereafter continental). Mean RWRs for each zone were derived using data from all dry spells in their respective climatic class. For these classes, we verified that composite daily values of incoming solar radiation, wind speed, and specific humidity were approximately constant on days 2-10 according to ERA-Interim.

Theoretical Illustration
To frame our observational analysis, we present a simple model of how SM influences LST during a rain-free period [F15]. In the transition from wet to dry soil, three drying stages can be distinguished. Stage I drying occurs when SM levels are above some critical value (S CRIT ) and ET is independent of SM. During stage II, SM constrains ET, while in stage III, SM is unavailable and ET becomes negligible .
Fluxes computed with a Penman-Monteith approach [Monteith, 1965] are coupled to a simple soil water budget model to simulate LST evolution during dry downs, as a function of initial SM. The Penman-Monteith equation can be rearranged to solve for LST as a function of evaporative resistance [F15]: where α is albedo, S d and L d are downwelling shortwave and longwave radiation, G is ground heat flux, ε is emissivity, σ is the Stefan-Boltzmann constant, ρ a and c p are, respectively, density and heat capacity of air, r a , is an aerodynamic resistance, λ is the latent heat of vaporization, q sat and q a are saturated and specific humidity, and Δ(T a ) is the slope of the temperature-vapor saturation curve. SM acts on ET through the surface resistance term, r s . We assume a linear reduction in surface conductance (1/r s ) for SM below S CRIT (see equation (7) in F15). Constant values of α = 0.25, S d = 450 W m À2 , L d = 300 W m À2 , G = 10 W m À2 , ε = 0.97, and q a = 0.0078 kg kg À1 are assumed for a 10 day dry spell. Drainage is considered to be zero. While F15 assumed a fixed T a , we prescribe here a constant rate of atmospheric warming over the dry spell of 0.25 K day À1 . We treat T a as a boundary condition to avoid the complexities of the coupled land-atmosphere system.

10.1002/2016GL068178
The RWR is calculated here as the slope of the linear regression of the daily values of LST-T a . Assuming that r a remains constant, RWR is proportional to changes in H during the dry spell. The modeled sensitivity of RWR to initial SM is illustrated in Figure 1, taking parameter values for short, shallow-rooted vegetation (r a = 40 s m À1 , S CRIT = 30 mm).
Evaporative fraction (EF) ranges from 1 to 0, from critical soil moisture (S CRIT ) to wilting point (S WILT ). The drying stages and evaporation regimes correspond to the short vegetation (green line). The shaded areas represent the sensitivity of RWR to ±0.1 K day À1 variation of atmospheric warming.
In the wet, energy-limited regime (initial SM > S CRIT ), negative RWRs prevail. Under these circumstances, the imposed atmospheric warming exceeds the increases in LST, indicating very effective surface cooling by evaporation. With decreasing initial SM (stages I to II), there is an increase of RWR, as ET becomes water stressed during the dry spell. In the transitional regime, positive RWRs prevail, indicating a gradual increase of H during dry downs. In this regime, the value of RWR increases with decreasing initial SM due to the nonlinearity of the coupled energy and water balance-small changes in SM force relatively large changes in LST when the soil is drier. However, the maximum RWR occurs for an initial SM above zero, because for very dry initial conditions, the accumulated evaporative losses over the first days of the dry spell take the soil into stage III drying. In fact, under the driest initial conditions RWR values are negative. In these circumstances, LST cannot keep up with the warming of the atmosphere, and SM-driven changes in surface fluxes are dominated by longwave emission. The sensitivity of RWR to the choice of atmospheric warming rate is shown in Figure 1; more rapid atmospheric warming, depicted by the lower edge of the shaded area, decreases values of RWR and vice versa.
The sensitivity of the relationship between RWR and initial SM to other key surface properties is also depicted in Figure 1. Two additional idealized surfaces are considered: bare soil (r a = 80 s m À1 , S CRIT = 20 mm) and tall deep-rooted vegetation (r a = 20 s m À1 , S CRIT = 40 mm). The lower aerodynamic resistance of tall, forest-like vegetation produces smaller RWRs. This simulation also includes a larger SM reservoir, consistent with deep roots. This produces weaker stage II variations of RWR with initial SM, as relative SM depletion is slower. Conversely, a smooth surface with shallow reservoir (bare soil) produces large values of RWR and a strong sensitivity to initial SM.

Results
Globally, there are large variations in the number of dry spell days and events according to our definition, due to the diversity of regional hydroclimates (Figure 2a). During dry spells, positive mean LST and ERA-Interim T a warming rates (not shown) were found for 97 and 93% of the sampled surface, respectively, of which 81 and 76% were significantly different from zero (p < 0.1). These high percentages demonstrate the reliability of the observations for characterizing dry spell behavior. The observed RWRs are presented at 1°resolution in Figure 2b. Positive RWRs were found for 73% of the grid boxes (of which 51% were significantly different from zero, p < 0.1). Since daily variations in solar radiation and wind speed are minimal during our dry spells, the observed surface warming can be associated with water stress; positive RWRs indicate an increase of H during dry spells. To test how robust these results are to our choice of T a data, we compared them with those derived from three alternative reanalyses. Each yielded very similar spatial structure, with spatial correlations (r 2 ) between ERA-Interim and JRA-55, MERRA, and CFSR of 0.71, 0.79, and 0.77, respectively.  [2000][2001][2002][2003][2004][2005][2006][2007][2008][2009][2010][2011][2012][2013][2014] and the corresponding rate of surface warming relative to overlying air (RWR, slope of (dashed line) the linear fit to LST-T a ) in six climatic zones; (b) RWRs per land cover type, stratified by 30 day antecedent precipitation in decile bins; circles are placed at the median precipitation of each bin; filled symbols denote significant differences (t test) between cover types; (c) Geographical distribution of Köppen-Geiger based climatic zones and dominant cover types (black denotes that other than bare, grass or forest is the dominant type) and (inset bar) number of observations per antecedent precipitation bin for each cover type.

10.1002/2016GL068178
Higher positive RWRs were observed for short vegetation (Figure 2c) compared to forest (Figure 2d), reflecting differences in aerodynamic properties and hydrological dry spell behavior. This is consistent with higher theoretical RWRs for short vegetation in the transitional regime ( Figure 1) and greater SM sensitivity for shallowrooted vegetation. Negative RWRs (where the surface warms up more slowly than the air) are evident in Figure 2b for two contrasting SM conditions: wet with a predominantly energy-limited regime and dry reaching stage III drying (see negative RWRs at the wet and dry end in Figure 1). On the wet side, negative RWRs correspond to areas with abundant water during rain-free periods. Clear examples (Figure 2b) are found in the extensively irrigated Indo-Gangetic plains, wetlands in northern Argentina and eastern China, and forest areas of the Amazon and Congo basins (Figure 2d). On the dry side (ET ≈ 0), the arid cold deserts of Central Asia present extensive negative RWRs. The model predicts this outcome for smooth, dry surfaces when synoptic warming is strong, and LST is high (summertime LST values around midday exceed 330 K in these regions). According to ERA-Interim, Saharan and Arabian dry spell atmospheric warming (~0.15 K d À1 ) is notably weaker than over the midlatitude deserts (~0.4 K d À1 ), and thus, RWRs are more positive. Additionally, if there are substantial dry spell increases in ground heat flux, RWR will be further decreased. Figure 3a shows the evolution of LST and T a composite anomalies from days 2 to 10 of all observed dry spells and their corresponding RWRs for the six climatic zones. Each zone displays a consistent increase of both T a and LST as the dry spell develops. The higher day-on-day increase of LST relative to T a results in positive RWRs, ranging from 0.01 K d À1 (arid cold) to 0.08 K d À1 (arid hot). For the arid cold zone, the negative signals observed in Central Asia deserts (Figures 2b and 2e) partially offset positive RWR in the western United States and elsewhere.
According to our simple model, the clear-sky surface warming signal during dry downs mainly depends on initial SM and surface properties. The model (Figure 1) shows low values of RWR for very wet and very dry conditions that can become negative when the air is warming rapidly. For transitional conditions, RWR increases as initial SM decreases in stage II drying and a maximum value is reached when the initial SM approaches stage III. To investigate the sensitivity of our observations to initial SM and land cover type, we stratified dry spells into deciles of 30 day antecedent rainfall and into three broad cover types for each climatic zone (Figure 3b). In principle, the sensitivity of RWR to antecedent rainfall (proxy for SM) can provide unique large-scale observational information about the dominant ET regime.
Considering all land covers first (Figure 3b, blue lines), a range of sensitivities of RWR to initial SM is apparent for all climatic zones. The full variation of RWR across all drying stages (theoretical "bell-shaped" curve, Figure 1) is only observed for continental and temperate humid zones. For these zones, RWRs increase as initial SM rises, up to intermediate conditions (deciles 5 and 6), before falling back with further increases in SM. In these zones, roughly half of the events fall within stages II to III drying. Values for the wetter deciles decrease but remain above zero, indicative of stage II drying. In the case of the drier climate zones (arid hot, arid cold, and temperate dry), only stage II to III drying (i.e., the dry part of the theoretical curve) is evident. In the arid hot and temperate dry zones there is a progressive increase of positive RWRs with wetness. This implies a stronger rate of change of H for higher levels of initial SM. In the driest cases, H reaches a maximum at some point during the dry spell (typically after 3-4 days, not shown), after which there is no further increase. In the arid cold zone, negative RWRs appear in the driest deciles, due to the negative values in central Asia (Figure 2e), and reflect faster warming of the air relative to surface. Finally, in contrast with the rest of the climatic zones, RWR in the equatorial domain showed little sensitivity to antecedent rainfall. This fact is likely linked to high SM availability and a deep root zone, although cloud contamination in the LST observations should be accounted for in tropical forest areas as a possible contributing factor.
The variation of RWR with antecedent rainfall for different land cover types is also shown in Figure 3b. RWRs per cover type were only computed if there were at least 1000 0.5°LST observations per day of dry spell per precipitation decile (inset plot in Figure 3c). Differences between land covers are the result of a combination of aerodynamic effects and hydrological dry spell behavior. Consistent with theory ( Figure 1), grasslands (higher r a ) exhibit overall higher RWRs than forests (lower r a ) (Figures 2c-2e). Interestingly, the shape of the RWR curve with antecedent rainfall varies between cover types (Figure 3b), suggesting differences in drying stage for similar SM levels at the onset of dry spells. For example, continental forests exhibit negative RWRs for the two wettest deciles after reaching maximum values for deciles 4 to 7; this suggests that forest areas which have received at least 90 mm rain in the previous month will be in stage I (unstressed) drying. On the other hand for grasslands in this climate zone, values of RWR are positive for even the wettest decile, implying stage II drying and modest Geophysical Research Letters 10.1002/2016GL068178 water stress. This difference in evaporative regime between forest and grassland is to be expected given the access of trees to a larger root zone soil reservoir . Similar, though weaker, forest-grassland differences can be detected for the temperate humid domain. For the remaining climate zones, Figure 3b provides only hints of expected forest-grassland differences in water stress.
Bare soil composites were examined for the two arid zones. We would expect bare soils to exhibit a stronger signal than vegetation, due to both aerodynamic and hydrological effects. Indeed, bare soil RWRs in the arid hot domain exhibit a steep increase with antecedent rainfall (Figure 3b). However, in the arid cold domain, while the bare soil and grassland curves are roughly parallel, the bare soil values are substantially lower, and mostly negative. While most bare soil observations in this climate zone originate in Central Asia, it is North America which contributes almost all of the forest observations. Comparison between cover types for the arid cold class should therefore be undertaken with care. Here we simply note that this is the only climate zone where values of RWR for forest exceed grassland (apart from the wettest decile). It should also be noted that grid boxes in dry areas classified as grassland or forest may have a substantial bare soil component.

Summary and Discussion
We have developed a unique observation-based diagnosis of SM control on land surface flux partition globally. To derive our diagnostic RWR, we compared how rapidly the land warms relative to the overlying atmosphere during dry spells of at least 10 days. Under clear skies, this is a proxy for the change in H as the soil dries out. We applied the method at 0.5°globally (excluding high latitudes). The sensitivity of RWR to initial SM content and surface properties was investigated for all the major nonpolar climate zones.
In our observed dry spell composites, we found that the land warms up in 97% of the grid boxes, and that in 73%, the land warms faster than the atmosphere. This positive RWR indicates water-stressed conditions and increasing H, at least during some of the dry spells. Theory indicates low values of RWR for the wettest and driest conditions, which can become negative when the atmosphere warms rapidly at the synoptic scale. Qualitatively consistent with this, we found the lowest RWRs for midlatitude deserts (at the dry end) and for tropical forests, irrigated land, and wetlands (at the wet end). Land cover is an important factor in the observed values of RWR. As expected, higher RWRs are found for shorter vegetation and bare soil compared to tall deep-rooted vegetation, due to differences in aerodynamic properties and hydrological dry spell behavior. Furthermore, the variation of RWR with antecedent rainfall provides information on which evaporation regime a particular region lies in climatologically. Consistent with expectations, ET in arid and dry climates remains strongly SM-constrained even following the wettest months. Temperate humid and continental zones on the other hand indicate decreasing SM stress with increasing antecedent rainfall roughly 50% of the time. Indeed for these two climatic zones, different drying stages can be observed depending on land cover type for a given value of antecedent rainfall. For instance, in our results forests in the continental climatic zone appear unstressed by SM during a 10 day dry spell provided the previous month saw at least~90 mm of rain. Conversely, RWR values indicate that grassland remains in a SM limited regime even with monthly antecedent rainfall of~150 mm.
Several sources of uncertainty and limitations for this observation-derived diagnostic should be noted. First, the reanalysis air temperate is not a purely observational product and relatively weak observational constraints in some regions (e.g., deserts) mean that model biases there may become important. Nevertheless, we found overall results robust to the choice of four independent reanalysis products. Second, this approach requires high-quality rainfall data sets to detect dry spells. We found systematic improvements (in terms of stronger warming rates) when using remotely sensed rather than reanalysis rainfall, particularly away from regions where precipitation is predominantly synoptically controlled. Moreover, RWR increased with the addition of a second and third observational data set. Third, to achieve global cover, the analysis was constrained to relatively short dry spells; longer dry spells would provide stronger SM forcing at the expense of poorer spatial coverage. Finally, we recognize that for partially vegetated surfaces, our LST-derived signal is affected by a combination of transpiration and soil evaporation. These processes operate on different time scales and cannot be disentangled here.
Noting these potential shortcomings, we highlight how well the observations, when averaged over enough cases, produce signals consistent with expectations. Given the importance of dry spells for heat waves, this approach should provide a valuable tool for assessing the role of SM in climate models. Other vegetationbased indices [Zscheischler et al., 2015] also help to distinguish different evaporative regimes and could Geophysical Research Letters 10.1002/2016GL068178 complement our analysis. An advantage of our method is its ability to be repeated with GCM simulations. Model RWR can be derived by carefully processing output to replicate the clear-sky sampling of LST during dry spells. This analysis could also potentially help to constrain future climate projections where change is linked to soil hydrology [Stegehuis et al., 2013].