Dynamic Behaviors of a Nonautonomous Discrete Predator-Prey System Incorporating a Prey Refuge and Holling Type II Functional Response

Abstract

A nonautonomous discrete predator-prey system incorporating a prey refuge and Holling type II functional response is studied in this paper. A set of sufficient conditions which guarantee the persistence and global stability of the system are obtained, respectively. Our results show that if refuge is large enough then predator species will be driven to extinction due to the lack of enough food. Two examples together with their numerical simulations show the feasibility of the main results.

1. Introduction

As was pointed out by Berryman [1], the dynamic relationship between predator and prey has long been and will continue to be one of the dominant themes in both ecology and mathematical ecology due to its universal existence and importance. Furthermore, the study of the consequences of the hiding behavior of the prey on the dynamics of predator-prey interactions can be recognized as a major issue in both applied mathematics and theoretical ecology [2]. In general, the effects of prey refuges on the population dynamics are very complex in nature, but for modeling purposes, it can be considered as constituted by two components [2]. The first one, which affects positively the growth of prey and negatively that of predators, comprises the reduction of prey mortality due to the decrease in predation success. The second one may be the tradeoffs and by-products of the hiding behavior of prey which could be advantageous or detrimental for all the interacting populations [3].

Sih [4] obtained a set of general conditions which ensure that the refuge use has a stabilizing effect on Lotka-Volterra-type predator-prey systems; he also examined the effect of the cost of refuge use in decreased prey feeding or reproductive rate. In [5], González-Olivares and Ramos-Jiliberto investigated the dynamic behaviors of predator-prey system incorporating Holling type II functional response and a constant refuge: 𝑥̇(𝑡)𝑥=𝑟𝑥1−𝐾−𝛽(1−𝑚)𝑥𝑦(),𝑦̇(𝑡)()1−𝑚𝑥+𝑎=−𝑑𝑦+𝑐𝛽1−𝑚𝑥𝑦(),1−𝑚𝑥+𝑎 (1.1) where 𝑥(𝑡),𝑦(𝑡) denote the densities of prey and predator population at any time 𝑡, respectively; 𝑐,𝑑,𝑘,𝑟,𝛽,𝑎,and𝑚 are positive constants; here 𝑟 is the intrinsic per capita growth rate of prey; 𝐾 is the prey environmental carrying capacity; 𝛽 is the maximal per capita consumption rate of predators; 𝑎 is the amount of prey needed to achieve one-half of 𝛽; 𝑐 is the conversion factor denoting the number of newly born predators for each captured prey; 𝑑 is the death rate of the predator; 𝑥𝑚 is the number of prey that refuge can protect at time 𝑡. Kar [6] also studies the dynamic behaviors of system (1.1). He obtained the conditions for the existence and stability of the equilibria and persistent criteria for the system. He also shows that the system admits a unique limit cycle when the positive equilibria is unstable. In these papers, all their finds indicate that the refuge influencing the dynamic behavior of predator-prey system greatly and increasing the amount of refuge could increase prey density and lead to population outbreaks. Kar [7] also studied the influence of harvesting on a system with prey refuge.

Some scholars argued that the nonautonomous case is more realistic, because many biological or environmental parameters do subject to fluctuate with time; thus more complex equations should be introduced. Already, many scholars [8–15] studied the dynamic behaviors of nonautonomous predator-prey system incorporating prey refuge. Recently, Xu and Jia [11] proposed and studied the nonautonomous predator-prey system incorporating prey refuge and Holling type II functional response, that is, 𝑥̇(𝑡)=(𝑎(𝑡)(𝑡)𝑥(𝑡𝑥(𝑡)−𝛽(𝑡(𝑡𝑥(𝑡)𝑦(𝑡)−𝑏)))(1−𝑚))(𝑡(𝑡𝑥(𝑡),𝑦̇(𝑡)=(𝑡)(𝑡)𝑦(𝑡)𝑦(𝑡)+𝑐(𝑡)𝛽(𝑡(𝑡𝑥(𝑡)𝑦(𝑡)1+𝛼)(1−𝑚))−𝛾−𝑑)(1−𝑚))(𝑡(𝑡𝑥(𝑡),1+𝛼)(1−𝑚)) (1.2) where 𝑥(𝑡) and 𝑦(𝑡) denote the density of prey and predator populations at time 𝑡, respectively; 𝑥(𝑡)𝑚(𝑡) denotes the number of prey that the refuge can protect at time 𝑡; 𝑎(𝑡),𝑏(𝑡),𝑐(𝑡),𝑑(𝑡),𝛼(𝑡),𝛽(𝑡),𝛾(𝑡),and𝑚(𝑡)(0≤𝑚(𝑡)<1) are nonnegative continuous function that have the upper and lower bounds.

