Hot-Moments of Soil CO 2 Efflux in a Water-Limited Grassland

The metabolic activity of water-limited ecosystems is strongly linked to the timing and magnitude of precipitation pulses that can trigger disproportionately high (i.e., hot-moments) ecosystem CO2 fluxes. We analyzed over 2-years of continuous measurements of soil CO2 efflux (Fs) under vegetation (Fsveg) and at bare soil (Fsbare) in a water-limited grassland. The continuous wavelet transform was used to: (a) describe the temporal variability of Fs; (b) test the performance of empirical models ranging in complexity; and (c) identify hot-moments of Fs. We used partial wavelet coherence (PWC) analysis to test the temporal correlation between Fs with temperature and soil moisture. The PWC analysis provided evidence that soil moisture overshadows the influence of soil temperature for Fs in this water limited ecosystem. Precipitation pulses triggered hot-moments that increased Fsveg (up to 9000%) and Fsbare (up to 17,000%) with respect to pre-pulse rates. Highly parameterized empirical models (using support vector machine (SVM) or an 8-day moving window) are good approaches for representing the daily temporal variability of Fs, but SVM is a promising approach to represent high temporal variability of Fs (i.e., hourly estimates). Our results have implications for the representation of hot-moments of ecosystem CO2 fluxes in these globally distributed ecosystems.


Introduction
Water-limited ecosystems cover over 30% of the land-area of the world and changes in their carbon dynamics could have important impacts on the global carbon budget [1,2].It is known that the metabolic activity of these ecosystems is strongly tied to the timing and magnitude of precipitation pulse (PP) events at global [3] and regional scales [4,5].Thus, the study of how water-limited ecosystems respond to different patterns of PPs has been an important research topic for more than 40 years [6], and multiple studies have highlighted the importance of these patterns on ecosystem responses [7][8][9][10].Despite these efforts, it is still unclear how different types of pulses trigger responses of ecosystem processes [11][12][13], and consequently it is challenging to represent these responses using modeling approaches across water-limited ecosystems [14][15][16].
Climate models project an increase in precipitation variability, including more extreme rainfall events followed by longer dry periods in water-limited ecosystems [17].These projections have motivated multiple experiments to better understand the role of changes in frequency and intensity of PPs on ecosystem processes [18][19][20][21].Studies have shown that changes in PPs substantially influence carbon dynamics including net primary production [9], gross primary production [22], net ecosystem exchange [23], ecosystem respiration [24], and soil CO 2 efflux (Fs, as a result of heterotrophic and autotrophic respiration within the soil) [18].Rapid changes in water availability after a PP could increase the entropy of water-limited ecosystems [25], especially after those initial PPs following a long dry period resulting in hot-moments with disproportionately high Fs [26].Such PPs could rapidly change soil CO 2 diffusion rates, increase photosynthetic substrate supplies for Fs, and enhance soil microbial and plant metabolism that ultimately increase Fs [8,18,27,28].The non-stationary nature of these events and the diversity of underlying mechanisms complicate measurement campaigns and modeling efforts [29].
There are several challenges for studying Fs hot-moments triggered by PPs.First, hot-moments are rare events and are constrained by unique biophysical conditions such as antecedent soil moisture [30], PP magnitude [31], and the available substrate supply for microbial respiration [32,33].Second, hot-moments are sporadic events resulting from a sharp input of a forcing variable (i.e., water input in water-limited ecosystems); therefore, continuous measurements are needed to accurately capture their succinct patterns and magnitudes [15,29,34].Third, semi-empirical functions based on temperature and moisture responses are commonly used to model Fs [35,36].These commonly used functions usually fail to represent Fs hot-moments because fluxes could increase >500% in just a few hours versus pre-pulse conditions [29,34].Considering the attention that water-limited ecosystems have gained for the global carbon cycle [1,2], there is a need for more information on the magnitudes and potential mechanisms of Fs in these globally distributed ecosystems.
The overarching goal of this study is to characterize the influence of PPs that could trigger Fs hot-moments in a water-limited grassland.We analyzed >2-years of continuous measurements of Fs under vegetation (Fs veg ) and at bare soil (Fs bare ).We postulate three interrelated hypotheses: (a) Not all PPs will trigger hot-moments, as hot-moments may depend on the intensity of the PP, pre-pulse soil moisture conditions (e.g., after a long dry period), and canopy metabolism (e.g., during growing season); (b) Soil moisture variability will have higher temporal correlation with Fs than soil temperature in this water-limited ecosystem; and (c) Hot-moments -as non-stationary events-are difficult to represent for semi-empirical models, but machine learning techniques could improve their representation as these methods can account for non-linear relationships.Here, we applied time series analysis to describe the temporal variability of Fs veg and Fs bare , and use three modeling approaches ranging in complexity from linear models to machine learning (i.e., support vector machine) to test their potential for representing temporal trends and hot-moments.The ultimate aim of this study is to enhance the discussion of measurement efforts and modeling approaches for Fs hot-moments in these globally distributed ecosystems.

Study Site
The Balsa Blanca site is located at 208 m.a.s.l and 6.3 km from the coast, in the Cabo de Gata Natural Park (Almería, Spain; N36 • 56 26.0 , W2 • 01 58.8 ).Following the Köppen classification, the site has a desert climate (Bwh) characterized by the Thermo-Mediterranean bioclimatic zone.The mean air temperature is 18 • C with mean annual precipitation of 200 mm year −1 .Bare soil, gravel and rock cover about 49% of the landscape.Vegetation is sparse with 60% cover, and is dominated by the perennial grass Macrochloa tenacissima (L.) Kunth with a mean height of 0.5 m.The rest of the vegetation is composed of Chamaerops humilis, Rhamnus lycoides, and Pistacia lentiscus.The soils are Mollic Leptosols with a sandy loam texture (61.0%sand, 22.8% silt, 19.5% clay), 1.9% organic carbon, 0.16% total nitrogen, 12.2 ratio C:N, 7.9 pH, 1.5% equivalent carbonates, and with a soil bulk density of 1.25 g cm −3 in the upper 30 cm of the soil.Analyzed data were collected between June 2011 and November 2013 (i.e., >2 years of measurements) at the study site.

