Stochastic spatial model of producerconsumer systems on the lattice
Abstract
The objective of this paper is to give a rigorous analysis of a stochastic spatial model of producerconsumer systems that has been recently introduced by Kang and the author to understand the role of space in ecological communities in which individuals compete for resources. Each point of the square lattice is occupied by an individual which is characterized by one of two possible types, and updates its type in continuous time at rate one. Each individual being thought of as a producer and consumer of resources, the new type at each update is chosen at random from a certain interaction neighborhood according to probabilities proportional to the ability of the neighbors to consume the resource produced by the individual to be updated. In addition to giving a complete qualitative picture of the phase diagram of the spatial model, our results indicate that the nonspatial deterministic meanfield approximation of the stochastic process fails to describe the behavior of the system in the presence of local interactions. In particular, we prove that, in the parameter region where the nonspatial model displays bistability, there is a dominant type that wins regardless of its initial density in the spatial model, and that the inclusion of space also translates into a significant reduction of the parameter region where both types coexist.
Stochastic spatial model of producerconsumer systems on the lattice
\runauthorN. Lanchier
\addressSchool of Mathematical and Statistical Sciences,
Arizona State University,
Tempe, AZ 85287, USA.
[class=AMS] \kwd[Primary ]60K35, 91A22
Interacting particle systems, voter model, Richardson model, threshold contact process.
1 Introduction
To understand the role of space in ecological communities in which individuals compete for resources, Kang and the author [10] recently introduced a stochastic spatial model of producerconsumer systems on the lattice. Based on their numerical simulations of this spatial model and their analytical results for its nonspatial deterministic meanfield approximation, they concluded that the inclusion of space drastically affects the outcome of such biological interactions. Their results for the models with two species are reminiscent of the ones found in [15] for the spatially explicit LotkaVolterra model and in [12] for a non Mendelian diploid model: in the parameter region where the nonspatial model displays bistability, there is a dominant type that wins regardless of its initial density in the spatial model, while the inclusion of space also translates into a significant reduction of the parameter region where both species coexist. In the presence of three species, other disagreements between the spatial and nonspatial models appear. The main purpose of this paper is to give rigorous proofs of some of the results stated in [10] for the stochastic spatial model with two species, which we now describe in detail.
Model description.
Each point of the dimensional square lattice is occupied by exactly one individual which is characterized by one of two possible types, say type 1 and type 2. Each individual is thought of as a producer and consumer of resources, and we described the dynamics based on four parameters that denote the ability of an individual of type to consume the resource produced by an individual of type , that we call resource . To avoid degenerate cases, we also assume that each resource can be consumed by at least one type: . Each individual dies at rate one, i.e., the type at each vertex is updated at rate one, with the new type being chosen at random from a certain neighborhood according to probabilities proportional to the ability of the neighbors to consume the resource produced by the individual to be updated. More precisely, the state space consists of all functions that map the square lattice into the set of types and the dynamics are described by the Markov generator
(1) 
where denotes the number of type neighbors
and where the configuration is obtained from the configuration by assigning type to vertex and leaving the type of all the other vertices unchanged. To avoid cumbersome notations, we will write later to indicate that vertices and are neighbors, i.e.,
The positive integer is referred to as the range of the interactions. The fraction in (1) means that the conditional probability that the new type is chosen to be type given that vertex is of type at the time of an update is equal to the overall ability of the neighbors of type to consume resource divided by the overall ability of all the neighbors to consume resource . Note that the assumption on the parameters, for , allows to define
and to rewrite the Markov generator (1) in the form
(2) 
indicating that our model reduces to a twoparameter system. Note also that if one of the original parameters in equation (1) is equal to zero then the denominator in
(3) 
may be equal to zero, in which case the numerator is equal to zero as well. When this happens, we have the following alternative and define (3) as follows:

. In this case, motivated by the fact that there is no neighbor of type , we define the probability (3) at which vertex becomes type to be zero.

. In this case, motivated by the fact that there are only neighbors of type , we define the probability (3) at which vertex becomes type to be one.
We point out that, in addition to be the most natural from a biological point of view, the two assumptions above also make the transition probability (3) continuous with respect to the parameters, which is a key to some of our proofs.
The meanfield approximation.
As pointed out in [10], one of the most interesting aspects of the spatial model is that it results in predictions that differ significantly from its nonspatial deterministic meanfield approximation, thus indicating that the spatial component plays a key role in these interactions. We recall that the meanfield model is obtained by assuming that the population is wellmixing, and refer to Durrett and Levin [8] for more details. This assumption results in a system of coupled differential equations for the densities of each type, and since the densities sum up to one, the meanfield model of (2) reduces to the onedimensional system
(4) 
where denotes the density of type individuals. Interestingly, the limiting behavior of the meanfield model depends on whether the parameters and are smaller or larger than one half, which has a biological interpretation in terms of altruism and selfishness. More precisely,

