Biogeosciences How do variations in the temporal distribution of rainfall events affect ecosystem fluxes in seasonally water-limited Northern Hemisphere shrublands and forests ?

Rainfall regimes became more extreme over the course of the 20th century, characterised by fewer and larger rainfall events. Such changes are expected to continue throughout the current century. The effect of changes in the temporal distribution of rainfall on ecosystem carbon fluxes is poorly understood, with most available information coming from experimental studies of grassland ecosystems. Here, continuous measurements of ecosystem carbon fluxes and precipitation from the worldwide FLUXNET network of eddy-covariance sites are exploited to investigate the effects of differences in rainfall distribution on the carbon balance of seasonally water-limited shrubland and forest sites. Once the strong dependence of ecosystem fluxes on total annual rainfall amount is accounted for, results show that sites with rainfall distributions characterised by fewer and larger rainfall events have significantly lower gross primary productivity, slightly lower ecosystem respiration and consequently a smaller net ecosystem productivity.


Introduction
During the 20th century, precipitation intensity and the frequency of large rainfall events increased (Easterling et al., 2000;Trenberth, 2011).There was a globally consistent increase in metrics of rainfall extremes (Tank and Konnen, 2003;Karl et al., 1996;Peterson et al., 2008;Alexander et al., 2006).Anthropogenic climate change is expected to intensify this shift towards fewer, larger rainfall events (Trenberth et al., 2003).For the 21st century, general circulation models predict little change in total rainfall, but an increase in the frequency of heavy rainfall events (Tebaldi et al., 2006).Many regional climate change projections anticipate increased seasonality in precipitation (Manabe and Wetherald, 1987), and consistently drier summers over entire regions (e.g., southern and central Europe: Rowell and Jones, 2006).
While the impacts of differences in total rainfall amount on terrestrial ecosystem productivity have been studied for a range of ecosystems and climate regions (Wu et al., 2011;Huxman et al., 2004b), the ecological implications of changes in the temporal structure of rainfall, i.e. a shift towards fewer and larger rainfall events, have received less attention; results have so far been inconclusive (Jentsch et al., 2007).Almost all experimental field studies manipulating rainfall amount or distribution have been conducted on annual plant communities or perennial grasslands (Knapp et al., 2002;Miranda et al., 2009), because of the difficulty of performing manipulation experiments on taller canopies (exceptions include Volder et al., 2010;Limousin et al., 2010;Misson et al., 2010).Despite the potential impact of changes in temporal patterns of rainfall distribution on shrublands and forests, data for these ecosystems is sparse.There are even few remote sensing studies that address this issue (Fang et al., 2005;Good and Caylor, 2011).
Studies on grassland ecosystems indicate that increases in rainfall variability may exert as strong a control on future biogeochemical cycles as changes in total rainfall, atmospheric warming or increases in atmospheric carbon dioxide concentrations (Knapp et al., 2002;Weltzin et al., 2003).It is clear that the "repackaging" of rainfall and its "translation" into soil moisture available for uptake by plants is strongly influenced by the temporal distribution of rainfall events (Heisler-White et al., 2008;Loik et al., 2004), but there is no clear consensus on the impact of varying rainfall distributions on the productivity of different ecosystems.From an ecohydrological point of view, the optimal distribution of rainfall event size and frequency for a given biogeochemical process is expected to vary with total rainfall quantity, and to depend on the details of the processes of interception, run-off, drainage and evapotranspiration (Rodríguez-Iturbe et al., 2001;Porporato et al., 2001;Laio et al., 2001).The effects on these processes of changes in the temporal pattern of delivery of rainfall, including feedbacks on transpiration and soil water balance via altered stomatal conductance, will depend on soil water retention properties and ecosystem type as well as total annual rainfall amounts, making prediction and modelling difficult.Knapp et al. (2008) propose a conceptual model, hypothesising that greater rainfall variability will increase soil water limitation in mesic systems but decrease stress in xeric systems, because increased soil water content variability will be unfavourable for normally unstressed mesic systems but favourable for predominantly stressed xeric systems.
The effects of changes in rainfall variability are likely to differ between ecosystems because of variations in canopy structure, rooting depth and the ability of some dominant species to tolerate water stress (Porporato et al., 2001).In particular, the effects of changes in rainfall variability on net ecosystem productivity (NEP) will depend on the magnitude and timing of changes of gross primary productivity (GPP) and ecosystem respiration (RE).There is evidence that the response timescales of microbial soil respiration and plant photosynthesis to rain pulses can be quite different (Williams et al., 2009;Huxman et al., 2004b;Scott et al., 2006); the relative timing as well as magnitudes of these changes will most likely vary between ecosystems (in this paper, we consider both RE and GPP to be positive, and use the term net ecosystem productivity, NEP, to refer to the difference NEP = GPP − RE ).
In this study, we focus on relationships between rainfall and ecosystem fluxes in seasonally dry shrublands and forests.These ecosystems are characterised by periods of pronounced summer drought, and the primary limitation on plant productivity is water availability, with water supply from precipitation significantly smaller than plant demand (as measured by potential evapotranspiration).Plants in these climatic regions have developed a range of adaptations to the intermittency of water supply and long periods of drought, including drought-deciduousness (e.g.some Quercus species), sclerophyllous leaf habits, stomatal regulation (Zavala and de la Parra, 2005;David et al., 2007;Rambal, 1992), and tap roots to exploit water reserves in deeper soil layers (Schenk and Jackson, 2002;Moreno et al., 2005;Rambal, 1984).Although not necessarily higher in diversity than grasslands at the level of species counts, these communities are physiologically and ecohydrologically more complex: woody perennials have more extensive root systems than grasses, allowing them to access water from different soil layers, characterised by different temporal responses to rainfall events.In addition, the combination of species of several different life habits opens up possibilities for complex inter-species interactions.These interactions and speciesdependent adaptations to water stress make it unlikely that a simple model of the effects of changes in rainfall distribution will serve for these mixed communities.
One approach to addressing these difficult questions is to exploit the continuous eddy covariance measurements of surface-atmosphere exchanges of carbon dioxide, water, and energy now being made at hundreds of research sites globally (Baldocchi, 2008).Flux tower data allow direct quantification of NEP and its decomposition into GPP and RE (Reichstein et al., 2005) and make it possible to analyse relationships between ecosystem fluxes and rainfall characteristics across ecosystem types and sites in a robust way.
Here, we use FLUXNET data to examine the relationship between temporal patterns of rainfall distribution and ecosystem fluxes in seasonally water-limited shrubland and forest ecosystems.We address four questions.
1. How does the distribution of rainfall event size influence ecosystem fluxes in shrublands and forests, independent of total rainfall amounts?
2. How does the influence of differences in rainfall distribution compare to the influence of differences in total rainfall amount?
3. Do drier and wetter shrublands and forests differ in sensitivity to these factors?eddy covariance and meteorological measurements at 30-min temporal resolution from 966 site-years at 253 sites.NEP flux data are gap-filled and partitioned into GPP and RE using common standardised algorithms (Moffat et al., 2007;Papale et al., 2006;Reichstein et al., 2005).From the available sites, we identified ecosystems experiencing a climate with a dry season by selecting sites based on their Köppen climate classification (Kottek et al., 2006).We included sites in the Northern Hemisphere with arid/semi-arid (BSh, BSk), Mediterranean (Csa, Csb) and sub-tropical (Cfa) climates.
From those sites, we selected shrubland, woody savanna, deciduous broadleaf, evergreen and mixed forest sites that had not been recently disturbed or heavily managed.We then eliminated site-years that did not contain a sufficient proportion of high quality data: specifically, we retained only those site-years where at least 80 % of the half-hourly data were either original or gap-filled with high confidence (as indicated by the quality flags provided in the La Thuile dataset).This site selection process yielded the 28 sites and 85 site-years displayed in Table 1.