Eddy Covariance
Net ecosystem exchange (NEE) was measured using an eddy covariance tower, and we use this information to describe the general metabolic activity (i.e., ecosystem acting as a net carbon source or sink) of the ecosystem during the study period.Previous studies have described in detail the instrumentation, data processing, quality assurance/quality control, and data gapfilling for the study site [37].Briefly, fluxes of CO 2 were estimated from fast-response instruments (10 Hz measurements) mounted atop a 3 m tower using an open-path infrared gas analyzer (IRGA, LI-7500, Licor; Lincoln, NE, USA).Winds and sonic temperature were measured by a three-axis sonic anemometer (CSAT-3, Campbell Scientific, Logan, UT, USA, hereafter CSI).Data QA/QC and post-processing was performed following standard procedures [38] and described previously for the study site [37].Storage fluxes were not calculated because the study site has short canopy with well mixing and it is assumed that the storage flux to be zero for a 24 h period.Gaps due to environmental conditions, instrument malfunction and nighttime low turbulence were filled using the marginal distribution sampling technique [39], replacing missing values using a time window of several adjacent days.Positive values of fluxes denote net release to the atmosphere, while negative values indicate net uptake by the ecosystem.

Soil CO 2 Efflux (Fs)
To estimate Fs we measured soil CO 2 concentrations within the soil and applied the gradient method [40].Because the vegetation cover is sparse we instrumented an area under vegetation (dominated by Macrochloa tenacissima (L.) Kunth) to estimate Fs veg , and another area with bare soil (situated at 50 cm away from plants) to estimate Fs bare .We consider that 50 cm away from Macrochloa tenacissima (L.) Kunth is enough to avoid most of the roots as a previous study demonstrated that most roots of Macrochloa tenacissima (L.) Kunth are localized right underneath its canopy [41].This approach has been followed at the study site by multiple studies to understand spatial variability of FS [42,43].All solid-state CO 2 sensors (GMM-222, 0-10,000 ppm, Vaisala, Inc., Vantaa, Finland) were installed with a soil temperature probe (107, CSI), and a water content reflectometer (CS616, CSI) at 5 cm depth (i.e., one sensor per depth at each location).All measurements were made every 30 s and stored as 5 min averages.
Soil CO 2 efflux (Fs) was calculated assuming that all transport is due to diffusion [40] as: where Fs is the upward gas flux (µmol CO 2 m −2 s −1 ), D S the soil CO 2 diffusion coefficient (m 2 s −1 ), ρ a the mean air molar density (µmol m −3 ), dχ c dz is the vertical CO 2 molar fraction gradient (ppm m −1 ).The CO 2 gradient was calculated using the difference between the mean air CO 2 molar fraction (obtained from the eddy covariance tower) and the value of each soil CO 2 sensor.The CO 2 molar fraction was corrected for variations in temperature and pressure.The mean air CO 2 molar fraction was obtained from days when the eddy covariance tower had been recently calibrated and used as a constant for all the studied period, due to the impossibility of using the data continuously due to soiling of the infrared gas analyzer (IRGA) lens.The calibration of this sensor was done monthly using a N 2 standard for zero (purity of 99.999%) and known CO 2 standards for span.We assumed a constant value for atmospheric CO 2 molar fraction, neglecting its small fluctuations that cause negligible errors Soil Syst.2018, 2, 0047 4 of 18 on the final fluxes [44].Due to the large difference in CO 2 molar fraction between the soil (when in some cases could be >1000 ppm) and the atmosphere we estimate that the systematic errors could be between 2% and 4% on the final fluxes.The soil CO 2 diffusion coefficient (D S ) was obtained as: where D a is the diffusion coefficient of the CO 2 in free air, calculated according to Jones [45], and ξ is the tortuosity.We calculated ξ following Moldrup et al. [46] as: where φ is soil porosity and θ is soil volumetric water content [47].In this manuscript positive values of Fs denote a CO 2 release to the atmosphere in consistency with NEE measurements.