type is said to be altruistic whenever since in this case the resources it produces are more beneficial for the other type than for its own type,

type is said to be selfish whenever since in this case the resources it produces are more beneficial for its own type than for the other type.
The meanfield model has two trivial equilibria given by and , respectively. By setting the righthand side of equation (4) equal to zero and assuming, in order to study the existence of nontrivial interior fixed points, that the product , we find
from which it follows that
is the unique fixed point that may differ from zero and one. To analyze the global stability of the fixed points, we observe that the sign of the derivative (4) is given by
(5) 
which leads to the following three possible regimes for the meanfield model.
1 – A selfish type always wins against an altruistic type.
Proof.
Without loss of generality, we may assume that indicating that type 1 is altruistic whereas type 2 is selfish. In this case, the first line of equation (5) implies that the derivative of is always negative when hence there is no fixed point in , and the trivial equilibrium is unstable while the trivial equilibrium is globally stable. ∎
2 – In the presence of two selfish types, the system is bistable.
Proof.
Since and we have . The second line of equation (5) implies that the derivative of has the same sign as so the interior fixed point is unstable and the boundary equilibria locally stable: the system converges to if the initial density of type 1 is strictly larger than , and to if the initial density of type 1 is strictly smaller than . ∎
3 – In the presence of two altruistic types, coexistence occurs.
Proof.
Since and the fixed point again belongs to as in the presence of two selfish types. However, the second line in equation (5) now implies that the derivative of has the same sign as so the interior fixed point is globally stable whereas the two boundary equilibria and are unstable. ∎
The stochastic spatial model.
In order to state our results for the spatial model, we start by defining some of the possible regimes the spin system (2) can exhibit.
Definition 1
– For the stochastic spatial model (2), we say that

type wins if starting from any configuration with infinitely many type ,

the system clusters if starting from any translation invariant configuration,

type survives if starting from any configuration with infinitely many type ,

both types coexist if they both survive.
As pointed out above, the analytical results for the spatial and the nonspatial models differ in many aspects, which reveals the importance of the spatial component in the interactions. These disagreements are more pronounced in low dimensions and emerge when both types are selfish or both types are altruistic, while in the presence of one selfish type and one altruistic type the analytical results for both the spatial and nonspatial models agree.
To begin with, we focus on the onedimensional nearest neighbor process. In this case, one obtains a complete picture of the phase diagram based on a simple analysis of the interfaces of the process which almost evolve according to a system of annihilating random walks, whereas the more challenging study in higher dimensions and for larger ranges of interactions relies on very different techniques. We point out that the onedimensional nearest neighbor process appears as a degenerate case in which coexistence is not possible while increasing the spatial dimension or the range of the interactions may result in coexistence of altruistic types.
Theorem 2
– Assume that . Then,

Except for the case , the process clusters.

