## Abstract

The analysis of soil water partitioning in seasonally dry climates necessarily requires careful consideration of the periodic climatic forcing at the intra-annual timescale in addition to daily scale variabilities. Here, we introduce three new extensions to a stochastic soil moisture model which yields seasonal evolution of soil moisture and relevant hydrological fluxes. These approximations allow seasonal climatic forcings (e.g. rainfall and potential evapotranspiration) to be fully resolved, extending the analysis of soil water partitioning to account explicitly for the seasonal amplitude and the phase difference between the climatic forcings. The results provide accurate descriptions of probabilistic soil moisture dynamics under seasonal climates without requiring extensive numerical simulations. We also find that the transfer of soil moisture between the wet to the dry season is responsible for hysteresis in the hydrological response, showing asymmetrical trajectories in the mean soil moisture and in the transient Budyko's curves during the ‘dry-down‘ versus the ‘rewetting‘ phases of the year. Furthermore, in some dry climates where rainfall and potential evapotranspiration are in-phase, annual evapotranspiration can be shown to increase because of inter-seasonal soil moisture transfer, highlighting the importance of soil water storage in the seasonal context.

## 1. Introduction

Seasonal variations in climatic inputs (in particular, rainfall and potential evapotranspiration) have garnered considerable attention in recent years as controlling factors for mean annual soil water partitioning [1,2], plant responses and adaptive strategies [3,4], regional vegetation distribution, carbon fluxes and primary productivity [5–8], with further implications for agriculture and land management [9]. These scientific emphases on the role of climate seasonality come at a time of discernible climate change. For instance, increase in the global mean temperature compounded with significant drying in some regions is likely to lead to increase in the frequency and the intensity of seasonal droughts [10]. Meanwhile, rainfall seasonality and its interannual variability have been observed to change in magnitude, timing and duration in the tropics [11]. Thus, it is important under these contexts to be able to evaluate the extent to which intra-annual variations in climate inputs may influence various ecohydrological processes.

Climate seasonality is a defining feature of many dry ecosystems, often characterized by a distinct non-uniformity in their timing of annual rainfall. This results in one or two wet seasons during which most of the annual rainfall occurs, separated by prolonged dry periods. Such seasonal rainfall variations are regularly found in tropical dry, monsoon and Mediterranean climates [12]. In terms of potential evapotranspiration, tropical dry and monsoon climates exhibit high year-round net radiation with little fluctuation between the seasons owing to their presence in the lower latitudes, whereas Mediterranean climates, which are typically found on the western continental margins of the mid-latitudes, experience larger seasonal variations. Mediterranean climates are also marked by out-of-phase cycles of rainfall and potential evapotranspiration which manifest as hot, dry summers and cool, wet winters; this is particularly challenging for plant life, because high productivity is limited to the beginning of the summer season before soil moisture depletes [13]. It is rare in the study of these climates (and their coexisting ecosystems) to reference commonalities shared with other seasonally dry regions, perhaps because of their collectively expansive geographical extent and disparate meteorological drivers. However, it is reasonable from an ecohydrological standpoint to treat them under a common framework considering their similarly large hydroclimatic contrast between the seasons.

Our approach focuses on the analysis of two hydrologically important quantities in seasonal climates: the mean soil moisture and the evapotranspiration ratio (defined as the proportion of total rainfall lost from the soil through evapotranspiration). Budyko [14] found that the long-term evapotranspiration ratios calculated over steady-state conditions for a large number of watersheds fall onto a single curve when plotted as a function of their dryness indices (mean annual potential evapotranspiration over rainfall), and since then, Budyko's curve has remained a staple in the analysis of hydrological partitioning 15–18]. Previously, the effect of climate seasonality on soil water balance has been studied analytically using simple models [16] and numerically for more complex models [19,20]. Others have adopted phenomenologically derived relationships based on boundary conditions imposed by water or energy balances [2,21]. Soil water storage has been found to play an important role in seasonal climates in reducing losses through leakage and deep percolation, though such losses can be increased when seasonal rainfall and potential evapotranspiration are out-of-phase, as in Mediterranean climates [1,16,21]. Still, the problem of quantifying the mean pathways within the water cycle remains a challenging one, requiring suitable characterization of random-like hydroclimatic forcings which are simultaneously embedded in a periodic seasonal cycle.

