On a Difference Equation with Exponentially Decreasing Nonlinearity

We establish a necessary and sufficient condition for global stability of the nonlinear discrete red blood cells survival model and demonstrate that local asymptotic stability implies global stability. Oscillation and solution bounds are investigated. We also show that, for different values of the parameters, the solution exhibits some time-varying dynamics, that is, if the system is moved in a direction away from stability by increasing the parameters , then it undergoes a series of bifurcations that leads to increasingly long periodic cycles and finally to deterministic chaos. We also study the chaotic behavior of the model with a constant positive perturbation and prove that, for large enough values of one of the parameters, the perturbed system is again stable.


Introduction
In 1950-1980 various biological and medical phenomena were for the first time described using differential equations, among them the blood cell production 1-3 .
The general approach was to describe the cell production dynamics as where N is the number of cells, p t, N is the cell production rate at time t under the condition that the population size is N, and d t, N is the mortality rate.It was a common belief that the mortality is proportional to the amount of blood cells in the circulation d t, N t μN t , where μ can be time dependent or not, and that higher per capita production rates correspond Discrete Dynamics in Nature and Society to smaller numbers of existing cells N.However, if there exists a constant equilibrium N, p and d are time independent and p N d N , while p N > d N for 0 < N < N, p N < d N for N > N, then the dynamics of 1.1 is very simple: all solutions with N 0 > 0 are positive and monotonically converge to the constant equilibrium.Real data suggest that such models poorly describe the oscillatory and chaotic behavior which frequently occurs in observations.As the next stage, it was suggested to introduce delay in the production term the maturation delay : in fact, cell production in the bone marrow takes some time before the cell is released into circulation.The resulting delay equations can demonstrate oscillatory and chaotic behavior; see, for example, 1, 2, 4 .
The delay differential equation was introduced by Wa żewska-Czy żewska and Lasota 3 to describe the survival of red blood cells in an animal.Here, as above, N t is the number of the mature red blood cells at time t, μ > 0 is the per capita death rate the probability of death for blood cells which currently circulate , and P > 0 and q > 0 define the red blood cell production function: P can be described as the limit production when the number of cells tends to zero, the decay of cell production for large cell number becomes faster with the growth of q, and delay τ is the time required to produce a red blood cell.Equation 1.2 is considered with a nonnegative initial function and a positive initial value which describe the number of cells in the past: N t ϕ t , ϕ t ≥ 0, −τ ≤ t < 0, N 0 N 0 > 0. 1.3 Evidently 1.2 with initial condition 1.3 has a unique solution which is positive for any t ≥ 0. Equation 1.2 has a positive equilibrium N which satisfies the equation After 1.2 was introduced in 3 , this equation and its nonautonomous modification where μ t > 0, P t > 0, q > 0, τ t ≥ 0, were intensively studied 4-8 , and some results were presented in the monograph 9 .Later on, several generalizations of 1.2 were investigated, for example, equations with impulses 10-13 and a distributed delay 14 .Discrete modifications of 1. where {x} x − x is the fractional part of the number, the same integration process as above leads to the delay difference equation Compared to ordinary and delay differential equations, difference equations have certain advantages, and their study can be motivated by the following reasons.
1 On the one hand, compared to ordinary differential equations, difference models have much richer and more complex dynamics, including oscillation, cycles, and transition to chaos.On the other hand, due to a simpler form, they allow a relatively simple analysis of stability, bifurcations, and some other dynamics-related issues.
Comparing two models, nondelay difference equation 1.10 and delay equation 1.12 , we notice that already 1.10 inherits some main properties of the dynamics of original delay differential equation 1.2 : existence of unstable positive solutions for some values of the parameters, sustainable oscillations, and chaotic behavior.
As mentioned above, there have already been some publications on the discrete Lasota-Wazewska model.However, the present paper is either different or completes the existing results from the following points of view.
1 The model considered in the present paper is relatively simple autonomous, without delay .However, for this equation we for the first time obtain strict solution bounds and sharp stability results, not just sufficient attractivity conditions.
2 The problems studied earlier for difference equations of this type were usually restricted to bounds, oscillation, and stability, and the analysis mainly referred to the range of parameters lower than those defining the first period-doubling bifurcation when the stable two-cycle appears.The present paper studies the whole range of parameters, including the chaotic areas.
3 For the first time, constant perturbation of the discrete Lasota-Wazewska model was introduced and the resulting period-halving bifurcations were investigated.
The paper is organized as follows.Section 2 presents known facts applied to 2.1 and contains some auxiliary results which will be applied later.In Section 3 we obtain lower and upper bounds for solutions of 2.1 , present a sharp nonoscillation condition, and demonstrate that the local asymptotic stability of the positive equilibrium of 2.1 is equivalent to its global asymptotic stability.In Section 4 we investigate bifurcations of the map described by 2.1 in two parameters: p and q.It is demonstrated that an addition of a small positive perturbation to the right hand side of 2.1 leads to stable behavior for q large enough but does not change the character of the bifurcation diagram in p.Finally, in Section 5 the results of the present paper are summarized and several open problems are outlined.

