# Spatial Dynamics of a Leslie–Gower Type Predator-Prey Model with Interval Parameters

Authors: Caiyun Wang; Min Guo; Wangsen Lan; Xiaoxin Xu
Journal: Discrete Dynamics in Nature and Society
 (2022)

---

## Abstract

Due to various imprecisions in nature, imprecise parameters in biological modeling should be taken into account. This paper studies the spatial dynamics of an imprecise prey-predator model of the Leslie–Gower type by presenting imprecise parameters as interval parameters. First, conditions of Turing instability are obtained via bifurcation analysis and interval-valued functions. Then, the effects of interval parameters on pattern selection are discussed via multiple-scale analysis. We discover that when all the parameters of the model are interval parameters, the value of the controlled parameter increases, and the range of the pattern selection domain expands as the value of the interval variable increases, i.e., both the controlled parameter and boundary of the pattern selection domain are interval numbers. Finally, under the effects of the interval parameters of diffusion and the prey’s conversion rate into biomass for the predator, the density of the prey decreases or increases, respectively, and the structure or the microstructure of the pattern of the model changes with the growing value of the interval variable. This paper provides a new perspective on the study of the spatial predator-prey model.


---

## Body

## 1. Introduction

The increasing interest in mathematical modeling of ecosystems has heightened the need for understanding the dynamic behaviors of species. Of particular interest are the spatial dynamics of the prey-predator model, one of which is the following mite model introduced by Wollkind et al. [1, 2], whose work was based on the following [3]:(1)dUdT=RU1−UK−PUV,dVdT=SV1−VHU..The predator’s and prey’s population densities are denoted byVT and UT, where R,S,K,H are all positive parameters. R is referred to as the intrinsic growth rate of the prey without predation. S is the intrinsic growth rate of the predator. The feature of this model is that the prey’s growth rate is logistic growth with capacity K. Yet, the predator’s growth rate is interpreted as “logistic growth” with “carrying capacity” HU, which is called the Leslie–Gower term [4, 5]. HU estimates the decrease in predator population density caused by the lack (per capita U/V) of its favorite prey, where H refers to the conversion rate into the predator’s biomass. Although the predator can convert to the substitute of its preferred resource, its population density will suffer due to the availability of its favorite food U. PU, termed as the functional response with the meaning of prey consumption rate of per predator at prey density U and predator density V, can be chosen as Holling type I–III functions [6–9]. Particularly, Holling type III function is more suitable for vertebrate species, i.e., PU=MU2/AU2+1, where M and A are referred to as half-saturation constant and capture rate. The process that the predators are able to change their predation rate and convert from less common to more common prey is the learning process. The learning process can make the foraging system more efficient because the predators can switch to alternative food resources when they know that their population fluctuates. Holling III functional response is for this situation [10].On the other hand, diffusion is one of the spatial motions in the prey-predator model, which depicts the random mobility of individual species. By incorporating spatial motion, model (1) becomes the spatial predator-prey model of the Leslie–Gower type. Although there are many papers about the spatial prey-predator model of the Leslie–Gower type, such as bifurcation analysis or global stability of steady states, and so on, we mainly discuss spatial patterns [11–15]. There have been plenty of studies about spatial patterns of this model [16–22]. For the spatial prey-predator model of the Leslie–Gower type, some researchers focused on the Holling–Tanner type, where PU is chosen as Holling type II response function. They mathematically analyzed the Turing instability of the model and numerically simulated different kinds of patterns [23, 24]. Some papers discussed specific effects on pattern formation of the model, such as cannibalism, delay effects, and Allee effect [25–27]. Other papers discussed the effects of the cross-diffusion [28–30].Yet, diffusion rate in above papers mentioned is constant. In fact, dispersal or diffusion behavior among species is not simple movement, but is a complex process which can be normally divided into three stages: emigration, movement between areas, and immigration. Moreover, dispersal rates in the process of three stages are condition-dependent and are density-dependent in particular [31]. The density-dependent dispersal rate can be either positive or negative. The common species with positive or negative density-dependent dispersal rate are insects, birds or mammals in the real word [32, 33]. Then, based on the work above, researchers have studied spatial dynamics of diffusive system with varying diffusion rate [34, 35]. However, coefficients in above papers are accurate, and all the spatial predator-prey models mentioned above fail to view the biological parameters as imprecise parameters. In reality, the impreciseness can happen because of incomplete information during the determination of the initial value, measurement, and data collection. After Zadeh’s pioneering work of introducing uncertainty into mathematically modeling [36], more and more researchers have incorporated impreciseness in the mathematical models via fuzzy sets [37–39]. As far as predator-prey models are concerned, although there are researches about the imprecise ODE model or imprecise diffusive model [40–44], few papers are about patterns of spatial model with interval diffusion rate. Motivated by the work [34, 35, 44], we will study Turing patterns of an imprecise prey-predator model of the Leslie–Gower type with Holling type III function in this paper. For simplicity, we set up the corresponding accurate model first:(2)∂U∂T=RU1−UK−MU2VAU2+1+D1∇2U,z,T∈Ω×0,∞,∂V∂T=SV1−VHU+D2∇2V,z,T∈Ω×0,∞,∂Uz,T∂n=∂Vz,T∂n=0,z,T∈∂Ω×0,∞,Uz,0=U0>0,z∈Ω,Vz,0=V0>0,z∈Ω,where Ω, a bounded domain, satisfies Ω⊂ℝ2, and n is the outward unit normal vector of the smooth bound of ∂Ω. The coefficients of diffusion of prey and predator are D1,D2.By introducing dimensionless variablesx=U/K,y=MKV/R, and t=RT, system (2) can be simplified as(3)∂x∂t=x1−x−x2yex2+1+d1∇2x,z,t∈Ω×0,∞,∂y∂t=yη1−γyx+d2∇2y,z,t∈Ω×0,∞,∂xz,t∂n=∂yz,t∂n=0,z,t∈∂Ω×0,∞,xz,0=u0>0,z∈Ω,yz,0=v0>0,z∈Ω,where η=S/R,e=AK2,γ=R/HMK2,d1=D1/R, and d2=D2/R.In Section2, interval number and interval-valued function are introduced, and the imprecise model corresponding to model (3) is set up. Turing instability and amplitude equations are analyzed mathematically in Section 3. Theoretical analysis is illustrated by numerical simulations in Section 4. At last, conclusions about the model are drawn in Section 5.