Empirical Modeling of Soil CO 2 Efflux (Fs)
We tested three hierarchical approaches to model hourly and daily means of Fs including a simple conditional linear approach, support vector machine (SVM) as a machine learning approach, and an 8-day moving window.First, we modeled Fs veg and Fs bare using soil moisture and soil temperature as forcing factors using a conditional approach (referred as Model 1 throughout the text).It is known that water availability is the main controlling factor in water-limited ecosystems, but the temperature dependence of Fs is only relevant when water is not a limiting factor [26,35,48].Thus, we used the following approach: If soil temperature > mean annual soil temperature (when the soil is likely to be hot and dry) If soil temperature < mean annual soil temperature (when it is likely to be cool and moist) log(Fs) = a+ b(Ts) + b(SM) + c(SM) 2 (5) where log(Fs) is either natural logarithm of Fs veg or Fs bare , SM is soil moisture in the bare soil or under vegetation, and Ts is soil temperature in the bare soil or under vegetation.This simple approach is likely to represent the temporal trends but unlikely to represent hot-moments of Fs.Second, we used a support vector machine (SVM) approach using soil temperature and soil moisture as predictor variables for Fs veg or Fs bare .This machine learning technique are supervised learning models/algorithms that analyze data for classification and regression analysis.To explore the nonlinear relations of the response of Fs we used a kernel function with SVM.The kernel function was set to be Gaussian and the kernel scale was set to 0.35.We propose that SVM is a flexible approach for exploring nonlinear relations for classification of data and testing its applicability for representing Fs pulses.Previous studies have described in detail the theory behind SVM [49,50] and the application of Gaussian kernels for SVM [51,52].We cross-validate the classifier using 10-fold cross-validation.
Third, we modeled Fs veg and Fs bare using an 8-day moving window (referred as Model 2 throughout the text).Thus, it allows for parameters to shift through time so that projections match the data as close as possible.For this approach, we only applied Equation (5) assuming that temperature and moisture are important within this 8-day moving window.In this case, different parameters for Equation (5) were fitted for each 8-day sliding window calculation.This moving window approach is conceptually analogous to the calculation of ecosystem respiration based on a moving window of nighttime NEE [39].The selection of an 8-day moving window was based on the fact that the effect of PPs on Fs usually last about 8 days at the study site [27] and this highly parameterized approach could better represent hot-moments of Fs.

Time Series Analyses
First, we used wavelet analysis to explore the spectral properties of the times series of Fs (Fs veg and Fs bare ) and data-model agreement (or disagreement) for Fs veg and Fs bare derived from the three empirical approaches described above.Wavelet analysis is a time series technique that has been widely applied in the geosciences [53], and recently utilized to analyze the temporal variability of ecosystem-scale fluxes [54-56], soil CO 2 effluxes [15,57], and to identify data-model agreement [57][58][59].This technique is used to quantify the spectral characteristics of time series that may be nonstationary and heteroscedastic.Specifically, we used the continuous wavelet transform because of its ability to produce a smooth picture in the frequency domain of a time series (e.g., soil CO 2 efflux, data-model residuals) and its suitability for visual interpretation.The ability to discern small intervals of scales (i.e., spectral resolution) depends on the choice of the mother wavelet function.For this, we used the widely used Morlet wavelet, a complex nonorthogonal wavelet with a good time and scale resolution that has been widely used for geophysical applications [53,60], and biometeorological measurements [57][58][59].We first analyzed time-series of measurements using a 1-h time step and then model residuals using a 1-day time step.
Second, to test the potential influence of soil temperature or soil moisture on Fs, we explored their partial temporal correlation with Fs veg or Fs bare .In other words, we explored the temporal correlation of Fs with soil temperature taking into account the influence of soil moisture and vice versa.We applied partial wavelet coherence analysis (PWC), as it can be interpreted as a technique similar to partial correlation that can identify significant temporal correlations between two different time series after eliminating the influence of a third one.Previous reports have described the PWC in detail for climate studies [61].The statistical significance (5% significance level) of common power between any two time series was assessed for PWC using 1000 Monte Carlo simulations of white noise time series [60].The time-series used for PWC were analyzed using a 1-h time step.

Description of Temporal Patterns
Daily soil temperature ranged from a minimum of 5 • C and a maximum of 38 • C; with an annual mean of 22 • C under vegetation and 23 • C for bare soil (Figure 1A).Between September 2011 and August 2012 (first hydrologic year) the study site received 216 mm of precipitation, and between September 2012 and August of 2013 (second hydrologic year) it received 221 mm (Figure 1B).There were no substantial differences in total precipitation between the hydrological years, but the distribution of the PPs resulted in different patterns of soil moisture along the years (Figure 1B).
Water pulses following the dry season (between September and October) substantially raised Fs veg and Fs bare (Figure 1C).The overall mean for Fs veg and Fs bare was 1.5 ± 1.4 µmol CO 2 m −2 s −1 but hot-moments were present with instantaneous values >20 µmol CO 2 m −2 s −1 (Figure 2A,B).We analyzed the cumulative sum of Fs based on hydrologic years (see above).Annual Fs veg was 480 and 819 gC m −2 , while annual Fs bare was 431 and 871 gC m −2 for the first and second hydrologic years, respectively.Positive NEE values were common during the dry season when the ecosystem acted as a net CO 2 source to the atmosphere (Figure 1D).
We explored the spectral characteristics of the time series of Fs veg and Fs bare (Figure 3).Results using wavelet analysis demonstrate that the 1-day period for Fs veg is more constant than for Fs bare (Figure 3A,B).We observed three distinct PPs that substantially influenced the spectral signature of Fs.First, the large precipitation pulses at the beginning of each growing season (i.e., September) resulted in discrete hot-moments influencing the periodicity for Fs veg (Figure 3A) and Fs bare (Figure 3B) between a time-period of 1-to 16-days (pulses Case I and II).Arguably, another "hot-moment" is observed during October of 2012, but for the purpose of this study (and to simplify the discussion) we consider it similar to a Case II and therefore is not further analyzed.Second, the large precipitation pulse during the growing season of year 2013 also influence the periodicity for Fs veg and Fs bare between 1-to 8-days time-periods (pulse Case III).

Influence of Soil Temperature and Soil Moisture on Fs
Using PWC analysis we found a lack of consistent temporal coherence between Fsveg or Fsbare with soil temperature when considering the effect of soil moisture (Figure 4A,B).In contrast, we found significant temporal coherence between Fsveg or Fsbare with soil moisture when considering the effect of soil temperature (Figure 4C,D).The temporal coherence was clearly influenced by PPs for Case I, II, and III with significant temporal correlations between 1-and 16-day periods.The seasonal influence of soil moisture on Fs is represented by the temporal correlation at scales > 128-days (Figure 4C,D).

