# Global Stability for a Predator-Prey Model with Dispersal among Patches

Authors: Yang Gao; Shengqiang Liu
Journal: Abstract and Applied Analysis
 (2014)

---

## Abstract

We investigate a predator-prey model with dispersal for both predator and prey amongn patches; our main purpose is to extend the global stability criteria by Li and Shuai (2010) on a predator-prey model with dispersal for prey among n patches. By using the method of constructing Lyapunov functions based on graph-theoretical approach for coupled systems, we derive sufficient conditions under which the positive coexistence equilibrium of this model is unique and globally asymptotically stable if it exists.


---

## Body

## 1. Introduction

In the literature of predator-prey population systems, both continuous reaction-diffusion systems and discrete patchy models are used to study the spatial heterogeneity [1, 2]; patchy models are often used to describe directed movement of population among niches or migration among habitats. It is naturally interesting problem to consider how the dispersal or migration of predator and prey influences the global dynamics of the interacting ecological system; thus patchy predator-prey model received lots of attentions [1, 3–6].Since the discrete patchy models usually involve high-dimensional system, it is rather mathematically challenging to study the uniqueness and stability of the positive equilibrium of the predator-prey patchy models, and the available global dynamics criteria in the literatures mainly focus on the special case of two-patch [3] or on the permanence and existence of periodic solutions [4–6].Recently, Li and Shuai [7] considered the following predator-prey model with dispersal for prey among n-patch:
(1)x˙i=xi(ri-bixi-eiyi)+∑j=1n‍dijx(xj-αijxxi),y˙i=yi(-γi-δiyi+εixi),i=1,…,n.
Here, xi, yi denote the densities of prey and predators on the patch i, respectively. The parameters ri, bi and γi, δi in the model are nonnegative constants. What is more, the parameters ei and εi in the model are positive constants. Constant dijx is the dispersal rate of the prey from patch j to patch i and constants αijx can be selected to represent different boundary conditions in the continuous diffusion case.In [7], the authors studied the global stability of the coexistence equilibrium of system (1), by considering (1) as a coupled n predator-prey submodels on networks. Using results from graph theory and a developed systematic approach that allows one to construct global Lyapunov functions for large-scale coupled systems from building blocks of individual vertex systems, Li and Shuai [7] obtain the following sharp results for (1).Proposition 1 (see [7, Theorem 6.1]).
Assume that(dijx)n×n is irreducible. If there exists k such that bk>0 or δk>0, then, whenever a positive equilibrium E* exists in (1), it is unique and globally asymptotically stable in the positive cone R2n+.Although well-improved results have been seen in the above work on dispersal predator-prey model, such models are not well studied yet in the sense that model (1) assumes no dispersal for predator, which is not realistic in many cases [1, 3]. Thus it is interesting for us to consider the global stability of the positive equilibrium for predator-prey model with dispersal for both predator and prey.Motivated by the above work in [7], in this paper we generalize model (1) into the following predator-prey model with dispersal for both predator and prey:
(2)x˙i=xi(ri-bixi-eiyi)+∑j=1n‍dijx(xj-αijxxi)y˙i=yi(-γi-δiyi+εixi)+∑j=1n‍dijy(yj-αijyyi),0000000000000000000000000000i=1,…,n.
Here, the parameters ri, bi, ei, γi, δi, and εi are defined the same as those in (1). The nonnegative constants dijy, αijy, and dijy are the dispersal rate of the predators from patch j to patch i, and αijy represents the different boundary conditions in the continuous diffusion case. Clearly, when dijy=0 for all i,j=1,…,n, model (2) directly reduces to (1); thus our model (2) directly extends model (1) in [7].The main purpose of this paper is to obtain the global stability for the coexistence equilibrium of (2). We will engage the techniques of constructing Lyapunov function based on graph-theory which were well developed by Li et al. in [7–9]; we refer to [10–12] for recent applications. Our study seems to be the first attempt in applying the network method for coupled network systems of differential equations to address the predator-prey system with dispersal for both predator and prey among patches. Networked method has been extensively investigated in the several fields. For example, multiagent systems can be seen as complicated network systems. A lot of researchers take their interest in flocking and consensus of the multiagent systems [13–17]. What is more, neural network systems can be seen as complicated network systems. Over the past few decades, various neural network models have been extensively investigated [18–20].A mathematical description of a network is a directed graph consisting of vertices and directed arcs connecting them. At each vertex, the local dynamics are given by a system of differential equations called the vertex system. The directed arcs indicate interconnections and interactions among vertex systems.A digraphG with n vertices for the system (2) can be constructed as follows. Each vertex represents a patch and (j,i)∈E(G) if and only if dijx,dijy>0. At each vertex of G, the vertex dynamics is described by a predator-prey system. The coupling among these predator-prey systems is provided by dispersal of predator and prey among patches.This paper is organized as follows. In the next section, we introduce preliminaries results on graph-theory based on coupled network models. In Section3, we obtain the main result of system (2). This is followed by a brief conclusion section.


