This is a preprint of: “Phase portraits of the Leslie-Grower system”, Jaume Llibre, Clàudia Valls, Acta Math. Sin. (Engl. Ser.), vol. 42B(5), 1734–1742, 2022. DOI: [10.1007/s10473-022-0502-4] PHASE PORTRAITS OF THE LESLIE-GOWER SYSTEMS JAUME LLIBRE1 AND CLAUDIA VALLS∗2 Abstract. In this paper we characterize the phase portraits of the Leslie-Gower model for competition among species. We give the com- plete description of their phase portraits in the Poincaré disc (i.e. in the compactification of R2 adding the circle S1 of the infinity) modulo topological equivalence. It is well-known that the equilibrium point of the Leslie-Gower model in the interior of the positive quadrant is a global attractor in this open quadrant, and in this paper we characterize where the orbits attracted by this equilibrium born. 1. Introduction and statement of the main results The dynamical relations between predators and prey is one of the funda- mental objects of study in population dynamics. This is mainly due to the fact that to its big number of applications [3] and because it allows a better understanding of the behavior of food chains or trophic networks (see for instance [14, 24]). The first predator-prey model was proposed by the italian mathematician Vito Volterra [1, 2] in its celebrate publish monograph in 1926 [25] where he described the model as a non-linear system of ordinary differential equations. This model coincided with a bidimensional model for biochemic interactions proposed earlier by the american physicist Alfred J. Lotka and this is the main reason why these type of models are called Lotka-Volterra models [2, 14, 15, 24]. The main feature of this model is that the unique positive equilibrium point is a center and all the orbits are closed concentric orbits around the equilibrium point. This implies that the population sizes of the predators and prey oscillate permanently around this point for any initial condition. This behavior of the solutions of the system was early questioned because this is not what happens in nature. After Volterra’s work, there were many attempts to solve the objec- tions of Volterra’s model. One of this attempts was Leslie–Gower model [9] whose main feature is that the predator growth equation is of logistic 2010 Mathematics Subject Classification. Primary 34A05. Secondary 34C05, 37C10. Key words and phrases. predator-prey models, Leslie–Gower system, Poincaré compacti- fication, global phase portraits. 1 2 J. LLIBRE AND C. VALLS type (see [1]) and so it considers implicitly that there exists competition or self-interference between the predators. The Leslie–Gower model has been studied from many points of view, see for instance the recent papers [5, 6, 7, 8, 10, 11, 12, 19, 20, 21, 22, 23, 26, 27, 28, 29, 30, 31, 32, 33]. As far as we know there are no results concerning the characterization of its global dynamics taking into account the orbits which come or scape at infinity. Thus, in this work we will study the global phase portraits of this model in the Poincaré disc. The Leslie–Gower model is described by the differential equations of the form  x  y  ẋ = r 1 − x − qxy, ẏ = s 1 − y, k nx where the dot means derivative with respect to the time t and all the param- eters r, k, q, s, n are positive and have the following meaning:r is the intrinsic rate of growth prey, k is the environmental bearing capacity for prey, q is the maximum consumption rate (maximum quantity of necessary prey that can be consumed for each predator in each unity of time), s is the intrinsic rate of growth predators and n is the quantity of prey with food for the predators. This system is not defined in (0, 0). Making the reparameterization dτ = nx dt, we obtain the equivalent system  nr  ẋ = nr − x − nqy x2 , ẏ = s(nx − y)y, k where the dot denotes derivative with respect to the new independent vari- able τ . Now doing the change (x, y, τ ) → (X, Y, T ) given by 1 1 (1) X= x, Y = y, T = knrτ, k kn we get (2) ẋ = (1 − x − Ay)x2 , ẏ = B(x − y)y, where the dot denotes derivative in the new time T and A = knq/r, B = s/r being A, B > 0 and we have renamed the variables (X, Y ) again as (x, y). Note that system (2) is a continuous extension of the original Leslie–Gower system. In fact we shall present the distinct global phase portraits of the Leslie– Gower system when its parameters A, B varies in R+ in the Poincaré disc. In this way we can describe the dynamics of their orbits which come or go to the infinity of R2 . For doing this we will use the Poincaré compactification. Roughly speaking the Poincaré compactification of a polynomial differ- ential system consists in extending this system to an analytic system on a closed disc D2 of radius one, whose interior is identified with R2 and its boundary, the circle S1 , plays the role of the infinity. This closed disc is called the Poincaré disc, because the technique for doing such an extension LESLIE-GROBER SYSTEMS 3 (a) 0 < B ≤ 1, S = 23, R = 6. (b) B > 1, S = 25, R = 8. Figure 1. Phase portraits of system (2) on the Poincaré disc. is due to Poincaré. For details on this compactification see [4, chapter 5] or the summary presented in subsection 2.1. Theorem 1. The phase portraits of system (2) in the Poincaré disc are totpologically equivalent to one of the 2 phase portraits presented in Figure 1. The proof of Theorem 1 is given in section 3. In Figure 1 S denotes the number of separatrices and R denotes the number of canonical regions. See for more details in subsection 2.1. In subsection 3.2 we will show that system (2) has three equilibria, pj for j = 1, 2, 3. In the next corollary we provide new information from the ecological Leslie-Gower predator-prey model. It is well-known that the equilibrium point p3 in the interior of the positive quadrant is a global attractor in this open quadrant, but it is was unknown where the orbits attracted by p3 born. Corollary 2. The following two staements hold. (a) If the parameter 0 < B ≤ 1, then all the orbits of the Leslie-Gower system (1) attracted by p3 come from the endpoint of the positive x-half-axis, except one orbit coming from the poinf p2 . (b) If the parameter B > 1, then all the orbits of the Leslie-Gower system (1) attracted by p3 come from either the endpoint of the positive x- half-axis, or from the endpoint of the positive y-half-axis, or from the origin, except one orbit coming from the point p2 . 4 J. LLIBRE AND C. VALLS The proof of Corollary 2 follows immediately from the phase portraits of Figure 1. 2. Preliminary results 2.1. Poincaré compactification. In order to classify the global dynamics of a polynomial differential system the first crucial step is to characterize their finite and infinite equilibrium points in the Poincaré compactification [18]. The second main step for determining the global dynamics in the Poincaré disc of a polynomial differential system is the characterization of their separatrices. For the polynomial differential systems in the Poincaré disc it is known that the separatrices are the infinite orbits, the finite equi- librium points, the separatrices of the hyperbolic sectors of the finite and infinite equilibrium points, and the limit cycles. If Σ denotes the set of all separatrices in the Poincaré disc D2 , Σ is a closed set and the components of D2 \ Σ are called the canonical regions. We denote by S and R the number of separatrices and canonical regions, respectively. We consider the set of all polynomial vector fields in R2 of the form (3) (ẋ, ẏ) = X(x, y) = (P (x, y), Q(x, y)), where P and Q are real polynomials in the variables x and y of degrees d1 and d2 , respectively. Take d = max{d1 , d2 }. Denote by Tp S2 be the tangent space to the 2-dimensional sphere S2 = {s = (s1 , s2 , s3 ) ∈ R3 : s21 + s22 + s23 = 1} at the point p. Assume that X is defined in the plane T(0,0,1) S2 = R2 . Consider the central projection f : T(0,0,1) S2 → S2 . This map defines two copies of X, one in the open northern hemisphere and the other in the open southern hemisphere. Denote by X ′ the vector field Df ◦ X defined on S2 except on its equator S1 = {s ∈ S2 : s3 = 0}. Clearly S1 is identified to the infinity of R2 . If X is a planar polynomial vector field of degree d, then p(X) is the only analytic extension of sd−1 ′ 3 X to S . The vector field p(X) is 2 called the Poincaré compactification of the vector field X, for more details see [4, chapter 5]. On the Poincaré sphere S2 we use the following six local charts, which are given by Ui = {s ∈ S2 : si > 0} and Vi = {s ∈ S2 : si < 0}, for i = 1, 2, 3, with the corresponding diffeomorphisms φi : Ui → R2 , ψi : Vi → R2 , defined by φi (s) = −ψi (s) = (sm /si , sn /si ) = (u, v) for m < n and m, n ̸= i. Thus (u, v) will play different roles in the distinct local charts. The expressions of the vector field p(X) are         d 1 u 1 u d+1 1 u (u̇, v̇) = v Q , − uP , , −v P , in U1 , v v v v v v LESLIE-GROBER SYSTEMS 5         u 1 u 1 u 1 (u̇, v̇) = vd P , − uQ , , −v d+1 Q , in U2 , v v v v v v (u̇, v̇) = (P (u, v), Q(u, v)) in U3 . We note that the expressions of the vector field p(X) in the local chart (Vi , ψi ) is equal to the expression in the local chart (Ui , ϕi ) multiplied by (−1)d−1 for i = 1, 2, 3. The orthogonal projection under π(y1 , y2 , y3 ) = (y1 , y2 ) of the closed northern hemisphere of S2 onto the plane s3 = 0 is a closed disc D2 of radius one centered at the origin of coordinates called the Poincaré disc. Since a copy of the vector field X on the plane R2 is in the open northern hemisphere of S2 , the interior of the Poincaré disc D2 is identified with R2 and the boundary of D2 , the equator S1 of S2 , is identified with the infinity of R2 . Consequently the phase portrait of the vector field X extended to the infinity corresponds to the projection of the phase portrait of the vector field p(X) on the Poincaré disc D2 . The equilibrium points of p(X) in the Poincaré disc lying on S1 are the infinite equilibrium points of the corresponding vector field X. The equilib- rium points of p(X) in the interior of the Poincaré disc, i.e. on S2 \ S1 , are the finite equilibrium points. We note that in the local charts U1 , U2 , V1 and V2 the infinite equilibrium points have their coordinate v = 0. For a polynomial vector field (3) if s ∈ S1 is an infinite equilibrium point, then −s ∈ S1 is another infinite equilibrium point. Thus the number of infinite equilibrium points is even and the local phase portrait of one is that of the other multiplied by (−1)d+1 . 2.2. Separatrix skeleton. Given a flow (D2 , ϕ) by the separatrix skeleton we mean the union of all the separatries of the flow together with one orbit from each one of the canonical regions. Let C1 and C2 be the separatrix skeletons of the flows (D2 , ϕ1 ) and (D2 , ϕ2 ) respectively. We say that C1 and C2 are topologically equivalent if there exists a homeomorphism h : D2 → D2 which sends orbits to orbits preserving or reversing the direction of all orbits. From Markus [13], Neumann [16] and Peixoto [17] it follows the next theorem which shows that is enough to describe the separatrix skeleton in order to determine the topological equivalence class of a differential system in the Poincaré disc D2 . Theorem 3 (Markus–Neumann–Peixoto Theorem). Assume that (D2 , ϕ1 ) and (D2 , ϕ2 ) are two continuous flows with only isolated equilibrium points. Then these flows are topologically equivalent if and only if their separatrix skeletons are equivalent. 6 J. LLIBRE AND C. VALLS 3. Proof of Theorem 1 We will state and prove several auxiliary results. As pointed out in the introduction we will work with system (2) instead of the original Leslie– Gower system. 3.1. Limit cycles. To state and prove the first of our results we recall the Bendixson-Dulac theorem: consider a planar continuously differentiable system (4) ẋ = P (x, y), ẏ = Q(x, y), defined in some open simply connected subset U ⊂ R2 and set X = (P, Q). Assume that there exists a continuously differentiable function F : U → R such that ∂(P F ) ∂(QF ) + ∂x ∂y does not change sign and vanishes only on a set of zero Lebesgue measure. Then system (4) has no periodic orbits contained in U . Proposition 4. System (2) has no limit cycles. Proof. By uniqueness of solutions it is clear that if system (2) has periodic orbits they do not intersect the coordinate axes. By making the change of variables x → ±x, y → ±y, if necessary, we can restrict our study to the first quadrant Ω = {(x, y) ∈ R2 : x > 0, y > 0} and prove that system (2) has no periodic orbits on it. We consider the function D : Ω → R, given by D(x, y) = 1/(x2 y). Note that it is continuously differentiable in Ω. Moreover, ∂(ẋD) ∂(ẏD) B 1 + = − 2 − < 0, ∂x ∂y x y because B > 0. Now we can apply the Bendixson-Dulac theorem and since Ω is simply connected the system has no limit cycles. This concludes the proof of the proposition. □ 3.2. Finite equilibrium points. Now we study the finite equilibrium points of system (2). We introduce the notation ∆ = (1 + B + AB)2 − 4(B + 2AB + A2 B). Proposition 5. System (2) has three singular points p1 = (0, 0), p2 = (1, 0) and p3 = (1, 1)/(A + 1). The singular p1 is formed by two hyperbolic and two parabolic sectors, the singular point p2 is a saddle and the singular point p3 is an stable node if ∆ ≥ 0 and a stable focus if ∆ < 0. LESLIE-GROBER SYSTEMS 7 In the proof of the proposition we provide the explicit position of the hyperbolic and parabolic sectors for p1 in function of the values of the pa- rameter B. Proof. It is clear that the singular points are p1 = (0, 0), p2 = (1, 0) and p3 = (1, 1)/(A + 1). The Jacobian matrix of the system is   −x(3x + 2Ay − 2) −Ax2 J= By B(x − 2y) The eigenvalues of the Jacobian matrix at the singular point p2 are −1, B. Since B > 0, the singular point p2 is a hyperbolic saddle.   0 0 The Jacobian matrix of the singular point p1 is and so p1 is 0 0 linearly zero. To obtain the local behavior of this point we will use the blow-up technique. To do so, let z = y/x and we get the equation ẋ = x2 (1 − x − Axz), ż = xz(B + x − Bz + Axz − 1). Rescaling by the independent variable (the time) by s = xT we obtain the system ẋ = x(1 − x − Axz), ż = z(B + x − Bz + Axz − 1), where the dot means derivative in the new variable s. The Jacobian matrix evaluated at (0, 0) has the eigenvalues 1 and B − 1. Therefore the point (0, 0) is a saddle if and only if B < 1, an unstable node if and only if B > 1 and a saddle-node if and only if B = 1. Moreover, the Jacobian matrix evaluated at (0, (B − 1)/B) has the eigenvalues 1 and −(B −1). Therefore the point (0, (B −1)/B) is a saddle if and only if B > 1, an unstable node if and only if B < 1 and a saddle-node if and only if B = 1. Doing the blow-down process we obtain that if B ∈ (0, 1] then the lo- cal phase portrait at the origin is the following: in the first quadrant (the positive quadrant) there is a hyperbolic sector and in the second quadrant there is an attractor parabolic sector. The boundary of this sector is on the negative x-half axis and on the positive y-half axis and all the other orbits of this parabolic sector arrive to the origin with slope (B − 1)/B. In the third quadrant we have another hyperbolic sector and, finally in the fourth quadrant we have a repeller parabolic sector. The boundary of this para- bolic sector is again the x and y axes and all the other orbits of this sector exit from the origin with slope (B − 1)/B. Furthermore, when B > 1 we get that topologically the local phase portrait of the origin is the same as the one for B ∈ (0, 1], but the separatrix of the sector which was on the positive x-half axis now is inside the first quadrant with slope (B − 1)/B, and the separatrix of the hyperbolic sector which was contained in the third quadrant now is in the interior of this quadrant with slope (B − 1)/B. 8 J. LLIBRE AND C. VALLS The eigenvalues of the Jacobian matrix at the singular point p3 is are √ 1 + B + AB ± ∆ λ± = − 2(A + 1)2 √ Since A, B > 0 we have that ∆ < |1 + B + AB| and so λ± < 0. Hence, the point p3 is locally a stable node if ∆ ≥ 0 and is locally a stable focus if ∆ < 0. This concludes the proof of the proposition. □ 3.3. Infinite equilibrium points. Proposition 6. System (2) has three pairs of infinite equilibrium points: the origin of the local chart U1 (denoted by q1 ) and its diametrally opposite point, q2 = (−1/A, 0) (and its diametrally opposite point) and the origin of the local chart U2 (denoted by q3 ) and its diametrally opposite point. The infinite equilibrium point q1 is an unstable node, q2 is a saddle-node and q3 is formed by two hyperbolic and two parabolic sectors. In the proof of the proposition we provide the explicit position of the hyperbolic and parabolic sectors for q3 . Proof. System (2) on the local chart U1 becomes u̇ = u(1 + Au + (B − 1)v − Buv), v̇ = v(1 + Au − v). The equilibrium points on v = 0 are the solutions of u(1 + Au) = 0 and so we have two equilibrium points q1 = (0, 0) and q2 = (−1/A, 0). Computing the eigenvalues of the Jacobian matrix at q1 we get that they are 1, 1 and so q1 is an unstable node. Computing the eigenvalues of the Jacobian matrix at q2 we get that they are (0, −1) so the point q2 is semi-hyperbolic. Using [4, Theorem 2.19] we conclude that q2 is a saddle-node. Indeed, translating the point q2 to the origin and making a rescaling of time s = −τ we get the new system −A + B + AB −A + 2B + AB u̇ = u + 2 v − Au2 − uv + Bu2 v, (5) A A v̇ = v(v − Au), where the dot means derivative in the new time s. In order to apply [4, Theorem 2.19]  we must  pass the linear part of system (5) to its Jordan −1 0 normal form . This is made doing the linear change of variables 0 0 A(B − 1) + B u=X− Y, v = Y. A2 LESLIE-GROBER SYSTEMS 9 In these new variables system (5) becomes B (A + 2)B(A(B − 1) + B) 2 Ẋ = X − AX 2 − XY + Y A A3 B(2B − A2 + 2AB − 2A) 2 B(A(B − 1) + B)2 3 − XY + Y A2 A4 (A + 1)B 2 Ẏ = −AXY + Y . A Applying [4, Theorem 2.19] we get that q2 is a saddle-node having the nodal sector in the local chart U1 and the two hyperbolic sectors in the local chart V1 . Now we compute system (2) on the local chart U2 and we obtain u̇ = u(−Au + Bv − u2 + (1 − B)uv), v̇ = B(1 − u)v 2 . We only need to study the origin of the local chart U2 which  we denote by  0 0 q3 . Note that the Jacobian matrix evaluated at q3 is and so the 0 0 point is linearly zero. We need to apply the blow-up technique. To do so, let w = v/u and we get the equation u̇ = u2 (−A − u + Bw + (1 − B)uw), ẇ = uw(A + u − uw). Reescaling by the independent variable (the time) by s = uτ we obtain the system u̇ = u(−A − u + Bw + (1 − B)uw), ẇ = Bw(A + u − uw), where the dot means derivative in the new variable s. The unique equilib- rium point over the straight line u = 0 is the origin which is a saddle because their eigenvalues are ±A. Doing the blow-down process, going back through the changes of variables and taking into account how is the flow on the in- variant axes, we obtain that the origin of U2 is formed by two hyperbolic sectors and two parabolic sectors. More precisely, in the first quadrant in the plane (u, v) we have a hyperbolic sector, the separatrices are on the axes. In the second quadrant we have a repeller parabolic sector whose separatrices are on the axes, in the third quadrant we have hyperbolic sector and, finally, in the fourth quadrant we have an attractor parabolic sector. □ 3.4. Global phase portraits. From subsection 2.2 and Proposition 4 in order to prove Theorem 1 we only need to determine the behaviour of the separatrices of the hyperbolic sectors of the finite and infinite singular points, i.e. where they born and where they die. Doing so we will have all the sepa- ratrix skeleton adding one orbit in each canonical region, and consequently we will have the global phase portraits of system (2). But taking into ac- count the two invariant straight lines and the local phase portraits at all the finite and infinite singular points, the place where born and die all the sep- aratrices of the hyperbolic sectors is determined in a unique way for every 10 J. LLIBRE AND C. VALLS one of the two cases in function of the parameter B either if 0 < B ≤ 1 or if B > 1. In this way we obtain the 2 phase portraits in the Poincaré disc of Figure 1. This completes the proof of Theorem 1. Acknowledgements The first author is supported by the Agencia Estatal de Investigación grant PID2019-104658GB-I00 and the H2020 European Research Council grant MSCA-RISE-2017-777911. The second author is partially supported by FCT/Portugal through UID/MAT/04459/2019. References [1] N. Bacaër, A short history of Mathematical Population Dynamics, Springer-Verlag, 2011. [2] A. Bazykin, Nonlinear Dynamics of interacting populations, World Scientific Publish- ing Co. Pte. Ltd., 1998. [3] A.A Berryman, A.P. Gutierrez and R. Arditi, Credible, parsimonious and use- ful predator-prey models. A reply to Abrams, Gleeson and Sarnelle, Ecology. 1995; 76:1980–1985. [4] F. Dumortier, J. Llibre and J.C. Artés, Qualitative Theory of Planar Differential Systems, Springer, 2006. [5] E. Gonzalez-Olivares and A. Rojas-Palma, Global stability in a modified Leslie-Gower type predation model assuming mutual interference among generalist predators, Math. Biosci. Eng. 17 (2020), 7708–7731. [6] Q. Han, L. Chen and D. Jiang, Periodic solution and stationary distribution for sto- chastic predator–prey model with modified Leslie–Gower and Holling type II schemes, Filomat 34 (2020), 1383–1402. [7] Z. Hou, Geometric method for global stability of discrete population models, Discrete Contin. Dyn. Syst. Ser. B 25 (2020), 3305–3334. [8] A.B Junior and N.A. Maidana, A modified Leslie–Gower predator–prey model with alternative food and selective predation of noninfected prey, Math. Methods Appl. Sci. 44 (2021), 3441–3467. [9] P.H. Leslie and J. C. Gower, The properties of a stochastic model for the predator–prey type of interaction between two species, Biometrika, 47 (1960), 219–234. [10] Y. Liu and J. Wei, Spatiotemporal dynamics of a modified Leslie–Gower model with weak Allee effect, Internat. J. Bifur. Chaos Appl. Sci. Engrg. 30 (2020), 2050169, 24 pp. [11] R. Ma, Y. Bai and F. Wang, Dynamical behavior analysis of a two–dimensional discrete predator–prey model with prey refuge and fear factor, J. Appl. Anal. Comput. 10 (2020), 1683–1697. [12] L. Ma and B. Liu, Dynamic Analysis and Optimal Control of a Fractional Order Singular Leslie-Gower Prey–Predator Model, Acta Math. Sci. Ser. B (Engl. Ed.) 40 (2020), 1525–1552. [13] L. Markus, Global structure of ordinary differential equations in the plane, Trans. Amer. Math. Soc. 76 (1954), 127–148. [14] R.M. May, Stability and complexity in model ecosystems (2nd edition), Princeton University Press, 2001. [15] J.D. Murray, Mathematical Biology, Springer–Verlag New-York, 1989. [16] D.A. Neumann, Classification of continuous flows on 2–manifolds, Proc. Amer. Math. Soc. 48 (1975), 73–81. LESLIE-GROBER SYSTEMS 11 [17] M. Peixoto, Dynamical Systems, Proceedings of a Symposium held at the University of Bahia, Acad. Press, New York, 1973, pp.389–420. [18] H. Poincaré, Mémoire sur les courbes définies par les équations différentielles, J. Math. 37 (1881), 375–422; Oeuvres de Henri Poincaré, vol. I, Gauthier-Villars, Paris, 1951, pp 3–84. [19] L. Puchuri, E. González-Olivares and A. Rojas-Palma, Multistability in a Leslie– Gower–type predation model with a rational nonmonotonic functional response and generalist predators. Comput. Math. Methods 2 (2020), e1070, 18 pp. [20] A. Singh, M. Preeti and P. Malik, Hopf bifurcation and chaos in a Leslie–Gower prey-predator model with discrete delays, Int. J. Biomath. 13 (2020), 2050048, 27 pp. [21] J. Su, Degenerate Hopf bifurcation in a Leslie–Gower predator–prey model with preda- tor harvest, Adv. Difference Equ. 2020, Paper No. 194, 19 pp. [22] V. Tiwari, J.P Tripathi, R.K. Upadhyay, K. Ranjit, Y.P. Wu, J.S. Wang, and G.Q. Sun, Predator–prey interaction system with mutually interfering predator: role of feedback control, Appl. Math. Model. 87 (2020), 222–244. [23] D. Tsvetkov and R. Angelova-Slavova, Positive periodic solutions for periodic predator–prey systems of Leslie–Gower or Holling–Tanner type, Nonlinear Stud. 27 (2020), 991–1002. [24] P. Turchin, Complex population dynamics. A theoretical/empirical synthesis, Mono- graphs in Population Biology 35, Princeton University Press, 2003. [25] V. Volterra, Variazioni e fluttuazioni del numero d’individui in specie animali con- viventi, Memorie della R. Accademia dei Lincei, S.VI, IT 1926; II:31–113. [26] X. Wang, Y. Tan, Y. Cai and W. Wang, Impact of the fear effect on the stability and bifurcation of a Leslie–Gower predator–prey model, Internat. J. Bifur. Chaos Appl. Sci. Engrg. 30 (2020), 2050210, 13 pp. [27] F. Wu, Propagation threshold in an integrodifference predator–prey system of Leslie– Gower type, J. Difference Equ. Appl. 27 (2021), 26–40. [28] X.P. Yan and C.H Zhang, Stability of a delayed diffusive predator–prey model with prey harvesting of Michaelis-Menten type, Appl. Math. Lett. 114 (2021), 106904, 6 pp. [29] P. Ye and D. Wu, Impacts of strong Allee effect and hunting cooperation for a Leslie– Gower predator–prey system, Chinese J. Phys. 68 (2020), 49–64. [30] H. Zhao and D. Wu, Point to point traveling wave and periodic traveling wave induced by Hopf bifurcation for a diffusive predator–prey system, Discrete Contin. Dyn. Syst. Ser. S 13 (2020), 3271–3284. [31] R. Zou and S. Guo, Dynamics of a Leslie–Gower predator–prey system with cross- diffusion, Electron. J. Qual. Theory Differ. Equ. 2020, Paper No. 65, 33 pp. [32] R. Zou and S. Guo, Dynamics of a diffusive Leslie–Gower predator–prey model in spatially heterogeneous environment, Discrete Contin. Dyn. Syst. Ser. B 25 (2020), 4189–4210. [33] W.Q. Zuo, Z.P. Ma, Zhan-Ping and Z.B. Cheng, Spatiotemporal dynamics induced by Michaelis–Menten type prey harvesting in a diffusive Leslie–Gower predator-prey model, Internat. J. Bifur. Chaos Appl. Sci. Engrg. 30 (2020), 2050204, 24 pp. 1 Departament de Matemàtiques, Universitat Autònoma de Barcelona, 08193 Bellaterra, Barcelona, Catalonia, Spain Email address: jllibre@mat.uab.cat 2 Departamento de Matemática, Instituto Superior Técnico, Universidade de Lisboa, Av. Rovisco Pais 1049–001, Lisboa, Portugal Email address: cvalls@math.tecnico.ulisboa.pt