2.2
By the biological interpretation x n is the red blood cells count , we refer to positive solutions of 2.1 only; in particular, we assume that the initial condition is positive: The following stability definitions follow the monograph 18 .
Definition 2.1.The point x * is a critical point of the difference equation x is a fixed point or an equilibrium point if f x x.The equilibrium point is locally stable if for any ε > 0 there exists δ > 0 such that x for any x 0 in the domain including x; everywhere in this paper we assume that x 0 ∈ 0, ∞ .
Our purpose is to establish a new global stability condition for 2.1 which is both necessary and sufficient and prove that local stability implies global stability.
We will also use the following definition.
Definition 2.2.The difference equation is persistent if there exists m > 0 such that any solution with x 0 > 0 satisfies lim sup x n ≥ m.

2.7
In the recent decades the problems associated with understanding the dynamical behavior of nonlinear difference equations have been receiving intensive attention.The asymptotic behavior of solutions of 2. The stability of fixed points is one of the most important issues in the study of population dynamics models.In fact, it is much easier to perform local stability analysis of a fixed point than its global stability analysis.For local stability of 2.1 we see that, if where then the fixed point x is locally stable.As a special case of the results established in the above mentioned papers for global stability, we have the following global stability conditions for the unique fixed point x of 2.1 .We note that β in C7 is the unique solution of the equation h 2 x x in 0, x , where h x 1 − α −1 pe −qx .1.
In general, local and global stability conditions for difference equations do not coincide and different methods should be used to investigate the conditions for global stability.On the one hand, it is possible for a model to be locally stable and fail to be globally stable when the peak maximum is too high and when the curve falls off too quickly.On the other hand, it is possible for a model to be globally stable and fail to be locally stable if the nonlinear function is not continuous.In fact, there is no obvious connection between local and global stability conditions, and there is still a gap between the most local and global stability tests.Cull 28,29 considered this problem and proved that local stability implies global stability for some discrete population models.The results in 28, 29 are based on the enveloping principle, where it is proved that the enveloping by a linear fractional function is sufficient for global stability.However, this proof is not easy for complicated models and the parameter in the enveloping function must be adjusted for each particular population model.Moreover, the map should be a unimodal function in the sense of 28 , which is not applicable to the model 2.1 that we will study in this paper.
Further we will apply the following result on the equivalence of local and global stability; its exact statement is taken from 25, Proposition 7 , see also 30, 31 .is negative for x ∈ a, b , x / x * , then x is globally asymptotically stable.