Annual rainfall amounts
We are interested in the impact of variations in the temporal distribution of rainfall, so need to control for the effect of overall variation in total rainfall, constructing models for the effect of variations in total rainfall on ecosystem fluxes, then using residuals of fluxes with respect to these models to analyse the effects of variations in rainfall distribution.We thus first examine the effect of variations in total annual rainfall amount (P tot , mm yr −1 ) on annual ecosystem fluxes (GPP, RE, and NEP, g C m −2 yr −1 ).We separate and exclude frozen precipitation from total annual liquid rainfall by counting only precipitation where the three-day moving mean screen air temperature is above freezing.It is not immediately clear what might be the best model to use to capture the variability due to total rainfall amount, since there is significant spread in the flux data.The most straightforward approach is to use a piecewise-constant model, with one flux value for "dry" sites and one for "wet" sites, using the mean flux values for each group of sites as unbiased estimators of the group flux values.The stratification into "dry" and "wet" sites is chosen to maximise the variance in fluxes (both GPP and RE) explained by the split (Table 2).Although this approach to modelling effects due to total rainfall amount is simplistic, by optimising the partition into "dry" and "wet" sites we are able to explain a greater proportion of the variance in fluxes than, for example, a more complicated exponential model fitted using a nonlinear least squares approach.The ecosystem flux residuals with respect to this simple piecewise-constant model show no linear correlation with the total annual precipitation amount (data not shown, but all slopes near zero and p-values greater than 0.4).
The only possible disadvantage to this optimisation approach is that it partitions the sites rather unevenly, with only 10 sites and 21 site-years in the "dry" group as compared to 18 sites and 64 site-years in the "wet" group.This is a consequence of the smaller variance of ecosystem fluxes for larger precipitation totals compared to smaller totals.From that perspective, the uneven distribution into "dry" and "wet" groups may be viewed as an advantage, since it eliminates a potential source of spurious within-group variability, although it does mean that sample sizes for comparing effect sizes in the "dry" group are rather small (Sect.2.3).An ancillary benefit of using a split into "dry" and "wet" sites that maximises the variance explained by the overall precipitation amount effect is that it removes a potentially arbitrary choice for a division into drier and wetter sites.

