2005-Oujda International Conference on Nonlinear Analysis. Electronic Journal of Differential Equations, Conference 14, 2006, pp. 241–248. ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu ftp ejde.math.txstate.edu (login: ftp)

PERIODIC SOLUTIONS FOR SMALL AND LARGE DELAYS IN A TUMOR-IMMUNE SYSTEM MODEL

RADOUANE YAFIA

Abstract. In this paper we study the Hopf bifurcation for the tumor-immune system model with one delay. This model is governed by a system of two differential equations with one delay. We show that the system may have periodic solutions for small and large delays for some critical value of the delay parameter via Hopf bifurcation theorem bifurcating from the non trivial steady state.

1. Introduction

In this paper, we consider a model that provides a description of tumor cells in competition with the immune system. This description is described by many authors, using ordinary and delayed differential equations to model the competition between immune system and tumor. In particular [26, 31, 32] other similar models can be found in the literature, see, [23, 35, 39] provide a description of the modelling, analysis and control of tumor immune system interaction.

Other authors use kinetic equations to model the competition between immune system and tumor. Although they give a complex description in comparison with other simplest models, they are, for example, needed to model the differences of virulence between viruses, see, [1, 2, 5, 6, 7, 15]. Several other fields of biology use kinetic equations, for instance [19] and [20] give a kinetic approach to describe population dynamics, [2] deals with the development of suitable general mathematical structures including a large variety of Boltzmann type models.

The reader interested in a more complete bibliography about the evolution of a cell, and the pertinent role that have cellular phenomena to direct the body towards the recovery or towards the illness, is addressed to [22, 27]. A detailed description of virus, antivirus, body dynamics can be found in the following references [10, 21, 34, 36]. The mathematical model with which we are dealing, was proposed in a recent paper by M. Galach [26]. In this paper the author developed a new simple model with one delay of tumor immune system competition, this idea is inspired from the paper of Kuznetsov and Taylor (1994) [32] and he recall some numerical results of Kuznetsov and Taylor in order to compare them with those obtained in his paper, see, [26].

2000 Mathematics Subject Classification. 34K18.

Key words and phrases. Tumor-Immune system competition; delayed differential equations;

stability; Hopf bifurcation; periodic solutions.

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

Published September 20, 2006.

2. Mathematical model

The mathematical model describing the tumor immune system competition is given by a system of two differential equations with one delay, see Galach [26],

\begin{matrix}\frac{dx}{dt} = σ + ωx(t − τ)y(t − τ) − δx \\ \frac{dy}{dt} = αy(1 − βy) − xy \end{matrix}  (2.1)

where the parameter ω describes the immune response to the appearance of the tumor cells and the constant τ is the time delay which the immune system needs to develop a suitable response after the recognition of non-self cells. Time delays in connection with the tumor growth also appear in Bodnar and Foryś [12] and [13], Byrne [14], Foryś and Kolev [24] and Foryś and Maciniak-Czochra [25]. For the meaning of the parameters α, β, δ and σ, see Kuznetsov and Taylor [32] and Kirschner and Panetta [31].

For τ = 0 system (2.1) becomes a system of ordinary differential equations:

\begin{matrix}\frac{dx}{dt} = σ + ωxy − δx \\ \frac{dy}{dt} = αy(1 − βy) − xy \end{matrix}  (2.2)

In [26], the author study the existence, uniqueness and nonnegativity of solutions and he show the nonexistence of nonnegative periodic solution of system (2.2), using the Dulac-Bendixon criteria, see [37]. The possible nonnegative steady states of system (2.2) and their stability are summarized in the Table 1; see also [26],

Table 1. Nonnegative steady states of system (2.2) and their stability

Region Conditions  P_{0}  P_{1}  P_{2}
1  ω > 0, αδ < σ  stable
2  ω > 0, αδ > σ  unstable  stable
3  {ω < 0, αδ > σ, \atop α(βδ − ω)^{2} + 4βωσ > 0}  unstable  stable
4  {ω < 0, αδ < σ, ω + βδ < 0, \atop α(βδ − ω)^{2} + 4βωσ > 0} stable unstable stable
5  ω < 0, α(βδ − ω)^{2} + 4βωσ > 0 stable