Though most dynamic behaviors of population models are based on the continuous models governed by differential equations, the discrete time models are more appropriate than the continuous ones when the size of the population is rarely small or the population has nonoverlapping generations [12]. It has been found that the dynamic behaviors of the discrete system is rather complex and contains more rich dynamics than the continuous ones [16]. Though the influence of prey refuge for continuous model has been extensively investigated, seldom did scholars investigated the influence of prey refuge for discrete predator-prey system. To the best of the authors’ knowledge, to this day, only Zhuang and Wen [17] studied the local property stability of the fixed points of the discrete Leslie-Gower predator-prey systems with and without Allee effect. In this paper, we study the corresponding discrete prey-predator system of (1.2): 𝑥()(𝑛)𝑎(𝑛)(𝑛)𝑥(𝑛)−𝛽(𝑛(𝑛𝑦(𝑛)𝑛+1=𝑥exp−𝑏)(1−𝑚))(𝑛(𝑛𝑥(𝑛),𝑦()(𝑛)(𝑛)(𝑛)𝑦(𝑛)+𝑐(𝑛)𝛽(𝑛(𝑛𝑥(𝑛)1+𝛼)(1−𝑚))𝑛+1=𝑦exp−𝛾−𝑑)(1−𝑚))(𝑛(𝑛𝑥(𝑛).1+𝛼)(1−𝑚)) (1.3) Here, we assume that 𝑎(𝑛),𝑏(𝑛),𝑐(𝑛),𝑑(𝑛),𝛼(𝑛),𝛽(𝑛),and𝑚(𝑛)(0≤𝑚(𝑛)<1),𝛾(𝑛) are all bounded nonnegative sequences. Noting that 𝑥(𝑛)=𝑥(0)exp^𝑛−1_𝑘=0𝑎(𝑘)(𝑘)𝑥(𝑘)−𝛽(𝑘(𝑘𝑦(𝑘)−𝑏)(1−𝑚))(𝑘(𝑘𝑥(𝑘),𝑦(𝑛)(0)1+𝛼)(1−𝑚))=𝑦exp^𝑛−1_𝑘=0(𝑘)(𝑘)𝑦(𝑘)+𝑐(𝑘)𝛽(𝑘(𝑘𝑥(𝑘)−𝛾−𝑑)(1−𝑚))(𝑘(𝑘𝑥(𝑘).1+𝛼)(1−𝑚)) (1.4)

For the point of view of biology, in the sequel, we assume that 𝑥(0)>0,𝑦(0)>0, then from (1.4), we know that the solutions of system (1.3) are positive. From now on, for any bounded sequence 𝑥(𝑛), 𝑥^𝑢=sup_𝑛∈𝑁𝑥(𝑛),𝑥_𝑙=inf_𝑛∈𝑁𝑥(𝑛). (1.5)

2. Permanence

We will investigate the persistent property of the system in this section.

Lemma 2.1 (see [12]). Assume that {𝑥(𝑘)}>0 and 𝑥(𝑘+1)≤𝑥(𝑘)exp{𝑎(𝑘)−𝑏(𝑘)𝑥(𝑘)} (2.1) for 𝑘∈𝑁, where 𝑎(𝑘) and 𝑏(𝑘) are nonnegative sequences bounded above and below by positive constants. Then limsup_𝑘→∞𝑥(𝑘)≤1𝑏_𝑙exp𝑎^𝑢−1. (2.2)

Lemma 2.2 (see [12]). Assume that {𝑥(𝑘)} satisfies 𝑥(𝑘+1)≥𝑥(𝑘)exp{𝑎(𝑘)−𝑏(𝑘)𝑥(𝑘)},𝑘≥𝑁₀, (2.3)limsup_𝑘→∞𝑥(𝑘)≤𝑥^∗, and 𝑥(𝑁₀)>0, where 𝑎(𝑘) and 𝑏(𝑘) are nonnegative sequences bounded above and below by positive constants and 𝑁₀∈𝑁. Then liminf_𝑘→∞𝑥(𝑘)𝑎≥min^𝑙𝑏^𝑢exp𝑎_𝑙−𝑏^𝑢𝑥^∗,𝑎_𝑙𝑏^𝑢. (2.4)

Theorem 2.3. Every positive solution (𝑥(𝑛),𝑦(𝑛)) of system (1.3) satisfies limsup_𝑛→∞𝑥(𝑛)≤𝑀₁,limsup_𝑛→∞𝑦(𝑛)≤𝑀₂. (2.5) Here 𝑀₁=(1/𝑏_𝑙)exp(𝑎^𝑢−1),𝑀₂=(1/𝑑_𝑙)exp(−𝛾_𝑙+(𝑐^𝑢𝛽^𝑢/𝛼_𝑙)−1).