Global Attractivity
In this section, we study qualitative properties of 2.1 and we establish a new global stability condition by applying Lemma 2.3.Before we state stability results, let us discuss the properties of the function We have f 0 p > 0. If pq ≤ 1 − α, then f x is increasing for any x > 0. If pq > 1 − α, then f x is decreasing for 0 < x < x * and f x is increasing for x > x * , where This means that 2.1 is not a population model in the sense of Cull's definition 28 , which describes positive unimodal functions increasing for x < x * and decreasing for x > x * .Thus the results obtained, for example, in 28, 29 cannot be applied to 2.1 .
then all positive solutions of 2.1 either exceed the equilibrium x for n ≥ 2 or are less than x, or coincide with x for any n ≥ 0. If pq ≤ 1 − α or 3.3 is satisfied, then all solutions converge monotonically to x.
Proof.Let pq ≤ 1 − α.Then f x is increasing for any x > 0, and all solutions converge monotonically to x thus, {x n } is increasing if x 0 < x, since x n < x implies x n < x n 1 < x, and is decreasing for x 0 > x .Further, let pq > 1 − α.Since f x is decreasing for 0 < x < x * and increasing for x > x * , then all positive solutions are oscillatory about the positive equilibrium whenever x < x * .In fact, if x < x * , then x n < x implies x n 1 > x, since f x > x for any x ∈ 0, x .There is also K > x such that f K x e.g., in Figure 2 b , K ≈ 3.9, where f K x ≈ 1.96 .If x n ≥ K then there exists j ∈ N such that either x n j < K, which implies x n j 1 < x, or all x n j x for j large enough.All solutions are nonoscillatory for x ≥ x * : they either exceed x beginning with the second iteration if at least one of the inequalities x 0 > x and f x 0 > x is satisfied or are less than x for any n, otherwise.The fact that x ≥ x * is equivalent to f x * ≥ x * , or completes the proof.
We note that in the nonoscillatory case we have monotone convergence to x; such monotonically convergent solutions for α 0.5, p 0.6, and q 2 are illustrated in Figure 1.The graph of the function 3.1 demonstrates that, unlike most unimodal functions in population dynamics models, f x decreases for small x and increases for large x, thus having a minimum, not a maximum.Here α 0.5, p 0.6, q 2, x * ≈ 0.438, x ≈ 0.47 > x * .The graph incorporates some iterations which converge monotonically, beginning with the second iteration to the positive equilibrium.1 with α 0.5, p 10, q 2 a and α 0.5, p 50, q 2 b .In the former case all positive solutions converge oscillatory to the positive equilibrium x ≈ 1.35, and in the latter case the solution converges to a stable two-cycle ≈ 2.86, 1.59.
Remark 3.2.Theorem 3.1 presents nonoscillation conditions for 2.1 .Moreover, the proof yields that these conditions are sharp: if they are not satisfied, all positive solutions which are not identically equal to x for n large enough are oscillatory about x.
Remark 3.3.In the proof of Theorem 3.1 we have applied the fact that will be used in future: if f x > 0 for x > 0, the map f has the only positive equilibrium x and is increasing on x * , ∞ , where x * < x, then all trajectories with x 0 > 0 are eventually monotone and converge to x.
Theorem 3.4.For any positive solution x n of 2.1 one has where f x and x are defined in 3.1 and 3.2 , respectively.
Proof.If either pq ≤ 1 − α or 3.3 is satisfied, then by Theorem 3.1 we have monotone convergence, so we will just consider the case when pq > 1 − α and 3.3 does not hold.Since x * is the absolute minimum point for f x , x ≥ 0, then all x n but probably x 0 satisfy x n ≥ f x * .We have x * > x, so the sequence x n with x n > x * is decreasing until f x n < x * .We note that, since f x is decreasing for x ∈ 0, x * , then max which completes the proof.
Remark 3.5.Let us note that in the case when we have a stable 2 k -cycle, sharper bounds than 3.6 can be obtained.Moreover, if q grows, then even for the deterministic chaos the attraction domain of the solution becomes narrow see Figure 5 a , for the bifurcation diagram .
then the positive equilibrium point x of 2.1 is globally stable, that is, the local stability of x implies its global stability.
Proof.We can omit the monotonic convergence case x * < x and assume that x ∈ 0, x * .Since f 0 p and f x is decreasing in 0, x * , then for any b ≥ max{p, x * } the function f maps the segment 0, b onto itself.The Schwarzian derivative of function f defined by 3.1 is negative, if x / x * ln pq/ 1 − α /q, so by Lemma 2.3 a locally asymptotically stable equilibrium x is also globally asymptotically stable.
To illustrate the main result in Theorem 3.6, we consider the following examples.