Temporal structure of rainfall
We next examine the effect on ecosystem carbon fluxes of differences in the temporal distribution of rainfall throughout the year.In order to compensate for differences in overall annual rainfall amount, we consider residuals of ecosystem fluxes with respect to our piecewise-constant flux/total rainfall model (Sect.2.2), denoted by GPP, RE and NEP, all measured in g C m −2 yr −1 .For NEP, we consider residuals with respect to the difference of our models for GPP and RE, in keeping with the definition of NEP above, i.e.NEP = GPP -RE.We hypothesise that at least some proportion of these flux residuals can be explained by differences in the distribution of rainfall through the year, in particular the degree to which the precipitation regime is characterised by "extreme" rainfall, i.e. large rainfall events separated by long periods of dryness.
We characterise the temporal distribution of rainfall events using two measures.The first is the annual mean daily rainfall intensity, i.e. the annual mean rainfall amount for days with rain, in mm d −1 , denoted by I .For a given total rainfall amount, rainfall regimes characterised by fewer and larger rainfall events will have larger I than regimes with lighter, more "continuous" rainfall.The second measure is the proportion of total annual rainfall amount due to days whose daily rainfall is greater than the 95th percentile of daily rainfall, denoted by R 95 .This statistic is commonly defined in terms of the climatological distribution of daily rainfall, (e.g.Alexander et al., 2006), but since we do not have detailed climatology information for our sites, we instead define an R 95 statistic based on each site-year's distribution of daily rainfall.The resulting statistic carries much the same information as the standard R 95 statistic, in that larger values of the statistic indicate rainfall less evenly distributed throughout the year, with larger and more widely spaced rainfall events.
The principal difference between our pragmatic singleyear definition of R 95 and the climatologically based one is that, in a long-term climatic time series, years are likely to exist without large rainfall events, meaning that the baseline precipitation distribution used for the calculation of the R 95 statistic would have a smaller 95th percentile (compared to Table 1.Study sites.Shown are FLUXNET site identifiers, site locations, years of data available and count of data site-years for each site, IGBP ecosystem classification (CSH: closed shrubland, DBF: deciduous broadleaf forest, EBF: evergreen broadleaf forest, ENF: evergreen needleleaf forest, MF: mixed forest, SAV: savanna, WSA: woody savanna), Köppen climate classification (Kottek et al., 2006; BSh: hot semiarid; BSk: cold semi-arid; Cfa: humid subtropical; Csa: hot-summer Mediterranean; Csb: warm-summer Mediterranean), site elevation (m), climatological annual rainfall (mm), climatological annual mean air temperature (T a , • C) and site description reference.

Site
Lat the distribution in a single year that does have large rainfall events).We would thus find larger R 95 values for the years considered in this case than we do using the "single year" R 95 statistic used here, since more of the rainfall in a year with large rainfall events would occur beyond the 95th percentile of the reference distribution.In general, it is impossible to say what impact this would have on the results of our study, since statistics for a long climatic reference period are not available to us for most of our sites.For some sites, it is possible that the years that we examine have an anomalously high number of large rainfall events compared to the long-term climate, while for other sites the years that we examine may have less large rainfall events than the long-term distribution.
In order to examine the effects that our choice of R 95 statistic might have, the best that we can do is to consider the impact of our definition of R 95 based on data from single years in the context of sites where we have multiple years of data.Figure 1 shows differences between R 95 values calculated in two different ways, first using all of the available years of data for each site together to calculate the 95th percentile of daily rainfall, then calculated using single siteyears of data.Although this approach gives only a rough indication of the likely effects of using a single year-based R 95 , we can gain some idea of whether the definition of R 95 that we use is reasonable.The results shown in Fig. 1 broadly conform to our expectations: there is some spread in the differences between the two calculation methods, but there is not a notable asymmetry, i.e. the R 95 values based on all site-years of data for each site are not very much different from those calculated using a single year of data in each case (the mean difference across all sites and years is −1.4 %).Although only a heuristic guide, we feel that this justifies the form of the R 95 statistic that we use in the remainder of the paper, calculated based on single site-years of data.This allows us to use a consistent statistic across all sites, independent of the number of years of data we have for the site.
Within each precipitation group defined in Table 2 and for all sites together, standardised residuals of ecosystem fluxes with respect to our piecewise-constant flux/total rainfall model were related to indices of rainfall distribution using linear regression.Rainfall indices were calculated from daily time series of rainfall amounts for each site-year, using a daily rainfall threshold of 5 mm to define days with rain.
There are two distinct potential types of variability in the flux residuals that we must address: inter-site spatial variability and intra-site temporal variability.We present three distinct analyses to examine these different types of variability, one combining inter-site spatial and intra-site temporal variability, one isolating inter-site spatial variability and the third isolating intra-site temporal variability at sites where sufficient data exists.
For the analysis of combined spatio-temporal variability, we consider standardised flux residuals, Z GPP , Z RE and Z NEP , and rainfall predictors, Z I and Z R 95 , where the standardised score for a variable X taking value x is defined as Here µ X and σ X are respectively the population mean and standard deviation of X, here calculated across values for all site-years.
To isolate the inter-site spatial variability, we consider the same ecosystem flux residuals and precipitation distribution predictors as for the combined spatio-temporal analysis, but standardise the data in a different way.We use mean values of all the site-years available for each site, then standardise these site-level values using the across-site mean and standard deviation.This isolates between-site gradients and removes within-site, year-to-year variation.We denote these site-level standardised variables as Z (S) X , for a variable X, one value per site.We then use these site-level standardised values to perform the same kind of regression analysis as for the combined spatio-temporal variability case.
We have enough sites (28, in total) that we can draw some conclusions concerning the inter-site spatial variability, but the situation regarding temporal variability within sites is www.biogeosciences.net/9/1007/2012/Biogeosciences, 9, 1007-1024, 2012 much less hopeful.Out of the 28 sites, we have only 8 with four or more years of data and 13 sites with only one or two years of data (Table 1).Considering just the sites with four or more years of data, we can calculate, for example, precipitation intensity calculated for each year at each site, and use the site-level mean across years and standard deviation across years to produce standardised values.In the analysis of intra-site temporal variability below, we denote these standardised values by Z (T ) X for a variable X, one value per site per year.The intention here is to remove the inter-site spatial variability and isolate intra-site anomalies.

