Stability analysis and optimal control on a multi-strain coinfection model with amplification and vaccination


                  In this paper, a multi-strain coinfection model with amplification (or mutation) is established to characterize the interaction between common strain and amplified strain, as well as vaccination. The basic reproduction number 
                        
                           
                              ℛ
                           
                           
                              0
                           
                        
                      is derived, from which the criteria on the existence and local (or global) stability of equilibria (including disease-free, dominant-strain and coexistence-strain) are established. By analyzing the effectiveness of vaccination, we find that a critical inoculation level could make the disease eliminate when 
                        
                           
                              
                                 ℛ
                              
                              
                                 0
                              
                           
                           <
                           1
                        
                     , while inefficient vaccines could cause backward bifurcation when 
                        
                           
                              
                                 ℛ
                              
                              
                                 0
                              
                           
                           <
                           1
                        
                     . Based on sensitivity analysis and realistic control policy, the optimal strategy of disease control is obtained. The theoretical results are illustrated by numerical simulation and clinical data of COVID-19 in Morocco.
               

Furthermore, mutations of pathogens are very common during the replication of DNA or RNA in multi-strain infectious diseases.Merker et al. [27] found that Beijing-lineage has been the most successful branch among all strains of tuberculosis mycobacterium over the past 200 years.Burgess [28] provided valuable strategies and analysis methods on how to study the functional consequences of disease-related mutation.Bonhoeffer et al. [17] pioneered the following SI epidemic model with multi-strain transmission and mutation, where () is the number of susceptible people at  time,   () is the number of infected people with the ith-strain at  time,  is the mutation rate,  is the natural mortality,   is the death rate,  is the birth rate of susceptible hosts,   is the infection rate of the ith-strain.They found that mutations and competition between different pathogens within host can influence the evolution of pathogen virulence.Castillo-Chavez and Feng [8] developed a TB model with drug-resistant and drug sensitive strains due to non-compliance of treatment.They found that there exists competition exclusion and coexistence between two strains in which antibiotic resistance will lead one strain to elimination.Feng et al. [9] proposed a two strains tuberculosis model with infection age to explore the influence of different incubation periods on disease transmission due to antibiotic resistance and incomplete treatment.
Recently, the phenomenon of amplification among multi-strain transmission diseases, i.e., mutation with acquired drug resistance, has https://doi.org/10.1016/j.rinp.2023.106556Received 28 March 2023; Received in revised form 14 May 2023; Accepted 15 May 2023 been a center issue [10,[29][30][31][32][33].Blower and Chou [10] proposed a multi-strain epidemic model to simulate the dynamic behavior of tuberculosis with amplification.They found that those areas with reducing wild susceptible strains tend to be the hotspots of MDR-TB.Kendall et al. [29] constructed a TB epidemic model with treatment-related amplification, they found that acquired resistance was the initial factor of transmission of drug-resistant strains due to low concentration levels of therapeutic drug or incorrect treatment [30][31][32][33].Kuddus et al. [30] proposed a two-strain epidemic model with amplification as follows where   is the effective contact rate between infected population   with drug-sensitive strain and susceptible population ,   is the cure rate of infected population with drug-sensitive strain,   is the death rate due to disease in infected population with drug-sensitive strain;   is the effective contact rate between infected population   with drug-resistant strain and susceptible population ,   is the cure rate for infected population with drug-resistant strain,   is the disease-related death rate of infected population with drug-resistant strain;  is the proportion of individuals amplified from drug-sensitive strain to drug-resistant strain, which occurs during treatment; () is the number of recovered people.The basic reproduction number was obtained, by which the existence and stability of three kinds of equilibria are discussed.They found that the contact rate between two strains has significant impact on the prevalence of disease, and lowquality treatment make it more possibly coexist for two strains.During the COVID-19 pandemic, some mutations in drug resistance have been reported [34].
Moreover, vaccination is one of the greatest achievements in public health prevention, which has completely eliminated the smallpox worldwide and the near extinction of polio [1].He et al. estimated that more than 1.5 million deaths are averted in 12 countries due to vaccination for COVID-19 [35].Vaccines could not guarantee entire protection against disease, there still exists potentially infection for those people who have been vaccinated, especially for multi-strain diseases, where some pathogens mutate even if the host develops antibodies [1], such as dengue and COVID-19.Pathogen adaptation to vaccination may lead to multiple forms of mutation in the pathogen acting on different phenotypic traits [36].Most pathogens can undergo rapid antigenic evolution that allows them to escape vaccine-induced immune responses (vaccine escape) [37][38][39].When the vaccine is not full effective, it cannot completely protect people from infection, which will make the pathogen experience double-selection of attacking target, namely the vaccinated host and the unvaccinated host, so as to mutate [40,41].McLeod and Gandon [36] considered two kinds of evolution come from by mutation: the escape of pathogens from vaccine-induced immune responses, and the adjustment of pathogen virulence.They elucidated the roles played by epistasis and recombination, and highlighted the different protective effects of vaccination.Obviously, the existence of multi-strain with amplification make vaccination more inefficient to confront disease transmission control.For example, COVID-19 has now emerged with some highly contagious mutant strains, such as Delta and Omicron variants [42], posing new challenges to control the outbreak of disease even though many kinds of vaccine have been used [11,12].
Based on the above considerations, according to the genetic variation of multi-strain amplification, vaccination impact, we propose a multi-strain SVIR epidemic model, where the whole infected population are divided into two compartments,   () denotes the number of infected population with common strain at time , for which the vaccine and drugs are effective;   () presents the number of infected people with amplified strain at time .() is the number of susceptible people at time ,  () is the number of vaccinated people at time .
For convenience, we assume that vaccine is completely effective to the common strain while ineffective to amplified strain.The effective contact rate between infected population with common strain and susceptible population is , where  is the contact rate between infected population and susceptible population,  ∈ [0, 1] is a parameter that measures population psychology or inhibition.(1 − )    indicates that the common strain develops amplification (or mutate) during transmission due to improper treatments, which further leads to amplification or even creation of new strains. is the recruitment rate,  is the natural mortality rate,  is the vaccination rate,   is the cure rate of infected population with common strain,   is the death rate due to disease infected population with common strain;   is the cure rate of infected population with amplified or mutant strain,   is the disease-related death rate of the infected population with amplified or mutant strain.
This paper is organized as follows, in the second section, some basic properties of a multi-strain epidemic model with amplification (or mutation) are obtained.In Section ''The existence and stability of equilibra'', the basic reproduction number ℛ 0 is established, the existence and local (or global) stability of disease-free equilibrium, dominant equilibrium and coexistence equilibrium are discussed respectively.Section ''Effects of vaccine effectiveness on multi-strain infectious disease'' analyzes effects of vaccine effectiveness for multi-strain infectious diseases, which have a critical level against strains.At the same time, unrealistic vaccines can cause backward bifurcation.In Section ''Sensitivity analysis and optimal control'', sensitivity analysis of ℛ 0 is carried out to determine the model parameters with the greatest influence on the prevalence of multiple strains infectious disease, and optimal control strategy is proposed.In Section ''Numerical simulation'', the theoretical results are illustrated by numerical simulation.Finally, a brief discussion and conclusion are given in Section ''Discussion and conclusion''.