3.10
Here p 10, q 2, and α 0.5.The fixed point of 3.10 is given by x ≈ 1.3484.One can see that conditions C1 -C9 are not satisfied, so these results cannot be applied to establish stability of 3.10 .Since x n 50e −2x n .

3.12
Here p 50, q 2, α 1/2.We note that for the positive fixed point x ≈ 1.965 of 3.12  4 demonstrates that p 50 is in the range of parameter p where 2.1 has a stable 2-cycle.
Remark 3.9.We also illustrate Theorem 3.6 by applying the Coppel result 32 which claims that the single positive equilibrium x of f x is globally stable whenever the equation x f f x has the only positive solution x.For 3.10 , by Theorem 3.6 the positive equilibrium point is globally stable.This is illustrated by the graph in Figure 3 a , where the curve of y f f x intersects y x only once.By Theorem 3.6 for 3.12 the positive equilibrium point is unstable.Figure 3 b , illustrates that there are more than one intersection points really, three: the middle one corresponds to the equilibrium, while the left and the right points form a stable two-cycle of y f f x and y x, so x is unstable.For p 10 there is one intersection point stability while for p 50 there are three intersection points the two-cycle is stable rather than the positive equilibrium .

Bifurcations and Chaos
Stability properties of 2.1 can also be illustrated using the bifurcation diagram.The critical point p c at which the unique positive equilibrium x loses its stability f x −1 can be easily computed as p c αx exp qx , where x 2 − α qα .

4.1
For p < p c , all positive solutions converge to the equilibrium point x see Figure 4 a .As p increases beyond p c , there is a series of period-doubling bifurcations leading to a deterministic chaos.This means that for parameters close to the stable region, this will be a stable two-cycle and if the system is moved in a direction away from stability, by increasing the parameters, then the dynamics becomes more complex and the system undergoes a series of bifurcations leading to increasingly longer periodic cycles and finally deterministic chaos.There are many models which demonstrate similar behavior period-doubling bifurcations leading to deterministic chaos , for example, the Ricker model 33 and the modification of the Hassel map 34 4.3 The bifurcation diagram for 2.1 , with α 0.5, q 2, as p is growing a , and for the perturbed model 4.6 with α 0.5, q 2, and λ 0.2 b .For p < p c 0.75e 3 ≈ 15.06 all trajectories of 2.1 converge to the equilibrium point x.As p increases beyond p c , there is a series of period-doubling bifurcations leading to a deterministic chaos.We see that the two bifurcation diagrams are similar, with bifurcations for larger p in 4.6 compared to 2.1 ; chaotic behavior is observed in the perturbed equation for any p, as well as in the original one.
Since in ecological systems chaotic dynamics is difficult to observe, it was suggested 35, 36 that perturbations can be a reason for this phenomenon.In fact, the perturbed Ricker model does not experience chaotic behavior for r large enough; see 37 for some more details and relevant references.With the growth of r the period-doubling route to chaos, which is characteristic to 4.4 , will break down, giving rise to distinctive period-halving bifurcations and a stable two-cycle 35, 36 for r large enough and also for perturbation λ < 1, since for λ > 1 there is a stable equilibrium .Not all perturbed models demonstrate the same phenomenon for r large.For example, the perturbed logistic model is chaotic for large r, whether λ 0 or not.
Let us consider the perturbed version of model 2.1 x n 1 1 − α x n pe −qx n λ, n ≥ 0.

4.6
There are two parameters p and q in the model 4.6 , and we will study bifurcation in both of them.
If we consider bifurcations of 2.1 in q, then there is again transition to deterministic chaos through the series of period-doubling bifurcations see Figure 5 a .The perturbed model 4.6 demonstrates a similar picture, up to a certain value of q; further, we observe The bifurcation diagram for 2.1 , with α 0.5, p 2, as q is growing a and for the perturbed model 4.6 with α 0.5, p 2, and λ 0.2 b .We see that there is a deterministic chaos in 2.1 for q large enough, while the perturbed model experiences chaos and then as q grows the break of chaos through period-halving bifurcations.
period-doubling reversals leading to a stable equilibrium.This type of behavior can be rigorously justified.
Proof.We recall that the right-hand side of 4.6 g x 1−α x n pe −qx n λ has a minimum point at x * , where Thus for fixed p, α, and λ, there exists q 0 > 0 such that x * 1 q ln pq 1 − α < λ for q > q 0 ; 4.8 hence we also have f x * > λ > x * .This means that g x is monotone increasing on the interval x * , ∞ including the equilibrium.By Remark 3.3 the only positive equilibrium of 4.6 is globally asymptotically stable; moreover, there is a monotone convergence to this equilibrium.
As far as bifurcations in q are considered, the perturbed equation 4.6 behaves like the perturbed version of the Ricker 4.2 and some other 4.3 population dynamics models: it demonstrates period-halving bifurcations and stable behavior for large q.If we follow bifurcations in p, then chaos is observed for any large p, similar to the perturbed logistic model 4.5 .

