About the Properties of a Modified Generalized Beverton-Holt Equation in Ecology Models

Abstract

This paper is devoted to the study of a generalized modified version of the well-known Beverton-Holt equation in ecology. The proposed model describes the population evolution of some species in a certain habitat driven by six parametrical sequences, namely, the intrinsic growth rate (associated with the reproduction capability), the degree of sympathy of the species with the habitat (described by a so-called environment carrying capacity), a penalty term to deal with overpopulation levels, the harvesting (fishing or hunting) regulatory quota, or related to use of pesticides when fighting damaging plagues, and the independent consumption which basically quantifies predation. The independent consumption is considered as a part of a more general additive disturbance which also potentially includes another extra additive disturbance term which might be attributed to net migration from or to the habitat or modeling measuring errors. Both potential contributions are included for generalization purposes in the proposed modified generalized Beverton-Holt equation. The properties of stability and boundedness of the solution sequences, equilibrium points of the stationary model, and the existence of oscillatory solution sequences are investigated. A numerical example for a population of aphids is investigated with the theoretical tools developed in the paper.

1. Introduction

The use of mathematical models in ecology is nowadays of important research interest, [1–10] because such models lead to a more precise study of the dynamics and interactions of populations among them and with the habitat. This paper is devoted to the study of a modified version of the well-known Beverton-Holt equation in ecology which describes the population evolution of some species in a certain habitat. The associate dynamics evolves subject to an intrinsic growth rate associated with the reproduction capability and the degree of sympathy of the species with the habitat described by a so-called environment carrying capacity, [2–6, 11, 12]. The standard Beverton-Holt equation becomes the so-called generalized Beverton-Holt equation [12], when extended by considering in the model two extra sequences, namely, the harvesting (namely, fishing or hunting) quota, or regulation through the use of pesticides and an additive disturbance sequence. Such a disturbance sequence consists of an independent consumption sequence, which describes a predation of the species though which it gains, develops, and grows plus a net migration sequence. The migration term plays a close role to that of the independent consumption in the population evolution and it is due to potential positive/negative migration towards/from the habitat under study. In general, the joint contribution of independent consumption plus positive or negative migration (resp., immigration and emigration) towards or from the habitat could also potentially include other extra effects of correction of contribution of population in the model. Some of those effects are, for instance, unexpected accidental mortality not fixed through properly adjusting the intrinsic growth rate, loss of population due to harvesting which has not been parameterized through the harvesting quota or measuring modeling errors in general. The study in this paper adds two more sequences to the generalized Beverton-Holt equation to conform a penalty term for eventual overpopulation at intermediate sampling points. This is the main novelty of the paper related to previous literature on the subject. The penalty term corrects the population evolution so that large amounts of population translate into the decrease of the number of individuals in the future. The motivation of the use of the penalty term is clear since it is well known that large amount of individuals in a habitat makes the overall population decrease by natural regulation, for instance, due to competition for food or difficulties for nesting. The properties of boundedness and stability of the solution, existence of equilibrium points of the limit stationary version of the equation, as well as the existence of oscillatory solutions are investigated. The following basic notation is being used through the manuscript.

Notation 1.1. 𝐍 and 𝐑 are the sets of natural and real numbers, respectively. 𝐍₀∶=𝐍∪{0};𝐑₊∶={𝑧∈𝐑∶𝑧>0};𝐑₀₊∶=𝐑₊∪{0}={𝑧∈𝐑∶𝑧≥0}.𝑘∶={1,2,…,𝑘};𝑘⁰∶=𝑘∪{0}={0,1,2,…,𝑘}. (1.1)𝑔^′(𝑥) denotes the first derivative with respect to 𝑥 of the differentiable real function𝑔(𝑥).
The symbols “⋀" and “∨" stand for logic conjunction and disjunction, respectively.
𝐼_𝑛 is the identity matrix of the 𝑛th order.

The population evolution is supposed to be governed by a modified Beverton-Holt equation, 𝑥_𝑘+1=𝛾_𝑘𝜇_𝑘𝐾_𝑘𝑥^𝑝𝑘_𝑘𝑒^{−𝛼𝑘𝑥𝑘}𝑥_𝑘−1𝐾_𝑘+𝜇_𝑘−1𝑥_𝑘−1+𝑑_𝑘;𝑘∈𝐍 (1.2) subject to initial conditions 𝑥₋₁>0,𝑥₀≥0, which is parameterized by the following real sequences.

𝛾_𝑘∶=1−𝑏_𝑘,𝑘∈𝐍, with 𝑏_𝑘∈[𝑏₁,𝑏₂] being the harvesting (fishing or hunting) quota subject to 0<𝑏₁≤𝑏₂≤1. Thus, 𝛾_𝑘∈[𝛾₁,𝛾₂],𝑘∈𝐍, is the harvesting quota excess subject to 0<𝛾₁∶=1−𝑏₂≤𝛾₂∶=1−𝑏₁≤1. If repopulation activity is incorporated to harvesting, then 𝑏₁<0 when repopulation dominates harvesting so that 𝛾₂>1.

𝜇_𝑘∈[𝜇₁,𝜇₂],𝑘∈𝐍, with ∞>𝜇₂≥𝜇₁>1 is the intrinsic growth rate which is associated to the population ability for reproduction.

𝐾_𝑘∈[𝐾₁,𝐾₂],𝑘∈𝐍, with ∞>𝐾₂≥𝐾₁>0, is the environment carrying capacity which is associated to the sympathetic degree, more or less favorable, of the habitat towards the species.

𝑑_𝑘∈[𝑑₁,𝑑₂],𝑘∈𝐍, with ∞>𝑑₂≥𝑑₁>−∞, is the external disturbance term which refers to joint effects of independent consumption plus migration contribution, namely, immigration (i.e., positive net migration from outside towards the habitat), null, or emigration (i.e., negative net migration from the habitat to outside) including atypical mortality not included in the model by the growth rate: the disturbance. The additive disturbance contribution has the generic form 𝑑_𝑘=𝑑_0𝑘+𝑑_𝑚𝑘 where 𝑑_0𝑘 is the independent consumption and 𝑑_𝑚𝑘 is the net migration contribution. The decomposition of 𝑑_𝑘=𝑑_0𝑘+𝑑_𝑚𝑘 into the two separated parts 𝑑_𝑚𝑘 and 𝑑_𝑚𝑘 is irrelevant for the modified generalized Beverton- Holt equation from an analysis point of view but it is quite relevant from a biological insight point of view, as pointed out in detail by one of the reviewers. The reason is that it allows the inclusion in the same evolution equation of two quite different effects, in nature, like predation and net migration.

𝑝_𝑘∈[𝑝₁,𝑝₂],𝛼_𝑘∈[𝛼₁,𝛼₂], with 𝛼_𝑘=0⇒𝑝_𝑘=0,𝑘∈𝐍, are sequences which conform a penalty term 𝑥^𝑝𝑘_𝑘𝑒^{−𝛼𝑘𝑥𝑘} to describe how the excess of population translates into a contributed tendency to its decrease at the next sampling time. The sequences are assumed to be subject to the constraints ∞>𝑝₂≥𝑝₁≥0,∞>𝛼₂≥𝛼₁≥0, and 𝛼₁=0⇒𝑝₁=0.

A particular case of (1) is obtained with 𝛾_𝑘=1 (no harvesting activity), 𝑝_𝑘=𝛼_𝑘=0 (no penalty correcting term for overpopulation control), 𝑑_𝑘=0 (zero jointly net independent consumption plus migration) resulting in the well-known standard Beverton-Holt equation, 𝑥_𝑘+1=𝜇_𝑘𝐾_𝑘𝑥_𝑘𝐾_𝑘+𝜇_𝑘−1𝑥_𝑘,𝑘∈𝐍, (1.3) which is a particular version of the Hassell model 𝑥_𝑘+1=𝜇_𝑘𝑥_𝑘/(1+(𝜇_𝑘−1)𝑥_𝑘/𝐾_𝑘)^𝑐𝑘,𝑘∈𝐍, for 𝑐_𝑘=1 and where a double sampling period is used related to (1.2). In other words, only the recruitment and storage measurements are relevant to describe the population evolution since there is not necessary the evaluation at intermediate samples of a penalty term associated with overpopulation. Another particular case of (1.2) is obtained with 𝑝_𝑘=𝛼_𝑘=0 resulting in the so-called generalized Beverton-Holt equation, 𝑥_𝑘+1=𝛾_𝑘𝜇_𝑘𝐾_𝑘𝑥_𝑘𝐾_𝑘+(𝜇_𝑘−1)𝑥_𝑘+𝑑_𝑘,𝑘∈𝐍, (1.4) also described equivalently as a two-stage evolution equation, 𝑥⁻_𝑘+1=𝜇_𝑘𝐾_𝑘𝑥_𝑘𝐾_𝑘+𝜇_𝑘−1𝑥_𝑘,𝑥_𝑘+1≡𝑥⁺_𝑘+1=𝛾_𝑘𝑥⁻_𝑘+1+𝑑_𝑘,𝑘∈𝐍, (1.5) which includes discontinuities at sampling times due to harvesting and disturbance activities associated with impulses in the corresponding logistic continuous-time differential equation. Finally, note that the modified Beverton-Holt equation leads to the Ricker model (which is also a limiting equation of the Hassell model and that of the Beverton- Holt one if 𝑐_𝑘=1): 𝑥_𝑘+1=𝑒^{𝜇𝑘(1−𝑥𝑘/𝐾𝑘)}𝑥_𝑘,𝑘∈𝐍, (1.6) if 𝑑_𝑘=0 and a harvesting quota 𝛾_𝑘=((𝐾_𝑘+𝜇_𝑘−1𝑥_𝑘−1)/𝐾_𝑘(𝜇_𝑘+𝑥_𝑘−1))𝑒^𝜇𝑘,𝑘∈𝐍, is applied provided that 𝐾_𝑘>𝑒^𝜇𝑘𝜇_𝑘−1/(1+𝜇_𝑘−1/𝑥_𝑘−1),𝑘∈𝐍, in order to ensure that 𝛾_𝑘<1. Thus, the environment carrying capacity satisfies liminf_{𝑥𝑘−1→∞}𝐾_𝑘>𝑒^𝜇𝑘𝜇_𝑘−1 and 𝐾_𝑘∈(𝑒^𝜇𝑘𝜇_𝑘−1−𝑓(𝑥_𝑘−1),∞),𝑘∈𝐍. The function 𝑓∶𝐑₊→𝐑₀₊ is defined by 𝑓(𝑥)=𝜇²𝑒^𝜇/(𝑥+𝜇) which is decreasing on 𝐑₊ with lim_𝑥→∞𝑓(𝑥)=0 with 𝐑₀₊∶={𝑧∈𝐑∶𝑧≥0} and 𝐑₊∶={𝑧∈𝐑∶𝑧>0}=𝐑₀₊⧵{0}. To investigate issues on the standard Beverton- Holt equation, it is useful to analyze the solution sequence of its inverse [11]. The mathematical properties of the modified generalized Beverton- Holt equation including boundedness properties, equilibrium points, and existence of oscillations are investigated. A numerical example is proposed for the study of the evolution of a population of aphids, a very resistive species of insects which causes considerable damage and reproduces according to several reproduction cycles per year [13, 14].