## 2. Prerequisites and Modeling

In this section, we introduce interval numbers, interval-valued function, and arithmetic operations about it. Then, we set up the imprecise predator-prey model of the Leslie–Gower type.


### 2.1. Definitions

Definition 1.
(see Reference [40]). An interval numberAis presented by a closed intervalaˇ,a^. LetRbe the set of real numbers and defineA=aˇ,a^=y|aˇ≤y≤a^,y∈R, wherea^,aˇare the upper and lower limit ofA. Each real numbera∈Rcan be denoted as interval numbera,a.
LetA=aˇ,a^ and B=bˇ,b^ be any two interval numbers, which define the following arithmetic operations:(i)
Addition:A+B=aˇ,a^+bˇ,b^=aˇ+bˇ,a^+b^(ii)
Subtraction:A−B=aˇ,a^−bˇ,b^=aˇ−bˇ,a^−b^(iii)
Scalar multiplication:αA=αaˇ,a^=αaˇ,αa^ifα>0αa^,αaˇifα<0(iv)
Multiplication:A⋅B=aˇ,a^⋅bˇ,b^=minaˇb^,a^bˇ,aˇb^,a^b^,maxaˇbˇ,a^bˇ,aˇb^,a^b^(v)
Division:A/B=aˇ,a^/bˇ,b^=aˇ,a^⋅1/b^,1/b⌣.Definition 2.
(see [40]). An interval-valued function can be defined ashq=a1−qbqfor any interval numbera,b, whereq∈0,1called the interval variable.


### 2.2. Imprecise Leslie–Gower Type Predator-Prey Model with Diffusion

Assume that parameters in system (3) are imprecise and presented as interval numbers. Let η˜,e˜,γ˜,d1˜, and d2˜ be the interval parameters corresponding to η,e,γ,d1, and d2. We have the imprecise predator-prey model of the Leslie–Gower type with diffusion.(4)∂x∂t=x1−x−x2ye˜x2+1+d˜1∇2x,z,t∈Ω×0,∞,∂y∂t=yη˜1−γ˜yx+d˜2∇2y,z,t∈Ω×0,∞,∂xz,t∂n=∂yz,t∂n=0,z,t∈∂Ω×0,∞xz,0=u0>0,z∈Ω,yz,0=v0>0,z∈Ω,where η˜=ηˇ,η^, e˜=eˇ,e^, γ˜=γˇ,γ^, d1˜=d1˜,d1^, and d2˜=d2ˇ,d2^. ηˇ>0,eˇ>0,γˇ>0,d1ˇ>0,d2ˇ>0. Apparently, the interval-valued function for the interval h˜=a,b is hq=a1−qbq, which is a continuous and strictly increasing function. Based on Theorem 1 in [40], the imprecise form of model (4) is given by(5)∂x∂t=x1−x−x2yeˇ1−qe^qx2+1+d1ˇ1−qd^1q∇2x,z,t∈Ω×0,∞,∂y∂t=yηˇ1−qη^q1−γˇ1−qγ^qyx+d2ˇ1−qd^2q∇2yz,t∈Ω×0,∞,∂xz,t∂n=∂yz,t∂n=0,z,t∈∂Ω×0,∞,xz,0=x0>0,z∈Ω,yz,0=y0>0,z∈Ω.


## 2.1. Definitions

Definition 1.
(see Reference [40]). An interval numberAis presented by a closed intervalaˇ,a^. LetRbe the set of real numbers and defineA=aˇ,a^=y|aˇ≤y≤a^,y∈R, wherea^,aˇare the upper and lower limit ofA. Each real numbera∈Rcan be denoted as interval numbera,a.
LetA=aˇ,a^ and B=bˇ,b^ be any two interval numbers, which define the following arithmetic operations:(i)
Addition:A+B=aˇ,a^+bˇ,b^=aˇ+bˇ,a^+b^(ii)
Subtraction:A−B=aˇ,a^−bˇ,b^=aˇ−bˇ,a^−b^(iii)
Scalar multiplication:αA=αaˇ,a^=αaˇ,αa^ifα>0αa^,αaˇifα<0(iv)
Multiplication:A⋅B=aˇ,a^⋅bˇ,b^=minaˇb^,a^bˇ,aˇb^,a^b^,maxaˇbˇ,a^bˇ,aˇb^,a^b^(v)
Division:A/B=aˇ,a^/bˇ,b^=aˇ,a^⋅1/b^,1/b⌣.Definition 2.
(see [40]). An interval-valued function can be defined ashq=a1−qbqfor any interval numbera,b, whereq∈0,1called the interval variable.


## 2.2. Imprecise Leslie–Gower Type Predator-Prey Model with Diffusion

Assume that parameters in system (3) are imprecise and presented as interval numbers. Let η˜,e˜,γ˜,d1˜, and d2˜ be the interval parameters corresponding to η,e,γ,d1, and d2. We have the imprecise predator-prey model of the Leslie–Gower type with diffusion.(4)∂x∂t=x1−x−x2ye˜x2+1+d˜1∇2x,z,t∈Ω×0,∞,∂y∂t=yη˜1−γ˜yx+d˜2∇2y,z,t∈Ω×0,∞,∂xz,t∂n=∂yz,t∂n=0,z,t∈∂Ω×0,∞xz,0=u0>0,z∈Ω,yz,0=v0>0,z∈Ω,where η˜=ηˇ,η^, e˜=eˇ,e^, γ˜=γˇ,γ^, d1˜=d1˜,d1^, and d2˜=d2ˇ,d2^. ηˇ>0,eˇ>0,γˇ>0,d1ˇ>0,d2ˇ>0. Apparently, the interval-valued function for the interval h˜=a,b is hq=a1−qbq, which is a continuous and strictly increasing function. Based on Theorem 1 in [40], the imprecise form of model (4) is given by(5)∂x∂t=x1−x−x2yeˇ1−qe^qx2+1+d1ˇ1−qd^1q∇2x,z,t∈Ω×0,∞,∂y∂t=yηˇ1−qη^q1−γˇ1−qγ^qyx+d2ˇ1−qd^2q∇2yz,t∈Ω×0,∞,∂xz,t∂n=∂yz,t∂n=0,z,t∈∂Ω×0,∞,xz,0=x0>0,z∈Ω,yz,0=y0>0,z∈Ω.


