key: cord-0932566-1la05cuu
authors: Liu, Chunhua; Kong, Lei
title: Dynamics of an HIV model with cytotoxic T-lymphocyte memory
date: 2020-10-17
journal: Adv Differ Equ
DOI: 10.1186/s13662-020-03035-8
sha: 7b605271b05a15bc83b115e36a61b7ddfc6a6608
doc_id: 932566
cord_uid: 1la05cuu

We consider a four-dimensional HIV model that includes healthy cells, infected cells, primary cytotoxic T-lymphocyte response (CTLp), and secondary cytotoxic T-lymphocyte response (CTLe). The CTL memory generation depends on CD4(+) T-cell help, and infection of CD4(+) T cells results in impaired T-cell help. We show that the system has up to five equilibria. By the Routh–Hurwitz theorem and central manifold theorem we obtain some sufficient conditions for the local stability, globally stability of the equilibria, and the bifurcations. We still discover the bistability case where in the system there may coexist two stable equilibria or a stable equilibrium together with a stable limit cycle. Several numerical analyses are carried out to illustrate the validity of our theoretical results.

Human immunodeficiency virus (HIV) and acquired immunodeficiency syndrome (AIDS) have become global health problems; there were estimated 38.0 million people living with HIV at the end of 2019 [1] . HIV can be transmitted via the exchange of a variety of body fluids from infected people, such as blood, breast milk, or semen and vaginal secretions, but it is known that current antiretroviral drugs cannot enucleate HIV from the body. Mathematical models have been formulated for various epidemiological diseases like novel coronavirus [2] [3] [4] or malaria [5] . Models of HIV as a chronic infectious disease have been investigated in many papers. Some models focus on the size of infection, starting from the population models [6] [7] [8] , whereas other focus on cell infection, starting from the virus models, which have become an important tool in both understanding HIV-1 infection in host and providing valuable insight into HIV pathogenesis. A well-known model for HIV infection is the following system [9, 10] : where variables x, y, and v represent the densities of the healthy cells, the infected cells, and the virus at time t, respectively. Here a mass action infection mechanism is adopted. The parameters β and λ stand for the infection and constant growth rates of the healthy cell, respectively. System (1.1) and its variations have been investigated in many papers [11] [12] [13] [14] [15] [16] . Stephen et al. [17] recently added the term ρTV C+V incorporating the homeostatic proliferation of T-cells, which leads to interesting dynamic results, such as bistability and Hopf bifurcation.

It is well-known that it takes a long period for an HIV to become an AIDS, and in the medical literature [18] , it was pointed out that the latent reservoir (i.e., latent infection) was the main obstacle to eradicate the virus. Therefore the four-dimensional mathematical model including the latent infection seems more reasonable [19, 20] . In recent years, latent cells were considered in many models, such as Beddington-DeAngelis function response with delay [20] , Crowley-Martin function response [21] , and general infection function with CTC and VTC transmission [22] .

To recover from a viral infection, the cytotoxic T lymphocyte (CTL), which can clear away the infected cells to prevent further viral replications, plays a particularly important role. In 1996, Nowak and Bangham [9] proposed the well-known model with immune response:

where the variable z represents the concentration of CTLs. Many authors have studied the infective models with different immune responses, such as lytic and nonlytic immune responses [23, 24] , cell-mediated immune mechanism or humoral immune mechanism [25] [26] [27] , delayed immune response with drug therapies [28] , and general CTL immune response with silent infected cell-to-cell spread [29] .

However, after a viral infection, the CTLs that are responsible for clearing away the infected cells become cytotoxic T-lymphocyte precursors (CTLp) and have receptors for detecting the virus from the previous infection [30] . Upon contacting with the virus during a subsequent infection, the precursors differentiate and become cytotoxic T-lymphocyte effectors (CTLe), and these cells are again responsible for clearing away the invading virus. Considering this infective mechanism, Wodarz et al. [31, 32] provided the following model with CTL response:

where the healthy cells x and the infected cells y are described similarly as in system (1.1). Instead of just one class of CTL response, the CTLp and CTLe are introduced. The CTLp and CTLe are represented by w and z. These precursors emerge at rate cyw and may become effectors at rate cqyw or cleared away naturally at rate bw. Similarly, the effectors are created at rate cqyw and cleared at rate hz.