For τ > 0, the existence and uniqueness of solutions of system (2.1) for every t > 0 are established in [26], using the results presented in Hale [29]. Based on the results of Bodnar [11], in [26] the author showed that: (1) If ω ≥ 0, these solutions are nonnegative for any nonnegative initial conditions (biologically realistic case).

(2) If ω < 0, there exist nonnegative initial condition such that the solution becomes negative in a finite time interval.

Our goal in this paper is to consider the case (1) when ω > 0, which is the most biologically meaningful one. We study the asymptotic behavior of the possible steady states P_{0} and P_{2} with respect to the delay τ. We establish that, the Hopf bifurcation may occur by using the delay as a parameter of bifurcation. We prove this result for small and large delays.

This paper is organized as follows. In section 3, we recall some results about the absolute and conditional stability of delay equations and the zeros of second order transcendental polynomials. In section 4, we investigate the results presented in section 3 to prove the stability of the possible steady states (trivial and non-trivial) of the delayed system (2.1). The main result of this paper is given in section 5 and section 6. Based on the Hopf bifurcation theorem, we show the occurrence of Hopf bifurcation for small and large delays.

3. Stability of delay equations and zeros of second order transcendental polynomials

In this section we recall some results on the stability of delay equations and on the zero of second order transcendental polynomials.

3.1. Absolute and conditional stability. Consider the following general nonlinear delay differential system

\frac{dx}{dt} = f(x(t), x(t − τ)),  (3.1)

where x ∈ ℝ^{n}, τ is constant, f : ℝ^{n} × ℂ^{n} → ℝ^{n} is smooth enough to guarantee the existence and uniqueness of solutions of (3.1) under the initial condition

x(θ) = φ(θ), θ ∈ [−τ, 0],  (3.2)

where C = C([−τ, 0], ℝ^{n}). Suppose f(x^{∗}, x^{∗}) = 0, that is x = x^{∗} is a steady state of system (3.1).

Definition 3.1. The steady state x = x^{∗} of system (3.1) is called absolutely stable (i.e., asymptotically stable independent of the delay τ) if it is asymptotically stable for all delays τ > 0. x = x^{∗} is called conditionally stable (i.e., asymptotically stable depending of the delay τ) if it is asymptotically stable for τ in some interval, but not necessarily for all delays τ > 0.

The linearized system of (3.1) at x = x^{∗} has the form

\frac{dX}{dt} = A_{0}X + A_{1}X(t − τ),  (3.3)

where X ∈ ℝ^{n}, A_{i} (i = 0, 1) is an n × n constant matrix. Then the characteristic equation associated with system (3.3) takes the form

det[λI − A_{0} − A_{1}e^{−λτ}].  (3.4)

The location of the roots of some transcendental equation (3.4) in it is general form has been studied by many authors, see Baptistini and Táboas [3], Bellman and Cooke [4], Boese [8], Brauer [9], Cooke and van den Driessche [18], Cooke and Grossman [17], Huang [30], Mahaffy [33], Ruan and Wei [38] an dthe references therein. The following result, which was proved by Chin [16], gives necessary and sufficient conditions for the absolute stability of system (3.3).

Lemma 3.1. The system (3.3) is absolutely stable if and only if (i) Re λI − A_{0} − A_{1}) < 0

(ii) det[iζ − A_{0} − A_{1}e^{−iζτ}] ≠ 0 for all ζ > 0.

3.2. A second degree transcendental polynomial. In this section, we state some results on the second degree transcendental polynomial (see Ruan and Wei [38]). For most system with discrete delay, the characteristic equation of the linearized system at a steady state is a second degree transcendental polynomial equation of the following form:

λ^{2} + pλ + r + (sλ + q)e^{−λτ} = 0  (3.5)

where p, r, q and s are real numbers. It is known that the steady state is asymptotically stable if all roots of the characteristic equation (3.5) have negative real parts.

Let define the following hypotheses:

(H1) p + s > 0.

(H2) q + r > 0.

(H3) (s^{2} − p^{2} + 2r < 0 and r^{2} − q^{2} > 0) or (s^{2} − p^{2} + 2r)^{2} < 4(r^{2} − q^{2}).

(H4) r^{2} − q^{2} < 0 or (s^{2} − p^{2} + 2r > 0 and (s^{2} − p^{2} + 2r)^{2} = 4(r^{2} − q^{2})).

