Continuous Solutions to a Viral Infection Model with General Incidence Rate, Discrete State-dependent Delay, Ctl and Antibody Immune Responses

A virus dynamics model with intracellular state-dependent delay and a general nonlinear infection rate functional response is studied. We consider the case of merely continuous solutions which is adequate to the discontinuous change of parameters due to, for example, drug administration. The Lyapunov functionals technique is used to analyze stability of an interior infection equilibrium which describes the case of both CTL and antibody immune responses activated.


Introduction
In our research we are interested in mathematical models of viral diseases.We notice that early models [16,18] contain three variables: susceptible host cells, infected cells and free virus.Such models do not take into account the immune responses.
Immune responses could be innate (nonspecific responses), and specific (adaptive responses).More on the basic immunological background see e.g.[32].We concentrate on specific, adaptive immune responses.Main of them are effector responses i.e, they directly fight the pathogen.The two effector responses are antibodies and CTL (cytotoxic T lymphocytes or killer T cells).Antibodies can attach to the pathogen and neutralize it while CTL attack infected cells.See also [37] and references therein.We notice that the anti-virus antibody detection is commonly used in the diagnostic laboratory.The relative balance of both types of adaptive immune response "can be a decisive factor that determines whether patients Email: rezounenko@yahoo.comare asymptomatic or whether pathology is observe" [31].These lead to introduction of two additional variables of both adaptive immune responses [31,32] (see also [35] and references therein).
We will study a generalization of the model (1.1) below which contains five variables: susceptible (noninfected) host cells T, infected cells T * , free virus V, a CTL response Y, and an antibody response A. In case of bilinear nonlinearities and one constant concentrated delay (see, for example, [34]) it has the following form (1.1) Here the dot over a function denotes the time derivative i.e., Ṫ(t) = dT (t)  dt , all the constants λ, d, k, δ, p, N, c, q, β, γ, g, b, ω are positive.As for the immune responses, the fourth equation describes the regulation of CTL response and pY(t)T * (t) (in the second equation) being the rate of killing of infected cells by lytic immune response.The fifth equation describes the regulation of antibody response and qA(t)V(t) (in the third equation) being the rate of virus neutralization by antibodies [31, p. 1744].In (1.1), h denotes the delay between the time the virus contacts a target cell and the time the cell becomes actively infected (starts to produce new virions).
In the above model (1.1), the standard bilinear incidence rate is used according to the principle of mass action.For more details and references on the models of infectious diseases with more general types of nonlinear incidence rates f (compare the first two equations in systems (1.1) and (1.2)) see e.g.[6,11] and our assumptions and examples below.In paper [36], following [10,29,34], authors assume that the infection rate of the virus dynamics models is given by the Beddington-DeAngelis functional response [1,3], f (T, V) = kTV 1+k 1 T+k 2 V , where k, k 1 ≥ 0, k 2 > 0 are constants.The Lyapunov asymptotic stability [13] of points of equilibrium is studied for the following model with constant concentrated delay It is evident that the constancy of the delay is an extra assumption which essentially simplifies the analysis, but is not motivated by the biological background of the model.It was a reason (see e.g.[14,30]) to discuss distributed delay models as an alternative to discrete constant delay ones.One could consider a time-dependent delay h(t) if some biologically motivated properties of h(t) are available.We propose an another approach.
Our first goal is to remove the restriction of the constancy of the delay and investigate the well-posedness and Lyapunov stability of the following virus infection model (1.3) with a general functional response f and state dependent delay.It appears that the analysis essentially differs from the constant delay case.To the best of our knowledge, such models have been considered for the first time in [23].It is well known that differential equations with state dependent delay are always nonlinear by its nature (see the review [9] for more details and discussion).
As usual in a delay system with (maximal) delay h > 0 [4,7,12], for a function We denote the space of continuous functions by C ≡ C([−h, 0]; R 5 ) equipped with the sup-norm.In the above notations, we use u(t) = (T(t), T * (t), V(t), Y(t), A(t)) and consider a continuous functional (state dependent delay) η : C → [0, h].The delay η is obviously bounded since it can not exceed the life span of the host target cells.Now we are ready to present the system under consideration with a general functional response f (T, V) satisfying natural assumptions presented below.
See also examples in Section 3.1.We notice that the term e −ωh in front of f (see in the second equation (1.3)), in fact, states that only a part of the cell population survived during the virus incubation period.Clearly, it should be less than 1.It is an assumption which is not too precise in nonlinear systems.It could be regarded as a coefficient belonging to (0, 1) and could be incorporated into the definition of the function f .We keep this coefficient in the form of e −ωh for the only reason to simplify for the reader the comparison of computations with the constant delay case.
It is well known that merely continuous solutions to differential equations with discrete state-dependent delay may be non-unique (see examples in [5]).There are two different ways to guarantee the uniqueness of solutions as well as the well-posedness.The first one is to restrict the set of initial functions to more smooth ones [9].This way was used for the viral model in [23].The second way is to restrict the class of state-dependent delays [19,21] and work with continuous initial functions and solutions.In the current note, in contrast to [23], we discuss the second way, which is more convenient for our second goal discussed below.
There is a number of papers on non-delayed and (constant) delay viral models which are concentrated on the local and/or global stability of stationary solutions (see e.g.[6,11,31,32,35]).In case of the global asymptotic stability of a nontrivial disease stationary solutions is proved, one should conclude that the virus will never be eradicated i.e. the disease is in the chronic stage.Such results are very important for diagnostic purposes.On the other hand, after the diagnosis of a viral disease is confirmed, the prime goal is to find a way to cure the patient.Since a medical research is quite expensive and could last for decades, the mathematical models proved to be important and efficient.Our second goal is to present a model and choose a proper space for solutions which could be appropriate for therapy, including drug administration.The main motivation is the situation (see e.g.[17,25]) when the drug effectiveness was decreased in a stepwise manner.In terms of (1.3), the parameter N could change its value in a discontinuous way (see equation (2) in [25, p. 920]).One could see that at time moment of discontinuity of (any) parameter, the solution is continuous, but not differentiable (cf. Figure 2-B in [25, p. 921] and also Fig. 1 in [17, p. 23]).To realize how frequent such a discontinuity of the time-derivative could appear, one should compare a virus generational time (which is h in our notations (1.3)) and a drug regimen of treatment.Taking as an example HIV, we found in [15] that "a total HIV generational time of 25 h in vitro and is much shorter than our 52 h estimate from in vivo delays".On the other hand, the standard treatment schedules are two times a day or once-a-day pills.It suggests that on any history time segment [t − h, t] one has one or more discontinuities of the time-derivative.Moreover, co-infections (which are not too rare) by other pathogens, consume resources of the immune system, which leads to changes in other parameters (in terms of (1.3) the parameters β and g could change as well).
In study of local stability of an equilibrium of a system one could also use the method of linearized stability.For state-dependent delay equations this method is available the C case [2,8] and in the C 1 case [9,28].For the continuous case we use the Lyapunov functions approach [13].
The paper is organized as follows.In Section 2, we discuss and choose a natural set of initial data and prove the existence and uniqueness of solutions.Next we prove that the set is invariant.Section 3 is devoted to the stability properties of a stationary solution.We study the stability of an interior equilibrium which describes the case when both CTL and antibody immune responses are activated.We believe this infection equilibrium is only biologically meaningful in the study of the disease.

