Abstract
In this paper, we study a model of oncolytic virus infection with two time delays, one of which is the time from the entry of viruses into tumor cells to start gene replication, and the other is the time from the entry of viruses into tumor cells to release new virus particles by infected tumor cells. In previous studies on oncolytic virus infection models, the infection rate was linear. Combined with the virus infection models, the saturated infection rate, is further considered to describe the dynamic evolution between viruses and tumor cells more objectively so as to further study the therapeutic effect of oncolytic viruses. This paper discusses the dynamics of the system under three conditions: (1) , (2) and , and (3) and , and proves the global stability and local stability of the virusfree equilibrium, the stability of the infection equilibrium, and the existence of Hopf bifurcation. Finally, the conclusions of the paper are verified by MATLAB numerical simulations.
1. Introduction
An oncolytic virus is a natural or genetically modified virus that can specifically infect and kill tumor cells without damaging healthy cells (see [1]). Oncolytic virotherapy is a very promising cancer treatment method. Initially, The Lancet reported that influenza viruses can cause tumors in patients to regress, resulting in the emerging concept of oncolytic viruses. Since then, more than 160 different oncolytic viruses have been in preclinical studies and clinical trials (see [2]). According to current development trends, oncolytic viruses can be positioned as the next major breakthrough in cancer treatment after checkpoint inhibitors in immunotherapy.
Many researchers in the fields of medicine, biology, and mathematics have been exploring and improving the related problems of oncolytic virotherapy, which has brought breakthroughs in the treatment of some malignant tumors, but there are still many problems to be solved at present (see [3]). For example, the proliferation of oncolytic viruses in host solid tumor cells is limited and the difference in tumor cell receptor expression limits the spread of the virus, which in turn affects the therapeutic effect. Immune-mediated pharmacokinetic therapy is relatively slow compared to drugs that directly kill tumor cells, so time delays need to be considered in clinical trial design and outcome evaluation. These problems have affected the antitumor effect of oncolytic viruses and failed to achieve the goal of systematic and reliable tumor elimination of oncolytic viruses. In the research process of oncolytic virotherapy, mathematical models can better understand the dynamic mechanism of tumor cells-virus-immune cells (Cancer-Virus-CTLs). It also provides an effective tool to predict the outcome of immunotherapy and optimize tumor treatment strategies.
In order to describe the dynamic changes of the immune response in viral therapy, a series of ordinary differential equation models have been constructed in the past decades (see [4–11]. Dingli et al. published a model based on population dynamics in the mathematical biosciences (see [12]), which focused on the basic elements of radiation viroid therapy and discussed and proved the stability of the equilibrium of complete cure, local treatment, and treatment failure. On the basis of D. Dingli et al., Kim et al. proposed a virus model that considered both virus-specific immunity and tumor-specific immunity (see [13]). In addition, since the life cycle of viruses contains many intracellular processes, the time from the entry of viruses into tumor cells to the start of gene replication and the time of the release of new virus particles by the infected tumor cells are considered to further establish mathematical models with time delays. Wang et al. (see [14]) considered the cyclic period of viruses in cells and established a time-delay differential equation. Kim et al. (see [15]) proposed a delayed ordinary differential equation model. Li et al. (see [16]) further considered the time delay of virus self-replication based on the research of Kim et al. proposed an ordinary differential equation model with two time delays to analyze the influence of the time delay on the stability of the equilibrium.
For most models of oncolytic virus infection, the infection rate of oncolytic viruses and uninfected tumor cells is bilinear, but the actual incidence is not strictly linear. Therefore, we further assume that the model has a saturated infection rate . In this paper, we study a virus infection model with a saturated infection rate .
The rest of this paper is organized as follows: the second section builds the model. The third section solves the equilibrium of the system. The fourth section discusses the stability of the virusfree equilibrium and the virus infection equilibrium, as well as the existence conditions of Hopf bifurcation of the virus infection equilibrium, and explains the dynamic properties under three conditions: (1) , (2) and , and (3) and . In the fifth section, numerical simulations are performed to verify the obtained results. Finally, the conclusions of this article are given.
2. Mathematical Model Formulation
On the basis of the research (see [16]), this chapter considers the impact of saturated infection rate and two time delays on the virus, and the following model is shown:where , and represent the number of uninfected tumor cells, infected tumor cells, oncolytic viruses, and CTL immune cells, respectively. represents the maximum growth rate of uninfected tumor cells. represents the average rate at which the infected tumor cells release viruses. and represent the rate of CTL immune cells regeneration. represents the maximum carrying capacity of the human environment to tumor cells. represents the delay from viruses entering tumor cells to the beginning of gene replication. indicates the time from the entry of viruses into tumor cells to the release of new virus particles by infected tumor cells. Obviously, there is . and represent the rate of CTL immune cells clearing uninfected tumor cells and infected tumor cells, respectively. , and indicate the natural mortality of infected tumor cells, oncolytic viruses, and immune cells, respectively. The rupture of tumor cells releases tumor antigens, thereby inducing a systemic antitumor immune response, which may enhance the cytolytic activity of viruses.
The initial condition of system (1) iswhere . . Note , define the norm , where is a continuous function, constitutes a Banach space.
3. The Existence of the Equilibrium
Through simple calculations, can be obtained from the first equation of system (1). The virusfree equilibrium of system (1) is , where and .
Let the infection equilibrium of system (1) be , which is the positive solution of the following equations. Based on the model, there are reasonable assumptions as follows: , and . All the above parameters are positive constants.
Through simple calculations, the following conclusions are obtained:where
When , we get . Substituting (4)–(6) into the first equation of (3), we obtain the quadratic equation of one variable with respect to .where
When , equation (8) has a unique positive root . Therefore, is called as the basic regeneration number of system (1), when and , system (1) has a unique positive equilibrium . is the number of basic regeneration, and is the number of immune response regeneration. and together ensure the existence of positive equilibrium .
4. Stability Analysis of the Equilibrium
In order to study the stability of the equilibrium, we consider the linear approximation equation of system (1) at the equilibrium .where ,
Letwhere is a fourth-order identity matrix.
4.1. Global and Local Stability of Virusfree Equilibrium
Theorem 1. If , the virusfree equilibrium is globally asymptotically stable for any and .
Proof. Assuming , we define the Lyapunov functionalFurther, calculating the derivative along with system (1)When , is obtained, therefore . It is obvious to obtain that if and only if . Assuming that M̃ is the largest invariant set contained in , it is easy to get M̃ = . According to the LaSalle’s invariant principle, the virusfree equilibrium is globally asymptotically stable.
Theorem 2. If , the virusfree equilibrium is locally asymptotically stable for any and . If , is unstable for any and .
Proof. According to equation (12), the characteristic equation of virusfree equilibrium is as follows:Through certain calculations, it follows thatConsidering the following two equations separately(i)Equation (17) can be written as . Because and , the roots of equation (15) have negative real parts for all and .(ii)To discuss the case of the roots of equation (18), first of all, we consider , so it is easy to get equation (18) as follows:When , it follows that , therefore, all the roots of (19) have negative real parts. When , is a single root of (18), and the remaining roots have negative real parts. When , , therefore, (19) has a root of positive real parts.
Secondly, we consider and . Assuming that for some and , is the roots of (18), substituting into (18) and separating the real and imaginary parts to get thatBased on , it follows thatwhereIf , because of the following inequalityEquation (21) has no positive root. Thus, when , (18) has no pure imaginary root. Consequently, it can be seen from (i) and (ii) that when , is locally asymptotically stable for any and . When , is continuous with respect to , and , respectively. For any , we get and . Therefore, there exists such that holds. In summary, when , the equilibrium is unstable for any and .
4.2. Local Stability Analysis of the Virus Infection Equilibrium
Suppose that the positive equilibrium of system (1) is , and the characteristic equation corresponding to the positive equilibrium of equation (12) is as follows:
which is equivalent to the following equation:where
The case and .
When and , the characteristic equation (25) becomeswhere
(H1) .
The following theorem is obtained by the Routh–Hurwitz criterion.
Theorem 3. Suppose and , if and condition (H1) holds, then the equilibrium is locally asymptotically stable.
4.2.1. The Case and
When and , the characteristic equation (25) is written as follows:where
Assuming that is a pure imaginary root of (29), substituting into (29) and separating the real and the imaginary parts to obtain that
According to , we can get thatwhere
Let , it follows from (32) that
Let be the positive roots of (32), then (29) has pure imaginary roots , where . From (31), the following solutions are obtained:where , . The distribution of the roots of (34) is discussed next, and the following conditions are given: (S1) (1) and . (2) and . It is sufficient if any one of the above conditions is satisfied. (S2) (1) and . (2) and . Meeting any one of the above conditions can be established. (S3) and . (S4) and . (S5) (1) and . (2) and . (3) and . (4) and .
Any one of the above conditions can be satisfied, where
Next, we discuss the relationship between the roots of the characteristic equation and the coefficient. According to conclusions (see [17]), the following lemma is obtained.
Lemma 1. For equation (34), the following conclusions are drawn:(i)If (S1) holds, then equation (34) has only one root of the positive real parts. Satisfying condition (1), the equation has one positive root and three negative roots. Satisfying condition (2), there is a positive root, a negative root, and a pair of conjugate virtual roots.(ii)If (S2) holds, then equation (34) has two roots with positive real parts. Satisfying condition (1), the equation has two positive roots and two negative roots. Satisfying condition (2), there are two positive roots and a pair of conjugate virtual roots.(iii)If (S3) holds, then equation (34) has three different positive real roots and one negative root.(iv)If (S4) holds, then equation (34) has four different positive real roots.(v)If (S5) holds, then equation (34) has no positive real root. That is, if (1) or (2) is satisfied, the equation has two pairs of different conjugate virtual roots. Satisfying condition (3), equation (34) has four different negative roots. Satisfying condition (4), there are two negative roots and a pair of conjugate virtual roots.
According to Lemma 1, the following lemma is obtained:
Lemma 2. For equation (29), the following conclusions are drawn:(i)If (S1) holds and , then equation (29) has a pair of pure imaginary roots , where and .(ii)If (S2) holds and , then equation (29) has two pairs of pure imaginary roots , where and .(iii)If (S3) holds and , then equation (29) has three pairs of pure imaginary roots , where and .(iv)If (S4) holds and , then equation (29) has four pairs of pure imaginary roots , where and .(v)If (S5) holds and , then equation (29) has no pure imaginary root.
Lemma 3. If is a solution to equation (32), then
Proof. Assuming that is the solution of equation (25), and both sides of the equation are differentiated from at the same time, it follows thatwhere , and denote the derivation of , respectively. According to , we obtain thatThrough further calculations, we get thatWhen , , it follows that , and then substitute into to get the following equation:Therefore,The proof is completed.
According to Lemmas 1–3, the following theorem is obtained.
Theorem 4. For equation (29), the following conclusions are established:(i)If (S1) is true, then .(ii)If (S2) is true, then and .(iii)If (S3) is true, then , (iv)If (S4) is true, then , .
Proof. (i)If (S1) is true, then equation (34) has one positive root and three negative roots or one positive root, one negative root, and a pair of conjugate imaginary roots, and is monotonically increasing at the positive root , so we get . Therefore, the conclusion is true.(ii)If (S2) is true, then equation (34) has two positive roots ( and ) and two negative roots or has two positive roots and a pair of conjugate imaginary roots. Function decreases monotonically at the positive root and increases monotonically at the positive root , so and are obtained. Therefore, the conclusion is true.(iii)If (S3) is true, equation (34) has three different positive real roots ( and ) and one negative root. Function increases monotonically at , decreases monotonically at , and increases monotonically at . So, we get , and . Therefore, the conclusion is true.(iv)If (S4) holds, equation (34) has four different positive real roots ( and ). Function is monotonically decreasing at and , and monotonically increasing at and . So we obtain , , and . Therefore, the conclusion is true.
According to the analysis in references (see [18–21]), Lemma 4 and Theorem 5 are obtained on the basis of Theorems 3 and 4.
Lemma 4. For equation (25), when and , the following conclusions are obtained:(i)If (H1) and (S1) are true, then for , all the roots have a negative real part. There is a pair of pure imaginary roots at . There is at least one root with a positive real part at .(ii)If (H1) and (S2) are true, then there exists a positive integer that makes . This indicates that the system goes through stability switches from stable state to an unstable state. All the roots have negative real parts at . There are two pairs of pure imaginary roots at or , and the other roots have negative real parts. There is at least one root with a positive real part at or .(iii)If (S2) is true and (H1) is not, then there exists a positive integer that makes . This means that the system goes through stability switches from an unstable state to a stable state. There is at least one root with a positive real part at or . There are two pairs of pure imaginary roots at or , and the other roots have negative real parts. All the roots have negative real parts at .(iv)If (H1) and (S3) are true, or (H1) and (S4) are true, then the system has at least one stability switch.(v)If (H1) and (S5) are true, then for any , all roots have negative real parts.
Theorem 5. Suppose and .(i)If (H1) and (S1) are true, then the equilibrium is locally asymptotically stable for any . When , system (1) appears Hopf bifurcation at the equilibrium .(ii)If (H1) and (S2) are true, when , the equilibrium is locally asymptotically stable. When or , the equilibrium is unstable. When or , system (1) appears Hopf bifurcation at .(iii)If (S2) is true and (H1) is not. When , the equilibrium is locally asymptotically stable. When or , is unstable. When or , Hopf bifurcation occurs at .(iv)If (H1) and (S3), or (H1) and (S4) are true, then there is a stability switch at . When satisfies , Hopf bifurcation occurs at .(v)If (H1) and (S5) are true, then is locally asymptotically stable for any .
4.2.2. The Case and
In order to discuss the possible Hopf bifurcation at and , we regard as a parameter of . According to the analyses and conclusions of (see [19,22]), system (1) is stable at and , there is , and the system reaches stability at . Therefore, the following theorem can be obtained:
Theorem 6. Let us assume that and .(i)If (H1) and (S1) are true, then for any , there exists , and the equilibrium is locally asymptotically stable at .(ii)If (H1)and (S2) are true, there exists positive integer , when , there exists , and the equilibrium is locally asymptotically stable at .(iii)If (H1) and (S5) hold, then for any , there exists , and the equilibrium is locally asymptotically stable at .
According to the analysis (see [23]), when is fixed, the characteristic equation of system (1) only depends on the change of time delay . Therefore, we fix as a constant and study the stability of the system when increases from zero. In this case, the characteristic equation (25) becomeswhere
Assuming that the characteristic equation has pure imaginary root , we substitute it into the characteristic (43) and then separate the real parts and imaginary parts to obtain that
According to , it follows that
Through further calculations, we get the following equation:where
In order to study the influence of two time delays on the stability of the system, the following five criteria are given according to reference (see [23]): (U1) For any , there is , that is, is not the root of the characteristic equation (25). (U2) If , then for any , inequality can be obtained. (U3) For any , holds. (U4) For any , the number of roots of the following function is limited. (U5) As long as has positive roots , then it is continuously differentiable with respect to .
Next, we prove that for any , where is the maximum value that enables the equilibrium to exist, each item in the characteristic (43) conforms to the above five criteria. Suppose that (U2) is true, and then we prove that the other four criteria are true.
By calculations, we obtain , . For the convenience of expression, we give conditions.
(G1) .
If (G1) holds, then (U1) is correct and has at least one positive root. Next, when , obviously, (U3) holds. Then, on the basis of (U2), we can conclude that (U4) is true. Finally, we verify (U5). Suppose that is a positive real root of and is a point in the domain of . Obviously, and exist and they are continuous in some regions of . In addition, there is . Therefore, according to the implicit function theorem, the verification of (U5) is true.
If is the set of roots of , for any , we obtain . Let us define , we get and .
Let us define the mapping
There exists a positive root of in . Obviously, is continuously differentiable with respect to . Therefore, according to [23], the following theorem can be obtained:
Theorem 7. Suppose that one of the following conditions is correct.(i)For any fixed , (H1) and (S1) hold.(ii)For any fixed , (H1) and (S2) hold.Then we obtain the following conclusions: if and , (H1) and (G1) are true, then has a pair of pure imaginary roots , where . For certain , when , if , then the root of the characteristic equation crosses the imaginary axis from left to right (from right to left) at , where Table 1
5. Numerical Simulations
In this section, the main conclusions drawn in the paper are verified through a series of numerical simulations, and the influence of two time delays on the stability of the system (1) is discussed. Some of the parameters of numerical simulations are selected with reference to [8, 16].
When , the following parameters are selected:Through calculation, the virusfree equilibrium and are obtained. Applying Theorem 1, we obtain that the virusfree equilibrium is globally asymptotically stable. Finally, Figure 1 verifies this conclusion.
(a)
(b)
When , the parameters we selected are as follows:
By simple calculations, we obtain that the positive equilibrium of system (1) is , , and . According to condition (i) of theorem 3, the equilibrium is locally asymptotically stable, and the numerical simulation is shown in Figure 2.
(a)
(b)
When and , (29) has no positive real part root. Moreover, , , and . So, (H1) and (S5) hold. According to (v) in theorem 5, the equilibrium is locally asymptotically stable. As is shown in Figure 3, the number of uninfected tumor cells eventually tends to be stable.
(a)
(b)
We set and , we still take the set of parameters above. By simple calculations, it can obtain . We can compute , , and . So, (H1) and (S5) hold, which satisfies the condition (iii) of theorem 6. Therefore, according to theorem 6, the equilibrium is locally asymptotically stable. As shown in Figure 4, the number of uninfected tumor cells gradually tends to stabilize Figure 5.
(a)
(b)
(a)
(b)
Then we discuss the partial Hopf bifurcation with a time delay. When and is used as the bifurcation parameter, the selected parameters are , , , , , , , , , , , and . When (), we obtain . According to condition (iii) of Theorem 6, the numerical simulation is shown in Figure 6. The number of uninfected tumor cells fluctuates periodically.
(a)
(b)
Next, we discuss the Hopf bifurcation with two time delays. is fixed and is taken as the bifurcation parameter. The parameters selected are , , , , , , , , , , , and . When is fixed and (), by calculations, it obtains , , , , and . As shown in Figure 7, the number of uninfected tumor cells tends to be ultimately stable. When is fixed and (), as shown in Figures 8(a) and 8(b), the number of uninfected tumor cells always fluctuates periodically. Figures 7 and8 demonstrate the condition (i) of Theorem 7, that there exists a Hopf bifurcation with two time delays.
(a)
(b)
(a)
(b)
The parameters selected are , , , , , , , , , , , and . It can obtain and . When , the positive equilibrium is unstable. With the increase of , the system switches from unstable to stable. The numerical simulations are shown in Figures 9 and 10, where Figure 10 is the bifurcation diagram showing the stability of positive equilibrium changing with . These images demonstrate Theorem 5.
(a)
(b)
(c)
The parameters selected are , , , , , , , , , , , , and . and are obtained. When , the positive equilibrium is stable. According to Theorem 5, with the increase of , the system appears to switch phenomenons from stable to unstable. The numerical simulations are shown in Figures 11 and 12, where Figure 12 is a bifurcation diagram showing the change of the stability of with . The following images prove Theorem 5.
(a)
(b)
(c)
(d)
6. Conclusions
In this paper, we study a model with two time delays of oncolytic virus infection with saturation incidence. The first delay is the time from oncolytic viruses entering tumor cells to gene replication. The second delay is the time from viruses entering tumor cells to infected tumor cells releasing new virus particles. In the previous studies of oncolytic virus infection models, the infection rate was linear. In this paper, we further consider the nonlinear infection rate. The global asymptotic stability and local asymptotic stability of the uninfected equilibrium, as well as the conditions for the local asymptotic stability of the virus-infected equilibrium, are discussed, and the effect of time delay on the stability of the equilibrium is also discussed. Through theoretical analysis and numerical simulation, it is obtained that the local stability of the uninfected equilibrium is independent of the size of the time delay, but the change in time delay can promote or destroy the stability of the virus infection equilibrium. That is, the time delay will destroy the stable state of tumor treatment, and even produce Hopf bifurcation and lead to periodic oscillation.
When only one time delay is considered, that is and , taking as the bifurcation parameter, the conditions for generating Hopf bifurcation are obtained. It is verified that the time delay from the virus entering the tumor cells to the tumor cells cracking to produce new virus particles makes the stable equilibrium of the system become an unstable equilibrium, that is, there is a critical value of . When , the virus-infected equilibrium is locally asymptotically stable; when , the virus-infected equilibrium is unstable. At , system (1) generates Hopf bifurcation, which indicates that the time delay of oncolytic virus infection will affect the stability of tumor treatment. Therefore, reasonable control of time delay plays a key role in effectively curbing the growth and recurrence of tumor cells.
When considering the two time delays of the system, that is and , taking as the parameter of and as the bifurcation parameter, the conditions for generating Hopf bifurcation are discussed. When both delays exist, there are critical values and . When , , the infection equilibrium is locally asymptotically stable, but when , the infection equilibrium becomes unstable, and Hopf bifurcation occurs at the critical value . When is fixed, the delay time of is small in a certain range, and the system is in a stable state. As the delay time increases, the system changes to an unstable state. In other words, when the time for tumor cell lysis to produce new virions is fixed, the faster the replication rate of oncolytic viruses is, the higher the concentration in the cell is, the more virions are released and the immune system is further activated. In this case, the system is in a stable state. When the amount of time delay increases, the replication rate of viruses is relatively slow, the number of tumor cells lysis and death due to infection with oncolytic virus is reduced, and the tumor immune system enters the escape stage in the tumor immune editing process. The tumor cells in the system are not effectively contained, and the system is in an unstable state. At the same time, the replication time of the oncolytic virus cannot be blindly accelerated, and the excessive intracellular concentration of the oncolytic virus cannot be pursued, because oncolytic viruses not only kill tumor cells, but also reduce the immune system response caused by tumor cell stimulation, which may also lead to the failure of tumor cell treatment. That is to say, too little concentration of oncolytic virus makes it difficult to effectively remove tumor cells, and a high concentration of oncolytic virus cannot fully enhance the immune system and may even have pathological effects. In the process of immune system interaction, different time delays can lead to different tumor states, so the selection of two time delays plays an important role in tumor cell therapy.
In modern medical research, genetic engineering techniques such as DNA recombination technology can change or obtain more efficient oncolytic viruses and delay the time of infection. By controlling infection delay, the efficiency of viral infection can be further improved, which provides an effective method for the treatment of oncolytic viruses and has important application value.
This paper explores the effect of discrete-time delays on the tumor treatment interaction system. In the future, we can extend the deterministic system to a stochastic system to more comprehensively and realistically study the stochastic dynamics of oncolytic virus infection models, which will enrich our understanding of the influence of biological delay in the process of virus infection and its interaction with immune response.
Data Availability
Previously reported data were used to support this study. These prior datasets are cited at relevant places within the text as references.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This work was supported by the National Natural Science Foundation of China (no. 11971055).