(H5) r^{2} − q^{2} > 0, s^{2} − p^{2} + 2r > 0 and (s^{2} − p^{2} + 2r)^{2} > 4(r^{2} − q^{2}).

Theorem 3.1 ([38]). Let τ_{j}^{±} (j = 0, 1, 2, …) defined by

τ_{j}^{±} = \frac{1}{ζ_{±}} arccos \{\frac{q(ζ_{±}^{2} − r) − psζ_{±}^{2}}{s^{2}ζ_{±}^{2} + q^{2}} + \frac{2jπ}{ζ_{±}} \}, j = 0, 1, 2, …

where ζ_{±} is given by

ζ_{±} = \frac{1}{2}(s^{2} − p^{2} + 2r) ± \frac{1}{2}[(s^{2} − p^{2} + 2r)^{2} − 4(r^{2} − q^{2})]^{\frac{1}{2}} .

(i) If (H1)-(H3) hold, then all roots of equation (3.5) have negative real parts for all τ ≥ 0. (ii) If (H1), (H2) and (H4) hold, then when τ ∈ [0, τ_{0}^{+}) all roots of equation (3.5) have negative real parts, when τ = τ_{0}^{+} equation (3.5) has a pair of purely imaginary roots ±iζ_{±}, and when τ > τ_{0}^{+} equation (3.5) has at least one root with positive real part. (iii) If (H1), (H2) and (H5) hold, then there is a positive integer k such that there are k switches from stability to instability to stability; that is, when τ ∈ [0, τ_{0}^{+}), (τ_{0}^{−}, τ_{0}^{+}), …, (τ_{k−1}^{−}, τ_{k}^{+}), all roots of equation (3.5) have negative real parts, and when τ ∈ (τ_{0}^{+}, τ_{0}^{−}), (τ_{1}^{+}, τ_{1}^{−}), …, (τ_{k−1}^{+}, τ_{k−1}^{−}), and τ > τ_{k}^{+}, equation (3.5) has at least one root with positive real part.

Remark 3.1. Theorem 3.1 was obtained by Cooke and Grossman [17] in analyzing a general second order equation with delayed friction and delayed restoring force. for other related work, see, Baptistini and Táboas [3], Bellman and Cooke [4], Boese [8], Brauer [9], Cooke and van den Driessche [18], Cooke and Grossman [17], Huang [30], Mahaffy [33], Ruan and Wei [38], etc.

4. Steady states and stability for positive delays

Consider the system (2.1), and suppose that ω > 0. From the table 1 (see, section 2), we distingue between two cases: αδ < σ and αδ > σ.

Case 1: ω > 0 and αδ < σ. The system (2.1) has a unique positive equilibrium P_{0} given by P_{0} = (\frac{σ}{δ} , 0) and the linearized system around P_{0} takes the form

\begin{matrix}\frac{dx}{dt} = ω \frac{σ}{δ} y(t − τ) − δx \\ \frac{dy}{dt} = (α − \frac{σ}{δ})y \end{matrix}  (4.1)

which leads to the characteristic equation

W (λ) = (λ + \frac{σ}{δ} − α)(λ + δ).  (4.2)

Then we have the following result.

Proposition 4.1. Under the hypotheses ω > 0 and αδ < σ, the equilibrium point P_{0} is absolutely stable.

Proof. From the characteristic equation (4.2) and lemma 3.1, it is easy to obtain the result (see, [28, 29]). □

Case 2: ω > 0 and αδ > σ. In this case, system (2.1) has two equilibrium points (see, Table 1) P_{0} = (\frac{σ}{δ} , 0) and P_{2} = (x_{2}, y_{2}) where

x_{2} = \frac{−α(βδ − ω) + \sqrt{Δ}}{2ω}, y_{2} = \frac{α(βδ + ω) − \sqrt{Δ}}{2αβω}

with Δ = α^{2}(βδ − ω)^{2} + 4αβσω. From the characteristic equation (4.2), we deduce the following result.

Proposition 4.2. Under the hypotheses ω > 0 and αδ > σ, the equilibrium point P_{0} is unstable for all positive time delay.

Proof. For the proof of this proposition, from the characteristic equation (4.2). It is obvious to check the result. □

In the next, we shall study the stability of the non-trivial equilibrium point P_{2}. Let u = x − x_{2} and v = y − y_{2}, by linearizing system (2.1) around the non-trivial equilibrium point P_{2}, we obtain the linear system