Basic properties
We equip the system (1.3) with an initial condition where We assume the nonlinear function Here T max , V max are defined in (2.2).
The first result is the following.Theorem 2.2.Let η : C → [0, h] (state dependent delay) and f be continuous functionals.Then (i) for any initial function ϕ ∈ C there exist continuous solutions to (1.3), (2.1); (ii) if additionally, η satisfies (H ign ) and f satisfies (H1 f ), then for any initial function ϕ ≡ (T 0 , T * 0 , V 0 , A 0 , Y 0 ) ∈ Ω C , the problem (1.3), (2.1) has a unique solution.The solution continuously depends on the initial function and satisfies Remark 2.3.In [23], the following assumption on the state-dependent delay has been used In the current work we do not need this restriction i.e. the delay could vanish.
Proof.(i) The existence of continuous solutions is guaranteed by the continuity of the righthand side of (1.3) and classical results on delay equations [4,7].
(ii) The proof follows the line of [23,Theorem 2].The main essential difference is that we could use the quasi-positivity property of the right-hand side of (1.3) (see e.g.[26, Theorem 2.1, p. 81]).We stress that in the case of state-dependent delay we cannot directly apply [26, Theorem 2.1, p. 81] because it relies on the Lipschitz property of the right-hand side of a system, which is not the case for (1.3).Instead, we use the corresponding extension to the state-dependent delay case [22] which relies on the assumption (H ign ).It essentially simplifies the proof (cf.[23,Theorem 2]).The upper bounds on all the five coordinates in (2.2) follow from an easy variant of the Gronwall's lemma, which is formulated for the simplicity as In the case b = +∞, for any ε > 0 there exists t ε ≥ a such that (t) ≤ c 1 c −1 2 + ε for all t ≥ t ε .It gives the invariance of the set Ω C (2.2).We do not repeat details here (one can check and find differences [23,Theorem 2])).The continuous dependence on the initial function follows from (H ign ) as in [19].
We could conclude that the problem (1.3), (2.1) is well-posed in Ω C ⊂ C in the sense of Hadamard.