Proof. Let (𝑥(𝑛),𝑦(𝑛)) be any positive solution of system (1.3). From the first equation of system (1.3) it follows that 𝑥(𝑛+1)≤𝑥(𝑛)exp[𝑎(𝑛)−𝑏(𝑛)𝑥(𝑛)]. (2.6) Applying Lemma 2.1 to (2.6) leads to limsup_𝑛→∞𝑥(𝑛)≤1𝑏_𝑙exp𝑎^𝑢−1=^def𝑀₁. (2.7) From the second equation of system (1.3), similarly to the analysis of (2.6)-(2.7), we can obtain limsup_𝑛→∞𝑦(𝑛)≤1𝑑_𝑙exp−𝛾_𝑙+𝑐^𝑢𝛽^𝑢𝛼_𝑙−1=^def𝑀₂. (2.8)
This ends the proof of Theorem 2.3.

Theorem 2.4. Assume that inequalities 𝑎_𝑙>𝛽^𝑢𝑀₂,𝛾^𝑢<𝑐_𝑙𝛽_𝑙1−𝑚^𝑢𝑚₁1+𝛼^𝑢1−𝑚_𝑙𝑀₁ 𝐻₁ hold. Let (𝑥(𝑛),𝑦(𝑛)) be any positive solution of system of (1.3), then liminf_𝑛→∞𝑥(𝑛)≥𝑚₁,liminf_𝑛→∞𝑦(𝑛)≥𝑚₂. (2.9) Here 𝑚₁=𝑎_𝑙−𝛽^𝑢𝑀₂𝑏^𝑢exp𝑎_𝑙−𝑏^𝑢𝑀₁−𝛽^𝑢𝑀₂,𝑚₂=−𝛾^𝑢+𝐵𝑑^𝑢exp−𝛾^𝑢+𝐵−𝑑^𝑢𝑀₂,𝑐𝐵=^𝑙𝛽_𝑙1−𝑚^𝑢𝑚₁1+𝛼^𝑢1−𝑚_𝑙𝑀₁. (2.10)

Proof. According to the first inequality of 𝐻₁, one could choose 𝜀>0 small enough, such that the inequality 𝑎_𝑙−𝛽^𝑢𝑀₂+𝜀>0. (2.11) holds. For the above 𝜀>0, according to Theorem 2.3, there exists an integer 𝑛^∗∈𝑁 such that for all 𝑛≥𝑛^∗, 𝑥(𝑛)≤𝑀₁+𝜀,𝑦(𝑛)≤𝑀₂+𝜀. (2.12) For 𝑛≥𝑛^∗, from (2.12) and the first equation of system (1.3), we have 𝑥()(𝑛)𝑎(𝑛)(𝑛)𝑥(𝑛)−𝛽(𝑛(𝑛𝑦(𝑛)𝑛+1=𝑥exp−𝑏)(1−𝑚))(𝑛(𝑛𝑥(𝑛)(𝑛)𝑎(𝑛)(𝑛)𝑥(𝑛)−𝛽(𝑛(𝑛𝑀1+𝛼)(1−𝑚))≥𝑥exp−𝑏)(1−𝑚))₂+𝜀(𝑛(𝑛𝑥(𝑛)(𝑛)𝑎(𝑛)(𝑛)𝑥(𝑛)(𝑛(𝑛𝑀1+𝛼)(1−𝑚))≥𝑥exp−𝑏−𝛽)(1−𝑚))₂(𝑛)𝑎+𝜀≥𝑥exp_𝑙−𝑏^𝑢𝑥(𝑛)−𝛽^𝑢(1−𝑚(𝑛))𝑀₂(𝑛)𝑎+𝜀≥𝑥exp_𝑙−𝑏^𝑢𝑥(𝑛)−𝛽^𝑢𝑀₂+𝜀. (2.13) As a direct corollary of Lemma 2.2, according to (2.7) and (2.13), one has liminf_𝑛→∞𝑥(𝑛)≥min𝐴_1𝜀,𝐴_2𝜀, (2.14) where 𝐴_1𝜀=𝑎_𝑙−𝛽^𝑢𝑀₂+𝜀𝑏^𝑢,𝐴_2𝜀=𝐴_1𝜀exp𝑎_𝑙−𝑏^𝑢𝑀₁+𝜀−𝛽^𝑢𝑀₂+𝜀. (2.15) Noting that 𝑀₁=exp𝑎^𝑢−1𝑏_𝑙≥𝑎^𝑢𝑏_𝑙≥𝑎_𝑙𝑏^𝑢, (2.16) and so 𝑎_𝑙−𝑏^𝑢𝑀₁≤0, consequently, for arbitrary 𝜀, 𝑎_𝑙−𝑏^𝑢𝑀₁+𝜀−𝛽^𝑢𝑀₂+𝜀≤0. (2.17) The above inequality leads to 𝐴_1𝜀≥𝐴_2𝜀. (2.18) Thus, liminf_𝑛→∞𝑥(𝑛)≥𝐴_2𝜀. (2.19) Letting 𝜀→∞, it follows that liminf_𝑛→∞𝑥(𝑛)≥𝐴₁=𝑎_𝑙−𝛽^𝑢𝑀₂𝑏^𝑢exp𝑎_𝑙−𝑏^𝑢𝑀₁−𝛽^𝑢𝑀=₂^def𝑚₁. (2.20) According to (2.7), (2.8), and (2.20), for any 𝜀>0, there exists 𝑛₁∈𝑁, such that for all 𝑛≥𝑛₁, 𝑚₁−𝜀≤𝑥(𝑛)≤𝑀₁+𝜀,𝑦(𝑛)≤𝑀₂+𝜀. (2.21) Similarly to the analysis of (2.13)–(2.20), by applying (2.21), from the second equation of system (1.3), it follows that liminf_𝑛→∞𝑦(𝑛)≥−𝛾^𝑢+𝐵𝑑^𝑢exp−𝛾^𝑢+𝐵−𝑑^𝑢𝑀=₂^def𝑚₂. (2.22) Here 𝐵=^def𝑐_𝑙𝛽_𝑙(1−𝑚^𝑢)𝑚₁/(1+𝛼^𝑢(1−𝑚_𝑙)𝑀₁).
This completes the proof of Theorem 2.4.