Our main contribution here is to extend a physically based, stochastic model of soil moisture [22] into a seasonal context. Previously, this model has been used to analyse the effect of climate, soil and vegetation on various quantities of interest, including the mean soil moisture, plant biomass and carbon uptake and storage [1,6,22,23], though the descriptions of climate seasonality have so far been kept piecewise-linear to allow for mathematical tractability. The novelty of the work presented here is in integrating the daily and seasonal variabilities consistently within a stochastic framework. This bypasses the ad hoc use of seasonally ‘averaged‘ forcings (e.g. from monthly data) as driving factors in the soil water balance, because most soil water dynamics occur at a much smaller (e.g. daily) timescale. In addition, our simple models bring out new insights into seasonal soil moisture without resorting to impractically large numerical simulations. We construct not only the long-term evapotranspiration ratio, but also transient departures from Budyko's curve at the intra-annual scale, showing hysteresis in the mean soil moisture that is mediated by seasonal soil water storage. We neglect the contribution of capillary rise from groundwater or deeper layers and consider only surface water-dependent systems where groundwater is either very deep or its amount negligible. We also neglect the contribution of lateral flow. These additional complexities would require substantial modifications to the existing model and are left for future studies. Here, we focus on cases where soil moisture dynamics are dominated by seasonally varying, stochastic rainfall.

In the following sections, we begin by introducing the stochastic mean soil moisture model with three approximations for mean leakage/run-off and the corresponding Budyko's formulation. Next, the model results are presented using climatic parameters typical of tropical dry and Mediterranean climates, with attention to the role of the phase difference between rainfall and potential evapotranspiration. Then, we highlight seasonal hysteretic behaviour in the transient evapotranspiration ratio and other hydrological terms as a result of seasonal soil water storage and trace their projection along Budyko's framework. Finally, the effect of climate seasonality on the annually averaged Budyko's curve is presented.

## 2. Modelling mean soil moisture and Budyko's curve at the seasonal level

### (a) Mean soil moisture dynamics

At the daily scale, assuming negligible horizontal redistribution via topographic effects, soil moisture at a point is recharged through intermittent rainfall pulses of random depths, and depleted through evapotranspiration, leakage and run-off. We neglect any contribution from groundwater and focus only on surface water-dependent systems. In what follows, we focus on the vertically averaged, plant available ‘effective‘ soil moisture, *x*, at the daily timescale. The range of *x* brackets the upper and lower limits of soil moisture available for plant uptake, with *x*=0 occurring at the plant wilting point *s*_{w}, and *x*=1 at *s*_{1}, which is the threshold above which all soil water is assumed to be immediately lost through leakage and run-off. The threshold *s*1 is physically related to soil properties and is typically situated between field capacity and complete soil saturation [22,24]. Thus, *x* is simply a standardized version of relative soil moisture *s*, defined by *x*=(*s*−*s*_{w})/(*s*_{1}−*s*_{w}). With the above assumptions, the effective soil moisture balance, vertically averaged over the rooting zone, is [17,22,25]
*w*_{0}=*n*Z_{r}(*s*_{1}−*s*_{w}) is the maximum plant-available soil water storage volume per unit ground area, *n* is the vertically averaged soil porosity and Z_{r} is the rooting depth (e.g. cm). The rate of change in the total volume of plant-available soil moisture (*w*_{0}(d*x*(*t*)/d*t*)) is governed by the rate of rainfall *R*(*t*) (e.g. cm per day), evapotranspiration ET(*x*(*t*),*t*) (e.g. cm per day) and leakage/run-off LQ(*x*(*t*),*t*) (e.g. cm per day). Rainfall is considered as a time-dependent stochastic process and, at the daily scale, idealized as a marked Poisson process that is non-homogeneous in time, with a time-dependent rate parameter λ(*t*) and the depth of rainfall drawn independently from an exponential distribution of mean *α*(*t*) [25]. The sources of seasonality may come from gradual changes in the mean frequency of rainfall and in the parameters that control evapotranspiration.