Biometeorological conditions
Rainfall distribution is strongly linked to other meteorological and hydrological variables.In particular, there are strong relationships with mean annual vapour pressure deficit, VPD (kPa), annual incoming solar radiation, R G (kJ yr −1 ), and soil water status.We present results below (Sect.3.4) that relate these quantities to the measures of rainfall intermittency we use.In order to provide a homogeneous measure of soil water storage across sites, including those where soil water content was not directly measured, we characterise soil moisture status using the mean daily relative soil water deficit (WD, %), defined as WD = 365 i=1 WD i /365, where WD i is the relative soil water deficit for day i of the year, WD i = 1 − (θ i − θ min )/(θ max − θ min ); here θ i is the mean daily soil water content (mm), and θ min and θ max are the annual minimum and maximum mean daily soil water content (mm), all taken from the upper layer of a simple bucket model driven by daily meteorological and radiation inputs (Reichstein et al., 2002).

Annual rainfall amounts
As expected, there is a strong relationship between total rainfall and both GPP and RE (Fig. 2a and b).Because the rate of increase of GPP with total rainfall amount is greater than the rate of increase of RE, NEP also increases with total rainfall amount (Fig. 2c).
We find that the relationship between ecosystem fluxes and total annual rainfall amounts is mostly simply treated using a simple piecewise-constant model that assigns a single flux value to all "dry" sites and a different value to all "wet" sites (modelled values are shown as horizontal lines in Fig. 2).Although a simplistic approach, this model, using an optimal partition of sites into "dry" and "wet", serves to explain a larger fraction of the variance in ecosystem flux values than either a simple linear regression or a exponential model (Table 3).

Temporal structure of rainfall
We first present an analysis of the joint spatio-temporal variability in our data.We consider regressions of standardised ecosystem flux residuals, Z GPP , Z RE and Z NEP , onto Z I (Fig. 3a-c and upper half of Table 4).Overall, negative correlations are observed between ecosystem flux residuals and rainfall intensity: regressions between standardised ecosystem flux residuals and Z I are significant at the 5 % level for GPP for all sites together and for "dry" sites and "wet" sites alone, and for RE, for all sites together and for "wet" sites.This appears to indicate that, after controlling for overall variations in annual rainfall amount, ecosystem fluxes are reduced at sites with rainfall regimes characterised by greater I .
The regressions against Z R 95 displayed in Fig. 3d-f and the lower half of Table 4 show similar behaviour to those in terms of Z I : there is a significant dependence of ecosystem flux residuals on the intermittency of the rainfall regime, here as measured by the proportion of total rainfall due to the upper tail of daily rainfall events.Sites with greater R 95 have both smaller GPP (significant at the 5 % level for all site-years together, and for "dry" sites and "wet" sites alone) and RE (significant at the 5 % level for all site-years together and for "wet" sites alone).The relative strengths of the GPP and RE residual dependence on R 95 are such that there is a significant linear regression between NEP residuals and R 95 (significant at the 5 % level for all site-years together, and for "dry" sites and "wet" sites alone).
We next seek to determine how much of the total variability exhibited in Fig. 3 and Table 4 is due to spatial variations between sites and how much is due to intra-site temporal variability.
Figure 4 and Table 5 show equivalent results to Fig. 3 and Table 4, but based on regressions between site-level standardised residuals, Z (S) X , thus isolating variability due to inter-site spatial effects.The results of this analysis are broadly similar to those of the combined spatio-temporal variability analysis.There are significant negative correlations between standardised precipitation intensity and R 95 and ecosystem flux residuals, indicating that sites whose rainfall regimes on average are characterised by fewer and larger rainfall events, for the same mean total rainfall amount, have less gross primary Each point represents a single site-year ("dry" sites as squares, "wet" as circles, colour-coded by the IGBP ecosystem type: grey -CSH, green -DBF, blue -EBF, black -ENF, magenta -MF, orange -SAV, red -WSA).Horizontal lines show the piecewise-constant model values for fluxes of the "dry" (black) and "wet" (grey) groups of sites; for NEP, the lines show the differences between the GPP and RE models.Model parameters and standard errors for the site-year fits are shown at the top of the GPP and RE panels.Rectangles cover the range of interannual variability for each site.
productivity and less ecosystem respiration, with the relative sizes of these effects being such that there is also a negative correlation between rainfall distribution measures and NEP.
The picture for intra-site temporal variability, as measured by regressions between per-site temporally standardised quantities Z (T ) X , is much less clear.Unfortunately, although there appear to be effects of the type that we would expect (smaller ecosystem fluxes for greater Z (T ) I or Z (T ) R 95 ), the lack of data renders most of these regressions very weak.Table 6 shows a few of the more significant regression results.In most cases, there is a negative correlation between changes in the rainfall distribution predictors and ecosystem fluxes, but the results are too weak to draw any substantive conclusions.