3. Global Stability

In this section, by developing the analysis technique of [18], we obtain the conditions which guarantee the global stability of the system (1.3).

Theorem 3.1. Assume that 𝐻₁ holds, assume further that 𝜆₁=max1−𝑀₁𝑝₁,1−𝑚₁𝑝₂+𝛽^𝑢1−𝑚_𝑙𝑀₂1+𝛼_𝑙(1−𝑚^𝑢)𝑚₁𝜆<1,₂=max1−𝑑^𝑢𝑀₂,1−𝑑_𝑙𝑚₂+𝑐^𝑢𝛽^𝑢1−𝑚_𝑙𝑀₁1+𝛼_𝑙(1−𝑚^𝑢)𝑚₁²<1. (3.1) Then for any two positive solutions (𝑥(𝑛),𝑦(𝑛)) and (𝑥(𝑛),𝑦(𝑛)) of system (1.3), one has lim_𝑛→∞(𝑥(𝑛)−𝑥(𝑛))=0,lim_𝑛→∞𝑦(𝑛)−𝑦(𝑛)=0. (3.2) Here, 𝑝=₁^def𝑏^𝑢−(𝛽_𝑙𝛼_𝑙𝑚₂(1−𝑚^𝑢)²/[1+𝛼^𝑢𝑀₁(1−𝑚_𝑙)]²),𝑝=₂^def𝑏_𝑙−(𝛽^𝑢𝛼^𝑢𝑀₂(1−𝑚_𝑙)²/[1+𝛼_𝑙𝑚₁(1−𝑚^𝑢)]²).

