State Variables at and Near Steady State.
Using Eqs. 21-23, with the transition functions specified in Eqs. 35-37 and R given by Eq. 3, we obtain (see SI-C for derivations):
\(\frac{\text{dN}}{\text{dt}}=\ \left[b_{0}-d_{0}\frac{(E}{E_{c})}\right][\frac{1.21N^{\frac{4}{3}}\operatorname{}{\left(\frac{1}{\beta}\right)\left(1+\frac{4N\ln(\frac{1}{\beta)}}{3E}\right)}}{E^{\frac{1}{3}}}-\frac{3N^{2}ln(\frac{1}{\beta)}}{2E}]+m_{0}\)(38)
and
\(\frac{\text{dE}}{\text{dt}}={[w}_{0}-d_{0}\frac{(E}{E_{c})}][\frac{2.42E^{\frac{2}{3}}N^{\frac{1}{3}}}{\ln^{\frac{2}{3}}\left(\frac{1}{\beta}\right)}-\frac{2.26E^{\frac{2}{3}}S^{\frac{1}{3}}}{\ln\left(\frac{1}{\beta}\right)}]\ -\frac{w_{10}E}{\operatorname{}\left(\frac{1}{\beta}\right)}+\ m_{0}.\)(39)
For the general case, with immigration and both speciation mechanisms operating, along with extinction, we have
\(\frac{\text{dS}}{\text{dt}}=m_{0}e^{-\mu S-\gamma}+\sigma_{1}\frac{\text{KS}}{K+S}+\sigma_{2}b_{0}\left(\frac{1.21N^{\frac{4}{3}}\operatorname{}\left(\frac{1}{\beta}\right)}{E^{\frac{1}{3}}}\left(1+\frac{4N\ln\left(\frac{1}{\beta}\right)}{3E}\right)-\frac{1.5N^{2}\ln\left(\frac{1}{\beta}\right)}{E}\right)-\frac{{1.35\ d}_{0}}{E_{c}}\frac{S^{4/3}E^{2/3}}{ln(1/\beta)}.\)(40)
Steady States. Setting the time derivatives in Eqs. 38-40 to zero, we obtain the following relationships among the static values of the state variables and the parameters that describe the dynamics. From Eq. 38:
\(\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ E=E_{c}\frac{b_{0}}{d_{0}}(1+\delta_{E})\)(41)
and from Eq. 39:
\(N=[\frac{0.41w_{10}}{w_{0}-d_{0}\frac{(E}{E_{c})}}]^{3}E\left(1-\delta_{N}\right)\text{.\ \ \ \ }\)(42)
The correction terms, \(\delta_{E}\ \)and \(\delta_{N},\ \)are of order\(\frac{{(m}_{0}}{b_{0})\left(\frac{E^{\frac{1}{3}}}{N^{\frac{4}{3}}}\right)}\)and \({(S/E)}^{\frac{1}{3}}\), respectively, which will generally be << 1.
From Eq. 40, for the immigration-only case (\(\sigma_{1}=\ \sigma_{2}\)= 0), we have
\(\text{\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ S}e^{3\mu\frac{S}{4}}=[\frac{0.41m_{0}ln(\frac{1}{\beta)}}{d_{0}(\frac{E}{E_{c})}}]^{\frac{3}{4}}E^{\frac{1}{4}}\text{\ \ }\)(43)
For the case \(m_{0}=\ \sigma_{2}=0,\) and S>> K ,
\(S=[\frac{{1.35\sigma}_{1}\ln(\frac{1}{\beta)}}{d_{0}(\frac{E}{E_{c})}}]^{3/4}E^{\frac{1}{4}}\)(44)
If S << K , then
\(S=[\frac{{1.35\sigma}_{1}\ln(\frac{1}{\beta)}}{d_{0}(\frac{E}{E_{c})}}]^{3}E\)(45)
Finally, for the case \(m_{0}=\sigma_{1}=0:\)
\(S=[\frac{{0.9\sigma}_{2}b_{0}}{d_{0}(\frac{E}{E_{c})}}]^{3/4}\ln(\frac{1}{\beta)}E\)(46)
The Steady-State Species-Area Relationship (SAR). We can derive nested SARs from Eqs 43-46 because in a nested design,Ec , E , and N scale linearly with area. In the immigration-only model, taking the logarithm of Eq. 43 gives:
\(\ S=\ \frac{\ln\left(E\right)}{3\mu}+\frac{\ln\left(m_{0}\right)}{\mu}-\frac{\operatorname{4ln}\left(S\right)}{3\mu}+\frac{\ln\left(\ln\left(\frac{1}{\beta}\right)\right)}{\mu}+\frac{\ln\left(\frac{0.41E_{c}}{d_{0}E}\right)}{\mu}.\)(47)
In a nested SAR design, Ec , E , andN scale linearly with area, and in general, E>> m 0. Hence Sscales logarithmically with area. There are ln(ln(area)) corrections arise from the third and fourth terms (see Eq. 6) on the right-hand side of Eq. 47 and the fifth term contributes a constant. Ifm 0 scales as a power of area, then the second term also gives a ln(area), which is generally smaller than the first term.
For sufficiently small S in Eq. 43, the exponential term can be ignored, and then the first term on the right hand side of Eq. 47 can be set equal to the third term. This results in S ~\(m_{0}^{\frac{3}{4}}E_{0}^{\frac{1}{4}}\ \sim\ \text{area}\ \)ifm 0 scales linearly with area, and S ~ area1/4 if m 0is scale-independent. This limit applies when S < 1/\(\mu\), which for a typical 50 ha tropical forest plot is approximately S< 50 (see SI-E).
The METE SAR was derived (Harte et al. 2009) from the MaxEnt-predicted species- abundance and the species-level spatial-occupancy distributions. In contrast, the DynaMETE prediction does not depend on a spatial distribution function, but it does depend on the choices made for the mechanisms that govern the transition functions, f ,h , and q and that led to Eq. 41. Nevertheless, METE and the static limit of DynaMETE in the migration-only model predict the same approximate functional form for the SAR: S ~ ln(area) with ln(ln(area)) corrections. Numerical simulations for a variety of choices of state variables show that, beyond this functional similarity, they predict nearly overlapping SARs (see for example Fig. 2).
In the first speciation model, from Eq. 44 with S>> K we get quarter-power scaling of species richness with E (up to a \(\ln\left(1/\beta\right)\)correction) and therefore with area. In that model, with S << K, and for any S in the second speciation model, S scales linearly with area, up to a\(\ln\left(1/\beta\right)\) ~ ln(N )-ln(S ) correction. The first model, with a saturation term K that is small compared to the steady state species richness, thus yields a more realistic species-area relationship.
A Productivity-biomass-diversity-abundance relationship at steady state. Metabolic scaling informs us that individual mass (m ) is related to individual metabolism by\(m(\varepsilon)=\ m(1)\varepsilon^{4/3}\). Summing over the steady-state structure function, total biomass, B , is then (following the methods in SI-C):
B \(=m(1)S\sum_{n,\varepsilon}{\varepsilon^{\frac{4}{3}}nR(n,\varepsilon}\left|S,N,E\right)=\ m(1)\frac{3.57E^{\frac{4}{3}}}{S^{\frac{1}{3}}ln(\frac{1}{\beta)}}\). (48)
This community mass-metabolism relationship thus involves species richness, and also total abundance via the ln(1/\(\beta)\ \)term. Interpreting the state variable E as total net productivityP of the community, Eq. 48 provides a relationship among productivity, biomass, species richness and abundance:
\({{P=0.385\ B}^{\frac{3}{4}}S}^{\frac{1}{4}}\ln^{\frac{3}{4}}(\frac{1}{\beta)}\)(49)
where B and P are measured in units such that m(\(\varepsilon=1\))=1. Eq. 49 does not depend on whether migration, speciation or a combination of both contributes to diversification because Eq. 40 was not used in its derivation. Nor does it depend on the forms of the transition functions and the rate constants, which will differ from habitat to habitat.
Noting the different scaling exponents in the contribution of biomass and species richness to P , and that ln(1/\(\beta)\ \) varies approximately as ln(N ) – ln(S ), the influence of biomass on productivity is considerably stronger than that of species richness, which in turn is stronger than that of abundance, which only enters via the ln(1/\(\beta)\) term. Empirical surveys of the productivity-biomass relationship (Ghedini et al. 2018; Jenkins 2015; Niklas 2007) are qualitatively consistent with these results but extensive analysis will be required to test Eq. 49.
Eq. 49 is an “ideal biodiversity law”, an analog of the ideal gas law that relates thermodynamic state variables. Because this equation was derived using the steady-state structure function in Eq. 3, it will likely no longer hold under non-steady-state conditions. Following disturbance, the full structure function (Eq. 19) is needed to derive the productivity-biomass-abundance-species richness relationship, and it will then depend on the details of the disturbance mechanism.
State variable dynamics near steady state. Time trajectories of the state variables near steady state follow from Eqs. 38-40, which were derived from the static structure function. We examine both the first order responses of the state variables to several kinds of disturbance, expressed by altered transition rate parameters, and the recovery to steady state from a depleted state.
For these dynamical simulations we need to specify the transition rate parameters. We choose a forest that resembles the 50 ha BCI tropical forest plot (Condit et al., 2000; Condit, 1998; Condit et al., 2019; Hubbell et al., 1999). Approximate transition rate constants for this site are given in Table 2 and the rationale for these values is given in SI-E.
Figures 3a–3d illustrate responses of the state variables, over 100 years, to different disturbances represented by changing the values of the rate parameters in the transition functions. Independent of the magnitude of the changes in these parameters, certain general patterns emerge. A decrease in the immigration rate constant, for example under habitat fragmentation that isolates an ecosystem from its meta-community, results in a linear decrease in S at smallt , as well as a weak, damped oscillatory response of N , and nearly undetectable change in E (Fig. 3a). An increase in the death rate (Fig. 3b) generates a slight decrease in S , an initial large decline in N and a weaker decline in E , followed by damped oscillatory behavior. A decrease in the ontogenic growth rate (Fig. 3c) generates a nearly indiscernible increase in S , a large damped oscillatory initial rise in N and weak damped oscillatory initial decrease in E.
Figs. 3d shows the effect on state variable trajectories of combining perturbations in migration, death and growth rates. We do not attempt here a detailed comparison to real data because of the first order approximation used to obtain these theoretical curves, but it is encouraging that the time trajectories of the BCI state variables over the period 1985-2015 (inset in Fig. 3d) also exhibit a steady decline inS , decline and then partial recovery in N , and weak variability in E .
If we fix the transition rate parameters at their undisturbed values (Table 2) in the immigration-only model, and initially reduce the state variables 20% from their steady state values, their return to steady state is shown in Fig. 4. Noteworthy is the monotonic recovery ofS , the large overshoot and then decline to steady state ofN , and a much weaker overshoot and decline of E .
If speciation is the driver of diversification, the pattern of recovery of species richness is markedly different. In the first speciation model, with K >> S , recovery ofS is sigmoidal and extremely slow, with 90% recovery taking approximately 104 years. In the second speciation model, and in the first with K << S ,S recovers to 90% of steady state in approximately 4000 years, and at an ever-slowing, rather than sigmoidal, rate. The recovery trajectories of N and E are nearly the same in both speciation models and very similar to that in the immigration-only case (Fig. 4).