Influence of Soil Temperature and Soil Moisture on Fs
Using PWC analysis we found a lack of consistent temporal coherence between Fs veg or Fs bare with soil temperature when considering the effect of soil moisture (Figure 4A,B).In contrast, we found significant temporal coherence between Fs veg or Fs bare with soil moisture when considering the effect of soil temperature (Figure 4C,D).The temporal coherence was clearly influenced by PPs for Case I, II, and III with significant temporal correlations between 1-and 16-day periods.The seasonal influence of soil moisture on Fs is represented by the temporal correlation at scales > 128-days (Figure 4C,D).

Model 1
Using the conditional approach for Equations ( 4) and ( 5) with daily Fs veg values resulted in an explained variance of 47%, with a RMSE of 0.49 µmol CO 2 m −2 s −1 .The continuous wavelet transform shows that errors in the residuals were focused on hot-moments during PPs with periodicities ranging between 2-and 16-days (Figure 5A,B).These underestimated hot-moments correspond to the PPs previously identified as Case I to III.This model has a larger probability to over represent Fs veg values of ~1.5 µmol CO 2 m −2 s −1 and under represent Fs veg values > 3.0 µmol CO 2 m −2 s −1 when compared to measurements (Figure 6A).Considering the full 2-years of measurements, the model underestimates the cumulative Fs veg flux by 12% mainly by over representing Fs veg values for year 2011 and underestimating Fs veg values of years 2012 and 2013 (Figure 6B).
Using PWC analysis we found a lack of consistent temporal coherence between Fsveg or Fsbare with soil temperature when considering the effect of soil moisture (Figure 4A,B).In contrast, we found significant temporal coherence between Fsveg or Fsbare with soil moisture when considering the effect of soil temperature (Figure 4C,D).The temporal coherence was clearly influenced by PPs for Case I, II, and III with significant temporal correlations between 1-and 16-day periods.The seasonal influence of soil moisture on Fs is represented by the temporal correlation at scales > 128-days (Figure 4C,D).Using this approach for Fs bare resulted in an explained variance of 35% and a RMSE of 0.67 µmol CO 2 m −2 s −1 .The continuous wavelet transform also shows that errors for Fs bare residuals were focused on hot-moments that correspond to the PPs previously identified as Case I to III (Figure 5G,H).This model also has a larger probability of overestimating Fs veg values of ~1.5 µmol CO 2 m −2 s −1 and underestimating Fs veg values > 2.5 µmol CO 2 m −2 s −1 when compared to measurements (Figure 6C).Considering the full 2-years of measurements, the model underestimates the cumulative Fs bare flux by 15% mainly by overestimating Fs bare values for year 2011 and underestimating Fs bare values of years 2012 and 2013 (Figure 6D).

Support Vector Machine (SVM)
Using SVM, we were able to explain 71% of the variance of Fs veg with a RMSE of 0.7, and 51% of the variance of Fs bare with a RMSE of 0.98 µmol CO 2 m −2 s −1 .The continuous wavelet transform shows that errors in the residuals were also focused on the previously identified Cases I to III for Fsveg (Figure 5C,D) and Fs bare (Figure 5I,J).This model also has a larger probability of overestimating Fs values of ~1.5 µmol CO2 m −2 s −1 and underestimating Fs veg values > 2.5 µmol CO2 m −2 s −1 when compared to measurements (Figure 6A,C).Considering the full 2-years of measurements, the model underestimates the cumulative Fs veg flux by 7% and underestimates Fs bare by 1% as a result of error cancelation.This approach overrepresents Fs values for year 2011 and underestimates Fs values of years 2012 and 2013 (Figure 6B,D).

Model 2
Using the 8-day moving window approach with Equation (5), we were able to explain 84% of the variance of Fs veg with a RMSE of 0.52, and 48% of the variance of Fs bare with a RMSE of 0.73 µmol CO 2 m −2 s −1 .The continuous wavelet transform shows that errors in the residuals were also focused on hot-moments during PPs with periodicities ranging between 2-and 16-days (Figure 5E,F).These underestimated hot-moments also correspond to the PPs previously identified as Case I to III.This model has a similar probability density function for Fs veg values when compared to measurements (Figure 5A), but slightly overestimates Fs bare values of ~1.5 µmol CO 2 m −2 s −1 and underestimates Fs bare values > 3.0 µmol CO 2 m −2 s −1 when compared to measurements (Figure 6A,C).Considering the full 2-years of measurements, the model underestimates the cumulative Fs veg flux by 2% and underestimates Fs bare by 3% as this approach closely followed the measurements (Figure 6B,D).

Responses to Different Precipitation Pulses
We selected three responses (i.e., Cases) to PPs based on the spectral properties of the time series of measurements and model residuals.These response Cases represent rapid discrete changes in the amplitude of soil moisture from a baseline (i.e., pre-pulse conditions), have substantial influence on the spectral properties of the time series (Figure 3), and are difficult to represent by the proposed modeling approaches (Figure 5).We identified that a Case I response was evident during September 2011 (i.e., days of the year 240 to 265), a Case II during September 2012 (i.e., days of the year 240 to 265), and a Case III during March 2013 (i.e., days of the year 110 to 135).
A Case I response followed a large PP (i.e., >20 mm) after the long dry season that sharply increased soil moisture (Figure 7A).This resulted in an increase of Fs veg to 11.7 µmol CO 2 m −2 s −1 representing a change of 5400% and of Fs bare to 8 µmol CO 2 m −2 s −1 (increase of 3885%) from pre-pulse conditions (Figure 7D).A Case II response followed a small PP (i.e., <5 mm) after the dry season that slightly increased soil moisture (Figure 7E).Despite this modest response in soil moisture Fs veg increased to 18 µmol CO 2 m −2 s −1 and Fs bare to 34.4 µmol CO 2 m −2 s −1 , representing a change of 9000% and 17,000% from pre-pulse conditions, respectively (Figure 7E).A Case III response was a large PP (i.e., >20 mm) during the moist conditions of the growing season that moderately increased soil moisture (Figure 7C) and increased Fs veg to 8.7 µmol CO 2 m −2 s −1 and Fs bare to 11.7 µmol CO 2 m −2 s −1 representing an increase of 954% and 1099% from pre-pulse conditions, respectively (Figure 7F).