Relative effect sizes
Having demonstrated that, after accounting for overall differences in rainfall amount, there remains a significant effect on ecosystem fluxes due to differences in rainfall distribution, most of which should apparently be ascribed to inter-site spatial variation, a natural question is to ask what proportion of the variability in ecosystem fluxes is due to differences in rainfall distribution compared to the portion of the variability due to overall rainfall differences.
We can most easily assess this difference by comparing ecosystem flux variations due to differences in rainfall distribution to the difference in fluxes assigned to the "dry" and "wet" precipitation groups by our piecewise-constant total rainfall amount model: These inter-group flux differences can then be compared to the 1SD values in Tables 4, 5 and 6, which show the changes in ecosystem fluxes associated with a one standard deviation change in the Z I or Z R 95 predictor variables.In each case, the appropriate standard deviation is used, i.e. for the combined spatial and temporal variability in Table 4, the standard deviation is across all site-years, for the inter-site spatial variability in Table 5, the standard deviation is across site mean values (means across the available data years for each site), and for Table 6, the standard deviation is calculated separately for each site across the available years of data for the site.The aim here is to provide an indication of the variability in ecosystem fluxes associated with typical variations in rainfall distribution: in this context, regression slopes can be thought of as the change in ecosystem flux values associated with a change of one standard deviation in the predictor variable.We use the appropriate ecosystem flux standard deviations to convert these slopes back into indicative changes in fluxes for comparison between the different types of variability that we consider.(d, e, f).Each point represents mean conditions for a single site ("dry" sites as squares, "wet" as circles, colour-coded by the IGBP ecosystem type: grey -CSH, green -DBF, blue -EBF, black -ENF, magenta -MF, orange -SAV, red -WSA).Lines show linear regressions significantly different from zero at the 5 % (thick lines) or 10 % (thin lines) level, solid lines for the "dry" (black) and "wet" (grey) groups separately and dashed lines for the whole data set.
Summary results are shown in Table 7.Although we do not claim that the approach taken here can give more than a rough idea of the relative effect sizes for ecosystem flux variations due to differences in overall rainfall amount compared to differences in rainfall distribution, it is clear from Table 7 that increases in I or R 95 have a consistently negative impact on ecosystem fluxes.This impact is smaller than the impact of differences in overall rainfall amount (as measured by the flux differences between the "dry" and "wet" groups), but is of a similar order of magnitude.For instance, for GPP, the impact of a one standard deviation change in R 95 on the inter-site (S) variation is −185 g C m −2 yr −1 , as compared to a total rainfall amount effect size of 851 g C m −2 yr −1 .For both GPP and RE, the effect sizes for changes in rainfall distribution are around 15-20 % of the size of those for differences in overall rainfall amount, indicating that differences in the distribution of rainfall, in terms of changes in the number and size of rainfall events and the lengths of periods of drought, can have a significant effect on ecosystem fluxes in the ecosystems that we consider here.The results in Table 7 are those for all sites; if we consider the "dry" sites or "wet" sites alone, the rainfall distribution effect sizes tend to be larger for dry sites than for wet sites (see column 1SD in Tables 4 and 5).For drier sites, as defined by annual total rainfall, variations in rainfall distribution thus contribute more strongly to variations in ecosystem carbon balance than at wetter sites, a result consistent with observations (Harper et al., 2005).