2. Boundedness of the Solutions

There is a wide variety of works concerning boundedness of the solutions, stability, investigation of equilibrium points, and existence of oscillatory solutions in discrete recursive sequences of several types (see, e.g., [15–23]). Parallel problems are also relevant in continuous-time differential systems [16, 24–26]. To investigate those problems in the context of the modified Beverton-Holt equation (1.2), define 𝑠_𝑘∶=𝑥_𝑘⁻¹,𝑘∈𝐍, being the solution sequence of the inverse of the modified Beverton-Holt equation provided that 𝑥_𝑘≠0,𝑘∈𝐍. Thus, one gets directly from (1.2) 𝑠⁻¹_𝑘+1=𝛾_𝑘𝜇_𝑘−1𝐾_𝑘−1𝜁_𝑘𝑠⁻¹_𝑘−1𝐾_𝑘−1+𝜇_𝑘−1−1𝑠⁻¹_𝑘−1+𝑑_𝑘=𝛾_𝑘𝜇_𝑘−1𝐾_𝑘−1𝜁_𝑘𝐾_𝑘−1𝑠_𝑘−1+𝜇_𝑘−1−1+𝑑_𝑘,𝑘∈𝐍, (2.1) where 𝜁_𝑘∶=𝑥^𝑝_𝑘𝑒^{−𝛼𝑥𝑘},𝑘∈𝐍. Then, 𝑠_𝑘+1=𝐾_𝑘−1+𝜇_𝑘−1−1𝑥_𝑘−1𝛾_𝑘𝜇_𝑘−1𝐾_𝑘−1𝜁_𝑘𝑥_𝑘−1+𝐾_𝑘−1+𝜇_𝑘−1−1𝑥_𝑘−1𝑑_𝑘=𝐾,𝑘∈𝐍_𝑘−1𝑠_𝑘−1+𝜇_𝑘−1−1𝛾_𝑘𝜇_𝑘−1𝐾_𝑘−1𝜁_𝑘+𝐾_𝑘−1+𝑠_𝑘−1+𝜇_𝑘−1−1𝑑_𝑘,𝑘∈𝐍 (2.2)(2.3) which can be rewritten as 𝑠_𝑘+1=𝐾_𝑘−1𝛾_𝑘𝜇_𝑘−1𝐾_𝑘−1𝜁_𝑘𝑠_𝑘−1+𝜇_𝑘−1−1𝛾_𝑘𝜇_𝑘−1𝐾_𝑘−1𝜁_𝑘−𝐾_𝑘−1𝑠_𝑘−1(𝐾_𝑘−1𝑠_𝑘−1+𝜇_𝑘−1−1)𝑑_𝑘𝛾_𝑘𝜇_𝑘−1𝐾_𝑘−1𝜁_𝑘𝛾_𝑘𝜇_𝑘−1𝐾_𝑘−1𝜁_𝑘+𝐾_𝑘−1𝑠_𝑘−1+𝜇_𝑘−1−1𝑑_𝑘,𝑘∈𝐍 (2.4) for the standard Beverton-Holt equation (𝑑_𝑘=𝛼_𝑘=𝑝_𝑘≡0 implying 𝜁_𝑘≡1 and 𝛾_𝑘≡1 and double sampling rate; see (1.3)), (2.2)–(2.4) becomes in particular [11] 𝑠_𝑘+1=𝜇_𝑘⁻¹𝑠_𝑘+1−𝜇_𝑘⁻¹𝐾_𝑘⁻¹,𝑘∈𝐍. (2.5)

The boundedness of the solutions for bounded initial conditions is a minimal requirement in order that a biological model is well posed. The subsequent result is concerned with boundedness and global Lyapunov's stability of the solution sequences of (1.2) and (2.3)–(2.4).

Theorem 2.1. The following properties hold.
(i)𝑑₁≤𝑥_𝑘<∞, for all 𝑘∈𝐍.(ii)If 𝑑₁>0, then {𝑥_𝑘}^∞₀ and {𝑠_𝑘}^∞₀ are positive bounded sequences. As a result, the modified Beverton-Holt equation (1.2) is globally Lyapunov's stable.(iii)If 𝑑₁≥0, then {𝑥_𝑘}^∞₀ and {𝑠_𝑘}^∞₀ are nonnegative sequences, the first one being bounded from above. As a result, the modified Beverton-Holt equation (1.2) is globally Lyapunov's stable. Furthermore, {𝑠_𝑘}^∞₀ is positive with limsup_𝑘→∞𝑠_𝑘≤∞.(iv)If 𝑑_𝑘≥0 and 𝑥_𝑘≥0, then 0≤𝑥_𝑘+1<∞, for all 𝑘∈𝐍₀.

Proof. (i) 𝑥_𝑘≥𝑑_𝑘≥𝑑₁, for all 𝑘∈𝐍 follows directly from (2.4). Also, 𝜁_𝑘=1 if 𝛼_𝑘=𝑝_𝑘=0. Otherwise, 𝛼_𝑘≠0 so that 𝜁_𝑘≤max_𝑘∈𝐍𝜁_𝑘∶𝛼_𝑘𝑝_𝑘≠0=max_𝑘∈𝐍𝜁𝑝_𝑘𝛼_𝑘∶𝛼_𝑘𝑝_𝑘≠0=max_𝑘∈𝐍𝑝_𝑘𝛼_𝑘𝑒^𝑝𝑘∶𝛼_𝑘𝑝_𝑘≠0≤𝑝₂𝛼₁𝑒^𝑝2 (2.6) since 𝑑𝜁(𝑥)/𝑑𝑥=0 for 𝛼≠0 only at its relative maximum 𝑥=𝑝/𝛼 so that 𝜁(𝑝/𝛼)=(𝑝/𝛼𝑒)^𝑝=max_{𝑥∈𝐑0+}𝜁(𝑥) since 𝜁(0)=0 and lim_𝑥→∞𝜁(𝑥)=0 implying that the relative maximum is also the absolute maximum. As a result, 𝜁_𝑘≤max(1,(𝑝₂/𝛼₁𝑒)^𝑝2)<∞,forall𝑘∈𝐍. Then, from (1.2) 𝑥_𝑘+1≤𝛾₂𝜇₂𝐾₂𝑥_𝑘−1𝐾₁+𝜇₁−1𝑥_𝑘−1max1,𝑝₂𝛼₁𝑒^𝑝2+𝑑₂≤𝛾₂𝜇₂𝐾₂max1,𝑝₂/𝛼₁𝑒^𝑝2𝜇₁+−1𝑑₂<∞, (2.7) for all 𝑘∈𝐍 and (i) has been proved.
If 𝑑₁>0, then (i) implies 0<𝑑₁≤𝑥_𝑘<∞ and 0<𝑠_𝑘≤𝑑₁⁻¹<∞, for all 𝑘∈𝐍 and (ii) has been proved.
If 𝑑₁≥0, then (i) implies 0≤𝑥_𝑘<∞ and 0<𝑠_𝑘≤∞, for all 𝑘∈𝐍 and, thus, liminf_𝑘→∞𝑥_𝑘≥0,limsup_𝑘→∞𝑥_𝑘<∞,liminf_𝑘→∞𝑠_𝑘>0, and limsup_𝑘→∞𝑠_𝑘≤∞. Property (iii) has been proved. Property (iv) follows directly from (1.2).

For formulation coherency, the population cannot be negative at any sample. This idea motivates the following axiom and simple-related assertions which will be then useful for some mathematical proofs concerning the case of the eventual negativity of the joint disturbance contribution.

Axiom 2.2. 𝑥_𝑘≥0,forall𝑘∈𝐍.

Assertion 2.3. 𝑥_𝑘=0⇒𝑑_𝑘≥0 for any 𝑘∈ℕ.

Assertion 2.4. 𝑑₂<0⇒𝑥_𝑘>0,forall𝑘∈𝐍.

Assertion 2.5. 𝑑_𝑘<0⇒𝜁_𝑘>0 for any 𝑘∈𝐍.

Assertion 2.6. 𝜁_𝑘=0 only if 𝑥_𝑘=0 and 𝑑_𝑘≥0 for any 𝑘∈𝐍.