Proof. Let 𝑥(𝑛)=𝑥(𝑛)exp(𝑢(𝑛)),𝑦(𝑛)=𝑦(𝑛)exp(𝑣(𝑛)). (3.3) then system (1.3) is equivalent to 𝑢()(𝑛)=𝑏(𝑛)−𝛽(𝑛)𝛼(𝑛(𝑛𝑛+1−𝑢)(1−𝑚))²𝑦(𝑛)𝛿(𝑛)(𝑛𝑥(𝑛)(𝑢(𝑛+𝛽(𝑛(𝑛1+(1+𝛿))1−exp)))(1−𝑚))𝛿(𝑛)𝑦(𝑛)(𝑣(𝑛,𝑣()(𝑛)𝑐(𝑛)𝛽(𝑛(𝑛1+1−exp))𝑛+1−𝑣=−)(1−𝑚))𝛿(𝑛)((𝑛𝑥(𝑛)(𝑢(𝑛(𝑛)𝑦(𝑛)(𝑣(𝑛.1+1+𝛿))1−exp))+𝑑1−exp)) (3.4) Here 𝛿(𝑛)=^def𝛼(𝑛)(1−𝑚(𝑛))𝑥(𝑛),𝛿(𝑛)=^def𝛼(𝑛)(1−𝑚(𝑛))𝑥(𝑛). By using the mean-value theorem, it follows that 𝑢()(𝑛𝑏(𝑛)−𝛽(𝑛)𝛼(𝑛(𝑛𝑛+1=𝑢)1−)(1−𝑚))²𝑦(𝑛)𝛿(𝑛)((𝑛𝑥(𝑛)𝜃1+1+𝛿))exp₁(𝑛)𝑢(𝑛−𝛽(𝑛(𝑛))(1−𝑚))𝛿(𝑛)𝑣(𝑛)𝑦(𝑛)𝜃1+exp₂(𝑛)𝑣(𝑛),𝑣()(𝑛)𝑐(𝑛)𝛽(𝑛(𝑛𝑥(𝑛)𝜃𝑛+1=𝑢)(1−𝑚))exp₄(𝑛)𝑢(𝑛)𝛿(𝑛)((𝑛(𝑛)(𝑛)𝑦(𝑛)𝜃1+1+𝛿))+𝑣1−𝑑exp₃(𝑛)𝑣(𝑛), (3.5) where 𝜃_𝑖(𝑛)∈(0,1)(𝑖=1,2,3,4). To complete the proof, it suffices to show that lim_𝑛→∞𝑢(𝑛)=lim_𝑛→∞𝑣(𝑛)=0. (3.6) In view of (3.1), we can choose 𝜀>0 small enough such that 𝜆^𝜀₁=max1−𝑀₁+𝜀𝑝^𝜀₁,1−𝑚₁−𝜀𝑝^𝜀₂+𝛽^𝑢1−𝑚_𝑙𝑀₂+𝜀1+𝛼_𝑙(1−𝑚^𝑢)𝑚₁𝜆−𝜀<1,_𝜀2=max1−𝑑^𝑢𝑀₂+𝜀,1−𝑑_𝑙𝑚₂+𝑐−𝜀^𝑢𝛽^𝑢1−𝑚_𝑙𝑀₁+𝜀1+𝛼_𝑙(1−𝑚^𝑢)𝑚₁−𝜀²<1. (3.7) Here, 𝑝=^𝜀₁^def𝑏^𝑢−(𝛽_𝑙𝛼_𝑙(𝑚₂−𝜀)(1−𝑚^𝑢)²/[1+𝛼^𝑢(𝑀₁+𝜀)(1−𝑚_𝑙)]²),𝑝=^𝜀₂^def𝑏_𝑙−(𝛽^𝑢𝛼^𝑢(𝑀₂+𝜀)(1−𝑚_𝑙)²/[1+𝛼_𝑙(𝑚₁−𝜀)(1−𝑚^𝑢)]²).
For the above 𝜀>0, according to Theorems 2.3 and 2.4, there exists a 𝑘^∗∈𝑁, such that 𝑚₁−𝜀≤𝑥(𝑛),𝑥(𝑛)≤𝑀₁+𝜀,𝑚₂−𝜀≤𝑦(𝑛),𝑦(𝑛)≤𝑀₂+𝜀, (3.8) for all 𝑛≥𝑘^∗.
Noticing that 𝜃_𝑖(𝑛)∈(0,1) (𝑖=1,2,3,4) implies that 𝑥(𝑛)exp(𝜃₁(𝑛)𝑢(𝑛)),  𝑥(𝑛)exp(𝜃₄(𝑛)𝑢(𝑛)) lie between 𝑥(𝑛) and 𝑥(𝑛), 𝑦(𝑛)exp(𝜃₂(𝑛)𝑣(𝑛)),  𝑦(𝑛)exp(𝜃₃(𝑛)𝑣(𝑛)) lie between 𝑦(𝑛) and 𝑦(𝑛). From (3.5), it follows that |𝑢(𝑛+1)|≤|𝑢(𝑛)|max1−𝑀₁+𝜀𝑝^𝜀₁,1−𝑚₁−𝜀𝑝^𝜀₂+|𝑣(𝑛𝛽)|^𝑢1−𝑚_𝑙𝑀₂+𝜀1+𝛼_𝑙(1−𝑚^𝑢)𝑚₁−𝜀≤𝜆^𝜀₁𝑢(𝑛,|𝑣(𝑛,|𝑣(≤|𝑣(𝑛max{|)|)|}𝑛+1)|)|max1−𝑑^𝑢𝑀₂+𝜀,1−𝑑_𝑙𝑚₂+|𝑢(𝑛𝑐−𝜀)|^𝑢𝛽^𝑢1−𝑚_𝑙𝑀₁+𝜀1+𝛼_𝑙(1−𝑚^𝑢)𝑚₁−𝜀²≤𝜆^𝜀₂max{|𝑢(𝑛)|,|𝑣(𝑛)|}. (3.9) Let 𝜆=max{𝜆^𝜀₁,𝜆^𝜀₂}, then 𝜆<1. In view of (3.9), we have max{|𝑢(𝑛+1)|,|𝑣(𝑛+1)|}≤𝜆max{|𝑢(𝑛)|,|𝑣(𝑛)|}≤𝜆^{𝑛−𝑘∗}max𝑢𝑘^∗,𝑣𝑘^∗. (3.10) Therefore (3.6) holds and the proof is complete.