## 3. Bifurcation Analysis

In this section, we analyze the existence of positive equilibrium of model (5) first and then find the conditions that ensure the appearance of the Turing pattern.


### 3.1. Analysis of Positive Equilibrium of Model (5)

The corresponding nondiffusion model of (5) is(6)dxdt=x1−x−x2yeˇ1−qe^qx2+1=fx,y,dydt=yηˇ1−qη^q1−γˇ1−qγ^qyx=gx,y.Obviously, model (6) has an equilibrium E0=1,0. We are interested in the interior equilibrium points E∗=x∗,y∗, where y∗=x∗/γˇ1−qγ^q, and x∗ satisfies(7)Hx≔x3+θ2x2+θx+θ0=0,where θ2=1/γˇ1−qγ^qeˇ1−qe^q−1,θ1=1/eˇ1−qe^q and θ^0=−1/eˇ1−qe^q<0. The numbers of positive equilibria of (6) are the same as the numbers of positive real roots of H in I0=0,1. According to [45], (7) has at most three roots and at least one root in the interval I0=0,1, so system (6) has at most three equilibria and at least one equilibrium.SupposeE∗=x∗,y∗ is any positive equilibrium of model (6), at which the Jacobian matrix is calculated by(8)J=∂f∂x∂f∂y∂g∂x∂g∂yx∗,y∗=fxfygxgy=−2eˇ1−qe^qx∗3−1+eˇ1−qe^qx∗2eˇ1−qe^qx∗2+1−x∗2eˇ1−qe^qx∗2+1ηˇ1−qη^qγˇ1−qγ^q−ηˇ1−qη^q.E∗=x∗,y∗ is stable when the following conditions hold:(9)trJ<0.And(10)detJ>0,where trJanddetJ are trace and determinant of matrix J given by(11)trJ=fu+gv=−2eˇ1−qe^qx∗3+eˇ1−qe^qx∗2−1eˇ1−qe^qx∗2+1−ηˇ1−qη^q,detJ=fxgy−fygx=ηˇ1−qη^qx∗eˇ1−qe^qx∗2+1H′x∗.


### 3.2. Analysis of Turing Instability of Model (5)

We obtain the conditions for Turing instability in this subsection. SupposeE∗=x∗,y∗ is one of the positive and steady equilibrium, whose characteristic polynomial is(12)λI−Jk=0,where the expression of Jk is(13)Jk=fx−d1ˇ1−qd1^qk2fygxgy−d2ˇ1−qd2^qk2,where k is the wavenumber. The eigenvalue of (13) can be obtained by the following equation:(14)λ2−trkλ+Δk=0,where(15)trk=fx+gy−d1ˇ1−qd1^q+d2ˇ1−qd2^qk2,Δk=a1k4−a2k2+fxgy−fygx,where(16)a1=d1ˇ1−qd1^qd2ˇ1−qd2^q,a2=d2ˇ1−qd2^qfx+d1ˇ1−qd1^qgy.Then, the eigenvalue can be expressed as(17)λk=trk±trk2−4Δk2.Reλk>0 guarantees Turing instability. Clearly trk>0, the necessary condition for Turing instability is Δk<0 for some value of k. That is,(18)a2>0.And(19)a22−4a1fxgy−fygx>0.As a result, if (9), (10), (18), and (19) hold, then Turing patterns of model (5) can appear.


### 3.3. Pattern Selection of Model (5)