Proof. Assertion 2.4 is equivalent to prove that its contrapositive proposition 𝑥_𝑘≤0forsome𝑘∈𝐍⇒𝑑₂≥0 is true. Assume that 𝑥_𝑘≤0 for some 𝑘∈𝐍. Then, Axiom 2.2 implies that 𝑥_𝑘=0 so that 𝜁_𝑘=0 and then 𝑥_𝑘+1=𝑑_𝑘≥0 from (1.2) and Axiom 2.2 so that Assertion 2.3 is proved. Also, 𝑥_𝑘≤0⇒𝑑₂≥max_𝑗∈𝐍(𝑑_𝑗)≥𝑑_𝑘≥0⇔𝑑₂<0⇒𝑥_𝑘>0,forall𝑘∈𝐍 what proves Assertion 2.4. Assertion 2.5 follows since Axiom 2.2⋀𝑑_𝑘<0⟹𝜁_𝑘=𝐾_𝑘−1+𝜇_𝑘−1−1𝑥_𝑘−1|𝑑_𝑘|𝛾_𝑘𝜇_𝑘−1𝐾_𝑘−1𝑥_𝑘−1⟺𝑥_𝑘+1≥≥0𝜇₁−1𝑑_𝑘𝛾₂𝜇₂𝐾₂+𝐾₁|𝑑_𝑘|𝛾₂𝜇₂𝐾𝑥_𝑘−1≥𝜇₁−1𝑑_𝑘𝛾₂𝜇₂𝐾₂>0forany𝑘∈𝐍. (2.8) Assertion 2.6 follows since for any 𝑘∈𝐍,𝜁_𝑘=0⇒𝑥_𝑘=0∨𝑥_𝑘=∞.𝑥_𝑘=∞ is impossible from Theorem 2.1(i) so that for any 𝑘∈𝐍,𝜁_𝑘=0⇒𝑥_𝑘=0. If 𝑥_𝑘=0, then 𝑥_𝑘+1≥0 from Axiom 2.2 so that 𝑑_𝑘≥0 which is Assertion 2.3.

The extinction of the standard Beverton-Holt equation is investigated in [27] when the intrinsic growth rate is less then unity. Note that eventual extinction is also admitted by Axiom 2.2 for intrinsic growth rate exceeding unity and negative disturbance contribution. This situation includes as particular cases that of local emigration plus zero-independent consumption and that of negative-independent consumption with no net migration contribution. However, the axiom prevents against eventual negative levels of population which is not physically possible.

Remark 2.7. Note that Assertion 2.5 together with an intermediate result in the proof of Theorem 2.1 implies 𝑑_𝑘<0⇒∞>max(1,(𝑝₂/𝛼₁𝑒)^𝑝2)≥𝜁_𝑘>0 for any 𝑘∈𝐍, 𝑑_𝑘≥0⟹∞>max1,𝑝₂𝛼₁𝑒^𝑝2≥𝜁_𝑘≥0forany𝑘∈𝐍. (2.9) Then, ∞>𝜁_𝑘≥0, for all 𝑘∈𝐍, and ∞>limsup_𝑘→∞𝜁_𝑘≥liminf_𝑘→∞𝜁_𝑘≥0. Also, if 𝑑₂<0 or if 𝑑_𝑘<0 (although a nonnegative upper-bound be assumed), for all 𝑘∈𝐍, then ∞>limsup_𝑘→∞𝜁_𝑘≥liminf_𝑘→∞𝜁_𝑘>0.

Now, Axiom 2.2 and Assertions 2.3–2.6 extend the results of Theorem 2.1 concerning positivity of the solutions and further stability results as follows.

Theorem 2.8. The following properties hold under Axiom 2.2.
(i){𝑥_𝑘}^∞₀ and {𝑠_𝑘}^∞₀ are nonnegative sequences which are bounded and positive and bounded from below, respectively, independently of the disturbance contribution. As a result, the modified Beverton-Holt equation (1.2) is globally Lyapunov's stable. Furthermore, if 𝑑_𝑘>0,forall𝑘∈𝐍, both sequences are positive and bounded.(ii)Define 𝜌_𝑘∶=1/𝛾_𝑘𝜇_𝑘−1𝜁_𝑘,𝛿_𝑘∶=𝐾_𝑘−1𝑠_𝑘−1+𝜇_𝑘−1−1=(𝐾_𝑘−1+(𝜇_𝑘−1−1)𝑥_𝑘−1)/𝑥_𝑘−1,𝜂_𝑘∶=𝑣_𝑘−𝜔_𝑘,𝑣_𝑘∶=(𝜇_𝑘−1−1)/(𝜌_𝑘⁻¹+𝛿_𝑘𝑑_𝑘), and 𝜔_𝑘∶=𝐾_𝑘−1𝑠_𝑘−1𝛿_𝑘𝑑_𝑘/𝜌_𝑘⁻¹(𝜌_𝑘⁻¹+𝛿_𝑘𝑑_𝑘). If 𝑑_𝑘>0,forall𝑘∈𝐍, then either 0<^𝑘_𝑖=0𝜌_𝑖<∞;∀𝑘∈𝐍𝑡ℎ𝑒𝑛0<liminf^𝑘_{𝑘→∞𝑖=0}𝜌_𝑖≤limsup^𝑘_{𝑘→∞𝑖=0}𝜌_𝑖∧<∞0≤^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝜂_𝑖<∞⟹0≤^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝑣_𝑖≤∞∧0≤^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝜔_𝑖≤∞,∀𝑘∈𝐍,𝑜𝑟lim^𝑘_{𝑘→∞𝑖=0}𝜌_𝑖=−lim_𝑘→∞^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝜂_𝑖⟹=∞lim_𝑘→∞^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝜔_𝑖=∞. (2.10)(iii)If 𝑑_𝑘≤0,forall𝑘∈𝐍, then either 0<^𝑘_𝑖=0𝜌_𝑖∧≤∞;∀𝑘∈𝐍0≤^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝜂_𝑖≤∞⟹0≤^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝑣_𝑖≤∞∧0≤^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝜔_𝑖≤∞,∀𝑘∈𝐍,𝑜𝑟lim^𝑘_{𝑘→∞𝑖=0}𝜌_𝑖=lim_𝑘→∞^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝜂_𝑖⟹=∞lim_𝑘→∞^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝜔_𝑖=∞. (2.11)

Proof. (i) Note from (1.2), Axiom 2.2, and Assertion 2.5 that, for 𝑑_𝑘<0, 𝜁_𝑘>0∧𝑑_𝑘=−𝑑_𝑘≤𝛾_𝑘−1𝜇_𝑘−1𝜇_𝑘−1𝜁_𝑘𝐾_𝑘−1+𝜇_𝑘−1−1𝑥_𝑘−1=𝛾_𝑘−1𝜇_𝑘−1𝜇_𝑘−1𝜁_𝑘𝑠_𝑘−1𝐾_𝑘−1𝑠_𝑘−1+𝜇_𝑘−1−1>0 (2.12) so that from (2.4) 𝑠_𝑘+1≥𝐾_𝑘−1𝛾_𝑘𝜇_𝑘−1𝐾_𝑘−1𝜁_𝑘𝑠_𝑘−1+𝜇_𝑘−1−1𝛾_𝑘𝜇_𝑘−1𝐾_𝑘−1𝜁_𝑘+𝐾_𝑘−1𝑠_𝑘−1𝐾_𝑘−1𝑠_𝑘−1+𝜇_𝑘−1−1𝑑_𝑘𝛾_𝑘𝜇_𝑘−1𝐾_𝑘−1𝜁_𝑘𝛾_𝑘𝜇_𝑘−1𝐾_𝑘−1𝜁_𝑘−𝐾_𝑘−1𝑠_𝑘−1+𝜇_𝑘−1−1𝑑_𝑘≥𝐾_𝑘−1𝛾_𝑘𝜇_𝑘−1𝐾_𝑘−1𝜁_𝑘𝑠_𝑘−1>0⟺0≤𝑥_𝑘+1<∞ (2.13) for any 𝑘∈𝐍₀ if 𝑑_𝑘<0. On the other hand, 𝑑_𝑘≥0⇒0≤𝑥_𝑘+1<∞ and 0<𝑠_𝑘+1≤∞ for any 𝑘∈𝐍₀ from Theorem 2.1(iv). Furthermore, 𝑑_𝑘>0⇒0<𝑥_𝑘+1<∞ and 0<𝑠_𝑘+1<∞ for any 𝑘∈𝐍₀ from Theorem 2.1(ii). As a result, Property (i) has been fully proved.
(ii) Recursive direct calculations from (2.4) yield taking (i) into account irrespective of the value of the disturbance contribution: 0<𝑠_𝑘=^𝑘_𝑖=0𝜌_𝑖𝑠₀+^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝜂_𝑖=^𝑘_𝑖=0𝜌_𝑖𝑠₀+^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝑣_𝑖−𝜔_𝑖≤∞,∀𝑘∈𝐍 (2.14) and, furthermore, if 𝑑_𝑘>0;forall𝑘∈𝐍 then 0<𝑠_𝑘=^𝑘_𝑖=0𝜌_𝑖𝑠₀+^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝜂_𝑖=^𝑘_𝑖=0𝜌_𝑖𝑠₀+^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝑣_𝑖−𝜔_𝑖<∞,∀𝑘∈𝐍 (2.15) since {𝑥_𝑘}^∞₀ is positive and bounded from Theorem 2.1. Thus, the first two constraints linked via logic conjunction follow directly from (2.4) follows directly from (2.14) since either ∏^𝑘_𝑖=0[𝜌_𝑖]<∞,forall𝑘∈𝐍, so that |∑^𝑘−1_𝑖=0∏^𝑘−1_𝑗=𝑖+1[𝜌_𝑗]𝜂_𝑖|<∞,forall𝑘∈𝐍, or lim_𝑘→∞∏^𝑘_𝑖=0[𝜌_𝑖]=−lim_𝑘→∞(∑^𝑘−1_𝑖=0∏^𝑘−1_𝑗=𝑖+1[𝜌_𝑗]𝜂_𝑖)=∞, or ∏^𝑘_𝑖=0[𝜌_𝑖]=−(∑^𝑘−1_𝑖=0∏^𝑘−1_𝑗=𝑖+1[𝜌_𝑗]𝜂_𝑖)=∞ for finite 𝑘∈𝐍. However, from the equivalent contrapositive logic proposition to Assertion 2.6, 𝑥_𝑘>0 (which holds from Theorem 2.1, since 𝑑_𝑘>0) ⇒𝜁_𝑘>0 which implies max(∏^𝑘_𝑖=0[𝜌_𝑖],|∑^𝑘−1_𝑖=0∏^𝑘−1_𝑗=𝑖+1[𝜌_𝑗]𝜂_𝑖|)<∞ from (2.15) and then ∏^𝑘_𝑖=0[𝜌_𝑖]=−(∑^𝑘−1_𝑖=0∏^𝑘−1_𝑗=𝑖+1[𝜌_𝑗]𝜂_𝑖)=∞ for finite 𝑘∈𝐍 is impossible. Property (ii) has been proved.
(iii) It is proved using similar techniques as those used in the proof of (ii) by noting that |𝑑_𝑘| is upper bounded with the bound of (2.12) for 𝑑_𝑘≤0 and 𝛿_𝑘>0,𝜌_𝑘>0,𝜁_𝑘≥0,∞≥𝜂_𝑘≥0,forall𝑘∈𝐍. The detailed proof is omitted.