In model (1.2) , there is no virus term whose population is assumed at a quasi-steady state, which is proportional to infected cells. Model (1.2) is completely analyzed by Bernard et al. [33] , who have found that the system transforms from one equilibrium to the next as the basic reproductive number R 0 increases. When R 0 increases further, they show that periodic solutions may arise from the third equilibrium via Hopf bifurcation.

In fact, another model with cytotoxic T-lymphocyte memory proposed in [31, 32] is given by

This model assumes that the target cells are CD4 + T cells; moreover, it includes the additional feature that expansion of the CTLp population is proportional to both antigen and the number of uninfected CD4 + T-cells capable of delivering T-cell help. The memory generation depends on CD4 + T-cell help, and infection of CD4 + T-cells results in impaired T-cell help. We also assume that differentiation into effector functions is independent of CD4 + T-cell help [31] . A detailed explanation of the model can be found in [31] . All the parameters are positive.

Dynamics of system (1.3) is numerically analyzed in [31, 32] . In this paper, we provide a rigorous analytical method of system (1.3), and the basic framework is as follows. In Sect. 2, we establish the well-posedness of the model including nonnegativity and boundedness of the solutions, the existence of equilibria, and local stability of the boundary equilibria. The local stability analysis of the positive equilibria and their bifurcations are presented in Sect. 3. Numerical illustrations are given in Sect. 4. Finally, we discuss both mathematical and biological perspectives of the findings in Sect. 5.

For mathematical simplicity, we do some rescallings in system (1.3) .

After changing back to the origin variables x, y, w, z, t, the scaled system is given by

where the horizontal lines on the heads of these parameters are removed, and the parameters d, a, c, q, b, h are replaced by d, a, c, q, b, h. Obviously, all the parameters are positive. The basic reproductive number of model (1.3) is R 00 = λβ ad , and for system (2.1), R 00 becomes R 0 = 1 ad .

Theorem 2.1 All solutions (x(t), y(t), w(t), z(t)) of system (2.1) are nonnegative for t > 0. Moreover, if x(0) ≥ 0, y(0) ≥ 0, w(0) ≥ 0, and z(0) ≥ 0, then all solutions of system (2.1) are ultimately bounded.

Proof By variation of constants we find the following solutions of (2.1):

which proves the nonnegativity of solutions of system (2.1). Note that the first equation of (2.1) implies dx dt ≤ 1dx. The solution is given by

Adding the first two equations of (2.1), we obtain

where d 1 = min{d, a}. It has the solution x + y ≤ (x(0) + y(0))e -d 1 t + 1 d 1 , which implies lim sup t→∞ (x(t) + y(t)) ≤ 1 d 1 , and thus x(t) and y(t) are bounded. Supposing that z is unbounded, by the second equation of (2.1) we have lim t→∞ y(t) = 0, which implies lim t→∞ w(t) = 0 from the third equation of (2.1). Then we get lim t→∞ z(t) = 0 from the fourth equation of (2.1), which contradicts with the unboundedness of z. Thus z must be bounded. Lastly, assume that w is unbounded. Based on the boundedness of z and the fourth equation of (2.1), we obtain lim t→∞ y(t) = 0, and from the third equation of (2.1) it follows that lim t→∞ w(t) = 0, which causes a contradiction. Hence w is bounded. The proof is complete. Theorem 2.1 shows that there exists a bounded positive invariant region ⊂ R 4 + for the system. Thus we concentrate on to discuss the dynamics. In fact, the infectionfree equilibrium E 0 = ( 1 d , 0, 0, 0) always exists, and there exists an infectious equilibrium without CTL, E 1 = (a, d(R 0 -1), 0, 0) if R 0 > 1. To find the infectious equilibrium with CTL, it suffices to solve the system If x 1 , x 2 are solutions of (2.4), then x 1 x 2 = q cd > 0 and x 1 + x 2 = -b-c-dq cd . Hence, to ensure the existence of a positive solution, we have b < c + dq.

The determinant of (2.4) is

