Spatiotemporal Dynamics of a Generalized Viral Infection Model with Distributed Delays and CTL Immune Response

In this paper, we propose and investigate a diffusive viral infection model with distributed delays and cytotoxic T lymphocyte (CTL) immune response. Also, both routes of infection that are virus-to-cell infection and cell-to-cell transmission are modeled by two general nonlinear incidence functions. The well-posedness of the proposed model is also proved by establishing the global existence, uniqueness, nonnegativity and boundedness of solutions. Moreover, the threshold parameters and the global asymptotic stability of equilibria are obtained. Furthermore, diffusive and delayed virus dynamics models presented in many previous studies are improved and generalized.


Introduction
During human infections with viruses such as human immunodeficiency virus (HIV), human T-cell leukemia virus (HTLV), hepatitis B virus (HBV) and hepatitis C virus (HCV), cytotoxic T lymphocyte (CTL) cells play a crucial role in antiviral defence by attacking and killing infected cells.So, modeling the role of CTL immune response in viral infection has attracted the attention of many researchers.In 1996, Nowak and Bangham [1] proposed a basic mathematical model by assuming that the infection process is bilinear and follows the principle of mass action.However, as a nonlinear relationship between parasite dose and infection rate has been frequently observed in experiments in [2,3], this bilinear incidence was replaced by Beddington-DeAngelis functional response in [4] and by a more general incidence function in [5].
In the above classical models that are formulated by ordinary differential equations (ODEs), the cell infection is instantaneous and only caused by contact with the free virus.In reality, there are two routes of infection and also time delays in cell infection and virus production.Motivated by these biological reasons, Li et al. [6] proposed a mathematical model formulated by delay differential equations (DDEs) to describe the global dynamics of HIV infection with CTL immune response.This delayed model is an extension of [1] that considers Holling type-II functional response and two kinds of discrete delays, one in cell infection and the other in virus production.Also, the authors of [7] improved the model of Nowak and Bangham [1] by introducing a discrete delay in cell infection and using a Crowley-Martin type incidence function.In 2016, Wang et al. [8] introduced an infinite distributed delay in cell infection in order to improve the basic model with CTL immune response [1], and they also considered both routes of infection, virus-to-cell infection and cell-to-cell transmission.Furthermore, a recent work presented in [9] studied the dynamical behavior of a viral infection model with two types of distributed time delays, CTL immune response and saturated incidence rates for both routes of infection.In this paper, we generalize all the ODE and DDE models presented in [1,[4][5][6][7][8][9] by proposing the following nonlinear system: where T(t), I(t), V(t) and Z(t) denote the densities of susceptible target cells, infected target cells, free virus particles and CTL cells at time t, respectively.The susceptible target cells are produced at constant λ, die at rate d and become infected by contact with free virus at rate f (T, I, V)V and by contact with infected cells at rate g(T, I)I.The parameters a and b are the death rates of infected cells and CTL cells.The parameter p represents the rate at which infected cells are killed by CTL cells, k is the production rate of free virus by an infected cell, and µ is the clearance rate of free virus.CTL cells expand in response to viral antigens derived from infected cells at rate cIZ.Further, we assume that the virus or infected cell contacts an uninfected cell at time t − τ and the cell becomes infected at time t, where τ is a random variable taken from a probability distribution f 1 (τ).The term e −α 1 τ represents the probability of surviving from time t − τ to time t, where α 1 is the death rate for infected but not yet virus-producing cells.In the same, we assume that the time necessary for the newly produced virions to become mature and infectious is a random variable with a probability distribution f 2 (τ).
The term e −α 2 τ denotes the probability of the immature virions surviving the delay period, where 1 α 2 is the average life time of an immature virus.Therefore, the integral ∞ 0 f 2 (τ)e −α 2 τ I(t − τ)dτ describes the mature viral particles produced at time t.The probability distribution functions f 1 (τ) and f 2 (τ) are assumed to satisfy f i (τ) ≥ 0 and As in [10,11], the incidence functions f (T, I, V) and g(T, I) for both routes of infection are continuously differentiable and satisfy the following hypotheses: (H 0 ) g(0,I) = 0, for all I ≥ 0; ∂g ∂T (T, I) ≥ 0 or g(T, I) is a strictly monotone increasing function with respect to T when f ≡ 0 and ∂g ∂I (T, I) ≤ 0, for all T ≥ 0 and I ≥ 0. (H 1 ) f (0, I, V) = 0, for all I ≥ 0 and V ≥ 0, (H 2 ) f (T, I, V) is a strictly monotone increasing function with respect to T or ∂ f ∂T (T, I, V) ≥ 0 when g(T, I) is a strictly monotone increasing function with respect to T , for any fixed I ≥ 0 and V ≥ 0, (H 3 ) f (T, I, V) is a monotone decreasing function with respect to I and V.
From a biological viewpoint, the above hypotheses are reasonable and consistent with reality.In fact, the first assumption (H 0 ) on the function g(T, I) means that the incidence rate by direct contact with infected cells is equal to zero if there are no susceptible cells.This incidence rate is increasing when the number of infected cells is constant and the number of susceptible cells increases.Also, it is decreasing when the number of susceptible cells is constant and the number of infected cells increases.Similarly, the second assumption (H 1 ) on the function f (T, I, V) means that the incidence rate by contact with free virus is equal to zero if there are no susceptible cells.By (H 2 ) and (H 3 ), this incidence rate is increasing when the numbers of infected cells and virus are constant and the number of susceptible cells increases.Also, it is decreasing when the number of susceptible is constant and the number of infected cells or free virus increases.Consequently, the more susceptible cells are, the more infectious events will occur.However, the higher the number of infected cells or the concentration of virus in the host is, the less infectious events will be [10,12,13].In addition, the functions f (T, I, V) and g(T, I) cover several types of incidence rates existing in the literature such as the classical bilinear incidence, standard incidence, Holling type-II functional response, Beddington-DeAngelis functional response, Crowley-Martin functional response and Hattaf-Yousfi functional response.
On the other hand, system (1) assumes that cells and viruses are well mixed, and ignores their mobility.Actually, viral propagation is a localized process [14] due to the fact that the virus is inherently unstable and the infection occurs mainly in lymphoid tissues.Also, the interaction between virus and the immune response tends to be local within the body of infected hosts [15].Further, cells are distributed in space and typically interact with the physical environment and other organisms in their spatial neighborhood [16].Therefore, it is more reasonable to study a reaction-diffusion version of system (1).So, the organization of this paper is as follows.In the next section, we present the reaction-diffusion version of (1) and some preliminary results.Section 3 is devoted to the global dynamics of the reaction-diffusion model.An application and some numerical simulations of our main results are presented in Section 4. Finally, the paper ends with mathematical and biological conclusions in the last section.

Model Formulation and Preliminaries
We first present a reaction-diffusion version of system (1) by taking into account the mobility of cells and viruses.Hence, system (1) becomes where T(x, t), I(x, t), V(x, t) and Z(x, t) are the densities of susceptible target cells, infected target cells, free virus particles and CTL cells at location x and time t, respectively.Here, we assume that the motion of the above four populations follows Fickian diffusion, meaning that the fluxes of these four populations are proportional to their concentration gradient and go from regions of high concentration to regions of low concentration, with the diffusion coefficients d T , d I , d V and d Z , respectively.is the Laplacian operator.The other parameters have the same biological meanings as those in system (1).
It is very important to note that our model (2) formulated by partial differential equations (PDEs) extends and generalizes many virus dynamics models existing in the literature.For instance, we obtain the diffused HBV infection model proposed by Wang et al. [17] when and g(T, I) = 0, where q > 0, β > 0 is a constant rate describing the infection process and δ(.) is the Dirac delta function.When and g(T, I) = 0, where 1 , 2 ≥ 0 are constants, we get the diffusive and delayed viral infection model with Beddington-DeAngelis functional response [18].Also, the diffusive and delayed viral infection model with Crowley-Martin functional response [19] is a special case of (2), it suffices to take and g(T, I) = 0.
Throughout this paper, we consider system (2) with initial conditions and zero-flux boundary conditions where Ω is a bounded domain in IR n with smooth boundary ∂Ω, and ∂ ∂ν indicates the outward normal derivative on ∂Ω.From the biological point of view, these conditions mean that the uninfected cells, infected cells, free virus particles and CTL cells do not move across the boundary ∂Ω.
We now study the well posedness of the PDE model ( 2) by establishing the global existence, uniqueness, nonnegativity and boundedness of solutions.To this end, we need some notations.Let X = C( Ω, IR 4 ) be the Banach space of continuous functions from Ω into IR 4 , and Moreover, we need the following lemma.Lemma 1.Let A, B and D be three constants with B = 0. Consider the following problem (5) Proof.Let ũ(t) be a solution to the ordinary differential equation Theorem 1.For any given initial condition φ ∈ C α satisfying (3), problem ( 2)-( 4) has a unique nonnegative solution.When the cells have the same diffusion coefficient (d T = d I = d Z ), this solution is defined on [0, +∞) and remains nonnegative and bounded for all t ≥ 0.
From the first equation of ( 2), we get By Lemma 1, we get This implies that T is bounded.Let The integral in G(x, t) is well-defined and differentiable with respect to t, due to T being bounded.Thus, where δ = min{a, b, d} and From Lemma 1, we have Thus, I and Z are bounded.It remains to prove that V is bounded.From the boundedness of I and ( 2)-( 4), we deduce that V satisfies the following system where This implies that V is bounded.From the above, we have proved that T(x, t), I(x, t), V(x, t) and Z(x, t) are bounded on Ω × [0, t max ).By the standard theory for semilinear parabolic systems [26], we deduce that t max = +∞.This completes the proof.
Clearly, system (2) has always one infection-free equilibrium E 0 (T 0 , 0, 0, 0), where T 0 = λ d , which represents the healthy state.Hence, we define the basic reproduction number for our PDE model as follows Biologically and as in [11,27], R 0 can be divided into parts as R 0 = R 01 + R 02 , where is the basic reproduction number corresponding to virus-to-cell infection mode, and R 02 = η 1 g( λ d , 0) a is the basic reproduction number corresponding to cell-to-cell transmission mode.
The other spatially uniform steady states of (2) satisfy the following system The last equation of (12) implies that Z = 0 or I = b c .Hence, we discuss two cases.
For the case when Z = 0, we get Then there is no biological equilibrium whenever and which implies that there exists a unique For the case when Z = 0, we have Since ] by If CTL immune response has not been established, we have cI 1 − b ≤ 0. So, we define the reproduction number for cellular immunity as follows where 1 b denotes the average life expectancy of CTL cells and I 1 is the number of infected cells at E 1 .
Hence, R Z 1 represents the average number of the CTL immune cells activated by infected cells.
Recapitulating the above discussions in the following theorem.