Remark 2.9. Note that the use of the contrapositive logic proposition to Assertion 2.6 is not feasible, as invoked in the last part of the proof of Theorem 2.8(ii), in order to prove some close result concerned with Theorem 2.8(iii). This follows since 𝑑_𝑘<0⇒𝜁_𝑘>0,forall𝑘∈𝐍 from Assertion 2.5.

3. Equilibrium Points

Assume that lim_𝑘→∞𝐾_𝑘=𝐾,lim_𝑘→∞𝜇_𝑘=𝜇, lim_𝑘→∞𝛾_𝑘=𝛾,lim_𝑘→∞𝛼_𝑘=𝛾,lim_𝑘→∞𝑝_𝑘=𝛾,lim_𝑘→∞𝑑_𝑘=𝑑. It is well known that both the standard Beveron-Holt equation, driven by two parameters, and the generalized one, driven by four parameters, have equilibrium points. In particular, the equilibrium point under nonextinction condition of the standard Beverton-Holt equation is the carrying capacity of the environment; see [2–6, 11, 12, 13, 27, 28]. The following result holds, which is concerned with the existence of equilibrium points in the limit equation (1.2), that is, the stationary solutions lim_𝑘→∞𝑥_𝑘=𝑥≥0 of (1.2) as the parameters converge to those finite limits.

Theorem 3.1. The following properties hold.
(i)If 𝑑>0, then it exists a unique equilibrium point 𝑥=𝑥 on 𝐑₀₊ which is, in addition, positive being subject either to 𝑥_𝑀>𝑥≥𝑥_𝑏, which requires 𝛼∈[0,(𝑝+1)/𝑑) as a necessary condition, or 𝑥_𝑀≤𝑥<𝑥_𝑏, which requires 𝛼≥(𝑝+1)/𝑑.(ii)If 𝑑<0, then an equilibrium point 𝑥=𝑥 on 𝐑₀₊ can exist only if the following two necessary conditions hold:(i)|𝑑|≤(1/𝐾)((𝑝+1)/𝛼𝑒)^𝑝+1;(ii)𝐾+(𝜇−1)|𝑑|<𝛾𝜇𝐾(𝑝+1)𝑥^𝑝𝑒^−𝛼𝑥−(𝛾𝜇𝐾𝛼𝑥^𝑝𝑒^−𝛼𝑥+2(𝜇−1))𝑥 for 𝑥∈𝑋_𝑀⊂[0,𝑥_𝑀). If there exists such an equilibrium point, then at most two distinct equilibrium points satisfying 0<𝑥₁≤𝑥_𝑀≤𝑥₂ can exist.(iii)If 𝑑=0, then 𝑥₁=0 is an equilibrium point of (1.2). Another equilibrium point might exist only if 𝐾[𝛾𝜇𝐾𝑝𝑥^𝑝𝑒^−𝛼𝑥(𝑝+1−𝛼𝑥)−1]>2(𝜇−1)𝑥 on some interval 𝑋_𝑀⊂[0,𝑥_𝑀)