Assuming in addition that , type 2 wins.
The exclusion of the case is simply due to the fact that, under this assumption, a vertex may flip to a new type if and only if all its neighbors already have this type, which implies that the process has infinitely many absorbing states. While in the absence of space one type wins regardless of the initial densities only in selfishaltruistic systems, Theorem 2 indicates that the most selfish/least altruistic type always wins in the onedimensional nearest neighbor case. In particular, coexistence no longer occurs even in altruistic systems. Although this case is pathological, it is also symptomatic of the general behavior of the spatial model for which the existence of a dominant type that always wins may occur even when both types are altruistic. The next result shows that, in any spatial dimension and for any dispersal range, selfishaltruistic interactions result in the same outcome in the spatial and nonspatial meanfield models.
Theorem 3
– Assume that . Then, type 2 wins.
Even though its rigorous proof is not straightforward, the result is predictable since, when one type is selfish and the other type is altruistic, regardless of the type of the vertex to be updated and the spatial configuration, each of the selfish neighbors is individually more likely to be selected to determine the new type than each of the altruistic neighbors. In contrast, in the case of selfishselfish interactions covered by Theorem 4 below, the spatial and nonspatial models disagree in all dimensions and for all dispersal ranges. Recall that in this case the meanfield model is bistable, indicating that the longterm behavior is determined by both the parameters and the initial densities. The next result shows on the contrary that, including local interactions, even if both types are selfish, enough asymmetry in the parameters implies that the most selfish type wins.
Theorem 4
– For all there is such that type 2 wins when .
Numerical simulations in two dimensions indicate more generally that, when type 2 is selfish, it wins whenever . In the neutral case, the obvious symmetry of the evolution rules implies that none of the types wins, but simulations indicate that the population again clusters, with boundaries getting sharper as the common value of increases. This suggests that, in the presence of two selfish types, the longterm behavior of the onedimensional nearest neighbor process occurs in two dimensions as well, a property that we believe is true in any spatial dimension and for any range of interactions. Looking finally at altruisticaltruistic interactions, recall that coexistence always occurs in the meanfield approximation whereas the least altruistic type always wins in the onedimensional nearest neighbor process. Our last two results show that the longterm behavior of the process in higher spatial dimensions and/or for larger dispersal ranges is intermediate between these two extreme behaviors.
Theorem 5
– For all , there is such that coexistence occurs when
In words, sufficiently strong mutual altruism translates into coexistence of both types. Numerical simulations of the twodimensional system suggest in addition that two altruistic types coexist in the neutral case therefore the coexistence region of the spatial model stretches up to the parameter point corresponding to voter model interactions. Interestingly, our last result shows that there exists a parameter region for which even altruistic types cannot coexist, thus indicating that the inclusion of local interactions translates into a reduction of the coexistence region in general and not only in the onedimensional nearest neighbor case.
Theorem 6
– There exists such that type 2 wins whenever
We note that, even though the inclusion of local interactions always translate into a reduction of the coexistence region, numerical simulations indicate that this reduction becomes less pronounced while increasing the spatial dimension or the range of the interactions. The results for the meanfield model and spatial model when are summarized in Figure 1. In the phase diagram of the spatial model, the dashed regions are the ones that are covered by our theorems, while the thick lines are the critical curves suggested by simulations. For a detailed analysis of the nonspatial model and additional results based on numerical simulations for the spatial model in the presence of three species, we refer the reader to the companion paper [10].
The rest of this paper is devoted to the proofs. To avoid cumbersome notations and messy calculations, we prove all our results, with the exception of Theorem 5, when the range of the interactions is equal to one. However, we point out that all the techniques used in the proofs easily extend to larger dispersal ranges. The reason for excluding Theorem 5 from this rule is that it is the only result that shows a disagreement between the onedimensional nearest neighbor case and the other cases, thus indicating that the range of the interactions plays a specific role in the presence of strong altruisticaltruistic interactions.
2 The onedimensional spatial model
This section is devoted to the proof of Theorem 2 which gives a complete picture of the longterm behavior of the onedimensional process with nearest neighbor interactions. To begin with, we study the evolution in the nonneutral case, and look at the process starting with only type 2 individuals in an interval of length . The following lemma shows that the survival probability of the individuals of type 1 decreases exponentially with when .
Lemma 7
– Assume that and . Then,
where .
Proof.
First, we observe that the process is attractive, which directly follows from the monotonicity of the probabilities in (3) with respect to the fractions of neighbors of each type. In particular, it suffices to prove the result when starting with exactly individuals of type 2. To understand how the number of type 2 vertices evolve, let and assume that the neighbors of have different types. Then, the rate at which vertex jumps from type 2 to type 1 is
(6) 
while the rate at which it jumps from type 1 to type 2 is
(7) 
In particular, as long as there are at least two individuals of type 2, the number of such individuals performs a random walk that jumps to the right at rate and to the left at rate since, due to onedimensional nearest neighbor interactions and the choice of the initial configuration, at all times the set of vertices of type 2 is an interval. Moreover, on the event that type 2 survives, the leftmost vertex of type 2 converges almost surely to while the rightmost vertex of type 2 converges almost surely to , therefore the probability to be estimated is equal to the probability that the asymmetric random walk described above tends to infinity. The result then reduces to a standard estimate for asymmetric random walks. More precisely, let denote the number of individuals of type 2 after updates of the system, and let
Then, using a firststep analysis, for all we have
where . Since , we deduce that and
This completes the proof. ∎
With Lemma 7 in hands, we can now establish the first part of Theorem 2 in the nonneutral case as well as the second part of the theorem.
Lemma 8
– Assume that and . Then, type 2 wins.
Proof.
Let be a positive integer. Then, starting from any initial configuration with infinitely many vertices of type 2, with probability one, there exist and such that
Using in addition that the process is attractive, and that the evolution rules are translation invariant in space and time, Lemma 7 implies that
Since this holds for all and since , the lemma follows. ∎
To complete the proof of Theorem 2, it remains to show that, when , the process clusters. Due to the particular geometry of the configurations in the onedimensional nearest neighbor case, the proof reduces to the analysis of an auxiliary process that we shall call the interface process: the spin system defined on the translated lattice by setting
In words, the process has a one at sites which are located between two vertices that have different types, and a zero at sites which are located between two vertices that have the same type, so the process keeps track of the interfaces of the spatial model.
Lemma 9
– Assume that and . Then, the process clusters.
Proof.
Thinking of each site as being occupied by a particle if and as empty otherwise, the idea of the proof is to establish almost sure extinction of this system of particles, which is equivalent to clustering of the spatial model. First, the proof of Lemma 7 indicates that a vertex of either type, say type , changes its type

