Journal of Theoretical Biology 363 (2014) 129–133
Contents lists available at ScienceDirect
Journal of Theoretical Biology journal homepage: www.elsevier.com/locate/yjtbi
Species distributions and area relationships C.J. Thompson a, T.E. Lee b,c,n, M.A. McCarthy b a b c
Department of Mathematics and Statistics, University of Melbourne, VIC 3010, Australia School of Botany, University of Melbourne, VIC 3010, Australia Biomathematics Unit, Department of Zoology, Faculty of Life Sciences, Tel-Aviv University, P.O. Box 39040, Tel-Aviv 69978, Israel
H I G H L I G H T S
Show a new, concise, derivation of Yule's equilibrium distribution. Derive the species–area relationship based on the evolution of species through mutations. Calculate the exponent of the species–area relationship for four different data sets. Confirm that the exponent term, based on mutation rates, agrees with previously observed values.
art ic l e i nf o
a b s t r a c t
Article history: Received 9 October 2013 Received in revised form 4 August 2014 Accepted 5 August 2014 Available online 20 August 2014
The well-known species–area relationship is one of many scaling laws, or allometries, in ecology and biology that have received much attention over the years. We present a new derivation of this relationship based on Yule's theory of evolution of species. Using definitions of mutation rates, our analysis yields species–area exponents that are in close agreement with previously observed values. & 2014 Elsevier Ltd. All rights reserved.
Keywords: Yule distribution Mutations Scaling exponents
1. Introduction The species–area relationship (SAR) S ¼ cAz ;
ð1Þ
dates back to at least 1921 when Arrhenius (1921) proposed Eq. (1) as an empirical formula for the number of species S of a given taxonomic group (or genus) found in a region of area A. Alternative SARS and their derivations, which are often based on assumed underlying species-abundance distributions (SADs), have been considered and debated over the years (May, 1975; May and Stumpf, 2000; McGuinness, 1984; Chisholm, 2007). While there has been no general consensus on the range of validity of particular SARs (for varying A), recent empirical and theoretical studies by Harte et al. (2009), O'Dwyer and Green (2010) and Storch et al. (2012) strongly suggest that Eq. (1) is not valid for ‘small A’, but that in certain situations, may be valid asymptotically n Corresponding author at: Mathematical Institute, Andrew Wiles Building, University of Oxford, Oxford, OX2 6GG. Tel.: þ 44 1865 611511. E-mail address: [email protected] (T.E. Lee).
http://dx.doi.org/10.1016/j.jtbi.2014.08.011 0022-5193/& 2014 Elsevier Ltd. All rights reserved.
for larger values of A. Our purpose here is not to discuss or comment on the large body of work that has been published on this problem (see May, 1975; May and Stumpf, 2000; McGuinness, 1984; Chisholm, 2007; Harte et al., 2009; O'Dwyer and Green, 2010; Storch et al., 2012 and references quoted therein), but rather to present a new derivation of Eq. (1) which is valid asymptotically for large A. Our derivation of Eq. (1) is based on ideas published over ninety years ago by Yule (1925) on the evolution of species through specific and generic mutations within and between genera of species. In the following section we present a stochastic dynamical systems model based on Yule's original ideas. In Section 3 we show that the species–area relation (1) follows from the asymptotic form of Yule's equilibrium distribution for a large number of species in a genus. Moreover, our asymptotic derivation of (1) provides an explicit expression for the exponent z in Eq. (1) in terms of Yule's species and generic mutation rates. In Section 4 we present some case studies as examples, showing close agreement between our formulae for z and the previous empirical studies based on Eq. (1). Our results are summarised and discussed in the final section.
130
C.J. Thompson et al. / Journal of Theoretical Biology 363 (2014) 129–133
2. The Yule distribution In his original work Yule (1925) proposed evolution of species and genera through two kinds of mutations which he described as follows:
Within any species, in any time interval, a ‘specific mutation’
can occur yielding a new species, but within the same genus as the parent. Within any genus, in any time interval, a ‘generic mutation’ can occur yielding a new species so different from the parent that it is placed in a new genus.
We assume (with Yule) that the rate s of specific mutations is a constant, and the same for all species within a genus; and that the rate g of generic mutations is a constant and same for all genera. Notice that at any time, the number of a species within a genera increases due to both s and the number of species within the genera. That is, the rich get richer (also known as ‘preferential attachment’). Yule (1925) was at pains to point out that these assumptions are simplistic and were designed solely to elicit quantitative comparisons with well-known facts. He provides the following example to demonstrate the last of his simplifying assumptions: if A, B, C, and D are the existing genera, and one of them throws a generic mutation, it is assumed that this will represent a new genus E. That is, the possibility that the new species may be classed under an existing genera is ignored. Combine this with the absence of a death process, and it becomes clear that the number of genera with few species could easily be overestimated. This overestimation has recently been addressed by the Birth–DeathMutation (BDM) Process (Maruvka et al., 2011). Despite these assumptions, we will see in Section 4 that Yule's comparisons with Willis' data proved to be successful. Not surprisingly, Yule's analysis (in the early 1920s) was long and complicated. Our main purpose in this section is to present a simple dynamical systems model for P n ðtÞ ¼ the probability that the genus has n species at time t:
ð2Þ
and to show that Yule's (equilibrium) distribution is a steady state solution of this system. Consider first the case n ¼1. In continuous time it follows from Yule's assumptions that dP 1 ðtÞ ¼ g ðs þ gÞP 1 ðtÞ; dt
ð3Þ
where the first (input) term on the right-hand side of (3) represents a generic mutation from another genus, and the second term represents losses from the single species in the (new) genus through specific and generic mutations. Similarly, for an established genus with at least two species, it follows from Yule's assumptions that dP n ðtÞ ¼ ðn 1ÞsP n 1 ðtÞ ðnsþ gÞP n ðtÞ; dt
n ¼ 2; 3; …:
ð4Þ
In Eq. (4) the first term arises from the n 1 possible specific mutations from n 1 species in the genus, and the second term represents losses from n species in the genus due to n specific mutations, and a constant loss rate g due to generic mutations. There are of course many possible dynamical models for Pn(t). We again note that Yule's scheme ignores death, and therefore extinction. As a consequence, and as mentioned by Yule (1925, p.38), ‘on our assumptions, the mean size of a genus after infinite time must itself be infinite’. Extinction could of course be included in our system, by adding appropriate P n þ 1 ðtÞ terms to the righthand side of (4) for example. We will return to these issues in the
final section. The important point to stress here is that Eqs. (3) and (4) although simplistic are biologically reasonable and logically consistent in the sense that on summation of Eqs. (3) and (4) we deduce that 1 d 1 ∑ P n ðtÞ ¼ g 1 ∑ P n ðtÞ ; ð5Þ dt n ¼ 1 n¼1 as required by conservation of probabilities, i.e., the initial condition 1
∑ P n ð0Þ ¼ 1
1
∑ P n ðtÞ ¼ 1;
implies
n¼1
n¼1
for all t 4 0. In the Appendix we present an exact solution of Eqs. (3) and (4) and show that in the limit t-1 Pn(t) approaches a globally stable equilibrium given by steady-state solution P nn obtained by setting the left hand sides of (3) and (4) to zero, i.e., P n1 ¼
g ¼ ð1 þ σ Þ 1 ; s þg
ð6Þ
and ðnσ þ 1ÞP nn ¼ ðn 1Þσ P nn 1 ;
n ¼ 2; 3; …:
ð7Þ
where s g
σ¼ ;
ð8Þ
is the ratio of the specific and generic mutation rates which we assume henceforth (with Yule) to be larger than unity. Iterating Eq. (7) we deduce (as shown in the Appendix) that P nn ¼
Γ ðnÞΓ ð1 þ 1=σ Þ ; σΓ ðn þ1 þ 1=σ Þ
n ¼ 1; 2; …;
ð9Þ
where Γ ðkÞ is the gamma function. Using Stirling's formula for Γ ðkÞ when k is large, we arrive at the asymptotic form (as shown in the Appendix) 1 1 ð1 þ 1=σ Þ as n-1: ð10Þ n P nn Γ 1 þ
σ
σ
The distribution P nn equation (9) was obtained by Yule in a rather lengthy and complicated derivation. This is not surprising since his work predated developments in stochastic processes and dynamical systems theory. His notion of equilibrium was also at odds with modern interpretations in terms of steady states. These shortcomings were in fact pointed out by Simon (1955) some thirty years after Yule's paper. Simon re-derived Yule's results, and Eqs. (9) and (10) in particular, in a more general and contemporary setting with interesting applications to word frequencies, city sizes, income distributions and frequencies of scientific publications. Simon (1955) correctly referred to (9) as the (equilibrium) Yule distribution, although that citation has been largely forgotten in recent times. We note in passing that although the P nn , from Eqs. (6) and (7), sum to unity as required, there is no reason for the asymptotic form values (10) to do likewise (see the Appendix). Nevertheless, the Yule distribution has some unusual and even ‘paradoxical’ properties (as stated by Yule, 1925, p. 38). In particular it follows from (9) and (10) that the (equilibrium) mean number of species in a genus, obtained by multiplying P nn , Eq. (9), by n and summing on n from 1 to 1, diverges when σ 4 1 by virtue of the asymptotic form (10) for the tail distribution P nn . In the following section we show that the asymptotic form of Yule's distribution (10) can be used to derive the species–area relationship (1) (for large A) with an explicit expression for the exponent z in terms of Yule's parameter σ defined in (8).
C.J. Thompson et al. / Journal of Theoretical Biology 363 (2014) 129–133
131
3. Derivation of the species–area relationship We consider a particular bioregion and assume that the region (with area A) can only support a maximum number N of species in a given genus. Clearly, this maximum occurs when there is one individual of each species in the genus in the region. While this situation is unstable in a biological sense (due to births and deaths of individuals), it is theoretically valid in a mathematical sense. Thus to calculate the expected number S of species in the region, it is important to note that the actual number of species n in the region is a random variable whose probability distribution must be truncated at the theoretical maximum possible value of N (as discussed in the Appendix). In our model, it is then reasonable to interpret N as the maximum possible number of individuals in a genus that could occupy the given region. With this interpretation, it is then not unreasonable to assume (May, 1975, Eqs. (1.1) and (1.5), for example) that N in the following equations is proportional to A. In the Appendix we show that the expected steady-state number of species (of a given genus) in the region is given asymptotically for large N by N
S ∑ nP nn ;
ð11Þ
n¼1
where P nn is the Yule distribution (9). In order to evaluate the asymptotic form for S on the right hand side of Eq. (11), we multiply (7) on both sides by n, and sum from n¼ 2 to n¼ N. The sum on the left hand side can be expressed as N
N
N
N
n¼2
n¼2
n¼1
n¼1
∑ n2 σ P nn þ ∑ nP nn ¼ ∑ n2 σ P nn ð1 þ σ ÞP n1 þ ∑ nP nn ;
N
N
n¼2
n¼1
ð13Þ
Equating (12) and (13) we see that the n2 terms cancel, and after some elementary algebra we obtain the identity N
∑ nP nn ≔
n¼1
1
σ 1
NðN þ1Þσ P nN 1 ;
ð14Þ
from (11) and (6). The identity (14) is valid for any N Z 1. For large N we deduce from (10), (11) and (14) that the expected number S of species is σ 1 1 S N 2 P nN Γ 1 þ N 1 1=σ as N-1: ð15Þ σ 1 σ 1 σ Finally, to arrive at the species–area relationship (1) we use the assumption that N is proportional to A, as discussed above, to obtain Eq. (1) valid asymptotically for large A with z ¼ 1 1=σ
(1925) are reproduced in Fig. 1 as a log–log plot, along with more recent data from Galapagos plants (Alan Tye, 2005), Yorkshire plants (Stace et al., 2003) and New Britain birds (World Wildlife Fund, 2006). As shown in Fig. 1, these plots reveal power laws of the Yule form (10). Values of σ derived from linear regression of the log–log plots of species-genera data, as shown in Fig. 1, were used to calculate the ‘inferred z’ values given in Table 1, evaluated from Eq. (16). The ‘observed z’ values given in Table 1 on the other hand used linear regression on log–log plots of species–area field data, assuming Eq. (1), to obtain estimates of z. The 95% confidence intervals (CI) for the inferred z values in Table 1 were obtained from regression lines for samples of (two or more) data points for each genus in Fig. 1.
ð12Þ
while the sum on the right hand side, after changing summation variables, can be expressed as ∑ nðn 1Þσ P nn 1 ¼ ∑ nðn þ 1Þσ P nn NðN þ 1Þσ P nN :
Fig. 1. Log–log plots of number of genera versus number of species per genus for four different taxa. Regression lines were fitted to the data. Slopes (of (1 þ 1/σ) from Eq. (10)) are respectively lizards: 1.69, G. plants: 1.68, Y. plants: 1.82 and birds: 1.87 with corresponding r2 values of 0.998, 0.990, 0.996 and 0.983.
ð16Þ
and c a constant which depends on σ. This should be considered an upper bound on the slope of the SAR because, in general, the maximum number of species per genus N will scale sub-linearly with area A.
4. Examples As suggested by the title of the original paper, Yule, 1925 had copious amounts of data to examine courtesy of Dr J.C. Willis FRS. Yules regression analysis of log–log plots (of genera frequency and numbers of species in a genus) of beetles, snakes and lizards in particular confirmed the asymptotic form of his equilibrium distribution (10). Willis's data for lizards extracted from Yule
5. Summary In this paper we have presented a new derivation of the wellknown species–area scaling relation (1). Our derivation is based on ideas published over ninety years ago by Yule (1925) on the evolution of species through specific and generic mutations of species and genera. Yule's derivation of the equilibrium species distribution (10) was long and complicated, and his interpretation of ‘equilibrium’ was at odds with modern developments in stochastic processes and dynamical systems theory. These shortcomings were partly overcome some thirty years after Yule's work by Simon (1955). Simon rederived Yule's equilibrium, or steadystate distribution, in a more general setting, with applications to word frequencies, city sizes, income distributions and frequencies of scientific publications. As noted in Section 2, the Yule distribution is somewhat paradoxical. Hence, power-law species–area relationships (1) may also be an exceptional case, as discussed in May (1975), Harte et al. (2009), O'Dwyer and Green (2010) and Storch et al. (2012). In this paper we have also re-derived Yule's results from a simple, recursive, dynamical systems model. We show in particular that Yule's equilibrium, or steady-state distribution, is in fact a globally attracting equilibrium state (i.e., irrespective of initial conditions) of our recursive, dynamical systems model (7). The details of this are provided in the Appendix. We consider a bioregion that can support a maximum number of species N in one genus, and derive an exact identity, Eq. (14), relating the expected number of species S in the region, and the Yule distribution P nN , Eq. (9). The assumption that N is proportional to the area A of the bioregion, which we discussed in some detail in Section 3, leads to the species–area relation (1) with an explicit expression for the exponent z (16) in terms of the ratio σ of Yule's specific and generic mutation rates. Note that the slope calculated from z is an upper bound since the maximum number of species per genus N will scale sub-linearly with area A.
132
C.J. Thompson et al. / Journal of Theoretical Biology 363 (2014) 129–133
Table 1 Values and confidence intervals (CI) for the exponent z in the species–area relationship (1) ‘inferred’ from (16) (z ¼1–1/σ) and the linear regression slopes given in Fig. 1, and previously ‘observed’ values obtained from species–area field data. Taxon
Location
Observed z
Inferred z
CI
Plants Plants
Yorkshire Galapagos
0.18 0.32
0.15–0.21 0.29–0.35
Lizards Birds
West Indies New Britain
0.21 (Usher (1973)) 0.31–0.33 (Preston, 1962; Hamilton et al., 1963; Johnson and Raven, 1973) 0.30 (Darlington, 1957) 0.18 (Diamond, 1975)
0.31 0.13
0.29–0.33 0.06–0.20
Finally we considered frequency data for the occurrence of species and genera for four different taxa (plants, birds and lizards). Log–log plots confirmed the general asymptotic form of Yule's equilibrium distribution (10) and provided accurate estimates, through linear regression, of Yule's ratio σ. Inferred values for z ¼ 1 1=σ for the four taxa given in Table 1 are in accord with previously published ‘observed z’ values (Usher, 1973; Preston, 1962; Hamilton et al., 1963; Johnson and Raven, 1973; Darlington, 1957; Diamond, 1975) derived from species–area field data. If one includes extinction in the Yule scheme, i.e., maintaining generic mutations (with rate g a 0), one also obtains the SAR equation (1) asymptotically for large A, but with σ in Eq. (16) replaced by σ 0 ¼ ðs μÞ=g, where μ is a specific extinction rate. The details of this analysis, and alternative ways of treating generic mutations, will be presented in subsequent publications.
the exponential terms exp½ ðlsþ gÞt, l ¼ 1; 2; …n. It follows that the steady state solution P nn Eq. (7) is in fact a globally attracting equilibrium state of the system of Eqs. (3) and (4). In order to obtain an asymptotic estimate for P nn when n is large (Eq. (10)), we iterate the recursion relation Eq. (7) for P nn , noting Eq. (A.3) for P n1 , to obtain P nn ¼
ðn 1Þσ ðn 2Þσ σ 1 1 ⋯ ¼ 1 þ nσ 1 þ ðn 1Þσ 1 þ 2σ 1 þ σ σ
Γ ðnÞΓ ð1 þ 1=σ Þ ; Γ ðn þ 1 þ 1=σ Þ
ðA:7Þ
where Γ ðkÞ is the gamma function defined for any (real part of) k 4 1 by Z 1 Γ ðk þ 1Þ ¼ xk e x dx ¼ kΓ ðkÞ: ðA:8Þ 0
Acknowledgments We are grateful to Robert May for his comment on a previous version of this paper. We also thank Ryan Chisholm and an anonymous referee for valuable suggestions and comments on revised versions of this paper.
Appendix A. Derivation of Eqs. (9)–(11) In order to solve the system of equations (3) and (4) for Pn(t), n ¼ 1; 2; … we begin with (3) for P 1 ðtÞ, i.e., dP 1 ¼ g ðs þ gÞP 1 ðtÞ; dt
ðA:1Þ
with general solution P 1 ðtÞ ¼
g þc1 e ðs þ gÞt sþg
ðA:2Þ
where c1 is a constant (determined by the initial value P 1 ð0Þ). Irrespective of the value of c1 we have Eq. (6), i.e., P n1 ¼
1 ¼ lim P 1 ðtÞ; 1 þ σ t-1
s g
σ¼ :
ðA:3Þ
Similarly the differential equation for P 2 ðtÞ (4) is dP 2 ¼ sP 1 ðtÞ ð2s þ gÞP 2 ðtÞ: dt
ðA:4Þ
Substituting (A.2) into (A.4) it is then straightforward to show that the general solution of (A.4) is n
P 2 ðtÞ ¼ P 2 þ c1 e
ðs þ gÞt
þ c2 e
ð2s þ gÞt
;
where c1 and c2 are constants and σ P n2 ¼ P n ¼ lim P 2 ðtÞ; 1 þ 2σ 1 t-1
For large (real) k we have the well-known Stirling's (asymptotic) formula (Abramowitz and Stegun, 1964), pffiffiffiffiffiffiffiffiffi Γ ðk þ 1Þ kk e k 2π k: ðA:9Þ In particular for large n (compared with 1/σ) pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Γ ðnÞ ðn 1Þn 1 e ðn 1Þ 2π ðn 1Þ pffiffiffiffiffiffiffiffiffi 1 ; nn 1 e n 2π n 1 þ O n and similarly 1 1 n þ 1=σ ðn þ 1=σ Þ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Γ nþ þ1 nþ e 2π ðn þ 1=σ Þ σ σ pffiffiffiffiffiffiffiffiffi 1 nn þ 1=σ e n 2π n 1 þ O ; n
ðA:10Þ
ðA:11Þ
where in deriving the second lines of Eqs. (A.10) and (A.11) we have used the fact that for large n x n 1þ ex ; ðA:12Þ n (with x ¼ 1 in Eq. (A.10) and x¼1/σ in Eq. (A.11)). Substituting (A.10) and (A.11) into Eq. (A.7) gives the asymptotic form Eq. (10) for the Yule distribution equation (A.7), and BYU. Finally, by summing (6) and (7), we deduce that for finite N N
∑ P nn ¼ 1 σ NP nN cN 1
n¼1
as N-1:
ðA:13Þ
ðA:5Þ
The truncated distribution (appropriate for a finite region) ( cN 1 P nn ; 1 r n r N P nn0 ¼ ðA:14Þ 0; n 4 N;
ðA:6Þ
is then properly normalised and the mean-number of species (in the region) is given by
is the steady-state solution of (A.4) in agreement with Eq. (7). Proceeding by induction it is then not difficult to show that the general solution of the differential equation (4) for Pn(t) n Z 3 is the sum of P nn , defined recursively by Eq. (7) and a linear combination of
1
N
n¼1
n¼1
S ∑ nP nn0 ∑ nP nn ; as in Eq. (11).
ðA:15Þ
C.J. Thompson et al. / Journal of Theoretical Biology 363 (2014) 129–133
References Abramowitz, M., Stegun, I.A., 1964. Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables. Data provided by Dr. Alan Tye, Head, Botany Department, Charles Darwin Research Station, Isla Santa Cruz, Galapagos, Ecuador. Arrhenius, O., 1921. Species and Area. J. Ecol. 9, 95–99. Chisholm, R.A., 2007. Sampling species abundance distributions: resolving the veil-line debate. J. Theor. Biol. 247, 600–607. Darlington, P.J., 1957. Zoogeography. John Wiley and Sons, New York. Diamond, J.M., 1975. Assembly of species communities. In: Cody, M.L., Diamond, J.M. (Eds.), Ecology and Evolution of Communities. Belknap/Harvard University Press, Cambridge, MA, pp. 342–444. Hamilton, T.H., Barth, R.H., Rubinoff, I., 1963. Species abundance: natural regulation of insular variation. Science 143, 1575–1577. Harte, J., Smith, A.B., Storch, D., 2009. Biodiversity scales from plots to biomes with a universal species–area curve. Ecol. Lett. 12, 789–797. Johnson, M.P., Raven, P.H., 1973. Species number and endemism: the Galapagos Archipelago revisited. Science 179, 893–895. Maruvka, Y.E., Kessler, D.A., Shnerb, N.M., 2011. The birth–death-mutation process: a new paradigm for fat tailed distributions. PloS One 6 (11), e26480. May, R.M., 1975. Patterns of species abundance and diversity. Ecol. Evol. Communities, 81–120.
133
May, R.M., Stumpf, M.P., 2000. Species–area relations in tropical forests. Science 290, 2084–2086. McGuinness, K.A., 1984. Species–area curves. Biol. Rev. 59, 423–440. O'Dwyer, J.P., Green, J.L., 2010. Field theory for biogeography: a spatially explicit model for predicting patterns of biodiversity. Ecol. Lett. 13, 87–95. Preston, F.W., 1962. The canonical distribution of commonness and rarity. Part I. Ecology 43, pp. 185–215 and 431–432. Simon, H.A., 1955. On a class of skew distribution functions. Biometrika 42, 425–440. Stace, C.A., Ellis, R.G., Kent, D.H., McCosh, D.J. (Eds.), 2003. Vice-County Census Catalogue of the Vascular Plants of Great Britain the Isle of Man and the Channel Islands. Botanical Society of the British Isles, London. Storch, D., Keil, P., Jetz, W., 2012. Universal species–area and endemics–area relationships at continental scales. Nature 488, 78–81. Usher, M.A.B., 1973. Biological Management and Conservation: Ecological Theory Applications and Planning. Chapman and Hall, London. World Wildlife Fund. WildFinder: Online Database of Species Distributions, ver. January 2006. 〈www.worldwildlife.org/WildFinder〉. Yule, G.U., 1925. A mathematical theory of evolution, based on the conclusions of Dr. JC Willis FRS. Philos. Trans. R. Soc. Lond. Ser. B, Contain. Pap. Biol. Charact. 213, 21–87.
Contents lists available at ScienceDirect
Journal of Theoretical Biology journal homepage: www.elsevier.com/locate/yjtbi
Species distributions and area relationships C.J. Thompson a, T.E. Lee b,c,n, M.A. McCarthy b a b c
Department of Mathematics and Statistics, University of Melbourne, VIC 3010, Australia School of Botany, University of Melbourne, VIC 3010, Australia Biomathematics Unit, Department of Zoology, Faculty of Life Sciences, Tel-Aviv University, P.O. Box 39040, Tel-Aviv 69978, Israel
H I G H L I G H T S
Show a new, concise, derivation of Yule's equilibrium distribution. Derive the species–area relationship based on the evolution of species through mutations. Calculate the exponent of the species–area relationship for four different data sets. Confirm that the exponent term, based on mutation rates, agrees with previously observed values.
art ic l e i nf o
a b s t r a c t
Article history: Received 9 October 2013 Received in revised form 4 August 2014 Accepted 5 August 2014 Available online 20 August 2014
The well-known species–area relationship is one of many scaling laws, or allometries, in ecology and biology that have received much attention over the years. We present a new derivation of this relationship based on Yule's theory of evolution of species. Using definitions of mutation rates, our analysis yields species–area exponents that are in close agreement with previously observed values. & 2014 Elsevier Ltd. All rights reserved.
Keywords: Yule distribution Mutations Scaling exponents
1. Introduction The species–area relationship (SAR) S ¼ cAz ;
ð1Þ
dates back to at least 1921 when Arrhenius (1921) proposed Eq. (1) as an empirical formula for the number of species S of a given taxonomic group (or genus) found in a region of area A. Alternative SARS and their derivations, which are often based on assumed underlying species-abundance distributions (SADs), have been considered and debated over the years (May, 1975; May and Stumpf, 2000; McGuinness, 1984; Chisholm, 2007). While there has been no general consensus on the range of validity of particular SARs (for varying A), recent empirical and theoretical studies by Harte et al. (2009), O'Dwyer and Green (2010) and Storch et al. (2012) strongly suggest that Eq. (1) is not valid for ‘small A’, but that in certain situations, may be valid asymptotically n Corresponding author at: Mathematical Institute, Andrew Wiles Building, University of Oxford, Oxford, OX2 6GG. Tel.: þ 44 1865 611511. E-mail address: [email protected] (T.E. Lee).
http://dx.doi.org/10.1016/j.jtbi.2014.08.011 0022-5193/& 2014 Elsevier Ltd. All rights reserved.
for larger values of A. Our purpose here is not to discuss or comment on the large body of work that has been published on this problem (see May, 1975; May and Stumpf, 2000; McGuinness, 1984; Chisholm, 2007; Harte et al., 2009; O'Dwyer and Green, 2010; Storch et al., 2012 and references quoted therein), but rather to present a new derivation of Eq. (1) which is valid asymptotically for large A. Our derivation of Eq. (1) is based on ideas published over ninety years ago by Yule (1925) on the evolution of species through specific and generic mutations within and between genera of species. In the following section we present a stochastic dynamical systems model based on Yule's original ideas. In Section 3 we show that the species–area relation (1) follows from the asymptotic form of Yule's equilibrium distribution for a large number of species in a genus. Moreover, our asymptotic derivation of (1) provides an explicit expression for the exponent z in Eq. (1) in terms of Yule's species and generic mutation rates. In Section 4 we present some case studies as examples, showing close agreement between our formulae for z and the previous empirical studies based on Eq. (1). Our results are summarised and discussed in the final section.
130
C.J. Thompson et al. / Journal of Theoretical Biology 363 (2014) 129–133
2. The Yule distribution In his original work Yule (1925) proposed evolution of species and genera through two kinds of mutations which he described as follows:
Within any species, in any time interval, a ‘specific mutation’
can occur yielding a new species, but within the same genus as the parent. Within any genus, in any time interval, a ‘generic mutation’ can occur yielding a new species so different from the parent that it is placed in a new genus.
We assume (with Yule) that the rate s of specific mutations is a constant, and the same for all species within a genus; and that the rate g of generic mutations is a constant and same for all genera. Notice that at any time, the number of a species within a genera increases due to both s and the number of species within the genera. That is, the rich get richer (also known as ‘preferential attachment’). Yule (1925) was at pains to point out that these assumptions are simplistic and were designed solely to elicit quantitative comparisons with well-known facts. He provides the following example to demonstrate the last of his simplifying assumptions: if A, B, C, and D are the existing genera, and one of them throws a generic mutation, it is assumed that this will represent a new genus E. That is, the possibility that the new species may be classed under an existing genera is ignored. Combine this with the absence of a death process, and it becomes clear that the number of genera with few species could easily be overestimated. This overestimation has recently been addressed by the Birth–DeathMutation (BDM) Process (Maruvka et al., 2011). Despite these assumptions, we will see in Section 4 that Yule's comparisons with Willis' data proved to be successful. Not surprisingly, Yule's analysis (in the early 1920s) was long and complicated. Our main purpose in this section is to present a simple dynamical systems model for P n ðtÞ ¼ the probability that the genus has n species at time t:
ð2Þ
and to show that Yule's (equilibrium) distribution is a steady state solution of this system. Consider first the case n ¼1. In continuous time it follows from Yule's assumptions that dP 1 ðtÞ ¼ g ðs þ gÞP 1 ðtÞ; dt
ð3Þ
where the first (input) term on the right-hand side of (3) represents a generic mutation from another genus, and the second term represents losses from the single species in the (new) genus through specific and generic mutations. Similarly, for an established genus with at least two species, it follows from Yule's assumptions that dP n ðtÞ ¼ ðn 1ÞsP n 1 ðtÞ ðnsþ gÞP n ðtÞ; dt
n ¼ 2; 3; …:
ð4Þ
In Eq. (4) the first term arises from the n 1 possible specific mutations from n 1 species in the genus, and the second term represents losses from n species in the genus due to n specific mutations, and a constant loss rate g due to generic mutations. There are of course many possible dynamical models for Pn(t). We again note that Yule's scheme ignores death, and therefore extinction. As a consequence, and as mentioned by Yule (1925, p.38), ‘on our assumptions, the mean size of a genus after infinite time must itself be infinite’. Extinction could of course be included in our system, by adding appropriate P n þ 1 ðtÞ terms to the righthand side of (4) for example. We will return to these issues in the
final section. The important point to stress here is that Eqs. (3) and (4) although simplistic are biologically reasonable and logically consistent in the sense that on summation of Eqs. (3) and (4) we deduce that 1 d 1 ∑ P n ðtÞ ¼ g 1 ∑ P n ðtÞ ; ð5Þ dt n ¼ 1 n¼1 as required by conservation of probabilities, i.e., the initial condition 1
∑ P n ð0Þ ¼ 1
1
∑ P n ðtÞ ¼ 1;
implies
n¼1
n¼1
for all t 4 0. In the Appendix we present an exact solution of Eqs. (3) and (4) and show that in the limit t-1 Pn(t) approaches a globally stable equilibrium given by steady-state solution P nn obtained by setting the left hand sides of (3) and (4) to zero, i.e., P n1 ¼
g ¼ ð1 þ σ Þ 1 ; s þg
ð6Þ
and ðnσ þ 1ÞP nn ¼ ðn 1Þσ P nn 1 ;
n ¼ 2; 3; …:
ð7Þ
where s g
σ¼ ;
ð8Þ
is the ratio of the specific and generic mutation rates which we assume henceforth (with Yule) to be larger than unity. Iterating Eq. (7) we deduce (as shown in the Appendix) that P nn ¼
Γ ðnÞΓ ð1 þ 1=σ Þ ; σΓ ðn þ1 þ 1=σ Þ
n ¼ 1; 2; …;
ð9Þ
where Γ ðkÞ is the gamma function. Using Stirling's formula for Γ ðkÞ when k is large, we arrive at the asymptotic form (as shown in the Appendix) 1 1 ð1 þ 1=σ Þ as n-1: ð10Þ n P nn Γ 1 þ
σ
σ
The distribution P nn equation (9) was obtained by Yule in a rather lengthy and complicated derivation. This is not surprising since his work predated developments in stochastic processes and dynamical systems theory. His notion of equilibrium was also at odds with modern interpretations in terms of steady states. These shortcomings were in fact pointed out by Simon (1955) some thirty years after Yule's paper. Simon re-derived Yule's results, and Eqs. (9) and (10) in particular, in a more general and contemporary setting with interesting applications to word frequencies, city sizes, income distributions and frequencies of scientific publications. Simon (1955) correctly referred to (9) as the (equilibrium) Yule distribution, although that citation has been largely forgotten in recent times. We note in passing that although the P nn , from Eqs. (6) and (7), sum to unity as required, there is no reason for the asymptotic form values (10) to do likewise (see the Appendix). Nevertheless, the Yule distribution has some unusual and even ‘paradoxical’ properties (as stated by Yule, 1925, p. 38). In particular it follows from (9) and (10) that the (equilibrium) mean number of species in a genus, obtained by multiplying P nn , Eq. (9), by n and summing on n from 1 to 1, diverges when σ 4 1 by virtue of the asymptotic form (10) for the tail distribution P nn . In the following section we show that the asymptotic form of Yule's distribution (10) can be used to derive the species–area relationship (1) (for large A) with an explicit expression for the exponent z in terms of Yule's parameter σ defined in (8).
C.J. Thompson et al. / Journal of Theoretical Biology 363 (2014) 129–133
131
3. Derivation of the species–area relationship We consider a particular bioregion and assume that the region (with area A) can only support a maximum number N of species in a given genus. Clearly, this maximum occurs when there is one individual of each species in the genus in the region. While this situation is unstable in a biological sense (due to births and deaths of individuals), it is theoretically valid in a mathematical sense. Thus to calculate the expected number S of species in the region, it is important to note that the actual number of species n in the region is a random variable whose probability distribution must be truncated at the theoretical maximum possible value of N (as discussed in the Appendix). In our model, it is then reasonable to interpret N as the maximum possible number of individuals in a genus that could occupy the given region. With this interpretation, it is then not unreasonable to assume (May, 1975, Eqs. (1.1) and (1.5), for example) that N in the following equations is proportional to A. In the Appendix we show that the expected steady-state number of species (of a given genus) in the region is given asymptotically for large N by N
S ∑ nP nn ;
ð11Þ
n¼1
where P nn is the Yule distribution (9). In order to evaluate the asymptotic form for S on the right hand side of Eq. (11), we multiply (7) on both sides by n, and sum from n¼ 2 to n¼ N. The sum on the left hand side can be expressed as N
N
N
N
n¼2
n¼2
n¼1
n¼1
∑ n2 σ P nn þ ∑ nP nn ¼ ∑ n2 σ P nn ð1 þ σ ÞP n1 þ ∑ nP nn ;
N
N
n¼2
n¼1
ð13Þ
Equating (12) and (13) we see that the n2 terms cancel, and after some elementary algebra we obtain the identity N
∑ nP nn ≔
n¼1
1
σ 1
NðN þ1Þσ P nN 1 ;
ð14Þ
from (11) and (6). The identity (14) is valid for any N Z 1. For large N we deduce from (10), (11) and (14) that the expected number S of species is σ 1 1 S N 2 P nN Γ 1 þ N 1 1=σ as N-1: ð15Þ σ 1 σ 1 σ Finally, to arrive at the species–area relationship (1) we use the assumption that N is proportional to A, as discussed above, to obtain Eq. (1) valid asymptotically for large A with z ¼ 1 1=σ
(1925) are reproduced in Fig. 1 as a log–log plot, along with more recent data from Galapagos plants (Alan Tye, 2005), Yorkshire plants (Stace et al., 2003) and New Britain birds (World Wildlife Fund, 2006). As shown in Fig. 1, these plots reveal power laws of the Yule form (10). Values of σ derived from linear regression of the log–log plots of species-genera data, as shown in Fig. 1, were used to calculate the ‘inferred z’ values given in Table 1, evaluated from Eq. (16). The ‘observed z’ values given in Table 1 on the other hand used linear regression on log–log plots of species–area field data, assuming Eq. (1), to obtain estimates of z. The 95% confidence intervals (CI) for the inferred z values in Table 1 were obtained from regression lines for samples of (two or more) data points for each genus in Fig. 1.
ð12Þ
while the sum on the right hand side, after changing summation variables, can be expressed as ∑ nðn 1Þσ P nn 1 ¼ ∑ nðn þ 1Þσ P nn NðN þ 1Þσ P nN :
Fig. 1. Log–log plots of number of genera versus number of species per genus for four different taxa. Regression lines were fitted to the data. Slopes (of (1 þ 1/σ) from Eq. (10)) are respectively lizards: 1.69, G. plants: 1.68, Y. plants: 1.82 and birds: 1.87 with corresponding r2 values of 0.998, 0.990, 0.996 and 0.983.
ð16Þ
and c a constant which depends on σ. This should be considered an upper bound on the slope of the SAR because, in general, the maximum number of species per genus N will scale sub-linearly with area A.
4. Examples As suggested by the title of the original paper, Yule, 1925 had copious amounts of data to examine courtesy of Dr J.C. Willis FRS. Yules regression analysis of log–log plots (of genera frequency and numbers of species in a genus) of beetles, snakes and lizards in particular confirmed the asymptotic form of his equilibrium distribution (10). Willis's data for lizards extracted from Yule
5. Summary In this paper we have presented a new derivation of the wellknown species–area scaling relation (1). Our derivation is based on ideas published over ninety years ago by Yule (1925) on the evolution of species through specific and generic mutations of species and genera. Yule's derivation of the equilibrium species distribution (10) was long and complicated, and his interpretation of ‘equilibrium’ was at odds with modern developments in stochastic processes and dynamical systems theory. These shortcomings were partly overcome some thirty years after Yule's work by Simon (1955). Simon rederived Yule's equilibrium, or steadystate distribution, in a more general setting, with applications to word frequencies, city sizes, income distributions and frequencies of scientific publications. As noted in Section 2, the Yule distribution is somewhat paradoxical. Hence, power-law species–area relationships (1) may also be an exceptional case, as discussed in May (1975), Harte et al. (2009), O'Dwyer and Green (2010) and Storch et al. (2012). In this paper we have also re-derived Yule's results from a simple, recursive, dynamical systems model. We show in particular that Yule's equilibrium, or steady-state distribution, is in fact a globally attracting equilibrium state (i.e., irrespective of initial conditions) of our recursive, dynamical systems model (7). The details of this are provided in the Appendix. We consider a bioregion that can support a maximum number of species N in one genus, and derive an exact identity, Eq. (14), relating the expected number of species S in the region, and the Yule distribution P nN , Eq. (9). The assumption that N is proportional to the area A of the bioregion, which we discussed in some detail in Section 3, leads to the species–area relation (1) with an explicit expression for the exponent z (16) in terms of the ratio σ of Yule's specific and generic mutation rates. Note that the slope calculated from z is an upper bound since the maximum number of species per genus N will scale sub-linearly with area A.
132
C.J. Thompson et al. / Journal of Theoretical Biology 363 (2014) 129–133
Table 1 Values and confidence intervals (CI) for the exponent z in the species–area relationship (1) ‘inferred’ from (16) (z ¼1–1/σ) and the linear regression slopes given in Fig. 1, and previously ‘observed’ values obtained from species–area field data. Taxon
Location
Observed z
Inferred z
CI
Plants Plants
Yorkshire Galapagos
0.18 0.32
0.15–0.21 0.29–0.35
Lizards Birds
West Indies New Britain
0.21 (Usher (1973)) 0.31–0.33 (Preston, 1962; Hamilton et al., 1963; Johnson and Raven, 1973) 0.30 (Darlington, 1957) 0.18 (Diamond, 1975)
0.31 0.13
0.29–0.33 0.06–0.20
Finally we considered frequency data for the occurrence of species and genera for four different taxa (plants, birds and lizards). Log–log plots confirmed the general asymptotic form of Yule's equilibrium distribution (10) and provided accurate estimates, through linear regression, of Yule's ratio σ. Inferred values for z ¼ 1 1=σ for the four taxa given in Table 1 are in accord with previously published ‘observed z’ values (Usher, 1973; Preston, 1962; Hamilton et al., 1963; Johnson and Raven, 1973; Darlington, 1957; Diamond, 1975) derived from species–area field data. If one includes extinction in the Yule scheme, i.e., maintaining generic mutations (with rate g a 0), one also obtains the SAR equation (1) asymptotically for large A, but with σ in Eq. (16) replaced by σ 0 ¼ ðs μÞ=g, where μ is a specific extinction rate. The details of this analysis, and alternative ways of treating generic mutations, will be presented in subsequent publications.
the exponential terms exp½ ðlsþ gÞt, l ¼ 1; 2; …n. It follows that the steady state solution P nn Eq. (7) is in fact a globally attracting equilibrium state of the system of Eqs. (3) and (4). In order to obtain an asymptotic estimate for P nn when n is large (Eq. (10)), we iterate the recursion relation Eq. (7) for P nn , noting Eq. (A.3) for P n1 , to obtain P nn ¼
ðn 1Þσ ðn 2Þσ σ 1 1 ⋯ ¼ 1 þ nσ 1 þ ðn 1Þσ 1 þ 2σ 1 þ σ σ
Γ ðnÞΓ ð1 þ 1=σ Þ ; Γ ðn þ 1 þ 1=σ Þ
ðA:7Þ
where Γ ðkÞ is the gamma function defined for any (real part of) k 4 1 by Z 1 Γ ðk þ 1Þ ¼ xk e x dx ¼ kΓ ðkÞ: ðA:8Þ 0
Acknowledgments We are grateful to Robert May for his comment on a previous version of this paper. We also thank Ryan Chisholm and an anonymous referee for valuable suggestions and comments on revised versions of this paper.
Appendix A. Derivation of Eqs. (9)–(11) In order to solve the system of equations (3) and (4) for Pn(t), n ¼ 1; 2; … we begin with (3) for P 1 ðtÞ, i.e., dP 1 ¼ g ðs þ gÞP 1 ðtÞ; dt
ðA:1Þ
with general solution P 1 ðtÞ ¼
g þc1 e ðs þ gÞt sþg
ðA:2Þ
where c1 is a constant (determined by the initial value P 1 ð0Þ). Irrespective of the value of c1 we have Eq. (6), i.e., P n1 ¼
1 ¼ lim P 1 ðtÞ; 1 þ σ t-1
s g
σ¼ :
ðA:3Þ
Similarly the differential equation for P 2 ðtÞ (4) is dP 2 ¼ sP 1 ðtÞ ð2s þ gÞP 2 ðtÞ: dt
ðA:4Þ
Substituting (A.2) into (A.4) it is then straightforward to show that the general solution of (A.4) is n
P 2 ðtÞ ¼ P 2 þ c1 e
ðs þ gÞt
þ c2 e
ð2s þ gÞt
;
where c1 and c2 are constants and σ P n2 ¼ P n ¼ lim P 2 ðtÞ; 1 þ 2σ 1 t-1
For large (real) k we have the well-known Stirling's (asymptotic) formula (Abramowitz and Stegun, 1964), pffiffiffiffiffiffiffiffiffi Γ ðk þ 1Þ kk e k 2π k: ðA:9Þ In particular for large n (compared with 1/σ) pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Γ ðnÞ ðn 1Þn 1 e ðn 1Þ 2π ðn 1Þ pffiffiffiffiffiffiffiffiffi 1 ; nn 1 e n 2π n 1 þ O n and similarly 1 1 n þ 1=σ ðn þ 1=σ Þ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Γ nþ þ1 nþ e 2π ðn þ 1=σ Þ σ σ pffiffiffiffiffiffiffiffiffi 1 nn þ 1=σ e n 2π n 1 þ O ; n
ðA:10Þ
ðA:11Þ
where in deriving the second lines of Eqs. (A.10) and (A.11) we have used the fact that for large n x n 1þ ex ; ðA:12Þ n (with x ¼ 1 in Eq. (A.10) and x¼1/σ in Eq. (A.11)). Substituting (A.10) and (A.11) into Eq. (A.7) gives the asymptotic form Eq. (10) for the Yule distribution equation (A.7), and BYU. Finally, by summing (6) and (7), we deduce that for finite N N
∑ P nn ¼ 1 σ NP nN cN 1
n¼1
as N-1:
ðA:13Þ
ðA:5Þ
The truncated distribution (appropriate for a finite region) ( cN 1 P nn ; 1 r n r N P nn0 ¼ ðA:14Þ 0; n 4 N;
ðA:6Þ
is then properly normalised and the mean-number of species (in the region) is given by
is the steady-state solution of (A.4) in agreement with Eq. (7). Proceeding by induction it is then not difficult to show that the general solution of the differential equation (4) for Pn(t) n Z 3 is the sum of P nn , defined recursively by Eq. (7) and a linear combination of
1
N
n¼1
n¼1
S ∑ nP nn0 ∑ nP nn ; as in Eq. (11).
ðA:15Þ
C.J. Thompson et al. / Journal of Theoretical Biology 363 (2014) 129–133
References Abramowitz, M., Stegun, I.A., 1964. Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables. Data provided by Dr. Alan Tye, Head, Botany Department, Charles Darwin Research Station, Isla Santa Cruz, Galapagos, Ecuador. Arrhenius, O., 1921. Species and Area. J. Ecol. 9, 95–99. Chisholm, R.A., 2007. Sampling species abundance distributions: resolving the veil-line debate. J. Theor. Biol. 247, 600–607. Darlington, P.J., 1957. Zoogeography. John Wiley and Sons, New York. Diamond, J.M., 1975. Assembly of species communities. In: Cody, M.L., Diamond, J.M. (Eds.), Ecology and Evolution of Communities. Belknap/Harvard University Press, Cambridge, MA, pp. 342–444. Hamilton, T.H., Barth, R.H., Rubinoff, I., 1963. Species abundance: natural regulation of insular variation. Science 143, 1575–1577. Harte, J., Smith, A.B., Storch, D., 2009. Biodiversity scales from plots to biomes with a universal species–area curve. Ecol. Lett. 12, 789–797. Johnson, M.P., Raven, P.H., 1973. Species number and endemism: the Galapagos Archipelago revisited. Science 179, 893–895. Maruvka, Y.E., Kessler, D.A., Shnerb, N.M., 2011. The birth–death-mutation process: a new paradigm for fat tailed distributions. PloS One 6 (11), e26480. May, R.M., 1975. Patterns of species abundance and diversity. Ecol. Evol. Communities, 81–120.
133
May, R.M., Stumpf, M.P., 2000. Species–area relations in tropical forests. Science 290, 2084–2086. McGuinness, K.A., 1984. Species–area curves. Biol. Rev. 59, 423–440. O'Dwyer, J.P., Green, J.L., 2010. Field theory for biogeography: a spatially explicit model for predicting patterns of biodiversity. Ecol. Lett. 13, 87–95. Preston, F.W., 1962. The canonical distribution of commonness and rarity. Part I. Ecology 43, pp. 185–215 and 431–432. Simon, H.A., 1955. On a class of skew distribution functions. Biometrika 42, 425–440. Stace, C.A., Ellis, R.G., Kent, D.H., McCosh, D.J. (Eds.), 2003. Vice-County Census Catalogue of the Vascular Plants of Great Britain the Isle of Man and the Channel Islands. Botanical Society of the British Isles, London. Storch, D., Keil, P., Jetz, W., 2012. Universal species–area and endemics–area relationships at continental scales. Nature 488, 78–81. Usher, M.A.B., 1973. Biological Management and Conservation: Ecological Theory Applications and Planning. Chapman and Hall, London. World Wildlife Fund. WildFinder: Online Database of Species Distributions, ver. January 2006. 〈www.worldwildlife.org/WildFinder〉. Yule, G.U., 1925. A mathematical theory of evolution, based on the conclusions of Dr. JC Willis FRS. Philos. Trans. R. Soc. Lond. Ser. B, Contain. Pap. Biol. Charact. 213, 21–87.
