Complexity and coexistence in a simple spatial model for arid savanna ecosystems

Tree–grass coexistence is broadly observed in tropical savannas. Recent studies indicate that, in arid savannas, such coexistence is stable and related to water availability. The role of different factors (from niche separation to demographic structure) has been explored. Nevertheless, spatial mechanisms of water–vegetation interactions have been rarely taken into account, despite their well-known importance for vegetation distribution. Here, we introduce a spatial model including tree and grass biomass dynamics, together with soil and surface water dynamics. We consider two water–vegetation feedbacks. Grasses increase water infiltration into the soil, while tree shadow limits evaporation, and both mechanisms increase soil water availability, leading to positive feedbacks. The infiltration feedback can also lead to spatial pattern formation. Despite the fact that trees and grasses compete for the same resource, namely water, we observe stable coexistence as a possible model outcome. The system displays a complex behavior, with multiple stable states and possible catastrophic shifts between states, e.g., patterned grassland, bare soil and forest. In our model, coexistence is always linked with multi-stability and spatial pattern formation, driven by grass infiltration feedback. Given such complex model solutions, we expect that, under real conditions, heterogeneities and disturbances, acting on the multi-stable states, may further foster coexistence.


Introduction
Herbaceous and woody plant coexistence is observed in many terrestrial ecosystems. Tropical savannas are the most common example. Savannas occupy about 20 % of the Earth land surface and they are observed in a large range of mean annual precipitation (e.g., in Africa, approximately between 100 and 1,500 mm, Lehmann et al. 2011). They are broadly defined as systems where woody and herbaceous vegetation coexist, competing mainly for the same resource, namely soil water. In the following we will use the word "tree" to indicate any kind of woody vegetation (including shrubs), and the word "grass" to indicate any grass or herbaceous species (House et al. 2003). Although widely debated in literature, the mechanisms behind their coexistence are not yet clear (e.g., Sankaran et al. 2004;Scholes and Walker 1993;Scholes and Archer 1997). Tree-grass coexistence is usually considered (sensu Chesson 2000) to be either stable, or unstable but stabilized by external disturbances (such as fires, grazing etc.). Most of the models explaining savannas as deterministic stable coexistence assume different root depths between trees and grasses, thus separating tree and grass niches (e.g., van Langevelde et al. 2003;Walker and Noy-Meir 1982;Walter 1971). Nevertheless, it has been repeatedly pointed out that this separation is not observed in all savannas (Higgins et al. 2000;House et al. 2003;Sankaran et al. 2004;Scholes and Archer 1997). Only few models can predict tree/grass coexistence without including resource separation, explaining coexistence as the deterministic outcome of internal dynamics (Baudena et al. 2010;Scheiter and Higgins 2007). Most of the other studies represent savannas as an unstable intermediate biome between forest and grassland, where coexistence of trees and grasses is observed because of disturbances or heterogeneities, such as fires, grazing, rainfall variability etc. (see, e.g., Beckage et al. 2009;D'Odorico et al. 2006;Gardner 2006;Higgins et al. 2000;Jeltsch et al. 2000). Recent developments suggest that the two theories, considering savannas as either deterministically stable or unstable equilibria maintained by disturbances, do not necessarily collide. In other words, depending on plant types and environmental conditions, both theories can be valid, and disturbances may simply enlarge the range where savannas are observed (Baudena et al. 2010;Higgins et al. 2010;Sankaran et al. 2005). In particular, in the driest savanna ecosystems, an upper limit to potential tree cover is observed, increasing with annual rainfall (Sankaran et al. 2005), indicating that tree-grass coexistence in this range is possibly stable and resource-driven. In these dry savannas, fires and herbivores still play an important role, but they are not necessary to maintain them (as shown also with fire exclusion experiments, e.g. Higgins et al. 2007). In the most humid ecosystems, potential tree cover is not resource limited, and thus fires prevent woody vegetation to close up, leading to disturbance-driven savannas (Hirota et al. 2011;Sankaran et al. 2005;Staver et al. 2011a, b).
Spatial mechanisms, including both vegetation and resource dynamics, have been rarely taken into account in the study of tree-grass dynamics. Mostly, space has been introduced as a source of heterogeneity, which maintains an otherwise unstable tree-grass coexistence (Jeltsch et al. 1998;Scholes and Walker 1993). Baudena et al. (2010) include space, describing the dynamics of vegetation cover, but the representation is implicit and no resource dynamics is taken into account.
Scale-dependent water-vegetation feedbacks have been shown to be of large importance in arid and semi-arid ecosystems, where water is the only limiting resource. They can determine spatial vegetation pattern formation, arising from the system dynamics without any underlying heterogeneity (e.g., Klausmeier 1999;von Hardenberg et al. 2001;Rietkerk et al. 2002;Sherratt 2005;Esteban and Fairen 2006;Meron 2011). These resource-concentration mechanisms, leading to short scale positive feedback and long range negative feedbacks, are related to the existence of bistability and catastrophic shifts in ecosystems (Rietkerk et al. 2004;Rietkerk and van de Koppel 2008). Considering the models that include explicit spatial representations, Pueyo et al. (2010) explore the role of different plant reproductive traits on the interplay between competition and facilitation among plants, but they do not investigate the possible outcomes of competition. Gilad et al. (2007) examine the interplay of different water-vegetation feedbacks in woody-herbaceous communities, evaluating how species interactions can vary from facilitation to competition depending on the relative feedback strength. Their model shows that stable coexistence of the two plant types, competing for the same water resource, is possible because of the role of spatial feedbacks. Nevertheless, their approach is not only resource-driven, because they represent biomass growth with a logistic term, assuming implicitly that other external limitations to plant growth exist, besides water, possibly explaining their finding of stable coexistence.
Given the importance of spatial resource-driven vegetation dynamics, here we study their effects for the possible outcomes of tree-grass competition for soil water in arid savannas. We analyze under which conditions spatial mechanisms may foster the stable coexistence of grasses and trees (mostly shrubs in the driest ecosystems). We introduce a modified version of the model by Pueyo et al. (2010), where grasses reproduce by vegetative propagation whereas trees disperse seeds, and we include some elements from the model of Gilad et al. (2007). The model represents water dynamics in space, and vegetation dynamics is dependent only on soil water content. Gilad et al. (2007) and Pueyo et al. (2010) assume that both plant types can increase water infiltration, and in Gilad et al. (2007) both plant types shade the soil as well. However, here, we consider that the two plant types have different water-vegetation feedbacks. Only grasses increase water infiltration into the soil (e.g. Rietkerk and van de Koppel 1997;Walker et al. 1981), while trees instead shade the soil and thus limit evaporation (see e.g. Baudena and Provenzale 2008;Gilad et al. 2007). The infiltration feedback can lead to vegetation pattern formation (see also Rietkerk et al. 2002). The tree shading feedback is positive (the larger the tree biomass, the larger the evaporation reduction), but it does not enhance water redistribution from bare to vegetated sites, and thus it does not lead to vegetation pattern formation (see also Gilad et al. 2007).
In the following, we analyze the role of these feedbacks on the possible outcomes of tree-grass competition in waterlimited ecosystems, examining the active role of spatial dynamics. We also explore the relation between pattern formation, tree-grass coexistence and alternative stable states.

Material and methods
We introduce a model with four partial differential equations (PDEs), which is a modified version of a spatially explicit integro-differential model introduced by Pueyo et al. (2010), including some elements from Gilad et al. (2007) (see below for details). The model describes the dynamics of tree density, b T (g m −2 ), grass density, b G (g m −2 ), soil water content, w, (mm) and surface water, h (mm): where t is time, r 2 ¼ @ 2 =@x 2 þ @ 2 =@y 2 , and (x, y) are space coordinates. Equation (1) represents the dynamics of tree biomass. Tree biomass grows linearly with the product of soil water and tree biomass itself (first term on the r.h.s.), c T is the conversion of water uptake into tree growth and γ T is the maximum specific tree water uptake. Trees decrease with a linear mortality (−m T b T ). In the model of Pueyo et al. (2010), the growth of both plant types with respect to water was instead represented by a Holling-type II functional response. Here, we simplify it to a linear term, a valid assumption when water availability is low. The third term on the r.h.s. represents the growth of tree biomass due to seed dispersal. S is an integral term that takes into account local seed production and seed dispersal: where Ω is the total domain, L is the mean dispersal distance and x ! À x ! 0 is the distance from the source x ! .
The seed dispersal kernel, exp À2 x only on the distance from the parent plant, i.e. seed dispersal is homogeneous (it does not depend on the plant position) and isotropic (it does not vary with direction; Mistro et al. 2005). The seed germination probability and growth depend linearly on soil water content w, while η is a parameter representing seed production (fecundity) and recruitment per water unit (see also Pueyo et al. 2008, for more details on seed dispersal modeling).
Grass biomass dynamics are represented in Eq. (2). As for the trees, grasses grow linearly with the product of soil water content and grass biomass, and they have a linear mortality term (respectively, first and second term on the r.h.s.). The coefficient c G is the conversion of water uptake into grass growth, while γ G is the grass maximum specific water uptake. Grasses have local vegetative propagation, represented by a diffusion term (r 2 b 2 G ) in the equation, where δ G is the diffusion coefficient of grass (HilleRisLambers et al. 2001).
In Eq.
(3), infiltration is the soil water input (first term on r.h.s.). Infiltration grows linearly with surface water h, and it grows with grass biomass saturating at large grass biomass values. In this term, α is the maximum infiltration rate, w 0 is water infiltration rate in the absence of grass plants and k is the saturation constant of water infiltration. The second term on the r.h.s. represents plant water uptake, which is linearly dependent on soil water content. Water uptake losses increase with tree and grass biomass, with γ T and γ G representing the maximum specific water uptake rate of trees and grasses respectively. Evaporation is the third term on the r.h.s., and it is linearly dependent on soil water content. As in Gilad et al. (2007), evaporation is reduced under trees because of tree shadow (σ is the coefficient of evaporation reduction, and b s is a biomass normalization coefficient). Soil water spreads within the soil with a diffusion constant δ w , see last term on the r.h.s. (Rietkerk et al. 2002). In the model of Pueyo et al. (2010), both water uptake and evaporation increase with soil water content saturating at high values. Here, we use a linear function, a common assumption in hydrology, especially in drier climates (see, e.g., Kim et al. 1996).
In Eq. (4), rainfall r represents the input of surface water, which diminishes when water infiltrate into the soil (second term on the r.h.s., same term as in Eq. 3). Surface water redistribution through surface runoff is modeled with the term δ h ∇ 2 h 2 , as in e.g. Gilad et al. (2007). This form can be obtained by describing runoff water as a shallow fluid layer on a flat surface (for the derivation, see Gilad et al. 2004, where δ h ∇ 2 h 2 is the only term left in case of flat surface). δ h is the coefficient of aboveground water redistribution.
We consider separately the effect of including infiltration feedback only (setting σ00 in Eq. 3), and the effect of adding the tree shading feedback (Eq. 3, σ01).
The model includes 18 parameters. Only 13 of them are free, while four of them are necessary to scale time, tree biomass, grass biomass, soil and surface water height (which both scale with the same parameter). The fifth superfluous parameter is σ, which we introduce only to turn on and off the shading feedback.
First, we analyze a non-spatial version of the model, excluding the diffusion terms in the equations of soil water, surface water, and grasses. In this case, infiltration simply becomes equal to rainfall. The integral S(b T ) in the seed dispersal term reduces to 2 pL 2 b T (the seeds are only produced locally). The equilibrium points of the non-spatial model are calculated analytically (see Appendix S1). To calculate the equilibrium states, we equal the time derivative of each variable to zero, i.e. db T If the equilibria are stable, the system returns to them after a small perturbation, whereas the system leaves an unstable equilibrium as soon as a small perturbation occurs. We calculate analytically the so-called "linear-stability" of the equilibrium points, which determines their stability in a small neighborhood of the equilibria. For this purpose, we calculate the eigenvalues of the Jacobian matrix associated to the system for each equilibrium point, which turns out to be stable if all of the real parts of the eigenvalues are negative (see Appendix S1 for the expression of the Jacobian matrix). The analysis of the non-spatial system equilibria helps understanding the role of the spatial mechanisms.
We solve the spatial model numerically, performing calculations along a grid of 64×64 cells with no-flow boundary conditions. The results we illustrate in the following do not vary substantially using different boundary conditions (e.g., periodic). Cell size is 2×2 m, and the time step is 0.002 day. The parameters chosen are reported in Table 1. These values are taken from earlier studies (Kletter et al. 2009;Pueyo et al. 2010;Rietkerk et al. 2002).
We start the simulations with homogeneous surface and soil water conditions (w t¼0 ¼ r=ρ e ; h t¼0 ¼ r=aw 0 ). Tree and grass biomass are randomly distributed on the plot with a large number of different configurations. For simplicity, we show here only four sets of initial conditions, but the final equilibria we find are the only qualitatively possible system configurations (as we have tested in a number of cases). The four sets of initial conditions used are: high tree and high grass (HTHG), with both trees and grasses covering randomly 30 % each of the surface (with initial density of 15 gm −2 for each pixel); high tree and low grass (HTLG), with trees covering 60 % of the plot and grasses covering 1 % at the initial stage (initial tree density, 15 gm −2 , initial grass density, 0.15 gm −2 ); high grass and low tree (LTHG), with grasses covering 60 % of the plot and trees covering 1 % at the initial stage (initial grass density, 15 gm −2 , initial tree density, 0.15 gm −2 ); low grass and low tree (LGLT), with grasses and trees covering 1 % each (initial density, 0.15 gm −2 ).
The results obtained are qualitatively robust and do not depend on fine-tuning of any particular parameter, but can be observed in general varying any of them. For illustration purposes, we choose to show the results in the parameter plane of two particular parameters (water/tree conversion rate c T , and rainfall r), but we have investigated also other parameter planes finding analogous results (not shown).
The analytical calculation of the value and linear stability of the equilibrium points of the non-spatial model, as well as the numerical integration of the full spatial PDEs, were implemented with a MATLAB 7.4 code.

Results
The non-spatial model displays three possible stable states: bare soil, homogeneous grass cover, and homogeneous tree cover. First, we consider only the grass infiltration feedback, Table 1 Model parameters (symbol, description, and value) These values are motivated by earlier studies (Gilad et al. 2007;Kletter et al. 2009;Pueyo et al. 2010;Rietkerk et al. 2002) Parameter Description Value Units c T Conversion of water uptake into tree growth 0-10 g m −2 mm −1 c G Conversion of water uptake into grass growth 7 g m −2 mm −1 γ T Maximum specific tree water uptake 0.1 m 2 g −1 day −1 γ G Maximum specific grass water uptake 0.1 m 2 g −1 day −1 m T Tree mortality rate 0. with trees not limiting evaporation through shading (σ00 in Eq. 3). The stable equilibrium points vary with changing values of water/tree conversion rate c T , and rainfall r (Fig. 1a). The other parameters are kept fixed (see Table 1 for parameter values). In this case, when rainfall is low, both plant functional groups do not survive. When rainfall is higher, if water-tree conversion rate c T is smaller than a certain value (c 1 ), only grasses survives, whereas above c 1 only trees are found. When we add shading feedback (σ01 in Eq. 3), a similar behavior is observed (Fig. 1b): bare soil at low rainfall; homogeneous grass cover at rainfall larger than a certain value r′ and c T lower than c 2 ; homogeneous tree cover for c T >c 2 at large enough rainfall. In addition, we observe bistability of homogeneous tree cover and bare soil in a region of parameter space with c T >c 2 and an intermediate rainfall rate (i.e., between the regions where the two states are separately stable). For such parameter values, if the initial conditions of tree cover are sufficiently high, the system displays a uniform tree cover, otherwise it ends up in a nonvegetated state.
Considering the full spatial model, we represent in Fig. 2 the different possible equilibria the system displays, as they vary with the water/tree conversion parameter c T and the rainfall rate r. This figure represents qualitatively where the different equilibria appear in parameter space. First, we analyze the case with infiltration feedback only (σ 00, Fig. 2a). At low c T values, the system at equilibrium displays bare soil at low rainfall and a uniform grass layer at high rainfall (as in the homogeneous model). For intermediate rainfall values, grass patterns are observed, and at lower rainfall values in their range of existence, they are bistable with bare soil. At high values of c T , the system goes from a bare soil state at low rainfall to a homogeneous tree-covered state at large rainfall, displaying the same behavior of the non-spatial model. At intermediate values of c T (slightly larger than c 1 ), the behavior of the system becomes more complex. At low rainfall, the system has only a bare soil equilibrium state. Increasing rainfall, we observe first bistability of bare soil and grass patterns (this case has parameter values in correspondence of part of the bare soil area in the non-spatial model, see Fig. 1a), then bistability of grass patterns with homogeneous tree cover (region γ in Fig. 2a), and for even larger rainfall (r>~0.9 mm day −1 ), we observe coexistence of trees and grasses in spatial patterns, bistable with homogeneous tree cover (region δ in Fig. 2a). The last two cases have parameter values corresponding to homogeneous tree cover in the non-spatial case (Fig. 1a).
We study also the spatial behavior of the system including the shading feedback (σ01 in Eq. 3), analyzing how the stable states varied when we change the values of c T and r (Fig. 2b). The system behavior is similar to what observed in the case with only infiltration feedback (Fig. 2a). Nevertheless, the system including both feedbacks displays a more complex behavior than the system including only infiltration feedback, showing even three states that are stable for the same parameter values (tri-stability). As c T gets just slightly larger than c 2 , we find as stable states, increasing rainfall: bare soil only; grass patterns bistable with bare soil (region α); tri-stability of homogeneous tree cover, grass patterns and bare soil (region ε); grass patterns bistable with homogeneous tree cover (region γ); tristability of homogeneous tree cover, grass patterns and tree-grass coexistence patterns (region φ); bistability of tree-grass coexistence patterns with homogeneous tree cover (region δ). Another difference with the case with only infiltration feedback is observed at large c T values, where at intermediate rainfall values both bare soil and homogeneous tree cover are stable (as in the non-spatial case, Fig. 1b).
We now study examples from the bistability areas in the spatial case with only infiltration feedback. In Fig. 3 Fig. 1 Non-spatial model, stable states in the parameter space of rainfall r and conversion of water uptake to tree growth, c T . a Only infiltration feedback is considered; b both infiltration and shading feedbacks are introduced. White areas represent bare soil; light gray, homogeneous grass cover; black, homogeneous tree cover; dark gray, bistability of bare soil and homogeneous tree cover. See Table 1 Table 1, with c T 06.1 gm −2 mm −1 ). For rainfall r00.6 mm day −1 (point I in Fig. 2a), if the initial conditions of both trees and grasses were high (HTHG), the system ends up in a grass pattern state, without any tree (Fig. 3a, top panels), while if the initial conditions of trees are high and grasses are low (HTLG) the system ends in a homogeneous treecovered state (Fig. 3a, bottom panels). When the system with the same parameter values gets larger rainfall (r0 0.9 mm day −1 , Fig. 2a, point II), we observe a state where trees and grasses coexist, both forming the same spatial pattern. The system reaches a state where trees and grasses coexist, displaying spatial patterns, when the initial conditions are high for both vegetation types (HTHG, top panels in Fig. 3b), while, if the initial conditions of grasses are low (HTLG), the system reaches the homogeneous state covered by tree only (bottom panels in Fig. 3b).