\begin{matrix}\frac{du}{dt} = ωx_{2}v(t − τ) − ωy_{2}u(t − τ) − δu \\ \frac{dv}{dt} = −y_{2}u + (α − 2αβy_{2} − x_{2})v \end{matrix}  (4.3)

The characteristic equation of equation (4.3) has the form

W (λ, τ) = λ^{2} + pλ + r + (sλ + q)e^{−λτ} = 0,  (4.4)

which is the same equation presented in section 3.2. Where p = δ + αβy_{2} > 0, r = δαβy_{2} > 0, s = −ωy_{2} < 0 and q = αωy_{2}(1 − 2βy_{2}) > 0. The stability of the equilibrium point P_{2} is a result of the localization of the roots of the equation

W (λ, τ) = 0.

Then we have the following theorem.

Theorem 4.1. Assume 0 < \frac{ω}{β} < α, αδ > σ, α > 0 and β > 0. Then P_{2} is conditionally stable.

Proof. From the expressions of p, q, s and r and the paper [26], we have p + s > 0 and q + r > 0. Therefore, the hypotheses (H1), (H2) are satisfied and the steady state P_{2} is asymptotically stable (see section 3.2) for τ = 0. Since r^{2} − q^{2} = α^{2}y_{2}^{2}(δ^{2}β^{2} − ω^{2}(1 − 2βy_{2})^{2}), the sign of r^{2} − q^{2} is deduced from the sign of (δβ − ω^{2}(1 − 2βy_{2})) = −2αω − \sqrt{Δ} which is negative. Therefore, r^{2} − q^{2} < 0

Now, we compute s^{2} − p^{2} + 2r. From the expressions of p, s and r, we have

s^{2} − p^{2} + 2r = (ω − αβ)(ω + αβ)y^{2} − δ^{2}.

As \frac{ω}{β} < α, we have

s^{2} − p^{2} + 2r < 0

and the hypothesis (H4) of section 3.2 is satisfied. From theorem 3.1 (ii), the equilibrium point P_{2} is conditionally stable and there exist τ_{l} such that: P_{2} is asymptotically stable for τ ∈ [0, τ_{l}) and unstable for τ > τ_{l}. For τ = τ_{l} the characteristic equation (4.4) has a pair of purely imaginary roots ±iζ, where

τ_{l} = \frac{1}{ζ_{l}} arccos \{\frac{q(ζ_{l}^{2} − r) − psζ_{l}^{2}}{s^{2}ζ_{l}^{2} + q^{2}} \},  (4.5)
ζ_{l} = \frac{1}{2}(s^{2} − p^{2} + 2r) + \frac{1}{2}[(s^{2} − p^{2} + 2r)^{2} − 4(r^{2} − q^{2})]^{1/2}.  (4.6)
□

In the next sections, we will study the occurrence of Hopf bifurcation for smaller and larger delays.

Notation: The index s is designed for small time delays and the index l is designed for large time delays.

Let z(t) = (u(t), v(t)) = (x(t), y(t)) − (x_{2}, y_{2}), then the system (2.1) is written as a functional differential equation (FDE) in C := C([−τ, 0], ℝ^{2}):

\frac{dz(t)}{dt} = L(τ)z_{t} + f(z_{t}, τ)  (4.7)

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

L(τ)φ = (\begin{matrix}ωy_{2}φ_{1}(−τ) + ωx_{2}φ_{2}(−τ) − δφ_{1}(0) \\ −y_{2}φ_{1}(0) + (α − 2αβy_{2} − x_{2})φ_{2}(0) \end{matrix})

and

f(φ, τ) = (\begin{matrix}σ + ωφ_{1}(−τ)φ_{2}(−τ) + ωx_{2}y_{2} − δx_{2} \\ −αβφ_{2}^{2}(0) + αy_{2} − αβy_{2}^{2} − φ_{1}(0)φ_{2}(0) − x_{2}y_{2} \end{matrix})

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

5. Hopf bifurcation occurrence for small delays

For small delays, let e^{−λτ} ≃ 1−λτ, then the characteristic equation (4.4) becomes

W_{0}(λ, τ) = (1 − sτ)λ^{2} + (p + s − qτ)λ + r + q = 0.  (5.1)