4. Extinction of Predator Species and Stability of Prey Species

In this section, by developing the analysis technique of [16], we show that under some suitable assumptions, the predator will be driven to extinction while prey will be globally attractive to a certain solution of a logistic equation.

We consider a discrete logistic equation: 𝑥(𝑛+1)=𝑥(𝑛)exp(𝑎(𝑛)−𝑏(𝑛)𝑥(𝑛)),𝑛∈𝑁. (4.1) For the above equation, we have the following lemma.

Lemma 4.1 (see [17]). For any positive solution 𝑥^∗(𝑛) of (4.1), one has 𝑚≤liminf_𝑛→∞𝑥^∗(𝑛)≤limsup_𝑛→∞𝑥^∗(𝑛)≤𝑀₁, (4.2) where 𝑚=(𝑎_𝑙/𝑏^𝑢)exp(𝑎_𝑙−𝑏^𝑢𝑀₁) and 𝑀₁ is defined by Theorem 2.3.

Theorem 4.2. Assume that the inequality 𝑚_𝑙𝛾>1−^𝑙𝑐^𝑢𝛽^𝑢𝑀₁ 𝐻₂ holds. Let (𝑥(𝑛),𝑦(𝑛)) be any positive solution of system (1.3), then 𝑦(𝑛)→0 as 𝑛→+∞.

Proof. 𝐻₂ is equivalent to the following inequality: −𝛾_𝑙+𝑐^𝑢𝛽^𝑢1−𝑚_𝑙𝑀₁<0. (4.3) From (4.3) we can choose positive constant 𝜀>0 small enough such that inequality −𝛾_𝑙+𝑐^𝑢𝛽^𝑢1−𝑚_𝑙𝑀₁+𝜀<0 (4.4) holds. Thus, there exists a 𝛿>0, −𝛾_𝑙+𝑐^𝑢𝛽^𝑢1−𝑚_𝑙𝑀₁+𝜀<−𝛿<0. (4.5) Let (𝑥(𝑛),𝑦(𝑛)) be any positive solution of system (1.3). For any 𝑞∈𝑁, according to the equation of system (1.3), we obtain 𝑦ln𝑞+1𝑦𝑞𝑞𝑞𝑦𝑞+𝑐𝑞𝛽𝑞𝑞𝑥𝑞=−𝛾−𝑑1−𝑚𝑞𝑞𝑥𝑞𝑞+𝑐𝑞𝛽𝑞𝑞𝑥𝑞1+𝛼1−𝑚≤−𝛾1−𝑚1+𝛼𝑞1−𝑚𝑞𝑥𝑞≤−𝛾_𝑙+𝑐^𝑢𝛽^𝑢1−𝑚𝑞𝑀₁+𝜀≤−𝛾_𝑙+𝑐^𝑢𝛽^𝑢1−𝑚_𝑙𝑀₁+𝜀<−𝛿<0. (4.6) Summating both sides of the above inequations from 0 to 𝑛−1, we obtain 𝑦(𝑛)ln𝑦(0)<−𝛿𝑛, (4.7) then 𝑦(𝑛)<𝑦(0)exp(−𝛿𝑛). (4.8) Theorem 2.3 implies that 𝑥(𝑛) are bounded eventually, which together with the above inequality shows that 𝑦(𝑛)→0, exponentially, as 𝑛→+∞. This completes the proof of Theorem 4.2.

Theorem 4.3. Assume that 𝑎_𝑙>𝛽^𝑢𝑀₂ and 𝐻₂ holds, also 𝑏^𝑢𝑏_𝑙exp𝑎^𝑢−1<2. 𝐻₃ Then for any positive solution (𝑥(𝑛),𝑦(𝑛)) of system (1.3) and any positive solution 𝑥^∗(𝑛) of system (4.1), one has lim_𝑛→∞𝑥(𝑛)−𝑥^∗(𝑛)=0,lim_𝑛→∞𝑦(𝑛)=0. (4.9)