we show two examples corresponding to points I and II in
In Fig. 4, we represent biomass and soil water along a transect, evolving in time at four different time shots. The first panel (a) shows the initialization with random tree and grass vegetation cover, and constant soil water. After 1 day we notice that the peaks of water follow the grass distribution (b). Grasses evolve over time and the soil water content keeps correlating with them (c), so that the trees are forced to grow where the grasses are (d).
When the spatial system includes shading feedback, tristability is observed, as we show in two examples ( Fig. 5; parameters correspond to Table 1, and c T 06.1 gm −2 mm −1 ). When rainfall r is equal to 0.55 mm day −1 (Fig. 2b, point I), high grass and tree cover (HTHG) lead to grass only with patterns (Fig. 5a, top panels), high trees and low grass cover (HTLG) lead to homogeneous tree-only cover (Fig. 5a, middle panels), while low grass and tree cover (LTLG) lead to bare soil equilibrium (Fig. 5a, bottom panels). Figure 5b shows the other tri-stability observed, corresponding to r0 0.75 mm day −1 (Fig. 2b, point II), between coexistence patterns (starting from HTHG initial conditions, top panels), grass patterns (starting with LTHG initial conditions, middle panels) and homogeneous tree cover (starting from HTLG initial conditions, bottom panels). Notice that the parameter set giving rise to such model outcome is not correspondent to a bistable area in the non-spatial model.