Note that c + dq -2 cdq ≥ 0. Combining this with the condition b < c + dq, we get that ≥ 0 only if 0 < b ≤ c + dq -2 cdq, and y > 0 and z > 0 imply x > q c and x > a by (2.3) . Set

The positive equilibria of model (2.1) are classified by the sign of . Let us consider three cases.

> 0, which is equivalent to 0 < b < c + dq -2 cdq. Here we consider three cases by the sign of f (a).

1. f (a) < 0. In this case, > 0 is obviously satisfied, and a > q c is equivalent to R 1 > 1, and thus R 1 < R 0 by (2.6). Combining this with R 1 > 1, we get that model (2.1) admits a unique positive equilibrium

2. f (a) = 0. In this case, R 0 = R 1 , which indicates that a is a solution of f (x) = 0, and the other root is x * = q acd . To get a positive equilibrium of model (2.1), x * = q acd > a and a > q c are required, which yields q > a 2 cd and ac > q. Combining this with > 0 (if and only if 0

3. f (a) > 0. Let us consider two cases.

(1) q c ≥ a, which is equivalent to R 1 ≤ 1. To ensure the positive equilibrium of model (2.1), the following conditions are required:

(2.7)

Note that f ( q c ) > 0 always holds, and > 0 if and only if 0

Combining this with c > dq, we have ad < 1 to ensure that the intersection is nonempty, and hence model

(2) q c < a, which is equivalent to R 1 > 1. Model (2.1) admits two positive equilibria if the following conditions are satisfied:

(2.8)

if and only if q < a 2 cd, and hence model (2.1) admits two positive equilibria if one of the following conditions holds:

= 0 if and only if b = c + dq -2 cdq, and model (2.1) admits a unique positive equilibrium if

Combining these two conditions with = 0 which is equivalent to b = c + dq -2 cdq, the following results are obtained:

Case III If < 0, then there is no positive equilibrium, because f (x) = 0 has no real root. The following theorem summarizes all positive equilibria of system (2.1).

(3) Model (2.1) admits a unique positive equilibrium E 20 = (x 20 , y 20 , w 20 , z 20 ) if b = c + dq -2 cdq, q > a 2 cd, and c > dq.

Next, we discuss the local stability and global stability of the boundary equilibria. The stability of the equilibria is based on the Jacobian matrix of model (2.1):

The characteristic roots at the equilibrium E 0 are given by

and

Hence the local stability of equilibrium E 0 is given by the following theorem.

Proof To show the global stability of equilibrium E 0 , we use the method of fluctuation lemma employed by Hirsch et al. [34] [35] [36] . We introduce some notations. For a continuous bounded function g : [0, ∞] → R, let

By Theorem 2.1 the solutions x(t), y(t), w(t), and z(t) are always nonnegative and bounded for any nonnegative initial conditions, and the limits lim sup t→∞ and lim inf t→∞ always exist for each of these solutions. By the fluctuation lemma there exists a sequence t n such that

(2.11)

By a similar argument, for the remaining three equations in model (2.1), we get

(2.14)

We claim that y ∞ = 0. Otherwise, y ∞ > 0; then it follows from (2.11) and (2.12) that

Therefore R 0 ≥ 1, which contradicts the condition R 0 < 1, and thus y ∞ = 0. And from equations (2.13) and (2.14) we get that w ∞ = 0 and z ∞ = 0. The conditions y ∞ = 0, w ∞ = 0, and z ∞ = 0 imply that y(t) → 0, w(t) → 0, and z(t) → 0 as t → ∞, respectively. Thus from the first equation of model (2.1) we have the asymptotic differential equationẋ = 1dx. By simple calculation we get the solution

d , which clearly shows that the solution x(t) → 1 d as t → ∞ by the theory of asymptotically autonomous systems. The proof is complete.

Letting (s + d + y 1 )(sx 1 + a) + x 1 y 1 = 0, we have s 2 + (d + y 1 )s + ay 1 = 0. Note that s 1 + s 2 = -(d + y 1 ) < 0 and s 1 s 2 = ay 1 > 0, which imply that s 1 and s 2 have negative real parts. Obviously, s 3 = -h has a negative real part. In order for all the roots of equation (2.15) to have negative real parts, it is required that

Theorem 2.5 Equilibrium E 1 is locally asymptotically stable for 1 < R 0 < R 1 or R 0 > 1 > R 1 and is unstable for

Remark Note that f (a) = ab (R 1 -1) (R 1 -R 0 ) = -as 4 , from which it follows that E 1 is locally asymptotically stable, whereas system (2.1) may have positive equilibrium under certain conditions by Theorem 2.2 (case I.3). Moreover, system (2.1) may have both stable positive equilibrium and stable boundary equilibrium E 1 .

By an easy calculation the characteristic equation of the positive equilibrium follows: It has a characteristic root s 1 = 0, and s 2 , s 3 , s 4 are determined by the equation where the horizontal lines on the heads of these letters are removed, and we still denote x, y, w, z by x, y, w, z. The third equation of model (3.7) has no linear term, and then the center manifold is a curve tangent to the w-axis.

To obtain an approximative expression of the center manifold, we set

where m 1 , n 1 , p 1 , m 2 , n 2 , p 2 are undetermined coefficients. It follows that

To find the unknown coefficients m 1 , m 2 , n 1 , n 2 , p 1 , p 2 , we substitute (3.7) and (3.8) into (3.9), compare the coefficients at w and w 2 , and obtain For the equilibria E 2-and E 2+ , we have the following properties. Proof Similarly to the analysis of equilibrium E 20 , we obtain α 4 

which is strictly increasing in (-∞, 1 d ). Suppose x 2+ < 1 d . Submitting x 2+ into this inequality and simplifying it, we have c + bdq > √ . Note that c + bdq > 0. Squaring both sides and simplifying, we obtain 4bc > 0, which implies x 2-< c+dq-b 2cd < x 2+ < 1 d . An easy calculation shows that g( c+dq-b 2cd ) = -4c by (3.3). Then g( c+dq-b 2cd ) < 0 when > 0. It follows that g(x 2-) < g( c+dq-b 2cd ) < 0, which yields det(J(E 2-)) = α 4 = y 2-w 2x 2-g(x 2-) < 0. Solving g(x) = 0, we obtain the roots x + = c+ √ bc

bc > 0, after squaring both sides and simplifying, we obtain (bcdq) 2 

x 2+ g(x 2+ ) > 0. This completes the proof.

Proof Since the determinant det(E 2-) < 0, there is at least one characteristic root that has no negative real part. Therefore E 2-is unstable.

Next, we discuss the stability of equilibria E 2+ . By (3.2) the characteristic equation at E 2+ is given by

Note that α 1 > 0, α 2 > 0, and α 3 > 0, and according to the proof of Propositiony 3.1, we know that α 4 > 0. The relevant Routh-Hurwitz determinants are

(3.12)

Note that 1 > 0 and the sign of 4 is the same as 3 . By the formulas of x 2+ , y 2+ , and z 2+ , 2 and 3 can be written more explicitly as

Next, we give the following lemma to show that if both 2 and 3 can become zero, then 3 will cross zero before 2 does. Remark The proof is similar to that in [33] , so here we omit it.

Thus, to consider the stability of the positive equilibrium E 2+ , we only need to consider the possibility of 3 = 0 (see [33] ). Note that x, y, z do not contain h and 3 = 0 is just a quadratic equation in terms of h. Let

(3.14)

Obviously, h * 1 < h * < h * 2 , and thus we have the following theorem.

(1) If C 3 < 0, then E 2+ is locally asymptotically stable when h ∈ (h * 2 , +∞) and is unstable when h ∈ (0, h * 2 ). 

Remarks 1. The proof is straightforward by considering the sign of the quadratic polynomial A 3 h 2 + B 3 h + C 3 , and the detailed proof is contained in Theorem 3.4, so we omit it.

2. If we use the computing method of paper [33] , then 2 and 3 can be written more explicitly as

).

