A Logistic-Harvest Model with Allee Effect under Multiplicative Noise
aa r X i v : . [ q - b i o . P E ] A ug A Logistic-Harvest Model with Allee E ff ect under Multiplicative Noise ✩ Almaz Tesfay a,b, ∗ , Daniel Tesfay a,b , James Brannan c , Jinqiao Duan d a School of Mathematics and Statistics & Center for Mathematical Sciences, Huazhong University of Science and Technology,Wuhan 430074, China b Department of Mathematics, Mekelle University, P.O.Box 231, Mekelle, Ethiopia c Department of Mathematical Sciences, Clemson University,Clemson, South Carolina 29634, USA d Department of Applied Mathematics, Illinois Institute of Technology, Chicago, IL 60616, USA
Abstract
This work is devoted to the study of a stochastic logistic growth model with and without the Allee e ff ect.Such a model describes the evolution of a population under environmental stochastic fluctuations and is inthe form of a stochastic di ff erential equation driven by multiplicative Gaussian noise. With the help of theassociated Fokker-Planck equation, we analyze the population extinction probability and the probability ofreaching a large population size before reaching a small one. We further study the impact of the harvest rate,noise intensity, and the Allee e ff ect on population evolution. The analysis and numerical experiments showthat if the noise intensity and harvest rate are small, the population grows exponentially, and upon reachingthe carrying capacity, the population size fluctuates around it. In the stochastic logistic-harvest modelwithout the Allee e ff ect, when noise intensity becomes small (or goes to zero), the stationary probabilitydensity becomes more acute and its maximum point approaches one. However, for large noise intensity andharvest rate, the population size fluctuates wildly and does not grow exponentially to the carrying capacity.So as far as biological meanings are concerned, we must catch at small values of noise intensity and harvestrate. Finally, we discuss the biological implications of our results. Keywords:
Stochastic dynamics; logistic growth model; threshold population; Fokker-Planck equation;harvesting factor; stochastic di ff erential equation. -Mathematics Subject Classification: 39A50, 45K05, 65N22.
1. Introduction
A group of individuals of the same species living in a limited place is called a population [24]. Thedynamical process of population growth and decline is a function of factors that are intrinsic to a populationand the environmental conditions.The well known logistic growth model describes the growth of population, followed by a reduction, andbound by the maximum population size (carrying capacity). This model is a nonlinear di ff erential equation dX t dt = rX t (cid:18) − X t K (cid:19) , X (0) = x , (1.1) ∗ Corresponding author
Email addresses: [email protected] (Almaz Tesfay ), [email protected] (Daniel Tesfay), [email protected] (James Brannan), [email protected] (Jinqiao Duan)
Preprint submitted to ... August 5, 2020 here r > X t is the population size at time t and K is the carrying capacity. Thismodel was first introduced by Verhust [12]. When X t is very small, the equation in (1.1) becomes dX t dt = rX t and dX t dt =
0, when X t nears the carrying capacity K .Equation (1.1) has a unique solution given by X t = K + Ae − rt , where A = ( Kx − t → ∞ .Allee e ff ect was studied widely in a biology book [13]. In this book, the authors cited many papersdealing with the Allee e ff ect. Allee [30] suggested that per capita birth rate declines at a low populationsize ( densities). In this case, the population may go to extinction. The logistic growth model with the Alleee ff ect is one of the most important models in mathematical ecology owing to its theoretical and practicalsignificance. An Allee e ff ect shows a non-negative association between reproduction and population size,and survival of individuals. There are two distinct variations of the Allee e ff ect. Namely, strong Allee e ff ectand weak Allee e ff ect. Strong Allee e ff ect introduces a population threshold [25] that the population mustexceed in order to grow, while the weak Allee e ff ect does not admit any threshold. For more details aboutthis model see [34, 27, 25] and the reference therein.The classic general logistic growth model with Allee e ff ect [32] is given by dX t dt = rX t (cid:18) X t S − (cid:19) (cid:18) − X t K (cid:19) , X (0) = x , (1.2)where X t is the population size at time t in a given area or place, r > K > S refers to the threshold population (Allee threshold) which is the minimumpopulation that is necessary for the species to survive with values 0 < S < K . Extinction occurs wheneverthe population decreases below the Allee threshold value S . Here the initial population size X must begreater than the threshold value S . Equation (1.2) has two stable equilibrium solutions at X ( t ) = X ( t ) = K , and an unstable equilibrium solution at X ( t ) = S .Based on the resources available to the system, the population should reach the carrying capacity K .If the initial population is below the critical threshold S , then it approaches extinction as time goes on.Thus the threshold population is useful to biologists in order to determine whether a given species shouldbe placed on the endangered list so that the survival of the species will then be given due attention andnecessary protection.Fishing has a lot of benefits to human beings and it has also a great impact on the socio-economic andinfrastructure development of a country. For example, it serves as food, generates income, and creates jobopportunities. Many scientists [3, 33, 40] devised strategies to prevent the extinction of renewable resourcessuch as fish by harvesting, and they agreed on the importance. See [15, 28 ? ] for further explanation onharvesting strategies.The logistic growth model, with and without the Allee e ff ect, and with harvesting has been used to studythe fishery farming [28]. Harvesting is an interesting research area in a population study. The most importantinput for the successful management of harvested populations is a sustainable strategy. Harvesting strategyshould not lead to instabilities or extinctions.In this paper, we focus on proportional harvesting which removes a fixed proportionality of individualseach time t ( year). In other words, if the population increases, the harvested also increases, and if thepopulation decreases the quantity harvested decreases.Now let us consider the mathematical model of the relative-rate harvesting on logistic growth model[12] in Eq. (1.3) and logistic growth with Allee e ff ect Eq. (1.4), respectively. dX t dt = rX t (cid:18) − X t K (cid:19) − λ X t , X (0) = x , (1.3)2nd dX t dt = rX t (cid:18) X t S − (cid:19) (cid:18) − X t K (cid:19) − λ X t , X (0) = x , (1.4)where again r is the population growth rate, K > S refers to the Alleethreshold, 0 < S < K and λ is harvest rate.Equilibria points of Eq. (1.3) lie at X t = , and X t = K (cid:16) − λ r (cid:17) for λ < r . The potential function V ( x )of Eq. (1.3) is given by V ( x ) = r h − x + K x i + λ x . For λ =
0, the function V ( x ) becomes the potentialfunction of equation (1.1).Model (1.4) has equilibria points at X t = , and at the solutions of λ = r (cid:16) X t S − (cid:17) (cid:16) − X t K (cid:17) . The maximumof the parabola is at X t = ( S + K )2 , where we have a saddle-node bifurcation at λ = r ( K − S ) S K . The potential function V ( x ) of Eq. (1.4) for λ , V ( x ) = r h x − ( S + K )3 S K x + S K x i + λ x . V ( x ) reduces to the potential function of equation (1.2)if λ = ? , 29] considered the deterministic model of logistic growth with andwithout Allee e ff ect under harvesting factor and studied the behavior of the deterministic model free of anystochastic element. Even though deterministic models are much easier to analyze than their correspondingstochastic models, they neglect of random influences on the growth process. stochastic di ff erential equationsmay be regarded as more adequate models for the development of a population. Since random events a ff ectpopulation dynamics.In our paper, we focus on both deterministic and stochastic model. Biological populations exhibit someform of stochastic behavior and that environmental noise should thus be an integral component of any dy-namic population model [25]. Population ecology deals with demographic and environmental stochasticity.In this work, we consider environmental stochasticity.Several factors a ff ect the environment population resides [42]. To model environmental e ff ects, onepossibility is to explicitly include additional variables, for example, chemical agents, food supply, rainfall,and average temperature into di ff erential equation (1.3) and (1.4). On the other hand, population systemsare often subject to environmental noise. Thus it is important to reveal how the noise a ff ects the populationsystems.According to Equation (12.20) in [25], a stochastic fishing model is given by a stochastic di ff erentialequation (SDE) dX t = ( H ( X t ) X t − λ X t ) dt + ǫ X t dB t , X (0) = x , (1.5)where H ( X t ) is natural growth rate of harvested population, and λ, ǫ are constants. The drift coe ffi cient anddi ff usion coe ffi cient of this SDE are f ( X t ) = H ( X t ) X t − λ X t and g ( X t ) = ǫ X t , respectively. This stochasticdi ff erential equation has a unique solution [25], and the solution is a homogenous di ff usion process.Here, we choose H ( X ) = rX (1 − XK ) and H ( X ) = rX ( XS − − XK ). Then by Eq. (1.5) the stochasticversion of the logistic-harvest model of (1.3) and (1.4), respectively are dX t = (cid:20) rX t (cid:18) − X t K (cid:19) − λ X t (cid:21) + ǫ X t dB t , X (0) = x , (1.6)and dX t = (cid:20) rX t (cid:18) X t S − (cid:19) (cid:18) − X t K (cid:19) − λ X t (cid:21) + ǫ X t dB t , X (0) = x , (1.7)where B t is a one-dimensional Brownian motion and ǫ is the Gaussian noise intensity with 0 < ǫ < . The objective of this work is to investigate the behavior of the logistic-harvest without or with Alleee ff ect, driven by multiplicative Gaussian noise. In other words, we will combine the theory of populationbiology with that of stochastic di ff erential equations. According to Drake and Lodge [17], there are three3tatistics most commonly used to evaluate the population helpful in studying stochastic population models.These quantities are the extinction probability, the first passage probability, and the mean time to extinction.In our study, we focus on the extinction probability.In this paper, we first review the deterministic logistic-harvest model with and without the Allee ef-fect, and then we investigate their stochastic counterpart. We further discuss the extinction probability ofthe stochastic models. To gain some insight into the logistic-harvest mechanism and consequently aboutthe underlying biological phenomenon, we apply the Euler-Maruyama scheme to approximate the samplesolution paths of the stochastic logistic-harvesting model. Finally, we present a short discussion on thecomparison between the deterministic models and stochastic models as parameter x , λ and ǫ vary.This paper is arranged as follows: After recalling basic facts about Brownian motion and stochasticdi ff erential equations in section 2, we review and discuss the behavior of the equilibrium solution of thedeterministic of the logistic-harvest model without the Allee e ff ect (1.3) and analyze its correspondingstochastic model (section 3). We drive the exact solution of model (1.6) and explain the e ff ect of the harvestrate λ , noise intensity ǫ and initial value x on the stationary density function of the Fokker-Plank equationfor the SDE in (1.6). In section 4, we review the deterministic logistic-harvest model with Allee e ff ect(1.4). We discuss the e ff ect of the harvest rate λ , noise intensity ǫ , and initial value x on the stationarydensity function of the Fokker-Plank equation for the SDE in (1.6). The Euler-Maruyama approximationis then used to approximate the solution of the stochastic model. In section 5, we summarize numericalexperiments to reveal the sample path behaviors of the deterministic and stochastic models. Finally, insection 6, we present a short conclusion about our findings.
2. Preliminaries
In this section, we recall some basic facts about Brownian motion and a stochastic di ff erential equations.Assume ( Ω , F , { F t } t > , P ) is a complete probability space with a filtration { F t } t > satisfying the usualconditions, i.e. { F t } t > is increasing and continuous while F contains all P − null sets. Brownian motion B t is an abstract of random walk process [20] defined on the filtered probability space ( Ω , F , { F t } t > , P ) whichsatisfies the following properties: • Stationary and normal increments: B t − B s , for s < t is normally distributed with mean is equal tozero and variance is equal to t − s , • Independence of increments: B t − B s , for s < t , is independent of the past, • Continuity of paths: B t is a continuous function of t , almost surely. • The process starts at origin: B = • Brownian motion is nowhere di ff erentiable, almost surely.Stochastic di ff erential equations [6] are often used in modeling biological phenomena, by taking the intrinsicrandom e ff ects into account. Intrinsic forcing induced SDE models are considered in population dynamics,epidemics, genetics, and oncogenesis.Consider a stochastic di ff erential equation driven by Gaussian noise dX t = f ( X t ) dt + g ( X t ) dB t , t ∈ (0 , ∞ ) . (2.1)If both drift f and noise intensity g satisfy a local Lipschitz condition, a growth condition or a priori estimateon the solution, then stochastic di ff erential equation (2.1) has a unique continuous solution X t on t ∈ (0 , ∞ ) . [1, 7, 16, 21]. 4 . Logistic-harvest model without Allee e ff ect Deterministic logistic-harvest model without Allee e ff ect Consider a population X t with dynamic according to the logistic growth model without Allee e ff ect. Theidea is how to guarantee maximum stable yield in a resource population harvested at rate λ X t members perunit time. The harvested population in model (1.3) can be written as: dX t dt = r X t − X t K ! , X (0) = x , (3.1)where r = r − λ and K = (1 − λ r ) K . We define F ( x ) = r x (cid:16) − xK (cid:17) .Equilibria points or constant solutions of (3.1), are X u = X s = K (1 − λ r ) , which is a non-trivial equilibrium point if λ < r . X u is unstable while X s is stable. For λ > r , i.e. if the harvesting e ff ort is very large, the population will die out. In this case X u = X s is anasymptotic growth value of the harvest population model. Since 1 < K for r > λ , this implies that theasymptotic values of harvesting population lower than the non-harvesting population; (See Fig 1).The function F ( x ) in model (3.1) is autonomous function, because it is independent of t and it is con-tinuously di ff erentiable ( class of C ). Thus it has a unique solution and its non-trivial solution to the initialvalue problem is [19 ? ] X t = X K X + ( K − X ) e − r t , X (0) = x , (3.2)The non-trivial solution (3.2) goes to the asymptotic value K as time goes to infinity, i.e., lim t →∞ X t = K ,for any X >
0. Hence X u = X makes dX t / dt > , whichfurther increases X t and the population rises towards K which is asymptotically stable. When x > K , dX t / dt < K . The function F ( x ) has maximum value at r K which isobtained by substituting X = K in Eq. (3.1). The deterministic model (3.1) can be written as dX t dt = − ∂ Vdx where V is the potential function defined by V ( x ) = − Z r x − xK ! dx = − − r x + r K x . The potential function has a local minimum corresponding to the stable equilibrium and a local maximum x = r − λ >
0. The system has only one stable equilibrium, so it iscalled monostable.For r − λ >
0, the population converges to the stable equilibrium X s = K (cid:16) − λ r (cid:17) , and the yield at thestable equilibrium, called the sustainable yield is ¯ K = λ X s = λ K (cid:16) − λ r (cid:17) . From this we can calculate thefishing e ff ort that maximize is λ MS Y = r which is called maximum sustainable yield (MSY) is ¯ λ = r K ,and the corresponding stable equilibrium is K max = K .From Figure 1b, we can observe that when λ = x < K , the phase point moves faster and fasteruntil it reaches K , and dXdt reaches its maximum value rK . While the phase point approaches to wardscarrying capacity K if K < x < K and x > K . 5 a) λ vary. (b) λ and initial value vary. Figure 1: The phase line and trajectories of dX t dt = rX t (1 − X t / K ) − λ X t . (a) As the value of harvesting e ff ort is su ffi ciently big ( λ > r ), the population extinction occurs. (b) The solution of model (1.6) for di ff erent value of λ and x . Here we can see that as λ increases, the population size X t goes to zero and X t has S − shape when x < K . While K < x < K and x > K , the populationsize approaches to K as t → ∞ . When λ , K ), the phase point moves fasterand faster until it reaches K , and dXdt reaches its maximum value r K . While if K < x < K and x > K ,the phase point goes to wards K .In a biological view, this tells us that the population initially growth faster and faster [5] and the graphof X t is concave up. But dXdt starts to decrease if the initial value passes half of carrying capacity K or halfof asymptotic value K . In this case, X t has concave down shape. For initial value below half of carryingcapacity K or half of asymptotic value K , X t has S -shaped; ( see Figure 1b ). Stochastic logistic-harvest model without Allee e ff ect We will consider stochastic perturbation of the logistic-harvest model without Allee e ff ect (1.6). dX t = (cid:20) rX t (cid:18) − X t K (cid:19) − λ X t (cid:21) dt + ǫ X t dB t , X = x . (3.3)Eq. (3.4) can be transformed into the form of the SDE as in our previous paper [39] and rewritten as dX t = ( r − λ ) X t − X t (cid:16) − λ r (cid:17) K dt + ǫ X t dB t , X = x . (3.4)Since this model has four parameters, we non-dimensionalize by rescaling the population size (variable) andtime. Then the new model (or SDE) will have fewer parameters. Because studying the qualitative behaviourof a SDE (or model) with many parameters is di ffi cult. Define Y = X t K (cid:16) − λ r (cid:17) = X t K , τ = ( r − λ ) . t = r t The new model becomes dY = Y (1 − Y ) d τ + ǫ YdB (cid:18) τ r − λ (cid:19) , Y = x K (cid:16) − λ r (cid:17) = y , (3.5)6r dY = Y (1 − Y ) d τ + ǫ √ r − λ YdB τ , Y = y , (3.6)where ǫ is a positive constant representing random growth e ff ects ( 0 < ǫ < r > λ . B τ is a Brownianmotion which has independent and stationary increments with stochastically continuous sample paths.The solution of the model 3.6 is a homogenous di ff usion process with the drift coe ffi cient µ ( t , y ) = y (1 − y ) and di ff usion term υ ( t , y ) = ǫ r − λ y . Finding the exact solution of the nonlinear SDE in (3.6) issimilar with [[41], Section 9.3]. Set a new variable Z = Y and apply Itˆo formula [11, 41]. Our goal is thatto reduce the nonlinear SDE in terms of Y in to a linear SDE in Z , which we then able to solve. Thus weget a new linear SDE dZ = " ǫ r − λ − ! Z + d τ − ǫ √ r − λ ZdB τ , (3.7) Z = n . According [36] and [[41], Theorem 9.4] , the solution of Eq. (3.7) is (a) λ = . ǫ = . x vary. (b) λ = . ǫ = . x vary. Figure 2: Sample solutions of dX t = h rX t (cid:16) − X t K (cid:17) − λ X t i dt + ǫ X t dB t . (a) λ = . ǫ = . x vary. (b) λ = . ǫ = . x vary. When x ∈ (0 , K ) and x > K , the population approaches its maximum population size.Parameters r = K = Z = ϕ τ Z + Z τ ϕ − s ds ! , τ ≥ , where ϕ τ = exp (cid:18) ǫ √ r − λ (cid:19) − ! τ − ǫ √ r − λ B τ ! . Since Y = Z , we obtain the unique, strong solution of equation (3.6) Y = y exp − (cid:18) ǫ √ r − λ (cid:19) ! τ + ǫ √ r − λ B τ !(cid:16) + y R t ϕ − s ds (cid:17) . (3.8)7rom equation (3.8), we observe that the solution exists for all τ > y >
0, then Y > √ < ǫ √ r − λ , then − (cid:18) ǫ √ r − λ (cid:19) ! τ + ǫ √ r − λ B τ = " − (cid:18) ǫ √ r − λ (cid:19) ! + ǫ √ r − λ B τ τ τ goes to −∞ as time τ → ∞ .According to the strong law of large numbers, we apply B τ τ = τ → ∞ . From this we have Y → τ → ∞ .When the value of Gaussian noise intensity ǫ is small, the solution in (3.8) become a solution of thedeterministic model in (1.3), i. e., lim ǫ → Y = + (cid:16) y − (cid:17) e − τ . The Euler-Maruyama method was implemented [10] in order to give an approximation for the samplepaths solution of the stochastic model. Some sample solution paths are plotted in Figure 2. We observe thatthe sample solution paths are positive.
Extinction probability
This subsection deals with the transition density function p ( y , τ ) for the process Y = { Y τ , τ > } whichsatisfies the following theorem. The stationary density gives important long time information about theprobabilistic behaviour of the solution of a given SDE. Theorem 1. ( Fokker-Plank equation (FPE)): [ Simon (2019) [36], Theorem 5.4]. The probability densityp ( x , t ) of the solution of the SDE in (3.6) solves the partial di ff erential equation ∂ p ∂τ = − ddy ( y (1 − y ) p ) + γ d dy ( y p ) , (3.9)where γ = ǫ √ r − λ and with initial condition p ( y s | y τ ) = δ ( y τ − y s ) for τ ≥ s . Proof : See Simon (2019) [36]. x S t a t i ona r y den s i t y p =0.75 =0.125=0.25=0.375=0.5 (a) λ = .
75 and ǫ vary. x S t a t i ona r y den s i t y p =0.2 =0.36=0.84=0/93==0.96 (b) ǫ = . λ vary. Figure 3: Stationary densities of model (1.6) for r =
1. (a) λ = . ǫ = . , . , . , .
5. (b) ǫ = . λ = . , . , . , .
8n our case, the density p satisfies the time-independent FPE. i.e. Eq. (1) is the second order di ff erentialequation as in Mackeri˘ c ius [41] ddy ( y (1 − y ) p ) − γ d dy ( y p ) = . Noting that p ( y ) ≥
0, for all y ∈ (0 , ∞ ), [31] and R ∞ p ( y ) dy =
1. For 1 ≥ γ , the stationary density p in(0 , ∞ ) is p ( y ) = M y γ − e − y γ (here M is the normalizing constant). In (0 , ∞ ), the function p is integrable if γ − > − γ <
2. Setting λ = λ = λ , γ ≥ ff usion process (SDE) in (3.6) has no stationary density. That means populationbecomes extinct, but we have a noise-induced transition for 0 < γ < √
2. In this case extinction can notoccurs; ( Figure 3). In fact lim y → p ( y ) = , > γ M = , = γ ∞ , γ < < γ The next step is to show how to find the maximum point y max of p ( y ). Using p ′ ( y ) = y max , so we have y max = − γ . When γ becomes small, y max =
1. In this case the stationary densitybecome more acute.
4. Logistic-harvest model with Allee e ff ect Deterministic logistic-harvest model with Allee e ff ect Now let’s nondimensionalize the Allee e ff ect model (1.4) dX t dt = rX t (cid:18) X t S − (cid:19) (cid:18) − X t K (cid:19) − λ X t , which helps to rescale variables such that the rescaled model has fewer parameters. Let’s rescale populationsize X t by expressing it relative to the carrying capacity K (scaling by S would work as well).Setting Y t = X t K and β = KS . The new di ff erential equation has the following form: dY t dt = rY t ( β Y t − − Y t ) − λ Y t , Y (0) = y , (4.1)where y = x K . Our new model has just three parameter, which makes the bifurcation analyses, computationof equilibria, etc. more transparent.It is clear that if λ =
0, then the logistic-harvesting model with Allee e ff ect in Eq. (1.4) reduces to thenon-harvesting logistic growth model with Allee e ff ect as given in Eq. (1.2).Equation (4.1) has one trivial equilibria point at Y t = Y t = β + ± √ ( β + − β (1 + λ/ r )2 β if ( β + − β (1 + λ/ r ) >
0. Thus model (4.1) has the following three equilibrium9
Phase line diagram Y Y (a) λ = . X -0.2-0.100.10.20.30.40.5 V ( x ) Y m (b) λ = . λ = .
2. (d) λ vary. Figure 4: (a) Phase line diagram of (1.4). (b) Potential function V of the model (1.4). (c) λ = .
2. In this cases, model (1.3) isalways positive while model (1.4) is negative when the population X t < S . (d) λ vary. Parameters r = S = K = points. Y = , Y = ( β + − p ( β + − β (1 + λ/ r )2 β , Y = ( β + + p ( β + − β (1 + λ/ r )2 β . (4.2) Y and Y are stable equilibria separated by unstable equilibrium Y . Set m = r ( ( β − β ) which is called thecritical point.Clearly, if we use λ > m in Eq. (4.2), no fixed point which shows Y = λ = m , there exists two equilibrium points, i.e. Y = Y m = β + β (unstable); ( See Fig. 4c).The non-trivial equilibrium point Y is an asymptotic growth value of the harvest model. Since Y < K for λ < m , this implies that the asymptotic values of harvesting fish population lower than the non-harvesting fish population; ( See Fig. 4d ). 10n Figure 4b, we plot the graph of the potential function V ( x ) for values of λ = .
15, defined by V ( x ) = − Z (cid:20) rx (cid:18) xS − (cid:19) (cid:18) − xK (cid:19) − λ x (cid:21) dx . In term of V ( x ), Eq. (1.4) can be written as: dxdt = − dVdx . The phase line diagram for Eq. (1.4) is shown in Figure 4a. Denoting the stable equilibrium point by Y and the unstable equilibrium point by Y , the separation between the two equilibrium points is Y − Y .If λ < m , Figure 4b shows the potential function V ( x ) has two local minima corresponding to the stableequilibrium Y and Y and one local maximum at Y which is an unstable equilibrium. The function V ( x )is called a double-well potential, because the two stable equilibrium Y and Y separated by an unstableequilibrium Y .From the biological point of view, it is meaningful to choose β >
1, and 0 < λ < r (cid:18) ( β + β − (cid:19) ( or0 < λ < m ), and the state Y represents to the population free state that means it is the state of populationextinction, in this case no population are present. The state Y implies the state of stable population, wherethe population density does not increase but stays at a constant level.The number of equilibrium points depends on the sign of m , where m = ( β + − β (1 + λ/ r ).( i ) If m >
0, then there are two fixed points; (Two non-trivial equilibrium points),( ii ) If m =
0, then there is only one point; (One non-trivial equilibrium points),( iii ) If m <
0, then there is no fixed point; (No non-trivial equilibrium points).When m < , i.e. λ > m the population will go extinct as t → ∞ . As far as biological meaning isconcerned, we must catch at a harvest rate λ < r (cid:18) ( β + β − (cid:19) . So in this case the model in (4.1) has twoequilibria, one stable Y and one unstable Y with Y < Y . In Figure 4c shows the phase line plots for Eq. (effort) harvestseparation Y - Y = 2.18 Figure 5: Blue: Harvest yield λ Y ( λ ) versus λ . Red: Separation Y ( λ ) − Y ( λ ) versus λ ). (4.1), dY t dt = rY t ( β Y t − − Y t ) − λ Y t for increasing λ . If λ is less than the critical point m , there are two11table equilibrium solutions and one unstable equilibrium solution. As λ increases beyond m , there is onestable equilibrium solution.Using r = . β = λ Y ( λ ) versus eort λ and separation Y ( λ ) − Y ( λ ) versus λ . Note that from Figure 5 we get maximum yield when λ = .
18. For the value of λ , the separation between the two equilibrium solutions is Y ( λ ) − Y ( λ ) = . Stochastic logistic-harvest model with Allee e ff ect We consider the dimensionaless stochastic perturbation of the logistic-harvest model with Allee e ff ect(1.7) by setting a new variable Y t = X t K . dY t = [ rY t ( β Y t − − Y t ) − λ Y t ] dt + ǫ Y t dB t , Y (0) = x K = y , (4.3)where β = KS > , K is the carrying capacity and S is the Allee parameter with 0 < S < K . B t is a Brownianmotion with stochastically continuous sample paths, as well as independent and stationary increments. Thestochastic perturbation of the logistic-harvest model with Allee e ff ort is discussed in [1, 7, 16, 21]. (a) λ = . ǫ =
0, and initial value x vary. (b) λ = . ǫ = .
02 and initial value x vary. Figure 6: Sample solutions of stochastic logistic-harvest model with Allee e ff ect. (a) the solution of (4.3) with λ = . ǫ = x vary. (b) λ = . ǫ = .
02 and x vary. Here it is clearly seen that when x < S , the populationextinct occurs. While for x > S and x > K , the population size approaches to its maximum size K . Parameters r = S = K = The non-trivial solution of the SDE in (4.3 can be found as follows [34, 27, 9]. Having in mind that S < Y (0) < K , let’s define C − function Z t : R + → R + as Z t = log( Y t ), and apply Itˆo formula to Z t , thesystem in (4.3) is converted to a SDE with additive noise (to remove any state or level-dependent noise fromthese trajectories): dZ t = f ( Z t ) dt + g ( Z t ) dB ( t ) , t ∈ (0 , T ) , (4.4)12here f ( Z t ) = r ( β e Z t − − e Z t ) − λ − ǫ , g ( Z t ) = ǫ. Note that f ( Z t ) = f ( Y t ) − ǫ =
0. This shows the equilibrium point of the deterministicterm of the additive noise system in (4.4) is a ff ected by the Gaussian noise intensity ǫ .Now we will show that Y t is the solution of the SDE in (4.3). Since Y t = e Z t , apply Itˆo formula to have dY t = de Z t = e Z t dZ t + e Z t ( dZ t ) = e Z t " r ( β e Z t − − e Z t ) − λ − ǫ ! dt − ǫ dB ( t ) + e Z t ǫ dt = Y t [ r ( β Y t − − Y t ) − λ ] dt + ǫ dB ( t )] . This solution is strong, continuous and positive, for S < Y (0) and 0 < S < K . The numerical simulation (solution paths) of the stochastic di ff erential equation in (4.3) is shown inFigure 6 with various initial values. To plot this we use the Euler-Maruyama method. As we can see inFigure 6, the sample path are positive and approaching to the carrying capacity Y when 0 < λ < . Whileit goes to extinction when λ ≥ . From Figure 6c we observe that all trajectories, except x = Y ( unstableequilibrium point) fall in to a potential pit ( x = x = Y ) the stable equilibrium points.The Euler-Maryuama approximation was used to approximate the solution of stochastic model. Di ff er-ent values of the constant in the drift coe ffi cient λ were applied.Next we prove that sample paths of X t of SDE (1.7) are uniformly continuous for a.e. t ≥
0. To showthis consider the following integral X t = X (0) + Z t f ( X ( s )) dt + Z t g ( X ( s )) dB s , (4.5)where f ( X s ) = rX s (cid:16) X s S − (cid:17) (cid:16) − X s K (cid:17) − λ X s , g ( X s ) = ǫ X s and 0 < S < X (0) < K . Suppose 0 < a < b < ∞ , b − a ≤
1, and p >
2. By applying the well known H ¨older inequality andmoment inequality for Itˆo integrals (4.5), we have E | X t − X s | p ≤ p − ( b − a ) p − Z ba E [ f ( X s )] p ds + p − p ( p − ! p Z ba E [ g ( X s ] p ds . (4.6)Since, E [ f ( X s )] p ≤ K ( K − S ) S r ! p + ( λ K ) p , E [ g ( X s )] p ≤ ( K ǫ ) p , The equation in (4.6) can be estimated by E | X t − X s | p ≤ Q ( b − a ) p , where Q = p − K p (( ( K − S ) S ) p r p + λ p + ( p ( p − ) p ǫ p ) . According to Kolmogorov-Centsov theorem on the continuity of a stochastic process [27], we knowthat almost every sample path of X t is locally but uniformly H ¨older-continuous with exponent 0 < γ < p − p .Therefore the SDE in (1.7) has uniformly continuous solution on t ≥
0. All solutions of this model goesto zero as t → ∞ . Since Y t = X t K , so the SDE in (4.3) has uniformly continuous solution on t ≥ t → ∞ . 13 .3. Extinction probability and first passage probability
This subsection explains where the probability that the population extinct will happens, and the proba-bility of reaching a large population size L before reaching a small one. Trajectories that start in the potentialwell on the right will eventually jump into the potential well on the left, even though it may take a very longtime. Once there, they rapidly move to the region around x = Y t , Eq. (4.3) can be approximated by dY t = − ( r + λ ) dt + ǫ Y t dB t . (4.7)Then the boundary value problem for the probability P ( y ) of exit at 0 before exit at L is12 ǫ y p ′′ − ( r + λ ) yp ′ = , p (0) = , p ( L ) = , (4.8)which has the solution p ( y ) = − y + ( r + λ ) /ǫ L + ( r + λ ) /ǫ Note that P ( y ) → L → ∞ , even though we are using small y approximations for values of y that are notsmall. We would obtain the same result even if we solved the probability of exit problem corresponding toEq. (4.3).According to Theorem 1, the Fokker-Planck equation corresponding to Eq. (4.3) is ∂ p ∂ t = ǫ ∂ ∂ y [ y p ] − ∂∂ y [ ry ( β y − − y ) − λ y ) p ] . (4.9)Since all solutions of Eq. (4.3) go to zero with probability one no matter what the starting point y , weexpect that the solution of Eq. (4.9) satisfieslim t →∞ p ( y , t / y ) , < y < ∞ . In other words, all populations eventually become extinct. However, it is reasonable, or realistic, values ofthe parameters, if y is in the potential well on the right in Figure 4b, it will take a very long time before thetrajectory jumps across the potential barrier into the potential well on the left. In this case it makes sense tolook at a quasi-stationary density [41], say q ( y ), that is obtained by solving12 ǫ ∂ ∂ y [ y p ] − ∂∂ y [ ry ( β y − − y ) − λ y ) p ] = Y m , ∞ ] where Y m is the location of the local maximum of the potential function, shown atthe vertical blank line in Figure 4b. The quasi-stationary density q ( y ), given by q ( y ) = y − − r + λ ) /ǫ e [ − r β y + r (1 + β ) y ] /ǫ R ∞ Y m y − − r + λ ) /ǫ e [ − r β y + r (1 + β ) y ] /ǫ dy is shown in Figure 7. 14 igure 7: The quasi-stationary density of model (4.3)
5. Numerical experiments
We summarize our numerical findings about the impact of parameters x , λ and ǫ on the solution of thedeterministic and stochastic models of logistic-harvest with and without Allee e ff ect.Here we apply the Euler-Maruyama (EM) method following [8] to Eq. (3.6). To apply this method inthe SDE (3.6) over time [0 , T ], we first need discretize the interval. For any positive n assume ∆ t = T / n ,and s j = jt , for j = , , ..., n . The numerical approximation to the solution X ( s j ) is denoted by X j . As in[8], the EM method has the following form: X j = X j − + f ( X j ) ∆ t + g ( X j )( B ( X j ) − B ( X j − ) , j = , , , .., n . (5.1) Numerical results and biological implications of logistic-harvest model without Allee e ff ect The phase line and trajectories of dX t dt = rX t (1 − X t / K ) − λ X t is plotted in Figure 1. Parameters r = K = ≤ λ ≤ r . In Fig. 1a as the value of harvesting e ff ort is su ffi ciently large ( overfishing ), thepopulation extinction occurs. Here the value of X u becomes small as λ increases. In Fig. 1b the solution ofmodel (1.6) for di ff erent values of λ and x . In this figure, we can observe that as λ increases, the populationsize X t decreases i.e. X t goes to zero and it has S − shape when x < K . While K < x < K and x > K ,the population size approaches to K as t → ∞ . When λ , x < K , the phase point moves faster andfaster until it reaches K , and dXdt reaches its maximum value r K . While if K < x < K and x > K , thephase point goes to wards K . The biological implications of this result tells us that the population initiallygrows faster and faster [5] and the graph of X t is concave up. But dXdt starts to decrease if x > K or x > K .In this case, the shape of X t is concave down.Figure 2 shows the numerical simulation of model dX t = [ rX t (1 − X t K ) − λ X t ] dt + ǫ X t dB t with fixedparameters r = K =
3. For λ = . ǫ = . x vary is plotted in Fig 2a. Fig.2b presents numerical simulation of stochastic logistic-harvest model without Allee e ff ect with λ = . ǫ = . ff erent values of x . We use f ( X j ) = X j [ r (1 − X j K ) − λ ] and g ( X j ) = ǫ X j in equation (5.1).15hen x ∈ (0 , K ) and x > K , the population approaches to its maximum size. In both deterministic andstochastic models the behaviour of the solution of the models are almost similar. In other words, for anypositive initial value x , X t goes to K as t → ∞ .The analysis of the stationary density of model (1.6) which varies under a proportional increase inGaussian noise ǫ and λ is drawn in Figure 3 for 0 < ǫ √ r − λ < r =
1. In this case we have noise-induced transition. In Fig. 3a the value of λ is 0.75 and ǫ = . , . , . , .
5. This tells us thatproportional increase in linear multiplicative noise can qualitatively change the behavior of the system.When ǫ becomes smaller and smaller then stationary density p ( y ) becomes more and more acute, and themaximum point y max of p ( y ) tends to one. From Fig. 3b we can see that the probabilistic qualitativebehaviour of the stationary densities is similar with Fig. 3a with fixed ǫ = . ff erent valueof harvest rate λ ( λ = . , . , . , . ǫ and λ ( if both go zero),then stationary density p ( y ) becomes more acute, and the maximum point of p ( y ) tends to y max =
1. If ǫ √ r − λ ≥ √
2, the population dynamic system in (1.6) has no stationary densities. This shows that thesolutions converge to zero (extinction). From this graph we know that the values of the extrema of stationarydensity depend on the noise intensity ǫ and harvest rate λ . Numerical results and biological implications of logistic-harvest model with Allee e ff ect We fixed the value of the parameters r = S = K =
3. Figure 4 plots about the phase line diagramof model (1.3) and (1.4), and potential function V of the model (1.4). Fig. 4a shows us deterministic model(1.4) has three equilibrium solution at X t = X t = Y (unstable) and X t = Y ( stable). In Fig. 4b,the area below Y m is an absorbing zone. For fixed λ = .
2, Fig. 4c shows model (1.3) is always positivewhile model (1.4) is negative when the population X t < S . The phase line diagram of model (1.4) is plottedin Fig. 4d for λ vary (or 0 < λ < m ) . We observe that the value of Y is smaller than the value of K , but Y > S .Figure 5 plots the harvest yield λ Y ( λ ) versus eort λ and separation Y ( λ ) − Y ( λ ) versus λ . From thisFigure, we obtain maximum yield when λ = . r = . β = λ , the separationbetween the two equilibrium solutions is Y ( λ ) − Y ( λ ) = . ff ect model is given in Figure 6. Fig.6a shows the solution of (4.3) with fixed value harvest rate λ = . ǫ = x vary. In Fig. 6b plots the numerical solution of (4.3) with noise ( ǫ = .
02) and λ = . x vary. Clearlyseen that the population extinct occurs when x < S . While for x > S and x > K , the population sizeapproaches to its maximum size K .The graph in Figure 7 presents the quasi-stationary density of model (4.3). Here p ( y ) → L → ∞ . Inother words the probability of reaching 0 (lower size) before reaching L (maximum size), when consideredas a function of initial population size, p ( y ) [2] has an inection point at deterministic unstable equilibrium Y . In this Figure inection point is Y m .
6. Conclusion
We have studied the logistic-harvest model with and without Allee e ff ect driven by multiplicative Gaus-sian noise. For the stochastic logistic-harvest model without Allee e ff ect we obtained exact solution, butfor stochastic logistic-harvest model with Allee e ff ect we proved the stability of the solution process. Weanalyzed the stationary density and the probability of reaching a large population size before reaching asmall one, for the stochastic models (1.6) and (1.7). 16ur numerical experiments demonstrated that the stochastic models in population growth are di ff erentfrom the deterministic models. The main result of our study is that the stochastic model under Gaussiannoise perturbation is asymptotically stable. This matches with an important result of the fishery theory.In the case of the logistic-harvest model without Allee e ff ect, when the harvesting rate λ is less than thegrowth rate r , we observe that there exists two equilibrium solutions X u = X s , which are less than thecarrying capacity of the population K . However, if the harvesting rate λ is greater than r (overfishing), thereis no fixed point at all and therefore no equilibrium solution.In the case of the logistic-harvest model with the Allee e ff ect, if the harvesting rate λ is equal to thecritical threshold m = r (cid:16) (1 + β )4 β − (cid:17) , we find that there exists only one nonzero equilibrium state. Thisequilibrium population is less than the carry capacity K . However, if the harvest rate λ is less than m , wehave two positive equilibrium solutions, both the stable and unstable equilibrium solutions are lower thanthe carrying capacity K and the unstable Y is less than the stable equilibrium solution Y . If the harvestingrate λ is greater than the critical threshold m , there is no fixed point at all and therefore no equilibriumsolution.As far as biological meaning is concerned, we have to catch at a harvest rate λ less than the growth rate r in the case of the logistic-harvest model(1.3), and less than the critical point in the case of the logistic-harvest model (1.4). Acknowledgments
This work was partly supported by the NSFC grants 11801192, 11771449 and 11531006.
ReferencesReferences [1] A. Friedman,
Stochastic Di ff erential Equations and Applications. Vol. 2, Academic Press San Diego, 1976.[2] B. Dennis, Allee e ff ects in stochastic populations, OIKOS (389401), 2002.[3] B. Dubey, P. Chandra, P. Sinha.A model for fishery resource with reserve area. Nonlinear Analysis: Real World Applications, ( 625-637), 2003.[4] B.Oksendal, Stochastic Di ff erential Equations: An Introduction with Applications . Springer, 2003.[5] C. Braumann, Introduction to Stochastic Di ff ential Equations with Application to Modellinf in Biology and Finance. Willey,2019.[6] C. Fima Kelebaner,
Introduction to stochastic calculus with applications , (Imperial College Press) 2nd Ed, London WC2H9HE, 2005.[7] D. Applebaum,
L´evy Processes and Stochastic Calculus , 2nd ed, Cambride Univ. Press, Cambridge, 2009.[8] D. Higham J., An algorithmic introduction to numerical simulation of stochastic di ff erential equations, SIAM Rev. (525-C546), 2001.[9] D. Jiang, N. Shi, Y. Zhao.Existance, uniquence and global stability of positive solutions to the food-limited population modelwith random perturbation.Math.Comp. Model (651-658), 2005.[10] D. Heath E. Platen. A Benchmark Approach to Quantitative Finance. SpringerVerlag, Berlin, 2006.[11] E. Allen, Modeling with Itˆo Stochastic Di ff erential Equations . ISBN:-13 978-1-4020-5952-0(HB), 2007.[12] F. Braue, C. C. Chavez, Mathematical models in population boilogy and epidemicology , 2nd Ed. Springer, 2011.[13] F. Courchamp, L.Berec, J. Gascoigne,
Allee E ff ects in Ecology and Conservation, Oxford Univ, 2008.[14] G. Grimmett R., D. Stirzaker R.,
Probability and Random Processes . 3rd Ed. Oxford: Oxford Univ. Press, 2006.[15] H. Kinfe, T. Zebebe, Mathematical modeling of Fish Resources Harvesting with Predator at Maximum Sustainable Yield.
Math. Theory Modeling . (2224-5804), 2016.[16] J. Duan An Introduction to Stochastic Dynamics . (Cambridge Univ. Press), 2015.[17] J. Drake M.,D. Lodge M., Allee e ff ects, propagule pressure and the probability of establishment: risk analysis for biologicalinvasions. Biological Invasions , (365-75), 2006.[18] J. Meiss D., Di ff erential Dynamical Systems , Society for Industrial and App. Math. Philadelphia,PA,USA, 2007.
19] J. Murray D.,
Mathematical Biology. An Introduction , Volume I, Springer Verlag, New York, 2002.[20] L. Arnold,
Random Dynamic Systems ,Springer,Germany, 2003.[21] L. Arnold,
Stochastic Di ff erential Equations: Theory and Applications , Wiley, New York, 1972.[22] L. Ricciardi M., Stochastic Population Theory: Birth and Death Processes. In: Hallam T. G., Levin S. A. (eds) MathematicalEcology. Biomathematics, vol 17. Springer, Berlin, Heidelberg, 1986.[23] M. Hasanbulli, P. Svitlana Rogovchenko, V.Yuriy Rogovchenko,Dynamics of a Single Species in a Fluctuating Environmentunder Periodic Yield Harvesting. J. Appl. Math.
Article ID 167671, 2013.[24] M. Panik J.,
Growth Curve Modeling: (Theory and Application) . ISBN:978-1-118-76404-6, 2014.[25] M. Panik J.,
Stochastic Di ff erential Equations: An Introduction with Applications in Population Dynamics Modeling . WilleyUSA, 2017.[26] M. Kot, Elements of Mathematical Ecology . Cambridge Univ. Press. New York, NY, 2011.[27] M. Krstic, M. Jovanovic, On stochastic population model with the Allee e ff ect, Math. and Computer Modelling (370-379),2010.[28] M. Rahmani H. Doust, The logistic modelling population; having harvesting factor. Yugoslav J. of Operations Research , (107-115) DOI:10.2298 / YJOR130515038R, 2015.[29] N. Britton F.,
Essential Mathematical Biology (Springer undergraduate) .Univ. of Bath, Claverton Down, Bath BA2 7AY, UK.mathematics series, 2003.[30] N. Saber Elaydi, J.Robert Sacker, Population model with Allee E ff ect: A new model, J.Bio.Dyam.Vol . , No. , 2009.[31] P. E. Kloeden, E.Platen, Numerical Solution of Stochastic Di ff erential equations. Springer Science & Business Media, 2013.[32] P. Amarasekare, Allee E ff ects in Metapopulation Dynamics.vol. , no. American naturalist , 1998.[33] P. Bhattacharya, S. Paul, K. S.Choudhury, Maximum Sustainable Yield Policy in PreyPredator System-A Study.2nd Interna-tional conference on Computing Communication and Sensor Network, 2013.[34] Q. Yang, D. Jiang, A note on asymptotic behaviors of stochastic population model with Allee e ff ect. App.MathModel . (4611-4619), 2011.[35] S. Lynch, Dynamical Systems with Applicatin using MATLAB , 2004 .[36] S. Simo , S. Arno,
Applied Stochastic Di ff erential Equations .Cambridge Univ. Press, 2019.[37] S. Pinheiro, A Stochastic Logistic Growth Model with Predation: An Overview of the Dynamics and Optimal Harvesting. In:Pinto A., Zilberman D. (eds) Modeling, Dynamics, Optimization and Bioeconomics III. Springer, (313-330), 2016.[38] S. Pinheiro, Optimal harvesting for a logistic growth model with predation and a constant elasticity of variance, Ann OperRes (461-480), 2018.[39] T. Almaz, T. Daniel,K. Anas, J. Brannan, Mean exit time and escape probability for the stochastic logistic growth model withmultiplicative α − stable L´evy noise. https: // arxiv.org / abs / ff ort Dynamics in a Prey-Predator Model with Harvesting. International Journal ofInformation and Systems Sciences, (318-332), 2010.[41] V. Mackeri˘ c ius, Introduction to Stochastic Analysis: Integral and Di ff erential Equation. WILEY, 2011.[42] X. Mao,E. Renshaw, G. Marion, Environmental Brownian noise suppresses explosions in population dynamic,
StochasticProcess.Appl . (95-110), 2002.(95-110), 2002.