Simple dynamics in non-monotone Kolmogorov systems

In this paper we analyse the global dynamical behaviour of some classical models in the plane. Informally speaking we prove that the folkloric criteria based on the relative positions of the nullclines for Lotka–Volterra systems are also valid in a wide class of discrete systems. The method of proof consists of dividing the plane into suitable positively invariant regions and applying the theory of translation arcs in a subtle manner. Our approach allows us to extend several results of the theory of monotone systems to nonmonotone systems. Applications in models with weak Allee effect, population models for pioneer-climax species, and predator–prey systems are given.


Introduction
The discrete equation is a popular modeling framework for analysing the dynamical behaviour of a single species. In (1.1), x n 0 is the population density in the n-th generation and g(x n ) 0 represents the density-dependent growth rate (or fitness function) from generation to generation. A common assumption in population dynamics is that g is decreasing. This means that the growth rate is mainly determined by negative density-dependent mechanisms such as intra-specific competition [5] or cannibalism. However, cooperative predation, resource defense, increased availability of mates and conspecific enhancement of reproduction are other biological mechanisms producing non-monotone growth rates, see [20]. For example, alders, big leaf maples, poplars and some pine trees thrive at low densities and decrease at high ones due to overcrowding effects and other ecosystem constraints [6,8,10,33].
In describing the interactions of k species, a natural extension of equation (1.1) is x i (n + 1) = x i (n)H i (x 1 (n), . . . , x k (n)), i = 1, . . . , k. (1. 2) The growth rate H i is typically of the form with g i : [0, +∞) −→ [0, +∞) a continuous function. In particular, if a ij > 0 and g i is decreasing for all i, j = 1, . . . , n, system (1.2) describes the evolution of k competing species. See [7,9,10,13,23] and the references therein for a detailed discussion on these models. Understading the dynamical behaviour of (1.2) is of critical importance from an applied point of view. There are several approaches for analysing this issue for competitive systems. For example, the convexity arguments given by Kon [22], the split Liapunov function method given by Baigent and Hou [2,18] or the theory of carrying simplex [31]. Recently, Hou [16,17] has provided a criterion of global attraction based on the relative position of the nullclines reminiscent to the classical results for the Lotka-Volterra system a ij x j ⎞ ⎠ , i = 1, . . . , k, (1.4) see [1] and the references therein. In contrast with the monotone case, the literature on non-monotone systems is relatively scarce. It is worth noting that system (1.2) can exhibit chaotic dynamics. In this paper we describe the global dynamical picture of model (1.2) when the functions g i are not necessarily decreasing. In particular, our criteria could be perceived as an extension of Hou's results to non-monotone systems. The method of proof is completely different from those papers mentioned above. First, we divide the phase space in suitable positively invariant regions and then we apply the theory of translation arcs [3,12,27] in a subtle manner. As emphasized in § 3, our conclusions are rather sharp.
The organization of the paper is as follows. In § 2, we introduce some notation and definitions. In § 3, we give the main results of the paper. Applications to population models for pioneer-climax species [7,9,10] and species with weak Allee effects [20] are discussed. To show the versatility of our results with different interactions, we discuss the dynamical behaviour of several classical predator-prey systems. We conclude the paper with a discussion on our findings.