Discussion
Our results support the paradigm that the distribution and intensity of PPs influence annual Fs in water-limited ecosystems.Relatively similar total precipitation but lower soil moisture variability resulted in 50% higher Fs emissions during the second hydrologic year.Lower soil moisture variability and prolonged rains (until June 2013) likely reduced water stress during the second hydrologic year.Higher ecosystem metabolic activity during the second year was observed by the consistent periodicity at the 1-day period for Fsveg and Fsbare (Figure 3).These results support observations that precipitation patterns which maintain higher soil moisture conditions for a longer time resulted in higher seasonal ecosystem metabolic rates [25,62].
Discrete PPs generate Fs responses that resulted in hot-moments of Fsveg and Fsbare (response Case I to III).The wavelet analyses demonstrate that these pulse responses had distinct spectral signatures We tested how the three modeling approaches were able to represent the different responses using a 1-h time-step.Overall, all models were able to represent Case I and III responses (Figure S1), but SVM was the model that better represented these responses (Figure 7H, I).None of the models were able to represent a Case II response.

Discussion
Our results support the paradigm that the distribution and intensity of PPs influence annual Fs in water-limited ecosystems.Relatively similar total precipitation but lower soil moisture variability resulted in 50% higher Fs emissions during the second hydrologic year.Lower soil moisture variability and prolonged rains (until June 2013) likely reduced water stress during the second hydrologic year.Higher ecosystem metabolic activity during the second year was observed by the consistent periodicity at the 1-day period for Fs veg and Fs bare (Figure 3).These results support observations that precipitation patterns which maintain higher soil moisture conditions for a longer time resulted in higher seasonal ecosystem metabolic rates [25,62].
Discrete PPs generate Fs responses that resulted in hot-moments of Fs veg and Fs bare (response Case I to III).The wavelet analyses demonstrate that these pulse responses had distinct spectral signatures localized within specific PPs (Figure 3).Not every PP resulted in an event with a distinct spectral signature in the time series of Fs demonstrate the uniqueness of the selected Cases (Figure 3).In other words, despite the fact that there were several PPs throughout the length of the study, we identify three distinct Cases that we interpret as hot-moments of Fs.We recognize that there is a hierarchy of PPs, but we bring attention to the use of time-series analysis and the interpretation of PPs and hot-moments by analyzing information of automated Fs measurements in the frequency-domain.Arguably, previous studies have likely underestimated hot-moments of Fs due to lack of continuous measurements as recent studies are demonstrating the importance of these transient but intense events [29,34,63].

Can We Model Daily Fs Using Temperature and Soil Moisture Information?
Our results show that the temporal influence of soil moisture on Fs overrides the influence of soil temperature on Fs in this water-limited ecosystem (Figure 4).This supports the fact that temperature is only relevant in these ecosystems when soil moisture is available for metabolic processes [26,35,48].The temporal influence of soil temperature in this water-limited ecosystem seems to be concentrated at the 1-day period during discrete days when water is available.In contrast, the influence of soil moisture has larger implications for the temporal variability of Fs at scales ranging from 2-to 16-days during the different responses (i.e., Case I to III).
Our results demonstrate the challenge that empirical modeling approaches have to represent hot-moments of Fs in water-limited ecosystems, where a Case II response was consistently challenging to represent by all approaches (Figure 5).A Case II response substantially increased the magnitude of Fs but the magnitude of the response was not proportional to the increase in soil moisture.Although this response is directly linked to an increase in soil moisture, the underlying biotic and abiotic mechanisms go beyond water availability as previously discussed for rewetting events [29,34].
The conditional approach (i.e., Model 1) was able to represent between 47% and 35% of the variability in daily Fs veg and Fs bare , respectively.This simple conditional approach can be easily applied across study sites and is interpretable, as it provides insights about temperature and soil moisture sensitivity of Fs (by comparing constants in the model).This approach is widely applicable at the daily-scale and parameters can be compared across sites and across site-years but has the largest uncertainty of all approaches.
The wealth of information from continuous measurements allows for the application of more complex models that can be highly parameterized such as the 8-day moving window approach (i.e., Model 2).This moving window approach follows the fact that water availability and precipitation pulses drive most of the metabolic pulses in this water-limited ecosystem [27].The use of an 8-day moving window was confirmed by the PWC showing that water pulses have a consistent temporal influence between 8-to 16-days regardless of their magnitude or antecedent moisture conditions.This approach significantly increased the proportion of explained variance to 84% for daily Fs veg and 48% to daily Fs bare .This approach could be applied as a gapfilling technique for time series of Fs [64], but has little interpretability and applicability beyond the parameterized time series because parameters are highly variable in order to match measurements for each moving window (Figure 6).Future applications could be done varying the size of the window to optimize for site-specific conditions and to improve gapfilling estimates.
Machine learning approaches are versatile and flexible enough to discover complicated nonlinear relationships and we argue that SVM has the potential for representing the non-stationary dynamics of Fs.This approach was able to increase the representation of daily Fs bare when compared with the approach of Model 2. Likely the use of more variables has the potential for improvement of predictions of hot-moments of Fs that are controlled by factors beyond soil moisture, but caution must be taken to avoid autocorrelation among variables and model over fitting [65].The opportunity for machine learning approaches to Fs is arguably starting with a few examples such as the use of random forests [66]; thus, it is expected that the use of these techniques will become more common in the near future.