Biometeorological conditions
Figure 5 shows mean biometeorological conditions for the two precipitation groups.Drier sites experience greater soil water deficit, atmospheric vapour pressure deficit and Table 7. Ecosystem flux effect sizes for total rainfall amount differences (P ) and differences in rainfall distribution measures (precipitation intensity I and R 95 ), for combined spatio-temporal (S + T ), spatial (S) and temporal (T ) regressions.All effect sizes are measured in g C m −2 yr −1 and regression-based values are for all sites together (for S + T and S), using the regression slope and the appropriate standard ecosystem flux standard deviation.Only effect sizes from regressions significant at the 5 % level are shown, and values in parentheses for the temporal regressions are indicative values from single sites.incoming radiation than wetter sites.These biometeorological variables were all significantly correlated with the R 95 statistic (Fig. 6).Mean daily relative soil water deficit, mean annual vapour pressure deficit and mean incoming radiation all increase in precipitation regimes characterised by fewer and larger rainfall events.

Rainfall variability and ecosystem fluxes
For the shrubland and forest sites considered here, rainfall regimes characterised by fewer and larger rainfall events, irrespective of total rainfall amount, have a negative impact on ecosystem fluxes.The effects of changes in rainfall distribution appear to be slightly greater at drier sites (total annual rainfall ≤725 mm).By examining the variability of ecosystem fluxes associated with typical spatial and temporal variations in total rainfall amount and rainfall distribution indices, we found that variations in rainfall distribution have a smaller but still considerable effect on ecosystem fluxes compared to variations in overall rainfall amount.
For drier sites, our results apparently contradict the conceptual model of Knapp et al. (2008).Their model implies that rainfall regimes characterised by fewer and larger rainfall events are expected to have a negative impact on NEP in mesic ecosystems, but a positive impact in xeric ecosystems.They argue that soil water fluctuations are amplified by decreasing precipitation frequency and increasing precipitation intensity, and that these larger fluctuations in soil water content will cause mesic systems to spend more time in conditions of water stress and xeric systems less (Knapp et al., 2008, Fig. 5).Previous observational studies have provided some evidence in support of this model in grassland ecosystems: Knapp et al. (2002) found that increasing precipitation interval decreased soil moisture and ANPP in a mesic grassland and Thomey et al. (2011) found results consistent with a slight modification of Knapp et al.'s model in an arid to semiarid grassland, although Heisler-White et al. (2008) observed the opposite effect in a semi-arid grassland.
It is likely that the discrepancy between these grassland studies and our results are due to differences in soil-plantwater interactions in the vegetation types considered.Grassland species are generally shallow-rooted and consequently sensitive to changes in water content in upper soil layers.
In dry climates, grasses have a short growing season (generally triggered by intense rainfall) and are highly sensitive to short rain pulses.Trees and shrubs, in contrast, rely more on year-round water availability.Knapp et al.'s model of soil water storage is a bucket model with a single shallow layer.In xeric ecosystems, such a model will necessarily lead to an increase in water content in the uppermost soil layer under rainfall regimes characterised by fewer and larger rainfall events, thus reducing water stress for shallow-rooted plants.Plants with deeper roots may respond to infrequent heavy rainfall events rather differently and may not benefit from short-term changes in soil moisture in the upper soil layers in the way predicted by Knapp et al.'s model.In fact, the response of vegetation to rain events is likely to be highly dependent on plant life form, which is itself a function of climate and soil conditions (Sala and Lauenroth, 1985;Schenk and Jackson, 2002;Ogle and Reynolds, 2004;Viola et al., 2008).In dry climates, Porporato et al. (2001) found that there is an intermediate rainfall frequency at which drought stress is minimum, and that this optimal frequency depends on plant and soil properties (see also Rodríguez-Iturbe et al., 2001 andLaio et al., 2001).An increase or decrease in precipitation frequency from this optimum at constant total rainfall will increase drought stress.
Given the number of different factors affecting productivity in these woody ecosystems and the complexity of the interactions between the factors, the prospects for introducing a conceptual model to serve as an alternative to the model of Knapp et al. (2008) appear poor.The Knapp et al. (2008) model is appropriate for arid and semi-arid grassland ecosystems where soil-water-plant interactions are relatively simple because of the restriction of the vegetation to shallow-rooted species of low stature.In woody ecosystems, soil-waterplant interactions are potentially much more complex, with multiple vegetation layers with different rooting depths and consequently different responses to changes in precipitation regime.
Much work has been done in physiological and gap modelling of seasonally dry forests in the Mediterranean -in particular, Zavala and de la Parra (2005) coupled a physiological model of a single Mediterranean tree species (Quercus ilex) with a stochastic model of rainfall variability and soil moisture based on the ideas of Rodríguez-Iturbe and Porporato (2004).They considered variations in rainfall regime relevant to our considerations here, but the results of these modelling efforts are difficult to interpret in a fashion that would allow us to offer a "schematic" model to compare with Knapp et al.'s grassland model, despite the fact that Zavala and de la Parra's model is an effort to build a more analytically tractable ecosystem model to help answer just the type of question we are addressing.A further impediment to presenting a conceptual model based on such work is that the applicability of the extensive work done on Mediterranean forests and savanna ecosystems to other regions of the world is uncertain.As far as our results for shrubland and forest ecosystems are concerned, a recent remote sensing study associating fractional woody cover with precipitation across Africa is of interest: Good and Caylor (2011) found that more frequent, less intense rainfall consistently led to a higher fraction of woody cover for a given amount of total rainfall, for a range of total rainfall amounts, a result more consistent with our results than with the hypothesis of Knapp et al. (2008).Although Good and Caylor's study is concerned with vegetation in Africa rather than the Northern Hemisphere regions we consider here, it is striking how consistent is the relationship between woody fraction and precipitation distribution across the whole range of African climates.
Of the factors that we have not considered here, perhaps the most important is variation in soil texture and type between sites.The response of soil respiration to rainfall events is strongly dependent on soil type (Inglima et al., 2009), and in many of the drier sites considered here, there may be additional geochemical carbon fluxes associated with the presence of calcareous substrates (Kowalski et al., 2008;Serrano-Ortiz et al., 2009).