Discussion and Open Problems
In the context of the model of red blood cells production, the mathematical properties of solutions can be interpreted in the following way.
Convergence to the equilibrium solution usually means healthy normal blood cell production, while sustainable oscillations correspond to a dynamic disease.In particular, if convergence to the equilibrium is monotone, then, in addition to normal blood cell production, parameter determination from data sampling is relatively easy compared to oscillatory cases.If the difference equation is considered as a numerical approximation of the continuous model, convergence to the equilibrium means that the numerical scheme is stable, in addition to the similar property of the continuous equation.
Deterministic chaos in the difference equation can mean one of the following two possibilities: first, blood cell production can be irregular which occurs in practice, but was usually modeled by the addition of noise, or more general stochastic terms, or by introducing randomized parameters in the model , or, second, the numerical discrete approximation can create irregularity.In the case of the chaotic behavior, the range of parameters which, as can be seen in Figures 4 and 5, can be separated from zero and bounded are indicative for comparison with the normal ideal blood cells production.A small positive constant perturbation can correspond to a drug administration or other external influence stimulating blood cells production at a constant level.
The results of the present paper can be summarized as follows.
1 The model studied in the present paper can be considered as a discretization of the Lasota-Wazewska equation 3 .For the difference equation, we have presented solution bounds and strict nonoscillation conditions.Moreover, we have demonstrated that local asymptotic stability of the unique positive equilibrium is equivalent to its global asymptotic stability.
2 Bifurcations were investigated in both parameters p and q involved in the equation.We also considered the linearly perturbed model.As p grows, the bifurcation diagram of the perturbed equation does not differ qualitatively from the nonperturbed model while demonstrating later bifurcations and transition to chaos than the original model .However, for large q the perturbed equation exhibits stable behavior.To the best of our knowledge, the chaotic behavior of 2.1 and 4.6 has never been studied before.
Finally, let us formulate some relevant open problems.
A Consider bifurcations and chaos in a one-parameter version of 4.6 when neither p nor q is fixed but there is a connecting relation p f q .For example, similar to Theorem 4.1, the global stability can be easily justified for p q k , k > 0 large enough and λ > 0 in 4.6 .Can a general theory be developed for some models similar to 2.1 which have a unique positive equilibrium and a unique critical point which is a minimum rather than a maximum and their linearly perturbed versions?Again, Theorem 4.1 demonstrates that, if the unique critical point tends to zero as the bifurcation parameter grows, then the perturbed model is stable for this parameter when it is large enough.Here we note that perturbed models with λ < 0 were investigated in 38 and encourage the study of models with negative constant perturbations which can describe immigration.B There are results for the analogue of 2.1 with periodic coefficients 27 .However, sharp oscillation and stability conditions are still an open problem.
C Strictly speaking, the Lasota-Wazewska equation 3 involves a delay, so its discretization should be 1 has been proposed as a research project by Kocić and Ladas 19, Project 4.6.1 .Since then several authors have expounded on various aspects of this model.For contributions, we refer the reader to the papers by Karakostas et al. 20 , Zheng et al.21 , Ma and Yu 22 , Meng and Yan 23 , Li and Cheng 17 , Ivanov 24 , Györi and Trofimchuk 25 , El-Morshedy and Liz 26 , Wang and Li 15 , and Saker 27 .