Proof. The equilibrium points of (1.2), if any, are the (real) values of 𝑥∈𝐑₀₊ such that 𝑓(𝑥)∶=(𝜇−1)𝑥²+𝐾−𝑑(𝜇−1)𝑥−𝑑𝐾=𝑔(𝑥)∶=𝛾𝜇𝐾𝜁(𝑥)𝑥. (3.1) Note that 𝑓∶𝐑→𝐑 is a convex parabola satisfying 𝑓(0)=−𝐾𝑑 whose zeros are 𝑥_𝑎,𝑏=(𝜇−1)𝑑−𝐾∓𝐾−(𝜇−1)𝑑²+4(𝜇−1)𝐾𝑑^1/2=−𝐾2(𝜇−1)𝑑𝜇−1<0 (3.2) since ((𝐾−(𝜇−1)𝑑)²+4(𝜇−1)𝐾𝑑)^1/2=𝐾+(𝜇−1)𝑑. Also, 𝑔∶𝐑→𝐑 satisfies 𝑔(0)=0,𝑔(∞)=lim_𝑘→∞𝑔(𝑥)=0 and has a relative maximum, which is also the absolute maximum in 𝐑₀₊) at 𝑥_𝑀=(𝑝+1)/𝛼 which is 𝑔_𝑀=((𝑝+1)/𝛼𝑒)^𝑝+1=max(𝑔(𝑥)∶𝑥∈𝐑₀₊). Also, 𝑔(𝑥) is increasing on 𝑥∈[0,𝑥_𝑀) and decreasing on 𝑥∈[𝑥_𝑀,∞). Thus, the equilibrium points 𝑥=𝑥, if any, related to 𝑥_𝑏 and 𝑥_𝑀 satisfy the subsequent constraints.
(a)If 𝑑>0, then 𝑓(0)=−𝐾𝑑<𝑔(0)=0. Thus, there is a unique equilibrium point 𝑥=𝑥 (i.e., 𝑓(𝑥)=𝑔(𝑥)) on 𝐑₀₊ which is, in addition, positive which is subject either to 𝑥_𝑀>𝑥≥𝑥_𝑏, which requires 𝛼∈[0,(𝑝+1)/𝑑) as a necessary condition, or 𝑥_𝑀≤𝑥<𝑥_𝑏, which requires 𝛼≥(𝑝+1)/𝑑.(b)If 𝑑<0, then 𝑓(0)=𝐾|𝑑|>𝑔(0)=0 and 𝑥_𝑎≤𝑥_𝑏<0. Thus, an equilibrium point 𝑥=𝑥 (i.e., 𝑓(𝑥)=𝑔(𝑥)) on 𝐑₀₊ can exist only if the following two necessary conditions hold:(i)|𝑑|≤(1/𝐾)((𝑝+1)/𝛼𝑒)^𝑝+1 implying 𝑓(0)≤𝑔_𝑀 since 𝑓(𝑥) is increasing on [𝑥_𝑏,∞) so that it is also increasing on 𝐑₀₊ since 𝐑₀₊⊂[𝑥_𝑏,∞);(ii)𝐾+(𝜇−1)|𝑑|<𝛾𝜇𝐾(𝑝+1)𝑥^𝑝𝑒^−𝛼𝑥−(𝛾𝜇𝐾𝛼𝑥^𝑝𝑒^−𝛼𝑥+2(𝜇−1))𝑥 for 𝑥∈𝑋_𝑀⊂[0,𝑥_𝑀) since 𝑓(0)>𝑔(0) and 𝑓(𝑥) being increasing on [0,∞) implies that solutions to 𝑓(𝑥)=𝑔(𝑥) might exist only if 𝑓(0)≤𝑔_𝑀 and the continuous function ℎ(𝑥)∶=𝑔(𝑥)−𝑓(𝑥) is monotonically strictly increasing (i.e., ℎ^′(𝑥)=𝑔^′(𝑥)−𝑓^′(𝑥)>0 on some interval 𝑋_𝑀⊂[0,𝑥_𝑀).
If there exists such an equilibrium point, then it can exist at most two distinct equilibrium points satisfying 0<𝑥₁≤𝑥_𝑀≤𝑥₂.
(c)If 𝑑=0, then 𝑓(0)=𝑔(0)=0 so that 𝑥₁ is an equilibrium point of (1.2). Another equilibrium point might exist only if 𝐾[𝛾𝜇𝐾𝑝𝑥^𝑝𝑒^−𝛼𝑥(𝑝+1−𝛼𝑥)−1]>2(𝜇−1)𝑥 in order that ℎ(𝑥)∶=𝑔(𝑥)−𝑓(𝑥) be monotonically strictly increasing on some interval 𝑋_𝑀⊂[0,𝑥_𝑀) on some interval 𝑋_𝑀⊂[0,𝑥_𝑀).

Remark 3.2. Particular cases of equilibrium points of interest are the following.
(i)The case 𝑝=𝛼=𝑑=0. Thus, the constraint (3.1) becomes 𝑓(𝑥)=(𝜇−1)𝑥²+𝐾𝑥=𝑔(𝑥)=𝛾𝜇𝐾𝑥 which yields an equilibrium point at 𝑥=𝐾(𝛾𝜇−1)/(𝜇−1), provided that 1>𝛾>𝜇⁻¹. If 𝛾=1, one obtains the equilibrium point of the standard Beverton-Holt equation 𝑥=𝐾.(ii)The case 𝑝=𝛼=0;𝑑>0. Thus, 𝑥=𝛾𝜇𝐾𝑥/(𝐾+(𝜇−1)𝑥) holds for 𝑥=𝑥₁=𝐾(𝛾𝜇−1)/(𝜇−1) provided that 1>𝛾>𝜇⁻¹, and which implies that lim_𝑘→∞𝑥⁻_𝑘=𝑥₁. Also, 𝑥=𝛾𝑥+𝑑 holds for 𝑥=𝑥₂=𝑑/(1−𝛾)>0 since 𝑑>0 and 1>𝛾>0, which implies that lim_𝑘→∞𝑥_𝑘=lim_𝑘→∞𝑥⁺_𝑘=𝑥₂. Then, 𝑥=𝑥₁=𝑥₂ is an equilibrium point of (1.2), provided that 1>𝛾>𝜇⁻¹, implying that lim_𝑘→∞𝑥_𝑘=lim_𝑘→∞𝑥⁺_𝑘=lim_𝑘→∞𝑥⁻_𝑘=𝑥 if 𝑥=𝐾(𝛾𝜇−1)/(𝜇−1)=𝑑/(1−𝛾) if 𝑑=(1−𝛾)(𝛾𝜇−1)𝐾/(𝜇−1). Lower and upper bounds for the equilibrium points of the limit stationary equation (1.2) are investigated in the subsequent result. The use of those bounds is important when the exact equilibrium point cannot be calculated by imprecise knowledge of the model parameterization or computational difficulties.

Theorem 3.3. The following properties hold.
(i)The equilibrium points 𝑥_1,2 of the limit stationary modified Beverton-Holt equation (1.2) satisfy the following properties: 0≤𝑑=𝑥₂₂≤𝑥₂≤𝑥₁₂≤𝑑𝐾𝜇−1 (3.3) provided that 𝐾>𝑑(𝜇−1) so that 𝑥₂=0 if and only if 𝑑=0. If 𝑑<0 and has a sufficiently small modulus, then 0>−|𝑑|=𝑥₂₂≤𝑥₂≤𝑥₁₂<0. Also, −𝑑𝐾𝜇−1≥𝑥₁₁≤𝑥₁≤𝑥₂₁𝐾=−𝑥(𝜇−1)<0𝑖𝑓𝑑≥0,min₁₁,𝑥₁₂≤𝑥₁−𝐾≤max(𝜇−1),−|𝑑|<0,𝑖𝑓𝑑<0. (3.4)(ii)Define the real functions 𝜓(𝑥)=𝛾𝐾𝜇𝑥/(𝐾+(𝜇−1)𝑥) and 𝜁(𝑥)=𝑥^𝑝𝑒^−𝛼𝑥. Then, any equilibrium point 𝑥 of the limit stationary Beverton-Holt equation (1.2) is locally asymptotically stable if the two constraints below jointly hold: 1+𝜁𝑥𝜓^′𝑥>0,𝜓𝑥𝜁^′𝑥+𝜁𝑥𝜓^′𝑥<1. (3.5)

Proof. Note from (3.1) that ℎ₂(𝑥)∶=(𝜇−1)𝑥²+𝐾−𝑑(𝜇−1)𝑥−𝑑𝐾≥ℎ(𝑥)∶=𝑓(𝑥)−𝑔(𝑥)=(𝜇−1)𝑥²+𝐾−𝑑(𝜇−1)−𝛾𝜇𝐾𝜁(𝑥)𝑥−𝑑𝐾≥ℎ₁(𝑥)∶=(𝜇−1)𝑥²+𝑝𝐾−𝑑(𝜇−1)−𝛾𝜇𝐾𝛼𝑒^𝑝𝑥−𝑑𝐾 (3.6) since 𝜁(𝑥)≤(𝑝/𝛼𝑒)^𝑝, for all 𝑥∈𝐑₀₊. From (3.1), the equilibrium points are the nonnegative real solutions to ℎ(𝑥)=0, if any. The zeros of the convex parabola ℎ₁(𝑥) are 𝑥₁₁=𝑑(𝜇−1)+𝛾𝜇𝐾(𝑝/𝛼𝑒)^𝑝−𝐾−𝑑(𝜇−1)+𝛾𝜇𝐾(𝑝/𝛼𝑒)^𝑝−𝐾²+4(𝜇−1)𝑑𝐾^1/22(𝜇−1)≤−𝑑𝐾𝜇−1≤0 (3.7) if 𝑑≥0,𝑥₁₁>0 if 𝑑<0 of sufficiently small modulus, and 𝑥₁₁<0 if 𝑑<0 of sufficiently large modulus 𝑥₁₂=𝑑(𝜇−1)+𝛾𝜇𝐾(𝑝/𝛼𝑒)^𝑝−𝐾+𝑑(𝜇−1)+𝛾𝜇𝐾(𝑝/𝛼𝑒)^𝑝−𝐾²+4(𝜇−1)𝑑𝐾^1/2≤2(𝜇−1)𝑑𝐾𝜇−1 (3.8) if 𝑑≥0 and 𝑥₁₂<0 if 𝑑<0 of sufficiently small modulus or complex if 𝑑<0 of sufficiently large modulus.
The zeros of the convex parabola ℎ₂(𝑥) are 𝑥₂₁=𝑑(𝜇−1)−𝐾−(𝑑(𝜇−1)+𝐾^21/2𝐾2(𝜇−1)=−𝑥(𝜇−1)<0,₂₂=𝑑(𝜇−1)−𝐾+𝑑(𝜇−1)+𝐾^21/22(𝜇−1)=𝑑. (3.9) In view of (3.6), since ℎ₁(𝑥),ℎ₂(𝑥), and ℎ(𝑥) are convex parabolas, subject to (3.6), the zeros 𝑥_1,2 ofℎ(𝑥), which are the equilibrium points of the limit stationary (1.2), satisfy 0≤𝑑=𝑥₂₂≤𝑥₂≤𝑥₁₂=𝑑(𝜇−1)+𝛾𝜇𝐾(𝑝/𝛼𝑒)^𝑝−𝐾+𝛾𝜇𝐾(𝑝/𝛼𝑒)^𝑝+𝑑(𝜇−1)−𝐾²+4(𝜇−1)𝑑𝐾^1/2≤2(𝜇−1)𝑑𝐾𝜇−1 (3.10) provided that 𝑑≥0 and 𝐾>𝑑(𝜇−1) so that 𝑥₂=0 if and only if 𝑑=0. If 𝑑<0 of sufficiently small modulus, then 0>−|𝑑|=𝑥₂₂≤𝑥₂≤𝑥₁₂=𝛾𝜇𝐾(𝑝/𝛼𝑒)^𝑝−|𝑑|(𝜇−1)−𝐾+𝛾𝜇𝐾(𝑝/𝛼𝑒)^𝑝−|𝑑|(𝜇−1)−𝐾²−4(𝜇−1)|𝑑|𝐾^1/2−2(𝜇−1)<0,𝑑𝐾𝜇−1≥𝑥₁₁≤𝑥₁≤𝑥₂₁𝐾=−𝑥(𝜇−1)<0if𝑑≥0,min₁₁,𝑥₁₂≤𝑥₁−𝐾≤max(𝜇−1),−|𝑑|<0if𝑑<0, (3.11) so that Property (i) has been proved. To prove Property (ii), first note that the modified limit Beverton-Holt equation (1.2) may be written more compactly as 𝑥_𝑘+1−𝑑=𝜓𝑥_𝑘−1𝜁𝑥_𝑘, (3.12) where 𝜓(𝑥)=𝛾𝐾𝜇𝑥𝐾+(𝜇−1)𝑥,𝜁(𝑥)=𝑥^𝑝𝑒^−𝛼𝑥, (3.13) which implies 𝜓^′(𝑥)=𝛾𝐾²𝜇𝐾+(𝜇−1)𝑥²;𝜁^′(𝑥)=𝑥^𝑝−1𝑒^−𝛼𝑥(𝑝−𝛼𝑥). (3.14) The equilibrium points 𝑥, if any, satisfy 𝑥−𝑑=lim_𝑘→∞𝑥_𝑘−𝑑=𝜓(𝑥)𝜁(𝑥). The linearized dynamics of (3.12) about the equilibrium points are [19] 𝑧_𝑘+1=𝜁𝑥𝜓^′𝑥𝑧_𝑘−1+𝜓𝑥𝜁^′𝑥𝑧_𝑘 (3.15) whose associated characteristic equation in the complex indeterminate 𝜆 is 𝜆²−𝜓𝑥𝜁^′𝑥𝜆−𝜁𝑥𝜓^′𝑥=0 (3.16) so that the equilibrium point is locally asymptotically stable if both zeros of (3.16) have modulus less than unity, namely, 1>−𝜁𝑥𝜓^′𝑥>𝜓𝑥𝜁^′𝑥−1 (3.17) or, equivalently, 1+𝜁𝑥𝜓^′𝑥>0;𝜓𝑥𝜁^′𝑥+𝜁𝑥𝜓^′𝑥<1. (3.18) For 𝑑≥0, since, 𝑥=𝑥₂≤𝑥₁₂, it follows that 𝜓(𝑥)≤𝛾𝜇𝑥≤𝛾𝜇𝑥₁₂≤𝛾𝜇𝑑𝐾/(𝜇−1) and also that 𝜓(𝑥)≤𝛾𝐾𝜇/(𝜇−1), for all 𝑥∈𝐑₀₊ directly from (3.13) so that 𝜓(𝑥)≤min𝛾𝐾𝜇(𝜇−1),𝛾𝜇𝑑𝐾,𝑑𝜇−1_𝑝𝑒^{−𝛼𝑑𝐾/(𝜇−1)}≤𝜁𝑥≤max𝜁(𝑥)∶𝑥∈𝐑₀₊=𝑝𝛼𝑒^𝑝 (3.19) if 𝑑≥0 from the constraint. From that constraint and 0≤𝑑=𝑥₂₂≤𝑥₂≤𝑥₁₂ and (3.13)-(3.14), 𝛾𝐾²𝜇𝐾+𝑑𝐾(𝜇−1)²≤𝜓^′𝑥𝐾≤𝛾𝜇min1,²(𝜇−1)²𝑑². (3.20)𝜁^′(𝑥) defined in (3.14) has a derivative (𝑑/𝑑𝑥)𝜁^′(𝑥)=𝑥^𝑝−2𝑒^−𝛼𝑥[(𝑝−1)𝑝−(𝑝+1)𝛼𝑥] which has a relative extremum at 𝑥=(𝑝−1)𝑝/(𝑝+1)𝛼. Thus, taking into account also (3.14) and 0≤𝑑=𝑥₂₂≤𝑥₂≤𝑥₁₂ for 𝑑≥0, one gets 𝑑^𝑝−1𝑒^{−𝛼𝑑𝐾/(𝜇−1)}𝑝−𝛼𝑑𝐾𝜇−1≤𝜁^′𝑥≤min𝑑𝐾𝜇−1^(𝑝−1)/2𝑒^−𝑑𝛼(𝑝−𝛼𝑑),2𝑝𝑝+1(𝑝−1)𝑝(𝑝+1)𝛼^𝑝−1𝑒^{(𝑝−1)𝑝/(𝑝+1)}. (3.21) Thus, from (3.18), sufficient conditions for the equilibrium point 𝑥=𝑥₂ to be locally stable if 𝑝≥𝛼𝑑 are from (3.19) to (3.21): 𝑑^𝑝𝑒^{−𝛼𝑑𝐾/(𝜇−1)}𝛾𝐾²𝜇𝐾+𝑑𝐾(𝜇−1)²>−1,min𝛾𝐾𝜇(𝜇−1),𝛾𝜇𝑑𝐾𝜇−1min𝑑𝐾𝜇−1^(𝑝−1)/2𝑒^−𝑑𝛼(𝑝−𝛼𝑑),2𝑝𝑝+1(𝑝−1)𝑝(𝑝+1)𝛼^𝑝−1𝑒^{(𝑝−1)𝑝/(𝑝+1)}𝑝+𝛾𝜇𝛼𝑒^𝑝𝐾min1,²(𝜇−1)²𝑑²<1. (3.22) Theorem 3.1(i)-(ii) establish that there is a unique equilibrium point on 𝐑₊ for 𝑑>0 of the limit stationary equation (1.2) and that at most two equilibrium points might exist on 𝐑₊ for 𝑑<0 being of small absolute value. Those results are obtained by investigating the zeros of the upper-bounding and lower-bounding parabolas to ℎ(𝑥). The subsequent result gives explicit conditions of the existence of two equilibrium points in such a situation provided that extra sufficiency-type conditions on the remaining limit parameters are fulfilled.

Theorem 3.4. The following properties hold.
(i)The limit equation (1.2) has two positive equilibrium points if 𝑑<0 and satisfies |𝑑|<((𝛾𝜇(𝑝/𝛼𝑒)^𝑝−1)𝐾−(𝜇−1)(𝑝/𝛼)/((𝜇−1)𝑝+𝛼𝐾))𝑝, provided that 1>𝛾>(((𝜇−1)𝑝+𝐾𝛼)/𝜇𝐾𝛼)(𝛼𝑒/𝑝)^𝑝,𝐾>max(0,(𝜇−1)𝑝/𝜇𝛼[(𝑝/𝛼𝑒)^𝑝−1/𝜇]).(ii)If 𝑑=0, then the limit equation (1.2) has a zero equilibrium point and another positive equilibrium points, and 𝐾>(𝜇−1)(𝑝/𝛼)/(𝛾𝜇(𝑝/𝛼𝑒)^𝑝−1).

Proof. (i) One gets from (3.6) for 𝑑<0 that ℎ(𝑥)=(𝜇−1)𝑥²+𝐾+(𝜇−1)|𝑑|−𝛾𝜇𝐾𝜁(𝑥)𝑥+𝐾|𝑑|. (3.23) For 𝑥=𝑝/𝛼>0 where 𝜁(𝑥) reaches its maximum value (𝑝/𝛼𝑒)^𝑝, and ℎ𝑝𝛼=𝑝(𝜇−1)𝛼𝑝+𝐾𝛼+𝑝(𝜇−1)𝛼+𝐾|𝑑|−𝛾𝜇𝐾𝑒^𝑝𝑝𝛼^𝑝+1<0, (3.24) if 1>𝛾>(((𝜇−1)𝑝+𝐾𝛼)/𝜇𝐾𝛼)(𝛼𝑒/𝑝)^𝑝 subject to (𝑝/𝛼𝑒)^𝑝>((𝜇−1)𝑝+𝐾𝛼)/𝜇𝐾𝛼=(𝜇−1)𝑝/𝜇𝐾𝛼+1/𝜇 for 𝐾>0⇔𝐾>max(0,(𝜇−1)𝑝/𝜇𝛼[(𝑝/𝛼𝑒)^𝑝−1/𝜇]). Since ℎ(0)=𝐾|𝑑|>0,ℎ(∞)=lim_𝑥→∞ℎ(𝑥)=∞,ℎ(𝑝/𝛼)=((𝜇−1)(𝑝/𝛼)+𝐾)(𝑝/𝛼)+((𝜇−1)(𝑝/𝛼)+𝐾)|𝑑|−(𝛾𝜇𝐾/𝑒^𝑝)(𝑝/𝛼)^𝑝+1<0 for 𝑥=𝑝/𝛼∈𝐑₀₊ and ℎ∶𝐑₀₊→𝐑 is a real continuous function, there exist 𝑥_1,2∈𝐑₊ with 𝑥₂>𝑥₁ such that ℎ(𝑥_𝑖)=0,𝑖=1,2 which are equilibrium points of the limit equation (1.2) and the proof is complete.
(ii) The proof is similar to that of (i) by noting that the particular condition |𝑑|=0<(((𝛾𝜇(𝑝/𝛼𝑒)^𝑝−1)𝐾−(𝜇−1)(𝑝/𝛼))/((𝜇−1)𝑝+𝛼𝐾))𝑝 is guaranteed by the given constraint on 𝐾.

The following result proves that if there is only one equilibrium point (𝑑>0) of the limit equation, then it is a globally stable attractor. Also, if there are two equilibrium points (𝑑≤0 with sufficiently small |𝑑| plus extra conditions on the parameters), then the smaller equilibrium point is locally unstable while the largest one is a globally stable attractor.

Theorem 3.5. The following properties hold.
(i)If 𝑑>0, then the unique positive equilibrium point is a global stable attractor.(ii)If 𝑑≤0 and of sufficiently small modulus, satisfying the constraints of Theorem 3.4, then there are two nonnegative distinct equilibrium points. The smallest one is locally unstable if 𝑝=0, and locally stable if 𝑝>0 while the largest one is a global stable attractor.

Proof. (i) Consider the initial value problem 𝑥_𝑖=𝑥_𝑖≥0(𝑖=0,1). Thus, one has for 𝑑>00<𝑑≤𝑥_𝑘+1=𝜇𝐾𝑒^{−𝛼𝑥𝑘}𝑥_𝑘−1𝐾+(𝜇−1)𝑥_𝑘−1+𝑑≤𝑥_𝑘+1≤𝜇𝐾𝑒^−𝛼𝑥_𝑘𝜇−1+𝑑≤limsup_𝑘→∞𝑥_𝑘≤𝜇𝐾𝜇−1+𝑑<∞. (3.25) Thus, 0<𝑑≤liminf_𝑘→∞𝑥_𝑘≤limsup_𝑘→∞𝑥_𝑘≤𝜇𝐾𝜇−1+𝑑<∞. (3.26) The equilibrium point is positive and unique from Theorem 3.1(i) and within [𝑑,𝑑𝐾/(𝜇−1)]⊂[𝑑,𝜇𝐾/(𝜇−1)+𝑑] from Theorem 3.3(i). Assume that the equilibrium point is locally unstable. Then, the stability constraints (3.18) do not jointly hold. Since (3.26) holds and the equilibrium is locally unstable, a stable limit oscillatory solution has to exist subject to (3.26). Then, there is 𝜃∈𝐍 (the oscillation period) such that 0<𝑑≤lim_𝑘→∞𝑥_𝑘+𝑘1=𝑥_𝑘+𝑘1≤𝜇𝐾/(𝜇−1)+𝑑<∞, for all 𝑘₁∈𝜃. Now, build the real vector 𝑥_𝑘⁽⁰⁾∶=(𝑥_𝑘,𝑥_𝑘−1,…,𝑥_{𝑘−𝜃+1})^𝑇→𝑥_𝑘⁽⁰⁾ as 𝑘→∞. The linearization about the limit oscillation (limit cycle) is defined by the dynamics 𝑧⁽⁰⁾_𝑘+1=𝐴𝑧_𝑘⁽⁰⁾ where 𝐴∈𝐑^𝜃×𝜃 is block-partitioned as 𝐴∶=[𝜁(𝑥)𝜓^′(𝑥)𝜓(𝑥)𝜁^′(𝑥)0⋯00𝐼_𝜃−1]. It is direct to see that 𝐴 has an eigenvalue with modulus greater than unity if and only if (3.16) has a zero with modulus greater than unity since they coincide. Therefore, the limit cycle is unstable which leads to a contradiction to the fact that it should exist and be stable. Then, the unique equilibrium point is locally stable and there is no stable limit cycle. Since the equilibrium point is stable and unique, it is also a global stable attractor and Property (i) has been proved.
(ii) If 𝑑=0, then 𝑥=0 is an equilibrium point. For any small perturbation population at sampled time 𝑘, which has to be positive if 𝑝=0,𝑥_𝑗>0 for 𝑗∈𝐽⧵𝑘 for some finite integer 𝐽>𝑘 from inspection of the limit equation (1.2). Then, the zero equilibrium is locally unstable. An alternative proof is the test of the stability conditions (3.18) for 𝑝=𝑥=0. The second condition |𝜓(0)𝜁^′(0)|+𝜁(0)𝜓^′(0)=𝜁(0)𝜓^′(0)=𝛾𝜇>1 is violated so that 𝑥=0 is locally unstable. For 𝑥=0 and 𝑝>0, both stability conditions are fulfilled as 1+𝜁(0)𝜓^′(0)=1>0, |𝜓(0)𝜁^′(0)|+𝜁(0)𝜓^′(0)=0<1 so that the zero equilibrium point is locally stable. The largest equilibrium point is a positive global stable attractor from Property (i). If 𝑑≤0 and satisfies the conditions of Theorem 3.4, then there are two positive equilibrium points. The largest one is a stable global attractor proved under similar considerations as those used to prove Property (i).

4. Existence of Oscillatory Solutions of the Nonstationary-Modified Beverton-Holt Equation

The main motivation of the study of this section is that very commonly the evolution equation of populations in biology problems follows oscillatory seasonal patterns according to reproduction needs, food supply from the environment, or metabolic cycles of the species [13]. In Theorem 2.8, it has been proved that 0<𝑠_𝑘=^𝑘_𝑖=0𝜌_𝑖𝑠₀+^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝜂_𝑖=^𝑘_𝑖=0𝜌_𝑖𝑠₀+^𝑘−1_𝑖=0^𝑘−1_𝑗=𝑖+1𝜌_𝑗𝑣_𝑖−𝜔_𝑖≤∞,∀𝑘∈𝐍, (4.1) where 𝜌_𝑖∶=1/𝛾_𝑖𝜇_𝑖−1𝜁_𝑖,𝜔_𝑖∶=𝐾_𝑖−1𝑠_𝑖−1𝛿_𝑖𝑑_𝑖/𝜌_𝑖⁻¹(𝜌_𝑖⁻¹+𝛿_𝑖𝑑_𝑖) and 𝛿_𝑖∶=𝐾_𝑖−1𝑠_𝑖−1+𝜇_𝑖−1−1, which implies and it is implied by 𝑥_𝑘∈[0,∞); for all 𝑘∈𝐍. Furthermore, if 𝑑_𝑘>0; for all 𝑘∈𝐍 then the above upper-bound is finite which implies and it is implied by 𝑥_𝑘∈[0,∞); for all 𝑘∈𝐍. It is of interest to investigate when the solution of the modified Beverton-Holt equation is oscillatory equivalent to the solution of its inverse to be oscillatory. It has to be pointed out that the generation of limit oscillatory solutions being independent of the initial conditions is very important in some engineering problems, like, for instance, all those associated with the synthesis of electronic oscillators (see, e.g., [28]). To that end in the context of the modified Beverton-Holt equation (1.2), define the discrete functions 𝜒_𝑥∶𝐍²₀→𝐑 and 𝜒_𝑥∶𝐍₀×𝐍²→𝐑 as follows: 𝜒_𝑥(𝑗,𝑘)∶=𝑥_𝑘−𝑥_𝑗,𝜒_𝑥(𝑘,𝑘+𝑗,𝑘+𝑗+ℓ)∶=𝜒_𝑥(𝑘,𝑘+𝑗)𝜒_𝑥(𝑘,𝑘+𝑗+ℓ)=𝑥_𝑘+𝑗−𝑥_𝑘𝑥_{𝑘+𝑗+ℓ}−𝑥_𝑘, (4.2) respectively. Thus, precise definitions of what is meant by oscillatory solution of (1.2) are given.

Definition 4.1. A particular solution of (1.2) is strongly oscillatory if for any 𝑘∈𝐍₀, there exist two finite natural numbers 𝑁_𝑖(𝑘);𝑖=1,2, depending on 𝑘, such that 𝜒_𝑥𝑘,𝑘+𝑁₁(𝑘),𝑘+𝑁₁(𝑘)+𝑁₂(𝑘)<0. (4.3)

Definition 4.2. A particular solution of (1.2) is weakly oscillatory if for any 𝑘∈[𝑘₁,𝑘₂]=𝑘₂⧵𝑘₁−1⊂𝐍₀ and some 𝑘₁∈𝐍₀ and 𝑘₂(≥𝑘₁)∈𝐍₀, there exist two finite natural numbers 𝑁_𝑖(𝑘);𝑖=1,2, depending on 𝑘, such that 𝜒_𝑥𝑘,𝑘+𝑁₁(𝑘),𝑘+𝑁₁(𝑘)+𝑁₂(𝑘)<0. (4.4) It turns out that any oscillation of the Beverton-Holt equation implies, and is implied by, an oscillation of the solution of its inverse. Definitions 4.1 and 4.2 imply two alternate changes of sign at some three samples everywhere for the solution sequence (strong oscillation) or within some finite interval (weak oscillation), respectively. For each interval [𝑘,𝑘+𝑁₁(𝑘),𝑘+𝑁₁(𝑘)+𝑁₂(𝑘)], only one oscillation test is needed. If 𝑥_𝑘=0, then either the population is extinguished and remains at the zero equilibrium point or it is recovering positive values at the next samples generated by immigration so that the oscillation existence text can be performed later on formally.

Assertion 4.3. 𝜒_𝑥(𝑘,𝑘+𝑁₁(𝑘),𝑘+𝑁₁(𝑘)+𝑁₂(𝑘))<0 for 𝑘∈𝐍₀,𝑁_𝑖(𝑘)∈𝐍(𝑖=1,2)⇔𝜒_𝑠(𝑘,𝑘+𝑁₁(𝑘),𝑘+𝑁₁(𝑘)+𝑁₂(𝑘))<0 for 𝑘∈𝐍₀,𝑁_𝑖(𝑘)∈𝐍(𝑖=1,2).

Proof. It is direct from 𝑠_𝑘=𝑥_𝑘⁻¹, for all 𝑘∈𝐍₀ if 𝑥_𝑘>0. If 𝑥_𝑘=0, either the population is extinguished for all samples and no test for oscillation applies or some population is recovered from immigration (positive disturbance) at future samples and the inverse Beverton-Holt equation can be defined.

Assertion 4.4. The subsequent technical result holds: 𝜒_𝑥𝑘,𝑘+𝑁₁(𝑘),𝑘₁<0∧𝜒_𝑥𝑘₁(𝑘),𝑘₁(𝑘)+𝑁₁𝑘₁(𝑘),𝑘₁(𝑘)+𝑁₁𝑘₁(𝑘)+𝑁₂𝑘₁(𝑘)<0 (4.5) with 𝑘₁(𝑘)∶=𝑘+𝑁₁(𝑘)+𝑁₂(𝑘) for some 𝑘∈𝐍₀,𝑁_𝑖(𝑘)∈𝐍,𝑁_𝑖(𝑘₁)∈𝐍(𝑖=1,2)⟺∃𝑁_𝑖(𝑘+𝑗)∈𝐍(𝑖=1,2)∶𝜒_𝑥𝑘+𝑗,𝑘+𝑁₁(𝑘+𝑗)+𝑁₂(𝑘+𝑗)<0, (4.6) for all 𝑗(≤𝑁₁(𝑘)+𝑁₂(𝑘))∈𝐍 such that 𝑥_𝑘+𝑗>0,𝑁_𝑖(𝑘)∈𝐍.

Proof. Note that sign𝑥_{𝑘+𝑁1(𝑘)}−𝑥_𝑘≠sign𝑥_{𝑘+𝑘1(𝑘)}−𝑥_𝑘,sign𝑥_{𝑘+𝑁1(𝑘)}−𝑥_𝑘≠sign𝑥_{𝑘+𝑘1(𝑘)}−𝑥_{𝑘+𝑁1(𝑘)},sign𝑥_{𝑘1(𝑘)+𝑁1(𝑘1(𝑘))}−𝑥_𝑘1(𝑘)≠sign𝑥_{𝑘1(𝑘)+𝑁1(𝑘1(𝑘))+𝑁2(𝑘1(𝑘))}−𝑥_{𝑘1(𝑘)+𝑁1(𝑘1(𝑘))},sign𝑥_{𝑘1(𝑘)+𝑁1(𝑘1(𝑘))}−𝑥_𝑘1(𝑘)≠sign𝑥_{𝑘1(𝑘)+𝑁1(𝑘1(𝑘))+𝑁2(𝑘1(𝑘))}−𝑥_{𝑘1(𝑘)+𝑁1(𝑘1(𝑘))}. (4.7) Thus, it is obvious that there exist 𝑁_𝑖(𝑘+𝑗)∈𝐍(𝑖=1,2)∶𝜒_𝑥(𝑘+𝑗,𝑘+𝑁₁(𝑘+𝑗)+𝑁₂(𝑘+𝑗))<0, for all 𝑗∈𝐍 such that 𝑥_𝑘+𝑗>0, for all 𝑗(≤𝑁₁(𝑘)+𝑁₂(𝑘))∈𝐍, since it always exist by construction at least one 𝑁₁(𝑘+𝑗)≤𝑁₁(𝑘₁(𝑘)) and at least one 𝑁₂(𝑘+𝑗)≤𝑁₂(𝑘₁(𝑘)), for all 𝑗∈𝐍 such that 𝑥_𝑘+𝑗>0. If 𝑥_𝑘+𝑗=0, for all 𝑗≠𝑘,𝑘+𝑁₁(𝑘),𝑘+𝑁₁(𝑘)+𝑁₂(𝑘), then there are always at least three sign changes in the interval [𝑘,𝑘+𝑁₁(𝑘)+𝑁₂(𝑘)] and the proof also follows.

Assertions 4.3-4.4 imply directly the following result which is concerned with the simplification of the test for strong oscillatory solutions.

Assertion 4.5. Assume that there exists a subset of infinite cardinal of the set of natural numbers 𝑆_𝑗∶={𝑘_𝑗∈𝐍₀∶1≤𝑘_𝑗+1−𝑘_𝑗<∞}⊂𝐍₀ such that 𝜒_𝑥𝑘,𝑘+𝑁₁(𝑘),𝑘+𝑁₁(𝑘)+𝑁₂(𝑘)<0∀𝑘∈𝑆_𝑗forsome𝑁_𝑖⟺(𝑘)∈𝐍(𝑖=1,2)𝜒_𝑠𝑘,𝑘+𝑁₁(𝑘),𝑘+𝑁₁(𝑘)+𝑁₂(𝑘)<0for𝑘∈𝐍₀,𝑁_𝑖(𝑘)∈𝐍(𝑖=1,2) (4.8) and some solution {𝑥_𝑘}^∞₁ generated by initial conditions 𝑥₋₁≥0,𝑥₀>0. Then, such a solution is strongly oscillatory. The converse is also true.

The main result of this section follows.

Theorem 4.6. A solution of the modified Beverton-Holt equation (1.2), subject to initial conditions 𝑥₋₁≥0,𝑥₀>0, is strongly oscillatory if and only if there exists an infinite sequence of triples {𝑘_𝑖,𝑁₁(𝑘_𝑖),𝑁₂(𝑘_𝑖)}_𝑖∈𝐍0∈𝐍³, where 𝑘_𝑖+1=𝑘_𝑖+𝑁₁(𝑘_𝑖)+𝑁₂(𝑘_𝑖)∈𝑆⊂𝐍 for some finite 𝑘₀∈𝐍₀, for all 𝑖∈𝐍₀ with 𝐍∋𝑁_𝑗(𝑘_𝑖)<∞(𝑗=1,2) such that 𝑠_𝑘𝑖𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)}𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)}>0 (4.9) and, furthermore, 𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)}<𝑠_𝑘𝑖𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)}−𝑠_𝑘𝑖𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)}−𝑠_𝑘𝑖=𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)}−𝑠_𝑘𝑖𝜅_𝑖−1⟺𝑥_{𝑘𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)}>𝑥_𝑘𝑖−𝑥_{𝑘𝑖+𝑁1(𝑘𝑖)}𝑥_𝑘𝑖𝑥_𝑘𝑖−𝑥_{𝑘𝑖+𝑁1(𝑘𝑖)}⟺𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)}<𝑠_𝑘𝑖<𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)}∨𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)}>𝑠_𝑘𝑖>𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)}⟺𝑥_{𝑘𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)}>𝑥_𝑘𝑖>𝑥_{𝑘𝑖+𝑁1(𝑘𝑖)}∨𝑥_{𝑘𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)}<𝑥_𝑘𝑖<𝑥_{𝑘𝑖+𝑁1(𝑘𝑖)}, (4.10) where 𝜅_𝑖∶=𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)}/𝑠_𝑘𝑖.