Given the stochastic nature of all variables presented in equation (2.1), it is useful to introduce some notation here to distinguish between ensemble averages and time averages. The ensemble average of a generic stochastic variable *p*(*u*,*t*) at time *t* is denoted by brackets as *p*(*u*,*t*) has an associated quasi-steady-state pdf *p*_{ss}(*u*,*t*), which is produced by applying the instantaneous conditions found at *t* constantly over an extended period of time, until *p*(*u*,*t*) has reached steady state, yielding *p*_{ss}(*u*,*t*). The resulting ensemble average produced from *p*_{ss}(*u*,*t*) is denoted by 〈*u*_{ss}(*t*)〉. In parallel, we also make use of temporal averages which are represented by overbars, defined as *t*_{0} and *T* are the initial time and period of interest. In the analyses that follow, temporal averages are taken of the ensemble average 〈*u*(*t*)〉 over a year (*T*=*T*_{year}), which is the natural period over which seasonal climatic variations occur. As such, once we neglect initial transients and consider only the seasonally periodic stochastic process, the temporal average of the ensemble over a year will be equivalent to the long-term temporal average; both will be designated by *x*(*t*) in grey, their associated pdf *p*(*x*,*t*) at six points in the year, their ensemble average 〈*x*(*t*)〉 as a bold line, and their long-term average as a dashed line.