In order to study the effects of interval parameters on pattern selection, we choosed1 as controlled number to deduce amplitude equations of (5) via multiple-scale analysis.Aroundd1=d1T, model (5) can be rewritten as(20)∂X∂t=LX+N,where(21)X=xy=∑j=13AjxAjyeikj⋅z+c.c…j=1,2,3,where Ajx and Ajy are the amplitudes of the pattern of x and y, c.c. denotes the complex conjugate item of ∑j=13AjxAjyeikj⋅z, and(22)L=fx+d1ˇ1−qd1^q∇2fygxgy+d2ˇ1−qd2^q∇2.And(23)N=16∂3f∂x3x3+12∂3f∂x2∂yx2y+∂2f∂xyxy+12∂2f∂x2x216∂3g∂x3x3+12∂3g∂x2∂yx2y+12∂3g∂y2∂xxy2+∂2g∂xyxy+12∂2g∂x2x2+12∂2g∂y2y2.The detailed expression of derivatives of all orders off,g to x,y can be seen (A.1)–(A.10) in Appendix.Near Turing bifurcation threshold,d1T, the bifurcation parameter d1ˇ1−qd1^q, the variable X, time t, and the nonlinear term N can be expanded as(24)d1T−d1ˇ1−qd1^q=εd11+ε2d22+…,t=T0+εT1+ε2T2,X=xy=εx1y1+ε2x2y2+ε3x3y3+oε4,N=ε2h2+ε3h3+oε4,where(25)h2=12∂2f∂x2x12+∂2f∂x∂yx1y112∂2g∂x2x12+∂2g∂x∂yx1y1+12∂2g∂y2y12,h3=16∂3f∂x3x13+12∂3f∂x2∂yx12x1+∂2f∂x2x2x1+∂2f∂x∂yx2y1+x1y216∂3g∂x3x13+12∂3g∂x2∂yx12y1+12∂3g∂y2∂xx1y12+∂2g∂x2x2x1+∂2g∂y2y2y1+∂2g∂x∂yx2y1+x1y2.The linear operatorL is dissembled via Taylor expansion as follows:(26)L=LT+d1T−d1ˇ1−qd1^qM,where(27)LT=fx+d1T∇2fygxgy+d2ˇ1−qd2^q∇2.And(28)M=−∇2000.The derivation of the amplitude is(29)∂A∂t=ε∂A∂T1+ε2∂A∂T2+⋯,where A=AjxAjy.Substitute (24)–(29) into (20) and collect ε, ε2, and ε3. Following linear systems can be given:(30)ε:LTx1y1=0,(31)ε2:LTx2y2=∂∂T1x1y1−d11Mx1y1−h2≜FxFy,(32)ε3:LTx3y3=∂∂T1x2y2+∂∂T2x1y1−d11Mx2y2−d22Mx1y1−h3≜HxHy,where(33)FxFy=∑j=13FxjFyjeikj⋅z+c.c.And(34)HxHy=∑j=13HxjHyjeikj⋅z+c.c.Solving (30), the solution is as follows:(35)x1y1=l1W1eik1⋅r+W2eik2⋅z+W3eik3⋅z+c.c.where l=fxd2ˇ1−qd2^2q−gyd1T/2gxd1T,kj=kT, c.c denotes the complex conjugate item of W1eik1⋅z+W2eik2⋅z+W3eik3⋅z and Wj is the amplitude of the pattern with the mode expikj⋅z under the first order perturbation, the form of which is determined by higher order term. Substituting equation (35) to (31), one can obtain the explicit expression of Fxj, and Fyjj=1,2,3 in (A.11) Appendix.From the Fredholm solvability condition [13], we know that if the vector function of the right side of (31) is orthogonal to the eigenvectors of the zero eigenvalue of the adjoint operator of LT, which can be denoted by LT+, there exists a nontrivial solution of (31). Substituting (35) into right side of (31), the eigenvector of the operator LT+ is given as follows:(36)1l2e−ikj⋅r+c.c.j=1,2,3,l2=2fyd1Tfxd2ˇ1−qd2^q−gyd1T.In the light of the Fredholm solvability condition, the orthogonality condition is(37)1,l2FxjFyj=0,j=1,2,3.Then, (37) yields(38)l+l2∂W1∂T1=d11kT2lW1+∂2f∂x2l2+2∂2f∂x∂yl+l2∂2g∂x2l2+2∂2g∂x∂yl+∂2g∂y2W¯2W¯3,l+l2∂W2∂T1=d11kT2lW2+∂2f∂x2l2+2∂2f∂x∂yl+l2∂2g∂x2l2+2∂2g∂x∂yl+∂2g∂y2W¯1W¯3,l+l2∂W3∂T1=d11kT2lW3+∂2f∂x2l2+2∂2f∂x∂yl+l2∂2g∂x2l2+2∂2g∂x∂yl+∂2g∂y2W¯1W¯2,where Wj¯ is the conjugation of Wjj=1,2,3.Let the form of the solution of (31) be(39)x2y2=X0Y0+∑j=13XjYjeikj⋅r+∑j=13XjjYijei2kj⋅r+X12Y12eik1−k2⋅r+X23Y23eik2−k3⋅r+X31Y31eik3−k1⋅r+c.c.where expressions of X0,Xj,Xjj,X12,X23Y0,Yj,Yjj,Y12,Y23, and X31Y31 are given as (A.12)–(A.18) in Appendix.Substituting (39) into (32) and collecting the coefficients of exp0,expikj⋅z,expi2kj⋅z, and expikj−ki⋅z, the expressions of Hxj, and Hyjj=1,2,3 can be given in (A.19) in Appendix.Employing the Fredholm solvability condition in (32) again, we have(40)l+l2∂Y1∂T1+∂W1∂T2=lkT2d11Y1+d22W1+EY¯2W¯3+Y¯3W¯2+D1W12+D2W22+W32W1,l+l2∂Y2∂T1+∂W2∂T2=lkT2d11Y2+d22W2+EY¯1W¯3+Y¯3W¯1+D1W22+D2W12+W32W2,l+l2∂Y3∂T1+∂W3∂T2=lkT2d11Y3+d22W3+EY¯1W¯2+Y¯2W¯1+D1W32+D2W12+W22W3,where Yj¯ is the conjugation of Yjj=1,2,3.Multiplying (38) and (40) by ε and ε2, in consideration of (29) and Aj=Ajx=lAiyj=1,2,3, the following amplitude equations are given:(41)τ0∂A1∂t=μA1+hA¯2A¯3−g1A12+g2A22+A32A1,τ0∂A2∂t=μA2+hA¯1A¯3−g1A22+g2A12+A32A2,τ0∂A3∂t=μA3+hA¯1A¯2−g1A32+g2A12+A22A3,where Aj¯ is the conjugation of Ajj=1,2,3, and(42)τ0=l+l2d1TlkT2,μ=d1T−d1ˇ1−qd1^qd1T,h=Ed1Tl2kT2,g1=−D1d1Tl3kT2,g2=−D2d1Tl3kT2,where the expression that E,D1, and D2 stand for can be seen (A.20)–(A.24) in Appendix.According to [13], we get the results of pattern selection in Table 1.Table 1
The domain of pattern selection.
μTypes of the patternStabilityμ2<μ<μ3Spot patternUnstableμ3<μ<μ4Mixture of spots and stripes patternStableμ>μ4Stripe patternStable


## 3.1. Analysis of Positive Equilibrium of Model (5)

The corresponding nondiffusion model of (5) is(6)dxdt=x1−x−x2yeˇ1−qe^qx2+1=fx,y,dydt=yηˇ1−qη^q1−γˇ1−qγ^qyx=gx,y.Obviously, model (6) has an equilibrium E0=1,0. We are interested in the interior equilibrium points E∗=x∗,y∗, where y∗=x∗/γˇ1−qγ^q, and x∗ satisfies(7)Hx≔x3+θ2x2+θx+θ0=0,where θ2=1/γˇ1−qγ^qeˇ1−qe^q−1,θ1=1/eˇ1−qe^q and θ^0=−1/eˇ1−qe^q<0. The numbers of positive equilibria of (6) are the same as the numbers of positive real roots of H in I0=0,1. According to [45], (7) has at most three roots and at least one root in the interval I0=0,1, so system (6) has at most three equilibria and at least one equilibrium.SupposeE∗=x∗,y∗ is any positive equilibrium of model (6), at which the Jacobian matrix is calculated by(8)J=∂f∂x∂f∂y∂g∂x∂g∂yx∗,y∗=fxfygxgy=−2eˇ1−qe^qx∗3−1+eˇ1−qe^qx∗2eˇ1−qe^qx∗2+1−x∗2eˇ1−qe^qx∗2+1ηˇ1−qη^qγˇ1−qγ^q−ηˇ1−qη^q.E∗=x∗,y∗ is stable when the following conditions hold:(9)trJ<0.And(10)detJ>0,where trJanddetJ are trace and determinant of matrix J given by(11)trJ=fu+gv=−2eˇ1−qe^qx∗3+eˇ1−qe^qx∗2−1eˇ1−qe^qx∗2+1−ηˇ1−qη^q,detJ=fxgy−fygx=ηˇ1−qη^qx∗eˇ1−qe^qx∗2+1H′x∗.