How Discrete PPs Influence Hot-Moments of Fs?
We identified three Cases of Fs to PPs where the Fs rate drastically changed (up to 17,000%) to affect the spectral properties of the time series and data-model agreement.Here, we describe the generalities of these responses, but we recognize that our results are based on the information within the available time series and longer records could show consistency or a larger diversity of responses.
Case I represents the ecosystem response after a large PP (>20 mm) following a long drought period, and is characteristic of the beginning of the rainy season in Mediterranean ecosystems [26,27].Case I is characterized by rapid water infiltration into the soil profile, which substantially increases soil water content in areas of bare soil and under vegetation (Figure 8A).Furthermore, the vegetation could provide a preferential pathway for water infiltration, as demonstrated by an increase of soil moisture with respect to the bare soil.In Case I we observed larger Fs veg than Fs bare with two possible explanations.First, we postulate that the autotrophic component of Fs (i.e., Fs veg minus Fs bare ) increased as plants rapidly start to use resources to activate the photosynthetic mechanism, therefore increasing their catabolism and autotrophic respiration [27].Second, it is likely that the heterotrophic component of Fs under vegetation has also increased because of higher substrate availability within the rhizosphere that is rapidly dissolved and available for microbes [28].Notably, the modeling approach by SVM was able to represent this response as Fs was highly correlated with soil moisture regardless of the underlying mechanisms of the heterotrophic and autotrophic components of Fs (Figure 7G).
Case II represents the ecosystem response after a small PP (<5 mm) following a long drought period.This distinction is important because in water-limited ecosystems, the beginning of the rainy season does not always start with a large (i.e., >20 mm) precipitation event [13,67,68].Case II is characterized by slow water infiltration that only permeates into soil surface layers because it is likely that vegetation intercepts most water from these small PPs and creates a shadow where soil moisture is lower (under vegetation) than at the bare soil (Figures 7B and 8B).Despite the small PP we observed a hot-moment for Fs veg (18 µmol CO 2 m 2 s −1 ) and Fs bare (>34 µmol CO 2 m 2 s −1 ) that increased Fs by 9000 and 17,000% from pre-pulse conditions, respectively.These responses are the highest reported for Fs for rewetting events [29] as a result of very low pre-pulse Fs rates (due to dry conditions during the non-growing season) and the sharp increase following the PP.Previous studies have demonstrated that Fs is spatially heterogeneous within a water-limited ecosystem [26] and that small PPs do not activate the metabolism at the plant-level [27].We postulate that increases in Fs veg and Fs bare are likely driven by an increase in the heterotrophic component of Fs because there is not enough available water to trigger plant catabolic metabolism to stimulate autotrophic respiration (Figure 8B).The lack of correlation between the size of the PP and the response of Fs make modeling of these responses Soil Syst.2018, 2, 0047 13 of 18 very challenging.None of the proposed approaches were able to represent the response as underlying variables (e.g., available substrate supply) may have influenced the sharp uncorrelated response to soil moisture.The implication of this modeling limitation is that low-to-medium Fs fluxes are usually overrepresented by models and high Fs fluxes are underestimated (Figure 6).Furthermore, these large responses could contribute up to 40% of total net CO 2 emissions during dry seasons [27], so accurate measurements and understanding of these hot-moments is needed to better understand local carbon dynamics.Case III represents the ecosystem response to a large precipitation event (~20 mm) during the cool and moist growing season (Figures 7C and 8C).We postulate that although the plant canopy can provide a preferential pathway for water infiltration into the soil (as seen in Case I), this season is characterized by less soil moisture limitation and a homogeneous distribution of soil water content across the soil profile (Figure 8C).Case III is characterized by a large change in soil water content with increases in Fs of about 1000%, but this sharp increase only represents <20% of the sharp response observed for Case I and II.Notably, the modeling approach by SVM was able to represent the response of Fsveg as it was highly correlated with soil moisture, but in a lesser degree the response of Fsbare which had a disproportionally low response to increased moisture in the soil profile (Figure 7I).We postulate that Case III is likely a result of an increase in heterotrophic and autotrophic activity as plant metabolism is active during the growing season [27].
We borrow the concept of the "bucket model" to describe how discrete PPs influence hot- Case III represents the ecosystem response to a large precipitation event (~20 mm) during the cool and moist growing season (Figures 7C and 8C).We postulate that although the plant canopy can provide a preferential pathway for water infiltration into the soil (as seen in Case I), this season is characterized by less soil moisture limitation and a homogeneous distribution of soil water content across the soil profile (Figure 8C).Case III is characterized by a large change in soil water content with increases in Fs of about 1000%, but this sharp increase only represents <20% of the sharp response observed for Case I and II.Notably, the modeling approach by SVM was able to represent the response of Fs veg as it was highly correlated with soil moisture, but in a lesser degree the response of Fs bare which had a disproportionally low response to increased moisture in the soil profile (Figure 7I).We postulate that Case III is likely a result of an increase in heterotrophic and autotrophic activity as plant metabolism is active during the growing season [27].
We borrow the concept of the "bucket model" to describe how discrete PPs influence hot-moments of Fs.The "bucket model" is a conceptual idea to describe the response of terrestrial ecosystems under different precipitation patterns [69].In this model, the "bucket" represents the uppermost soil layers with maximum root density and is characterized by upper and lower water stress thresholds (i.e., a stress gradient) [69].We propose that an analogous conceptual approach can be used to explain the hot-moments of Fs in our study.Therefore, "the soil water bucket" represents the uppermost soil layers within maximum root density and is characterized by upper and lower water stress thresholds that represent a gradient of ecosystem stress (Figure 8D-F).In Case I the ecosystem is stressed as a consequence of the long dry season, but the PP is able to increase soil moisture that moves the ecosystem to a state of less stress for a prolonged period of time until the following precipitation event (Figure 8D); this response is highly correlated with the temporal patterns of soil moisture.In Case II the ecosystem is also stressed as a consequence of the long dry season, but soil moisture is not substantially increased, resulting in a short period of reduced stress followed by a sharp return to a previous stressed state (Figure 8E).Finally, in Case III the ecosystem is less stressed and subsequent changes in soil moisture create variability under a stress gradient (Figure 8F).This conceptual idea aims to highlight that variability in PPs defines water-limited ecosystems, pushes these ecosystems towards less-stressed conditions, and could trigger a diverse response of hot-moments of Fs.
Climate models indicate a future with altered precipitation patterns where extreme PPs could be followed by long dry periods in water-limited ecosystems.Thus, it is critical to understand how diverse PPs influence ecosystem processes under different metabolic states.Automated measurements of Fs provide the opportunity to capture high-temporal resolution of ecosystem responses, providing information on disproportionately high (i.e., hot-moments) CO 2 fluxes following precipitation events [34].Although CO 2 fluxes in water-limited ecosystems are low compared to mesic ecosystems, their sensitivity to changes in precipitation pattern is high and consequently influences their annual net fluxes.Our results support the application of a machine learning approach (i.e., SVM) based on information of soil moisture and temperature to represent Fs, but we recognize that machine learning is parameterized with available data and consequently is not process based.Hot-moments of Fs as a result of a small precipitation event (<5 mm; Case II pulse) following a long drought period appear to be the most challenging events to represent and support the need of continuous measurements to capture the effects of this discrete but sharp response.We demonstrate that soil moisture has high temporal correlation with Fs and overshadows the influence of soil temperature in this water-limited ecosystem.Finally, because the variability of carbon dynamics in water-limited ecosystems influences the global carbon cycle, it is essential to quantify the responses of non-stationary ecosystem CO 2 fluxes to transient (and potentially extreme) precipitation events.