at rate one if none of its two nearest neighbors is of type ,
This induces the following dynamics for the interface process: a particle at site jumps to each of its empty nearest neighbors at rate and two particles distance one apart annihilate each other at rate one. More formally, the Markov generator is given by
where configuration is obtained from by exchanging the contents of vertices and while configuration is obtained by killing the particles at and if they exist. Since particles can only annihilate, the probability that a given site is occupied by a particle decreases in time so it has a limit when time goes to infinity. Since in addition the evolution rules of the process are translation invariant, this limit does not depend on the site under consideration:
In other respects, given any two particles alive at time , at some random time almost surely finite, either one of these particles is killed due to a collision with a third particle or both particles annihilate each other. The latter follows immediately from the fact that onedimensional random walks are recurrent while each time two particles are distance one apart, there is a positive probability that they annihilate at their next jump when . This implies that the limit must equal zero, from which we conclude that
for all . This completes the proof. ∎
3 Altruisticselfish interactions
This section is devoted to the proof of Theorem 3 which states that, in the presence of altruisticselfish interactions, the selfish type always wins. The proof is divided into four steps. The first step is to show that, under the assumptions of the theorem, the set of type 2 dominates its counterpart in a process that we shall call perturbation of the voter model. In particular, in order to establish the theorem, it suffices to prove that type 2 wins for this new process. The reason for introducing a perturbation of the voter model is that, contrary to process (2), it can be studied using duality techniques, and the second step of the proof is to describe its dual process in detail, while the third step is to construct selected dual paths that are key to proving that type 2 survives. Finally, the fourth step combines these selected dual paths to a block construction in order to prove that not only type 2 survives but also type 1 goes extinct. Before giving the details of the proof, we note that the third step will be used again in the proof of Theorem 6 while the fourth step in the proof of both Theorems 4 and 6, but these two steps are detailed in this section only. We also point out that Theorem 3 can be proved without the use of a block construction, but since it is needed in the proof of Theorems 4 and 6, we follow the same approach for all three theorems.
Coupling with a voter model perturbation.
We first observe that, under the assumptions of the theorem, there exists a constant fixed from now on such that
(8) 
in which case the set of type 2 for the spatial model (2) dominates stochastically its counterpart in a certain perturbation of the voter model that we denote later by . The dynamics of this voter model perturbation depend on a single parameter:
(9) 
and can be described as follows. As in the original spatial model, the type at each vertex is updated at rate one, but the new type is now chosen according to the following rules:

with probability , the new type at vertex is chosen uniformly at random from the set of the nearest neighbors,

with the residual probability , the new type at vertex is chosen to be type 1 if all the nearest neighbors are of type 1, and type 2 otherwise.
More formally, the dynamics are described by the Markov generator
(10) 
where configuration is obtained from by assigning type to vertex and leaving the type of all the other vertices unchanged. Then, we have the following result.
Lemma 10
– Assume that and (8) holds. Then,
Proof.
This follows from certain inequalities between the transition rates. When all the neighbors of vertex are of the same type at the time of an update, then vertex becomes of this type with probability one in both processes. Also, given that is of type 2 and has at least one neighbor of each type, the rate at which it becomes 1 for process (2) is bounded from above by
(11) 
which is the rate at which it becomes 1 for process (10). Similarly, given that is of type 1 and has at least one neighbor of each type, the rate at which it becomes 2 for process (2) is
(12) 
which is the rate at which it becomes 2 for process (10). To see that the last inequality is indeed true, we introduce the functions
and observe that we have the equalities
Since in addition is concave on whereas is linear, the last inequality in (12) follows. The inequalities on the transition rates in (11) and (12) give the lemma. ∎
Duality with branching coalescing random walks.
In view of Lemma 10, standard coupling arguments – see Liggett [13], Section II.1, for details about coupling techniques – imply that, to prove the theorem, it suffices to prove that type 2 wins for the process (10). Recall that the reason for introducing this process is that, contrary to the original spatial model, it can be studied using duality techniques. First, we construct the process (10) graphically as follows:

