2004-Fez conference on Differential Equations and Mechanics Electronic Journal of Differential Equations, Conference 11, 2004, pp. 167–173. ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu ftp ejde.math.txstate.edu (login: ftp)

STABILITY AND HOPF BIFURCATION IN A HAEMATOPOIETIC STEM CELLS MODEL

HAMAD TALIBI ALAOUI, RADOUANE YAFIA

Abstract. We consider the Haematopoietic Stem Cells (HSC) Model with one delay, studied by Mackey [4, 5] and Andersen and Mackey [1]. There are two possible stationary states in the model. One of them is trivial, the second E^{∗}(τ), depending on the delay, may be non-trivial . This paper investigates the stability of the non trivial state as well as the occurrence of the Hopf bifurcation depending on time delay. We prove the existence and uniqueness of a critical values τ_{0} and \overline{τ} of the delay such that E^{∗}(τ) is asymptotically stable for τ < τ_{0} and unstable for τ_{0} < τ < \overline{τ}. We show that E^{∗}(τ_{0}) is a Hopf bifurcation critical point for an approachable model.

1. Introduction

The population of haematopoietic stem cells (HSC) give rise to all of the different elements of the blood: the white blood cells, red blood cells, and platelets, which may be either actively proliferating or in a resting phase. After entering the proliferating phase, a cell is committed to undergo cell division at a fixed time τ later. The generation time τ is assumed to consist of four phases, G_{1} the pre-synthesis phase, S the DNA synthesis phase, G_{2} the post-synthesis phase and M the mitotic phase. Just after the division, both daughter cells go into the resting phase called G_{0}-phase. Once in this phase, they can either return to the proliferating phase and complete the cycle or die before ending the cycle.

The dynamics of the (HSC) are governed by the coupled differential delay equation (see [1, 4, 5, 6]):

\begin{matrix}\frac{dN}{dt} = −δN − β(N)N + 2e^{−γτ} β(N_{τ})N_{τ} \\ \frac{dP}{dt} = −γP + β(N)N − e^{−γτ} β(N_{τ})N_{τ} \end{matrix}  (1.1)

where β is a monotone decreasing function of N which has the explicit form of a Hill function,

β(N) = β_{0} \frac{θ^{n}}{θ^{n} + N^{n}} .  (1.2)

2000 Mathematics Subject Classification. 34K18.

Key words and phrases. Haematopoietic stem cells model; delayed differential equations;

Hopf bifurcation; periodic solutions.

\begin{matrix}c \\ ◯ \end{matrix}2004 Texas State University - San Marcos.

Published October 15, 2004.

The symbols in equation (1.1) have the following interpretation. N is the number of cells in non-proliferative phase, N_{τ} = N(t − τ), P the number of cycling proliferating cells, γ the rate of cells loss from proliferative phase, δ the rate of cells loss from non-proliferative phase, τ the time spent in the proliferative phase, β the feedback function, rate of recruitment from non-proliferative phase, β_{0} the maximum recruitment rate, and θ and n the control shape of the feedback function.

2. Stability without delay τ = 0

For τ = 0 the equation (1.1) reads to

\begin{matrix}\frac{dN}{dt} = −δN + β(N)N \\ \frac{dP}{dt} = −γP \end{matrix}  (2.1)

Theorem 2.1. Assume δ ∈ (0, β_{0}]. The system (2.1) has a positive equilibrium (N^{∗}, 0) = (β^{−1}(δ), 0) which is asymptotically stable. The trivial one (0, 0) is unstable.

Proof. The characteristic equation of the linearized equation of (2.1) around E^{∗} = (N^{∗}, 0), has two roots given by λ_{1} = −δ + α′(N^{∗}) and λ_{2} = −γ, where

α(N) = β(N)N  (2.2)

and α′(N) its derivative. Since β is a decreasing function, E^{∗} is asymptotically stable. For the trivial equilibrium, the roots of the characteristic equation of the linearized equation of (2.1) around (0, 0) are λ_{1} = −δ + α′(0) and λ_{2} = −γ. Since α′(0) = β_{0} > δ, (0, 0) is unstable. □

3. Stability for positive delay