Proof. Since 𝐻₂ holds, it follows from Theorem 4.2 that lim_𝑛→∞𝑦(𝑛)=0. (4.10) To prove lim_𝑛→∞(𝑥(𝑛)−𝑥^∗(𝑛))=0, let 𝑥(𝑛)=𝑥^∗(𝑛)exp(𝑢(𝑛)), (4.11) then from the first equation of system (1.3) and (4.11), 𝑢(𝑛+1)=𝑢(𝑛)−𝑏(𝑛)𝑥^∗(𝑛)(𝑢(𝑛−𝛽(𝑛(𝑛𝑦(𝑛)exp))−1)(1−𝑚))(𝑛(𝑛𝑥(𝑛).1+𝛼)(1−𝑚)) (4.12) Using the mean-value Theorem, one has exp(𝑢(𝑛)−1)=exp(𝜃(𝑛)𝑢(𝑛))𝑢(𝑛),𝜃(𝑛)∈(0,1). (4.13) Then the first equation of system (1.3) is equivalent to 𝑢(𝑛+1)=𝑢(𝑛)1−𝑏(𝑛)𝑥^∗(𝑛)(𝜃(𝑛)𝑢(𝑛−𝛽(𝑛(𝑛𝑦(𝑛)exp)))(1−𝑚))(𝑛(𝑛𝑥(𝑛).1+𝛼)(1−𝑚)) (4.14) To complete the proof, it suffices to show that lim_𝑛→∞𝑢(𝑛)=0. (4.15) We first assume that 𝜆=max1−𝑏^𝑢𝑀₁,1−𝑏_𝑙𝑚₁<1, (4.16) then we can choose positive constant 𝜀>0 small enough such that 𝜆_𝜀=max1−𝑏^𝑢𝑀₁+𝜀,1−𝑏_𝑙𝑚₁−𝜀<1. (4.17) For the above 𝜀, according to Theorems 2.3 and 2.4, Lemma 4.1, and (4.10), there exists an integer 𝑛₂∈𝑁 such that 𝑚₁−𝜀≤𝑥(𝑛)≤𝑀₁+𝜀,𝑚−𝜀≤𝑥^∗(𝑛)≤𝑀₁+𝜀,𝑦(𝑛)≤𝜀,𝑛≥𝑛₂. (4.18) Noting that 𝑚₁≤𝑚, then 𝑚₁−𝜀≤𝑥(𝑛),𝑥^∗(𝑛)≤𝑀₁+𝜀,𝑦(𝑛)≤𝜀,𝑛≥𝑛₂. (4.19) It follows from (4.19) that 𝛽(𝑛)(1−𝑚(𝑛))(𝑛(𝑛𝑥(𝑛)≤𝛽1+𝛼)(1−𝑚))^𝑢1−𝑚_𝑙1+𝛼_𝑙(1−𝑚^𝑢)𝑚₁−𝜀^def=𝑀_𝜀,𝑛≥𝑛₂. (4.20) Noting that 𝜃(𝑛)∈(0,1), it implies that 𝑥^∗(𝑛)exp(𝜃(𝑛)𝑢(𝑛)) lies between 𝑥^∗(𝑛) and 𝑥(𝑛). From (4.14), (4.17)–(4.20), we get |𝑢(𝑛+1)|≤|𝑢(𝑛)|max1−𝑏^𝑢𝑀₁+𝜀,1−𝑏_𝑙𝑚₁+𝛽−𝜀^𝑢1−𝑚_𝑙1+𝛼_𝑙(1−𝑚^𝑢)𝑚₁𝜀−𝜀=𝜆_𝜀|𝑢(𝑛)|+𝑀_𝜀𝜀,𝑛≥𝑛₂. (4.21) This implies that |𝑢(𝑛)|≤𝜆^𝑛−𝑛2_𝜀𝑢𝑛₂+1−𝜆^𝑛−𝑛2_𝜀1−𝜆_𝜀𝑀_𝜀𝜀,𝑛≥𝑛₂. (4.22) Since 𝜆_𝜀<1 and 𝜀 is arbitrary small, we obtain lim_𝑛→∞𝑢(𝑛)=0; it means that (4.15) holds when 𝜆<1.
Note that 1−𝑏^𝑢𝑀₁≤1−𝑏_𝑙𝑚₁<1. (4.23) Thus, 𝜆<1 is equivalent to 1−𝑏^𝑢𝑀₁>−1, (4.24) or 𝑏^𝑢𝑀₁=𝑏^𝑢𝑏_𝑙exp𝑎^𝑢−1<2. (4.25)
Now, we can conclude that (4.15) is satisfied as 𝐻₃ holds, and so lim_𝑛→∞(𝑥(𝑛)−𝑥^∗(𝑛))=0. This completes the proof of Theorem 4.3.

5. Examples and Numeric Simulations

In this section, we will give two examples to show the feasibility of our results.