Since the equilibrium point P_{2} is asymptotically stable for τ = 0, by Rouche’s theorem, there exist τ_{s} such that P_{2} asymptotically stable for τ < τ_{s} and unstable for τ > τ_{s}, where τ_{s} is the value for which the characteristic equation (5.1) has a pair of purely imaginary roots. Let λ = iζ, then W_{0}(iζ, τ) = 0 if and only if

\begin{matrix}(1 − sτ)ζ^{2} − q − r = 0, \\ p + s − qτ = 0 . \end{matrix}  (5.2)

Then, from the second equation of (5.2), we have τ_{s} = \frac{p+s}{q} and ζ_{s} = \sqrt{\frac{q(q+r)}{q−s(p+s)}}.

We deduce the following result of stability of the non-trivial equilibrium point P_{2} for small delays.

Theorem 5.1. Assume 0 < ω/β < α, αδ > σ, α > 0 and β > 0. Then, there exists τ_{s} such that: P_{2} is asymptotically stable for τ ∈ [0, τ_{s}) and unstable for τ > τ_{s}. For τ = τ_{s} the characteristic equation (5.1) has a pair of purely imaginary roots ±iζ_{s},

where τ_{s} = \frac{p+s}{q} and ζ_{s} = \sqrt{\frac{q(q+r)}{q−s(p+s)}}.

Proof. Since q−s(p+s) = q+ωδy_{2}+ωy_{2}^{2}(αβ−ω) and q > 0 and from the hypothesis αβ > ω, we have p + s > 0 and q − s(p + s) > 0. Then the quantities of τ_{s} = \frac{p+s}{q}

and ζ_{s} = \sqrt{\frac{q(q+r)}{q−s(p+s)}} are well defined. □

Now, we apply the Hopf bifurcation theorem, see [28], to show the existence of a non-trivial periodic solutions of system (4.7) bifurcating from the non trivial steady state P_{2}. 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.

Next we state the main result of this paper for small delays.

Theorem 5.2. Assume 0 < ω/β < α, αδ > σ, α > 0 and β > 0. There exists ε_{s} > 0 such that, for each 0 ≤ ε < ε_{s}, equation (4.7) has a family of periodic solutions p_{s}(ε) with period T_{s} = T_{s}(ε), for the parameter values τ = τ(ε) such that

p_{s}(0) = P_{2}, T_{s}(0) = \frac{2π}{ζ_{s}} and τ(0) = τ_{s}, where τ_{s} = \frac{p+s}{q} and ζ_{s} = \sqrt{\frac{q(q+r)}{q−s(p+s)}} are given in equation (5.2).

Proof. We apply the Hopf bifurcation theorem introduced in [28]. From the expression of f in (4.7), we have

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

From equation (5.2) and theorem 5.1, the characteristic equation (5.1) has a pair of simple imaginary roots λ_{s} = iζ_{s} and \overline{λ}_{s} = −iζ_{s} at τ = τ_{s}.

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

λ′(τ) = −\frac{∂W_{0}(λ, τ)/∂τ}{∂W_{0}(λ, τ)/∂λ},  (5.3)

for τ in a neighborhood of τ_{s}. Letting, λ(τ) = p(τ) + iq(τ), from (5.3) we have

p′(τ)_{/τ=τ_{s}} = \frac{q}{2(1 − sτ_{s})} .

From the hypothesis 0 < ω/β < α, we conclude that

p′(τ)_{/τ=τ_{s}} > 0,

which completes the proof. □

6. Hopf bifurcation occurrence for large delays

For large delays τ, let λ = κ + iζ. According to the Hopf bifurcation theorem [28], we come to the main result of this paper for large time delays.

Theorem 6.1. Assume 0 < ω/β < α, αδ > σ, α > 0, β > 0 and δ close to 0. There exists ε_{1} > 0 such that, for each 0 ≤ ε < ε_{1}, equation (4.7) has a family of periodic solutions p_{l}(ε) with period T_{l} = T_{l}(ε), for the parameter values τ = τ(ε) such that p_{l}(0) = P_{2}, T_{l}(0) = \frac{2π}{ζ_{l}} and τ(0) = τ_{l}, where τ_{l} and ζ_{l} are given respectively in equations (4.5) and (4.6).

