Continuous models for single species

Malthus model

In Nature there are many groups of individuals of various species occupying a specific area at a specific time. A population is a group of individuals (all members of a single species) who live together in the same habitat and are likely to interbreed. Each population has a unique physical distribution in time and space. It may contain individuals of different ages and its size (density) is likely to change over time, growing or shrinking according to the reproductive success of its members. The study of population dynamics focuses on these changes – how, when, and why they occur. Understanding of population dynamics is useful for interpreting survey data, predicting pest outbreaks, and evaluating the effectiveness of control tactics. We will present a series of population models, each of which is applicable to different environmental circumstances. We start with the simplest models: one population. Next, we can complicate models by adding influence of various factors from environment.

Mathematical modeling is based on properties of functions and their derivatives. The quantity that determines the change rate of the function N(t) is its derivative \(\displaystyle\frac{dN(t)}{dt}\). If \(\displaystyle\frac{dN(t)}{dt}> 0\), the function \(N(t)\) increases if the argument \(t\) increases. On contrary, if \(\displaystyle\frac{dN(t)}{dt}<0\), the function \(N(t)\) decreases if \(t\), increases. It is sufficient information to start the modeling.

The first step of modeling: why the number N(t) of individuals in the population changes in time? One can consider three primary processes: birth of the individual, death of somebody or migration: somobdy leaves the population, somebody comes in. In mathematical language it means that the rate of changes is determined by three elements:

\[\frac{dN(t)}{dt}= +\text{gain (birth)} - \text{loss (death)} \pm \text{migration}\]
  1. The left side describes the rate of changes of population : \(\displaystyle\frac{dN(t)}{dt}\)

  2. the sign + : the growth of population caused by birth of new individuals is an increassing function \(f(N)\) of the number \(N\) of individuals in the population: Much more individuals are born in the population of 1 mln individuals than in the population of 50 thousands of individuals (it is a rule and in many cases it can be accepted but there are cases when it is not a true assumption). The simplest function is the linear function:

\[f(N) = a N\]

where the coefficient \(a \gt 0\) is the birth rate of the population.

  1. the sign - : the loss of population caused by death of individuals is an increasing function \(g(N)\) of \(N\): In the population 1 mln more individuals pass away than in the population of 50 thousands. Again, the simples function is a linear one,

\[g(N) = b N\]

where the coefficient \(b \gt 0\) jest the death rate.

At this stage of modleing, we neglect the migration processes. From the above considerations it follows that the simplest model which describe time changes in population can take the form

\[\frac{dN(t)}{dt}= a N(t) - b N(t) = (a-b) N(t) = k N(t), \quad \quad k=a-b\]

It is a linear ordinary differential equation of the first order. From the theory of differential equations we know that the initial condition has to be imposed:

\[N(0) = N_0\]

It is a simple differential equation which can be solved by the method of variable separations:

\[N(t)= N_0 e^{kt}\]

From this solution one can infer that:

  1. if the growth rate is the same as the death rate, i.e., \(k=0\) or \(a=b\), then \(N(t) = N(0)\) and the population is stable with fixed number of individuals.

  2. if the birth rate is greater than the death rate \(k\gt 0\) or \(a \gt b\), then the population grows: \(N(t) \gt N(0)\).

  3. if the birth rate is less than the death : \(k\lt 0\) czyli \(a \lt b\), then the population is smaller and smaller and can dissappear: \(N(t) \lt N(0)\).

The above model was proposed by Malthus in 1798. The main feature of this model is exponential growth (for a > b) or exponential decay (for a < b) of population. There are examples for which this model is quite good but only in some limited interval of time (e.g. the growth of bacteria colony under favorable conditions)

Applying Sage for analysis of the Malthus model

Analytical solution

In SAGE it is easily to obtain an analytical form of the solution of this differential equation.

Experiment with Sage!

Run the code cell below and obtain analytical solution of the Mathus model.