Positivity and boundedness of solutions
For convenience, let   =  +   +   ,   =  +   +   .For the model (3), we have the following initial condition By the theorem of the existence of solution of ordinary differential equation [1], it is easy to obtain that model (3) Considering the first case, for every  ∈ [0, t), we have We further have we can get ( t) > 0, this is a contradiction, so the case (1) is not true.

The existence and stability of equilibra
The disease-free equilibrium Let   () =   () = 0, we can get the following disease-free equilibrium of model ( 5) Obviously, the disease-free equilibrium always exists.We calculate the basic reproduction number ℛ 0 for model (5) by the next generation matrix method [43].
Consider   () and   () compartments, It is easy to get ) .

The dominant equilibrium of amplified strain
Consider the dominant equilibrium of common strain, i.e.,   ≠ 0,   = 0.According to the last equation of model ( 5),   (1 − )  = 0, therefore   = 0.This is a contradiction, therefore the dominant equilibrium of common strain of model ( 5) does not exist.
} , then the dominant equilibrium  * of model ( 5) is locally asymptotically stable.
Proof.The Jacobian matrix of model (5) at We have the following characteristic equation where it can be known that  3 and  4 have negative real parts from the Routh-Hurwitz criterion [44].In conclusion, if } , the dominant equilibrium  * is locally asymptotically stable.□