Secondary factors
In general, we observe that rainfall regimes characterised by fewer and larger rainfall events, are correlated with increases in soil water deficit (Fig. 6a), atmospheric vapour pressure deficit (Fig. 6b), and incoming solar radiation (Fig. 6c).The effects on atmospheric vapour pressure deficit and incoming radiation both arise because less frequent rainfall corresponds to longer drought periods.
It thus appears that under such rainfall conditions, not only is soil water availability lower, but atmospheric water demand is greater.These hydrological factors are more limiting in drier climates, so we expect them to have a greater effect on ecosystem processes for drier sites than for wetter sites (Fig. 5a and b).Indeed, there is evidence that high vapour pressure deficit plays a strong role in controlling stomatal conductance and photosynthetic efficiency in conditions of soil water limitation (Aires et al., 2008;Siam et al., 2008).
Water limitation in drier climates also prevents plants from making effective use of the generally higher radiation levels seen in conditions of infrequent intermittent rainfall.At wetter sites, biological activity is more likely to be limited by other constraints such as soil nutrient or light availability (Huxman et al., 2004a).In addition, excess light producing photo-oxidative damage may contribute a supplementary stress factor for vegetation at drier sites (Martínez-Ferri et al., 2000).
Because of its direct and overriding influence on plant function, most ecohydrological studies have focused on the effect of rainfall variability on ecosystems via its control on soil water content (Eamus et al., 2006;Eagleson, 2002;Rodríguez-Iturbe and Porporato, 2004).However, atmospheric water demand and incoming solar radiation are potentially important secondary factors that should be also taken into account.Further experimental and modelling studies will be needed to disentangle the effect of these different factors.

Differential ecosystem flux responses
Our results show that the influence of variations in rainfall distribution is generally stronger for GPP than RE, leading to www.biogeosciences.net/9/1007/2012/Biogeosciences, 9, 1007-1024, 2012 a decrease in NEP for rainfall regimes with fewer and larger rainfall events.Studies in a mesic grassland have previously observed a decrease in both leaf photosynthetic carbon gain and soil respiration as a result of more intermittent rainfall patterns (Knapp et al., 2002;Fay et al., 2002;Harper et al., 2005), while Fay et al. (2008) reported increases in leaf-level photosynthesis and decreases in soil respiration due to increased rainfall event size for levels of total rainfall ranging from 400 to 1000 mm yr −1 .These studies did not provide information concerning the relative sensitivity of these individual processes to variations in rainfall pattern, so it is not possible to reconstruct the dependence of overall ecosystem carbon flux on precipitation frequency for these sites.
It is likely that the greater sensitivity of GPP compared to RE to changes in precipitation frequency is due to a combination of contributions from the factors described in Section 4.2, i.e.GPP is affected by concurrent increases in soil water stress, atmospheric vapour pressure deficit and excess radiation.RE is indirectly affected by these factors, through its dependence on GPP (Migliavacca et al., 2011), but is more directly responsive to changes in soil moisture.Moreover, increases in air temperature associated with decreased precipitation frequency may even provide a positive influence on RE, partially compensating for the decrease in RE due to decreased soil water content (Davidson et al., 1998).Other studies in arid and semi-arid ecosystems have observed differing responses of GPP and RE to rainfall pulse size (Huxman et al., 2004b;Potts et al., 2006;Arneth et al., 2006;Williams et al., 2009).The structure of respiration response has been attributed to the high sensitivity of soil microbes to rainfall events following drought conditions (Inglima et al., 2009;Misson et al., 2006;Xu et al., 2004;Lee et al., 2004;Jenerette et al., 2008), although the relationship between pulse size and duration of active soil respiration seems to saturate at moderate event sizes (Huxman et al., 2004b;Sponseller, 2007), perhaps due to exhaustion of labile carbon sources.The response of vascular plant photosynthetic activity to rainfall events is generally of longer duration than the response of microbial respiration, with both the magnitude of the response (Ignace et al., 2007;Chen et al., 2009;Williams et al., 2009;Scott et al., 2006) and the duration of physiological activity (Huxman et al., 2004b;Williams et al., 2009) increasing following larger rainfall pulses.
The differential responsiveness of respiration and photosynthesis to discrete rainfall events seen in the studies reported above leads to variations in NEP with changes in precipitation intensity and distribution independent of total rainfall amount, moderated by other factors, such as temperature, light availability, initial soil water content or canopy conditions (Schwinning and Sala, 2004).Arneth et al. (2006) found remarkable plasticity in canopy photosynthetic parameters in response to intermittent dry periods during the rainy season in a semi-arid woodland.Ultimately, any decoupling between GPP and RE can only be temporary because, at larger temporal and spatial scales, respiration fluxes are controlled by substrate supply (Campbell et al., 2004;Janssens et al., 2001;Reichstein et al., 2003;Misson et al., 2007).