Obviously, A 22 > 0, B 22 > 0, C 22 > 0, A 33 > 0, B 33 > 0, and C 33 > 0, which indicates that 22 > 0 and 33 > 0 as long as h > b. As in [33] , the infectious equilibrium E 2+ with CTL response is always stable if the death rate of the CTLe is higer than that of the CTLp. Proof As in [37] , a four-dimensional model has a Hopf bifurcation if 1 > 0, 2 > 0, and 3 = 0. Obviously, 1 > 0 always holds by (3.12) , and according to Lemma 3.1, 2 > 0 as 3 crosses zero for some change of parameters. Therefore a Hopf bifurcation occurs as 3 = 0, which is equivalent to finding the roots of the quadratic polynomial equation (3.13) and (3.14) . Hence

exists, which indicates that E 2+ is stable for h ∈ (0, h * 2 ) and E 2+ is unstable for h ∈ (h * 2 , +∞), whereas 3 = 0 if h = h * 2 . Hence by [37] a Hopf bifurcation occurs at h = h * 2 . However, if C 3 > 0 and B 3 < 0, then 3 > 0 as h ∈ (0, h * 1 ) ∪ (h * 2 , +∞), which implies that E 2+ is stable, whereas 3 < 0 as h ∈ (h * 1 , h * 2 ), which ensures that E 2+ is unstable, and 3 = 0 as h = h * 2 or h = h * 1 . Hence by [37] two Hopf bifurcations occur at h = h * 2 and h = h * 1 . If C 3 = 0 and B 3 < 0, then 3 > 0 as h ∈ (-B 3 A 3 , +∞), which implies that E 2+ is stable, whereas 3 < 0 as h ∈ (0, -B 3 A 3 ), which ensures that E 2+ is unstable, and 3 = 0 as h = -B 3 A 3 . Thus by [37] 