Lemma 2 . 3 .
Let the function g : a, b → a, b , where g ∈ C 3 a, b , have at most one critical point x * in a, b .If a unique fixed point x ∈ a, b of g is locally asymptotically stable and the Schwarzian derivative

Figure 1 :
Figure 1:The graph of the function 3.1 demonstrates that, unlike most unimodal functions in population dynamics models, f x decreases for small x and increases for large x, thus having a minimum, not a maximum.Here α 0.5, p 0.6, q 2, x * ≈ 0.438, x ≈ 0.47 > x * .The graph incorporates some iterations which converge monotonically, beginning with the second iteration to the positive equilibrium.

Figure 2 :
Figure2: Some iterations for 2.1 with α 0.5, p 10, q 2 a and α 0.5, p 50, q 2 b .In the former case all positive solutions converge oscillatory to the positive equilibrium x ≈ 1.35, and in the latter case the solution converges to a stable two-cycle ≈ 2.86, 1.59.

3 . 11 thenExample 3 . 8 .
by Theorem 3.6 the equilibrium x is globally stable.For illustration see Figure 2 a .Consider the model x n 1 1 2

Figure 3 :
Figure 3:The unique equilibrium point is stable if the equation has no two-cycle, or the curve y f f x intersects y x only once.Here f x is defined in 3.1 , α 0.5, q 2, where p 10 a and p 50 b .For p 10 there is one intersection point stability while for p 50 there are three intersection points the two-cycle is stable rather than the positive equilibrium .

Figure 4 :
Figure 4:The bifurcation diagram for 2.1 , with α 0.5, q 2, as p is growing a , and for the perturbed model 4.6 with α 0.5, q 2, and λ 0.2 b .For p < p c 0.75e 3 ≈ 15.06 all trajectories of 2.1 converge to the equilibrium point x.As p increases beyond p c , there is a series of period-doubling bifurcations leading to a deterministic chaos.We see that the two bifurcation diagrams are similar, with bifurcations for larger p in 4.6 compared to 2.1 ; chaotic behavior is observed in the perturbed equation for any p, as well as in the original one.

Figure 5 :
Figure 5:The bifurcation diagram for 2.1 , with α 0.5, p 2, as q is growing a and for the perturbed model 4.6 with α 0.5, p 2, and λ 0.2 b .We see that there is a deterministic chaos in 2.1 for q large enough, while the perturbed model experiences chaos and then as q grows the break of chaos through period-halving bifurcations.

x n 1 1 − 2 and its nondelay version x n 1 1 −E The difference equation x n 1 1 − α x n px n ln K x n 5. 4 can
α x n pe −qx n−k , n ≥ 0. 5.1 Prove that the local stability of the positive equilibrium of 5.1 implies its global stability.D As was mentioned in the beginning, the difference equation 2.1 was considered as a numerical discretization of the Lasota-Wazewska equation.A similar model for the Mackey-Glass equation was studied in 39 , where the right-hand side is a unimodal function.If we consider the model corresponding to the Nicholson blowflies equation x n 1 1 − α x n px n − k e −qx n−k 5.α x n px n e −qx n , 5.3 then the right-hand side of 5.3 , generally, has more than one critical point.However, it still has the only positive equilibrium, and it was demonstrated in 40 that local asymptotic stability implies global asymptotic stability of 5.3 .Establish necessary and sufficient oscillation and global stability conditions for 5.2 .Explore bifurcations and chaos for the perturbed version of 5.3 .be considered as a discretization of the Gompertz differential equation.Obtain conditions on α, p, and K which guarantee that all solutions with 0 < x 0 < M for some M > 0 remain positive.Deduce local stability conditions for the unique positive equilibrium and demonstrate that they imply global stability.Explore the delay version of this model x n 1 1 − α x n px n ln K x n − k 5.5 and deduce global stability results.

Table 1
One can see that the global stability conditions C1 -C9 are different from the local stability condition 2.8 .The natural question now is the following: can we obtain a new global stability condition for 2.1 which improves the above conditions?Our aim in this paper is to give an affirmative answer to this question and prove that the local stability condition 2.8 is also a global stability condition for 2.1 .This means that the local stability implies the global stability for 2.1 Table