Conclusions
We have shown that more intermittent rainfall regimes, characterised by fewer and larger rainfall events, can have a strong negative effect on both GPP and NEP of woody ecosystems, particularly in drier climates, independent of total rainfall amount.Future amplification of the hydrological cycle caused by global warming may thus pose a threat to the productivity and sustainability of shrubland and forest ecosystems in these climates.Increased rainfall variability will add a supplementary constraint to the more commonly known limitation imposed by an expected future decrease in total annual or seasonal rainfall in these regions (Gao and Giorgi, 2008;Giorgi and Lionello, 2008).
Here, we consider only annual quantities for a subset of climate conditions and ecosystem types.FLUXNET sites archive continuous time series data collected at a temporal resolution of 30 min, a level of detail that should permit more sophisticated analysis of the questions addressed here, including a comprehensive cross-site examination of the phasing of GPP, RE and NEP changes following rainfall pulses.

Fig. 1 .
Fig. 1.Differences between R 95 statistics calculated based on a daily rainfall distribution derived from all precipitation data available for each site and R 95 values calculated for each site year separately.Points show values for individual site-years; vertical lines highlight the range of differences for each site.Sites are colour-coded by the IGBP ecosystem type: grey -CSH, green -DBF, blue -EBF, black -ENF, magenta -MF, orange -SAV, red -WSA.

Fig. 2 .
Fig. 2.Relationship between total annual rainfall and (a) GPP, (b) RE, and (c) NEP.Each point represents a single site-year ("dry" sites as squares, "wet" as circles, colour-coded by the IGBP ecosystem type: grey -CSH, green -DBF, blue -EBF, black -ENF, magenta -MF, orange -SAV, red -WSA).Horizontal lines show the piecewise-constant model values for fluxes of the "dry" (black) and "wet" (grey) groups of sites; for NEP, the lines show the differences between the GPP and RE models.Model parameters and standard errors for the site-year fits are shown at the top of the GPP and RE panels.Rectangles cover the range of interannual variability for each site.

Fig. 3 .Fig. 4 .
Fig. 3. annual ecosystem flux residuals (Z GPP , Z RE and Z NEP ) versus Z I (a, b, c) and Z R 95(d, e, f).Each point represents a single site-year ("dry" sites as squares, "wet" as circles, colour-coded by the IGBP ecosystem type: grey -CSH, green -DBF, blue -EBF, black -ENF, magenta -MF, orange -SAV, red -WSA).Lines show linear regressions significantly different from zero at the 5 % (thick lines) or 10 % (thin lines) level, solid lines for the "dry" (black) and "wet" (grey) groups separately and dashed lines for the whole data set.Rectangles cover the range of interannual variability for each site.

Table 2 .
Precipitation groups defined by a partition into "dry" and "wet" sites that optimises the flux variance explained by a piecewise-constant model of the effect of differences in total annual rainfall amount across all 28 sites.

Table 3 .
R 2 values for different models relating ecosystem fluxes to total annual rainfall amounts (all three models contain the same number of parameters, so these values are directly comparable between models).

Table 4 .
Slopes, intercepts (both with standard errors in parentheses) and R 2 values for linear regressions between standardised ecosystem flux residuals and precipitation variables, for all site-years and for precipitation groups defined in Table2.Regressions significant at the 5 % level are indicated in bold.Since all values are standardised, all values in the table are unitless, except for column 1SD , which shows the change in ecosystem flux (g C m −2 yr −1 ) associated with a one standard deviation change in the predictor variable, Z I or Z R 95 .

Table 5 .
Slopes, intercepts (both with standard errors in parentheses) and R 2 values for linear regressions between site-level standardised ecosystem flux residuals and precipitation variables, for all sites and for precipitation groups defined in Table2.Regressions significant at the 5 % level are indicated in bold.Since all values are standardised, all values in the table are unitless, except for column 1SD , which shows the change in ecosystem flux (g C m −2 yr −1 ) associated with a one standard deviation change in the predictor variable, Z I or Z R 95 .