Remark Since the expressions of x, y, w, z are very complicated, it is too complicate to directly discuss the sign of 4 = B 2 3 -4A 3 C 3 . In the next section, by numerical calculations we will demonstrate that the three cases 4 > 0, 4 = 0, and 4 < 0 are possible.

In this section, we demonstrate the theoretical results by numerical simulations. For convenience, we will work on the scaled model (2.1) instead of the original model (1.3) .

The values of λ used for simulations in [32] are λ = 1 and λ = 10, and as a bifurcation in [33] , this indicates that the values of λ have a fairly large variation; therefore we take λ = 0.5657. Note that β = 0.5 and β = 0.001 were used in [32] , whereas β = 3 400 was used in [33] , and hence choosing β = 0.0707 is reasonable. The parameter q ∈ [0, 1] represents the decomposition rate; here we choose q = 0.7071. The same parameter values as in [32] are taken, and appropriate value for h is chosen:

According to the relationship between the parameters of the original model (1.3) and the scaled model (2.1), the parameter values in the scaled model (2.1) are chosen as follows:

Choose a as the bifurcation parameter. With these parameter values, we obtain

To get R 0 < 1, a > 2 is required, which in turn yields that the infection-free equilibrium E 0 = ( 1 d , 0, 0, 0) = (2, 0, 0, 0) is globally stable by Theorem 2.4. For a = 3, system (2.1) has a globally stable infection-free equilibrium E 0 , which is shown in Fig. 1 . As shown in Fig. 1 , the healthy and infected cells decrease and converge to the equilibrium E 0 directly, whereas CTLe and CTLp firstly increase, then decrease, and finally converge to zero, which revels that the infected cells CTLp and CTLe die out directly after a brief fluctuation.

As a decreases to the critical value a c 1 . = 2, R 0 increases and passes the threshold 1, which reveals that E 0 becomes unstable. However, the infectious equilibrium E 1 = (a, d(R 0 -1), 0, 0) = (a, 0.5( 2 a -1), 0, 0) without CTL appears. As stated in Theorem 2.5,

2 ) > 0, which implies that 1 < R 0 < R 1 is equivalent to 1.707 ≈ 2+ √ 2 2 < a < 2 or 1 4 < a < 2- According to Theorem 2.2, system (2.1) admits a unique positive equilibrium

2 ) (Theorem 2.2 (2)(i),2(ii)), then system (2.1) admits two positive equilibria E 2+ = (x 2+ , y 2+ , w 2+ , z 2+ ) and

which is unstable once it exists.

Next, we will prove Theorem 3. it is locally stable, as shown in Fig. 3 .