Stationary solutions
We look for nontrivial disease stationary solutions (1.3).Consider the system with u(t) = u(t − η(u t )) = u and denote the coordinates of a stationary solution by Since the stationary solutions of (1.3) do not depend on the type of delay (state-dependent or constant) we have (see e.g.[23,36]) Our case differs from [23,36] by a more general class of nonlinearities f .The last two equations in (2. where T is the unique positive root of (2.4).
We notice that, from biological point of view, (H2), (H3) are standard assumptions on reproduction numbers, which are given here in a short form.We could summarize the above calculations in the following Proposition 2.5 (cf.[23,Lemma 7]).Let assumptions (H2) and (H3) be satisfied and f satisfy (H1 f ).Then the system (2.3) has a unique solution ( T, T * , V, Y, A) (the unique stationary solution of (1.3)).All the coordinates are positive, T is the unique positive root of the equation (2.4) and coordinates satisfy (2.5) We will use these equations connecting the coordinates of the stationary solution in our study of the stability properties (cf.(2.3)).

Stability properties
In study of differential equations with nonnegative variables the function v(x) = x − 1 − ln x : (0, +∞) → R + plays an important role in construction of Lyapunov functionals.One can see that v(x) ≥ 0 and v(x) = 0 if and only if x = 1.The derivative equals v(x) = 1 − 1 x , which is evidently negative for x ∈ (0, 1) and positive for x > 1.The graph of v explains the use of the composition v x x 0 in the study of the stability properties of an equilibrium x 0 .Another important property is the following estimate To check it, one simply observes that all three functions vanish at x = 1 and As before, we denote u(t) = (T(t), T * (t), V(t), Y(t), A(t)) and ϕ = ( T, T * , V, Y, A) the stationary solution of (1.3).
We assume f satisfies in a neighborhood of ( T, V) Remark 3.1.It is easy to see the clear geometrical meaning of the assumption (H2 f ).First, we mention that the first coordinates of f in (H2 f ) are equal.For any fixed first coordinate we could consider f T (V) ≡ f (T, V) and its graph.In the case of differentiable functions, (H2 f ) implies 0 < d dV f T (V) < f T ( V)/ V.In our study we do not assume differentiability of f .More discussion and examples are in Section 3.1.
The following assumption on the state-dependent functional η is based on the property (H ign ).Consider an arbitrary ϕ ∈ C and its arbitrary extension ϕ ext (s), s ∈ [−h, η ign ] with constant η ign > 0 defined in (H ign ).Due to the property (H ign ) we could define an auxiliary function η ϕ (t) ≡ η(ϕ ext t ) for t ∈ [0, η ign ].Since both η and ϕ are continuous we see that η ϕ ∈ C[0, η ign ].We are interested in the (right) derivative of η ϕ at zero and its properties.Now we are ready to formulate our next local assumption on η.
(H2 η ) There is a µ-neighborhood of the stationary point ϕ such that (for any ϕ ∈ C satisfying ϕ − ϕ C < µ) the following two properties hold

Remark 3.2. It is easy to see, by the definition of η
We also mention that the existence of η + (ϕ) ∈ R does not require the differentiability of ϕ (see Section 3.1 for examples).
Our result is the following.
Proof.Let us introduce the following Lyapunov functional with state-dependent delay along a solution of (1.3) A particular case of the constant delay functional and the Beddington-DeAngelis functional response f has been considered in [23,36].The main difference is in the state-dependence of the lower bound of the last integral in (3.3).In [23], a particular case of such a functional was studied along continuously differentiable solutions.
Let us compute the time derivative of the last integral in (3.3) along a solution of (1.3) We see the main difference with the constant-delay case in the appearance of the term Remark 3.4.In case (investigated in [23]) of both solutions and state-dependent delay η are continuously differentiable, for any , where [(Dη)(u t )](•) is the Fréchet derivative of η at point u t .Hence, (for a solution in µ-neighborhood of the stationary solution ϕ) the estimate guarantees the property (3.2) due to the boundedness of (Dη)(ψ) L(C;R) as µ → 0 (here ψ − ϕ C < µ).Now our solutions are merely continuous, so we cannot use the above arguments.Instead, we use (H2 η ).
We use the same notations as in [23,36] to simplify for the reader the comparison of the computations.We also remind that the state-dependence is present in both the system and the Lyapunov functional and the class of nonlinear functions f is wider.
We have along a continuous solution Here we used the last equality in (2.5) and notation S sdd (t) defined in (3.4).Opening parenthesis, grouping similar terms and canceling some of them, we obtain To save the space we omit long computations where we intensively used equations (2.5), for example, . Next, we add To save the space, let us denote the sum {...} above as R 1 (t) i.e., Now we add Considering the first four terms in [. ..] above we suggest to split the logarithm as follows Substitution of (3.7) into (3.6)implies As before we used the function v(x) = x − 1 − ln x to save the space.Next, we can rewrite the first term in (3.8) as The last inequality due to the monotonicity of f (see assumption (H1 f )).Now we transform R 1 (t), defined in (3.5).Calculations give It is clear that assumption (H2 f ) gives R 1 (t) ≥ 0 in a neighborhood of the stationary solution.Moreover it implies the existence of constants We substitute (3.9) into (3.8) to get where and S sdd (t) is defined in (3.4).One can see, using v(x) ≥ 0, (3.9) and (3.11) that D sdd1 (t) ≥ 0.
Remark 3.5.It is easy to check that D sdd1 (t) = 0 if and only if It follows from the property v(x) = 0 if and only if x = 1 and also (3.11).
Our goal is to prove that there is a neighborhood of u ∈ C, where d dt U sdd (t) < 0 (except the stationary point u).We notice that D sdd1 (t) ≥ 0, while the sign of S sdd (t) is undefined.We plan to show that there is a neighborhood of the stationary point, where |S sdd (t)| < D sdd (t).
Here L ξ,µ is the Lipschitz constant of ξ.Hence, the delay (3.18) satisfies all the needed conditions.
It is easy to see that more general delay terms could be used.For example, with a differentiable ρ : R → [0, h].The example (3.18) is a particular case with ρ(s) ≡ s and κ(s) ≡ 1.

2.
One can check that the Beddington-DeAngelis functional response [1,3] of the form f (T, V) = kTV 1+k 1 T+k 2 V , with k, k 1 ≥ 0, k 2 > 0 satisfies (H2 f ) globally.We also mention that the Beddington-DeAngelis functional response includes as a special case (k 1 = 0) the saturated incidence rate f (T, V) = kTV 1+k 2 V .In our study we need the property (H2 f ) in a small neighborhood of ( T, V) only.

3.
Another example of the nonlinearity is the Crowley-Martin incidence rate f (T, V) = kTV (1+k 1 T)(1+k 2 V) , with k ≥ 0, k 1 , k 2 > 0 (see e.g.[33]).Remark 3.7.In this article we propose a rather general framework for state-dependent delay viral models.The assumptions on the state-dependent delay guarantee the well-posedness and local stability of stationary solutions.Since a mathematical model is always a simplification of real life processes, many biologically important factors are to be reflected in more complex systems.The life cycle of particular cells and their interaction with viral particles could essentially differ from one organ to another.Further assumptions on the delay functional could naturally appear when studying a particular viral infection with biological characteristics of a target organ, its cells and type of the virus.

4 )
3) imply T * = γ β , V = b g .This and the third equation give A = Nδγg−βcb βqb .The positivity of A holds provided the constants in the system (1.3) satisfy (H2) Nδγg > βcb.Substitution of the value V = b g into the first equation of (2.3) gives λ − d T = f T, Since f (•, b g ) is strictly increasing in the first coordinate, continuous and f (0, b g ) = 0, it is easy to see that the last equation (2.4) has the unique solution T ∈ (0, b g ).This unique positive root is denoted below by T. The first two equations in (2.3) give (remind that T * is already known) Y = λ−d T−e ωh δ T * e ωh p T * .The positivity of Y follows from the following assumption