Next, we introduce the invasion reproduction number [1]
where  * denotes the number of secondary infections, i.e, the number of new born infection via an infected individual with common strain; denotes the number of secondary infections that one infected individual with common strain could create during its lifetime.Mathematically, the invasion number ℛ 2 1 gives a threshold for the stability of the dominant equilibrium of amplified strain.And biologically, the invasion number gives the number of secondary infections one individual infected with common strain can produce in which amplified strain is at equilibrium during its lifetime.
As a result of Theorem 4, we can obtain the next corollary: 5) is globally asymptotically stable.
Proof.Consider the third equation of the model ( 5), Similarly, by the fourth equation of the model ( 5), we obtain further, H. Wu et al.
By Eqs. ( 6) and ( 7), it can be seen And because   () and   () are bounded, then we have , the hyperplane   = 0 attracts all the solutions of model (5).

The coexistence equilibrium
Suppose   ≠ 0 and   ≠ 0, from model ( 5) it is easy to get that That is, if the coexistence equilibrium exists, ℛ 0 > 1; when  †  ,  †  > 0, then  † > 0. The existence of coexistence equilibrium is equivalent to the existence of the positive root of   .In addition, if the positive root of   is satisfied )   ((+)+) , then  †  > 0.
In the case of  3 < 0, we can get the similar result.Therefore, we can have the following proposition.1.

Proposition 1. As shown in Table
(1) If the condition in case 1 is satisfied, there is no positive root of Eq. ( 8); (2) If one of the conditions in case 2 is satisfied, there is only one positive root of Eq. ( 8); (3) If one of the conditions in case 3 is satisfied, there are no or two positive roots of Eq. ( 8); (4) If the condition in case 4 is satisfied, there are three positive roots or one positive root of Eq. (8).
Remark 2. For parameters  3 ,  2 and  1 , the similar results can be obtained when at least one of them is zero.Because of the harsh conditions of this degenerate case, we assume  3 ≠ 0,  2 ≠ 0 and  1 ≠ 0. Remark 3. The Proposition 1 only can give sufficient condition of the positive root of   , however, from Remark 1, we see that the existence of the coexistence equilibrium  † = ( † ,  † ,  †  ,  †  ) of model ( 5) requires additional condition: that is The Jacobian matrix of model (5) at We have where From Routh-Hurwitz criterion [44], it can be known that if then the roots of Eq. (10)

Effects of vaccine effectiveness on multi-strain infectious disease
It is worth noting that the vaccines discussed in this paper are assumed to be effective on common strain but ineffective on amplified strain, which is crucial on control the out break of multi-strain disease.Let 1 −  1 be the effectiveness of the vaccine against common strain, 1 −  2 be the effectiveness of the vaccine against amplified strain. 1 is the temporary recovery rate of infected population with common strain,  2 is the temporary recovery rate of infected population with amplified strain,  1 is the proportion of infected population with common strain that return to the susceptible population,  2 is the proportion of infected population with amplified strain that return to the susceptible population.Other parameters and variables are defined in accordance with model (3), we propose the flow chart as follow Fig. 1.
We find that there is a critical vaccination threshold, which could control the outbreak of disease if it is effective against the strains, while backward bifurcation will occur if it is ineffective.

