Deterministic approximations for stochastic processes in population biology

https://doi.org/10.1016/S0167-739X(00)00067-4Get rights and content

Abstract

Differential equations are frequently used as deterministic models for stochastic processes. They approximate the expectations of the modeled stochastic processes, and are presumably good approximations when the state space (population) is sufficiently large. Although experimental simulations often support this approach, the justification of the use of deterministic equations for describing stochastic processes is not obvious. Some examples show that the deterministic approximation can deviate considerably from the exact result even for large populations.

In this paper, we study two model examples: (a) A one sex pair formation process occurring in a closed population, with a nonlinear transition rate that results from assuming the mass action law for the pairing rule. This process has a finite state space. (b) A one sex pair formation process obtained from the first scenario by adding the effects of immigration and death. This process has an infinite state space.

For both scenarios we compute analytically the equilibrium behavior of the stochastic and the deterministic models. We use these computations to compare the equilibrium expectations (in the infinite case) and the asymptotic expansions of the equilibrium expectations (in the finite state space model) of the stochastic processes with their deterministic counterparts. We use this comparison to study the appropriate representation of the mass action type reactions, and to study the quality of the deterministic model as an approximation for the stochastic scenario.

Introduction

Many biological processes have a stochastic nature and can be described by stochastic models. Investigation of such stochastic models is often difficult and therefore deterministic differential equations, for which many analytical tools are available, are frequently used instead. When the stochastic transition rates are nonlinear, the deterministic models become merely approximations that hope to approximate the stochastic process if the discussed population is sufficiently large. In general, it is difficult to determine quantitatively to what extent the deterministic models approximate the underlying stochastic scenario. To this end, one first needs to attribute an appropriate interpretation to the deterministic equations. To illustrate, we mention the classical example given by the Lotka–Volterra (LV) predator–prey model for the demographics of a predator population (P) and a prey population (R) that interact through predation. LV equations describe the temporal changes in the P and R populations via the two coupled ODEs Ṙ=aR−cRF; Ḟ=dRF−bF, where the positive parameters a, b, c, and d denote, respectively, the natural growth rate of R, the natural death rate of P, the death rate of R due to unfortunate interactions with P, and the reproduction rate of P due to successful predation on R. The quadratic reaction terms cRF and dRF are known as the “mass action law”. The stochastically varying population sizes (P and R) attain integer values, unlike the solution of the differential equations model. Hence, the LV equations are to be interpreted as a model for the time evolution of the expected number of prey/predator individuals, which is not necessarily an integer. This view is reasonable only if one can prove that the LV equations provide a good approximation for these expectations when both populations are sufficiently large and t→∞.

The relations between stochastic processes and their deterministic approximations have been studied in the literature. We mention here a few citations: some surveys in the context of population dynamics can be found in the books of Bartholomew [1], Pielou [11], and Renshaw [12], and general theoretical results can be found in the work of Kurtz [7].

Stochastic and deterministic descriptions are often considered as “equivalent” modeling strategies (see e.g. [12]), where equivalent means that the deterministic description is a good approximation for the expected value of the corresponding stochastic process when the population is large. Numerical simulations of particular examples (e.g. [12]) often support this view. However, simulations do not guarantee a good fit between the stochastic scenario and its deterministic approximation, and the coagulation–fragmentation process (CFP) is a good example for illustrating this difficulty.

A coagulation–fragmentation process describes the stochastic time evolution of a population of N particles distributed into groups which coagulate and fragment at rates that depend on the sizes of the interacting groups. CFPs and their application to various fields have been studied extensively (see [5] for an overview and references). The expected cluster size distributions were studied by means of the deterministic conservation lawċj=12k=1j−1(Rj−kkcj−kck−Qjkcj)−k=1(Rjkcjck−Qj+kkcj+k),j=1,2,…,where cj=cj(t) denotes the expected number (density) of j-particle clusters and Rjk=Rkj and Qjk=Qjjk(j>k) are the reaction coefficients. Rjk and Qjk represent the coagulation rate of clusters of size j and k into clusters of size j+k and the fragmentation rate of clusters of size j into fragments of sizes k and jk, respectively. These equations are named after Smoluchowski who formulated them in 1916 and 1917 for the pure coagulation case only.

The derivation of Smoluchowski’s equations is heuristic [5], and their value as an approximation for the expectation of the stochastic process is not proved. Gueron [5] showed how equations similar to (1.1) can be derived systematically as a deterministic approximation for cj(t), by neglecting all the effects of the correlations. In the “corrected” form of Smoluchowski’s equations, the exact combinatorial count used for describing the mass action law for the coagulation rate of equal size clusters leads to terms of the form 12cj(cj−1) instead of the homogeneous terms of the form cjcj that appear in (1.1).

CFPs have also been modeled with integro-differential equations of the form∂f(x,t)∂t=−120xφ(x,y)f(x,t)dy−0f(x,t)f(z,t)ψ(x,z)dz+120xf(y,t)f(x−y,t)ψ(y,x−y)dy+xf(y,t)ψ(y,x)dy,where x∈[0,∞), f(x,t) is the time dependent expected group size distribution, φ(x, y) is the rate at which groups of size y are generated from the fragmentation of groups of size x, and ψ(x,y) is the rate of fusion of groups of size x with groups of size y. Eq. (1.2) was proposed by Schumann in 1940 [13] for the pure coagulation process and by Melzak [10] for a CFP (see [5] for details).