The mean soil moisture balance corresponding to equation (2.1) can be written in its normalized form as
*γ*(*t*)=*w*_{0}/*α*(*t*).

When considering averages over a large area of heterogeneous soil and vegetation, evapotranspiration may be assumed to depend linearly on *x*, taking a value of 0 at *x*=0 to *x*=1 22]. Thus, the evapotranspiration term on the right-hand side of equation (2.3) can be simplified as
*t*. Equation (2.3), now approximating evapotranspiration, can be reduced to
*x*_{t}〉.

Because it will be used later, we note here that the pdf of *x* and its ensemble mean are already known under steady-state conditions for constant parameters, in the form of a truncated gamma distribution with shape parameter *a* and rate parameter *b* [22],
*Γ*(⋅,⋅) indicates a truncated gamma function [26].

### (b) Approximating the average leakage/run-off (LQ)

The complication of solving equation (2.5) under fully seasonal conditions arises from the nonlinear time dependency in each of the parameters, and our lack of information on the full pdf *p*_{t}(*x*) at every time point, which is required to quantify the leakage/run-off term, which in turn governs the evolution of mean soil moisture. Here, we present, in order of increasing complexity, four treatments to the LQ term in equation (2.6). These approximations reduce equation (2.5) into an ordinary differential equation (ODE) for 〈*x*_{t}〉.

We start by considering cases where 〈LQ_{t}〉 can be effectively neglected. Such an assumption is justified when the mean rainfall depth *α*_{t} is low relative to the soil root depth *w*_{0} or if rainfall is infrequent, corresponding to small *λ*_{t}. Such is often the case in extremely dry climates or during the dry season. Equation (2.5) subsequently becomes
*x*_{0} if all parameters are considered to be constant (e.g. within the same growing season), resulting in an analytical solution of

Indeed, the full solution to the stochastic differential equation (2.1), *p*_{t}(*x*), is available when constant parameters are assumed for this case [13]. However, because these assumptions do not place an upper bound on the value of 〈*x*_{t}〉, they also result in significant overestimation of 〈*x*_{t}〉 when there is leakage or run-off. Other implications for soil moisture and plant water stress are extensively discussed in reference [13] under constant climate conditions for the growing season. The next three cases are newer improvements on equation (2.5) that approximate 〈LQ_{t}〉 using various assumptions for *p*_{t}(*x*).

#### (i) Quasi-steady-state approximation

The first treatment approximates the instantaneous pdf of soil moisture by its quasi-steady-state pdf as determined by the corresponding environmental parameters,
*a*_{t}=*λ*_{t}/*k*_{t} and *b*_{t}=*γ*_{t}. That is, the parameters found at each instant *t*, i.e. *λ*_{t}, *k*_{t}, and *γ*_{t}, are ‘frozen‘ as constants and applied to a parallel homogeneous stochastic process until it has reached steady state, represented by *p*_{t,ss}(*x*). This is similar to assuming that the inhomogeneous process can instantaneously reach steady state at every point in time. Examples of *p*_{t,ss}(*x*) are shown in figure 2 at different times of the year.

Hence, we substitute *p*_{t,ss}(*x*) in place of *p*_{t}(*x*) in the formula for mean leakage, shown in equation (2.6) as _{t}〉 is approximated by the steady-state leakage value 〈LQ_{t}〉_{ss}, calculated using instantaneous environmental values at *t* and given by

Equation (2.5) can now be written by applying equation (2.7) for 〈*x*_{t}〉_{ss} as
_{t}′ is strictly positive. Additionally, equation (2.10) can be cast into an equivalent form of
*x*_{t}〉 and its quasi-steady-state value 〈*x*_{t}〉_{ss}. The approximation of *p*_{t}(*x*) by *p*_{t,ss}(*x*) is most appropriate when soil moisture can quickly adjust to changes in its environment, which is the case when soil water storage is small. Figure 2 also shows that the difference between true steady-state soil moisture and *p*_{t,ss}(*x*) diminishes during persistently wet or dry periods. During the dry-down and wetting-up periods, the quasi-steady-state approximation predicts soil moisture values that are, respectively, lower and higher than true steady-state values because it does not capture the effect of soil moisture transfer over time. Thus, this treatment of the leakage term is expected to work well at small *γ*_{t} values and in conditions where the environmental parameters do not change quickly. Because 〈*x*_{t}〉_{ss} can be calculated using environmental parameters only and do not depend on the state variable 〈*x*_{t}〉, this approximation results in an inhomogeneous, linear, time-dependent ODE that can be formally solved (see, for example, page 14 of [27]).

#### (ii) Negligible fluctuations approximation

The next treatment assumes that the soil moisture pdf itself is concentrated entirely on its mean value, represented by a Dirac delta function at 〈*x*_{t}〉, i.e.
*x*_{t}〉 in figure 2. The substitution of the Dirac delta function into the pdf in equation (2.5) results in 〈*e*^{−γ(1−xt)}〉=*e*^{−γ(1−〈xt〉)} in the leakage term, and correspondingly, the new soil moisture equation becomes

Because *e*^{−γ(1−x)} for any value of *x* is a convex function, by Jensen's inequality, 〈*e*^{−γ(1−xt)}〉≥*e*^{−γ(1−〈xt〉)}. While we thus know that the approximated mean leakage must be smaller or equal to the actual mean leakage at any given time point, the cumulative result of this estimation is more difficult to predict, because the evolution of 〈*x*_{t}〉 necessarily depends also on the history of its seasonal variations. The difference between 〈*e*^{−γ(1−xt)}〉 and *e*^{−γ(1−〈xt〉)} decreases as the soil moisture pdf becomes naturally concentrated around its mean; this happens as the true value of soil moisture approaches its upper or lower bounds (for example, as during the wet period in figure 2). Thus, we may expect this approximation to perform better in conditions that allow soil moisture to remain perennially close to its bounds, such as in extremely wet or arid environments where soil moisture is continuously saturated or otherwise dry. Unlike the previous approximations, equation (2.12) is a nonlinear ODE. While an analytical solution is unavailable, its closed form allows numerical solutions to be easily found using simple ODE solvers.

#### (iii) Self-consistent truncated gamma approximation

The last treatment approximates the instantaneous soil moisture pdf *p*_{t}(*x*) by a truncated gamma distribution with parameters consistent with the evolution of the mean soil moisture. The shape parameter of the truncated gamma distribution—which is the pdf of *x* under steady-state conditions (equation (2.7))—is constantly adjusted, so that the resulting pdf will have its mean centred on the current value of 〈*x*_{t}〉. This combines the advantages of the quasi-steady-state approximation to better capture the overall shape of *p*_{t}(*x*) with the negligible fluctuations approximation's ability to track its mean.

The mean of the truncated gamma distribution with given shape parameter *a* and rate parameter *b* can be calculated using equation (2.7). To match it to the mean soil moisture value at any given time, we set the instantaneous rate parameter to *b*_{t}=*γ*_{t} (similar to the quasi-steady-state model) and obtain an implicit function *a*_{t}=*f*^{−1}(〈*x*_{t}〉;*γ*_{t}) for the instantaneous shape parameter *a*_{t} in terms of 〈*x*_{t}〉 and *γ*_{t}, where
*b*_{t}=*γ*_{t} and *a*_{t}=*f*^{−1}(〈*x*_{t}〉;*γ*_{t}), where *f*^{−1} is the inverse of equation (2.13), has its mean equal to 〈*x*_{t}〉 (figure 2). The shape of the inverse function *f*^{−1} is illustrated in figure 3 for *γ*_{t}=5.5.

Substituting the formula for *f*^{−1}(〈*x*_{t}〉;*γ*_{t}) in place of *a*_{t} ensures that the ODE, as a function for 〈*x*_{t}〉, can be solved using simple numerical solvers in much the same way as equation (2.12). The additional complexity of this model comes only from having an implicit function embedded within the ODE. In fact, this ODE is similar to the ODE derived from the quasi-steady-state approximation (equation (2.10)), where the rainfall and leakage/run-off process can again be considered as a new censored rainfall process with an equivalent, strictly positive rainfall frequency, *a*_{t}=*f*^{−1}(〈*x*_{t}〉). We have chosen to solve for *a*_{t} as function of 〈*x*_{t}〉 while keeping *b*_{t}=*γ*_{t}, because the run-off generation mechanism acts as a censoring of the rainfall process; in a censored process, the size of the marks (e.g. rainfall depths) can be drawn from an exponential distribution with the same normalized mean of *γ*_{t} as the original process, owing to the memoryless property of exponential distributions [28]. Furthermore, the parameter *a*_{t} for a given *b*_{t} is directly proportional to the mode of the truncated gamma distribution, suggesting that the mode might be better used to summarize the effects of *p*_{t}(*x*) instead of its mean.

The difference between this approximation and the quasi-steady-state approximation is that the leakage/run-off term is no longer predetermined by environmental conditions but rather needs to be solved implicitly using current values of 〈*x*_{t}〉. This modification produces considerable improvements in the approximation to the true values (figure 2 and appendix A). While relatively more computationally costly because of the inversion required for *f*^{−1}, it still reduces the problem to that of solving a nonlinear ODE and eliminates the need to simulate the full random processes associated with stochastic rainfall realizations.

### (c) Budyko's formulation

In addition to the mean soil moisture evolution over a year, we are also interested in the partitioning of rainfall into evapotranspiration and leakage/run-off. A synthetic representation of the evapotranspiration ratio, defined as the ratio of evapotranspiration over rainfall, is provided through Budyko's curve as a function of the dryness index, defined as the ratio of potential evapotranspiration over rainfall [14]. While Budyko first introduced this framework under steady-state conditions using long term means, where

For example, the instantaneous partitioning terms using time-dependent parameters are

The annually averaged curve using seasonal parameters that change over time can be constructed by using the time average of the instantaneous dryness index and the evapotranspiration ratio, i.e.
*T*_{year} is the length of the year. While *t* is contained in the integral limits, the resulting quantifies do not depend on *t*, because the integrations are conducted over their natural period of variation, namely over a year. For that same reason, in the absence of inter-annual variabilities which is the case considered here,

## 3. Results

### (a) Comparison with stochastic simulations

The three new approximations of the LQ term described in the previous section—quasi-steady-state, negligible fluctuations and self-consistent truncated gamma— allow complete flexibility in the climatic inputs during modelling of seasonal soil moisture dynamics. Their results are compared here with numerical simulations using parameters typical of tropical dry and Mediterranean climates. We start each run with an arbitrary initial condition *x*_{0}=0.5. At each time step, changes in the effective soil moisture *x* is determined by inputs from rainfall and output from ET and LQ as determined by the time-dependent climate parameters. Rainfall is generated according to a non-homogeneous, marked Poisson process with mean frequency λ_{t} and mean intensity *γ*_{t}. We adopt sinusoidal forms to describe rainfall and potential evapotranspiration inputs λ_{t} and *k*_{t}, e.g.
*v*_{t} is a stand-in variable for λ_{t} or *k*_{t}, *μ*_{v} is the annual mean, *A*_{v} is the amplitude of seasonal variation, *ϕ*_{v} is the phase and the period *ω* is set to a year. If rainfall depth at a time exceeds soil water storage (i.e. 1−*x*, with 1 being the upper bound of effective soil moisture), soil moisture is filled to its upper bound at *x*=1 and the rest of rainfall is lost to leakage/run-off. Loss from evapotranspiration occurs as a linear function of *x*. The mean soil moisture 〈*x*_{t}〉 is calculated over many iterations, with the initial transient state discarded.

Results from a representative tropical dry and Mediterranean climates are shown in figure 4. For Mediterranean climates, specific climate parameters are set to (*μ*_{λ},*A*_{λ})=(0.3,0.2) and (*μ*_{k},*A*_{k})=(0.03,0.02), with the phase difference *ϕ*_{k}−*ϕ*_{λ}=180° to reflect a hot, dry summer and a cool, wet winter (top panel, figure 4). For the tropical dry climate, potential evapotranspiration is set to high year-round with mild fluctuations, e.g. (*μ*_{k},*A*_{k})=(0.06,0.02), with more dramatic changes in rainfall between the seasons, (*μ*_{λ},*A*_{λ})=(0.6,0.575). Both climate types are set to a constant soil storage index of *γ*_{t}=5.5, consistent with a globally averaged effective rooting depth [22].

The bottom panels of figure 4 show comparisons of the simulation and the model outputs for mean soil moisture 〈*x*_{t}〉, water partitioning components 〈ET_{t}〉 and 〈LQ_{t}〉, and an instantaneous Budyko-type partitioning, e.g. 〈ET_{t}〉/〈*R*_{t}〉 versus 〈*D*_{t}〉, for every time point of the year. As expected, the self-consistent truncated gamma approximation followed the trajectory of the simulations most closely for all terms of comparison. The other two models do well enough in capturing the time evolution of the mean leakage/run-off and mean evapotranspiration as well as the general arc of seasonal mean soil moisture evolution. For example in the Mediterranean climate, all three models were able to show the surge in evapotranspiration at the beginning and end of the summer season. Further analyses in appendix A show that when considering annually averaged values, all three models are also able to capture the annual evapotranspiration ratio

In addition, all three models also demonstrate hysteretic behaviour, wherein the same climate condition may not result in the same hydrological response owing to the effect of transient soil water storage. This can be seen in the ‘loop‘ of the transient Budyko curves, as well as a similar loop in the relationship between 〈*x*_{t}〉 and 〈LQ_{t}〉 (not shown). For any given dryness index 〈*D*_{t}〉, the corresponding 〈ET_{t}〉/〈*R*_{t}〉 value can fall within two domains—one along the upper ‘dry-down‘ trajectory and one along the lower ‘rewetting‘ trajectory. The reason more evapotranspiration can occur on average during the ‘dry-down‘ trajectory is that stored soil moisture can be carried over in greater amount from a wet period than from a dry period, accentuating the role of soil water storage.

### (b) Transient departures from Budyko's curve

The previously described hysteresis of mean soil moisture and evapotranspiration ratio is further explored in this section. Here, we expand the projected views in the bottom panels from figure 4 with an additional dimension for time and show that increasing the phase difference between rainfall and potential evapotranspiration will further increase the asymmetry in the curves between the wet and dry season.

Figure 5 shows expanded views of the transient Budyko's curves produced from the self-consistent truncated gamma approximation, parametrized by the same mean and amplitude for both rainfall and potential evapotranspiration inputs as in the Mediterranean (figure 5*a*) and tropical dry (figure 5*b*) climates shown previously in figure 4. In other words, the difference between each transient curves are attributed only to the phase difference between rainfall and potential evapotranspiration, with the thick black line showing a phase difference of 0° and the thick grey line showing 180°. For comparison, the classical steady-state Budyko's curve and the annually averaged Budyko's curves parametrized with seasonal climatic inputs are shown, respectively, as a thin black line and two dashed lines, black and grey, for a greater range of the dryness index. While the instantaneous dryness index and evapotranspiration ratio for the transient curves might fluctuate over a range of values over the year, each transient curve collapses to only a single point placed on the annually averaged curves.

Figure 5 shows that in both climate types, increase in the phase difference results in a larger, more asymmetrical hysteresis loop between the two seasons. Because this is due to the transfer of stored soil water from the wet season to the dry season, it implies that seasonal soil water storage becomes more important as rainfall and potential evapotranspiration become more out-of-phase. In Mediterranean climates, there is a sizable part of the year during which the instantaneous evapotranspiration ratio is above 1, where the high evapotranspiration is sustained not only from rainfall inputs but also from residual soil water carried over from the wet season. This ‘water-deficit‘ portion of the curve diminishes as the climatic inputs become more in-phase, when as a consequence evapotranspiration occurs mostly when the rainfall input is also high. In tropical dry climates, even though the climatic inputs are in-phase, there is still a relatively large excursion of the instantaneous evapotranspiration ratio above 1, because rainfall is so low during the dry season that even a small amount of stored water from the wet season will allow evapotranspiration to surpass the rainfall input during the dry season.

### (c) Annual Budyko's curve

The annually averaged Budyko's curves calculated from seasonal climatic inputs are shown in figure 6 as a function of seasonal rainfall amplitudes, with the rainfall signal both in-phase (black) and out-of-phase (grey) with potential evapotranspiration, with a phase difference of 0° and 180°, respectively. The evapotranspiration ratios *x*_{t}〉 have reached periodic steady state over the year. The particular curves on which the Mediterranean and tropical dry climates shown in figure 4 reside are highlighted as dashed curves. In fact, these dashed curves are exactly the same as the dashed seasonal curves shown in figure 5 of the same colour. The reason that the graphs do not extend to all _{t} at any given point cannot fall below 0.

It can be seen that the annually averaged evapotranspiration ratio is lower for an out-of-phase climate than an in-phase climate, regardless of seasonal rainfall amplitude. In general, increasing seasonal rainfall amplitude decreases the annual evapotranspiration ratio. However, in some cases (left panel, at high

## 4. Discussion and conclusion

To facilitate the analyses of hydrological processes under seasonally varying climates, both at the annual and intra-annual timescales, we introduced three new approximations to a stochastic soil moisture model which were used to construct relatively parsimonious ODEs that closely described the seasonal trajectories of mean soil moisture and its associated pdfs. This links daily stochastic soil moisture dynamics to its seasonal patterns in a direct and consistent way. It is also a considerable advantage over having to simulate stochastically rainfall input at every time step and taking the ensemble mean over many realizations, and moreover is amenable to analytical formulations.

The resulting mean soil moisture trajectories are used to investigate the effect of transient soil water storage in seasonal climates. The mean soil moisture and instantaneous evapotranspiration ratio both exhibit hysteresis because of seasonal soil moisture storage, whereby evapotranspiration occurs during the dry season using soil moisture carried over from the wet season as supplements to minimal rainfall inputs. This results in a significant portion of the year during which the instantaneous evapotranspiration ratio stays above 1, which can be further prolonged as rainfall and potential evapotranspiration become more out-of-phase. In some dry environments, the annual evapotranspiration ratio can increase under seasonal climates precisely, because soil water storage allows a more ‘efficient‘ mode of evapotranspiration by transferring soil moisture between seasons.

The main results brought out by our models are consistent with those from other empirical studies. The importance of monthly soil moisture carryover has previously been pointed out by Jothityangkoon & Sivapalan [29], who, by incorporating seasonality and soil water transfer between storm events, were able to achieve a much better match between observed and predicted interannual variability of water balance for selected catchments in the United States, Australia and New Zealand. Likewise, inclusion of a monthly carryover factor in Gerrits *et al.* [30] significantly improved predictions of annual water balance in semiarid areas in Africa. While these studies highlight the role of seasonal soil water storage on annually integrated quantities, we have explicitly related soil water storage to their seasonal, hysteretic dynamics. Petrie & Brunsell [31] used the seasonality of potential evapotranspiration to induce such hysteresis in the presence of soil water storage, but we have extended the analysis here to account also for the complex interplay between seasonal magnitude and timing of rainfall and potential evapotranspiration (via their amplitude and phase difference). Our models also predict that increase in the phase shift between atmospheric demand and water supply (towards more Mediterranean-type climates) decreases the annual evapotranspiration ratio. This has been shown through first-order modelling exercises [1,16] and has recently been corroborated with observations across a global network of flux towers [32].

While our approach captures well the hydrological influence of seasonal climate variability, the simultaneous contributions from seasonal variation in soil water storage (as a result of between-season changes in rooting depths or vegetation cover), groundwater and interannual variability in the seasons, which are all known to affect annual soil water partitioning, remain to be explored. Vico *et al.* [4] recently investigated plant water and carbon uptake based on different physiological and phenological strategies in seasonally dry ecosystems using a coupled model of soil water (including deeper, more persistent storage) and plant carbon balances. The versatility of the models introduced here can allow for similar future extensions to account for the role of groundwater dynamics and vegetation water uptake strategies on seasonal and annual water partitioning.

## Funding statement

X.F. acknowledges funding from the NSF Graduate Research Fellowship Programme. A.P. acknowledges NSF grant nos. CBET 1033467, EAR 1331846, EAR 1316258, FESD1338694, as well as the US DOE through the Office of Biological and Environmental Research, Terrestrial Carbon Processes programme (de-sc0006967), the Agriculture and Food Research Initiative from the USDA National Institute of Food and Agriculture (2011-67003-30222).

## Acknowledgements

The authors gratefully acknowledge helpful discussions and advice from Valerie Isham and David R. Cox. The authors also thank Tal Svoray, Gabriel G. Katul and Shmuel Assouline, for organizing the workshop on Ecohydrology of Semiarid Environments: Confronting Mathematical Models with Ecosystem Complexity, held at Ben-Gurion University of Negev in Beer-Sheva, Israel in May 2013, which stimulated many interesting and valuable discussions, as well as helpful suggestions made by two anonymous reviewers.

## Appendix A. Model performances

The accuracy of the three models introduced in this study is checked against simulation results under more general climate conditions. Figure 7 shows differences in the annual evapotranspiration ratio *μ*_{λ},*A*_{λ})=(0.2,0.1); (ii) seasonal, (*μ*_{λ},*A*_{λ})=(0.5,0.5) and (iii) wet, (*μ*_{λ},*A*_{λ})=(0.9,0.1), with three constant soil storage indices *γ*_{t}=3,5.5,30, and five phase difference regimes between rainfall and potential evapotranspiration (while potential evapotranspiration itself remains the same at (*μ*_{k},*A*_{k})=(0.03,0.01) for all cases). The three clustered bars represent results from each model, with the leftmost for the quasi-steady-state approximation, the middle for the negligible fluctuations approximation and the rightmost for the self-consistent truncated gamma approximation. As can be seen, all models are able to capture simulation results for *γ*_{t}=5.5, which is adopted for our tropical dry and Mediterranean examples, the quasi-steady-state model can capture

We also evaluate the ability of the models to follow intra-annual variations in mean soil moisture, which influences their ability to reproduce the transient trajectories of figure 5. The modelled values are evaluated using their pattern correlation (PC) and centred RMSE over an annual cycle when compared with simulated values. These measures are related to each other through *σ*_{f} and *σ*_{r} are the standard deviations, respectively, of the modelled and simulated values. This relationship can be concisely summarized in the Taylor diagrams [33] of figure 8. The triangle, diamond and circular markers now represent, respectively, the quasi-steady-state model, the negligible fluctuations model and the self-consistent truncated gamma model. Each is subjected to the same climate and soil conditions as in figure 7, except the phase difference between rainfall and potential evapotranspiration is restricted to 0° (in-phase) and 180° (out-of-phase).

We find in figure 8 that the self-consistent truncated gamma model vastly outperforms the other two models, especially in highly seasonal climates. In the Taylor diagrams at the top, the self-consistent truncated gamma model markers (circles) have the shortest distance to the simulation marker (star), which indicates the lowest values of RMSE. Here, all values are normalized with respect to the standard deviation of the simulation [33]. The bottom panels show the RMSE for all models using non-normalized results; thus, the values represent how much modelled mean soil moisture deviates from the simulated results in absolute terms. The self-consistent truncated gamma model shows a remarkable level of fidelity to simulation results under all conditions considered. The RMSEs for the other two models tend to grow with increasing soil depth, with especially large errors in highly seasonal climates attributable to their tendency to exaggerate the variations in mean soil moisture over each season, resulting in a larger standard deviation compared with simulations. These plots collectively show that we can be confident in our use of the self-consistent truncated gamma model to produce transient trajectories in figure 5.

## Appendix B. Nonlinear effects of climate seasonality

Failure to consider seasonal variations in the climate can lead to significant underestimation of the proportion of the annual rainfall lost through leakage/run-off, and thus an overestimation in the annual evapotranspiration ratio. We illustrate this point by comparing

Figure 9 considers the effects of nonlinearity in the soil moisture model resulting from explicitly considering seasonal rainfall parameters. The standard sinusoidal form (3.1) is used for both rainfall and potential evapotranspiration, with (*μ*_{k},*A*_{k})=(0.03,0.01). The amplitude of rainfall seasonality, *A*_{λ}, varies along with the phase shift Δ*ϕ*=*ϕ*_{k}−*ϕ*_{λ}, which changes from 0° (in-phase) to 180° (out-of-phase). The grey and black meshes show the results, respectively, for a deeper soil of *γ*_{t}=30 and a typical soil depth of *γ*_{t}=5.5. It is worthwhile to stress that the results in each panel are derived from the same mean climatic parameters.

Errors in *μ*_{λ}=0.3, errors in *γ*_{t}=5.5. This is probably caused by the generation of a substantial amount of leakage through to a relative lack of soil water storage. On the other hand, for high rainfall of *μ*_{λ}=0.7, increasing rainfall amplitude *A*_{λ} may actually decrease the calculated error, especially for deeper soils when the climate inputs are in-phase. However, this trend quickly reverses, with the errors increasing again once rainfall amplitude increases past a certain level.

Altogether, these observations show that assessing hydrological variables containing inherently nonlinear behaviour, such as the case with

- Received August 17, 2014.
- Accepted November 18, 2014.

- © 2014 The Author(s) Published by the Royal Society. All rights reserved.