The effectiveness of vaccine against mutated or amplified strain
Now we assume that the vaccine is partially effective against infected individuals with mutated or amplified strain, the validity is expressed as 1 − .We have the following model as the special case of model (5).  ) .
After calculation, we have , further, it is obtained that Consider the limit state of inoculating , If ℛ 0 < 1, that means the vaccine is sufficiently effective against amplified strain, there is a critical level of inoculation  * , such that By calculation, we have And lim →∞ ℛ * 0 = 0, which is consistent the hypothesis that the vaccination is completely effective against common strain, and there is a critical level of inoculation  # , such that ℛ * 0 ( # ) = 1, namely, By calculation, we have there is a critical threshold value of vaccination  for the multi-strain disease when the vaccination is sufficiently effective against amplified strain.This has important implications for multi-strain diseases, e.g., COVID-19, for which vaccines should be as effective as possible against all strains.

Backward bifurcation caused by inefficient vaccine
In this paper, we assume that the vaccination is completely effective against common strain and inefficient against amplified strain, but in reality, some common strains are often resistant to the vaccination due to mutation.Vaccination effectiveness is expressed as 1 −  1 .For simplicity of calculation, let  = 0 in model (11), and then the epidemic model with vaccination and common strain is obtained, Similar to the discussion in Section ''The existence and stability of equilibria'', the disease-free equilibrium and basic reproduction number of model (12) are , 0, 0 Let   ≠ 0, we can get that .

Sensitivity analysis and optimal control
Sensitivity analysis is conducted to discuss the influence of parameters on dynamic behavior of model (5).We focus on local sensitivity analysis of basic reproduction number ℛ 0 by changing one input parameter, while fixing all other parameters to their baseline values.
The sensitivity index  ℛ 0  of ℛ 0 to parameter  [47] is as follows which are shown in Fig. 2. The positive sensitivity index indicates that ℛ 0 or ℛ 0 increases with some parameters; in contrast, the negative sensitivity indicator indicates that ℛ 0 or ℛ 0 decreases with parameters.Obviously, from Fig. 2, we can see that ℛ 0 and ℛ 0 are positively correlated with , and ℛ 0 is negatively correlated with .
To sum up, we get that (, ) is convex for any (, ).And  is compact, the solution of model ( 3) is bounded.By Filippov-Cesari Existence Theorem [1], there exists an optimal solution.□ Take the Lagrange function as and Hamiltonian function be where   represents the right end of -th state variable in model ( 14),   is the adjoint variable.Next, we use the Pontryagin's minimum principle [1] to seek the optimal control of model (14).