## 3.2. Analysis of Turing Instability of Model (5)

We obtain the conditions for Turing instability in this subsection. SupposeE∗=x∗,y∗ is one of the positive and steady equilibrium, whose characteristic polynomial is(12)λI−Jk=0,where the expression of Jk is(13)Jk=fx−d1ˇ1−qd1^qk2fygxgy−d2ˇ1−qd2^qk2,where k is the wavenumber. The eigenvalue of (13) can be obtained by the following equation:(14)λ2−trkλ+Δk=0,where(15)trk=fx+gy−d1ˇ1−qd1^q+d2ˇ1−qd2^qk2,Δk=a1k4−a2k2+fxgy−fygx,where(16)a1=d1ˇ1−qd1^qd2ˇ1−qd2^q,a2=d2ˇ1−qd2^qfx+d1ˇ1−qd1^qgy.Then, the eigenvalue can be expressed as(17)λk=trk±trk2−4Δk2.Reλk>0 guarantees Turing instability. Clearly trk>0, the necessary condition for Turing instability is Δk<0 for some value of k. That is,(18)a2>0.And(19)a22−4a1fxgy−fygx>0.As a result, if (9), (10), (18), and (19) hold, then Turing patterns of model (5) can appear.


## 3.3. Pattern Selection of Model (5)

In order to study the effects of interval parameters on pattern selection, we choosed1 as controlled number to deduce amplitude equations of (5) via multiple-scale analysis.Aroundd1=d1T, model (5) can be rewritten as(20)∂X∂t=LX+N,where(21)X=xy=∑j=13AjxAjyeikj⋅z+c.c…j=1,2,3,where Ajx and Ajy are the amplitudes of the pattern of x and y, c.c. denotes the complex conjugate item of ∑j=13AjxAjyeikj⋅z, and(22)L=fx+d1ˇ1−qd1^q∇2fygxgy+d2ˇ1−qd2^q∇2.And(23)N=16∂3f∂x3x3+12∂3f∂x2∂yx2y+∂2f∂xyxy+12∂2f∂x2x216∂3g∂x3x3+12∂3g∂x2∂yx2y+12∂3g∂y2∂xxy2+∂2g∂xyxy+12∂2g∂x2x2+12∂2g∂y2y2.The detailed expression of derivatives of all orders off,g to x,y can be seen (A.1)–(A.10) in Appendix.Near Turing bifurcation threshold,d1T, the bifurcation parameter d1ˇ1−qd1^q, the variable X, time t, and the nonlinear term N can be expanded as(24)d1T−d1ˇ1−qd1^q=εd11+ε2d22+…,t=T0+εT1+ε2T2,X=xy=εx1y1+ε2x2y2+ε3x3y3+oε4,N=ε2h2+ε3h3+oε4,where(25)h2=12∂2f∂x2x12+∂2f∂x∂yx1y112∂2g∂x2x12+∂2g∂x∂yx1y1+12∂2g∂y2y12,h3=16∂3f∂x3x13+12∂3f∂x2∂yx12x1+∂2f∂x2x2x1+∂2f∂x∂yx2y1+x1y216∂3g∂x3x13+12∂3g∂x2∂yx12y1+12∂3g∂y2∂xx1y12+∂2g∂x2x2x1+∂2g∂y2y2y1+∂2g∂x∂yx2y1+x1y2.The linear operatorL is dissembled via Taylor expansion as follows:(26)L=LT+d1T−d1ˇ1−qd1^qM,where(27)LT=fx+d1T∇2fygxgy+d2ˇ1−qd2^q∇2.And(28)M=−∇2000.The derivation of the amplitude is(29)∂A∂t=ε∂A∂T1+ε2∂A∂T2+⋯,where A=AjxAjy.Substitute (24)–(29) into (20) and collect ε, ε2, and ε3. Following linear systems can be given:(30)ε:LTx1y1=0,(31)ε2:LTx2y2=∂∂T1x1y1−d11Mx1y1−h2≜FxFy,(32)ε3:LTx3y3=∂∂T1x2y2+∂∂T2x1y1−d11Mx2y2−d22Mx1y1−h3≜HxHy,where(33)FxFy=∑j=13FxjFyjeikj⋅z+c.c.And(34)HxHy=∑j=13HxjHyjeikj⋅z+c.c.Solving (30), the solution is as follows:(35)x1y1=l1W1eik1⋅r+W2eik2⋅z+W3eik3⋅z+c.c.where l=fxd2ˇ1−qd2^2q−gyd1T/2gxd1T,kj=kT, c.c denotes the complex conjugate item of W1eik1⋅z+W2eik2⋅z+W3eik3⋅z and Wj is the amplitude of the pattern with the mode expikj⋅z under the first order perturbation, the form of which is determined by higher order term. Substituting equation (35) to (31), one can obtain the explicit expression of Fxj, and Fyjj=1,2,3 in (A.11) Appendix.From the Fredholm solvability condition [13], we know that if the vector function of the right side of (31) is orthogonal to the eigenvectors of the zero eigenvalue of the adjoint operator of LT, which can be denoted by LT+, there exists a nontrivial solution of (31). Substituting (35) into right side of (31), the eigenvector of the operator LT+ is given as follows:(36)1l2e−ikj⋅r+c.c.j=1,2,3,l2=2fyd1Tfxd2ˇ1−qd2^q−gyd1T.In the light of the Fredholm solvability condition, the orthogonality condition is(37)1,l2FxjFyj=0,j=1,2,3.Then, (37) yields(38)l+l2∂W1∂T1=d11kT2lW1+∂2f∂x2l2+2∂2f∂x∂yl+l2∂2g∂x2l2+2∂2g∂x∂yl+∂2g∂y2W¯2W¯3,l+l2∂W2∂T1=d11kT2lW2+∂2f∂x2l2+2∂2f∂x∂yl+l2∂2g∂x2l2+2∂2g∂x∂yl+∂2g∂y2W¯1W¯3,l+l2∂W3∂T1=d11kT2lW3+∂2f∂x2l2+2∂2f∂x∂yl+l2∂2g∂x2l2+2∂2g∂x∂yl+∂2g∂y2W¯1W¯2,where Wj¯ is the conjugation of Wjj=1,2,3.Let the form of the solution of (31) be(39)x2y2=X0Y0+∑j=13XjYjeikj⋅r+∑j=13XjjYijei2kj⋅r+X12Y12eik1−k2⋅r+X23Y23eik2−k3⋅r+X31Y31eik3−k1⋅r+c.c.where expressions of X0,Xj,Xjj,X12,X23Y0,Yj,Yjj,Y12,Y23, and X31Y31 are given as (A.12)–(A.18) in Appendix.Substituting (39) into (32) and collecting the coefficients of exp0,expikj⋅z,expi2kj⋅z, and expikj−ki⋅z, the expressions of Hxj, and Hyjj=1,2,3 can be given in (A.19) in Appendix.Employing the Fredholm solvability condition in (32) again, we have(40)l+l2∂Y1∂T1+∂W1∂T2=lkT2d11Y1+d22W1+EY¯2W¯3+Y¯3W¯2+D1W12+D2W22+W32W1,l+l2∂Y2∂T1+∂W2∂T2=lkT2d11Y2+d22W2+EY¯1W¯3+Y¯3W¯1+D1W22+D2W12+W32W2,l+l2∂Y3∂T1+∂W3∂T2=lkT2d11Y3+d22W3+EY¯1W¯2+Y¯2W¯1+D1W32+D2W12+W22W3,where Yj¯ is the conjugation of Yjj=1,2,3.Multiplying (38) and (40) by ε and ε2, in consideration of (29) and Aj=Ajx=lAiyj=1,2,3, the following amplitude equations are given:(41)τ0∂A1∂t=μA1+hA¯2A¯3−g1A12+g2A22+A32A1,τ0∂A2∂t=μA2+hA¯1A¯3−g1A22+g2A12+A32A2,τ0∂A3∂t=μA3+hA¯1A¯2−g1A32+g2A12+A22A3,where Aj¯ is the conjugation of Ajj=1,2,3, and(42)τ0=l+l2d1TlkT2,μ=d1T−d1ˇ1−qd1^qd1T,h=Ed1Tl2kT2,g1=−D1d1Tl3kT2,g2=−D2d1Tl3kT2,where the expression that E,D1, and D2 stand for can be seen (A.20)–(A.24) in Appendix.According to [13], we get the results of pattern selection in Table 1.Table 1
The domain of pattern selection.
μTypes of the patternStabilityμ2<μ<μ3Spot patternUnstableμ3<μ<μ4Mixture of spots and stripes patternStableμ>μ4Stripe patternStable