Proof. Expand 𝜒_𝑠(𝑘_𝑖,𝑘_𝑖+𝑁₁(𝑘_𝑖),𝑘_𝑖+𝑁₁(𝑘_𝑖)+𝑁₂(𝑘_𝑖)) using Definition 4.1 to obtain 𝜒_𝑠𝑘_𝑖,𝑘_𝑖+𝑁₁𝑘_𝑖,𝑘_𝑖+𝑁₁𝑘_𝑖+𝑁₂𝑘_𝑖=𝑠_𝑘𝑖𝑠_𝑘𝑖−𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)}−𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)}+𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)}𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)}<0. (4.11)
Then, the proof follows directly from Assertions 4.4-(4.5) and direct computations since (4.11) is equivalent to the four equivalent propositions for the existence of strong oscillatory solution.

Remark 4.7. Note that the statement after the disjunction symbol in the third proposition in Theorem 4.6 for existence of strong oscillation may be expanded as follows: 1𝛾_{𝑘𝑖+𝑁1(𝑘𝑖)−1}𝜇_{𝑘𝑖+𝑁1(𝑘𝑖)−2}𝜁_{𝑘𝑖+𝑁1(𝑘𝑖)−2}^{𝑘𝑖+𝑁1(𝑘𝑖)−2}_{𝑘=𝑘𝑖}𝜌_𝑗+^{𝑘𝑖+𝑁1(𝑘𝑖)−1𝑘}_𝑗=0^{𝑖+𝑁1(𝑘𝑖)−1}_ℓ=𝑗+1𝜌_ℓ𝜂_𝑗𝑠_𝑘𝑖>1>𝑠_{𝑘𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)}=1𝛾_{𝑘𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)−1}𝜇_{𝑘𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)−2}𝜁_{𝑘𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)−2}=^{𝑘𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)−2}_{𝑘=𝑘𝑖}𝜌_𝑗+^{𝑘𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)−1𝑘}_𝑗=0^{𝑖+𝑁1(𝑘𝑖)+𝑁2(𝑘𝑖)−1}_ℓ=𝑗+1𝜌_ℓ𝜂_𝑗𝑠_𝑘𝑖 (4.12) since 𝜁(𝑥) varies from zero to infinity, it is always possible to have strong oscillatory solutions by accomplishing with the inequalities with appropriate parameterizations of combinations of the harvesting quota and intrinsic growth rate at certain samples sufficiently far away from each former test for oscillation (see also Assertion 4.4). Similar considerations apply for the reverse inequalities guaranteeing also strong oscillatory solutions conditions for weak oscillatory solutions are similar but they only apply on some finite interval.