Proof. As in the previous section, we apply the Hopf bifurcation theorem introduced in [28]. From the expression of f in (4.7), we have,

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

From equation (4.4) and theorem 4.1, the characteristic equation (4.4) has a pair of simple imaginary roots λ_{l} = iζ_{l} and \overline{λ}_{l} = −iζ_{l} at τ = τ_{l}. From equation (4.4), W (λ_{l}, τ_{l}) = 0 and \frac{∂}{∂λ}W (λ_{l}, τ_{l}) = 2iζ_{l} + p + (s − τ(isζ_{l} + q))e^{−iζ_{l}τ_{l}} ≠ 0. According to the implicit function theorem, there exists a complex function λ = λ(τ) defined in a neighborhood of τ_{l}, such that λ(τ_{l}) = λ_{l} and W (λ(τ), τ) = 0 and

λ′(τ) = −\frac{∂W (λ, τ)/∂τ}{∂W (λ, τ)/∂λ},  (6.1)

for τ in a neighborhood of τ_{l}. Then

λ′(τ) = −λ\frac{sλ^{3} + (s^{2}p + q)λ^{2} + (sr + pq)λ + qr}{τsλ^{3} + (s + τ(sp + q))λ^{2} + (2q + τ(sr + pq))λ + pq − sr + qr} (6.2)

From equation (6.2) we have

\begin{matrix}κ′(τ)|_{τ=τ_{l}} = \\ ζ_{l}^{2} \frac{s^{2}ζ_{l}^{4} + (sqr(τ − 1) + 2q^{2})ζ_{l}^{2} + sr^{2}(q − sr) + pq^{2}(p + r) − qr(2q + τ(sr + pq))}{A^{2} + B^{2}}, \end{matrix}

(6.3)

where

A = −(s + τ(sp + q))ζ_{l}^{2} + pq − sr + qr,
B = −τsζ_{l}^{2} + (2q + τ(sr + pq))ζ_{l}.

From the expression of r, when δ is close to 0, then r is very small. From equation (6.3), we conclude that,

κ′(τ)_{/τ=τ_{l}} > 0

Therefore, the transversality condition is verified, which completes the proof. □

References

[1] J. Adam and N. Bellomo, A survey of models on tumor immune systems dynamics, Birkhäuser, Boston (1996).

[2] L. Arlotti, N. Bellomo and E. De Angelis; Generalized kinetic (Boltzman) models: Mathematical structures and application, Math. Models Meth. Appl. Sci., Vol. 12, pp. 567-592 [3] M. Baptistini and P. Táboas, On the stability of some exponential polynomials, J. Math.

Anal. Appl., Vol. 205, pp. 259-272 (1997).

[4] R. Bellman and K. L. Cooke, Differential-difference equations, Academic Press, New York (1963).

[5] N. Bellomo and L. Preziosi, Modeling and mathematical problems related to tumor immune system interactions, Math. Comp. Modelling, Vol. 31, pp. 413-452 (2000).

[6] Eds. Bellomo and M. Pulvirenti, Modeling in applied sciences: A kinetic theory approch, Birkhäuser (2000).

[7] Eds. Bellomo and M. Pulvirenti, Special issue on the modeling in applied sciences by methods of transport and kinetic theory, Math. Comp. Modelling, Vol. 12, pp. 909-990 (2002).

[8] F. G. Boese, Stability criteria for second-order dynamical systems involving several time delays, SIAM J. Math. Anal., Vol. 26, pp. 1306-1330 (1995).

[9] F. Brauer, Absolute stability in delay equations, J. Differential Equations, Vol. 69, pp. 185-191 (1987).

[10] R. Bürger, The mathematical theory of selection, recombination and mutation, Wiley (2000).

[11] M. Bodnar, The nonnegativity of solutions of delay differential equations, Appl. Math. Lett., Vol. 10, No. 6, pp. 91-95 (2000).

[12] M. Bodnar and U. Foryś, Behavior of solutions to Marchuk’s model depending on a time delay, Int. Math. Comput. Sci., Vol. 10, No. 1, pp. 97-112 (2000).

[13] M. Bodnar and U. Foryś, Periodic dynamics in the model of immune system, Appl. Math., Vol. 27, No. 1, pp. 113-126 (2000).