19 Figure 1 .
Figure 1.Time series of daily mean of soil temperature (A), daily sum of precipitation, daily mean of soil moisture (B), daily mean of soil CO2 efflux (Fs) nearby vegetation (Fsveg) and in bare soil (Fsbare; C), and daily mean of net ecosystem exchange (NEE; D).Dashed lines in panels A-C represent variables measured under vegetation and solid lines represent variables measured in bare soil.Gray areas in panels B and C represent precipitation pulses selected as response cases (i.e., Case I, Case II, Case III; see results section).Solid circles in panel B represent precipitation.Letters in the x-axis represents months of the year from July 2011 to November 2013.

Figure 2 .
Figure 2. Time series of half-hour values of soil CO2 efflux under vegetation (Fsveg; A) and in bare soil (Fsbare; B).Letters in the x-axis represents months of the year from July 2011 to November 2013; y-axis represent each hour of the day.

Figure 1 . 19 Figure 1 .
Figure 1.Time series of daily mean of soil temperature (A), daily sum of precipitation, daily mean of soil moisture (B), daily mean of soil CO 2 efflux (Fs) nearby vegetation (Fs veg ) and in bare soil (Fs bare ; C), and daily mean of net ecosystem exchange (NEE; D).Dashed lines in panels A-C represent variables measured under vegetation and solid lines represent variables measured in bare soil.Gray areas in panels B and C represent precipitation pulses selected as response cases (i.e., Case I, Case II, Case III; see results section).Solid circles in panel B represent precipitation.Letters in the x-axis represents months of the year from July 2011 to November 2013.

Figure 2 .
Figure 2. Time series of half-hour values of soil CO2 efflux under vegetation (Fsveg; A) and in bare soil (Fsbare; B).Letters in the x-axis represents months of the year from July 2011 to November 2013; y-axis represent each hour of the day.

Figure 2 .
Figure 2. Time series of half-hour values of soil CO 2 efflux under vegetation (Fs veg ; A) and in bare soil (Fs bare ; B).Letters in the x-axis represents months of the year from July 2011 to November 2013; y-axis represent each hour of the day.

Figure 3 .
Figure 3. Continuous wavelet power spectra using the continuous wavelet transform of the time series of soil respiration under vegetation (Fsveg; A) and in bare soil (Fsbare; B).The color codes for power values are from yellow (low values) to dark red (high values).Black contour lines represent the 5% significance level and thin black line indicates the cone of influence that delimits the region not influenced by edge effects.Letters in the x-axis represents months of the year from July 2011 to November 2013.Arrows represent "hot-moments" of Fs identified for this study.

Figure 3 .
Figure 3. Continuous wavelet power spectra using the continuous wavelet transform of the time series of soil respiration under vegetation (Fs veg ; A) and in bare soil (Fs bare ; B).The color codes for power values are from yellow (low values) to dark red (high values).Black contour lines represent the 5% significance level and thin black line indicates the cone of influence that delimits the region not influenced by edge effects.Letters in the x-axis represents months of the year from July 2011 to November 2013.Arrows represent "hot-moments" of Fs identified for this study.