Experiment with Sage!

Run the code cell below and obtain plot solutions of Mathus model for various parameters. The plot depicts the growth, lack of changes and decay of population dependening on value of the parameter k=a-b.

Experiment with Sage!

Rerun the above example with:

  • different positive the growth rates \(k=a-b>0\)

  • different negative the growth rates \(k=a-b<0\)

How does the evolution of population behave?

The Verhulst (logistic) model

Let us assume that in the Malthus model we take into account only birth process and neglect the death process, i.e. \(b=0\). Then the population grows exponentially:

\[N(t)= N_0 e^{at}\]

and it is visualized in the above figure for the case \(k=a >0\). Such a fast rate of growth can be observed in some specific systems and only in a limited interval of time. Generally a rapid growth can lead to an imbalance in nature. Moreover, food resources are limited in the finite regions where the population lives. It means that the birth rate \(a\) starts to slow down and self-regulation is observed. The first model which take into account this effect was proposed by Verhulst in 1838. In this model, the coefficient \(a\) is not fixed but it depends on a number of individuals \(N\) in the population:

\[a \to a(N)\]

and the evolution equation assumes the form:

\[\frac{dN}{dt} = a(N) \; N, \quad N(0)=N_0\]

The dependence of the function \(a(N)\) upon \(N\) should have the following feature: If the population increases then the growth rate should decrease because the more difficult is access to food and the amount of food available is smaller if \(N\) grows. Of course there are many functions which fulfil this property: all decreasing function of \(N\). Verhulst proposed the following function:

\[a(N)= r \left[1- \frac{N}{K}\right]\]

where \(r \gt 0\) is a parameter whic has similar interpretation as the parameter \(a\) in the Malthus model (it characterizes the growth rate). The parameter \(K \gt 0\) characterizes resources of food and is named carrying capacity of the environment. Let us note that the constant \(K\) is in the quotient and \(N/K\) is dimensionless. It means that \(K=N_s\) is a characteristic number of population. If \(N \gt K\) then \(a(N) \lt 0\) and population decreases. In turn, if \(N \lt K\) then \(a(N) \gt 0\) and population grows.

In what a way K influences of the growth rate? We note that when K increases (more food) then N/K decrease and in consequence 1-N/K increase. It means that a(N) increases or the growth rate increases and the population grows. Hence the conclusion is:

Note

When K increases the growth rate increases

From the above consideration we obtain the evolution equation in the form

\[\frac{dN}{dt} = r \left[1- \frac{N}{K}\right] N, \quad \quad N(0)=N_0\]

Thie equation is called the Verhulst equation or the logistic equation. There are two parameters: r and K. There is the question: how many relevant parameters are in the Verhuls model. To answer this question we should scale this equation to the dimensionless form.To this aim, let us define the relative number of individuals in the population:

\[x= \frac{N}{K}\]

The dimensionless time is defined as

\[s=r t\]

In the new variables the Verhulst equation takes the form:

\[\frac{dx}{ds} = f(x) = x[1-x], \quad x = x(s), \quad \quad x(0) = x_0 = \frac{1}{K} N(0)\]

We see that in this equation, there are no parameteres. We can conclude that:

Important

Relevant properties of the system do not depend on values of r oraz K. The dependence is only quantitative, not qualitative.

Worked example - stationary states of Verhulst model

We want to find stationary states of the system and determine their stability.

  1. we have to solve the algebraic equation \(f(x) =0\), i.e. \(x(1-x)=0\),

