International Institute for Applied Systems Analysis Schlossplatz 1
A-2361 Laxenburg, Austria
Tel: +43 2236 807 342 Fax: +43 2236 71313 E-mail: publications@iiasa.ac.at Web: www.iiasa.ac.at
Interim Reports on work of the International Institute for Applied Systems Analysis receive only limited review. Views or opinions expressed herein do not necessarily represent those of the Institute, its National Member Organizations, or other organizations supporting the work.
Interim Report IR-05-085
Self-thinning and Community Persistence in a Simple Size- structured Dynamical Model of Plant Growth
F. Dercole (dercole@elet.polimi.it) K. Niklas (kjn2@cornell.edu) R. Rand (rhr2@cornell.edu)
Approved by Ulf Dieckmann
Program Leader, Evolution & Ecology
December 2005
Contents
Abstract... 1
Introduction ... 1
Mathematical Model... 3
Community persistence and equilibrium size density distributions ... 7
Examples ... 9
Transient behavior ... 10
Discussion... 11
Conclusions ... 16
Appendix I: Stability of the trivial solution... 18
Appendix II: Steady state equilibrium... 19
References ... 20
Digital Object Identifier (DOI):
10.1007/s00285-005-0322-x
Mathematical Biology
F. Dercole·K. Niklas·R. Rand
Self-thinning and community persistence in a simple size-structured dynamical model of plant growth
Received: 29 September 2004 / Revised version: 2 February 2005 / Published online: 2 May 2005 – cSpringer-Verlag 2005
Abstract. This paper presents a size-structured dynamical model of plant growth. The model takes the form of a partial differential-integral equation and includes the effects of self- shading by leaves. Closed form solutions are presented for the equilibrium size density dis- tribution. Analytic conditions are derived for community persistence, and the self-thinning exponent is obtained as a function of species characteristics and environmental conditions.
1. Introduction
Most mathematical models agree that migration, birth, and death are the three most important processes dictating the growth and dynamics of populations or entire communities (Chew and Chew, 1965; Harper, 1977; Nagano, 1978; Hara, 1984a,b, 1985; Takada and Iwasa, 1986; Petersen, 1988; West et al., 1989; Petersen et al., 1990; Weiner, 1990; Kohyama, 1991, 1992). However, the search for general prin- ciples of plant population dynamics is unlikely to produce robust generalizations, unless specific habitat and life forms are taken into account. Therefore, population dynamic models of plant growth should be specifically addressed to certain habi- tat/forms type (e.g. trees in forest, shrubs in desert, annuals in early successional communities, herbaceous perennials in grasslands) and should account for both the spatial and size structures of the population (Leslie, 1945; Lefkovitch, 1965;
Schaffer and Leigh, 1976; Antonovics and Levin, 1980; Pacala and Silander, 1985;
Metz and Diekmann, 1986; Crawley, 1990; Cz´ar´an and Bartha, 1992; De Angelis and Gross, 1992; Bishir and Namkoong, 1992; Hara, 1992; Hara and Wyszomirski, 1994; Durrett and Levin., 1994; Takada and Hara, 1994; Diekmann et al., 1998;
Cushing, 1998; Diekmann et al., 2001; Bolker et al., 2003).
In this regard, population size structure is the most apparent and the most appli- cable to applied biological sciences such as agronomy and forestry. For example, world-wide data sets describing the size density distribution of tree-dominated F. Dercole: Department of Electronics and Information, Politecnico di Milano, Milano, Italy.
e-mail:dercole@elet.polimi.it
K. Niklas: Department of Plant Biology, Cornell University, Ithaca, NY USA.
e-mail:kjn2@cornell.edu
R. Rand: Department of Theoretical and Applied Mechanics, Cornell University, Ithaca, NY USA. e-mail:rhr2@cornell.edu
Key words or phrases: Allometric relationships – Basal stem diameter – Dynamical models – Persistence conditions – Plant size – Self-thinning – Size density distribution – Size structure
1
communities are available and have been used to calibrate allometric relation- ships, namely how trees attributes such as number per unit area, standing leaf biomass, basal stem diameter, height, change with variation in tree size (Yoda et al., 1963; Enquist and Niklas, 2001). Theoretical allometric studies based on mechan- ical (McMahon, 1973) and biophysical (Niklas, 1994; West et al., 1997, 1999b,a;
Enquist et al., 1998, 1999; Niklas et al., 2003) principles, result in one of plant ecology’s most robust empirical generalizations: the so-called self-thinning rule. It predicts that the number of trees per unit area scales as the -2 power of basal stem diameter or as the -3/4 power of tree biomass, which also implies, by assuming that tree biomass is proportional to tree volume (approximately basal stem diameter times height), that tree height scales as the 2/3 power of basal stem diameter. In contrast to such theoretical studies, statistical analysis of real data sets confirms the allometric power law relationships but yields values of the power law exponent which may differ from the theoretically derived values. E.g. (Niklas et al., 2003) found that the number of trees per unit area scales as the -1.75 power of basal stem diameter based on a world-wide data set for the size density distributions of 226 communities consisting of woody plant species.
In this paper, we develop a non-spatial, size-structured continuum model of plant growth, without focusing on a particular species, but with emphasis on a dense tree-dominated forest. We intentionally sacrifice biological realism for mathemati- cal tractability, but we derive a closed form solution for the equilibrium size density distribution and analytical conditions for communities persistence. Our model is based on the following features: (i) all individuals larger than a designated size class contribute new individuals to the smallest size class (in the form of seeds), (ii) individuals continue to grow in size (and thus pass into progressively larger size classes) until they reach a species-specific maximum size class beyond which they do not grow owing to death, (iii) individuals in each size class (other than the smallest class) shade those in the smaller size classes, and (iv) individuals in all size classes die as a result of stochastic abiotic or biotic processes. Intraspe- cific competition for water and soil nutrients and soil nutrient recycling through decomposition of dead plants are omitted from the model. This is a simplification of reality but is expected to be valid when water and soil nutrients are available at nonlimiting quantity and concentrations. Herbivores effects are also omitted. How- ever, there is general agreement among plant ecologists that asymmetric intra- and inter-specific competition for light is the driving force in the dynamics of most plant communities living in non-desert environments (Crawley, 1990) and that removal of herbivores often has little effect on plant population dynamics and persistence (Crawley, 1989).
The model permits us to assess how the dynamics of the size density distri- bution of a population depends on its species-specific growth and seeding rates, shading capacity, and mortality. The model allows us to simulate the dynamics of a population from its inception (germination of broadcast seeds) to maturity as well as to evaluate the equilibrium condition, i.e., the structure of a stable size density distribution and its persistence under the foregoing parameter perturbation. Despite its simplicity the model shows a remarkable agreement with the self-thinning rule
as tested by prior works using a world-wide data set for the size density distributions of 226 communities consisting of woody plant species (Enquist and Niklas, 2001).
The rest of the paper is organized as follows. Section 2 presents the mathemat- ical model. We first consider a discrete size structure and we later take the limit to an infinite number of size classes to obtain a continuum model. We feel that this approach is more biologically transparent than directly deriving the continuum model. Section 2.1 reports on the analytical conditions for the persistence of the plant community and the closed form solutions for the equilibrium size density distributions (mathematical details are relegated to Appendices I and II), which are pictured in some examples in Section 2.2. Section 2.3 shows the transient behav- ior of the model reporting some examples obtained via numerical simulation. The discussion of the mathematical results and their biological interpretations in terms of community persistence and self-thinning are given in Section 3, while a criti- cal revisitation of the model and some remarks on possible future works close the paper.
2. Mathematical model
We consider the size density distribution of a single species withM+1 discrete size bins, where the size variable is plant basal stem diameterD. Emphasis is placed on basal stem diameter because it is the easiest to measure and most often used size variable in the field. We note, however, that the model is equally appropriate and robust if the ‘single species’ is considered to reflect an ‘average species’ for each of theM+1 size bins. We now consider the discrete version of the model, but we will later extend the model by taking the limit as M approaches infinity to obtain a continuum model, which is more tractable mathematically. For direct (but equivalent) approaches to the continuum model see Webb (1985); Metz and Diekmann (1986); Cushing (1998); Diekmann et al. (1998, 2001).
For the discrete version, the variables describing the state space areN0,N1,..., NM, each corresponding to a measure of the plants which lie in size binsD0,D1,..., DM, respectively. The equations of motion are:
dN0
dt = −γ N0−λN0(N1+N2+ · · · +NM)−µN0+σ (N1+N2+ · · · +NM) (1) dNi
dt =γ (Ni−1−Ni)−λNi(Ni+1+Ni+2+ · · · +NM)−µNi (1≤i≤M) (2) whereγ is a scaling coefficient for the average rate of growth in size (basal stem diameter) for all the individuals in size binDi. This parameter dictates the rate at which individuals ‘migrate’ from one size bin to the next larger bin size as a consequence of increasing size. The parameterλdenotes asymmetric competition for light, i.e. the negative effects of larger individuals on the growth of smaller counterparts in terms of light capture. Larger (taller and wider) individuals will necessarily shade smaller individuals growing in their shadow as a result of their elevated and broader leaf canopies. The mortality coefficientµreflects the effects of
a variety of factors that have the potential to cause death (e.g., pathogens, herbivory, and local mechanical perturbation such as wind-throw). This parameter reduces the number of individuals in any given size bin. Note that asymmetric competition in model (1), (2) is modeled as an extra mortality (Hara, 1985). Other authors have interpreted asymmetric competition as a limitation to growth, i.e. as a reduction of the rate of ‘migration’ to the next class (Takada and Iwasa, 1986; Hara, 1992).
While it is obvious that shading has both the effects of increasing mortality and lim- iting growth, the two effects appear as different terms in the equations. Our choice here has been made for convenience, in order to obtain closed form equilibrium solutions to the continuum model (see Appendix II). The seeding coefficientσ is a measure of the fecundity of the individuals comprising a population (or that of the
‘average species’ reflecting a mixed-species community), i.e., the number of viable seeds produced per individual per unit time. All individuals (except those in the first class) produce seeds until they die, since no evidence for a threshold plant size below which reproduction is impossible has been documented (Rees and Crawley, 1989). Finally, note that the model proscribes a maximum size. Individuals reach- ing this size survive for a time (they grow in size within the largest bin size) but they eventually die and leave this size bin.
Briefly then,γis a growth coefficient,λis a shading coefficient,µis a mortality coefficient, andσis a seeding coefficient. See Fig. 1.
We have taken the coefficientsγ,λ,µandσ as constants. In the case ofλ,µ andσ this has been done for convenience. A more general (but complex) model
Fig. 1. Schematic of a simple size density distribution withNinumber of plants in each of three size categories (defined in terms of basal stem diameterDi). Four parameters are used to model the dynamics ofNi: the seeding coefficientσ, which characterizes the number of seeds produced per individual per unit time; the shading coefficientλ, which characterizes the extent to which larger plants shade smaller plants; the growth coefficientγ, which char- acterizes the rate at which plants grow and thus ‘migrate’ from one size category to the next (e.g. fromD0toD1); and the mortality coefficientµ, which characterizes the loss of plants from each size category as the result of stochastic abiotic or biotic processes (e.g. wind throw and disease, respectively)
would allow these parameters to depend oni, i.e. to vary as a function of basal stem diameter. In contrast, the assumption thatγ is constant is predicated on empirical observations (Chew and Chew, 1965; Niklas, 1997). Specifically, 1/γ is the time required for a plant to grow from one size bin to the next, i.e.,γ is the diameter growth rate. Allometric theory shows that the height of a vascular plantHscales as a power law of the basal stem diameter, i.e.H ≈Dh(Enquist and Niklas, 2001), from which it follows that plant volume and plant mass are proportional toD2+h. Taking a time derivative we see that biomass growth rate is proportional toD1+h times diameter growth rateγ. Moreover, it has also been shown that the biomass growth rate at the level of an individual plant is proportional to the standing leaf biomass, which scales asH D≈D1+h(Enquist and Niklas, 2002). Therefore, the diameter growth rateγ is constant across all size bins.
To transform the discrete model (1), (2) into a continuum model, we replace the discrete model’s bin measureNi by a densityn(D, t)in the continuum model (in the following we will omit the time dependence when not necessary). The two are related by the equation:
Ni = Di+h
Di
n(D)dD (3)
wherehis the bin width (which is typically constant in the available size density distributions for real plants). Eq. (3) is approximated for smallhby the relation
Ni =n(Di) h (4)
We now nondimensionalize the bin size variableDso thatDmax=1 by defining D˜ = D
Dmax
and then dropping the∼for convenience.
Note that
∂n
∂D ≈n(Di)−n(Di−1)
h ≈Ni −Ni−1
h2 (5)
and 1
D
n dD ≈(n(Di+1)+n(Di+2)+ · · · +n(DM))h
≈(Ni+1+Ni+2+ · · · +NM) (6) Dividing Eq. (2) byhgives
1 h
dNi dt =
h2(Ni−1−Ni)−λNi
h (Ni+1+Ni+2+ · · · +NM)−µNi h where=γ h. Note that ash→0,γ → ∞such that the productγ his constant.
Using Eqs. (4), (5), (6) and taking the limit ash→0, we obtain:
∂n
∂t = −∂n
∂D −λn
1
D
n dD−µn (7)
We next divide Eq. (1) byh:
1 h
dN0
dt = −
h2N0−λN0
h (N1+ · · · +NM)−µN0
h +σ
h(N1+ · · · +NM) Again using Eqs. (4), (5), (6), we obtain the approximate intermediate result:
∂n
∂t = − hn−λn
1
0 n dD−µn+σ h
1
0 n dD, D=0
where we have taken the smallest size class asD=0 for convenience. Multiplying byhand taking the limit ash→0, we obtain the boundary condition:
n=σ
1
0
n dD, D=0 (8)
The continuum model consists of Eq. (7) with the boundary condition (8). Eqs. (7), (8) are a particular version of the Sharpe and Lotka (1911)-McKendrick (1926)-von Foerster (1959) equation
∂n
∂t = − ∂
∂D(G(D, t)n(D, t))−M(D, t)n(D, t) with boundary condition
G(D, t)n(D, t)=R(t), D=0
which has been used by many authors in plant ecology (Nagano, 1978; Takada and Iwasa, 1986; Kohyama, 1992) as well as in various other fields of biology (Trucco, 1965; Shinko and Streifer, 1967; Levin and Paine, 1974; Paine and Levin, 1981;
Kirkpatrick, 1984). Likewise, the Sharpe-Lotka-McKendrick-von Foerster equa- tion is a particular version of the more general Fokker-Plank equation (also known as diffusion or Kolmogorov forward equation), which describes the moment dynam- ics of a probability distribution and has been used by Hara (1984b) and many others to incorporate the spatial and stochastic fluctuations of the size density distribution.
The original Sharpe-Lotka-McKendrick-von Foerster equation is linear. Nonlinear age- and size-structured population models have been first proposed by Gurtin and MacCamy (1974) and Hoppensteadt (1974) (see also Webb, 1985; Metz and Diekmann, 1986; Cushing, 1998; Diekmann et al., 2001). Traditionally, the source of nonlinearity is the density dependence of some of the characteristic population rates (birth, death, migration) on the total population size
NT =
1
0
n dD (9)
(see Eq. (6)). Such nonlinearities pose integral conditions on the solution of the boundary-value problem, which typically prevent any possibility of deriving closed form equilibrium solutions. Conversely, the nonlinearity present in our model is due to asymmetric competition, often called one-side competition because the char- acteristic population rate, mortality due to shading, is only affected by one side of
the size distribution, with respect to the shaded plant, i.e. by bigger plants. This results in the partial differential-integral Eq. (7), whose equilibria can be solved analytically (see Appendix II). A similar model, in which the appearance of such a nonlinearity arises from hierarchical competition, is described in Kinzig et al.
(1999) and inspired the present work.
Although model (7), (8) involves four parameters,λ,µandσ, the number of parameters can be reduced to two by defining the following normalized quantities:
tˆ=t, nˆ = λ
n, µˆ = µ
, σˆ =σ
(10)
In terms of these quantities, Eqs. (7), (8) become:
∂nˆ
∂tˆ = −∂nˆ
∂D − ˆn
1
D
ˆ
n dD− ˆµnˆ (11) ˆ
n= ˆσ
1
0 n dD,ˆ D=0 (12)
2.1. Community persistence and equilibrium size density distributions
Note that the model (11), (12) exhibits the trivial solutionnˆ ≡0, which leads to the local extinction of the species being modeled. We show in Appendix I that the trivial solutionnˆ ≡0 is stable if and only if
ˆ
σ < µˆ
1−e− ˆµ (13)
and unstable otherwise. Eq. (13) indicates that the plant species will not survive unless the seeding coefficientσˆ is large enough relative to the mortality coefficient
ˆ
µ. See Fig. 2.
The equations governing the equilibrium size density distribution (the steady state equilibrium) in the model (11), (12) are:
dnˆ dD+ ˆn
1
D
ˆ
n dD+ ˆµnˆ=0 (14) ˆ
n= ˆσ
1
0
ˆ
n dD, D=0 (15)
In Appendix II we show that a nontrivial steady state solution exists in the region of theµ-ˆ σˆ parameter plane in which the trivial solution is unstable, see Figure 2.
The analytical expression of the steady state solution is different in the two regions A and B in Fig. 3, which are separated by the following curve:
ˆ σ = 2
2− ˆµ (16)
In region B, the steady state solution takes the form:
ˆ
n(D)= 4k2(µˆ2−4k2)
4k2− ˆµ2+(4k2+ ˆµ2)cosh 2k(1−D)−4kµˆsinh 2k(1−D) (17)
Fig. 2. Stability of the trivial solutionnˆ ≡0. See Eq. (13)
wherekis a parameter related toµˆ andσˆ by the equation:
ˆ
µσˆsinh2k−2σ kˆ sinhkcoshk+2k2=0 (18) and in region A, we have the steady state solution:
ˆ
n(D)= 4c2(µˆ2+4c2)
4c2+ ˆµ2+(4c2− ˆµ2)cos 2c(1−D)−4cµˆsin 2c(1−D) (19) where:
ˆ
µσˆsin2c−2σ cˆ sinc cosc+2c2=0 (20) On the curve (16) separating regions A and B, the steady state solution takes the form:
ˆ
n(D)= 2µˆ2
((1−D)µˆ −2)2 (21) See Appendix II for details.
Fig. 3. The curves given in Eqs. (13), (16) divide theµˆ-σˆ plane into 3 regions (here des- ignated by A, B, C). In region C, the only steady state is the trivial solutionnˆ ≡0, which is stable in region C and unstable in regions A and B. The representation of the nontrivial steady state solution (17) is valid in region B, while the representation (19) is valid in region A. On the curve (16) separating regions A and B, the representation (21) is valid
2.2. Examples
As an example of the steady state solution (17) in region B, we take ˆ
µ=3, σˆ =4 ⇒ k=1.46227 (from Eq. (18)) The corresponding steady state solution (17) is plotted in Fig. 4.
As an example of the steady state solution (19) in region A, we take ˆ
µ=1, σˆ =3 ⇒ c=0.6011 (from Eq. (20)) The corresponding steady state solution (19) is plotted in Fig. 4.
As an example of the steady state solution (21) on the boundary (16) between regions A and B, we take
ˆ
µ=1.5, σˆ =4
The corresponding steady state solution (21) is plotted in Fig. 4.
These solutions have been checked by numerically integrating the partial differ- ential-integral equation (11) with the boundary condition (12). These numerical integrations show that these nontrivial steady states are stable because they were approached astˆ→ ∞for all initial conditions tried.
Fig. 4. Steady state equilibrium solutions. The lower curve is Eq. (17) forµˆ = 3,σˆ =4 which lies in region B. The middle curve is Eq. (19) forµˆ =1,σˆ =3 which lies in region A. The upper curve is Eq. (21) forµˆ =1.5,σˆ =4 which lies on the boundary (16) between regions A and B. Note that while the analytical solutions (17), (19) and (21) have different forms, they all have the same general appearance when plotted
2.3. Transient behavior
The transient behavior exhibited by the model is shown here by numerically inte- grating the partial differential-integral equation (11) with the boundary condition (12) for a variety of parameters and initial conditions and displaying plots of nˆ versusDfor a time sequence. Fig. 5 corresponds to the parametersµˆ =3,σˆ =4, which lies in region B (see Fig. 3), for the initial conditionn(0)ˆ =0.1(1−D).
This scenario corresponds to an environment in which there is initially a very small population compared to the steady state. The growth is monotonic.
Figure 6 corresponds to the same parameters but a different initial condition, i.e.,n(ˆ 0)=12/(1+D)10. In this example, we have an initial population that is comparable to that of the steady state, consisting of mainly smaller sized plants.
The result is that we see a ‘bump’ in then(D)ˆ curve, which propagates to the right as time increases. This bump corresponds to growth of the plants which were initially present and has been observed in real plant populations (Sprugel, 1976). We also see slower growth occurring due to seeds generated by the bump and other plants.
The final time shown,tˆ=9 is very close to the steady state (which has already been shown as the bottom curve in Fig. 4). Comparing the initial distribution with
Fig. 5. Transient behavior obtained by numerical integration of Eqs. (11), (12). Parameters areµˆ = 3,σˆ = 4 and the initial condition isn(ˆ 0)= 0.1(1−D). The lowest curve in this figure corresponds totˆ=0, and each curve above it corresponds in sequence to times {1,2,· · ·,9}. This scenario corresponds to an environment in which there is initially a very small population compared to the steady state
the steady state, we see that the measure of smaller plants has decreased while that of larger plants has increased.
Figure 7 corresponds to the same initial condition as Fig. 6, and the same value of the mortality coefficientµˆ =3, but to a smaller value of the seeding coefficient
ˆ
σ =2. This system lies in region C of Fig. 3 in which the only steady state is the trivial solutionnˆ ≡0. We see the same bump as in Fig. 6, again corresponding to the growth of the plants which were initially present. However, here the seeding rate is not large enough to sustain a nontrivial steady state, and the population is seen to be decreasing to zero as time gets larger.
3. Discussion
An interesting feature of the model is the absence of a nontrivial steady state in the region C of theµ-ˆ σˆ parameter plane. This ‘tipping point’ phenomenon is to be con- trasted with another plausible scenario in which the total population size at steady state gradually gets smaller as the seeding coefficientσˆ is reduced to zero. Let us
Fig. 6. Transient behavior obtained by numerical integration of Eqs. (11), (12). Parameters areµˆ =3,σˆ =4 and the initial condition isn(ˆ 0)=12/(1+D)10. Each curve corresponds in sequence to times{0,0.1,0.2,· · ·,0.9}. This scenario corresponds to an environment in which there is an initial population which is comparable in measure to the steady state. The bump seen propagating to the right as time increases corresponds to the growth of plants which were initially present. We also see slower growth occurring due to seeds generated by the bump and other plants
designate byNˆT the normalized total steady state population size (see Eq. (9)):
NˆT =
1
0
ˆ
n dD (22)
We show in Appendix II how the steady state solutions obtained above can be used to get expressions forNˆT. For example, in region B we obtain
NˆT = (µˆ2−4k2)sinhk
2kcoshk− ˆµsinhk (23)
By combining Eq. (23) with Eq. (18), we may find the locus of points in theµ-ˆ σˆ parameter plane which correspond to a fixed steady state population sizeNˆT. The result is shown in Fig. 8 forNˆT = 2,4,6,8 and 10. By inspection of Fig. 8 we obtain the following approximate relationship betweenNˆT and parametersµˆ and
ˆ
σ, valid forµ >ˆ 3:
NˆT ≈2(σˆ − ˆµ) (24)
Fig. 7. Transient behavior obtained by numerical integration of Eqs. (11), (12). Parameters areµˆ =3,σˆ =2 and the initial condition isn(ˆ 0)=12/(1+D)10. Each curve corresponds in sequence to times{0,0.1,0.2,· · ·,0.9}. This scenario is similar to that of Fig. 6 except here the seeding rateσˆ is too small to sustain a nontrivial steady state, and the population is seen to be decreasing to zero as time gets larger
Figure 8 illustrates that although nontrivial steady states occur at all points in regions A and B, the steady state population size gets larger as we get further away from the boundary of region C, see Eq. (13).
It is of interest to convert Eq. (24) from normalized quantities back into physical quantities, (see Eqs. (10), (9)), i.e.
NT ≈2
σ−µ
λ
valid for µ >3 (25)
This shows that although increases in both the mortality coefficientµand the shad- ing coefficientλtend to decrease the total steady state population size, they do so in a markedly different manner. An increase inλreducesNT in a continuous fashion, whereas an increase inµhas a tipping point at approximatelyµ=σsuch that any larger values ofµcauseNT to be zero.
The model presented here mimics many of the features of real plant size den- sity distributions. For example, as described in the introduction, prior work using a world-wide data set for the size density distributions of 226 communities consisting of woody plant species (Enquist and Niklas, 2001) indicates that these distributions are well approximated by log-log linear regression curves, taking the general form
Fig. 8. Normalized total steady state population sizeNˆT. Each curve corresponds to the locus of points on whichNˆT =2,4,6,8,or 10 respectively. Cf. Figs. 2, 3
Ni =βDi−α (26)
whereNiis the number of individuals in the basal stem diameter size binDi,βis an allometric constant, often called the Y-intercept of the log-log linear regression curve, and −αis the slope of the regression curve. Similarly, our model can be used to generate log-log size density distributions for all simulated plant popula- tions reaching an equilibrium condition by approximating the distribution equilib- ria (17), (19) by a regression curve of the form of Eq. (26). The result is shown in Figs. 9, 10 forαandβ, respectively.
As was mentioned in the introduction, both theoretical and experimental stud- ies have shown thatαtypically takes on values around 2. E.g., across 226 plant communities, the mean±SE of the slopes of the log-log linear regression curves was found to be−1.75±0.03. A particular value ofαcorresponds in our model to a relationship betweenµˆ andσˆ which plots as a curve in theµ-ˆ σˆ plane. See Fig. 11. Note that the model predicts that higher mortalityµˆ accompanies lower seeding coefficientσˆ.
Our model indicates that the shading coefficientλreduces the total number of individuals in a population, but that this parameter has no effect on the slope−αof a log-log linear regression of the numbers of individuals in size bins versus plant size of a population’s size density distribution.
Fig. 9. This figure shows the result of least squares fitting a straight line to the lognˆ vs.
logDversion of Eq. (17), yielding an approximation of the form lognˆ =logβ−αlogD ornˆ =βD−α. Here we displayαas a function ofµˆ andσˆ. Points in region C are displayed as havingα=0. Cf. Fig. 3
Our model shows that the shape of the size density distribution will change over time as a population grows in size (number of individuals) such that the distribution is skewed initially to the left (toward smaller plant sizes) as a result of seeding but subsequently becomes less skewed (as individuals increase in their average size and migrate to larger bin sizes). This prediction is consistent with the behavior of real plant populations and mixed-species communities. It also sheds light on the supposition that the numerical value of the slope of log-log linear size density distributions decreases as a population or community ages or fails to experience global physical or biological disturbance (Niklas et al., 2003). This supposition is predicated on the fact that larger individuals, which take time to grow to their large size, shade their smaller counterparts and thus directly or indirectly reduce the numbers of smaller individuals in a population or a community.
Fig. 10. This figure shows the result of least squares fitting a straight line to the lognˆ vs.
logDversion of Eq. (17), yielding an approximation of the form lognˆ =logβ−αlogD ornˆ =βD−α. Here we displayβas a function ofµˆ andσˆ. Points in region C are displayed as havingβ=0. Cf. Fig. 3
4. Conclusions
Our model is extremely simple in terms of using only four parameters to generate plant size density distributions, and each of these parameters (,µ,λandσ) math- ematically subsumes a number of biological phenomena, each of which taken in isolation may profoundly influence the shape of real size distributions. Likewise, we have adopted a number of simplifying assumptions. For example, we assume nonlimiting water abundance and soil-nutrient concentrations, and that the mortal- ity, shading and seeding coefficientsµ,λandσ are invariant across the different size categories in each size distribution. Finally, we note that our model treats a population of conspecifics (i.e. a population composed of many plants of the same species) rather than a community composed of many species, each characterized by its own specific biological features (reflected by the numerical values of,µ,λand σ). Future refinements of the model are required to remedy these simplifications.
Fig. 11. This figure shows curves in theµˆ-σˆ plane which correspond to constantαvalues in the regression relationnˆ = βD−α. These curves are level sets of the surface in Fig. 9.
Theαvalues chosen, 1.75 and 2, correspond respectively to experimentally obtained and theoretically derived values
Thus, for example, we may assume thatλscales like a power law ofD, as the capac- ity to shade smaller individuals increases as a function of absolute size. Another refinement involves simulating the size density distributions of mixed-species com- munities to evaluate whether their mathematical behavior accords with the behavior of real communities.
Despite its apparent simplicity, our model obtains a general phenomenology for size distributions that is consistent with empirical observations on the size density distributions of vascular plant communities world-wide. This congruity suggests that the model is biologically robust. For example, our model draws sharp attention to the influence of the relative numerical values of the seeding coefficient and the mortality coefficient,σandµ, respectively. Model behavior for whichσsufficiently exceedsµpredicts stable equilibria for size density distributions, whereas distribu- tions that eventually collapse are obtained whenσis sufficiently small compared to µ(see Eq. (13) or (25)). This result is biologically logical and consistent with obser- vations on the distribution of conspecifics in mixed-species communities, which show that species failing to produce sufficient seeds eventually disappear from such communities.
The model as presented here appears to mimic the general phenomenology observed for real plants and offers a mathematically explicit and heuristic approach to dissecting the complex affects of growth in size, fecundity, mortality, and the capacity to shade smaller individuals on the general shape of a size distribution.
In particular, the model offers a closed form expression for the equilibrium size density distribution, from which explicit algebraic conditions for persistence have been derived. Moreover, when the results of the model are approximated by a log- log regression curve, a given slope (such as -2) has been shown to correspond a specific relationship between mortality and fecundity.
Appendix I: Stability of the trivial solution
In order to investigate the linear stability of trivial solution, we linearize Eqs.
(11), (12) aboutnˆ=0:
∂nˆ
∂t = −∂nˆ
∂D − ˆµnˆ ˆ
n= ˆσ
1
0
ˆ
n dD, D=0
To determine stability, we setn(D, t)ˆ =ertf (D), which gives:
f+(µˆ+r)f =0 (27)
f = ˆσ
1
0
f dD, D=0 (28)
where by we mean ∂/∂D. The transition from stable to unstable occurs when r=0. Using this condition in Eq. (27) gives
f =exp(− ˆµD) (29)
Inserting (29) into the boundary condition (28), we get
1= ˆσ
1−e− ˆµ ˆ µ
Thus, the trivial solutionnˆ ≡0 is stable if and only if ˆ
σ < µˆ
1−e− ˆµ (13)
and unstable otherwise.
Appendix II: Steady state equilibrium
The equations governing the steady state equilibrium in the model (11), (12) are:
dnˆ dD+ ˆn
1
D
ˆ
n dD+ ˆµnˆ =0 (14)
ˆ n= ˆσ
1
0
ˆ
n dD, D=0 (15)
here reported for convenience. These equations may be solved by setting z(D)=
1
D
ˆ
n dD ⇒ n(D)ˆ = −dz
dD (30)
which gives the following ordinary differential equation:
z+zz+ ˆµz=0 (31) with the boundary conditions
z=0 when D=1 (32)
z+ ˆσ z=0 when D=0 (33) Equation (31) can be solved to satisfy the boundary condition (32):
z= (µˆ2−4k2)sinh(1−D)k
2kcosh(1−D)k− ˆµsinh(1−D)k (34) wherekis an arbitrary constant. Substituting (34) into the boundary condition (33) gives the relation (18) between k andµˆ andσˆ. The corresponding steady state expression (17) is obtained by differentiating Eq. (34), while Eqs. (19), (20) are obtained from (17), (18) by settingk=ic, respectively.
Each of the two representations (17), (19) is valid in a distinct region of the ˆ
µ-σˆ parameter plane. To find these regions, we note that the two representations are separated by the valuek =0. The curve separating the two regions is there- fore obtained by settingk =0 in Eq. (18), which gives Eq. (16). Straightforward computations show that the corresponding steady state solution takes the limiting form of Eq. (21) ask→0.
The Eqs. (19), (20) simplify in the important case of no mortality,µˆ =0:
ˆ
n(D)=2k2sec2k(1−D), µˆ =0 wherekis related toσˆ by:
ˆ
σ =2kcsc 2k, µˆ =0
The normalized total steady state population sizeNˆT (see Eq. (22)) is related to the quantityz(D)defined in Eq. (30):
NˆT =z(0)
In region B, for example, we may use Eq. (34) to easily obtain expression (23).
References
Antonovics, J., Levin, D.A.: The ecological and genetic consequences of density-dependent regulation in plants. Ann. Rev. Ecol. Syst. 11, 411–452 (1980)
Bishir, J., Namkoong, G.: Density-dependent dynamics in size structured plant populations.
J. Theor. Biol. 154, 163–188 (1992)
Bolker, B.M., Pacala, S.W., Neuhauser, C.: Spatial dynamics in model plant communities:
What do we really know? The American Naturalist 162, 135–148 (2003)
Chew, R.M., Chew, A.E.: The primary productivity of a desert shrub (Larrea tridentata).
Ecological Monographs 35, 355–375 (1965)
Crawley, M.J.: Insect herbivores and plant population dynamics. Ann. Rev. Entomol. 34, 531–564 (1989)
Crawley, M.J.: The population dynamics of plants. Proceedings of the Royal Society of London B 330, 125–140 (1990)
Cushing, J.M.: An Introduction to Structured Population Dynamics. SIAM, Philadelphia, cbms-nsf regional conference series in applied mathematics, 71 edition, 1998
Cz´ar´an, T., Bartha, S.: Spatiotemporal dynamic models of plant populations and communi- ties. Trends in Ecology & Evolution 7, 38–42 (1992)
De Angelis, D.L., Gross, L.J., eds.: Individual-Based Models and Approaches in Ecology:
Populations, Communities and Ecosystems. Chapman & Hall, New York, 1992 Diekmann, O., Gyllenberg, M., Metz, J.A.J., Thieme, H.R.: On the formulation and analysis
of general deterministic structured population models I. Linear theory. J. Math. Biol. 36, 349–388 (1998)
Diekmann, O., Gyllenberg, M., Huang, H., Kirkilionis, M., Metz, J. A.J., Thieme, H.R.: On the formulation and analysis of general deterministic structured population models II.
Nonlinear theory. J. Math. Biol. 43, 157–189 (2001)
Durrett, R., Levin., S.A.: The importance of being discrete (and spatial). Theoretical Popu- lation Biology 46, 363–394 (1994)
Enquist, B.J., Niklas, K.J.: Invariant scaling relations across tree-dominated communities.
Nature 410, 655–660 (2001)
Enquist, B.J., Niklas, K.J.: Global allocation rules for patterns of biomass partitioning in seed plants. Science 295, 1517–1520 (2002)
Enquist, B.J., Brown, J.H., West, G.B.: Allometric scaling of plant energetics and population density. Nature 395, 163–165 (1998)
Enquist, B.J., West, G.B., Charnov, E.L., Brown, J.H.: Allometric scaling of production and life- history variation in vascular plants. Nature 401, 907–911 (1999)
Gurtin, M.E., MacCamy, R.C.: Nonlinear age-dependent population dynamics. Archive for Rational Mechanics and Analysis 54, 281–300 (1974)
Hara, T.: Dynamics of stand structure in plant monocultures. J. Theor. Biol. 110, 223–239 (1984a)
Hara, T.: A stochastic model and the moment dynamics of the growth and size distribution in plant populations. J. Theor. Biol. 109, 173–190 (1984b)
Hara, T.: A model for mortality in a self-thinning plant population. Annals of Botany 55, 667–674 (1985)
Hara, T.: Effects of the mode of competition of stationary size distribution in plant popula- tions. Annals of Botany 69, 509–513 (1992)
Hara, T., Wyszomirski, T. Competitive asymmetry reduces spatial effects on size-structure dynamics in plant populations. Annals of Botany 73, 285–297 (1994)
Harper, J.L.: Population Biology of Plants. Academic Press, New York, 1977
Hoppensteadt, F.: An age dependent epidemic model. J. Franklin Institute Journal 297, 325–
333 (1974)
Kinzig, A.P., Levin, S.A., Dushoff, J., Pacala, S. Limiting similarity, species packing, and system stability for hierarchical competition-colonization models. The American Natu- ralist 153, 371–383 (1999)
Kirkpatrick, M.: Demographic models based on size, not age, for organisms with indeter- minate growth. Ecology 65, 1874–1884 (1984)
Kohyama, T.: Simulating stationary size distribution of trees in rain forest. Annals of Botany 68, 173–180 (1991)
Kohyama, T.: Size-structured multi-species model of rain forest trees. Functional Ecology 6, 206–212 (1992)
Lefkovitch, L.P.: The study of population growth in organisms grouped by stages. Biometrics 21, 1–18 (1965)
Leslie, P.H.: On the use of matrices in certain population mathematics. Biometrika 33, 183–
212 (1945)
Levin, S.A., Paine, R.T.: Disturbance, patch formation and community structure. Proceed- ings of the National Academy of Science USA 71, 2744–2747 (1974)
McKendrick: Applications of mathematics to medical problems. Proceedings of the Edin- burgh Mathematical Society 44, 98–130 (1926)
McMahon, T.: Size and shape in biology. Science 179, 1201–1204 (1973)
Metz, J.A.J., Diekmann, O.: The Dynamics of Physiologically Structured Populations.
Springer Verlag, Berlin, lecture notes in biomathematics, 68 edition, 1986
Nagano, M.: Dynamics of stand development. In: T. Kira, Y. Ono, T. Osokawa (eds.), Bio- logical Production in a Warm-Temperate Evergreen Oak Forest of Japan, volume 18, University of Tokyo Press, JIBP Synthesis edition, 1978, pp. 21–31
Niklas, K.J.: Plant Allometry. The University of Chicago Press, Chicago, 1994
Niklas, K.J.: Size- and age-dependent variations in the properties of sap- and heartwood in black locust (Robinia pseudoacacia L.). Annals of Botany 79, 473–478 (1997)
Niklas, K., Midgley, J.J., Rand, R.H.: Tree size frequency distributions, plant density, age and community disturbance. Ecology Letters 6, 405–411 (2003)
Pacala, S.W., Silander, J.A.: Neighborhood models of plant-population dynamics. 1. Single- species models of annuals. The American Naturalist 125, 385–411 (1985)
Paine, R.T., Levin, S.A.: Intertidal landscapes: Disturbance and the dynamics of pattern.
Ecological Monographs 51, 145–178 (1981)
Petersen, T.D.: Effects of interference from Calamagrostis Rubescens on size distributions in stands of Pinus Ponderosa. J. Political Economy 25, 265–272 (1988)
Petersen, T.D., Ning, Z., Newton, M.: Dynamics of size structure in seedling stands to Frax- inus Mandshurica in Northeast China. Annals of Botany 66, 255–263 (1990)
Rees, M., Crawley, M.J.: Growth, reproduction and population dynamics. Functional Ecol- ogy 3, 645–653 (1989)
Schaffer, W.M., Leigh, E.G.: The prospective role of mathematical theory in plant ecology.
Systematic Botany 1, 209–232 (1976)
Sharpe, F.R., Lotka, A.J.: A problem in age distribution. Philosophical Magazine 21, 435–438 (1911)
Shinko, J.W., Streifer, W.: A new model for age-size structure of a population. Ecology 48, 910–918 (1967)
Sprugel, D.G.: Dynamic structure of wave-generated Abies balsamea forests in the north- eastern United States. J. Ecol. 64, 889–911 (1976)
Takada, T., Hara, T.: The relationship between the transition matrix model and the diffusion model. J. Math. Biol. 32, 789–807 (1994)
Takada, T., Iwasa, Y.: Size distribution dynamics of plants with interaction by shading. Eco- logical Modelling 33, 173–184 (1986)
Trucco, E.: Mathematical models for cellular systems the von Foerster equation. Part I.
Bulletin of Mathematical Biophysics 27, 285–304 (1965)
von Foerster, H.: Some remarks on changing populations. In: The Kinetics of Cellular Pro- liferation, Grune and Stratton, 1959, pp. 382–407
Webb, G.F.: Theory of Nonlinear Age-Dependent population dynamics. Marcel Dekker, New York, 1985
Weiner, J.: Asymmetric competition in plant populations. Trends in Ecology & Evolution 5, 360–364 (1990)
West, G.B., Brown, J.H., Enquist, B.J.: A general model for the origin of allometric scaling laws in biology. Science 276, 122–126 (1997)
West, G.B., Brown, J.H., Enquist, B.J.: The fourth dimension of life: Fractal geometry and allometric scaling of organisms. Science 284, 1677–1679 (1999a)
West, G.B., Brown, J.H., Enquist, B.J.: A general model for the structure and allometry of plant vascular systems. Nature 400, 664–667 (1999b)
West, P.W., Jackett, D.R., Borough, C.J. Competitive processes in a monoculture of Pinus Raidata D. Don. Oecologia 81, 57–61 (1989)
Yoda, K., Kira, T., Ogawa, H., Hozumi, K.: Self-thinning in overcrowded pure stands under cultivated and natural conditions. Journal of Biology Osaka City University 14, 107–129 (1963)