Figure 4 .
Figure 4. Partial wavelet coherence analysis (PWC) to test the influence of soil temperature (controlling for soil moisture) on soil CO 2 efflux in areas under vegetation (Fs veg ; A) or in bare soil (F sbare ; B).PWC analysis to test the influence of soil moisture (controlling for soil temperature) on CO 2 efflux in areas under vegetation (Fs veg ; C) or in bare soil (Fs bare ; D).The color codes for temporal correlation are from blue (low values) to yellow (high values).Yellow areas within black contour lines represent the 5% significance level for power values (i.e., high temporal correlation); the thin black line indicates the cone of influence that delimits the region not influenced by edge effects.Gray areas in panels C and D represent precipitation pulses analyzed to characterize the response cases (i.e., Case I, Case II, Case III).Letters in the x-axis represents months of the year from July 2011 to November 2013.Soil Syst.2018, 2, x FOR PEER REVIEW 9 of 19

Figure 5 .
Figure 5. Analysis of model residuals of soil CO2 efflux (Fs) under vegetation (Fsveg; A-F) or in bare soil (Fsbare; G-L).Residuals using Equations (2) and (3) (Model 1), Support Vector Machine (SVM) or an 8-day window approach (Model 2) for Fsveg (A,C,E) and Fsbare (G,I,K).Continuous wavelet power spectra using the continuous wavelet transform for residuals of each model for Fsveg (B,D,F) and Fsbare (H,J,L).The color codes for power values are from yellow (high values) to red (low values).Yellow areas within black contour lines represent the 5% significance level for power values (i.e., high value of residuals); the thin black line indicates the cone of influence that delimits the region not influenced by edge effects.

Figure 5 .
Figure 5. Analysis of model residuals of soil CO 2 efflux (Fs) under vegetation (Fs veg ; A-F) or in bare soil (Fs bare ; G-L).Residuals using Equations (2) and (3) (Model 1), Support Vector Machine (SVM) or an 8-day window approach (Model 2) for Fs veg (A,C,E) and Fs bare (G,I,K).Continuous wavelet power spectra using the continuous wavelet transform for residuals of each model for Fs veg (B,D,F) and Fs bare (H,J,L).The color codes for power values are from yellow (high values) to red (low values).Yellow areas within black contour lines represent the 5% significance level for power values (i.e., high value of residuals); the thin black line indicates the cone of influence that delimits the region not influenced by edge effects.

Figure 6 .Figure 6 .
Figure 6.Probability density functions of soil CO2 efflux (Fs) values for measurements and each model under vegetation (Fsveg; A) and bare soil (Fsbare; C).Cumulative sums of soil CO2 fluxes for measurements and each model under vegetation (Fsveg; B) and bare soil (Fsbare; D).Letters in the x-Figure 6. Probability density functions of soil CO 2 efflux (Fs) values for measurements and each model under vegetation (Fs veg ; A) and bare soil (Fs bare ; C).Cumulative sums of soil CO 2 fluxes for measurements and each model under vegetation (Fs veg ; B) and bare soil (Fs bare ; D).Letters in the x-axis represents months of the year from July 2011 to November 2013.For a description of the models see methods section.SVM = support vector machine.

Figure 7 .
Figure 7. Pulse responses of soil moisture (A-C), Fs under vegetation (Fsveg) or in bare soil (Fsbare; D-F), and model output using support vector machine (SVM) for Fsveg or Fsbare (G-I).Responses are defined as Case I-III, and are identified in Figures 1 and 3. Solid lines represent variables measured under vegetation and dashed lines represent variables measured in bare soil.Fs means soil CO2 efflux.The numbers in the x-axis represent the day of the year (DOY) for years 2011 to 2013 for each panel.

Figure 7 .
Figure 7. Pulse responses of soil moisture (A-C), Fs under vegetation (Fs veg ) or in bare soil (Fs bare ; D-F), and model output using support vector machine (SVM) for Fs veg or Fs bare (G-I).Responses are defined as Case I-III, and are identified in Figures 1 and 3. Solid lines represent variables measured under vegetation and dashed lines represent variables measured in bare soil.Fs means soil CO 2 efflux.The numbers in the x-axis represent the day of the year (DOY) for years 2011 to 2013 for each panel.

19 Figure 8 .
Figure 8. Summary of selected Fs responses to precipitation pulses (A-C), and conceptual depiction of the contingent responses of the ecosystem to precipitation pulses (D-F).The rectangles in D-F represent a "soil water bucket", where the "y" axis represents a stress gradient and the dotted line represents a reference point for comparison among panels.Selected responses are defined as Case I-III, and are identified in Figures 1, 3 and 7. Case I results in less stressful conditions for several days, Case II is characterized by a short period of stress relief, and Case III is in between fluctuations of less water stress.Fsveg = soil CO2 efflux in areas under vegetation; Fsbare = soil CO2 from areas without vegetation (i.e., bare soil); SM = soil moisture, where the number of "+" signs denotes higher soil moisture either at bare soil or under vegetation.

Figure 8 .
Figure 8. Summary of selected Fs responses to precipitation pulses (A-C), and conceptual depiction of the contingent responses of the ecosystem to precipitation pulses (D-F).The rectangles in D-F represent a "soil water bucket", where the "y" axis represents a stress gradient and the dotted line represents a reference point for comparison among panels.Selected responses are defined as Case I-III, and are identified in Figures 1, 3 and 7. Case I results in less stressful conditions for several days, Case II is characterized by a short period of stress relief, and Case III is in between fluctuations of less water stress.Fs veg = soil CO 2 efflux in areas under vegetation; Fs bare = soil CO 2 from areas without vegetation (i.e., bare soil); SM = soil moisture, where the number of "+" signs denotes higher soil moisture either at bare soil or under vegetation.