5. Numerical Example

In this section, the evolution of a population of aphids is investigated under the modified generalized Beverton-Holt equation through a numerical tested example. A brief empirical description of those species follows. Aphids, or plant lice, are small, soft-bodied, pear-shaped insects which are commonly found on nearly all indoor and outdoor plants, as well as vegetables, field crops, and fruit trees. Most of them are about 1/10 inch long. They feed on plants by piercing them with syringe-like mouths pars and sucking the sap. Their diet is rich in carbohydrates and deficient in amino acids. Some of these amino acids cannot be synthesized by the insect but are supplied by the intracellular symbiont bacteria Buchnera aphidicola. Such a symbiont lives inside huge host cells (bacteriocytes) of which there are about 60–80 per individual and are transmitted to eggs through generations. Aphids have unusual and complex life cycles which allow them to build up huge levels of population in very short periods of time. Furthermore, they overwinter as fertilized eggs. Nymphs which hatch from these eggs become wingless females (stem mothers) which reproduce without mating holding new amounts of eggs in their bodies until hatching born alive youths. This pattern continues for as long as conditions are favorable. When the days get shorter in the fall and there are cooler temperatures, a generation appears which includes both males and females which are the production of fertilized eggs which overwinter [14]. As a short summary of the above description, it can be emphasized that the populations grow very fast in short periods of time during the same year and they also decrease fast to reach low survival thresholds but not up to extinction.