Normalizing the delay τ by the time scaling t → \frac{t}{τ}, effecting the change of variables u(t) = N(tτ) and v(t) = P (tτ), the system (1.1) is transformed into

\begin{matrix}u̇(t) = τ[−δu(t) − α(u(t)) + 2e^{−γτ} α(u(t − 1))] \\ v̇(t) = τ[−γv(t) + α(u(t)) − e^{−γτ} α(u(t − 1))] \end{matrix}  (3.1)

where α is given by equation (2.2). Let (H0) δ < \frac{β_{0}}{2} and denote by \overline{τ} = \frac{1}{γ} ln (_{\frac{2}{1+ \frac{2δ}{β_{0}}}} ).

Note that (H0) implies that for each 0 < τ < \overline{τ}, α′(u^{∗}) < 0 and β_{0}(2e^{−γτ} − 1) > δ and system (3.1) has a unique positive equilibrium E^{∗}(τ) = (u^{∗}(τ), v^{∗}(τ)) with

u^{∗}(τ) = θ(\frac{β_{0}(2e^{−γτ} − 1) − δ}{δ} )^{1/n}, v^{∗}(τ) = \frac{δu^{∗}}{γ} (\frac{1 − e^{−γτ}}{2e^{−γτ} − 1})

and the characteristic equation of the linearized equation associated with (3.1) around E^{∗}(τ) is

W (λ, τ) = (λ + τγ)(λ − τa(τ) − τb(τ)e^{−λ}) = 0,  (3.2)

with a(τ) = −(δ + α′(u^{∗})) and b(τ) = 2e^{−γτ} α′(u^{∗}) and

α′(u^{∗}) = \frac{δ}{β_{0}(2e^{−γτ} − 1)^{2}} [β_{0}(1 − n)(2e^{−γτ} − 1) + nδ] .

Since τγ > 0, the stability of the positive equilibrium E^{∗}(τ) follows from the study of roots of the equation

Δ(λ, τ) = λ − τa(τ) − τb(τ)e^{−λ} = 0  (3.3)

corresponding to the characteristic equation associated to the first equation in (3.1). To obtain the switch of stability of E^{∗}(τ), one needs to find the imaginary root of equation (3.3). Let λ = iζ, then Δ(iζ, τ) = 0 if and only if

ζ = arccos ( − \frac{a(τ)}{b(τ)}) ∈ (0, π) for 0 ≤ |\frac{a(τ)}{b(τ)}| ≤ 1 and

τ\sqrt{b^{2}(τ) − a^{2}(τ)} = arccos(−\frac{a(τ)}{b(τ)}) for 0 ≤ |\frac{a(τ)}{b(τ)}| < 1.

(3.4)

Let (H1) a(τ) < 0 and |b(τ)| < −a(τ) for all τ > 0. (H2) τa(τ) < 1, and |a(τ)| < |b(τ)| for all τ > 0.

Theorem 3.1. Under assumption (H0), we have:

(1) The trivial equilibrium (0, 0) is unstable for 0 < τ < \overline{τ}. (2)

(i) If a and b satisfy (H1), then E^{∗}(τ) is asymptotically stable for 0 < τ < \overline{τ}.