Hence, we find two stationary states \(x_1=0\) and \(x_2=1\)

  1. to determine their stability we have to calculate the derivative: \(\lambda = f'(x) = 1- 2 x\)

and find the exponent \(\lambda_1= f'(x_1) = 1 \gt 0\) (i.e. it is unstable), \(\lambda_2= f'(x_2) =-1 \lt 0\) (asymptotically stable)

Because there is only one stable stationaryy state \(x=1\), all solutions of the Verhulst equation \(x(s)\) with the initial condition \(x_0 \gt 0\) tends to this state \(x=1\).

If the initial condition is \(x_0=0\) the solution is \(x(s)=0\). Howver, an arbitrary small perturbation throws out the system from this state and the next the system starts to evlove to the stable state \(x=1\).

The verhulst equation can be solved analytically by the methos of the variable separation in the following way:

\[\frac{dx}{x(1-x)} = ds \quad \quad \mbox{or equivalently} \quad \quad \left[ \frac{1}{x}-\frac{1}{x-1}\right] dx = ds\]

Alternatively we can use Sage to obtain indefinite integral:

Next we integrate both sides:

\[\int_{x_0}^{ x(s)} \left[\frac{1}{x} -\frac{1}{x-1}\right] dx = \int_0^{ s} ds\]

The final result is

\[x(s) = \frac{x_0 e^s}{1+ x_0(e^s -1)} = \frac{x_0}{x_0 + e^{-s}(1-x_0)}\]

If \(x_0=0\) then we see that \(x(s)=0\). If \(x_0 \gt 0\) then it is seen that \(x(s)\) tends to the stationary state \(x_2=1\).

Now, we want to return to the “old” variables \(N=N(t)\) and time \(t\). We insert \(x=N/K\) and \(s=rt\) to the solution for \(x(t)\) and obtain

\[N(t) = \frac{K N_0}{N_0 +(K-N_0)e^{-rt}}\]

The stable state is \(x = N/K=1\), i.e., \(N=K\). So, \(K\) is a characteristic number of individuals which is found in the stationary state. In the Verhuls model, the evolution is not exponential but the saturation effect is observed: at the fixed carrying capacity of the environment the population is stable at he value \(K\), i.e. it is equal to the paramtere \(K\) in the model. Therefore in many case it seems to be more realistic than the Malthus model.

In the Figure we depict 3 characteristic solution in dependence on initial conditions:

\((A) \quad N_0 \in(0, K/2), \quad \quad N(t) \quad\mbox{is an increasing function of time and tends to the saturated values wartości K}\)

\((B) \quad N_0 \in [K/2, K), \quad \quad N(t) \quad\mbox{is also an increasing to K}\)

\((C) \quad N_0 \gt K, \quad \quad \quad \quad N(t) \quad \mbox{is a decreasing function of time and tends to the value K}\)

In the case (A), the solution has a shape of the capital letter S and therefore is sometimes called the sigmoid function, which is popular in the problems of artificial intelligence and neuron networks.

Experiment with Sage!

A solution of the Velhulst equation can be obtained using the system of computer algebra. Sage can have some problems to get an explicit result (i.e., \(x(t)=...\)). However, if we apply the simplification - full_simplify() then the result is as below.

Time evolution of the population in the Verhulst model

We want to plot for the time evolution of the population for selected initial conditions. If can be easily done with Sage.

Experiment with Sage!

Run code below to obtain time evolution of the Verhulst model with different initial conditions. How the future of the population depend on the current population?

Generalizations and modifications of the Verhulst model

The model with predation (the Ludwig-Jones-Holling model)

In the Verhulst model, the death process can be implicitly contained in the parameter \(r \gt 0\): If \(r\) is smaller, the grow process is smaller because the death process increases. It corresponds to the regime when in the Malthus model \(a \gt b\). In other words the growth rate is greater than the death rate. The Verhulst model can be generalized to the case of presence of predators which kill individuals in the Verhulst model. Because in the presence of predators, a number of preys decreases and according to the rules of modelling, we should add the negative term to the Verhulst equation. So, we postulate the following equation:

\[\frac{dN}{dt} = r \left[1- \frac{N}{K}\right] N -F(N), \quad \quad N(0)=N_0\]

The function \(F(N)\) describes the influence of predators and should exhibit the following properties:

  1. \(F(N=0) =0\) - it means that when there are no preys then there is no influence of predators.

  2. For large values of N, the function \(F(N)\) should saturate, i.e. for \(N\to \infty\) the finction \(F(N)\) tends to a constant finite value. It means that the predator can consume only a limited amount of preys.

If population is too small the predator preferes to change the region of hunting. In this modeling one can apply the Hill function exploited in biochemistry, namely,

\[F(N)= F_n(N)= \frac {BN^n}{A^n + N^n}, \quad n\gt 0\]

Below we present the shape of this function for 3 values of the exponents \(n=2, 4, 6.\) The case \(n=2\) was applied to description of ecological systems, in particular the spruce budworm-forest interaction in the eastern North America (D. Ludwig, D. D. Jones and C. S. Holling, “Qualitative Analysis of Insect Outbreak Systems: The Spruce Budworm and Forest”, The Journal of Animal Ecology, Vol. 47, p. 315, 1978). This model is determined by the equation

\[\frac{dN}{dt} = r \left[1- \frac{N}{K}\right] \, N - \frac {BN^2}{A^2 + N^2}, \quad n \gt 0\]

This model and its modifications have been applied for description of various phenomena and processes. It is the first model which exhibits non-trivial properties. With this respect, the Malthus and Verhulst models can be treated as trivial for analysis.

Experiment with Sage!

Investigate how the The Hill function depens on value the the expenent. For this purpose plot the function for exponent \(n=2,3,4\). Try yourself also different values.

Example of a model which describes Allee effect

Simplifying, we can say that the Allee effect (1931) occur in small or sparse populations (where group size is important for both reproduction and survival) that can result in critical population thresholds below which populations crash to extinction. There are a variety of mechanisms that can create Allee effects, including mate finding, predation, environmental modification, social interactions, reproductive mechanisms, including fertilization efficiency in broadcast spawners, pollen limitation, sperm limitation, reproductive facilitation by conspecifics, and female choice. According to Allee, it should be a threshold of population below which the population will die out. An example of modeling of such an effect is the modified Verhulst equation

\[\frac{dN}{dt} = G(N) = r (N- N_c) \left[1- \frac{N}{K}\right] \;N\]

where \(0 \lt N_c \lt K\) is the next stationary state from the Allee prediction. It is an example of modification. One can construct many other models but this is the simplest. We observe that in this model there are 3 stationary states:

\[\begin{split}\begin{aligned} G(N)= r (N- N_c) \left[1- \frac{N}{K}\right] \;N =\frac{r}{K} (N- N_c) (K- N) \;N = 0, \\ \mbox{ hence} \quad N_1=0, \quad N_2=N_c, \quad N_3=K \end{aligned}\end{split}\]

Stability analysis:

\[G'(N)= \frac{r}{K} \left[(K-N) N -(N-N_c) N + (N-N_c) (K-N)\right]\]
\[\begin{split}\begin{aligned} \lambda_1= G'(N_1) =G'(0) = -K N_c \lt 0 \quad \mbox{(stable state)}, \\ \lambda_2= G'(N_2) =G'(N_c) = (K-N_c)N_c \gt 0 \quad \mbox{(unstable state)}, \\ \lambda_3= G'(N_3) =G'(K) = -K(K- N_c) \lt 0 \quad \mbox{(stable state)} \end{aligned}\end{split}\]

There are 2 stable stationary states: \(N=0\) oraz \(N=K\). The stationary state \(N=N_c\) is unstable. The stationary state \(N=K\) is the same as in the Verhulst model and the statble state \(N=0\) describes the Alee effect: the population can dissapear when initially it is not enough large.

Time evolution of population in the Allee model in dependence of initial conditions (for \(K=1\) and \(N_c=0.5\))

The green and red curves present change on population for sufficiently large initial population (like in the Verhulst model). The blue curve presents the Allee effect: if initially the population is not large, it diminishes in time and finally ceases to exist.