for each oriented edge with

we let be a Poisson point process with intensity ,

we draw an arrow from vertex to vertex at times to indicate that the individual at mimics the individual at .


for each vertex

we let be a Poisson point process with intensity ,

we draw a set of arrows from each to vertex at times to indicate that becomes of type 1 if all its neighbors are of type 1, and of type 2 otherwise.

We refer the reader to the lefthand side of Figure 2 for an illustration of the graphical representation. The type at any spacetime point can be deduced from the graphical representation and the configuration of the system at earlier times by going backwards in time. Declare that a dual path from spacetime point down to exists, which we write , whenever there are sequences of times and vertices
such that the following two conditions hold:

for , there is an arrow from to at time and

for , there is no arrow that points at the segments .
The dynamics imply that vertex is of type 1 at time if and only if
(13) 
Note that the setvalued process
defines a system of branching coalescing random walks where particles jump at rate to one of their nearest neighbors chosen uniformly at random and split at rate into particles which are sent to each of their neighbors. In addition, whenever two particles are located on the same vertex they coalesce. See the lefthand side of Figure 2 for an illustration where the system of random walks is represented in thick lines. This process is reminiscent, though not identical, of the dual process of the biased voter model [4] and we also refer to [1, 11] for similar dual processes.
Selected dual paths and random walk estimates.
To bound the probability that a given spacetime point is of type 1, we now construct a dual path that keeps track of a specific particle in the dual process. For each , we define the hyperplane
as well as the deterministic times where . Recall that the constant is defined in equation (9) above. Also, we let denote the Euclidean distance between a vertex and its orthogonal projection on a set . The dual path starts at and is defined recursively from the graphical representation as follows. For all , define
For all , we set while at time we have the alternative:

If for some then .

If and for some , then

If and then
Note that there is a unique that satisfies the second condition above so the dual path is welldefined. In words, the dual path travels backwards in time down the graphical representation following the arrows from tip to tail. When a set of arrows is encountered, the process jumps in the direction that makes it closer to the hyperplane between time and time while after the last time the process jumps uniformly at random in all directions. This, together with the properties of the graphical representation introduced above, implies that the selected dual path always jumps at rate one. At each jump, with probability , the target site is chosen uniformly at random from the set of the nearest neighbors whereas, with probability , the process moves so as to decrease its distance to until time when a similar mechanism operates in the direction of the second axis, and so on. In order to later compare the process (10) with oriented site percolation, we start by collecting important properties of the selected dual path. These properties are based on random walk estimates and are given in the following three lemmas.
Lemma 11
– There exist and such that, for all ,
for all sufficiently large.
Proof.
First, we observe that for all we have
from which we deduce that, for all ,
(14) 
In particular, recalling the definition of and , we obtain
The lemma then follows from large deviation estimates for the Poisson distribution and the Binomial distribution together with standard coupling arguments. ∎
Lemma 12
– For any constant and all ,
for all sufficiently large.
Proof.
Denote by the discretetime onedimensional simple symmetric random walk and introduce the discrete random variable that counts the number of times the th coordinate of the process jumps by time . Since the th coordinate evolves according to the continuoustime simple symmetric random walk run at rate after time , we have
for all . Taking , large deviation estimates for the Poisson distribution imply that the first term on the righthand side can be bounded by
for suitable constants and , and all sufficiently large. To estimate the second term, we use the reflection principle and Chebyshev’s inequality. For all , we have
where the function is the moment generating function of . Since
when is small enough, taking , we can conclude that
for all sufficiently large. This completes the proof. ∎
Lemma 13
– Let and . Then, for all large,
Proof.
Since is stochastically smaller than the rate one simple symmetric random walk on the set of nonnegative integers with a reflecting boundary at zero and starting from the same initial state, the proof of Lemma 12 directly implies that
(15) 
for all large and all . Moreover, by Lemmas 11 and 12,
(16) 
for all large and all . Finally, combining (15) and (16), we conclude that