(ii) If a and b satisfy (H2), n is sufficiently large and γ is close enough to 0, there exists a unique τ_{0} in ]0, \overline{τ}[ such that E^{∗}(τ) is asymptotically stable for τ ∈]0, τ_{0}[ and unstable for τ ∈ (τ_{0}, \overline{τ}).

Proof. (1) The characteristic equation of the linearized equation associated to (3.1) around (0, 0) is

λ + τ(δ + β_{0}) − 2τe^{−γτ} β_{0}e^{−λ} = 0  (3.5)

From (H0), we have β_{0}(2e^{−γτ} − 1) > δ, thus (3.5) has a real root which is positive. Then (0, 0) is unstable.

(2) part (i): Let λ = μ + iν be a root of equation Δ(λ, τ) = 0 for 0 < τ < \overline{τ}. We have

\begin{matrix}μ − τa(τ) − τb(τ)e^{−μ} cos(ν) = 0 \\ ν + τb(τ)e^{−μ} sin(ν) = 0 \end{matrix}  (3.6)

If there exists a root μ_{0} ≥ 0 of (3.6), then −a(τ) ≤ b(τ)e^{−μ_{0}} cos(ν). Since −1 ≤ cos(ν) ≤ 1 and 0 < e^{−μ_{0}} < 1 and b(τ) < 0 for 0 < τ < \overline{τ}, we have b(τ) ≤ a(τ), which contradicts the assumption (H1). So for all 0 < τ < \overline{τ}, the roots of the equation (3.3) have negative real parts, and therefore E^{∗}(τ) is asymptotically stable. □

For the proof of the stability in (2) part (ii), we need the following lemmas.

Lemma 3.2 (Hale 1993 [2]). All roots of the equation (z + c)e^{z} + d = 0, where c and d are real, have negative real parts if and only if: (i) c > −1, (ii) c + d > 0, and (iii) \sqrt{d^{2} − c^{2}} < ζ, where ζ is the root of ζ = −c tan ζ, 0 < ζ < π, if c ≠ 0 and ζ = \frac{π}{2} if c = 0.

Lemma 3.3. Under hypotheses (H0) and (H2), for n sufficiently large and γ close enough to 0, there exists a unique solution τ_{0} of the second equation of (3.4) in ]0, \overline{τ}[, such that iζ_{0} is a purely imaginary root of equation (3.3), with ζ_{0} = arccos(−\frac{a(τ_{0})}{b(τ_{0})}). Furthermore, the following inequalities hold

τ\sqrt{b^{2}(τ) − a^{2}(τ)} < arccos(−\frac{a(τ)}{b(τ)}) for τ ∈ (0, τ_{0})

τ\sqrt{b^{2}(τ) − a^{2}(τ)} > arccos(−\frac{a(τ)}{b(τ)}) for τ ∈ (τ_{0}, \overline{τ})

(3.7)

Lemma 3.4. Let f : (0, π) → ℝ be defined by f(x) = α tan x, α < 1 and α ≠ 0. Then, f has a unique fixed point ζ ∈ (0, π), such that: For 0 < α < 1, f(x) < x if x ∈ (0, ζ) ∪ (\frac{π}{2} , π) and f(x) > x if x ∈ (ζ, \frac{π}{2} ); and for α < 0, f(x) < x if x ∈ (0, \frac{π}{2} ) ∪ (ζ, π) and f(x) > x if x ∈ (\frac{π}{2} , ζ).

Proof of (2) part (ii) of theorem 3.1. We only have to verify the three conditions (i), (ii) and (iii) of lemma 3.2. The assertions (i) and (ii) follow from (H2) with c = −τa(τ) and d = −τb(τ).

For condition (iii), let τ ∈ (0, τ_{0}) and f(ζ) = τa(τ) tan ζ. From the first equation of (3.7) we have: If a(τ) = 0, the first inequality of (3.7) becomes −τb(τ) < \frac{π}{2} , and (iii) is satisfied. If 0 < τa(τ) < 1 or a(τ) < 0, since

f( arccos(− \frac{a(τ)}{b(τ)})) = τ\sqrt{b(τ)^{2} − a(τ)^{2}},

the first equation of (3.7) implies that

f( arccos(− \frac{a(τ)}{b(τ)})) < arccos(− \frac{a(τ)}{b(τ)}),

with arccos(−\frac{a(τ)}{b(τ)}) ∈ (0, π). From lemma 3.4 and the graph of f, if ζ is the fixed point of f in (0, π), we have,

f( arccos(− \frac{a(τ)}{b(τ)})) < ζ,  (3.8)

that is \sqrt{(τb(τ))^{2} − (τa(τ))^{2}} < ζ, which leads to the desired assertion. This complete the stability of E^{∗}(τ) for 0 < τ < τ_{0}.

To prove the unstability of E^{∗}(τ) in (2) part (ii), for τ_{0} < τ < \overline{τ}, we will show that the characteristic equation (3.3) has at least one root with positive real part. Let τ_{0} < τ < \overline{τ}. If all the roots of the characteristic equation (3.3) have negative real parts, the properties (i), (ii) and (iii) of lemma 3.2 are satisfied. From the second equation of (3.7) and from (3.8) we have

f( arccos(− \frac{a(τ)}{b(τ)})) > arccos(−\frac{a(τ)}{b(τ)})

f( arccos(− \frac{a(τ)}{b(τ)})) < \overline{ζ}

Henceforth, from lemma 3.4 and the graph of f, we have

arccos(− \frac{a(τ)}{b(τ)}) < \overline{ζ}, and arccos(− \frac{a(τ)}{b(τ)}) > \overline{ζ}

which is impossible.

Now, suppose that there is one root with zero real part with all the remaining roots having negative real parts. From (3.4) and lemma 3.3 we deduce that τ = τ_{0}, which contradicts the assumption τ > τ_{0}. Then E^{∗}(τ) is unstable for τ_{0} < τ < \overline{τ}

□

Proof of Lemma 3.3. In view of (H0) and (H2), to find a root of second equation of (3.4) is equivalent to find a root of the equation

τ = −\frac{arccos(−\frac{a(τ)}{b(τ)})}{b(τ) sin(arccos(−\frac{a(τ)}{b(τ)}))}.  (3.9)

Let y(τ) = arccos(−\frac{a(τ)}{b(τ)}), and F (τ) = −\frac{y(τ)}{b(τ) sin(y(τ))}. Besides, in the hypotheses (H0) and (H2), F is continously differentiable on τ_{0} ∈ [0, \overline{τ}]. As F (0) > 0, for sufficiently large n and F (\overline{τ}) < \overline{τ} for γ close enough to 0, then there exits at least one solution τ_{0} of equation (3.9) in ]0, \overline{τ}[. Now, for the uniqueness of τ_{0}, let g(τ) = τ − F (τ), then

g′(τ) = 1 − \frac{y′(τ)b(τ) sin(y(τ)) − y(τ)b′(τ) sin(y(τ))}{(b(τ) sin(y(τ)))^{2}} − \frac{y(τ)b(τ) cos((τ))y′(τ)}{(b(τ) sin(y(τ)))^{2}}

where

y′(τ) = −\sqrt{1 − (\frac{a(τ)}{b(τ)})^{2}} \frac{a′(τ)b(τ) − a(τ)b′(τ)}{b^{2}(τ)}.

Since lim_{γ→0} \frac{d}{dτ} α′(u^{∗}) = 0, from (3.2), we have

lim_{γ→0}b′(τ) = 0 and lim_{γ→0}a′(τ) = 0.

Then lim_{γ→0}g′(τ) = 1 > 0, for 0 ≤ τ ≤ \overline{τ}. Since g′ > 0 and g is an increasing function on the interval ]0, \overline{τ}[ for γ close enough to 0, τ_{0} is unique in ]0, \overline{τ}[. By the continuity property of F, we have F (τ) > τ for τ ∈]0, τ_{0}[ and F (τ) < τ for τ ∈]τ_{0}, \overline{τ}[. □

4. Hopf Bifurcation Occurrence

Below, we will show that the following system has a Hopf bifurcation at τ = τ_{0},

\begin{matrix}\frac{dN}{dt} = −δN − β(N)N + 2e^{−γτ_{0}}β(N_{τ})N_{τ} \\ \frac{dP}{dt} = −γP + β(N)N − e^{−γτ_{0}}β(N_{τ})N_{τ} \end{matrix}  (4.1)

This system is equivalent to

\begin{matrix}u̇(t) = τ[−δu(t) − α(u(t)) + 2e^{−γτ_{0}}α(u(t − 1))] \\ v̇(t) = τ[−γv(t) + α(u(t)) − e^{−γτ_{0}}α(u(t − 1))] \end{matrix}  (4.2)

with u(t) = N(tτ) and v(t) = P (tτ). System (4.2) has a unique positive equilibrium E^{∗} = (u^{∗}, v^{∗}) = (u^{∗}(τ_{0}), v^{∗}(τ_{0})), for all τ > 0.

By the translation z(t) = (u(t), v(t)) − (u^{∗}, v^{∗}), system (4.2) is written as a functional differential equation (FDE) in C := C([−1, 0], ℝ^{2}):

ż(t) = L(τ)z_{t} + f_{0}(z_{t}, τ)  (4.3)

where L(τ) : C → ℝ^{2} is a linear operator and f_{0} : C ×ℝ → ℝ^{2} are given respectively by

L(τ)φ = τ (\begin{matrix}−(δ + α′(u^{∗}))φ_{1}(0) + 2e^{−γτ_{0}}α′(u^{∗})φ_{1}(−1) \\ −γφ_{2}(0) + α′(u^{∗})φ_{1}(0) − e^{−γτ_{0}}α′(u^{∗})φ_{1}(−1) \end{matrix})

f_{0}(φ, τ) = τ (\begin{matrix}−α(φ_{1}(0) + u^{∗}) + α′(u^{∗})φ_{1}(0) − 2e^{−γτ_{0}}α′(u^{∗})φ_{1}(−1) \\ +2e^{−γτ_{0}}α(φ_{1}(−1) + u^{∗}) − δu^{∗} \\ α(φ_{1}(0) + u^{∗}) − α′(u^{∗})φ_{1}(0) − e^{−γτ_{0}}α(φ_{1}(−1) + u^{∗}) \\ +e^{−γτ_{0}}α′(u^{∗})φ_{1}(−1) − γv^{∗}. \end{matrix})

for φ = (φ_{1}, φ_{2}) ∈ C.

Now, we apply the Hopf bifurcation theorem, see [2], to show the existence of a non-trivial periodic solution to (4.2) bifurcating from the non trivial equilibrium E^{∗}. We use the delay as a parameter of bifurcation. Therefore, the periodicity is a result of changing the type of stability, from stationary solution to limit cycle. Let

(H3) a(τ_{0}) < \frac{1}{\overline{τ}} and |a(τ)| < |b(τ)| , for 0 < τ < \overline{τ}.

Theorem 4.1. Under hypotheses (H0) and (H3) if n is sufficiently large and γ is close enough to 0, then, for τ ∈]0, τ_{0}[, E^{∗} is asymptotically stable; it is unstable for τ ∈]τ_{0}, \overline{τ}[, where τ_{0} is stated in lemma 3.3.

The proof of the above theorem follows the same procedure as that the proof of theorem 3.1 (2) (ii). Therefore, we omit it.

Theorem 4.2. Assume (H0) and (H3) hold, n is sufficiently large and γ is sufficiently small. There exists ε_{0} > 0 such that, for each 0 ≤ ε < ε_{0}, equation (4.2) has a family of periodic solutions p(ε) with period T = T (ε), for the parameter values τ = τ(ε) such that p(0) = E^{∗}, T (0) = \frac{2π}{ζ_{0}} and τ(0) = τ_{0}, where τ_{0} stated in lemma

3.3 and ζ_{0} = arccos ( − \frac{a(τ_{0})}{b(τ_{0})}).

Proof. We apply the Hopf bifurcation theorem introduced in [2]. From the expression of f_{0} in (4.3), we have

f_{0}(0, τ) = 0  and  \frac{∂f_{0}(0, τ)}{∂φ} = 0, for all τ > 0

The linearized equation associated to (4.2) around E^{∗} has the following characteristic equation:

Δ_{0}(λ, τ) = λ − τa(τ_{0}) − τb(τ_{0})e^{−λ} = 0,  (4.4)

Firstly, let λ = iζ. From (3.4) and lemma 3.3, we have

Δ_{0}(iζ, τ) = 0 ⟺ ζ_{0} = arccos ( − \frac{a(τ_{0})}{b(τ_{0})}) and τ = τ_{0}

where τ_{0} is unique in (0, \overline{τ}). Thus, the characteristic equation (4.4) has a pair of simple imaginary roots λ_{0} = iζ_{0} and \overline{λ}_{0} = −iζ_{0} at τ = τ_{0}.

Lastly, we need to verify the transversality condition. From (4.4), Δ_{0}(λ_{0}, τ_{0}) = 0 and \frac{∂}{∂λ}Δ_{0}(λ_{0}, τ_{0}) = 1 − τ_{0}a(τ_{0}) + λ_{0} ≠ 0. According to the implicit function theorem, there exists a complex function λ = λ(τ) defined in a neighborhood of τ_{0}, such that λ(τ_{0}) = λ_{0} and Δ_{0}(λ(τ), τ) = 0 and

λ′(τ) = −\frac{∂Δ_{0}(λ, τ)/∂τ}{∂Δ_{0}(λ, τ)/∂λ},  (4.5)

for τ in a neighborhood of τ_{0}. Let λ(τ) = p(τ) + iq(τ). From (4.5) we have

p′(τ)_{/τ=τ_{0}} = \frac{τ_{0}(b^{2}(τ_{0}) − a^{2}(τ_{0}))}{(1 + τ_{0}b(τ_{0}) cos ζ_{0})^{2} + (τ_{0}b(τ_{0}) sin ζ_{0})^{2}}

From (H3), we conclude that p′(τ)_{/τ=τ_{0}} > 0. □

5. Discussions

It’s known (Mackey (1997) [5]) that when taking γ as a bifurcation parameter and allowing γ to increase, a supercritical Hopf bifurcation of (1.1) is followed by an inverse Hopf bifurcation. Considering the delay τ as a parameter of bifurcation makes the study of bifurcation more complicated. In [1] the following conditions of stability of the non-trivial steady state of (1.1)

were proposed (Hayes (1950) [3]) |\frac{a(τ)}{b(τ)}| > 1 or |\frac{a(τ)}{b(τ)}| ≤ 1 and τ < \frac{arccos(− \frac{a(τ)}{b(τ)} )}{\sqrt{b(τ)^{2}−a(τ)^{2}}}

where 0 < τ < \frac{1}{γ} ln(_{\frac{2}{1+ \frac{δ}{β_{0}}}} ), δ < β_{0}.

In sections 2 and 3 of this paper it’s shown that if the loss rate γ from proliferating cells is smaller and the control shape n is large, then the steady state E^{∗}(τ) may be stable for τ = 0 and hence it’s stable for 0 < τ < τ_{0} and unstable for τ_{0} < τ < \overline{τ}, where \overline{τ} = \frac{1}{γ} ln(_{\frac{2}{1+ \frac{2δ}{β_{0}}}} ), 2δ < β_{0}. But at τ = τ_{0} we cannot give any result of stability

of E^{∗}(τ_{0}), because the dependance of E^{∗}(τ) on the delay τ, which makes the study of the Hopf bifurcation more difficult.

In the rest of the paper to study the Hopf bifurcation around the critical value τ = τ_{0}, we propose the approchable model (4.1) of (1.1). Then E^{∗}(τ_{0}) is the unique non-trivial steady state of (4.1) for all 0 < τ < \overline{τ}, which is stable for 0 < τ < τ_{0} and unstable for τ_{0} < τ < \overline{τ} and the Hopf bifurcation occures at τ = τ_{0}.

The results proposed in this paper should hopefully improve the understanding of the qualitative properties of the description delivered by model (1.1). So far we have now a description of stability properties and Hopf bifurcation with a detailed analysis of the influence of delays terms.

References

[1] L. K. Andersen and M. C. Mackey, Resonance in Periodic Chemotherapy: A Case Study of Acute Myelogenous Leukemia. J. theor. Biol. (2001) 209, 113-130.

[2] J. K.Hale and S. M. Verduyn Lunel, Introduction to functional Differential equations. Springer-Verlag, New-York, (1993).

[3] N. D. Hayes, Roots of the transcendental equation associated with a certain differencedifferential equation. J. London Math. Soc. 25 226-232 (1950).

[4] M. C. Mackey, Unified Hypothesis for the Origin of Aplastic Anemia and Periodic Hematopoiesis. Blood 51, 5 (1978).

[5] M. C. Mackey, Mathematical Models of Haematopoietic Cell Replication and Control. In: The Art of Mathematical Modelling: Case Studies in Ecology, Physiology and Biofluids (Othmer, H. G., Adler, F. R., Lewis, M. A. and Dallon, J. C., eds), pp. 149-178. New York: Prentice-Hall (1997).

[6] M. C. Mackey, Cell Kenitec Status of Haematopoietic Stem Cells. cell prolif. (2001), 34, 71-83.

Hamad Talibi Alaoui Université Chouaib Doukkali Faculté des Sciences, Département de Mathématiques et Informatique, B.P. 20, El Jadida, Morocco

   E-mail address: talibi@math.net

Radouane Yafia Université Chouaib Doukkali Faculté des Sciences, Département de Mathématiques et Informatique, B. P. 20, El Jadida, Morocco

   E-mail address: yafia radouane@hotmail.com