Discussion
We studied the effect of spatial dynamics and water-vegetation feedbacks on tree-grass coexistence in arid savannas. Previous works suggest that, in savannas with low water availability, tree-grass coexistence is a stable equilibrium, due to internal ecosystem dynamics. This coexistence has been explained introducing different mechanisms, but usually the role of water-vegetation feedbacks in space is not Fig. 2 Spatial model, conceptual scheme of the stable states in the parameter space of rainfall r and conversion of water uptake to tree growth, c T . a Only infiltration feedback is considered; b both infiltration and shading feedbacks are introduced. White areas represent bare soil; light gray, homogeneous grass cover; black, homogeneous tree cover; dark gray, bistability of bare soil and homogeneous tree cover; yellow (!), bistability of bare soil and grass pattern; orange ("), grass patterns; green (+), bistability of grass patterns and homogeneous trees; light blue (%), bistability of tree-grass coexistence patterns and homogeneous trees; red (ε), tri-stability of bare soil, homogeneous tree cover and grass patterns; dark blue (φ), tri-stability of homogeneous tree cover, grass patterns and tree-grass coexistence patterns. Continuous lines in the two plots represent the limits of the stability areas in the non-spatial case (Fig. 1). The white dots in panels a and b indicate the points in parameter space chosen for the runs represented in Figs. 3 and 5, respectively. See Table 1 for other parameter values included. Gilad et al. (2007) underline for the first time the role of spatial water-vegetation feedbacks in leading to tree-grass coexistence, although their approach is not totally resource-driven.
In the model we presented, vegetation dynamics depended only on water. Grasses had the ability to increase water infiltration into the soil, and this water-vegetation feedback could induce spatial vegetation instabilities, leading to pattern formation. Stable tree-grass coexistence was one of the possible outcomes of the internal dynamics. Interestingly, coexistence was always observed in connection with spatial pattern formation and multi-stability. At present, we cannot state whether multi-stability is necessary to observe coexistence in our model. Nevertheless, coexistence always occurred within patterns. Grass-water infiltration feedback was the mechanism fostering it. Grasses could self-organize in space, forming spatial patterns that were reflected into soil water distribution. Tree dynamics were influenced by such soil patterns, and when rainfall was abundant enough, trees could coexist with grasses, in the same spatial pattern. Without infiltration feedback (as studied with the non-spatial model), the system would only have homogeneous solutions, and no coexistence would be observed. At low resource (rainfall) level, the infiltration feedback led to grass surviving in patterns and outcompeting trees, for parameter values where grasses would not survive if they did not increase infiltration. At higher resource levels, some of the water accumulated by grasses below them became available for trees as well, and thus trees could persist together with grasses in their same patterns. This model outcome is in line with observations reporting that increased water infiltration due to grass presence is an important mechanism for savanna dynamics (Walker et al. 1981). Facilitation of trees by grasses, especially at the stage of seedlings establishment, has been observed in arid environments (see e.g. Anthelme and Michalet 2009;Cipriotti and Aguiar 2010;Maestre et al. 2003), and it is explained mainly with the increased water availability due to grass infiltration enhancement. Besides, in the model, tree shading feedback increased water availability below tree canopies, and thus allowed for tree coexistence with grasses at lower rainfall level than without tree Fig. 3 Model with only infiltration feedback. Different equilibria are reached starting from different initial conditions. The plots display the spatial distribution of biomass density (increasing from yellow to green, with yellow corresponding to bare soil). a Case with rainfall r00.6 mm day −1 , corresponding to point I in Fig. 2a. Top panels: starting with high tree and grass initial conditions (HTHG, top left panels, t00 years), after 10 years the system displays grass patterns (top right panels); initializing the system with high tree cover and low grass cover (HTLG, bottom left panels), the system reaches an equilibrium where a homogeneous tree layer is present (t010 years, bottom right panels). b Case with rainfall r00.9 mm day −1 , corresponding to point II in Fig. 2a. Top panels: starting with high tree and grass initial conditions (HTHG, top left panels, t00 years), after 10 years the system displays tree-grass coexistence patterns (top right panels); initializing the system with high tree cover and low grass cover (HTLG, bottom left panels), the system reaches an equilibrium where a homogeneous tree layer is present (t010 years, bottom right panels). In every run c T is set to 6.1 gm −2 mm −1 , see Table 1 for other parameter  values shading. Tree shading feedback did not give rise to spatial effects, because soil water diffusion counteracted the shading effect, decreasing the gradient between vegetated and bare soil. Therefore, this feedback did not lead to pattern formation and, alone, could not drive coexistence (not shown).
Despite its relative simplicity, the system we introduced displayed rather complex features. Infiltration feedbacks led to spatial self-organization. As observed in the non-spatial version of the model, shading feedback added complexity to the system, since homogeneous tree cover became bistable with bare soil because of the shading feedback. Complexity was even more evident when we observed the combined effect of both feedbacks in space, where multi-stability was observed. The concept of alternative stable states has been pointed out to be very important for ecosystem dynamics in general (e.g., Rietkerk et al. 2004;Scheffer et al. 2009), and for understanding savanna dynamics (e.g., van Langevelde et al. 2003;Staver et al. 2011a;Walker et al. 1981). Many critical thresholds were present in our system, where "catastrophic shifts" could occur, changing dramatically the ecosystems, e.g., going from a state with grass only, in patterns, to a homogeneous state with trees, or to a bare soil.
Coexistence of shrubs and grasses within the same patterns has been observed in arid savannas (e.g., Cipriotti and Aguiar 2010;Meyer et al. 2008), as well as in other arid ecosystems (e.g., Maestre et al. 2003;Pugnaire and Luque 2001;Shachak et al. 2008), usually classified as drylands or steppe, rather than savannas.
The infiltration feedback led to stable coexistence of woody and herbaceous vegetation in arid savannas. This did not rule out the role of disturbances. Because of the complexity of the model solutions, the wide occurrence of bistability, and even tri-stability, we speculate that the occurrence of disturbances (such as fires, grazing, etc.) or a heterogeneous environment could increase noticeably the range of parameter values where coexistence may be observed, as discussed elsewhere (e.g., Baudena et al. 2010). At the same time, stochastic disturbances can enlarge the range where spatial patterns occur (Butler and Goldenfeld 2009;D'Odorico et al. 2007), most likely expanding also the range where coexistence  Table 1; c T 06.5 gm −2 mm −1 and r00.9 mm day −1 patterns could be observed in our model. The complexity of the outcomes at the interface between homogeneous grass cover and homogeneous woody cover, possibly indicate that in reality it is not surprising to observe savanna so widely in arid areas. Coexistence may not appear only as a stable equilibrium, but external disturbances acting on the multi-stable equilibria may drive it, possibly leading to intermediate states where both plant types could survive.
The model we introduced is very simplified, and clearly has several limitations. We did not include topography; rainfall intermittency and heterogeneity; the effects of fires, grazers, and browsers. These elements could possibly further enhance coexistence, and it would be very interesting to investigate their effects, but their inclusion was beyond the purpose of this paper. We did not include light competition, because in general water was the only limiting resource in savannas. We did not represent explicitly the root system, Fig. 5 Model with both infiltration and shading feedbacks. Different equilibria are reached starting from different initial conditions. The plots display the spatial distribution of biomass density (increasing from yellow to green, with yellow corresponding to bare soil). a Case with rainfall r00.55 mm day −1 , corresponding to point I in Fig. 2b. Top panels: starting with high tree and grass initial conditions (HTHG, top left panels, t00 years), after 10 years the system displays grass patterns (top right panels); initializing the system with high tree cover and low grass cover (HTLG, middle left panels, t00 years), the system reaches an equilibrium where a homogeneous tree layer is present (t010 years, middle right panels); initializing the system with low tree cover and low grass cover (LTLG, bottom left panels, t00 years), the system reaches an equilibrium where bare soil only is present (t010 years, bottom right panels). b Case with rainfall r00.75 mm day −1 , corresponding to point II in Fig. 2b. Top panels starting with high tree and grass initial conditions (HTHG, top left panels, time t00 years), after 10 years the system displays tree-grass coexistence patterns (top right panels); initializing the system with low tree cover and high grass cover (LTHG, middle left panels), grass patterns are observed (t0 10 years, middle right panels); when the system is initialized with high tree low grass cover (HTLG, bottom left panel, t00 years), a homogeneous tree layer is observed (bottom right panels, t010 years). In every run c T is set to 6.1 gm −2 mm −1 , see Table 1 for other parameter values thus neglecting competition between roots and shoots , and at the same time also neglecting the so-called root-water feedback, which is another important water-vegetation feedback. The larger the plant biomass, the further the root system develops, reaching more water. Gilad et al. (2007) show that this feedback is relevant for arid ecosystem dynamics, and may influence the way tree and grass patterns overlap. Finally, we did not include separately tree demographic stages, although we acknowledge that the main interactions (such as competition or facilitation) occur between grasses and tree seedlings (e.g. Baudena et al. 2010;Sankaran et al. 2004).
Despite these limitations, we demonstrated how vegetation-water feedbacks and spatial interactions could lead to stable coexistence of woody and herbaceous vegetation in arid savannas (as expected, e.g. Sankaran et al. 2005), within overlapping patterns (as observed, e.g. Cipriotti and Aguiar 2010;Meyer et al. 2008;Scholes and Archer 1997).
Our simple model, including two vegetation types competing for one resource, soil water, displayed a truly complex behavior, including pattern solutions and multi-stable states, emerging as a consequence of the vegetation feedbacks included (grass infiltration and tree shading feedbacks). These results underline clearly that such ecosystems can be sensitive to environmental (e.g., climatic) variations, and sudden catastrophic changes could occur. Therefore, such spatial dynamics must be taken into account when considering ecosystem structure and resilience.