Mathematical framework
The Euclidean disk with centre at p = (p 1 , p 2 ) ∈ R 2 and radius R > 0 is denoted by with (x, y) = x 2 + y 2 . A subset of the plane homeomorphic to a (closed) Euclidean disk is called a topological disk. A map h : V ⊂ R 2 −→ V which is continuous and injective is called an embedding. Notice that if h(V ) = V , then h is a homeomorphism. In this section, we work, without further mention, with embeddings defined on topological disks. We recall that h : D −→ D is an orientation preserving embedding if it has degree one, that is, for all x ∈ D, then h is an orientation preserving embedding. Now we give a practical criterion to deduce when a map is an orientation preserving embedding. We have employed the notation ∂D to refer to the boundary of D in R 2 .
Assume the following conditions: Given h : D −→ D an embedding, the ω-limit set of a point p ∈ D is defined as We say that h has trivial dynamics if for all x ∈ D, In the previous theorem, index(h, q) denotes the usual topological index of h at q. We mention that theorem 2.1 in [26] deals with homeomorphisms. However, the same proof works for embeddings.
Notice that G([0, +∞) 2 ) ⊂ [0, +∞) 2 . A common assumption for the single-species model (1.1), already stated by Ricker [28] and by Moran [25], is that there is a positive equilibrium p (the carrying capacity) so that The previous condition can be adapted in model (3.1) as follows: Biologically, condition (P) means that the density of population of species i increases (resp. decreases) in the next generation when the weighted total biomass of both species is below (resp. above) a threshold.
In this section, we will make the assumption Condition (S) implies that G is one-to-one on R. The biological meaning of this is simple. If we take two different initial data, the densities of population are different each other in any future generation. To avoid technical problems, we also suppose: The system has, at most, four equilibria namely (0, 0), (r 1 , 0), (0, r 2 ) and Obviously, the last equilibrium is located at the intersection of L 1 and L 2 . In our analysis, we exclude the case L 1 = L 2 (α = β = 1 and r 1 = r 2 ). We prove that the relative position of L 1 and L 2 completely determines the dynamical behaviour of (3.1). This is a well-known result when the system is monotone. Our contribution will be to show that it is also true for non-monotone systems. We stress that if (x, y) ∈ (0, +∞) 2 \R, then Proposition 3.1. Assume that (P) and (S) are satisfied. Then, G(R) ⊂ R and is an orientation preserving embedding.
We observe that A 1 + αy r 1 and y + βA 1 r 2 . Then, by condition (P), Now, it is clear that Analogously, we can prove that We also observe that condition (S) ensures that G 1 (x, 0) = xg 1 (x) and G 2 (0, y) = yg 2 (y) are strictly increasing (they are locally injective and have two fixed points).
Collecting the above information, we deduce that The conclusion now follows from proposition 2.1 immediately.
Assume that (P) and (S) hold. Then, one of the following properties is satisfied: (ii) There exists n 0 ∈ N so that (x n , y n ) ∈ R for all n n 0 .
Proof. Assume that (x 0 , y 0 ) ∈ (0, +∞) 2 . If (x n , y n ) ∈ R for all n ∈ N, condition (P) ensures that {x n } and {y n } are strictly decreasing, see (3.3). Then, (i) holds. If there is n 0 ∈ N so that (x n0 , y n0 ) ∈ R, we deduce that (x n , y n ) ∈ R for all n n 0 . Notice that G(R) ⊂ R by the previous proposition. The case x 0 = 0 (resp. y 0 = 0) is treated in an analogous manner. Observe that in such a case x n = 0 (resp. y n = 0) for all n ∈ N and y n (resp. x n ) is decreasing provided y n > A 2 (resp. x n > A 1 ).
Remark 3.3. Lemma 3.2 says that it is enough for analysing the dynamical behaviour of (3.1) in R. We repeatedly use this fact in the subsequent results.
The following result describes the behaviour of (3.1) on the axes.
Proof. We study the equilibrium r 1 . The other case is analogous and we omit the details. By remark 3.3, we can restrict our attention on [0, We remark that the unique fixed points of ϕ are 0 and r 1 , and the fact ϕ( Now we are in a position to give the main result of this section. (ii) (0, r 2 ) is a global attractor in (0, +∞) 2 provided r2 β > r 1 and r1 α r 2 .
Proof. First we realize that the origin is always a local repeller in (0, +∞) 2 . Indeed, by condition (P), there is a neighbourhood U of the origin in R so that g 1 (x + αy) > 1 and g 2 (βx + y) > 1 for all (x, y) ∈ U with x = 0 and y = 0. This implies that G 1 (x, y) > x and G 2 (x, y) > y for all (x, y) ∈ U with x = 0 and y = 0. Since G((0, +∞) 2 ) ⊂ (0, +∞) 2 , it is clear that the origin is a local repeller in (0, +∞) 2 . We also stress that by remark 3.3, it is enough to study the dynamical behaviour in R. Now we are ready to prove the theorem.
Next we prove that q = (0, 0) and q = (0, r 2 ). To see this, we check that both fixed points are local repellers in (0, +∞) 2 . At the beginning, we have already mentioned this fact for the origin. On the other hand, the eigenvalues of the linearized system at (0, r 2 ) are where the associated eigenvectors are (0, 1) and (w 1 , w 2 ) with w 1 = 0 respectively. Since r 1 > αr 2 and (P), we conclude that g 1 (αr 2 ) > 1. On the other hand, we have already seen in the proof of proposition 3.1 that φ(y) = yg 2 (y) is an increasing function. Moreover, r 2 > 0 is a global attractor of φ by lemma 3.4. Using these two facts together with (H), we conclude that φ (r 2 ) = 1 + r 2 g 2 (r 2 ) ∈ [0, 1). Observe that by (S), φ (r 2 ) = 0. Now, it is clear that (0, r 2 ) is a hyperbolic saddle point. In particular, it is a local repeller in (0, +∞) 2 . The proof of (ii) is analogous and we omit the details.

1, we have that G(A) ⊂ A and G| A : A −→ A
is an orientation preserving embedding. A critical fact is that Thus, if we apply theorem 2.2 to G| A , we conclude that for each p ∈ A, there exists q ∈ F ix(G| A ) so that On the other hand, we notice that r2 β > r 1 , r1 α > r 2 and (P) imply that g 2 (βr 1 ) > 1 and g 1 (αr 2 ) > 1.
Repeating the argument of the proof of (i), we can prove (r 1 , 0) and (0, r 2 ) are local saddle points. Recall that the origin is always a local repeller in (0, +∞) 2 .
Consequently, for all p ∈ A ∩ (0, +∞) 2 , Finally, we consider a sequence {(x n , y n )} obtained from (3.1) so that (x 0 , y 0 ) ∈ R\A with x 0 > 0 and y 0 > 0. By the same argument as that in lemma 3.2, one of the following facts holds: h1 (x n , y n ) tends to a fixed point of G.
h2 (x n , y n ) ∈ A for all n n 0 with n 0 large enough. If h1 holds, then the fixed point has to be (x * , y * ) because the other fixed points are local repellers in (0, +∞) 2 . If h2 holds, then we apply the above argument. The proof of (iii) is now completed.
Example 3.6. Models with weak Allee effect.
The predation by a generalist predator with a saturating functional response is a common mechanism associated with the presence of weak Allee effects [20]. A natural extension in the plane of the single species models with this Allee effect is x n+1 = x n g(x n + αy n ) y n+1 = y n g(y n + βx n ), (3.7) where α, β > 0, with r > a > 0 and b > 0, see lemma 1.1 in [20].
Although g is not always monotone (see Fig. 1), we can apply theorem 3.5. Some particular choices of parameters are r = 0.5, a = 0.45 and b = 3 together with α = 0.5 and β = 0.6 (i); α = 1.1 and β = 0.6 (iii); α = 1.1 and β = 1.4 (iv). In general, condition (S) is difficult to verify in (3.7) because det DG(x, y) has a very complex expression. Notice that detDG(x, y) can be negative for some points (x, y) ∈ R and for some choice of the parameters. Actually, system (3.7) can exhibit chaotic dynamics.
As a direct application of theorem 3.5, we obtain that (0.1, 0) is a global attractor in case 1, there is an interior fixed point (x * , y * ) that is a global attractor in case 2 and there is trivial dynamics with (0.1, 0) and (0, 0.1) as local attractors in case 3. As mentioned in the title of the example, (3.8) is a variant of the classical model discussed in [24]. We mention that this type of growth rates also appear in the evolution of climax species, see [7,9,10]. As shown in [7], there is an asymptotically stable 2-cycle for this model so that the two species can coexist. This indicates that the conclusions in theorem 3.5 is not true for the model (3.9) by noticing that r2 β > r 1 and r1 α r 2 .

Extinction in planar systems beyond (3.1)
Predator-prey models are prototypes of Kolmogorov systems in which the results of § 3 are not directly applicable. However, some ideas developed in theorem 3.5 also work in these models. This shows the versatility of the mathematical framework given in § 2.
Condition (PP1) means that the intra-specific competition is contest (see [15]). (PP2) indicates that the growth rate of the predator is the result of the conjunction of two biological facts: the intra-specific competition and the consumption of the prey. (PP3) has an analogous meaning for the prey. (LG) says that each species in isolation has logistic growth rate with carrying capacities p * and q * respectively. (PP4) says that the predator density decreases in the next generation when it is above a suitable threshold. We say that system (4.1) is permanent if there are two constants ε, M > 0 so that ε lim inf x n lim sup x n M ε lim inf y n lim sup y n M for all sequence {(x n , y n )} of (4.1) with initial condition (x 0 , y 0 ) ∈ (0, +∞) 2 . Informally speaking, the notion of permanence excludes the extinction of some species in the system. To apply the classical theory of permanent systems we need to recall the notion of absorbing set. We say that a positively invariant set R ⊂ [0, +∞) 2 is an absorbing set for (4.1) if for all (x 0 , y 0 ) ∈ [0, +∞) 2 , there is n 0 ∈ N so that F n (x 0 , y 0 ) ∈ R for all n n 0 .
The following lemma is an immediate consequence of lemma 2.1 in [19] by considering the average Liapunov function V (x, y) = x μ1 y μ2 . See also theorem 3.1 in [29]. Next we give the main result of this section.

Discussion
It is well known that the relative position of the nullclines determines the dynamical behaviour of the classical Lotka-Volterra model (5.1) In this paper we have proved that the same results remain true in a broad family of discrete systems, namely x n+1 = x n g 1 (x n + αy n ) y n+1 = y n g 2 (y n + βx n ). and G(x, y) = (xg 1 (x + αy), yg 2 (y + βx)).
The role of these conditions is critical. Notice that (P) is a necessary condition to maintain the dynamical behaviour of (5.1) in (5.2). On the other hand, if we drop (S), as discussed in § 3, new phenomena emerge in (5.2) in comparison with (5.1). For example, the presence of 2-cycles or chaotic dynamics. It is worth stressing that condition (P) encompasses functions that are not monotone. Particularly, we can describe the dynamical behaviour of species subject to Allee effects. Other marked examples are the population models for pioneer-climax species that appears when g 1 is decreasing and g 2 is one-humped, see [14,21,33,34]. For a long time, the dominant topic in these models was mainly the exclusion of the pioneer species, ignoring other dynamical patterns. In contrast with this point of view, theorem 3.5 describes all the possible dynamical patterns. In particular, we study the exclusion for the climax species, which has also been analysed recently by Gilbertson and Kot [11].
There are many models in population dynamics that are not of the form (5.2), i.e., the growth rates are not a scalar function composed by a linear combination of the densities of the species. Nevertheless, many results are valid when the map associated with the model satisfies condition (S). This is the case of most predator-prey systems with a generalist predator. For these models, we have proved that the absence of coexistence states leads to the exclusion of the prey. A direct consequence of this is that the presence of any oscillatory behaviour in (0, +∞) 2 implies the existence of a coexistence state.