As can be seen in Fig. 3 (a) and Fig. 3 2 )). Intersecting with the existence condition of equilibrium E 2+ and B 3 < 0, we obtain 4 < 0 for a ∈ (a * , 2+ √ 2 2 ) and 4 = 0 for a = a * , whereas 4 > 0 for a ∈ (0, 1 4 ) ∪ ( 1 4 , 2-

2 , a * ), which implies that the two roots of 3 = 0 are (by (3.13)) [38] , we obtain that the first Lyapunov coefficient is -0.07004781 as h = h * 1 ≈ 0.029245, whereas another first Lyapunov coefficient is -0.001386494 as h = h * 2 ≈ 0.304727, and thus the two Hopf bifurcations are supercritical, and the limit cycles are stable. The phrase diagrams of system are shown in Fig. 4 .

As can be seen in Figs. 4(a), 4(b), and 4(c), a stable limit cycle appears, and stable periodic solutions bifurcate from it.

We must point out that when a ∈ (0, 1 4 ) ∪ ( 1 4 , 2-√ 2 2 ), the infectious equilibrium E 1 without CTL is stable, whereas the infectious equilibrium E 2+ with CTL exists. To display this case, let a = 0.27 and h = 0.4, which indicates that both E 1 and E 2+ are stable. By direct calculation we obtain E 1 ≈ (0.27, 3.2037, 0, 0) and E 2+ ≈ (1.7071, 0.0858, 1.4371, 6.7009), and we draw the diagrams with different initial values.

As shown in Figs. 5(a) and 5(b), there are two stable equilibria E 1 and E 2+ , that is, bistability occurs, and the infected cells converge to one of them depending on the initial values.

Besides, when equilibrium E 2+ is unstable, then a stable limit cycle occurs, and E 1 is still stable. To display this phenomenon directly, letting a = 0.1 and h = 0.3, we draw the diagrams with different initial values.

As shown in Figs. 6(a) and 6(b), a stable equilibrium E 1 , an unstable equilibrium E 2+ , and a stable limit cycle, which is bifurcated from E 2+ , appear, and the infected cells converge to one of them depending on the initial values.

For the case of C 3 < 0, in Theorems 3.3 and 3.4, we give another set of parameters d = 3.6 = b, q = 1, c = 28.8, h = 0.8. We do not discuss it now.

This paper studies an HIV model proposed by Wodarz et al. [31, 32] to describe the interaction between healthy cells and infected cells as well as primary and secondary immune response. Compared with [33] , in this model, we assume that the production of primary immune response is not only connected with infected cells but also with healthy cells. We also assume that virus at its steady state is proportional to infected cells. The structure of equilibria is analyzed in [31, 32] . But for a higher-dimensional system, stability and bifurcation analysis is important and complex for the full range of possibilities. Because of adding the healthy cells to the produced CTLp, the model shows rich dynamic behavior on stability and bifurcations. It is interesting that this model displays the bistability phenomenon, that is, two stable equilibria E 1 and E 2+ or a stable equilibrium E 1 and a stable limit cycle, which is bifurcated from the unstable equilibrium E 2+ . Which one is stable not only depends on relationship of parameters but also depends on the initial values of cells. As shown in Figs. 5(a) and 5(b), high initial virus load close to E 1 leads to the convergence to E 1 , which means that CTL memory fails to establish. Low initial virus load close to E 2+ leads to the convergence to E 2+ , CTL memory successfully establishes, and the virus load first increases and then decreases to stay at a low level. Therefore we must increase dosage to inhibit the virus replication in a brief period and help the immune response establish. As shown in Fig. 6(a) , if the initial healthy cells and virus load are the same, and high initial CTL account means that the CTL response establishes, virus load may decrease at primary process, which was considered as the CTL clear away some virus, but over time, it oscillates and cannot be completely eradicated. This phenomenon can be viewed as an individual having a chronic disease that may flare up from time to time, and it is a long struggle between virus and immune response. Initial values with little CTL lead to high virus load. As can been in Fig. 6(c) , high or low initial healthy cells do not change the development of disease. However, high initial healthy cells first lead to decreasing the virus load to a very low state, which is almost clear away the virus, which means that the healthy cells play an important role in clearing away the virus under certain circumstances.