[14] H. M. Byrne, The effect of time delay on the dynamics of avascular tumor growth, Math. Biosci., Vol. 144, No. 2, pp. 83-117 (1997).

[15] M. A. J. Chaplain Eds., Special issue on mathematical models for the growth, development and treatment of tumors, Math. Models Meth. Appl. Sci., Vol. 9 (1999).

[16] Y. S. Chin, Unconditional stability of systems with time lags, Acta Math. Sinica, Vol. 1, pp. 125-142 (1960).

[17] K. L. Cooke and Z. Grossman, Discrete delay, distributed delay and stability switches, J. Math. Anal. Appl., Vol. 86, pp. 592-627 (1982).

[18] K. L. Cooke and P. van den Driessche, On zeros some transcendental equations, Funkcialaj Ekvacioj, Vol. 29, pp. 77-90 (1986).

[19] L. Desvillettes and C. Prvots, Modelling in population dynamics through kinetic-like equations, preprint n. 99/19 of the University of Orlands, dpartement de mathmatiques.

[20] L. Desvillettes, C. Prévots and R. Ferrieres, Infinite dimensional reaction-diffusion for population dynamics, (2003).

[21] O. Dieckmann and J. P. Heesterbeek, Mathematical epidemiology of infectious diseases, Wiley, New York (2000).

[22] G. Forni, R. Fao, A. Santoni and L. Frati Eds., Cytokine induced tumor immunogeneticity, Academic Press, New York (1994).

[23] U. Foryś, Marchuk’s model of immune system dynamics with application to tumor growth, J. Theor. Med., Vol. 4, No.1, pp. 85-93 (2002).

[24] U. Foryś and M. Kolev, Time delays in proliferation and apoptosis for solid avascular tumor, Prep. Institute of Applied Mathematics and Machanics, No. RW 02-10 (110), Warsaw University (2002).

[25] U. Foryś and A. Marciniak-Czochra, Delay logistic equation with diffusion, Proc. 8-th Nat. Conf. Application of Mathematics in Biology and Medicine, Lajs, pp. 37-42 (2002).

[26] M. Galach, Dynamics of the tumor-immune system competition the effect of time delay, Int. J. Appl. Comput. Sci., Vol. 13, No. 3, pp. 395-406 (2003).

[27] L. Greller, F. Tobin and G. Poste, Tumor hetereogenity and progression: Conceptual foundation for modeling, Invasionand Metastasis, Vol. 16, pp. 177-208 (1996).

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

[29] J. K. Hale, Theory of functional differential equations, New York, Springer-Verlag (1997).

[30] W. Huang, Algebraic criteria on the stability of the zero solutions of the second order delay differential equations, J. Anhui University, pp. 1-7 (1985).

[31] D. Kirschner and J. C. Panetta, Modeling immunotherapy of the tumor-immune interaction, J. Math. Biol., Vol. 37, No. 3, pp. 235-252 (1998).

[32] V. A. Kuznetsov and M. A. Taylor, Nonlinear dynamics of immunogenic tumors: Parameter estimation and global bifurcation analysis, Bull. Math. Biol., Vol. 56, No. 2, pp. 295-321 (1994).

[33] J. M. Mahaffy, A test for stability of linear differential delay equations, Quart. Appl. Math., Vol. 40, pp. 193-202 (1982).

[34] R. M. May and M. A. Nowak, Virus dynamics (mathematical principles of immunology and virology), Oxford Univ. Press (2000).

[35] H. Mayer, K. S. Zanker and U. der Heiden, A basic mathematical model of the immune response, Chaos, Vol. 5, No. 1, pp. 155-161 (1995).

[36] A. S. Perelson and G. Weisbuch, Immunology for physists, Rev. modern phys. Vol. 69, pp. 1219-1267 (1997).

[37] L. Perko, Differential equations and dynamical systems, New York: Springer (1991).

[38] S. Ruan and J. Wei, On the zeros of transcendental functions with applications to stability of delay differential equations with two delays, Dynamics of Continuous, Discrete and Impulsive Systems Series A: Mathematical Analysis 10, pp. 863-874 (2003).

[39] J. Waniewski and P. Zhivkov, A simple methematical model for tumor-immune system interactions, Proc. 8-th Nat. Conf. Application of Mathematics in Biology and Medicine, Lajs, pp. 149-154 (2002).

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@math.net