Numerical simulation
In this section, we illustrate the theoretical results through numerical simulation.We summarize parts of the numerical simulations of the model (3) dynamics as shown in Table 2.
It can be seen from Fig. 5, the disease-free equilibrium  0 and dominance equilibrium  * are unstable, and the solutions starting near the disease-free equilibrium  0 and the dominant equilibrium  * both tend to the coexistence equilibrium  † .
(d) Optimal control Let  = 10.0,  = 0.01,  = 0.06,  = 0.02,   = 0.008989,   = 0.245,  = 0.1,   = 0.1,   = 0.1,  = 0.0001.The weights of three measures are  1 = 1,  2 = 9,  3 = 0.01, respectively.Fig. 6 shows that the optimal control strategy of  1 (),  2 (),  3 () is true.Fig. 7 depicts the change of each compartments over time with or without optimal control.We focus on the number of   () with common strain, and   () with amplified strain.From Fig. 8, we can see that both the number of   () and   () decrease significantly after control, by which we can see that   () tends to 0 under the control measures due to sufficient awareness of the common strain, media coverage can improve the degree of self-protection against   infection, and the vaccine is fully effective against the transmission of multi-strain disease.The amount of   () decreases significantly under the control measures, but did not approach 0. Isolation is the most effective method for those infected individuals with amplified strain.The amplified strain do not become extinct because the weight corresponding to isolation measures  is much larger than the other two measures.We can see that increasing awareness of strains and vaccination are effective weapons to control transmission of multi-strain disease.At the same time, protective measures such as isolation are essential, especially when a new disease emerges or mutated and resistant amplified strain emerge.
(e) Backward bifurcation Let  = 1,  = 0.01,  1 = 0.1,  = 1,   = 0.00108989,   = 0.002.As shown in Fig. 9(a), when  0 < 1, there is an endemic equilibrium.If  0 <   0 = 0.6970, the model has only a globally asymptotically stable disease-free equilibrium; if   0 <  0 < 1, the model has a locally asymptotically stable disease-free equilibrium, an unstable endemic equilibrium, and a locally asymptotically stable endemic equilibrium; if  0 > 1, the model has an unstable disease-free equilibrium and a locally asymptotically stable endemic equilibrium.From Fig. 9(b), we can see that the solutions from different initial values tend to the disease-free equilibrium  0 , where  0 = 0.43001.Fig. 9(c) shows when  0 = 0.96753, the model has a disease-free equilibrium  0 , two endemic equilibria  1 and  2 .From Fig. 9(d), we can see when  0 = 1.075, the model has a disease-free equilibrium  0 , an endemic equilibrium  1 , the solutions starting around  0 and  1 all approach  1 , which means that if  0 < min { 1,   0 } , the disease will tend to eliminate.Mathematically, the disease-free equilibrium is not globally stable when  0 < 1, and biologically, the infection can persist even if the reproduction number is less than one.The phenomenon of backward bifurcation   has serious consequences for disease control, therefore, when  > 0, disease control becomes more difficult.If  < 0 (Form Remark 4, when the vaccination is highly effective, i.e.,  1 = 0, then  = −   < 0), then the disease will die out as long as the basic reproduction number  0 is less than 1 by some measures.Of course, we should consider the optimal control in combination with our control model.As we all know, the COVID-19 is a multi-strain infectious disease.The daily number of existing cases in Morocco during the period from Dec. 28th, 2021 to Jan. 16th, 2022 is selected.Fig. 10 shows that the temporal evolution of existing cases of infection with the parameters shown in Table 3.There is a good agreement between the multi-strain epidemic model (4) and the observed data in Morocco.

Discussion and conclusion
In this paper a multi-strain coinfection model with amplification is established to simulate the interaction between common and amplified strains.The basic reproduction number ℛ 0 = max { ℛ 0 , ℛ 0 }  is derived, the existence and uniqueness of disease-free equilibrium are obtained, when ℛ 0 < 1 it is globally asymptotically stable, when ℛ 0 > 1 it is unstable; the dominant equilibrium of common strain does not exist; If )   ((+)+) .
From Fig. 4, we find that the disease-free equilibrium is unstable, if the coexistence equilibrium does not exist, the dominant equilibrium would be globally asymptotically stable if In another case, from Fig. 5, we can find that the disease-free equilibrium and dominant equilibrium are unstable, and the coexistence equilibrium is globally asymptotic stability.An interesting open problem is how to give precise threshold criteria for global asymptotic −  −   − )  +  −1 ,  = 0, 1, … , ,

Case Condition 1 𝑏 3 >Fig. 1 .
Fig. 1.A flow chart of the spread of multiple infectious strains with vaccine.

H.
Wu et al.

Fig. 6 .
Fig. 6.The value of control measures over time.

H
.Wu et al.

Fig. 7 .
Fig. 7. (a) All the solution sequences of the model, with control.(b) All the solution sequences of the model, without control.

Fig. 8 .
Fig. 8. (a) The number of   (), with and without control.(b) The number of   (), with and without control.(c) Total number of infected, that is the sum of   () and   (), with and without control.(d) Forward bifurcation of model (12), solid lines represent stability, dotted lines represent instability.
(f) Application to COVID-19: fitting the daily number of existing COVID-19 cases in Morocco

Fig. 10 .
Fig. 10.A daily fit of existing confirmed cases in Morocco.