## 4. Numerical Simulations

In this part, we numerically simulated patterns of system (5) in 2-dimensional and 200×200 space area, employing Euler-forward finite difference method and zero-flux boundary condition. We chose space step size as Δh=1 and time step size as Δt=0.05 with time interval 0,8000. For the sake of similar patterns of both prey and predator, the prey’s pattern was chosen to investigate in detail, and the simulations will not stop until the patterns reach their steady state. Additionally, we varied the interval variable q to show how the interval parameter would affect the spatial dynamics of system (5).


### 4.1. Effects of Interval Parameters on Controlled Parameter and Pattern Selection

In this section, we selected all parameters as interval parameters (Table2) and studied the effects of interval parameters on controlled parameter and pattern selection.Table 2
Parameters set in subsection4.1.
γˇ,γ^eˇ,e^d2ˇ,d2^ηˇ,η^[0.007384231991, 0.007434351518][523, 525][3.3, 4][0.7, 0.75]It can be seen from Table3 that the values of controlled parameter d1T increase from 0.2860078311 to 0.3264191705 as the values of interval variable q increase. According to the definition of interval number, it means that when parameters of system (5) are interval numbers, the controlled parameter is also an interval number with the value [0.2860078311, 0.3264191705]. By the same way, the boundary μ3 and μ4 of the domain of pattern selection (Table 1), where the controlled parameter lies, are also interval numbers. Moreover, the domains, i.e., μ2,μ3 and μ3,μ4, of pattern selection expand as an interval variable q increases. Figures 1–3 show the patterns that appear in corresponding pattern selection domain when q=0.3, q=0.6 and q=0.9. All the patterns in Figures 1–3 are in agreement with the theoretical analysis in Table 1, i.e., spot patterns (Figures 1(a), 2(a), and 3(a)) appear when μ lies in μ2,μ3, and mixture of spots and stripes (Figures 1(b), 2(b), and 3(b)) emerge when μ lies in μ2,μ3, and stripe patterns ((Figures 1(c), 2(c), and 2(c)) arise when μ>μ4. From this, we conclude that interval parameter of system (5) makes the boundary μ3 and μ4 be an interval number or changes the range of pattern selection domain μ2,μ3 and μ3,μ4, but it does not change the types of patterns in the corresponding pattern selection domain.Table 3
Parameters calculated in pattern selection corresponding to different values of interval variableq.
qd1Tg1g2hμ1μ2μ3μ400.286007831161.14349619130.9295721−1.936013253−0.00290101909900.047057555480.19488176630.30.297864947156.75265547123.6049461−2.283089563−0.00428712187700.066191084920.27654362640.60.309955239952.38574633116.2542829−2.620820927−0.00602741266000.088209203560.27654362640.90.322267718948.04954243108.8918077−2.949632089−0.00818213342300.11293126800.372171993310.326419170546.61205020106.4372849−3.057325591−0.00900551227800.12173424490.4817919226Figure 1
Pattern selection whenq=0.3: (a) μ2<μ=0.05<μ3, (b) μ3<μ=0.26<μ4, and (c) μ4<μ=0.38.
(a)(b)(c)Figure 2
Pattern selection whenq=0.6: (a) μ2<μ=0.08<μ3, (b) μ3<μ=0.28<μ4, and (c) μ4<μ=0.4.
(a)(b)(c)Figure 3
Pattern selection whenq=0.9: (a) μ2<μ=0.1<μ3, (b) μ3<μ=0.22<μ4, and (c) μ4<μ=0.5.
(a)(b)(c)


### 4.2. Effects of Interval Parameterd1 on Patterns of System (5)

The diffusion coefficient is an important parameter with the meaning of random motion, so we figured out what influences the impreciseness of diffusion would exert on system (5) by choosing values of parameters in Table 4, where we selected d1 as interval parameters and others as accurate.Table 4
Parameters set in subsection4.2.
d1ˇ,d1^eˇ,e^d2ˇ,d2^γˇ,γ^ηˇ,η^[0.54, 0.93]5304.10.0019964347240.65It can be seen from Figure4 that the microstructure of spot pattern of the prey changes with the increasing of the interval variable q: the density of the cold spot pattern gradually becomes sparse. Specifically, there are fewer and fewer isolated areas formed by the low density of the prey. To differentiate the patterns in Figure 4 more clearly, we plotted time-series of prey’s density at spatial location (21, 21) in Figure 5 and spatial mean density of the prey with respect to interval variable q in Figure 6. From Figure 5, when t>200, the density of prey lies in (0.4, 0.45) with q=0 in Figure 5(a), (0.35, 0.4) with q=0.4 in Figure 5(b), (0.3, 0.35) with q=0.6 in Figure 5(c), (0.24, 0.26) with q=0.8 in Figure 5(d), and (0.2, 0.21) with q=1 in Figure 5(e). One can see that the density of the prey decreases along with the increasing of q, which is in accordance with results of Figure 6. From above analysis, it means that the interval parameter d1 not only affects the microstructure of pattern of the prey but also the density of the prey.Figure 4
Patterns transitions induced by interval parameterd1. (a) q=0, (b) q=0.4, (c) q=0.6, (d) q=0.8, and (e) q=1.
(a)(b)(c)(d)(e)Figure 5
Time-series of prey’s density at spatial location21,21. (a) q=0, (b) q=0.4, (c) q=0.6, (d) q=0.8, and (e) q=1.
(a)(b)(c)(d)(e)Figure 6
Spatial mean density of the prey with respect to interval variableq.


### 4.3. Effects of Interval Parameterγ on Patterns of System (5)

The parameterγ denotes the reliance of predators on prey, which is an intrinsic factor of predator-prey systems of the Leslie–Gower type, so we studied the effects of the interval parameter γ on model (5). For this purpose, we chose γ as the interval parameter and others as accurate, the values of which were provided in Table 5. Then, we simulated patterns of prey numerically.Table 5
Parameters set in subsection4.3.
d1ˇ,d1^eˇ,e^d2ˇ,d2^γˇ,γ^ηˇ,η^0.243604.8[0.00221,0.004]0.7From Figure7, patterns of the prey transit orderly with the increasing of the interval variable q: hot spot pattern (Figure 7(a)), mixture of hot spot and stripe pattern (Figure 7(b)), stripe pattern (Figure 7(c)), mixture of cold spot and stripe pattern (Figure 7(d)), cold spot pattern (Figure 7(e)) with the increasing of the interval variable of q. It is easy to see the differences of patterns in Figure 7. Compared with Figure 4, Figures 8 and 9 show that the changing of prey’s density show an opposite trend, i.e., the density of the prey increases as the values of the interval variable q grows. We conclude that the interval parameter γ has a remarkable effect on types of pattern and the density of the prey.Figure 7
Patterns transitions induced by interval parameterγ. (a) q=0, (b) q=0.3, (c) q=0.7, (d) q=0.85, and (e) q=1.
(a)(b)(c)(d)(e)Figure 8
Time-series of prey’s density at spatial location61,11: (a) q=0 and (b) q=1.
(a)(b)Figure 9
Mean density of the prey with respect to interval variableq.


## 4.1. Effects of Interval Parameters on Controlled Parameter and Pattern Selection

In this section, we selected all parameters as interval parameters (Table2) and studied the effects of interval parameters on controlled parameter and pattern selection.Table 2
Parameters set in subsection4.1.
γˇ,γ^eˇ,e^d2ˇ,d2^ηˇ,η^[0.007384231991, 0.007434351518][523, 525][3.3, 4][0.7, 0.75]It can be seen from Table3 that the values of controlled parameter d1T increase from 0.2860078311 to 0.3264191705 as the values of interval variable q increase. According to the definition of interval number, it means that when parameters of system (5) are interval numbers, the controlled parameter is also an interval number with the value [0.2860078311, 0.3264191705]. By the same way, the boundary μ3 and μ4 of the domain of pattern selection (Table 1), where the controlled parameter lies, are also interval numbers. Moreover, the domains, i.e., μ2,μ3 and μ3,μ4, of pattern selection expand as an interval variable q increases. Figures 1–3 show the patterns that appear in corresponding pattern selection domain when q=0.3, q=0.6 and q=0.9. All the patterns in Figures 1–3 are in agreement with the theoretical analysis in Table 1, i.e., spot patterns (Figures 1(a), 2(a), and 3(a)) appear when μ lies in μ2,μ3, and mixture of spots and stripes (Figures 1(b), 2(b), and 3(b)) emerge when μ lies in μ2,μ3, and stripe patterns ((Figures 1(c), 2(c), and 2(c)) arise when μ>μ4. From this, we conclude that interval parameter of system (5) makes the boundary μ3 and μ4 be an interval number or changes the range of pattern selection domain μ2,μ3 and μ3,μ4, but it does not change the types of patterns in the corresponding pattern selection domain.Table 3
Parameters calculated in pattern selection corresponding to different values of interval variableq.
qd1Tg1g2hμ1μ2μ3μ400.286007831161.14349619130.9295721−1.936013253−0.00290101909900.047057555480.19488176630.30.297864947156.75265547123.6049461−2.283089563−0.00428712187700.066191084920.27654362640.60.309955239952.38574633116.2542829−2.620820927−0.00602741266000.088209203560.27654362640.90.322267718948.04954243108.8918077−2.949632089−0.00818213342300.11293126800.372171993310.326419170546.61205020106.4372849−3.057325591−0.00900551227800.12173424490.4817919226Figure 1
Pattern selection whenq=0.3: (a) μ2<μ=0.05<μ3, (b) μ3<μ=0.26<μ4, and (c) μ4<μ=0.38.
(a)(b)(c)Figure 2
Pattern selection whenq=0.6: (a) μ2<μ=0.08<μ3, (b) μ3<μ=0.28<μ4, and (c) μ4<μ=0.4.
(a)(b)(c)Figure 3
Pattern selection whenq=0.9: (a) μ2<μ=0.1<μ3, (b) μ3<μ=0.22<μ4, and (c) μ4<μ=0.5.
(a)(b)(c)


## 4.2. Effects of Interval Parameterd1 on Patterns of System (5)

The diffusion coefficient is an important parameter with the meaning of random motion, so we figured out what influences the impreciseness of diffusion would exert on system (5) by choosing values of parameters in Table 4, where we selected d1 as interval parameters and others as accurate.Table 4
Parameters set in subsection4.2.
d1ˇ,d1^eˇ,e^d2ˇ,d2^γˇ,γ^ηˇ,η^[0.54, 0.93]5304.10.0019964347240.65It can be seen from Figure4 that the microstructure of spot pattern of the prey changes with the increasing of the interval variable q: the density of the cold spot pattern gradually becomes sparse. Specifically, there are fewer and fewer isolated areas formed by the low density of the prey. To differentiate the patterns in Figure 4 more clearly, we plotted time-series of prey’s density at spatial location (21, 21) in Figure 5 and spatial mean density of the prey with respect to interval variable q in Figure 6. From Figure 5, when t>200, the density of prey lies in (0.4, 0.45) with q=0 in Figure 5(a), (0.35, 0.4) with q=0.4 in Figure 5(b), (0.3, 0.35) with q=0.6 in Figure 5(c), (0.24, 0.26) with q=0.8 in Figure 5(d), and (0.2, 0.21) with q=1 in Figure 5(e). One can see that the density of the prey decreases along with the increasing of q, which is in accordance with results of Figure 6. From above analysis, it means that the interval parameter d1 not only affects the microstructure of pattern of the prey but also the density of the prey.Figure 4
Patterns transitions induced by interval parameterd1. (a) q=0, (b) q=0.4, (c) q=0.6, (d) q=0.8, and (e) q=1.
(a)(b)(c)(d)(e)Figure 5
Time-series of prey’s density at spatial location21,21. (a) q=0, (b) q=0.4, (c) q=0.6, (d) q=0.8, and (e) q=1.
(a)(b)(c)(d)(e)Figure 6
Spatial mean density of the prey with respect to interval variableq.


## 4.3. Effects of Interval Parameterγ on Patterns of System (5)

The parameterγ denotes the reliance of predators on prey, which is an intrinsic factor of predator-prey systems of the Leslie–Gower type, so we studied the effects of the interval parameter γ on model (5). For this purpose, we chose γ as the interval parameter and others as accurate, the values of which were provided in Table 5. Then, we simulated patterns of prey numerically.Table 5
Parameters set in subsection4.3.
d1ˇ,d1^eˇ,e^d2ˇ,d2^γˇ,γ^ηˇ,η^0.243604.8[0.00221,0.004]0.7From Figure7, patterns of the prey transit orderly with the increasing of the interval variable q: hot spot pattern (Figure 7(a)), mixture of hot spot and stripe pattern (Figure 7(b)), stripe pattern (Figure 7(c)), mixture of cold spot and stripe pattern (Figure 7(d)), cold spot pattern (Figure 7(e)) with the increasing of the interval variable of q. It is easy to see the differences of patterns in Figure 7. Compared with Figure 4, Figures 8 and 9 show that the changing of prey’s density show an opposite trend, i.e., the density of the prey increases as the values of the interval variable q grows. We conclude that the interval parameter γ has a remarkable effect on types of pattern and the density of the prey.Figure 7
Patterns transitions induced by interval parameterγ. (a) q=0, (b) q=0.3, (c) q=0.7, (d) q=0.85, and (e) q=1.
(a)(b)(c)(d)(e)Figure 8
Time-series of prey’s density at spatial location61,11: (a) q=0 and (b) q=1.
(a)(b)Figure 9
Mean density of the prey with respect to interval variableq.


## 5. Discussion and Conclusions

The spatial dynamics of a prey-predator model with an interval parameter were initially studied in this paper, compared with the former work [23–30]. First, conditions of Turing instability of the model with an interval parameter were gained via an interval-valued function and mathematical analysis. Then, the effects of all or some interval parameters on the spatial dynamics of system (5) were discussed. Specifically, via multiple-scale analysis and numerical simulation, pattern selection and pattern transitions for system (5) under interval parameters were obtained. It is well known that multiple-scale analysis is effective when the values of parameters are very close to the controlled parameter. According to our results, whatever the sensitive parameter, the controlled parameter is also an interval number. The results in Section 4.2 and 4.3 show that the effects of different interval parameters on system (5) are different. In this paper, the impreciseness is presented by the interval parameter, so we conclude that the impreciseness of the parameter does affect the behavior of the entire spatial biological system. Based on the conclusion, our hypothesis proposed in the introduction that researchers should incorporate impreciseness when they make mathematical modeling is reasonable.It should be noted that we choosed1 as a controlled number to study the effects of interval parameters on pattern selection. Other parameters can also be chosen as controlled parameters. Moreover, the interval variable q describes any position of the value of an interval-valued function on the closed interval that presents the interval number. According to the definition of an interval-valued function, the bigger value of q means the value of the interval-valued function lies closer to the upper limit of the interval number. Our results reveal the tendency of the spatial dynamics in the process of the values of an interval-valued function increasing to the upper limit of the interval number (see Tables 6 and 7). This study has concentrated on the effects of interval number on the spatial dynamics of a spatial prey-predator model. Despite its preliminary nature, imprecise parameters make the model closer to the reality of a lack of accurate information and provide a new perspective on studying spatial predator-prey models. In future work, it intrigues us to think about the effects of impreciseness on the spatial epidemic model.Table 6
Different values ofd1 corresponding to q.
q00.40.60.81d10.540.67116520790.74825094340.83419025240.93Table 7
Different values ofγ corresponding to q.
q00.30.70.851γ0.00210.0026405390170.0033478013180.0036593995780.004



---

*Source: 2483688-2022-12-12.xml*