Global Stability
Regarding the global stability of the infection-free equilibrium E 0 , we have the following theorem.
Theorem 3. The infection-free equilibrium E 0 of system ( 2) is globally asymptotically stable when R 0 ≤ 1.
Proof.Based on the method proposed in [28], we construct the Lyapunov functional for system (2) at E 0 as follows For convenience, we let ϕ = ϕ(x, t) and ϕ τ = ϕ(x, t − τ) for any ϕ ∈ {u, w, v, z}.The time derivative of L 0 along the solution of system (2) satisfies From (7) and by applying Lemma 1, we get lim sup t→∞ T(x, t) ≤ T 0 .This implies that all omega limit points satisfy T(x, t) ≤ T 0 .Hence, it is sufficient to consider solutions for which T(x, t) ≤ T 0 .From the explicit formula of R 0 given in (11) and (H 0 )-(H 3 ), we get Hence, R 0 ≤ 1 ensures dL 0 dt ≤ 0. In addition, it can be shown that the largest compact invariant set in {(T, I, V, Z)| dL 0 dt = 0} is the singleton {E 0 }.Therefore, it follows from LaSalle's invariance principle [29] that E 0 is globally asymptotically stable when R 0 ≤ 1.
For the global stability of the two infection steady states E i of system (2), we suppose that R 0 > 1 and the incidence functions f and g satisfy for each infection equilibrium E i the following further hypothesis Therefore, we have the following result.
Theorem 4. Assume R 0 > 1 and (H 4 ) holds for each E i .
(i) The infection equilibrium without cellular immunity E 1 of system ( 2) is globally asymptotically stable The infection equilibrium with cellular immunity E 2 of system ( 2) is globally asymptotically stable if R z 1 > 1.
Calculating the time derivative of L 1 along the solution of system (2), we obtain Since the function f (T, I, V) is strictly monotonically increasing with respect to T, we have for It follows from (H 4 ) that Since H(ξ) ≥ 0 and R Z 1 ≤ 1, we have For (ii), we construct the Lyapunov functional as follows Calculating the time derivative of L 2 along the solution of system (2) and using and kη 2 I 2 = µV 2 , we have From ( 16)-( 18), we have dL 2 dt ≤ 0. Further, it is not hard to see that the largest invariant set in {(T, I, V, Z)| dL 2 dt = 0} is {E 2 }.By LaSalle's invariance principle, we deduce that E 2 is globally asymptotically stable.This ends the proof of Theorem 4.
Remark 1.The hypothesis (H 4 ) comes from (17) and (18).This hypothesis is a sufficient condition for that the time derivatives of the Lyapunov functionals L 1 and L 2 to be non-negative.When cell-to-cell mode is ignored (i.e., g ≡ 0), the assumption (H 4 ) can be reduced to which is verified by many types of the incidence rate including the bilinear incidence, the saturation incidence, the Beddington-DeAnglis functional response, the Crowley-Martin functional response and the Hattaf-Yousfi functional response.
In 2017, Xu et al. [30] proposed a PDE model with two discrete delays, cell-to-cell transmission and CTL immune response.They considered the spatial diffusion only in virus.This model is given by where the functions f and g satisfy the following properties: Choose the functions f (T, I, V) and g(T, I) as follows and g(T, I) = Clearly, f (T, I, V)V = β 1 T f (V) and g(T, I)I = β 2 T g(I) for all T, I, V ≥ 0. Based on f (V) V ≤ 0, g(I) I ≤ 0 and the last inequality of Lemma 3.1 in [30], it is not hard to prove that the above incidence functions f (T, I, V) and g(T, I) satisfy the five hypotheses (H 0 ) − (H 4 ).Therefore, the model and results investigated in [30] are extended and generalized.