, have motivated a tremendous amount of mathematical research. Some research has focused on the existence and uniqueness of positive solutions, on the conservation of the total mass, and on the existence of phase transitions (the appearance of an infinitely large cluster, known also as gelation). As Gueron [5] pointed out, these are mathematical properties or artifacts that result from the deterministic approximation, i.e. they are not corresponding properties of the associated stochastic CFP with a finite state space. This discrepancy may be due entirely to the fact that , were not derived as a limit when N→∞ of a finite state space stochastic CFP, and hence, it is not known how well they may approximate the expected group size distribution when N is large. Durrette et al. [4] computed, analytically, the steady state group size distribution for the class of reversible CFPs and studied the asymptotic behavior of the associated group size distribution for large N. They proved that even in the limit as N→∞, there are cases where the deterministic model given by Eq. (1.2) may deviate considerably from the exact solution.

The failure of the deterministic model to approximate the stochastic scenario in the example of CFPs, motivates our interest in the comparison between the stochastic processes and their deterministic counterparts in other cases, and in the appropriate deterministic model representation for the mass action law. In particular, we ask whether or not one should use the form that is derived from an exact combinatorial count or use homogeneous terms that are related to the stochastic scenario only heuristically. To study these questions we deal here with the pair formation process. This process has various applications in the study of population dynamics and the transmission of sexual diseases (see e.g. [3], [8], [9]). We investigate here two types of a one sex pair formation process with mass action kinetics: when no immigration and death are considered, the stochastic process has a finite state space, whereas including these factors makes the state space infinite. For the deterministic model we use the two variants for representing the mass action law, as mentioned above. Our goal is to compare the equilibrium expectations (in the infinite case) and the asymptotic expansions of the equilibrium expectations (in the finite state space model) of the stochastic processes to their deterministic counterparts.

Section snippets

A pair formation process on a finite state space

Consider a one sex stochastic pair formation process in a population of M individuals. All individuals remain in the system, and the pair formation process only redistributes the sub-populations of the couples and singles. The divorce rate, i.e. the rate at which couples break into two single individuals, is linear: it is proportional to the momentary number of pairs with a constant of proportionality denoted by γ. For the pairing process we assume perfect mixing in the population and apply the

A pair formation process on an infinite state space

The previous example dealt with a stochastic process on a finite state space. We now add immigration and death elements to the same model and end up with a stochastic process on an infinite state space.

Consider a one sex pair formation process in a population where singles join the system at the constant rate α, and the death rate of singles is linear with a constant of proportionality λ (couples do not die). The linear divorce rate is γ times the momentary number of pairs. The marriage rate is

Discussion

The justification of explicit omission of stochasticity in order to obtain deterministic model equations, is an important theoretical modeling question. In general, the deterministic models are merely an approximation. In order to pursue the view that deterministic models capture the main behavior of the dynamics of the modeled population, one must show that when the population is sufficiently large, the deterministic solution approaches the stochastic one at least as t→∞, i.e. at equilibrium.

Acknowledgements

This research was supported by the Institute for Mathematics and its Applications (IMA) with funds provided by the National Science Foundation. Part of the research was completed during my stay as a Long Term Visitor at the IMA, during the year devoted to Mathematics in Biology. I thank Prof. Carlos Castillo-Chavez, who was also a Long Term Visitor at the IMA, for many insightful discussions and for his help in improving the presentation.

Shay Gueron received his PhD in applied mathematics from the Technion–I.I.T. in 1999. He spent two years (1991–1993) at Cornell and Princeton, working in the field of mathematical biology. He joined the Department of Mathematics at the Technion in 1993. In 1999, he moved to the Department of Mathematics of the University of Haifa, as an Associate Professor. His research interests are in mathematical biology, bio-fluid dynamics (ciliary motion), population dynamics, stochastic processes,

References (13)

  • D.J. Bartholomew, Stochastic Models for Social Processes, Wiley, New York,...
  • C.M. Bender, S.A. Orszag, Advanced Mathematical Methods for Scientists and Engineers, McGraw-Hill, New York,...
  • K. Dietz et al.

    Epidemiological models for sexually transmitted diseases

    J. Math. Biol.

    (1988)
  • R. Durrett et al.

    Equilibrium behavior of the reversible coagulation–fragmentation processes

    J. Theoret. Probab.

    (1999)
  • S. Gueron

    The Steady-state distributions of coagulation–fragmentation processes

    J. Math. Biol.

    (1998)
  • S. Gueron, Asymptotic expansions for birth–death Markov chains with applications to population dynamics,...
There are more references available in the full text version of this article.

Cited by (3)

  • Dynamics of a nutrient-phytoplankton model with random phytoplankton mortality

    2020, Journal of Theoretical Biology
    Citation Excerpt :

    However, the deterministic approximation of stochastic processes can deviate considerably from the exact result, even for large populations such as the coagulation-fragmentation process and one sex stochastic pair formation process (Gueron, 2003). Many studies have been devoted to expounding the relations between stochastic processes and their deterministic approximations (see Gueron (2003) and the reference cited therein). Generally, it is difficult to determine quantitatively to what extent the deterministic models approximate the underlying stochastic scenario, and it is important in theoretical modeling to justify the explicit omission of stochasticity to achieve deterministic models.

  • Prospects of a computational origin of life endeavor

    2004, Origins of Life and Evolution of the Biosphere

Shay Gueron received his PhD in applied mathematics from the Technion–I.I.T. in 1999. He spent two years (1991–1993) at Cornell and Princeton, working in the field of mathematical biology. He joined the Department of Mathematics at the Technion in 1993. In 1999, he moved to the Department of Mathematics of the University of Haifa, as an Associate Professor. His research interests are in mathematical biology, bio-fluid dynamics (ciliary motion), population dynamics, stochastic processes, computational mathematics and numerical techniques. He is also interested in mathematical enrichment of high achievers, mathematical competitions, problem solving and problem creation, and has been the Israeli Team Leader for the International Mathematical Olympiad, since 1994.

View full text