The interaction of virus and the host immune system is a complicated and long process as HIV has a long latent period, and the disease cannot cure completely. We have shown rich dynamic patterns, but the model considered here is just a simple one. It is easy to improve and expand the model. For example, we can add the virus equation in model (1.3) (see a more detailed description in [32] ), we may consider delay in this model as in [39] , drug treatment [13] , or latent cells [40] mentioned in the introduction. Such modifications should more precisely react the reality and give us more advice in understanding the infection process, which leads to a more challenging mathematical analysis.

Hunaids: unaids data

Modeling the impact of non-pharmaceutical interventions on the dynamics of novel coronavirus with optimal control analysis with a case study

Modeling the dynamics of novel coronavirus (2019-ncov) with fractional derivative

Modelling the spread of covid-19 with new fractal-fractional operators: can the lockdown save mankind before vaccination?

Analysis of fractal-fractional malaria transmission model

Stability analysis and numerical solutions of fractional order hiv/aids model

A fractional order HIV-TB coinfection model with nonsingular Mittag-Leffler law

Fractional model of HIV transmission with awareness effect

Population dynamics of immune responses to persistent viruses

Global properties of basic virus dynamics models

Virus dynamics and drug therapy

Virus dynamics: a global analysis

Mathematical analysis of an HIV model with impulsive antiretroviral drug doses

Global analysis of HIV-1 dynamics with hill type infection rate and intracellular delay

Bifurcations of a mathematical model for HIV dynamics

Dynamics of an HIV infection model with two infection routes and evolutionary competition between two viral strains

Bistable dynamics and Hopf bifurcation in a refined model of early stage HIV infection

Presence of an inducible HIV-1 latent reservoir during highly active antiretroviral therapy

Modeling latently infected cell activation: viral and latent reservoir persistence, and viral blips in HIV-infected patients on potent therapy

Global stability of a delayed virus model with latent infection and Beddington-DeAngelis infection function

HIV-1 infection dynamics and optimal control with Crowley-Martin function response

Modeling the dynamics of viral infections in presence of latently infected cells

Viral infection model with periodic lytic immune response

Global stability in a viral infection model with lytic and nonlytic immune responses

An HIV infection model based on a vectored immunoprophylaxis experiment

A chronic viral infection model with immune impairment

Modeling the cell-to-cell transmission dynamics of viral infection under the exposure of non-cytolytic cure

Hopf bifurcation analysis of nonlinear HIV infection model and the effect of delayed immune response with drug therapies

Stability of a general CTL-mediated immunity hiv infection model with silent infected cell-to-cell spread

Killer Cell Dynamics: Mathematical and Computational Approaches to Immunology

Dynamics of cytotoxic T-lymphocyte exhaustion

A new theory of cytotoxic T-lymphocyte memory: implications for HIV treatment

Bifurcation analysis in a model of cytotoxic T-lymphocyte response to viral infections

Differential equation models of some parasitic infections: methods for the study of asymptotic behavior

Dynamics of an HIV-1 therapy model of fighting a virus with another virus

Dynamics of an HIV-1 infection model with cell mediated immunity

Closed-form conditions of bifurcation points for general differential equations

Matcont: a Matlab package for numerical bifurcation analysis of ODEs

A mathematical model of cell-to-cell spread of HIV-1 that includes a time delay

Asymmetric division of activated latently infected cells may explain the decay kinetics of the HIV-1 latent reservoir and intermittent viral blips

Our deepest gratitude goes to the editor, referees, and professor Changrong Zhu for their suggestions, which have helped to substantially improve this paper.

The research is supported by the Science and Technology Department of Yangtze Normal University (NO.2016XJQN05).

Data sharing not applicable to this paper as no datasets were generated or analyzed during the current study.

Now we verify the Hopf bifurcation by Theorem 3.4. The formulas given by (3.14) are It follows that the existence conditions of positive equilibrium E 2+ (a ∈ (0, 12 )) also directly guarantee that A 3 > 0 and C 3 > 0 (Theorem 3.4 (2)). Since B 3 < 0 (Theorem 3.4 (2) 

The authors declare that they have no competing interests.Authors' contributions CL found the model and carried out the numerical simulation, and CL and LK did the mathematical analysis and wrote the paper together. All authors read and approved the final manuscript. 

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.Received: 19 June 2020 Accepted: 6 October 2020