1 Introduction
Banded vegetation patterns exist widely in the world's ecosystems [1–6]. Such patterns comprise alternating bands of vegetation and bare ground, and are a characteristic feature of landscapes in many arid and semiarid areas [7–9]. These patterns have been widely studied through field investigations and theoretical analysis due to their widespread occurrence and special features [10].
A large number of hypotheses aim to comprehend the self-organization of these patterns [8,10–14]. One important hypothesis addresses the nonlinear feedback between plant biomass and water resources of runoff–runon systems [15,16]. Rainfall on bare soil bands barely infiltrates, but runs downhill into vegetated bands and accumulates. Plant biomass can be greatly magnified in vegetated areas due to uptake of soil water by plant roots.
Many dynamic models based on feedback mechanisms of plant biomass and water resources have been established [8,10,11,15–20]. Two are regarded as core spatial vegetation-water models. The Klausmeier model addresses the formation of banded vegetation patterns in semiarid regions driven by feedback between biomass and water infiltration [8]. The second, proposed by HilleRisLambers et al. [15] and Rietkerk et al. [11], describes the water budget in detail as both soil and surface water. This model produces vegetation patterns of spots, labyrinths, gaps, and regular bands, dependent on rainfall and slope gradient [11]. The two models and their modified versions provide significant understanding of core mechanistic processes of vegetation pattern formation [21].
In the published literature on banded vegetation patterns, only a few researchers have considered the importance of micro-geomorphic processes [22–27]. Processes of soil erosion, sediment transport, and sediment deposition exert great influence on vegetation growth. A substantial amount of sediment is intercepted, trapped, and accumulated at vegetation patches [28,29]. Because the topsoil layer varies with soil erosion and sediment deposition, vegetation growth is greatly affected. Field observations and laboratory experiments show that soil erosion and/or sediment deposition may promote the development of banded vegetation patterns on hillslopes [1,28,30–35].
Considering the influences of soil variation aids understanding of the formation of banded vegetation patterns. For vegetation growth on hillslopes, soil factors dominate the water resource. First, mineral nutrition and water resources necessary for vegetation growth are stored in the soil layer. Second, soil provides a stable environment for vegetation growth. Disturbance of the soil by geomorphic processes leads to variations of water resources and nutrients for plants, affecting vegetation growth [36]. Thornes [37] and Zhang [38] proposed theoretical frameworks regarding the ecology of erosion and eroded ecosystems. From this perspective, focusing on influences of soil variation is more important than considering water resource effects on banded vegetation patterns. Moreover, soil thickness can be directly measured, reflecting variations of microtopography in banded patterns.
The natural world indeed offers examples of the formation of vegetation patterns resulting from soil variation. As observed by Gallart et al. [32], “grassed stairs” and terracette patterns are common features in the central Pyrenees. Generally, the terracettes are restricted to an altitudinal range between 1700 and 2700 m. The total plant cover of the terracette surfaces is often between 30 and 35%, and the main building species are bunchgrasses such as Festuca eskia, Festuca gautieri, Sesleria cerulea, etc. The growth habit of bunchgrasses results in a slow radial growth rate of the clump, and therefore the bunch reaches a high plant density in the clumps, while leaving areas of bare soil among them. On flat or gentle slopes, the bunch growth habit produces a mosaic of more or less circular random phases. As the gradient becomes steeper, bunches become roughly semicircular and tend to connect together, giving a pattern of continuous treads and risers following topographic contours.
Based on the ecological observations of Gallart et al. [32], substantial sediments are transported by overland flow and settle in vegetated areas intercepted by vegetation. Under negative influences of the accumulated sediments, small-scale vegetation bands appear on hillslopes. Thus, the interactions between sediment deposition and vegetation growth promote self-organization of banded vegetation patterns on hillslopes. However, little literature exists on this topic. A dynamic model of sediment deposition and vegetation growth on hillslopes is still needed for quantitative analysis of patterns that result from such interactions.
In this study, a spatial dynamic model is established, describing dynamics of sediment deposition and vegetation growth, based on the field observations of Gallart et al. [32]. The model aims to theoretically investigate the mechanism for the spontaneous emergence of banded patterns and stepped topography in the above case of the central Pyrenees as well as other cases with similar interactions and environmental conditions. Considering the characteristics of vegetation growth and sediment deposition in the above case, the spatial model is parameterized and utilized to simulate and analyze the formation process of banded vegetation patterns. However, it should be noted that the Pyrenees sites studied here are characterized by more humid climatic conditions than the arid and semiarid regions studied in much of the literature. Therefore, the model developed here is site-specific in that the feedbacks between sediment deposition and vegetation growth are likely more important than those between water redistribution and vegetation growth in the emergence of banded vegetation.
In developing the dynamic model of sediment deposition and vegetation growth, Thornes’ thought, which directly integrates geomorphic and ecological processes [37,39–41], is accepted. Via Turing analysis and numerical simulations on the model, the process and mechanism of formation of banded vegetation patterns are addressed. In this research, banded vegetation patterns describe a stable state of spatially heterogeneous vegetation, which shifts under influences exerted by sediment accumulation on vegetation growth.
2 Development of the model for formation of banded vegetation patterns
The dynamic model describes the formation of banded vegetation patterns on hillslopes, which is mainly influenced by the interactions between sediment deposition and vegetation growth.
2.1 Interactions between sediment deposition and vegetation growth
Although many studies use the feedback mechanisms between water and biomass to interpret vegetation pattern formation, vegetation pattern formation can result from the interactions between geomorphic processes and vegetation growth [23,27,37,40,42–44]. Hoffman et al. [45] concluded that the “current mathematical models of pattern formation in drylands typically consider only an infiltration-vegetation feedback and root augmentation growth as driving mechanisms, yet patterns of sedimentation and erosion on the soil surface have a strong effect on hydrological processes, stressing the need to introduce a soil-vegetation and an annual-shrub feedback.”
Gallart et al. [32] observed that banded vegetation patterns in high-mountain terracettes were common features on steep slopes between 1700 and 2700 asl in the central Pyrenees. As described by Gallart et al. [32], terracettes “show a convex downslope riser of bunchgrasses retaining a tread of almost bare soil upslope.” The authors understood the formation of such patterns as the result of the interactions between bunchgrass growth and geomorphic processes. Following the approach of Gallart et al. [32], the present study focuses on the mechanism of interactions between sediment deposition and vegetation growth and its effects on the formation of banded vegetation patterns in ecosystems.
As described in Gallart et al. [32], “most terracettes show coarse deposits which cannot be transported by runoff.” In this situation, two aspects are suggested: on the one hand, sediment deposition takes place along with other sediment processes including soil erosion and sediment transfer; on the other hand, sediment deposition prevails over other sediment processes in the formation of banded vegetation patterns. To quantitatively describe the effects of such sediment movement on formation of banded patterns; in this study, we consider a net sediment deposition, which is the result of sediment deposition subtracting the effects of soil erosion and sediment transfer. In the following, “sediment deposition” as used in this paper always means “net sediment deposition”. Net sediment deposition in nature is manifested by actual thickness of the sediment layer.
As a natural micro-geomorphic process, sediment deposition on hillslopes is promoted by vegetation cover [46,47], agreeing with tread retention in the terracettes of Pyrenees [32]. Persistent accumulation of sediments can negatively influence vegetation growth. Gallart et al. [32] showed that sediment deposition restricted vegetation growth and led to pronounced mortality and lower seedling establishment. The negative influence of sediment deposition is reflected in two aspects. First, sediment deposition affects the properties and structure of soil, leading to deterioration of the environment in which vegetation grows [48]. Second, the sediment layer formed by deposition of sand or even coarser sediments hardly retains water, which infiltrates quickly and runs away as subsurface flow [32]. Lack of water resource in the sediment/soil layer results in slow growth or mortality of vegetation.
In this study, the sediment deposition process is considered as the main dynamic process that promotes formation of banded patterns. The feasibility of this proposition is addressed:
- • As shown by Gallart, infiltration rates on the deposited area are as high as 700 mm/h, and water quickly exfiltrates at the feet of vegetated bands. Water redistribution is therefore determined by the property of sediments;
- • As the sediment deposits increase, the surface where plants establish rises. Since the water runs as subsurface flow, the root systems of plants need to change to reach the water. The sediment deposition thus influences vegetation establishment and how vegetation absorbs the water resource.
More importantly, sediment deposition changes the runoff–runon mechanism [25]. Deposit of coarse sediments causes a small difference in infiltration rate between bare and vegetated areas. Sediment deposition is more informative than water redistribution in describing the pattern formation.
Sediment deposits act as substitutes of the soil layer for vegetation establishment. The contents of mineral materials and organic matter vary between sediment deposits. Since sediment deposition results in a change in nutrient distribution, sediment deposition is considered as a more comprehensive index than water resources.
Sediment deposition, instead of water redistribution, is thus considered as the main dynamic process. Sediment redistribution is dynamic, affects vegetation growth, and promotes the formation of banded patterns. In describing the formation of banded patterns, Gallart et al. [32] did not even consider water redistribution. Thornes [37] and Burg et al. [41] also emphasized effects of sediment movement on vegetation growth, without considering water redistribution.
The relationship between vegetation biomass and deposited sediment is interactive. Two assumptions of the dynamic model are proposed:
- • Sediment deposition is the main dynamic process that affects vegetation pattern formation. The influence of water resource is considered as a subfactor included in the sediment deposition;
- • Sediments are transported by un-concentrated runoff. Banded vegetation cannot appear in landscapes with incised rills and gullies, where flow concentration precludes generation of sheet flow [26].
Interactions between sediment deposition and vegetation growth are briefly emphasized:
- • Vegetation promotes sediment deposition. Higher plant biomass brings more sediment deposits;
- • Accumulated sediment harms vegetation growth in that it deters its development into a homogeneous cover.
Therefore, a dynamic model that integrates vegetation growth and sediment deposition is established to study the self-organization of the banded vegetation patterns. In fact, direct integration of ecological and geomorphic processes is supported in Thornes’ work, which provides a framework for the ecology of erosion from a holistic perspective [37,40]. This model is applied to theoretically explain the self-organization of banded patterns and stepped topography in the case of the central Pyrenees as well as other cases with similar interactions and environmental conditions.
2.2 Model establishment
The dynamic model has two state variables for vegetation growth and sediment deposition. Thickness of the deposited sediment layer, S, is used for net sediment deposition, and ground plant biomass, V, is used for vegetation growth. The change of S is related to natural deposition and promotion of deposition by plant biomass. The change in V is determined by the intrinsic growth of vegetation and negative effect of accumulated sediments on vegetation.
The spatial distribution of S and V on planar hillslopes is considered. With the introduction of spatial terms, the model describes rates of change of S and V in both space and time. The dynamics of the two state variables are expressed in the following equations, after which a detailed interpretation of the equations is provided.
$$\frac{\partial S\left(x,y,t\right)}{\partial t}={A}_{0}\left(1-\frac{S\left(x,y,t\right)}{{S}_{\mathrm{m}}}\right)+bV\left(x,y,t\right)+u\frac{\partial S\left(x,y,t\right)}{\partial x}\text{,}$$ | (1a) |
$$\frac{\partial V\left(x,y,t\right)}{\partial t}=hV\left(x,y,t\right)\left(1-\frac{V\left(x,y,t\right)}{{V}_{\mathrm{m}}}\right)-pS\left(x,y,t\right)+d\left(\frac{{\partial}^{2}}{\partial {x}^{2}}+\frac{{\partial}^{2}}{\partial {y}^{2}}\right)V\left(x,y,t\right)\text{.}$$ | (1b) |
In Eq. (1a), A_{0}(1–S/S_{m}) describes the net sediment deposition process. A sediment resource term A_{0} represents the rate of net sediment deposition per unit of time. The parameter S_{m} is the maximal thickness of the sediment layer on bare ground. Since the rise of the sediment layer changes the routing of runoff and reduces the rate of sediment deposition, A_{0}S/S_{m} encompasses the self-damping of the sediment deposition process. Based on the assumption that the increase in sediment deposition caused by vegetation is simply proportional to plant biomass, bV is proposed. The parameter b is dependent on how plants retard the flow and trap sediments [29]. Sediment movement is kept with the flow and can be modeled by the convection term u∂S/∂x, where u is determined by the slope gradient. This convection term describes the deposited sediment flux along the downslope direction.
In Eq. (1b), a classic logistic growth term hV (1 − V/V_{m}) is used to describe plant biomass growth. h is the intrinsic vegetation growth rate and V_{m} is the maximum ground plant biomass given the ratio of actual to potential evapotranspiration [37]. The negative effect of deposited sediments on vegetation is modeled by a simple linear term pS. The lateral growth of vegetation is modeled by a diffusion term $d\left({\partial}^{2}V/\partial {x}^{2}+{\partial}^{2}V/\partial {y}^{2}\right)$, where d is a diffusion coefficient for plant biomass [8,11,15]. Table 1 provides interpretations and units of all the variables and parameters shown in Eq. (1).
Symbol interpretation, symbol units, and reference values of parameters.
Symbol | Interpretation | Units | Reference value | Source |
S | Thickness of deposited sediment layer | cm | – | – |
V | Ground plant biomass | kg·m^{−2} | – | – |
t | Time | Day | – | – |
x, y | Space | dm | – | – |
A _{0} | Maximal rate of sediment deposition at S = 0 | cm·day^{−1} | 0.01 | Descheemaeker et al. [60], Mabit et al. [61] |
S _{m} | Maximal thickness of deposited sediment layer on bare ground | cm | 20 | Descheemaeker et al. [60], Mabit et al. [61] |
b | Increased thickness of sediment layer per unit biomass | m^{2}·cm·kg^{−1}·day^{−1} | 0.005∼0.11 | – |
h | Maximal vegetation growth rate at V = 0 | day^{−1} | 0.005 | Klausmeier [8], Gilad et al. [16] |
V _{m} | Maximal ground plant biomass | kg·m^{−2} | 2 | Thornes [39] |
p | Decreased plant biomass per unit deposited sediments | kg·m^{−2}·cm^{−1}·day^{−1} | 0.00005 | – |
u | Speed of deposited sediment flux along the downslope | dm·day^{−1} | 0.01 | – |
d | Plant dispersal rate | dm^{2}·day^{−1} | 0.01 | Gilad et al. [16] |
In Eq. (1), we apply simple linear terms bV and pS as specific functions to describe the interactions between vegetation and sediment deposition. It should be noted that many other functional forms satisfying the above assumptions can be used to describe the interactions and the results are not sensitive to the exact functional forms. With the employment of the model described by eq. (1), we now can investigate the formation of banded vegetation patterns generated on hillslopes from interactions between sediment deposition and vegetation growth. Theoretical analysis proves that this model can generate heterogeneous vegetation patterns and determine conditions for pattern formation. Numerical simulations are performed to form banded vegetation patterns.
3 Conditions for Turing bifurcation and banded pattern formation
For ecological significance, dynamics of the two state variables, S and V, should be restricted in the region S ≥ 0, V ≥ 0. Also, all the parameters in the model such as A_{0}, S_{m}, b, u, h, V_{m}, p, d, should be positive.
The first step in analyzing how the model generates banded vegetation patterns is to determine stationary states of the non-spatial system in which the space derivatives of eq. (1) are set to zero. By solving the equations of the two nullclines, ∂S/∂t = 0, ∂V/∂t = 0, i.e.
$${A}_{0}\left(1-\frac{S}{{S}_{\mathrm{m}}}\right)+bV=0\text{,}$$ | (2a) |
$$hV\left(1-\frac{V}{{V}_{\mathrm{m}}}\right)-pS=0\text{.}$$ | (2b) |
Three cases of stationary states corresponding to different parameter conditions can be described as follows:
- • when Δ < 0, no stationary state exists for vegetation and deposition;
- • when Δ = 0 and hA_{0}/pbS_{m} > 1, there exists a unique positive stationary state denoted by E_{c} (S_{c}, V_{c}), where S_{c} = S_{m}(1 + bV_{c}/A_{0}) and V_{c} = V_{m}(1–pbS_{m}/hA_{0})/2;
- • when Δ > 0 and hA_{0}/pbS_{m} > 1, two positive stationary states exist, denoted by E_{1}(S_{1},V_{1}) and E_{2}(S_{2},V_{2}), where S_{1,2} = S_{m}(1 + bV_{1,2}/A_{0}) and ${V}_{\mathrm{1,2}}={V}_{m}\left(1-pb{S}_{\mathrm{m}}/h{A}_{0}\pm \sqrt{\Delta /{h}^{2}}\right)/2$ (here set V_{1} > V_{2}).
In these three cases, the expression of Δ is described by
$$\Delta ={\left(\frac{pb{S}_{\mathrm{m}}}{{A}_{}}-h\right)}^{2}-\frac{4hp{S}_{\mathrm{m}}}{{V}_{\mathrm{m}}}\text{.}$$ | (3) |
By the method of the Jacobian matrix, the local stability of stationary states in different cases is analyzed. The Jacobian matrix associated with the non-spatial system can be expressed by
$$J\left(S,V\right)=\left[\begin{array}{cc}-\frac{{A}_{0}}{{S}_{\mathrm{m}}}& b\\ -p& h\left(1-\frac{2V}{{V}_{\mathrm{m}}}\right)\end{array}\right]\text{.}$$ | (4) |
The expressions of the stationary states are substituted into Eq. (4) and the two eigenvalues of the Jacobian matrix are calculated. The signs of the two eigenvalues λ_{1} and λ_{2} determine local stability of corresponding stationary states. Through these steps, the local stability of various stationary states is determined:
- • The two eigenvalues of E_{c} satisfy λ_{1}λ_{2} = 0. At this stationary state, the dimension of the non-spatial system can be reduced to one. In this case, E_{c} is a saddle-node point;
- • For stationary state E_{2}, the two eigenvalues satisfy λ_{1}λ_{2} < 0. The positive and negative eigenvalues reveal that E_{2} is a saddle point;
- • For stationary state E_{1}, the two eigenvalues satisfy λ_{1}λ_{2} > 0. There are three cases for the local stability of E_{1}. Since ${\lambda}_{1}+{\lambda}_{2}=pb{S}_{\text{m}}/{A}_{0}-{A}_{0}/{S}_{\text{m}}-\sqrt{\Delta}$, when $pb{S}_{\text{m}}/{A}_{0}>{A}_{0}/{S}_{\text{m}}+\sqrt{\Delta}$, E_{1} is an unstable node or unstable focus. When $pb{S}_{\text{m}}/{A}_{0}<{A}_{0}/{S}_{\text{m}}+\sqrt{\Delta}$, E_{1} is a stable node or stable focus. When $pb{S}_{\text{m}}/{A}_{0}={A}_{0}/{S}_{\text{m}}+\sqrt{\Delta}$, the situation is complicated because the non-spatial system may go through a Hopf bifurcation. Therefore, the local stability of E_{1} cannot be determined without given specific values of the parameters.
Analyzing spatial system (1) can occur after obtaining the stable stationary state, which demonstrates that the vegetation can survive in the environments of deposited sediments. Moreover, the stable stationary state E_{1} implies that system (1) is also stable under spatially homogeneous perturbations. However, when spatially heterogeneous perturbations at this stationary state exist, the homogeneous stationary state may be unstable. Such a property determines the occurrence of Turing's instability and the formation of banded patterns of system (1).
Starting from the spatially stable homogeneous stationary state, a small but spatially heterogeneous amount of plants and sediment is added or removed. If these perturbations diverge with time, the system (1) will enter a new spatial state, implying Turing instability [15]. The perturbation equations at a spatially homogeneous stable stationary state E_{1} are expressed as:
$$S\left(x,y,t\right)={S}_{1}+s\left(x,y,t\right)\text{,}$$ | (5a) |
$$V\left(x,y,t\right)={V}_{1}+v\left(x,y,t\right)\text{.}$$ | (5b) |
Substituting the perturbed equations into the system (1) gives
$$\frac{\partial s\left(x,y,t\right)}{\partial t}=-\frac{{A}_{0}}{{S}_{\mathrm{m}}}s\left(x,y,t\right)+bv\left(x,y,t\right)+u\frac{\partial s\left(x,y,t\right)}{\partial x}\text{,}$$ | (6a) |
$$\frac{\partial v\left(x,y,t\right)}{\partial t}=-ps\left(x,y,t\right)+h\left(1-\frac{2{V}_{1}}{{V}_{\mathrm{m}}}\right)v\left(x,y,t\right)+d\left(\frac{{\partial}^{2}v\left(x,y,t\right)}{\partial {x}^{2}}+\frac{{\partial}^{2}v\left(x,y,t\right)}{\partial {y}^{2}}\right)\text{.}$$ | (6b) |
Expanding the small space-time heterogeneous perturbation in Eqs. (5) around the stable stationary state in Fourier space gives the following form of the perturbations:
$$s\left(x,y,t\right)=\overline{s}\left(t\right)\phantom{\rule{0.25em}{0ex}}{\mathrm{e}}^{\mathrm{i}\left({l}_{1}x+{l}_{2}y\right)}\text{,}$$ | (7a) |
$$v\left(x,y,t\right)=\overline{v}\left(t\right)\phantom{\rule{0.25em}{0ex}}{\mathrm{e}}^{\mathrm{i}\left({l}_{1}x+{l}_{2}x\right)}\text{,}$$ | (7b) |
$$\frac{\text{d}\overline{s}\left(t\right)}{\text{d}t}=\left(-\frac{{A}_{0}}{{S}_{\mathrm{m}}}+\text{i}u{\text{l}}_{1}\right)\overline{s}\left(t\right)+b\overline{v}\left(t\right)\text{,}$$ | (8a) |
$$\frac{\text{d}\overline{v}\left(t\right)}{\text{d}t}=-p\overline{s}\left(t\right)+\left(h\left(1-\frac{2{V}_{1}}{{V}_{\mathrm{m}}}\right)-d\left({l}_{1}^{2}+{l}_{2}^{2}\right)\right)\overline{v}\left(t\right)\text{.}$$ | (8b) |
The Jacobian matrix of the ordinary differential eq. (8) is
$${J}_{2}=\left[\begin{array}{cc}-\frac{{A}_{0}}{{S}_{\mathrm{m}}}+\text{i}u{l}_{1}& b\\ -p& h\left(1-\frac{2{V}_{1}}{{V}_{\mathrm{m}}}\right)-d\left({l}_{1}^{2}+{l}_{2}^{2}\right)\end{array}\right]\text{.}$$ | (9) |
Since the eigenvalues of matrix (9) determines the occurrence of Turing instability [15], the following characteristic equation can be calculated:
$${\left(\lambda \left({l}_{1},{l}_{2}\right)\right)}^{2}-\left(\gamma +n\right)\lambda \left({l}_{1},{l}_{2}\right)+\gamma \eta +bp=0\text{,}$$ | (10) |
$$\gamma =-\frac{{A}_{0}}{{S}_{\mathrm{m}}}+\text{i}u{l}_{1}\text{,}$$ |
$$\eta =h\left(1-\frac{2{V}_{1}}{{V}_{\mathrm{m}}}\right)-d\left({l}_{1}^{2}+{l}_{2}^{2}\right)\text{.}$$ |
The dispersion relation is obtained:
$$\lambda \left({l}_{1},{l}_{2}\right)=\frac{1}{2}\left(\gamma +\eta \pm \sqrt{\epsilon +\text{i}\theta}\right)\text{,}$$ | (11) |
$$\epsilon ={\left(-\frac{{A}_{0}}{{S}_{\mathrm{m}}}-\eta \right)}^{2}-{u}^{2}{l}_{1}^{2}-4bp\text{,}$$ |
$$\theta =2u{l}_{1}\left(-\frac{{A}_{0}}{{S}_{\mathrm{m}}}-\eta \right)\text{.}$$ |
Straightforward manipulation yields:
$$\mathrm{Re}\left(\lambda \left({l}_{1},{l}_{2}\right)\right)=\frac{1}{2}\left(-\frac{{A}_{0}}{{S}_{m}}+\eta +\frac{j}{\sqrt{2}}\sqrt{\sqrt{{\epsilon}^{2}+{\theta}^{2}+\epsilon}}\right)\text{,}$$ | (12) |
and
$$\mathrm{Im}\left(\lambda \left({l}_{1},{l}_{2}\right)\right)=\frac{1}{2}\left(u{l}_{1}+\frac{j}{\sqrt{2}}\text{sign}\left(\theta \right)\sqrt{\sqrt{{\epsilon}^{2}+{\theta}^{2}}-\epsilon}\right)\text{,}$$ | (13) |
As described above, Turing's instability and patterned vegetation need two conditions to be satisfied. First, system (1) is stable to spatially homogeneous perturbations, i.e. the non-spatial system has a locally stable stationary state. E_{1} satisfies this condition when $pn{S}_{\text{m}}/{A}_{0}<{A}_{0}/{S}_{m}+\sqrt{\Delta}$. Second, system (1) is unstable to spatially heterogeneous perturbations, i.e. the perturbed system (6) (or (8)) is unstable. At least one Re (λ(l_{1}, l_{2})) must be larger than zero, i.e.
$$-\frac{{A}_{0}}{{S}_{\mathrm{m}}}+\eta +\frac{1}{\sqrt{2}}\sqrt{\sqrt{{\epsilon}^{2}+{\theta}^{2}}+\epsilon >0}\text{.}$$ | (14) |
By direct calculation, this condition for the occurrence of Turing's instability and banded patterns is:
$$u>\frac{{A}_{0}-{S}_{\text{m}}\eta}{{S}_{\mathrm{m}}{l}_{1}}\sqrt{-\frac{bp{S}_{\text{m}}}{{A}_{0\eta}}}-1\text{.}$$ | (15) |
4 Numerical analyses
The results provided in Section 3 reveal that the spatial model of vegetation growth and sediment deposition yields banded vegetation patterns with proper parameters. However, the dynamical behaviors of the spatial model are difficult to obtain theoretically. It is thus necessary to perform numerical analyses, the most reliable approach with satisfactory accuracy and efficiency in comprehensively investigating characteristics of the model with feasible parameters.
Numerical analyses mainly focus on the stable stationary state and banded vegetation patterns. Three aspects of numerical simulations are performed:
- • cases of stationary states in the phase space and possible bifurcations of the stationary states for different parameters;
- • extensive simulations on parameter conditions of emerging banded vegetation patterns based on theoretical analysis;
- • display of banded patterns with given parameters and their possible analysis.
Table 1 provides parameter reference values for numerical simulations. For various situations in numerical simulations, the parameter values may be changed around these reference values. The reference values of some parameters (A_{0}, S_{m}, h, V_{m}, d) are literature-based, while values of other parameters (b, p, h) are mathematically determined according to the conditions for the occurrence of Turing patterns. Using these values, numerical simulations are performed analogously to real cases.
Stability analysis on the non-spatial model proves the existence of stationary states. These states imply stable or unstable establishment of vegetation covers in depositional areas. Moreover, the stability of stationary states varies with parameter changes. Fig. 1 shows six possible cases of stationary states with different groups of parameters in phase diagrams. The stationary state E_{1} can either be attracting (the stable node in Fig. 1a and the stable focus in Fig. 1b) or repelling (weakly repelling equilibrium in Fig. 1c, unstable focus in Fig. 1d, and the unstable node in Fig. 1e). Simultaneously, stationary state E_{2} is always a saddle in the five cases represented by Fig. 1a–e. For Fig. 1f, E_{1} and E_{2} merge together into the saddle-node E_{c}.
The dynamic behavior of stationary state E_{1} presents gradual changes among the six cases, as the value of any one parameter continuously varies. Changes suggest two types of bifurcations that appear in the interactions between sediment deposition and vegetation growth described in Eq. (1). The saddle-node bifurcation occurs at the transformation between E_{c} and E_{1,2}. When the stationary state E_{1} transforms from stable focus into unstable focus, the Hopf bifurcation generates a weakly repelling equilibrium. The Hopf bifurcation indicates a sudden shift between two different dynamic outcomes, namely, that vegetation will either stably exist or die in depositional environments.
The bifurcation diagrams for parameters A_{0}, b, h, and p are plotted to better comprehend them. As shown in each graph of Fig. 2, two branches of the stationary states, E_{1} and E_{2}, are present. Due to the Hopf bifurcation, E_{1} transforms between stable states and unstable states. Meanwhile, as a result of saddle-node bifurcation, saddle-branch E_{2} is always at an unstable state.
The saddle-node bifurcation provides parameter values for the appearance of coexisting stationary states between vegetation and deposited sediments. In Fig. 2, the parameter values at which saddle-node bifurcation occurs are A_{0s} = 8.82×10^{−3}, b_{s} = 1.84×10^{−2}, h_{s} = 4.68×10^{−3}, and p_{s} = 5.34×10^{−5}. More importantly, the Hopf bifurcation determines critical values at which a sudden shift occurs between the vegetated and bare ground states. The critical values are A_{0H} = 9.90×10^{−3}, b_{H} = 1.64×10^{−2}, h_{H} = 4.96×10^{−3}, and p_{H} = 5.04×10^{−5}. Only when A_{0} > A_{0H}, b < b_{H}, h > h_{H} or p < p_{H}, can a stable vegetation cover exist and banded patterns be observed.
Emergence of banded vegetation patterns requires Turing instability. According to the analysis provided in Section 3, the parameters must satisfy $\mathrm{Re}\left(\lambda \left(\mathrm{0,0}\right)\right)<0$ and $\mathrm{Re}\left(\lambda \left({l}_{1},{l}_{2}\right)\right)>0$ when l_{1} and/or l_{2} are not equal to zero. From eq. (12), satisfying the latter condition requires 2 V_{1} < V_{m}. Therefore, as shown in Fig. 2, a value range for each parameter is simulated. The shaded area in each graph represents this value range, describing a necessary combination of A_{0}, b, h, p, S_{m} and V_{m} for Turing's instability. The parameter ranges described in Fig. 2 are A_{0H} < A_{0} < 1.08×10^{−2}, 1.50×10^{−2} < b < b_{H}, h_{H} < h < 5.24×10^{−3}, and 4.72×10^{−5} < p < p_{H}.
Considering the necessary conditions for Turing's pattern (i.e. Δ > 0, $h\left(1-2{V}_{1}/{V}_{\text{m}}\right)-{A}_{0}/{S}_{\text{m}}<0$, and 2 V_{1} < V_{m}), two cases of region diagrams are plotted for possible vegetation patterns for the non-spatial system (Fig. 3). In each region diagram, three parameter areas, homogeneous vegetation, no vegetation, and banded vegetation, are separated. Homogeneous vegetation implies that the S–V system exhibits a state of uniform vegetation biomass density. No vegetation is defined by the absence of a stationary state or an unstable one. Turing instability and banded vegetation occur in the parameter area enclosed by C1 and C2 or C3 and C2. From combination of α, β and δ, the value range of any needed parameter(s) for vegetation pattern formation can be easily calculated.
The ranges of parameters u and d for emerging banded vegetation patterns demand simulation of the dispersion relation. As shown in Fig. 4, Re (λ(l_{1}, l_{2})) first increases and then drops as wavenumber l_{1} increases. When the vertex point of the curve Re (λ(l_{1}, l_{2})) is above zero, i.e. u > u_{c} or d < d_{c}, the spatial model generates banded vegetation patterns. The threshold values u_{c} and d_{c} imply that the slope gradient and the plant dispersal rate determine pattern formation.
With the parameters selected, numerical simulations are provided through discretizing partial differential equations of system (1) in two-dimensional space. The first-order upwinding finite difference scheme is applied to discretize the convection term and finite difference approximation to the diffusion term. An explicit Euler method is applied for the time integration with a time step size Δt = 0.05[day] and space step size Δh = 2[dm] [49]. The space and time scale is averaged for the Euler method. Spatial vegetation patterns are shown in a rectangular domain including 100 × 100 grid elements, representing 20 [m] × 20 [m]. Periodic boundary conditions are employed in the pattern simulations [8]. Initial conditions are given as perturbing 9% of the grid elements in stable homogeneous vegetation and letting biomass in these grid elements equal zero.
Fig. 5 shows the evolution of the banded vegetation pattern at 6, 55, 165, and 494 years, with given parameters. In Fig. 5a, b, and c, the vegetated patches grow laterally to gradually generate parallel vegetation bands. These vegetation bands are perpendicular to the flow direction. Fig. 5d shows a completely developed banded vegetation pattern. Moreover, this banded pattern oscillates in time and migrates in the upslope direction.
Fig. 6 shows how the change of sediment deposition process influences vegetation pattern features. As shown in Fig. 6a–b, the rate of sediment deposition doubles. Accordingly, the wavelength of the banded pattern decreases, but the number of vegetation bands increases. Fig. 6c–d show that decreasing S_{m} leads to a similar consequence. Results in Fig. 6 reveal that sediment deposition acceleration results in decreasing wavelength in banded vegetation patterns. Moreover, there is a time lag before the stable pattern is reached.
A phenomenon of the banded pattern formation is shown in Fig. 7. The parameters presented in this case imply that the stationary state E_{1} is unstable (see Fig. 1d) and that the Turing instability cannot appear. However, the banded vegetation pattern emerges from the dynamic influences of system (1) on the boundary S = 0 and V = 0. Since S and V must be positive for practical ecological sense, S or V is set to zero when each value becomes negative in the simulations. This results in self-organization of a banded vegetation pattern on the slopes. Moreover, the formation of banded patterns is accelerated, with the lag time in their formation reduced to about 439 years.
When the patterns of S and V are observed along the downslope direction, waves of plant biomass and sediment layer are found. Fig. 8a shows the alternating bands of vegetation and bare ground along the downslope direction. The horizontal line represents the plant biomass at the stable stationary state. Plant biomass at vegetation bands is much higher than at the stable stationary state. Asynchronism between the peak of vegetation bands and the bottom of sediment bands reveals the lag in the effect of interactions between sediment deposition and vegetation growth.
To show the banded vegetation pattern on hillslopes, Fig. 8b is plotted. For analogy to the real pattern in the central Pyrenees, we chose a particular hillslope gradient (25%), as recorded in Gallart et al. [32] to plot Fig. 8b. Fig. 8b shows similar characteristics with the vegetation patterns described in literature. First, alternation of vegetated areas and bare areas is demonstrated. Second, the deposited sediments trapped in vegetated areas lead to the formation of treads on hillslopes. With the demonstration in Fig. 8b, the coevolution of vegetation patterns and landforms is quantitatively described in this study.
5 Discussion
Vegetation pattern formation on hillslopes mainly results from the balance of positive and negative feedback between resources for plants and plant biomass [10,11,16,50,51]. Since erosion and sediment deposition disturb the soil layer, affecting resources stored for plants, the spatial distribution of plant biomass is greatly influenced. This also leads to the formation of banded vegetation patterns [28,32,52].
Field investigations show that natural formation of banded vegetation patterns in the central Pyrenees is related to sedimentation on hillslopes and interactions between deposition and vegetation growth [32]. We determined that banded vegetation patterns self-organize when model parameters assume proper values (see the conditions for pattern formation described in Section 3). Using feasible parameter values provided in the literature (see Table 1), numerical simulations are performed analogously to real cases, and the main characteristics of banded pattern formation are quantitatively described.
Banded vegetation patterns were also discovered on semi-natural and abandoned land in the upper Guadalentín basin in southeastern Spain [34]. Cammeraat and Imeson found a strong interrelationship between plant development, sediment deposition, and soil surface crusting, leading to the development of fine-scale banded vegetation structures in degraded lands, dominated by Plantago albicans. At much lower slope positions, sediment deposition clearly becomes more and more important, and blankets of silt-sized material are deposited there. A thick silty crust is present, and the land shows bare and strongly crusted soils alternating with bands of P. albicans, developed in a direction parallel to the contours, at a more or less regular interval. The small bands of plants act as a sediment trap during the rare periods of overland flow, when silty sediments are transported downslope. Due to the accumulation of mainly fine-soil material and some organic debris, micro-topographical steps are also developed, coinciding with the bands of P. albicans. Directly under the small steps, erosion takes place during overland flow, creating small steep cuts in the surface with a height of about 2 cm, exposing coherent layered silts. Cammeraat and Imeson [34] suggested that the vegetation banding and resultant stepped topography were closely related to dynamics of sediment deposition. However, patterns described by Gallart et al. [32] and Cammeraat and Imeson [34] resulted from sediment deposits of different properties.
The field survey cases as described in Gallart et al. [32] and Cammeraat and Imeson [34] provide evidence of banded vegetation patterns that we can model. As can be seen from Figures 5 through 7, banded patterns of vegetation are completely developed under the conditions considered. Development of vegetation bands can be quantitatively described with this model. Employing the model, we found the following ecologically significant results:
- • Decrease/increase of the sediment layer leads to increase/reduction of biomass in vegetation bands. Peaks of vegetation bands correspond to the thinnest parts of sediment layer. However, the lagging effect results in asynchronism between vegetation and sediment bands;
- • Plants at boundaries of upslope vegetation bands suffer from negative effects of deposited sediments. This ultimately restricts the lateral growth of vegetation patches, and restriction results in the formation of parallel vegetation bands.
Considering ecosystem features, the simulated wavelengths and the real patterns are compared. Wavelengths of vegetation bands in simulated banded vegetation patterns range from 1.6 m to 2.5 m. The range observed by Gallart et al. [32] was from 0.2 m to 5 m, while Cammeraat and Imeson [34] reported wavelengths of 1.2 m in average. Thus, our simulated results are comparable.
Emergence and migration of vegetation bands require unidirectional runoff flow and sloped terrains [8,9,11,13,53]. Vegetation band formation is not directly dependent on slope gradient. Since sediment deposition is the main dynamic process in pattern formation, blocking of runoff flow and interception of sediments are important factors.
Field investigators have reported sediment deposition as being important for vegetation pattern formation [34,52,54,55]. These patterns also are formed by different mechanisms [7,56,57]. Bryan and Oostwoud Wijdenes [31], and Bryan and Brun [33] postulated that these bands result from positive effects of accumulated sediments in depositional zones. Interference factors such as overgrazing and fire lead to fragmentation of vegetation bands in nature.
The spatial vegetation-sediment model considers the effectiveness of sediment fluxes in the formation of banded vegetation. In contrast, the Klausmeier model considers the striped vegetation patterns to be generated merely by interactions between vegetation growth and water resources [8]. Sediment fluxes are more effective in forming banded patterns than water fluxes because water must be restored in soil before absorption. Redistribution of water in forming banded vegetation patterns is regarded as a consequence of disturbing soil layers by sediment deposition.
The important role of sediment fluxes in determining vegetation pattern formation is widely recognized, and many studies have contributed to the literature on this topic [25,28,45,52,58]. However, still few models are documented in the literature for predicting the coevolution of vegetation patterns and hillslope landforms [26,32,45,59]. Based on previous studies, the research reported here makes improvements. Directly integrating sediment flux and vegetation growth in a dynamic model allows the coevolution of vegetation patterns and hillslope landforms be calculated explicitly [32,37,41]. By comparing our work with that of Gallart et al., the effects of parameter variations on the coevolution and the characteristics of vegetation patterns and landforms can be quantitatively determined.
6 Conclusion
Formation of banded vegetation patterns on hillslopes, as an effect of sediment deposition, is investigated theoretically and numerically. With the consideration of the interactions between sediment deposition and vegetation growth, a relevant and novel spatially dynamic model is thus developed. To the best of our knowledge, such a model describing banded pattern formation on hillslopes due to sediment deposition has not yet appeared in the literature.
In this study, the analysis of Turing mechanisms has led to determination of the parameter conditions for the occurrence of banded vegetation patterns. Numerical simulations based on the model established simulate the self-organization process of the patterns under different conditions. Vegetation bands are found perpendicular to the direction of sediment flux. Also, with the model established, it is evident that treads on hillslopes emerge with the evolution of vegetation bands. The coevolution of banded vegetation patterns and landforms due to the interactions between sediment deposition and vegetation growth, as described quantitatively with the model, agrees with the field observations reported in the literature.
Disclosure of interest
The authors declare that they have no competing interest.
Acknowledgments
The authors would like to acknowledge with great gratitude the support of National Major Science and Technology Program for Water Pollution Control and Treatment (No. 2017ZX07101-002), the Fundamental Research Funds for the Central Universities (No. JB2017069), the China Scholarship Council (No. 201206730024), and the Chinese Natural Science Foundation (Project 39560023 to H.Z.). The authors also would like to thank the anonymous referees for their carefulness, helpful comments and suggestions, which greatly improved this paper.