Population dynamics and mathematical models of animal population growth are an important area of ecology, for both understanding how animal populations change through time, as well as managing populations to achieve natural resource objectives. Populations can be classified more or less into one of three categories, depending on their attributes and human values:
Notably, some species may at various times or places fall into more than one of these categories. For example, white-tailed deer in parts of North America occurred in low numbers in much of the early 20th century and were nearly extirpated in some locales (3). Populations of this species have, because of a combination of land use changes and management actions (e.g., harvest regulation, reintroduction) have recovered to allow sustainable harvest (1) throughout most of the range, and are in many places "overabundant", with consequent negative impacts on agriculture, forestry, and human safety (2).
Population models, suitably informed by data, can be useful tools for achieving whatever management goals are appropriate. Furthermore, the basic principles of modeling are the same whether we are talking about sustainable harvest, control, or conservation, although the specific details and emphases of the models will differ.
There are many ways to construct population models in R (or other platforms including spreadsheets) and there is no one right way. Fundamentally the approaches boil down to
We will start out with simple growth models that involve neither stochastic effects nor density feedback, and then create more complex and realistic models with random effects,density, or both. We will also briefly consider age structure and models based on direct modeling of stochastic processes.
Exponential growth model
The exponential growth model assumes that the population starts out at time t=0 at some initial N(0), and growth at constant exponential rate r. If r>0 this results in theoretically unlimited growth over time, if r<0 eventual asymptotic decline to 0, and if r=0 the population stays where it is, N0. The fundamental equation of growth is
N(t+1)= N(t) exp(r)
This can be solve to find a value for N at any future t, given N0 and r as
N(t)= N(0) exp(rt)
There are thus 2 basic ways to simulate growth: directly (by the 2nd equation) or iteratively (by the first). The 2 give equivalent results (this can also be proven mathematically). Both approaches are illustrated in the attached code, which you can run and compare. You can also change the values of N0, r, and T (number of years simulated) in any sensible (no negative N or time please!) fashion you wish. See this script file.
Adding stochastic effects
As a (small) concession to realism, suppose that we assume that population growth rates are not constant every year, but instead allow rates to vary stochastically, as they would under environmental influences. A simple way to do this is to specify a mean (average) growth rate e.g. r.mean=0.1, a standard deviation e.g. r.sd= 0.2, and then draw annual values of r from a Normal distribution with this mean and sd. The attached code does this, using the iterative form of the growth model. This is a sensible way to proceed since the stochastic influence needs to propagate through the values of N over time, and I don't really see how to do that in a 1-step solution.
I have added a couple of other modifications to this model that might make it more useful in some situations, particularly involving small populations. First, I've rounded N to integer values, and second I've assigned an extinction threshold of N=2. In combination these features avoid nonsensical situations such as N is a small fractional value, or N=1, both equivalent to extinction in sexually reproducing organisms. Otherwise it could be possible for a population to 'recover' numerically from what should be extinction. The script file contains these additional changes to the exponential model.
Adding density feedback-- The logistic growth model
The above population models assume that there is no upper limit to population growth, and no density feedback (slower growth rate as populations increase). While this might be realistic for some species, for many species there is ample evidence that
For additional realism, I have code 2 versions of the logistic model of density dependent growth. This model involves the parameters
The first version is the exact, continuous model based on integration of differential equations, and the second is the discrete-time approximation for the same model (see Williams et al. 2002 Chapter 8 for more). The second will give results that are numerically very close but not the same as the first, as you can see by running each under the same conditions. The second approach tends to be used more in simulation, because it is easier to add features such as stochastic growth (as in the second example). The script file contains both the continous and discrete- time models, as well as a modification allowing for stochastic variation in rmax over time.
Modeling stochastic processes directly
The above models are fine for many applications, particularly where were dealing with large populations. They have some undesirable features however, such as allowance of inadmissible values such as noninteger values of N. In addition, they somewhat "black box" the actual processes that are occurring in nature, which are fundamentally stochastic birth and death processes. An alternative approach which is very flexible is to model these processes directly, for instance F(t) as a Poisson process based on current N(t) and some average recruitment rate, and D(t) as a binomial process with some average rate. Looked at this way we start with some current N(t) and get
F(t) ~ Poisson( f*N(t))
D(t) ~Binomial (N(t), d)
where f and d are respectively per-capita average recruitment and mortality probability. Then we get next year's N as
N(t+1) = N(t) +F(t) -D(t)
This approach works well at simulating small populations and totally avoids the non-integer and zeros problems. See the attached code for some simple examples.
Incorporating age effects
It is well known that populations with age structure may grow at variable rates, even without stochastic effects or density feedbacks, due to instability in the age structure. Entire volumes have been written on models involving age (and size ) structure, e.g., by Hal Caswell.; see the Wikipedia summary and references therein. We briefly cover age structured population growth in Chapter 3 of Conroy and Carroll (2009); a more in-depth coverage is provided by Williams et al. (2002).
Briefly, the idea is that the population is now stratified into a vector of abundance, with the number in each element representing number alive in an age class at time t. If for example we had 3 possible age classes (birth year, one year , and two year) we might have
N (t)' = [ 100, 50, 25]
If the population is experiencing constant (but age-specific) survival and recruitment, we describe growth by multiplying the current population vector by a projection matrix A containing age- specific recruitment and survival rates. Examples of this are provided in the accomanying R code. Then growth rate is very succinctly described by
N (t+1)= A N (t)
We then try to find a point a which both the proportions of N in each age class as well as overall growth become constant, known as the stable age distribution. As with simpler population growth models there is a iterative solution to finding this point: we can start with some initial N (0) and iteratively find the point where stability occurs. Alternatively, we can use eigenanalysis (see also) to find a direct, analytical solution. Both approaches are illustrated in the attached code.
Next: Applications