Abstract
In this paper, an infection model with delay and general incidence function is formulated and analyzed. Theoretical results reveal that positive equilibrium may lose its stability, and Hopf bifurcation occurs when choosing delay as the bifurcation parameter. The direction of Hopf bifurcation and the stability of the periodic solutions are also discussed. Furthermore, to illustrate the numerous changes in the local stability and instability of the positive equilibrium, we conduct numerical simulations by using four different types of functional incidence, i.e., bilinear incidence, saturation incidence, Beddington–DeAngelis response, and Hattaf–Yousfi response. Rich dynamics of the model, such as Hopf bifurcations and chaotic solutions, are presented numerically.
1. Introduction
Mathematical modeling has been proven to be valuable in exploring mechanisms and dynamical behaviors of the viral infection process. The analysis of these models can provide insights into developing effective control strategies for infections and evaluating antiviral therapies. Based on the basic virus infection model introduced by Nowak et al. [1], various models have been formulated by many authors to describe the dynamics of virus population qualitatively and quantitatively, and lots of interesting phenomena have been discussed, such as the global stability of the equilibria, bifurcations, periodic oscillations, limit cycles, and the effects of time delays.
For models of virus infections, it is observed that the time delays introduced into them can predict the viral infection process better when compared to models without delays, such as the time between viral entry into a target cell and the production of subsequent viral particles, the time necessary for the newly produced virus to become mature and then infectious particles, and that needed to activate the immune response, all of which cannot be ignored when describing the interaction of them [2, 3]. In the field of virus dynamics, many models with discrete or distributed delays have been studied, showing that rich dynamics including stability switches, Hopf bifurcation, and chaotic oscillations can be generated by the mechanisms of time delays. Specifically, stability switches can occur due to a delayed immune response and exponentially decayed delay-dependent parameters [4–9], which leads to the switches in the stability of equilibria because of the change in the value of parameters related with these delays.
It is known that the function forms of the incidence rate of the infection have a crucial role in the modeling of the virus dynamics, which are important in determining qualitative behavior of the proposed models and in giving a reasonable description of the dynamics. Huang et al. [10] and Zhang and Xu [11] formulated a virus dynamical model with Beddington–DeAngelis (BD) infection rate , where x represents the concentration of uninfected cells and υ represents that of virus, with being constants. This function is similar to the well-known Holling type II functional response, and the term in the denominator reflects the mutual interference between viruses. In [12], Zhou and Cui used a Crowley–Martin (CM) function response of the form , with , which are constants. It can be seen that in BD and CM functional responses, when , these two functions are simplified to the Holling type II functional response. And when , they express a saturation response. Moreover, when , they are the mass action process (or Holling type I functional response). Other more general incidence rates have also been proposed. McCluskey and Yang used function as incidence rate under some biologically motivated assumptions to study the global stability of a virus model [13]. Hattaf et al. in [14, 15] used an incidence rate that covers a variety of incidence functions such as BD response when and CM response when if , which will be called Hattaf–Yousfi response in the following, with y being infected cells.
In this paper, we will study the influence of general incidence function and time delay on the dynamical behaviors of an infection model. Motivated by the works of [7, 13], we will focus on the dynamics of a delayed model that incorporates the immune response in antiviral defense, time delay in activating the immune response for the virus, and two general incidence functions for the transmission of virus-to-cell and cell-to-cell, respectively. The main aim of this paper is to study stability switches of the positive equilibrium, which are generated by the mechanism of time delay, implying the existence of changes in local stability and the instability of this equilibrium. By using the characteristic equation with delay-dependent parameters, normal theory, and center manifold theorem, we can indicate the existence of pure imaginary roots and then Hopf bifurcation. Furthermore, to verify the theoretical results, numerical simulations are performed for four different forms of incidence functions, including bilinear incidence, saturation incidence, Beddington–DeAngelis response, and Hattaf–Yousfi functional response, respectively.
The organization of the paper is as follows. In the next section, we introduce the virus dynamic model. The boundedness and nonnegativity of solutions of the model as well as the existence and uniqueness of equilibria are discussed in Section 3. In Section 4, both the stabilities of two equilibria and the conditions for the existence of Hopf bifurcation are studied. Furthermore, numerical simulations are presented in Section 5. A brief summary is given in Section 6, and the properties of the Hopf bifurcation solutions have been provided in Appendix.
2. Model Formulation
To account for the possible effect of the latent period in viral infection, we introduce a constant time delay τ to represent the time needed to activate the immune response. In this model, two transmissions of virus-to-cell and cell-to-cell have been both taken into consideration, and we hope to employ the function of incidence rate given by more general form. Consequently, the recruitment of infected cells is given by two functions, i.e., and , which were used in [16–18], here u, , and denote the concentrations of uninfected cells, infected cells, and free virus particles. z is used to denote the concentration of CTL cells. Assume that the uninfected cells are produced at a constant rate s. , , , and represent the death rates of uninfected cells, infected cells, free virus, and CTL cells. p denotes the killing rate of infected cells by CTL cells. k and c are the production rates of free viral particles and the effector cells, respectively. Then, we construct the following model:with initial conditionswhere is the Banach space of continuous functions mapping the interval into , with (), and all constant parameters are positive.
Function () is continuously differentiable and is assumed to satisfy the following conditions:
These fundamental assumptions are biologically motivated, and it is easy to check that class of functions (and ) satisfying these hypotheses include incidence functions such as (mass incidence), (saturation incidence), (Beddington–DeAngelis response), and (Hattaf–Yousfi response).
3. Preliminaries and Equilibria
We denote the state space of (1) by , equipped with the sup-norm . For any , the existence and uniqueness of the solutionof model (1) with initial condition (2) follow from the standard theory of functional differential equations [19]. Furthermore, X is positively invariant for (1). In the following, we discuss the boundedness [20] and nonnegativity of solutions of model (1) as well as the existence of equilibria.
Theorem 1. The solutions of model (1) with initial condition (2) are nonnegative and ultimately uniformly bounded for all .
Proof. For , if there exists such that for and , then it is obvious that , which implies that for sufficiently small , for , contradicting with for all . Thus, it is valid that for all . In order to show , , and for all , assume that there exists such that = 0 for the first time. Firstly, if , then from the second equation of (1); we haveand then there exists small enough such that for , which contradicts with the facts that for and . Thus for all . For and , sinceit is obvious that or contradicts with the above expression about and at . Therefore, the solution with initial condition (2) is nonnegative for all .
To prove the boundedness of the solutions, letting , we havewhere , then , implying that and are ultimately bounded. For , it follows from the third equation of (1) thatwhich leads to . Similarly, we have . Therefore, the solutions of system (1) are ultimately uniformly bounded.
System (1) always has an infection-free equilibrium , where . In the following, we analyze the existence and uniqueness of the positive equilibrium .
For convenience, we denoteDefine the basic reproductive number asFrom biological viewpoint, can be divided into two parts, with measuring the average number of secondary infected generation caused by an existing free virus and being that caused by an infected cell, which give the basic reproduction number corresponding to virus-to-cell and cell-to-cell infections, respectively.
It is easy to obtain and . By the first and second equations in (1), we havewhich gives . In order to get the infected steady state, it is necessary that , and then, we must have , with and .
From the equation of , it follows thatDenote the left side of (11) as , and then the positive equilibrium of model (1) are given by for . Noting that andfor when , and it is sufficient to show that if there is a positive equilibrium. In fact, when , we haveand then is positive for sufficiently small ; therefore, there exists an infected steady state .
Now we examine the derivative of at . Noting that , and thenand thus, we havewheredenoting the sum of negative terms, and thus function G is strictly decreasing at each of its zeros. Suppose there exists more than one infected steady state, then there must exist an equilibrium such that , which contradicts with the above fact. Therefore, there is only one positive equilibrium when the basic reproductive number , and it is easy to verify that there exists no infection equilibrium when . Then the following result is obtained.
Theorem 2. For system (1), when , there exists only the infection-free equilibrium . If , there is a unique infection equilibrium as well as .
When , by using the persistence theory in [21], it can be shown that system (1) is uniformly persistent.
4. Stability and Hopf Bifurcation
In this section, we will study the stability of the two equilibria and the conditions for the existence of Hopf bifurcation.
4.1. Stability of
Firstly, we linearize model (1) at the steady states to discuss its local stability. The Jacobian matrix leads to the following characteristic equation:and thus, the linearization of system (1) at can be expressed by
It is clear that (21) has two negative real roots and . The stability of is determined completely by using the following equation:
For the case , the roots of (22) have only negative real parts. Then the infection-free equilibrium is locally asymptotically stable. Otherwise, when , equation (22) has at least one root with positive real part, which implies that is unstable.
Theorem 3. If , then the infection-free equilibrium of model (1) is globally asymptotically stable, and it is unstable when .
Proof. For a continuous and bounded function , we defineThen for any solution of (1), we haveBy using the fluctuation lemma in [22], we know that there is a sequence with such thatSubstituting the sequence into the first equation of (1) and taking the limit givesand then we have . A similar argument to the second and third equations of model (1) yieldsandBy the assumption of (4), combining equality (28) into (27) giveswhich leads to . Noticing that is nonnegative since it is the supremum of the function , then there are possible cases of or . If , then we have , which gives a contradiction. Therefore, is valid. From (28), we have . When considering the fourth equation about in model (1), we obtain . As for , we can find that by applying the limit theory to the first equation in system (1), and thus, the proof is completed.
4.2. Stability of
For convenience, we introduce the following notations:
From the assumption (5), it follows that .
In the following calculation, we will use the following equilibrium equations:
Then, we have
Linearized model (1) at the infection equilibrium and the characteristic equation can be expressed aswhere and with
It can be seen that and are all positive. When the delay , and are independent with τ and the characteristic equation (33) can be reduced to
By using the Routh–Hurwitz criterion, we know that all solutions of (35) have negative real parts if and only if the following conditions are satisfied:
Let and assume
Theorem 4. If (37) holds and , then the infection equilibrium is globally asymptotically stable when .
Proof. Construct a Lyapunov functional asand then calculating the time derivation of along system (1) leads towith . By using (4), (37), and the property of function , we have and , which implies that and it can verified that is the largest compact invariant set where . Therefore, is globally asymptotically stable by using LaSalle’s invariance principle.
However, if the delay , it may destabilize the infected steady state and lead to Hopf bifurcation, which will be discussed in the following subsection.
4.3. Hopf Bifurcation
It is shown that all roots of (33) locate in the left side of the imaginary axis if and (36) is satisfied, which leads to local stability of the infected equilibrium . When τ increases from 0 to variant positive delays, it is possible that the roots of (35) pass through the imaginary axis and enter the right side in the complex plane. Therefore, the stability switch may occur, and it is necessary to study the transcendental equation (33) when , which is a critical value under small perturbation. In the following analysis, we will study the occurrence of any possible stability switching in this case.
Substituting into (33) and separating the real and imaginary parts give
Squaring and adding the two equalities lead towhere
Note that equation has at least one positive root when , i.e., . And let , then the characteristic equation (33) has a purely imaginary root which is equivalent to has a positive real root z. Take the transformation ; then for , it is equivalent to with , , and it has roots of , correspondingly, (); thus we have the following result about the roots of [23].
Lemma 1. Denoting , for the roots of ,(i)If , then there exists at least one positive root(ii)If and , then there are positive roots if and only if and (iii)If and , then there are positive roots if and only if there exists at least a positive such that Now assume that has four roots with ; then by the relation of , we have the fact of . From the equalities in (40), we can express and asThen for , we can obtain the following expression of delay :
Lemma 2. For the roots of characteristic equation (35), suppose conditions in (36) hold, then there are two possibilities:(i)Equation (35) has only roots with negative real parts for if any one of the following conditions ()() holds(ii)Equation (35) has only roots with negative real parts for if the conditions ()() are not satisfied, with () () , , and () , , there is a positive such that From this lemma, we can see that the infection equilibrium is asymptotically stable for all if the conditions ()() are not satisfied. Otherwise, if one of these three conditions is satisfied, then the infection equilibrium is asymptotically stable for , and a Hopf bifurcation can occur at this equilibrium when with being the critical value of .
Theorem 5. Suppose (36) holds, when (correspondingly for some ), then the characteristic equation (35) admits a pair of simple conjugate pure imaginary roots and , which crosses the imaginary axis from left to right (from right to left) if , where
Proof. Differentiating the characteristic equation (33) with respect to delay τ and arranging it can give the expression of asMeanwhile, noticing the factand then some computation leads towhich gives the desired result.
When model (1) undergoes stability switch at for , the conditions and theorem for the direction and stability of Hopf bifurcation can be discussed and proven by the normal theory and center manifold theorem in [24], which are presented in Appendix.
5. Numerical Simulation
In model (1), the incidence rates for transmission of virus-to-cell and cell-to-cell are taken as the general form. By choosing four specific types of functions, i.e., bilinear incidence, saturation incidence, Beddington–DeAngelis response, and Hattaf–Yousfi response, we conduct numerical simulations to investigate the complex dynamics that model (1) can have, when taking time delay τ as the bifurcation parameter. The values of parameters, , , , and , remain the same in the following different cases. Case 1. Consider the functions of and , then the incidence rates have the bilinear forms of and . We take the other parameters as , , , , and . Figure 1 presents the different results for varying producing rate of uninfected cells, plotted with (Figure 1(a)) and (Figure 1(b)), which shows the effect of parameter s on dynamical behavior of the model in this case. In this figure, we can see that there exists periodic solutions bifurcated from the infection equilibrium. The maximal and minimal values of are denoted by the two curves, and the line implies the local stability of . It is clearly shown that the value of increases with that of s. Figure 2 shows that the infection equilibrium is asymptotically stable when and a periodic solution exists when , which corresponds to the bifurcation diagram in Figure 1(b). Case 2. The type of saturation incidence is used in this case, i.e., and ; then the dynamical behaviors of the model are simulated in Figures 3 and 4, with the parameters taken as , , and . In Figures 3(a), we set , , , , and , and in Figure 3(b), we set , , , , and . By the set of parameters, Figure 3(a) shows that the time delay does cause a bifurcation at . The line is given by , and the infected steady state is locally asymptotically stable when . More complicated dynamics for model (1) are caused by the delay perturbation in Figure 3(b). A period-three solution occurs for and , which are displayed clearly by time series of for the time interval and their corresponding phase portraits of and plotted in Figure 4, respectively. Case 3. Beddington–DeAngelis response is used in this case; here, and . By this type of incidence function, chaotic motions also occur when increasing τ from 0 to 50 with , , , , , and . When altering parameter values of and , it can lead to the different bifurcation diagrams, which can be seen in Figures 5 and 6. We take , and , , respectively, and then the simulation results in Figure 5 show the effect of these two parameters on the system. From Figure 5, we can see that model (1) with these parameter values can experience complicated dynamics. Furthermore, for in Figure 5(b), the solution of the model is presented by its orbit, plot of time series, and phase portrait correspondingly in Figure 6. For , there exists a Hopf bifurcation, and when it increases to , chaotic motion occurs, shown in Figures 6(a) and 6(b). The plot of time series of and phase portrait of and are given to illustrate the dynamical behaviors when in Figures 6(c) and 6(d), showing that there exists a period-five solution. Case 4. When Hattaf–Yousfi response is used in this case, the functions take the form of and , respectively. If the parameter values are chosen as , , , , , , , , and , only Hopf bifurcation exists and no chaotic motions occur when , shown in Figure 7(a). For a different set of values, let , , , , , , , , and , then the resulting bifurcation diagram in Figure 7(b) shows that system (1) can have solution of period three. Figure 8 gives the time series of and phase portrait of and for to illustrate the solution of system (1), which corresponds to Figure 7(b).
(a)
(b)
(a)
(b)
(c)
(d)
(a)
(b)
(a)
(b)
(c)
(d)
(a)
(b)
(a)
(b)
(c)
(d)
(a)
(b)
(a)
(b)
6. Discussion
To model the mechanisms of infectious diseases mathematically and explore dynamical behaviors of infection processes, some types of functional incidence have been used in [10, 12–15], which play an important role in determining qualitative behaviors of the proposed models and in giving reasonable descriptions of the dynamics. General incidence rate has been introduced into many epidemic models with an aim to include different situations as much as possible. In our model, we incorporate the effect of time delay needed to activate the immune response for the virus and two general incidence functions for the transmission of virus-to-cell and cell-to-cell. The basic reproductive number, which works as a crucial threshold and determines the dynamics of the model, has been defined as the sum of two parts related with infections of virus-to-cell and cell-to-cell, respectively.
In this paper, we are concerned with the stabilities of two equilibria and the existence of Hopf bifurcation through which the positive equilibrium loses its stability and periodic solutions occur. In the analysis, the time delay is chosen as the bifurcation parameter, which can destabilize the positive equilibrium when it increases. By using the characteristic equation with delay-dependent parameters, normal theory, and center manifold theorem, the existence of pure imaginary roots is verified and then the system experiences Hopf bifurcation. To illustrate the dynamical behavior of stability switches, simulations are conducted to show the process numerically when taking general incidence function as four specific types of bilinear incidence, saturation incidence, Beddington–DeAngelis response, and Hattaf–Yousfi response. In each case, the dynamical behavior is simulated when increasing τ from 0 to 50, and two sets of parameters are taken to compare the bifurcation results. Furthermore, for some fixed values of τ, we also present the plot of time series, phase portrait, or solution orbit to show the complicated dynamics. The effect of delay perturbation and different forms of incidence can be seen in the figures.
It should be noted that when other factors are taken into consideration for more realistic mechanism, such as the delays describing the intracellular latency for virus-to-cell infection and cell-to-cell infection, or varying producing rate of the uninfected cells instead of constant, the model can be extended reasonably and the dynamical analysis will become more challenging than the present one, which we may discuss in future work.
Appendix
Theorem 6. For system (1),(i)The direction of Hopf bifurcation is determined by the sign of , i.e., it is a supercritical bifurcation when and a subcritical bifurcation when .(ii)The stability of the bifurcated periodic solution is determined by , i.e., the periodic solution is stable when and unstable when .(iii)The periodic of bifurcated periodic solutions is determined by , i.e., the period increases when and decreases when , with
Proof. In order to rewrite model (1) as a functional differential equation in , we take the following transformations for variables and time scale:Then model (1) becomes the following form:where , , ,withandwhere andBy Riesz representation theorem, there exists a matrix components in , which is a bounded variation function such that for . As for (A.3), has the expression ofwhere δ denotes the Dirac delta function. For , we defineandThen, model (1) is equivalent toFurthermore, for , defineand a bilinear inner productwhere . and are adjoint operators with eigenvalues , and then we discuss the eigenvector of and corresponding to these two eigenvalues, respectively. First, for , if is the eigenvector corresponding to , then it leads to , i.e.,Denoting the above equation as and being the element of in row and column, then we haveSimilarly, for , if is the eigenvector of corresponding to , we can obtainwhere is the element of matrix B in row and column andSince , we can chooseIn the following, by using the notations in [24], we compute the center manifold when . Correspondingly, assuming be the solution of (A.6), and define and asThen on the center manifold , we haveandwith z and being local coordinates for in the direction and . Thus,When considering real solution of (A.6), we can obtainLetOn the other side, note that . Then by using the expression of in (A.4), we haveSubstituting the expressions of and into (A.12) and comparing the coefficients of , , , and lead towhereNext, in order to determine the value of , we compute and briefly. According to the results in [24], we havewhere is a constant vector. Furthermore, it can be deduced that satisfieswherewithBy using Cramer’s rule, we havewhere is the matrix that the column in is substituted by . Thus , , and then can be determined. Furthermore, the values of , , , and can be computed by using (A.1).
Data Availability
No data were used to support the article.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
This research was supported by the National Natural Science Foundation of China (nos. 11501443 and 11701445) and Natural Science Basic Research Plan in Shaanxi Province of China (no. 2019JM338).