A common mathematical model for the population evolution of aphids is given in [13] as follows: 𝑎_𝑘+1=𝑓𝑟(1−𝑚)𝑎_𝑘=𝑓𝑟(1−𝑚)^𝑘+1𝑎₀, (5.1) where 𝑎_𝑘 is the number of adult female aphids in the 𝑘th generation, 𝑝_𝑘=𝑓𝑎_𝑘 is the number of progeny in the 𝑘th generation, 𝑚 is the fractional mortality of the young aphids, 𝑓 is the number of progeny per female aphid and 𝑟 is the ratio of female aphids to total adult aphids. The above parameters are all nonnegative by nature. The model might be generalized to extend the above parameters to be varying sequence. From standard stability results the model is globally asymptotically stable (and also exponentially stable) if and only if 𝑔∶=𝑓𝑟(1−𝑚)<1, globally stable with the population being constant if 𝑔=1 and unstable (with the solution diverging at exponential rate) if 𝑔>1. It turns out that such a description does not fix properly the above empirical description from biology knowledge. Exhaustive simulation via the modified generalized Beverton-Holt equation shows that the nontrivial solutions are always bounded and positive as expected. They have equilibrium points which depend on parameters. For instance, for zero harvesting quota, zero disturbance contribution, and very small or zero powers 𝑝 of the penalty term, the solutions reach an equilibrium which lies below that reached by the standard Beverton-Holt model, defined by the carrying capacity. However, there are several cycles of alternate fast increase/decrease of the populations during the transients. This model property is useful to describe the above biological behavior during one year of population evolution. If the parameters are readjusted, the solution changes accordingly while keeping a similar shape. For instance, the upper limit is smaller as the penalty term increases which also produces more abrupt cycles of increase/decrease of the population.

Figure 1 describes the population evolution for constant parameterizations 𝜇=2,𝐾=200,𝑑=0,𝑝=0.01,𝛼=0.002, and 𝛾=0.8. It is seen that the upper limit is smaller than the level of 200 corresponding to a standard Beverton-Hold equation, governed by 𝜇 and 𝐾, and also smaller than the level of 160, corresponding to a (unmodified) generalized Beverton-Holt equation, governed by 𝜇,𝐾,𝛾, and 𝑑=0. The reason is that a penalty term is introduced for the overpopulation of aphids. The penalty term translates into a certain level of population control since the population individuals compete in their habitat. The transient shows alternate cycles of increase/decrease of population levels. If 𝛼=0.001, with the remaining parameters being identical, then the upper limit increases as expected since both the standard and unmodified generalized Beverton-Holt equations are parameterized by a carrying capacity of 200. The population evolution is displayed in Figure 2 where the existence of the upper limit is obvious and a zoom for the first part of transient exhibiting the alternate cycles is given in Figure 3. Finally, Figure 4 shows the sample-to-sample ratio between the population evolution of the modified generalized Beverton-Holt equation to the standard one over a transient of 30 samples for zero harvesting quota.

Figure 1 

𝑑=0 and 𝛼=0.001.

Figure 2 

𝜇=2,𝐾=200,𝑑=0,𝑝=0.01,𝛼=0.002, and 𝛾=0.8.

Figure 3 

𝜇=2,𝐾=200,𝑑=0,𝑝=0.01,𝛼=0.001, and 𝛾=0.8.

Figure 4 

Relative comparison with the standard Beverton-Holt equation.

On the other hand, Figures 5 and 6 display the population evolution under moderate disturbance contribution consisting of joint independent consumption plus net migration fixed to 1/10 of the initial population with harvesting quotas of 20% and 60%, respectively. Those quotas can be interpreted in the plagues context as the use of pesticides to fight the plague. In the first case, the population grows slightly for each cycle, but exhaustive simulation for more samples demonstrates that the population remains bounded to smaller values than the carrying capacity, as expected. In the second case, it is shown that the levels of population decrease at each successive life cycle.

Figure 5 

𝜇=2,𝐾=200,𝑑=0,𝑝=0.01,𝛼=0.001,.

Figure 6 

𝛾=0.8.

6. Conclusions

This paper has been devoted to the study of stability and equilibrium point of a proposed modified generalized Beverton-Holt equation in ecology. Such an equation is governed, in general, by six parametrical sequences, namely,

(i)the intrinsic growth rate and the environment carrying capacity which define the standard Beverton-Holt equation,(ii)the harvesting quota and the disturbance contribution which, in addition, to the two former ones parameterize the generalized Beverton-Holt equation,(iii)a penalty term for the eventual overpopulation which includes, in general, the product of a potential term and a decreasing exponential one of the previous population. It has been proved that the proposed model is Lyapunov stable and possesses stable equilibrium points. It has also been proved that even under constant parameterizations the model can exhibit an oscillatory behavior to describe the alternative cycles of increasing/decreasing levels of population evolution. In this sense, it can be useful to correct the foreseen evolution through the standard Beverton-Holt equation which, under nonextinction conditions obtained from a carrying capacity greater than unity, leads to the population convergence to a finite set point limit defined by the environment carrying capacity. It can also be useful to introduce a correction of the behavior foreseen by its limiting equation, which is the Ricker model, which leads also to a limit set point being typically smaller than that associated with the standard Beveerton- Holt equation. The modified generalized Beverton-Holt equation is also shown to be useful to describe during its transient alternate cycles of the population levels for certain populations of insects since the solution is bounded and extinction-free under standard parameterization conditions so that it exhibits an oscillatory behavior within finite positive lower and upper bounds. In this context and in order to show its usefulness the proposed model has been applied to study a population of aphids which exhibit alternate cycles of increase and decrease of population and which state in latent levels of populations against very adverse conditions of their habitat.

Acknowledgments

The author is very grateful to MCYT by its partial support of this work through Grant no. DPI2006-00714 and to the Basque Government by its support of this work via Research Grants: Research Groups nos. IT-269-07, SAIOTEK SPE07UN04, and SAIOTEK SPED06UN10. He is also very grateful to the reviewers by their interesting comments which have helped to improve the content of the original manuscript.