## 2. Preliminaries

In this section, we will list some definitions and Theorems that we will use in the later sections.A directed graph or digraphG=(V,E) contains a set V={1,2,…,n} of vertices and a set E of arcs (i,j) leading from initial vertex i to terminal vertex j. A subgraph H of G is said to be spanning if H and G have the same vertex set. A digraph G is weighted if each arc (j,i) is assigned a positive weight. aij>0 if and only if there exists an arc from vertex j to vertex i in G.The weightw(H) of a subgraph H is the product of the weights on all its arcs. A directed path P in G is a subgraph with distinct vertices i1,i2,…,im such that its set of arcs is {(ik,ik+1):k=1,2,…,m}. If im=i1, we call P a directed cycle.A connected subgraphT is a tree if it contains no cycles, directed or undirected.A treeT is rooted at vertex i, called the root, if i is not a terminal vertex of any arcs, and each of the remaining vertices is a terminal vertex of exactly one arc. A subgraph Q is unicyclic if it is a disjoint union of rooted trees whose roots form a directed cycle.Given a weighted digraphG with n vertices, define the weight matrix A=(aij)n×n whose entry aij equals the weight of arc (j,i) if it exists, and 0 otherwise. For our purpose, we denote a weighted digraph as (G,A). A digraph G is strongly connected if for any pair of distinct vertices, there exists a directed path from one to the other. A weighted digraph (G,A) is strongly connected if and only if the weight matrix A is irreducible.The Laplacian matrix of(G,A) is denoted by L. Let ci denote the cofactor of the ith diagonal element of L. The following results are listed as follows from [7].Proposition 2 (see [7]).
Assumen≥2. Then
(3)ci=∑T∈Ti‍w(T),
where Ti is the set of all spanning trees T of (G,A) that are rooted at vertex i, and w(T) is the weight of T. In particular, if (G,A) is strongly connected, then ci>0 for 1≤i≤n.Theorem 3 (see [7]).
Assumen≥2. Let ci be given in Proposition 2. Then the following identity holds:
(4)∑i,j=1n‍ciaijFij(xi,xj)=∑Q∈Q‍w(Q)∑(s,r)∈E(CQ)‍Frs(xr,xs),
where Fij(xi,xj), 1≤i, j≤n, are arbitrary functions, Q is the set of all spanning unicyclic graphs of (G,A), w(Q) is the weight of Q, and CQ denotes the directed cycle of Q.Given a network represented by digraphG with n vertices, n≥2, a coupled system can be built on G by assigning each vertex its own internal dynamics and then coupling these vertex dynamics based on directed arcs in G. Assume that each vertex dynamics is described by a system of differential equations
(5)ui′=fi(t,ui),
where ui∈Rmi and fi:R×Rmi→Rmi. Let gij:R×Rmi×Rmj→Rmi represent the influence of vertex j on vertex i, and let gij≡0 if there exists no arc from j to i in G. Then we obtain the following coupled system on graph G:
(6)ui′=fi(t,ui)+∑j=1n‍gij(t,ui,uj),i=1,2,…,n.
Here functions fi,   gij are such that initial-value problems have unique solutions.We assume that each vertex system has a globally stable equilibrium and possesses a global Lyapunov functionVi.Theorem 4 (see [7]).
Assume that the following assumptions are satisfied.(1)
There exist functionsVi(t,ui), Fij(t,ui,uj) and constants aij≥0 such that
(7)V˙i(t,ui)≤∑i,j=1n‍aijFij(t,ui,uj),t>0,ui∈Di.(2)
Along each directed cycleC of the weighted digraph (G,A), A=(aij),
(8)∑(s,r)∈E(C)‍Frs(t,ur,us)≤0.(3)
Constantsci are given by the cofactor of the ith diagonal element of L.
Then the functionV(t,u)=∑i=1n‍ciVi(t,ui) satisfies V˙(t,u)≤0 for t>0,  u∈D; namely, V is a Lyapunov function for the system (6).