Example 5.1. Consider the following system: 𝑥()(𝑛)((𝑛𝑥(𝑛)−()𝑦(𝑛)𝑛+1=𝑥exp0.7−0.9+0.1cos))0.11−0.4()𝑥(𝑛),𝑦()(𝑛)(𝑛)+()𝑥(𝑛)1+0.71−0.4𝑛+1=𝑦exp−0.01−0.1𝑦0.7×0.11−0.4()𝑥(𝑛).1+0.71−0.4 (5.1) One could easily see that 𝑎_𝑙−𝛽^𝑢𝑀₂≈0.2975>0, (𝑐_𝑙𝛽_𝑙(1−𝑚^𝑢)𝑚₁/(1+𝛼^𝑢(1−𝑚_𝑙)𝑀₁))−𝛾^𝑢≈0.002>0, then, condition 𝐻₁ is satisfied. According to Theorem 2.3, system (1.3) is permanent. Numerical simulation (see Figure 1) indicates the permanence of system (5.1).

Figure 1 

Dynamic behaviors of system (5.1) with the initial values (0.14,0.3), (1.0,0.24), (0.6,0.5), and (0.4,0.8).

Example 5.2. Consider the following system: 𝑥()(𝑛)((𝑛𝑥(𝑛)−()𝑦(𝑛)𝑛+1=𝑥exp0.7−0.9+0.1cos))0.11−0.85()𝑥(𝑛),𝑦()(𝑛)(𝑛)+()𝑥(𝑛)1+0.71−0.85𝑛+1=𝑦exp−0.01−0.1𝑦0.7×0.11−0.85()𝑥(𝑛).1+0.71−0.85 (5.2) We could easily see that 1−(𝛾_𝑙/𝑐^𝑢𝛽^𝑢𝑀₁)=0.8457<𝑚_𝑙=0.85,𝑎_𝑙−𝛽^𝑢𝑀₂≈0.2975>0, (𝑏^𝑢/𝑏_𝑙)exp(𝑎^𝑢−1)≈0.926<2. Clearly, conditions of Theorems 4.2 and 4.3 are satisfied. And so, lim_𝑛→∞(𝑥(𝑛)−𝑥^∗(𝑛))=0 and lim_𝑛→∞𝑦(𝑛)=0, where 𝑥^∗(𝑛) is any positive solution of system (4.1). Figure 2 shows the dynamic behaviors of system (5.2).

Figure 2 

Dynamic behaviors of system (5.2) with the initial values (0.1,0.3),(0.27,0.48),(0.6,0.18), and (0.8,1.0).

6. Discussion

We proposed a nonautonomous discrete predator-prey system incorporating a prey refuge and Holling type II functional responses. It is well known that prey species makes use of refuges to decrease predation risk and refuge plays an important role on the dynamic behaviors of predator-prey populations. For system (1.3), we showed that the predator and prey will be coexistent in a globally stable state under some suitable conditions. However, in Section 4, we found that if the refuge is enough large, the predator species will be driven to extinction due to the fewer chances of predation. Obviously, increasing the amount of refuge can increase prey densities and lead to population outbreaks; such kind of finding is consistent with the continuous ones as shown by Kar [6]. In [11], Xu and Jia studied the continuous system (1.2). Sufficient conditions which guarantee the persistence and global stability of positive periodic solution of the system are obtained. Comparing the results of [11] with ours, we found that the conditions which guarantee the persistence of continuous system were similar to the discrete. However, for the conditions that guarantee the global stability of system, the discrete system is more complicated than that of the continuous system. Maybe the reason lies in that for the discrete population dynamics, the constructing and computing of Lyapunov function are relatively complicated than the continuous ones. Unlike the work of [11], we argued that it is an important topic to study the extinction of the species; since more and more species are driven to extinction with the development of modern society, this motivated us to study the extinction of the predator species.

At the end of the paper we would like to mention that one of the referees pointed out that “the nonautonomous character of the model is introduced to simulate the time dependent fluctuating properties of the environment. But it seems to me that a more realistic description of this dependence should be done in terms of stochastic variables rather than deterministic ones.” Indeed, recently, many excellent works concerned with the continuous population model with stochastic variables had been done, see [19–21] and the references cited therein. However, to the best of the authors’ knowledge, to this day, for discrete population dynamics, no similar work has been done. We leave this problem to future research.

Acknowledgments

The authors are grateful to anonymous referees for their excellent suggestions, which greatly improve the presentation of the paper. Also, research was supported by the Natural Science Foundation of Fujian Province (2011J01007), the Foundation of Fujian Education Bureau (JB09001), and the Technology Innovation Platform Project of Fujian Province (2009J1007).