Application and Numerical Simulations
In this section, we first apply our main results obtained in this study to the following model: where β 1 and β 2 denote, respectively, the virus-to-cell infection rate and the cell-to-cell transmission rate.The non-negative constants 1 and 2 measure the saturation effect.The other state variables and parameters have the same biological meanings as in models ( 1) and (2).Notice that system (21) extends the DDE model presented in [9] by introducing the spacial diffusion in both cells and viruses.Also, this system is a particular case of (2) with f (T, . As before, we consider system (21) with initial conditions and Neumann boundary conditions It is easy to check the first four hypotheses (H 0 )-(H 3 ).For the fifth hypothesis, we have Thus, the last hypothesis (H 4 ) is verified.By applying Theorems 3 and 4, we obtain the following result.1.If R 0 ≤ 1, then the infection-free equilibrium E 0 of system ( 21) is globally asymptotically stable.2. If R 0 > 1, then system (21) has two infection equilibria that are: (i) the infection equilibrium without cellular immunity E 1 that is globally asymptotically stable if R Z 1 ≤ 1; (ii) the infection equilibrium with cellular immunity E 2 that is globally asymptotically stable if R Z 1 > 1.
For the numerical simulations, we choose f i (τ) = γ i e −γ i τ for i = 1, 2. Clearly, ∞ 0 γ i e −γ i τ dτ = 1.Also, we consider the following new variables: Then the variables T, Y, I, U, V and Z satisfy the following system: The threshold parameters R 0 and R Z 1 for (24) are given by ( 11) and ( 15) with η 1 = γ 1 α 1 + γ 1 and . For the simplicity of numerical illustrations, we consider one-dimensional bounded spatial domain Ω = [0, 10] with d T = d I = d Z = 0.01 and d V = 0.02.Also, we consider β 2 and c as free parameters.All other parameter values are mentioned in Table 1.When β 2 = 1.5 × 10 −5 and c = 0.02, we have R 0 = 0.8539.By the first result given in Corollary 1, the infection-free equilibrium E 0 (719.4245,0, 0, 0) is globally asymptotically stable.This means that the virus is cleared, the infection dies out and the patient will be completely cured (see Figure 1).

Conclusions
In this article, we have proposed and investigated a generalized viral infection model with two infinite distributed delays, CTL immune response and spatial diffusion in both cells and virus.Also, the proposed model incorporated the classical virus-to-cell infection and the direct cell-to-cell transmission.Both routes of infection are modeled by two general incidence functions.Under some assumptions on these incidence functions, we have shown that the global dynamics of the model is completely determined by two threshold parameters that are the basic reproduction number R 0 and the reproduction numbers for cellular immunity R Z 1 .From the viewpoint of biology, we have proved that when R 0 ≤ 1 the infection-free equilibrium is globally asymptotically stable, which means that the virus is cleared and the infection dies out.Whereas, the virus persists in the host if R 0 > 1 and two steady states appear, one without cellular immunity which is globally asymptotically stable if R Z 1 ≤ 1 and the other with cellular immunity which is globally asymptotically stable if R Z 1 > 1.Hence, the activation of the CTL immune response is unable to eliminate the virus in vivo, but plays a fundamental role in the reduction of virus particles and infected cells.This last biological result can be easily deduced by comparing the components of virus particles and infected cells before and after the activation of cellular immunity.Since R 0 and R Z 1 have no relation to the diffusion coefficients d T , d I , d V and d Z , we conclude that the diffusion of cells and virus has no effect on the global stability of the three steady states of our PDE model with Neumann homogeneous boundary conditions.On the other hand, we have extended the models with ODEs [1,4,5], with DDEs [6][7][8][9] and with PDEs [17][18][19].Moreover, the more recent works presented in [30,31] are improved and generalized.

Table 1 .
List of parameters and their values used in numerical simulations.