## 3. Main Results

In this section, the stability for the positive equilibrium of then-patch predator-prey model (2) is considered. We regard (2) as a coupled system on a network. Using a Lyapunov function for the n-patch predator-prey model with dispersal and Theorem 4 of Section 2, we will establish that a positive equilibrium of the n-patch predator-prey model (2) with dispersal is globally asymptotically stable in R+2n as long as it exists.First of all, we will give a lemma for the system (2).Lemma 5.
The setR+2n is the positive invariant set for the system (2).The next Theorem gives the globally asymptotically stable condition for the positive equilibrium of the system (2).Theorem 6.
Assume that a positive equilibriumE*=(x1*,y1*,x2*,y2*,…,xn*,yn*) exists for the system (2) and the following assumptions hold.(1)
Dispersal matrixes(dijx)n×n, (dijy)n×n are irreducible; moreover there exists k such that bk>0 or δk>0.(2)
There exists nonnegative constantλ such that λ·dijxεixj*=dijyeiyj* for 1≤i, j≤n, or dijxεixj*=λ·dijyeiyj* for 1≤i, j≤n.
Then, the positive equilibriumE* is unique and globally asymptotically stable in R+2n.Proof.
Let(9)Zi1(xi,yi)=ri-bixi-eiyi,Zi2(xi,yi)=-γi-δiyi+εixi.
In the sequel, we have
(10)Zi1(xi*,yi*)=-1xi*∑j=1n‍dijx(xj*-αijxxi*),Zi2(xi*,yi*)=-1yi*∑j=1n‍dijy(yj*-αijyyi*).
Set Lyapunov functions as
(11)Vi(xi,yi)=εi(xi-xi*-xi*ln⁡xixi*)+ei(yi-yi-yi*ln⁡yiyi*).
Direct differentiating Vi along the system (2), we have
(12)V˙i(xi,yi)=εi(xi-xi*)[Zi1(xi,yi)-Zi1(xi*,yi*)]+ei(yi-yi*)[Zi2(xi,yi)-Zi2(xi*,yi*)]+εi(xi-xi*)Zi1(xi*,yi*)+εi(xi-xi*)xi*∑j=1n‍dijx(xj-αijxxi)+ei(yi-yi*)Zi2(xi*,yi*)+ei(yi-yi*)yi*∑j=1n‍dijy(yj-αijyyi)=εi(xi-xi*)[Zi1(xi,yi)-Zi1(xi*,yi*)]+ei(yi-yi*)[Zi2(xi,yi)-Zi2(xi*,yi*)]+∑j=1n‍dijxεixj*Fijx(xi,xj)+∑j=1n‍dijyeiyj*Fijy(yi,yj)=-εibi(xi-xi*)2-εi(xi-xi*)ei(yi-yi*)-eiδi(yi-yi*)2+ei(yi-yi*)εi(xi-xi*)+∑j=1n‍dijxεixj*Fijx(xi,xj)+∑j=1n‍dijyeiyj*Fijy(yi,yj)=-εibi(xi-xi*)2-eiδi(yi-yi*)2+∑j=1n‍dijxεixj*Fijx(xi,xj)+∑j=1n‍dijyeiyj*Fijy(yi,yj),
where
(13)Fijx(xi,xj)=xjxj*-xixi*+1-xi*xjxixj*,Fijy(yi,yj)=yjyj*-yiyi*+1-yi*yjyiyj*.
Set aijx=dijxεixj*, bijy=dijyeiyj*,A=(aijx)n×n, and B=(bijy)n×n. One has
(14)Gix(xi)=-xixi*+ln⁡xixi*,Giy(yi)=-yiyi*+ln⁡yiyi*.
Next, we have two cases to consider.
Case  I.
d
i
j
x
ε
i
x
j
*
=
λ
·
d
i
j
y
e
i
y
j
* for 1≤i, j≤n.
Case  II.
λ
·
d
i
j
x
ε
i
x
j
*
=
d
i
j
y
e
i
y
j
* for 1≤i, j≤n.
For Case I, from the fact thataijx=dijxεixj* and bijy=dijyeiyj*, we obtain that aijx=λbijy; thus A=λ·B. Then we obtain that
(15)V˙i(xi,yi)≤-εibi(xi-xi*)2-eiδi(yi-yi*)2+∑j=1n‍aijx(Gix(xi)-Gjx(xj))+∑j=1n‍aijx(1-xi*xjxixj*+ln⁡xi*xjxixj*)+∑j=1n‍bijy(Giy(yi)-Gjy(yj))+∑j=1n‍bijy(1-yi*yjyiyj*+ln⁡yi*yjyiyj*)≤-εibi(xi-xi*)2-eiδi(yi-yi*)2+λ∑j=1n‍bijy(Gix(xi)-Gjx(xj))+λ∑j=1n‍bijy(1-xi*xjxixj*+ln⁡xi*xjxixj*)+∑j=1n‍bijy(Giy(yi)-Gjy(yj))+∑j=1n‍bijy(1-yi*yjyiyj*+ln⁡yi*yjyiyj*).
Let ciy denote the cofactor of the ith diagonal element of the matrix B. From the irreducible character of matrix B, we have ciy>0.
Furthermore, set Lyapunov functions as(16)V(x,y)=V(x1,y1,…,xn,yn)=∑i=1n‍ciyVix(xi)+∑i=1n‍ciyViy(yi).
Then differentiating V along the solution of the system (2), we obtain that
(17)V˙(x,y)≤-∑i=1n‍ciyεibi(xi-xi*)2-∑i=1n‍ciyeiδi(yi-yi*)2+∑i,j=1n‍λbijyciy(Gix(xi)-Gjx(xj))+∑i,j=1n‍bijyciy(Giy(yi)-Gjy(yj)).
Let G represent the directed graph associated with matrix B. Then G has vertices 1,2,…,n with a directed arc (k,j) from k to j if and only if bkjy≠0. Then E(G) is the set of all directed arcs of G. By Kirchhoff’s Matrix-Tree Theorem (see Proposition 2) we know that υk=Ckk can be expressed as a sum of weights of all directed spanning subtrees T of G that are rooted at vertex k. Thus, each term in υkakj is the weight ω(Q) of a unicyclic subgraph Q of G obtained from such a tree T by adding a directed arc (k,j) from the root k to vertex j. Because the arc (k,j) is a part of the unique cycle CQ of Q and that the same unicyclic graph Q can be formed when each arc of CQ is added to a corresponding rooted tree T, then the double sum can be expressed as a sum over all unicyclic subgraphs Q containing vertices 1,2,…,n. Therefore, following from the irreducible character of matrix B and Theorem 2.3 in [7], we obtain
(18)∑i,j=1n‍λbijyciy(Gix(xi)-Gjx(xj))=0,∑i,j=1n‍bijyciy(Giy(yi)-Gjy(yj))=0.
Combining with the fact that 1-a+ln⁡a≤0, therefore we have
(19)V˙(x,y)≤0.
When we consider V˙(x,y)=0, by condition 1, there exists k∈N+ such that
(20)(xk-xk*)2=0or(yk-yk*)2=0.
It means that xk=xk*  or yk=yk*.
Ifi and k can be connected with an arc from k to i in G, then we have aiky>0 and biky>0. Furthermore,
(21)1-xi*xkxixk*+ln⁡xi*xkxixk*=0,1-yi*ykyiyk*+ln⁡yi*ykyiyk*=0.
Because of 1-a+ln⁡a≤0 and 1-a+ln⁡a=0,⇔a=0. we deduce that
(22)xixi*=xkxk*,yiyi*=ykyk*.
From xk=xk*, or yk=yk*, we obtain that xi=xi* and yi/yi*=yk/yk* or yi=yi* and xi/xi*=xk/xk*.
By condition 1 and the definition of matrixesA, B, we get that B are irreducible. By strong connectivity of G, there exists a directed path P from any i to k. Then we have that, for any i=1,2,…,n, there must be
(23)xi=xi*,yiyi*=μ,μ≥0,
or for any i=1,2,…,n, there must be
(24)yi=yi*,xixi*=μ,μ≥0.
Next, we will prove that the largest compact invariant subset of {(x,y)∣V˙(x,y)=0} is the singleton {E*}.
We only consider the case that(25)xi=xi*,yiyi*=μ,i=1,2,…,n,μ≥0.
The case that
(26)yi=yi*,xixi*=μ,i=1,2,…,n,μ≥0
is similar to this case. So we omit it.
Ifμ=0, we have yi=0 for any i=1,2,…,n, and then we have
(27)xi*(ri-bixi*)+∑j=1n‍dijx(xj*-αijxxi*)=0,
which contradicts to the fact that
(28)xi*(ri-bixi*-eiyi*)+∑j=1n‍dijx(xj*-αijxxi*)=0.
Ifμ>0 and μ≠1, we have yi=μyi* for any i=1,2,…,n, and then we have
(29)xi*(ri-bixi*-eiμyi*)+∑j=1n‍dijx(xj*-αijxxi*)=0,
which also contradicts to the fact that
(30)xi*(ri-bixi*-eiyi*)+∑j=1n‍dijx(xj*-αijxxi*)=0.
Therefore, we obtain that μ=1, which means
(31)xi=xi*,yi=yi*,i=1,2,…,n.
Namely, we get that the largest compact invariant subset of {(x,y)∣V˙(x,y)=0} is the singleton {E*}. Therefore, by the LaSalle Invariance Principle ([21]), E* is globally asymptotically stable in R+2n.
With the similar arguments to the Case I, we can prove thatE* is globally asymptotically stable in R+2n for Case II. This completes the proof.Remark 7.
Theorem6 is applicable to model (1): consider model (2) with dijy=0, i,j=1,…,n, and let λ=0; thus Theorem 6 directly reduces to Proposition 1 by Li and Shuai [7] for (1).By Theorem6 and similar arguments to Remark 7, we directly have the following global stability theorem for the predator-prey model with discrete dispersal of predator among patches.Corollary 8.
Consider the model(32)x˙i=xi(ri-bixi-eiyi),y˙i=yi(-γi-δiyi+εixi)+∑j=1n‍dijy(yj-αijyyi),000000000000000000000000i0000i=1,…,n.
Assume that the matrix (dijy)n×n is irreducible. If there exists k such that bk>0 or δk>0; then, whenever a positive equilibrium E* exists in (32), it is unique and globally asymptotically stable in the positive cone R+2n.


## 4. Discussion

In this paper, we generalize the model of then-patch predator-prey model of [7] to the general model (2) that both the prey and the predator have dispersal among n-patches. Based on the network method for coupled systems of differential equations developed in [7–9], we prove that the positive equilibrium of (2) is globally asymptotically stable given some conditions on the coupling (see Theorem 6). Our main theorem generalizes Theorem 6.1 in [7] and our results also cover the other case of (2) in that only the predators disperse among patches.Biologically, our result of Theorem6 implies that if predator-prey system is dispersing among strongly connected patches (which is equivalent to the irreducibility of the dispersal matrixes of predator and prey) and if the system is permanent (which guarantees the existence of positive equilibrium), then the numbers of both predators and prey in each patches will eventually be stable at some corresponding positive values given the well-coupled dispersal (condition 2 of Theorem 6).We remark that our Theorem6 requires the extra condition 2 for the coupling dispersal coefficients and that the global dynamics for the coexistence equilibrium of (2) without condition 2 of Theorem 6 are still unclear. It remains an interesting future problem for the patchy dispersal predator-prey model.



---

*Source: 176493-2014-03-12.xml*