Loading [MathJax]/jax/output/SVG/jax.js
Research article

Real-time PCR assay for detection of Staphylococcus aureus, Panton-Valentine Leucocidin and Methicillin Resistance directly from clinical samples

  • Rapid detection of Methicillin Resistant Staphylococcus aureus (MRSA) is an important concern for both treatment and implementation of infection control policies. The present study provides an ‘in house’ real-time PCR assay to detect directly nuc, pvl, and mecA genes. The assay is able to perform identification of MRSA, Methicillin-Sensitive S. aureus, Methicillin-Resistant Coagulase Negative Staphylococci and the Panton-Valentine leukocidin virulence gene from rectal and pharyngeal swab samples in a screening context. We found an analytical sensitivity of this current Triplex PCR assay of 514 CFU/mL. Analytical specificity was tested with different Gram-positive and Gram-negative species and yielded no false-positive PCR signal. The sensitivity and specificity of the Triplex Real Time PCR were both 100% for these targets when compared with the culture and conventional methods. This assay is readily adaptable for routine use in a microbiology laboratory, as it will enable the implementation of timely and properly guided therapy and infection control strategies.

    Citation: Liliana Galia, Marco Ligozzi, Anna Bertoncelli, Annarita Mazzariol. Real-time PCR assay for detection of Staphylococcus aureus, Panton-Valentine Leucocidin and Methicillin Resistance directly from clinical samples[J]. AIMS Microbiology, 2019, 5(2): 138-146. doi: 10.3934/microbiol.2019.2.138

    Related Papers:

    [1] Xiao Liang, Taiyue Qi, Zhiyi Jin, Wangping Qian . Hybrid support vector machine optimization model for inversion of tunnel transient electromagnetic method. Mathematical Biosciences and Engineering, 2020, 17(4): 3998-4017. doi: 10.3934/mbe.2020221
    [2] Wangping Qian, Yikang Xu, Zhenyuan Gu, Ziru Xiang, Zongyang Li, Jinchao Liu, Jiujiang Wu . Characterisation of transient electromagnetic signals during fixed interference sources in tunnel structure. Mathematical Biosciences and Engineering, 2021, 18(5): 6907-6925. doi: 10.3934/mbe.2021343
    [3] Yikang Xu, Zhaohua Sun, Wei Gu, Wangping Qian, Qiangru Shen, Jian Gong . Three-dimensional inversion analysis of transient electromagnetic response signals of water-bearing abnormal bodies in tunnels based on numerical characteristic parameters. Mathematical Biosciences and Engineering, 2023, 20(1): 1106-1121. doi: 10.3934/mbe.2023051
    [4] Naigong Yu, Hongzheng Li, Qiao Xu . A full-flow inspection method based on machine vision to detect wafer surface defects. Mathematical Biosciences and Engineering, 2023, 20(7): 11821-11846. doi: 10.3934/mbe.2023526
    [5] Xian Fu, Xiao Yang, Ningning Zhang, RuoGu Zhang, Zhuzhu Zhang, Aoqun Jin, Ruiwen Ye, Huiling Zhang . Bearing surface defect detection based on improved convolutional neural network. Mathematical Biosciences and Engineering, 2023, 20(7): 12341-12359. doi: 10.3934/mbe.2023549
    [6] Jianzhong Peng, Wei Zhu, Qiaokang Liang, Zhengwei Li, Maoying Lu, Wei Sun, Yaonan Wang . Defect detection in code characters with complex backgrounds based on BBE. Mathematical Biosciences and Engineering, 2021, 18(4): 3755-3780. doi: 10.3934/mbe.2021189
    [7] Qingfu Li, Huade Zhou, Hua Zhang . Durability evaluation of highway tunnel lining structure based on matter element extension-simple correlation function method-cloud model: A case study. Mathematical Biosciences and Engineering, 2021, 18(4): 4027-4054. doi: 10.3934/mbe.2021202
    [8] Venkat Anil Adibhatla, Huan-Chuang Chih, Chi-Chang Hsu, Joseph Cheng, Maysam F. Abbod, Jiann-Shing Shieh . Applying deep learning to defect detection in printed circuit boards via a newest model of you-only-look-once. Mathematical Biosciences and Engineering, 2021, 18(4): 4411-4428. doi: 10.3934/mbe.2021223
    [9] Bo Wu, Weixing Qiu, Wei Huang, Guowang Meng, Jingsong Huang, Shixiang Xu . Dynamic risk evaluation method for collapse disasters of drill-and-blast tunnels: a case study. Mathematical Biosciences and Engineering, 2022, 19(1): 309-330. doi: 10.3934/mbe.2022016
    [10] Xinglong Yin, Lei Liu, Huaxiao Liu, Qi Wu . Heterogeneous cross-project defect prediction with multiple source projects based on transfer learning. Mathematical Biosciences and Engineering, 2020, 17(2): 1020-1040. doi: 10.3934/mbe.2020054
  • Rapid detection of Methicillin Resistant Staphylococcus aureus (MRSA) is an important concern for both treatment and implementation of infection control policies. The present study provides an ‘in house’ real-time PCR assay to detect directly nuc, pvl, and mecA genes. The assay is able to perform identification of MRSA, Methicillin-Sensitive S. aureus, Methicillin-Resistant Coagulase Negative Staphylococci and the Panton-Valentine leukocidin virulence gene from rectal and pharyngeal swab samples in a screening context. We found an analytical sensitivity of this current Triplex PCR assay of 514 CFU/mL. Analytical specificity was tested with different Gram-positive and Gram-negative species and yielded no false-positive PCR signal. The sensitivity and specificity of the Triplex Real Time PCR were both 100% for these targets when compared with the culture and conventional methods. This assay is readily adaptable for routine use in a microbiology laboratory, as it will enable the implementation of timely and properly guided therapy and infection control strategies.


    The COVID-19 has been the most serious outbreak worldwide since 2000, and it is currently in the prevention and control phase. The novel coronavirus is in lineage B of the genus $ \beta $-coronavirus of the coronavirus family, in which severe acute respiratory syndrome-related coronavirus (SARS-CoV) and Middle East respiratory syndrome-related coronavirus (MERS-CoV) are also included. These pathogens are enveloped positive-sense RNA viruses that are widespread in mammals, including humans, and they can destroy the respiratory systems of humans and cause severe acute respiratory syndrome [1,2,3]. Despite global vaccination efforts, the epidemic situation is still critical for two main reasons. One reason is that the COVID-19 vaccination proportion remains low in some countries, such as low-income countries in Africa [4], leading to a lack of global immune barriers. Another issue that cannot be ignored is that the antibodies produced in the body do not fully defeat variable strains, and the effectiveness of antibodies decreases as the virus mutates [5]. According to data published by the Israeli Ministry of Health, the protective effect of the Pfizer COVID-19 vaccine decreases after 6 months of vaccination. Based on available statistics, it has been confirmed that some people who have been vaccinated against COVID-19 have been diagnosed, and this number is increasing. Therefore, there has been a global consensus to implement the booster vaccinations. By receiving the booster vaccine, a strong immune response can be rapidly produced in the body, and the level of neutralizing antibodies can significantly increase, thereby improving the overall immunity level and achieving outbreak control.

    In 1927, Kermack and McKendrick first used dynamic research methods to study the spread of infectious diseases [6]. Since then, an increasing number of researchers have engaged in mathematical modeling of infectious diseases considering their transmission mechanisms and control measures [7,8,9,10]. Vaccination and even booster shots were added to the model to study the dynamic behaviors of the corresponding infectious diseases [11,12,13]. In particular, since the outbreak of COVID-19, researchers worldwide have been very concerned, and corresponding infectious disease models have been developed to simulate epidemic trends to provide policy recommendations [14,15,16,17,18,19]. The models are not limited to ordinary differential equations (ODEs), and the use of partial differential equations (PDEs) is also a good way to study such problems [20,21,22].

    During the early stage of the outbreak, the treatment process was difficult because people had never been exposed to this disease. Although we have better research on treatments and have continued to learn about the disease, studies have shown that COVID-19 drugs have side effects on our bodies [23,24]. Therefore, vaccination is considered an effective strategy to prevent COVID-19 and control the outbreak. Thus, as vaccines gradually enter clinical trials and continue to spread worldwide, an increasing number of scholars have added vaccination to their models[25,26]. Given the current situation, most people have been vaccinated, but the mutations may lead to the failure of original antibodies. To maintain sustained immunity, it is feasible and meaningful to vaccinate with booster shots [27].

    Adding the necessary "time delay" to the model tends to make it more realistic. Sometimes, the model is sufficient to describe facts and phenomena accurately even with only one time delay [28,29,30,31]. However, sometimes, one time delay is not enough to describe complex phenomena fully, and double time delays or even multiple time delays are necessary [32,33,34,35,36,37,38]. When we add the effect of booster vaccination to the epidemic model, there are two periods of time that we cannot ignore: the time of antibody failure and the time of booster vaccination. Therefore, we introduce two time delays into the model to describe the problem better.

    For differential equation models, stability and Hopf bifurcation analysis have always been popular research topics. Mahata et al. [39] studied a fractional-order dynamical system of Susceptible-Exposed-Infected-Recovered-Vaccinated ($ SEIRV $) in the infectious population. The stability of all equilibria was analyzed with respect to the delay parameter, and threshold values of delay were found, beyond which the system exhibited Hopf bifurcation and the solutions were no longer periodic. Zhang et al. [40] studied an epidemic model involving nonlinear birth in population and vertical transmission. They discussed the existence of Hopf bifurcation and obtained the stability and direction of Hopf bifurcation by using normal form theory and center manifold theorem. References [41,42,43,45,46,47] discussed different dynamic properties of Hopf bifurcation in different systems. We hope that we can also analyze the dynamic properties after establishing the model and use the results to make reasonable recommendations for vaccination.

    The motivations of our study are as follows. First, although COVID-19 vaccines have shown good protective efficacy and great safety in both clinical trials and real-life applications, there are some cases of infection in people who have been vaccinated, suggesting that vaccine-induced immunity is not permanent. Therefore, booster vaccination is necessary and we need to take this into account in our model. Second, there is a time delay between booster vaccination and initial vaccination, which has a profound impact on the development of the epidemic. Different booster vaccination time will cause different phenomena in the model, which is essential to understand the development of epidemic and to formulate the effective policies. Third, when is the booster vaccination most beneficial for the development of the epidemic when taking into account the economy, national power and medical level? What measures should we take to minimize the impact of the epidemic on human life? Finally, COVID-19 exhibits a constantly and dynamically evolving mutational landscape, with relatively abundant genetic diversity and a high evolutionary capability over time[48,49,50]. As the virus continues to mutate, what effect will this have on the time of booster vaccination?

    Through the study of this paper, we have a more accurate grasp of when to carry out booster vaccination. By simulating the current epidemic state, a reasonable time range for booster vaccination is proposed, which can ensure that the epidemic can reach a stable state instead of getting out of control. This is meaningful to the formulation of epidemic prevention and control policies.

    This paper is organized as follows. In Section 2, we establish a delayed differential system with two time delays to study the impact of vaccination on the epidemic and determine the optimal time for the booster vaccination. In Section 3, we analyze the conditions for the existence and stability of equilibria and the existence of Hopf bifurcation. Next, the normal form of Hopf bifurcation is deduced by using the multiple time scales method in Section 4. In Section 5, numerical simulations are carried out to verify the theoretical analysis based on the parameters through statistical analysis. Finally, we draw conclusions about the optimal time range of booster vaccination by combining all aspects and provide some reasonable suggestions in Section 6.

    In the traditional model of infectious disease, susceptible people ($ S $) can be infected through contact with a source of infection, which converts them into infected people ($ I $). With treatment, those infected people can be cured and then become recovered people ($ R $) with antibodies. Considering the characteristics of COVID-19, vaccination has become the global trend, and vaccinated people have become a large group that cannot be ignored. Therefore, we need to represent them as a separate compartment denoted as vaccinated people ($ V $), which is a more accurate representation of the current outbreak. However, it is important to note that vaccinated people ($ V $) are those who have completed the vaccination but have not fully developed antibodies. This group of people is still at risk of being infected. For those who have fully developed antibodies, we classify them as recovered people ($ R $) since the process of producing antibodies is similar to the mechanism of cure after infection. This can be understood as recovering from the risk of infection. Combined with the above considerations, we choose the Susceptible-Vaccinated-Infective-Recovered ($ SVIR $) model for our study of COVID-19 with booster vaccination.

    Similar to the traditional $ SIR $ model, susceptible people ($ S $) will be transformed into infected people ($ I $) when they effectively contact $ I $. Infected people ($ I $) will be transformed into recovered people ($ R $) after healing through treatment. As schematically shown in Figure 1, after vaccination, susceptible people ($ S $) can be transformed into vaccinated people ($ V $). However, not all vaccinated people are able to complete the immune response and further produce antibodies. For those who are able to produce the immune response ($ mV $), we also classify them as recovered people ($ R $) since they share the same characteristics as those who are cured after infection. By contrast, those who cannot produce antibodies due to individual variability ($ (1-m)V $) are clearly at risk of infection with an infection rate of $ k $. Therefore, the term to transfer from $ V $ to $ I $ is $ k(1-m)VI $. In the case of COVID-19, even though antibodies have been produced in the body, this immunity is temporary. After a period of time, recovered people ($ R $) are faced with losing immunity, thus returning to susceptible people. We call this "antibody failure". However, if the booster vaccination is given before antibody failure, there is a return to vaccinated people ($ V $) and thus a chance to reproduce the immune response.

    Figure 1.  Schematic diagram of the $ SVIR $ model for COVID-19.

    In our model, we consider antibody failure in two main ways, partly because antibody levels decrease significantly over time and partly because the efficiency of existing antibodies against mutated strains is greatly reduced due to virus mutations. Both of these conditions increase the risk of infection. At the time of booster vaccination, the antibody level decreases compared to before and mutated strains are produced. Therefore, we believe that the current antibodies are much less efficient at the time of booster vaccination. For different mutated strains, booster shots are developed for their specific characteristics. Due to individual differences, it may occur that some people produce antibodies completely with the first vaccination but do not complete the immune response after booster vaccination, increasing the risk of infection. The specific definitions of variables and parameters are given in Table 1, and the rates are all expressed as the average over a year.

    Table 1.  Descriptions of parameters and variables.
    Symbol Descriptions Unit
    $ S $ Number of suspected people $ 10^5 $ persons
    $ V $ Number of vaccinated people $ 10^5 $ persons
    $ I $ Number of infected people $ 10^5 $ persons
    $ R $ Number of recovered people $ 10^5 $ persons
    $ \Lambda $ Population nature growth $ 10^5 $ persons
    $ d $ Population natural mortality rate
    $ c $ COVID-19 mortality rate of $ I $
    $ \beta $ Infection rate from $ S $ to $ I $
    $ k $ Proportional coefficient of transmission
    $ \alpha $ Vaccination rate
    $ m $ Transmission rate from $ V $ to $ R $
    $ \rho $ Booster vaccination rate
    $ \mu $ Vaccine failure rate
    $ u $ Cure rate

     | Show Table
    DownLoad: CSV

    Since we want to know when the booster vaccination is necessary and when it is most effective, we let $ \tau_2 $ represent the time delay between $ R $ and $ V $ for booster vaccination and simulate different vaccination time by changing the size of $ \tau_2 $. However, in the case of COVID-19 vaccines, the vaccines are nonpermanent, that is, antibodies produced by the vaccines have a problem of failure. The failure time is not negligible relative to $ \tau_2 $, so we use $ \tau_1 $ to express the time for antibody failure. In order to improve herd immunity and control the epidemic better, booster vaccination usually occurs before antibody failure, so we assume $ \tau_2 < \tau_1 $.

    Therefore, the $ SVIR $ epidemic model is constructed to explore the optimal booster vaccination time. The equation is presented as follows:

    $ {dSdt=Λ+μR(tτ1)αSβSIdS,dVdt=αS+ρR(tτ2)mVdVk(1m)VI,dIdt=βSIdIuIcI+k(1m)VI,dRdt=mV+uIμR(tτ1)ρR(tτ2)dR, $ (2.1)

    where $ \tau_1 $ is the time delay of antibody failure, $ \tau_2 $ is the time delay of booster vaccination and the specific definitions of variables and parameters are given in Table 1.

    For convenience of calculation, we let $ \gamma\triangleq k(1-m) $ in the following text.

    The initial condition of system (2.1) is given by $ S (\theta) = {\phi _1}(\theta), V(\theta) = {\phi _2}(\theta), I(\theta) = {\phi _3}(\theta), R(\theta) = {\phi _4}(\theta) $, $ \theta\in[-\tau, 0] $, $ \tau = \max{\{\tau_1, \tau_2\}} $, with $ \phi = [\phi _1, \phi _2, \phi _3, \phi _4]\in \mathbb{C} $ such that $ \phi _i\ge 0, i = 1, 2, 3, 4 $ for $ \theta\in[-\tau, 0] $ where $ \mathbb{C} $ denotes the Banach space of continuous functions mapping the interval $ [-\tau, 0] $ into $ \mathbb{R}_{ + 0}^4 = \{(S, V, I, R)|S > 0, V > 0, I > 0, R > 0\} $.

    When $ t\in[0, \tau] $, the solution of system (2.1) is as follows:

    $ S(t)=et0(α+d+βI)dξ[S(0)+t0(Λ+μR(ητ1))eη0(α+d+βI)dξdη],V(t)=et0(m+d+γV)dξ[A(0)+t0(αS(η)+ρR(ητ2))eη0(m+d+γV)dξdη],I(t)=et0(βSγV+d+u+c)dξI(0),R(t)=e(d)t[R(0)+t0(mV(η)+uI(η)μR(ητ1)ρR(ητ2))edηdη]. $ (2.2)

    Obviously, $ I(t) > 0 $. In the current state, the open policy leads to an increase in the number of infected people. The mutation of the virus makes many people choose booster vaccination when they face body failure. Therefore, the numbers of $ V $ and $ I $ might be higher. Although $ R(t) $ is large, $ R(t-\tau_1) $ and $ R(t-\tau_2) $ will not be large compared to the current values since the time delays are in years as a unit. Therefore, it is reasonable to assume that $ mV+uI-\mu R(t-\tau_1)- \rho R(t-\tau_2) > 0 $. Thus, $ R(t) > 0 $ based on the expression of $ R(t) $ in Eq (2.2). Furthermore, we can get $ R(t-\tau_1) > 0, R(t-\tau_2) > 0 $ when $ t\in[0, \tau] $. Thus, $ S(t) > 0, V(t) > 0 $.

    Overall, when $ t\in[0, \tau] $, all the solutions of system (2.1) are positive after we consider the realistic situation and make some assumptions. Through the similar approach, it can be shown that all solutions of system (2.1) are also positive on $ [\tau, 2\tau] $ under some restrictions. This conclusion can be extended to the interval $ [n\tau, (n+1)\tau]\; (n\in N) $. Therefore, for any $ t > 0 $, all solutions of system (2.1) are positive with the restrictions.

    Although the positivity of our model is conditional, our analysis of the epidemic shows that the present state satisfies the conditions for the existence of positivity, therefore our model is reasonable and meaningful.

    Let $ N(t) = S(t)+V(t)+I(t)+R(t) $, and $ N(t) $ represents the total number of people at time $ t $. Adding four equations of system (2.1), we obtain $ N'(t) = \Lambda-dN(t)-cI(t) $. According to the discussion on the positivity of solutions, we can find that $ I(t) > 0 $. Thus, $ N'(t) = \Lambda-dN(t)-cI(t)\leq \Lambda-dN(t) $. Based on the comparison theorem, we can obtain:

    $ 0 < N(t)\leq\mathrm{e}^{-dt} \left[N(0)+\Lambda\int_{0}^{t}\mathrm{e}^{d\eta}d\eta\right] = \mathrm{e}^{-dt}[N(0)-\frac{\Lambda}{d}]+ \frac{\Lambda }{d }. $

    Therefore, $ \mathop {\lim }\limits_{t \to \infty } \text{sup}\; N(t) \leq \frac{\Lambda }{d } $ and the solutions $ S(t) $, $ A(t) $, $ I(t) $, $ R(t) $ of system (2.1) are bounded when $ t > 0 $.

    The boundedness of the solutions of system (2.1) must hold. Its practical significance is that the number of people in each cabin can be stabilized within a specific range regardless of the current epidemic situation, which is optimistic for controlling the epidemic.

    Under the positivity and boundedness of solutions discussed above, we mark the feasible region of the system (2.1) as:

    $ Ω={S(t),A(t),I(t),R(t)|N(t)=S(t)+I(t)+Q(t)+R(t)Λd,S(t)0,A(t)0,I(t)0,R(t)0,mV+uIμR(tτ1)ρR(tτ2)>0.} $ (2.3)

    Next, we will analyze system (2.1). It has three equilibria,

    $ E1=(S1,V1,I1,R1),E2=(S2,V2,I2,R2),E3=(S3,V3,I3,R3), $ (3.1)

    with

    $ S1=Λ[m(μ+d)+d(u+d+ρ)]Π>0,V1=(μ+ρ+d)αΛΠ>0,I1=0,R1=αmΛΠ>0,S2=d+u+cβγβV2,V2=A+A24γpB2γp,I2=pV2+q,R2=mV2+uI2μ+ρ+d,S3=d+u+cβγβV3,V3=AA24γpB2γp,I3=pV3+q,R3=mV3+uI3μ+ρ+d,Π=αmμ+d(α+d)(m+u+d+ρ)+mdμ,q=[Λβd(d+u+c)](μ+ρ+d)βu(μ+ρ)+β(d+u+c)(μ+ρ+d),p=β(μ+ρ)mβ(m+d)(μ+ρ+d)+dγ(μ+ρ+d)βu(μ+ρ)+β(d+u+c)(μ+ρ+d),A=αγβρmμ+ρ+dρupμ+ρ+d+m+d+γq,B=α(d+u+c)βρuqμ+ρ+d. $

    Then, we show the following assumptions:

    $ {\bf (H1)}\; \; V_2^* = \frac{-A+\sqrt{A^2-4\gamma pB}}{2rp} > 0, S_2^* = \frac{d+u+c}{\beta}-\frac{\gamma}{\beta}V_2^* > 0, I_2^* = pV_2^*+q > 0, \\ {\bf (H2)}\; \; V_3^* = \frac{-A-\sqrt{A^2-4\gamma pB}}{2rp} > 0, S_3^* = \frac{d+u+c}{\beta}-\frac{\gamma}{\beta}V_3^* > 0, I_3^* = pV_3^*+q > 0,$

    where $ A, B, p $ are given in Eq (3.1). We can easily find that there is always an equilibrium $ E_1 $. If equilibrium $ E_1 $ is stable, there are no infected people in the final stable state, that is, the virus can finally be eliminated, which is called the disease-free equilibrium. Then, equilibrium $ E_2 $ makes sense under (H1) and equilibrium $ E_3 $ is meaningful under (H2). If equilibrium $ E_2 $ or $ E_3 $ is stable, the final state is always accompanied by some infected people, which means that the virus will always coexist with humans and we have no way to eliminate it. We call it the endemic equilibrium.

    Without loss of generality, we assume that the system (2.1) has three equilibria $ E_k\; (k = 1, 2, 3) $. Then, we transfer the equilibrium $ E_k $ into the original point: $ \tilde{S} = S-S_{k}^{*} $, $ \tilde{V} = V-V_{k}^{*} $, $ \tilde{I} = I-I_{k}^{*} $, $ \tilde{R} = R-R_{k}^{*} $ and re-denote $ \tilde{S}, \tilde{V}, \tilde{I}, \tilde{R} $ as $ S, V, I, R $. We obtain the following model:

    $ {dSdt=Λ+μR(tτ1)αSβSIβSkIβSIkdS,dVdt=αS+ρR(tτ2)mVdVγVIγVkIγVIk,dIdt=βSI+βSkI+βSIkdIuIcI+γVI+γVkI+γVIk,dRdt=mV+uIμR(tτ1)ρR(tτ2)dR. $ (3.2)

    For the equilibrium $ E_1 = (S_1^*, V_1^*, I_1^*, R_1^*) $, we can obtain the characteristic equation of system (3.2) as follows:

    $ (λβS1γV1+d+u+c)(λ+d)[(λ+m+d)(λ+α+d)+(λ+m+α+d)μeλτ1+(λ+α+d)ρeλτ2]=0. $ (3.3)

    Similar to the equilibrium $ E_1 $, we can also obtain the characteristic equation at the equilibrium $ E_k\; (k = 2, 3) $:

    $ λ4+ϱ1,kλ3+ς2,kλ2+ϱ3,kλ+ϱ4,k+μeλτ1(λ3+σ1,kλ2+σ2,kλ+σ3,k)+ρeλτ2(λ3+ς1,kλ2+ς2,kλ+ς3,k)=0, $ (3.4)

    where

    $ ϱ1,k=d+Φ2,k,ϱ2,k=dΦ2,k+Φ1,k,ϱ3,k=dΦ1,k+Φ0,k,ϱ4,k=dΦ0,k,σ1,k=Φ2,k,σ2,k=Φ1,k+Ψ1,k,σ3,k=Φ0,k+Ψ0,k,ς1,k=Φ2,km,ς2,k=Φ1,kΩ1,k,ς3=Φ0,kΩ0,k,Φ2,k=m+d+γIkβSkγVk+d+u+c+α+d+βIk,Φ1,k=(m+d+γIkβSkγVk+d+u+c)(α+d+βIk)+(m+d+γIk)(βSkγVk+d+u+c)+γ2IkVk+β2SkIk,Φ0,k=(α+d+βIk)[(m+d+γIk)(βSkγVk+d+u+c)+γ2IkVk]+βSk[αγIkβIk(m+d+γIk)],Ψ1,k=αmuβIk,Ψ0,k=m[α(βSkγVk+d+u+c)+βγVkIk]u[αγIk+βIk(m+d+γIk)],Ω1,k=m(α+d+βIkβSkγVk+d+u+c)+uγIk,Ω0,k=m[(α+d+βIk)(βSkγVk+d+u+c)+β2SkIk]+uγIk(α+d+βIk),k=2,3. $

    When $ \tau_1 = \tau_2 = 0 $, (3.3) becomes:

    $ (λβS1γV1+d+u+c)(λ+d)[λ2+(m+α+2d+μ+ρ)λ+(m+d)(α+d)+μ(m+α+d)+ρ(α+d)]=0. $ (3.5)

    Clearly, Equation (3.5) has the roots $ \lambda_1 = \beta S_1^*+\gamma V_1^*-(d+u+c) $, $ \lambda_2 = -d $, and the remaining roots $ \lambda_3, \lambda_4 $ are given by the following equation:

    $ λ2+(m+α+2d+μ+ρ)λ+(m+d)(α+d)+μ(m+α+d)+ρ(α+d)=0. $ (3.6)

    According to Vieta's theorem, we can obtain $ \lambda_3+\lambda_4 = -(m+\alpha+2d+\mu+\rho) < 0 $ and $ \lambda_3\lambda_4 = (m+d)(\alpha+d)+\mu(m+\alpha+d)+\rho(\alpha+d) > 0 $, thus, $ \mathrm{Re}(\lambda_3) < 0 $ and $ \mathrm{Re}(\lambda_4) < 0 $. Then, we show the following assumption:

    $ {\bf (H3)}\; \; \beta S_1^*+\gamma V_1^*-(d+u+c) < 0.$

    Under (H3), we can find that all the characteristic roots of Eq (3.5) are negative, so the equilibrium $ E_1 $ is locally asymptotically stable when $ \tau_1 = \tau_2 = 0 $. If (H3) does not hold, $ \lambda_1 = \beta S_1^*+\gamma V_1^*-(d+u+c) > 0 $, the equilibrium $ E_1 $ is unstable when $ \tau_1 = \tau_2 = 0 $.

    For the equilibrium $ E_k\; (k = 2, 3) $, Equation (3.4) is transformed to the following form when $ \tau_1 = \tau_2 = 0 $:

    $ λ4+a1,kλ3+a2,kλ2+a3,kλ+a4,k=0, $ (3.7)

    where $ a_{1, k} = \varrho_{1, k}+\mu+\rho, \; a_{2, k} = \varrho_{2, k}+\mu\sigma_{1, k}+\rho\varsigma_{1, k}, \; a_{3, k} = \varrho_{3, k}+\mu\sigma_{2, k}+\rho\varsigma_{2, k}, \; a_{4, k} = \varrho_{4, k}+\mu\sigma_{3, k}+\rho\varsigma_{3, k} $ with $ \varrho_{1, k} $, $ \varrho_{2, k} $, $ \varrho_{3, k} $, $ \varrho_{4, k} $, $ \sigma_{1, k} $, $ \sigma_{2, k} $, $ \sigma_{3, k} $, $ \varsigma_{1, k} $, $ \varsigma_{2, k} $, $ \varsigma_{3, k} $ given in Eq (3.4). According to Routh-Hurwitz stability criterion, we propose the following hypothesis:

    $ {\bf (H4)} \; \; a_{1, k}a_{2, k}-a_{3, k} > 0, a_{3, k}\left(a_{1, k}a_{2, k}-a_{3, k}\right) > a_{1, k}^{2}a_{4, k}, a_{4, k} > 0. $

    If (H4) is satisfied, all eigenvalues of Eq (3.7) have negative real parts, the equilibrium $ E_k\; (k = 2, 3) $ of model (2.1) is locally asymptotically stable when $ \tau_1 = \tau_2 = 0 $.

    Theorem 3.1. For system (2.1) with $ \tau_1 = 0, \tau_2 = 0 $, the stability results for equilibrium $ E_k\; (k = 1, 2, 3) $ are given as follows.

    1) Equilibrium $ E_1 $ always exists. If (H3) holds, equilibrium $ E_1 $ is locally asymptotically stable. On the contrary, if (H3) does not hold, equilibrium $ E_1 $ is unstable;

    2) When the assumption (H1) or (H2) holds, the equilibrium $ E_2 $ or $ E_3 $ is meaningful respectively. Further, if (H4) holds, equilibrium $ E_2 $ or $ E_3 $ is locally asymptotically stable. If (H4) does not hold, equilibrium $ E_2 $ and $ E_3 $ are unstable.

    For the equilibrium $ E_1 $, when $ \tau_2 = 0, \tau_1 > 0 $, the characteristic equation of system (3.3) becomes:

    $ (λβS1γV1+d+u+c)(λ+d)[(λ+m+d+ρ)(λ+α+d)+(λ+m+α+d)μeλτ1]=0. $ (3.8)

    If (H3) holds, $ \lambda_1 = \beta S_1^*+\gamma V_1^*-(d+u+c) < 0 $, $ \lambda_2 = -d < 0 $. We only need to consider the remaining part of Eq (3.8):

    $ (λ+m+d+ρ)(λ+α+d)+(λ+m+α+d)μeλτ1=0. $ (3.9)

    Substituting $ \lambda = \mathrm{i}\omega\; (\omega > 0) $ into Eq (3.9) and separating the real and imaginary parts, we obtain:

    $ {μωsin(ωτ1)+μ(m+α+d)cos(ωτ1)=ω2(m+d+ρ)(α+d),μ(m+α+d)sin(ωτ1)μωcos(ωτ1)=ω(α+m+2d+ρ). $ (3.10)

    Thus,

    $ {sin(ωτ1)=ω3(m+α+d)(α+d)ω+(α+m+2d+ρ)(m+α+d)ωμ(m+α+d)2+μω2,cos(ωτ1)=ω2(d+ρ)(m+α+d)(m+d+ρ)(α+d)μ(m+α+d)2+μω2. $ (3.11)

    Adding the square of the two equations in Eq (3.11) and letting $ z = \omega^2 $, we obtain:

    $ h(z)=z2+κ1z+κ0, $ (3.12)

    where $ \kappa_1 = -2(m+d+\rho)(\alpha+d)+(\alpha+m+2d+\rho)^2-\mu^2, \kappa_0 = (m+d+\rho)^2(\alpha+d)^2-mu^2(m+\alpha+d)^2. $ Therefore, we show the following assumptions:

    $(H5)κ0<0,(H6)κ214κ0>0,κ1<0,κ0>0.$

    If (H5) and (H6) do not hold, the equilibrium $ E_1 $ is always locally asymptotically stable for any $ \tau_1 > 0 $. If (H5) holds, Equation (3.12) has the unique positive root $ z_1 $, $ h'(z_1) > 0 $. Under (H6), Equation (3.12) has two positive roots: $ z_2 $ and $ z_3 $. Assuming $ z_2 > z_3 $, we get $ h'(z_2) > 0, \; h'(z_3) < 0 $. By substituting $ \omega_{11, l} = \sqrt{z_l}\; (l = 1, 2, 3) $ into Eq (3.11), we can obtain:

    $ τ(j)11,l={1ω11,l[arccos(P11,l)+2jπ],           Q11,l0,1ω11,l[2πarccos(P11,l)+2jπ],     Q11,l<0, $ (3.13)

    where $ Q_{11, l} = {\sin(\omega\tau_1)}|_{\omega = {\omega _{11, l}}, \tau = \tau _{11, l}^{(j)}}, {P_{11, l}} = \cos(\omega\tau _{1}) |_{\omega = {\omega _{11, l}}, \tau = \tau _{11, l}^{(j)}}, $ $ \sin(\omega\tau_1) $ and $ \cos(\omega\tau _1) $ are given in Eq (3.11).

    Then the following transversality condition holds:

    $ \left. \mathrm{Re}\left(\frac{\mathrm{d}\lambda }{\mathrm{d}\tau_1 }\right)^{ - 1} \right|_{\tau_1 = \tau _{11, l}^{\left( j \right)}} = \left.\mathrm {Re}\left( \frac{\mathrm{d}\tau_1}{\mathrm{d}\lambda} \right)\right|_{\tau_1 = \tau _{11, l}^{\left( j \right)}} = \frac{h'\left( \omega_{11, l}^{2} \right)}{\mu^2\left[ \omega_{11, l}^2+\left(m+\alpha+d\right)^2 \right]}\neq0\; \left(j = 0, 1, 2, \cdots\right).$

    For the equilibrium $ E_k\; (k = 2, 3) $, Equation (3.4) becomes the following form when $ \tau_1 > 0 $ and $ \tau_2 = 0 $.

    $ λ4+υ1,kλ3+υ2,kλ2+υ3,kλ+υ4,k+μeλτ2(λ3+σ1,kλ2+σ2,kλ+σ3,k)=0, $ (3.14)

    where $ \upsilon_{1, k} = \varrho_{1, k}+\rho, \upsilon_{2, k} = \varrho_{2, k}+\rho\varsigma_{1, k}, \upsilon_{3, k} = \varrho_{3, k}+\rho\varsigma_{2, k}, \upsilon_{4, k} = \varrho_{4, k}+\rho\varsigma_{3, k} $ with $ \varrho_{1, k} $, $ \varrho_{2, k} $, $ \varrho_{3, k} $, $ \varrho_{4, k} $, $ \sigma_{1, k} $, $ \sigma_{2, k} $, $ \sigma_{3, k} $, $ \varsigma_{1, k} $, $ \varsigma_{2, k} $, $ \varsigma_{3, k} $ given in Eq (3.4).

    Assuming that $ \lambda = \mathrm{i}\omega \left(\omega > 0\right) $ is a pure imaginary root of Eq (3.14), we substitute it into Eq (3.14) and separate the real and imaginary parts, and we have:

    $ {ω4υ2,kω2+υ4,k=μ(σ1,kω2+σ3,k)cos(ωτ1)μ(ω3+σ2,kω)sin(ωτ1),υ1,kω3+υ3,kω=μ(σ1,kω2+σ3,k)sin(ωτ1)μ(ω3+σ2,kω)cos(ωτ1). $ (3.15)

    Thus,

    $ {sin(ωτ1)=(υ1,kω3+υ3,kω)(σ1,kω2+σ3,k)(ω4υ2,kω2+υ4,k)(ω3+σ2,kω)μ(σ1,kω2+σ3,k)2+μ(ω3+σ2,kω)2,cos(ωτ1)=(ω4υ2,kω2+υ4,k)(σ1,kω2+σ3,k)+(υ1,kω3+υ3,kω)(ω3+σ2,kω)μ(σ1,kω2+σ3,k)2+μ(ω3+σ2,kω)2. $ (3.16)

    Adding the square of the two equations in Eq (3.16) and letting $ z = \omega^2 $, we get:

    $ h(z)=z4+˜κ1,kz3+˜κ2,kz2+˜κ3,kz+˜κ4,k, $ (3.17)

    where $ \tilde{\kappa}_{1, k} = -2\upsilon_{2, k}+\upsilon_{1, k}^2-\mu^2, \; \tilde{\kappa}_{2, k} = \upsilon_{2, k}^2+2\upsilon_{4, k}-2\upsilon_{1, k}\upsilon_{3, k}-\mu^2\left(\sigma_{1, k}^2-2\sigma_{2, k}\right), $$\tilde{\kappa}_{3, k} = -2\upsilon_{2, k}\upsilon_{4, k}+\upsilon_{3, k}^2-\mu^2\left(-2\sigma_{1, k}\sigma_{3, k}+\sigma_{2, k}^2\right), \; \tilde{\kappa}_{4, k} = \upsilon_{4, k}^2-\mu^2\sigma_{3, k}^2 $ with $ \sigma_{1, k}, \sigma_{2, k}, \sigma_{3, k} $ given in Eq (3.4) and $ \upsilon_{1, k}, \upsilon_{2, k}, \upsilon_{3, k}, \upsilon_{4, k} $ given in Eq (3.14). We hypothesize that Eq (3.17) has $ l \left(l = 1, 2, 3, 4\right) $ positive roots marked as $ z_1 > z_2 > \dots > z_{l} $, $ \left(l = 1, 2, 3, 4\right) $, so $ h'(z_l)\neq0 $. Substituting $ \omega_{1k, l} = \sqrt{z_{l}} $ into Eq (3.16), we can get the expression of $ \tau_{1k, l}^{(j)} $:

    $ τ(j)1k,l={1ω1k,l[arccos(P1k,l)+2jπ],           Q1k,l0,1ω1k,l[2πarccos(P1k,l)+2jπ],     Q1k,l<0, $ (3.18)

    where $ Q_{1k, l} = {\sin(\omega\tau_1)}|_{\omega = {\omega _{1k, l}}, \tau = \tau _{1k, l}^{(j)}}, {P_{1k, l}} = \cos(\omega\tau _{1}) |_{\omega = {\omega _{1k, l}}, \tau = \tau _{1k, l}^{(j)}}, $ $ \sin(\omega\tau_1) $ and $ \cos(\omega\tau _1) $ are given in Eq (3.16).

    Thus, we obtain the transversality condition:

    $ \left. \mathrm{Re}\left(\frac{\mathrm{d}\lambda }{\mathrm{d}\tau_1}\right)^{ - 1} \right|_{\tau_1 = \tau _{1k, l}^{(j)}} = \left.\mathrm{ Re}\left( \frac{d\tau_1}{d\lambda} \right)\right|_{\tau_1 = \tau_{1k, l}^{\left(j\right)}} = \frac{h'\left( \omega_{1k, l}^{2} \right)}{\mu ^2\left[ \left( -\sigma_{1, k}\omega_{1k, l}^{2}+\sigma_{3, k} \right)^2+\left( -\omega_{1k, l}^{3}+\sigma_{2, k}\omega_{1k, l} \right)^2 \right]}\neq0\; \left(j = 0, 1, 2, \cdots\right).$

    Theorem 3.2. For system (2.1) with $ \tau_2 = 0, \tau_1 > 0 $, the following conclusions hold when the stability conditions of Theorem 3.1 hold.

    1) Equilibrium $ E_1 $

    (a) If (H5) and (H6) do not hold, the equilibrium $ E_1 $ is locally asymptotically stable when $ \tau_1 > 0 $.

    (b) If (H5) holds, $ h(z) $ has only one positive root $ z_1 $, the system (2.1) undergoes Hopf bifurcation at $ E_1 $ when $ \tau_1 = \tau_{11, 1}^{(0)} $, where $ \tau_{11, 1}^{(0)} $ is given by Eq (3.13). If $ \tau_1 \in [0, \tau_{11, 1}^{(0)}) $, the equilibrium $ E_1 $ is locally asymptotically stable. When $ \tau_1 \in [\tau_{11, 1}^{(0)}, +\infty) $, the equilibrium $ E_1 $ is unstable.

    (c) If (H6) holds, $ h(z) $ has two positive roots $ z_2 $ and $ z_3 $, and system (2.1) undergoes Hopf bifurcation at $ E_1 $ when $ \tau_1 = \tau_{11, 2}^{(0)}, \tau_{11, 3}^{(0)} $, where $ \tau_{11, 2}^{(0)} $ and $ \tau_{11, 3}^{(0)} $ are given by Eq (3.18). Supposing $ z_2 > z_3 $, we get $ h'(z_2) > 0, h'(z_3) < 0 $. Then, $ \exists\; m \in N $, which can make $ 0 < \tau _{11, 2}^{(0)} < \tau _{11, 3}^{(0)} < \tau _{11, 2}^{(1)} < \tau _{11, 3}^{(1)} < \cdots < \tau _{11, 2}^{(m)} < \tau _{11, 2}^{(m+1)} $. When $ \tau_1 \in(0, \tau _{11, 2}^{(0)})\cup \bigcup_{i = 1}^{m}{(\tau _{11, 3}^{(i-1)}, \tau _{11, 2}^{(i)})} $, the equilibrium $ E_1 $ of the model (2.1) is locally asymptotically stable. When $ \tau_1 \in \bigcup_{i = 0}^{m-1}{(\tau _{11, 2}^{(i)}, \tau _{11, 3}^{(i)})\cup(\tau_{11, 2}^{(m)}, +\infty)} $, the equilibrium $ E_1 $ is unstable.

    2) Equilibrium $ E_k\; (k = 2, 3) $

    (a) If $ h(z) $ has no positive root, the equilibrium $ E_k $ is locally asymptotically stable when $ \tau_1 > 0 $.

    (b) If $ h(z) $ only has one positive root $ z_1 $, the system (2.1) undergoes Hopf bifurcation at $ E_k $ when $ \tau_1 = \tau_{1k, 1}^{(j)} $ and $ h'(z_1) > 0 $. We get $ \forall\; 0 < \tau_1 < \tau_{1k, 1}^{(0)} $, the equilibrium $ E_k $ is asymptotically stable, and when $ \tau_1 > \tau_{1k, 1}^{(0)} $, the equilibrium $ E_k $ is unstable.

    (c) If $ h(z) $ has two positive roots $ z_1 $, $ z_2 $, the system (2.1) undergoes Hopf bifurcation at $ E_k $ when $ \tau_1 = \tau_{1k, 1}^{(j)} $ and $ \tau_1 = \tau_{1k, 2}^{(j)} $. Assuming $ z_2 < z_1 $, we can get $ h'(z_1) > 0, h'(z_2) < 0 $. Thus, assuming $ \tau_{1k, 1}^{(0)} < \tau_{1k, 2}^{(0)} $, there exists $ m $, which makes $ 0 < \tau_{1k, 1}^{(0)} < \tau_{1k, 2}^{(0)} < \tau_{1k, 1}^{(1)} < \tau_{1k, 2}^{(1)} < \dots < \tau_{1k, 1}^{(m)} < \tau_{1k, 1}^{(m+1)} $. When $ \tau_1\in[0, \tau _{1k, 1}^{(0)}) \cup \bigcup_{i = 1}^m{(\tau _{1k, 2}^{\left(i-1 \right)}, \tau _{1k, 1}^{\left(i \right)})} $, the equilibrium $ E_k $ is locally asymptotically stable. When $ \tau_1\in\bigcup_{i = 0}^{m-1}{(\tau_{1k, 1}^{(i)}, \tau_{1k, 2}^{(i)})}\cup (\tau_{1k, 1}^{(m)}, +\infty) $, the equilibrium $ E_k $ is unstable.

    (d) If $ h(z) $ has three positive roots $ z_1 $, $ z_2 $, $ z_3 $, the system (2.1) undergoes Hopf bifurcation at $ E_k $ when $ \tau_1 = \tau_{1k, l}^{(j)} (l = 1, 2, 3) $. We assume $ z_3 < z_2 < z_1 $, so we have $ h'(z_1) > 0, h'(z_2) < 0, h'(z_3) > 0 $. Similar to the analysis of (c), the equilibrium $ E_k $ switches between stability and instability with increasing $ \tau_1 $. Finally, the equilibrium $ E_k $ is unstable.

    (e) If $ h(z) $ has four positive roots $ z_1 $, $ z_2 $, $ z_3 $ and $ z_4 $, the system (2.1) undergoes Hopf bifurcation at $ E_k $ when $ \tau = \tau_{1k, l}^{(j)} (l = 1, 2, 3, 4) $. Assuming that $ z_4 < z_3 < z_2 < z_1 $, we can obtain $ h'(z_1) > 0, h'(z_2) < 0, h'(z_3) > 0, h'(z_4) < 0 $. Similar to the analysis of (c), the equilibrium $ E_k $ switches between stability and instability with increasing $ \tau_1 $. Ultimately, the equilibrium $ E_k $ is unstable.

    In this subsection, we suppose $ \tau_1 = \tau_1^*\in(0, \tau_{1k, 0}), \tau_{1k, 0} = \underset{l, j}{\min}\; \tau_{1k, l}^{(j)}\; (k = 1, 2, 3;l = 1, 2, 3, 4;j = 1, 2, \cdots) $ which means that the antibody failure time is $ \tau_1^* $. Then, we choose $ \tau_2 $ as the bifurcation parameter, which is the booster vaccination time. Similar to the situation of $ \tau_1 > 0, \tau_2 = 0 $, we have the following conclusions. The calculation details are presented in Appendix A.

    Theorem 3.3. For system (2.1) with $ \tau_2 > 0, \tau_1 > 0 $, when the stability conditions of Theorems 3.1 and 2 hold, the following conclusions hold.

    1) For equilibrium $ E_1 $, under (H7) and (H8), the equilibrium $ E_1 $ of system (2.1) is locally asymptotically stable when $ \tau_2\in[0, \tau_{21, 0}) $ for the chosen $ \tau_1^* $.

    2) For equilibrium $ E_k\; (k = 2, 3) $, under (H9) and (H10), the equilibrium $ E_2 $ or $ E_3 $ of system (2.1) is locally asymptotically stable when $ \tau_2\in[0, \tau_{2k, 0}) $ for the chosen $ \tau_1^* $.

    In this section, we derive the normal form of Hopf bifurcation about system (2.1) by using multiple time scales method. When $ \tau_1 = \tau_1^*, \tau_2 = \tau_{2k, 0} \; (k = 1, 2, 3) $, we assume that the characteristic equation (3.3) or (3.4) has a pair of pure imaginary roots $ \lambda = \pm\mathrm{i}\omega^* $ at which system (2.1) undergoes Hopf bifurcation at equilibrium $ E_k = (S_k^*, V_k^*, I_k^*, R_k^*)\; (k = 1, 2, 3) $. For convenience, we let $ S^*\triangleq S_k^* $, $ V^*\triangleq V_k^* $, $ I^*\triangleq I_k^* $, $ R^*\triangleq R_k^* $. Since we are more concerned about the effect of booster vaccination time on the epidemic, we select $ \tau_2 $ as the bifurcation parameter. We can obtain the normal form as follows. The calculation details can be referred to Appendix B.

    $ ˙G=MτεG+HG2ˉG, $

    where $ G $ is the coordinate projected on the center manifold and $ M, H $ are given in Eqs (B9) and (B14).

    Let $ G = r\text{e}^{\text{i}\theta} $ and substitute it into Eq (B15), we can obtain the normal form of Hopf bifurcation in polar coordinates:

    $ {˙r=Re(M)τεr+Re(H)r3,˙θ=Im(M)τε+Im(H)r2. $ (4.1)

    Then, we have the theorem as follows.

    Theorem 4.1. If $ \frac{\mathrm{Re}\left(M \right) \tau _{\varepsilon}}{\mathrm{Re}\left(H \right) \tau _c} < 0 $ holds, the system (2.1) exists periodic solution around the equilibrium.

    1) If $ \mathrm{Re}(M){{\tau }_{\varepsilon }} < 0 $, the periodic solution of system (2.1) is unstable.

    2) If $ \mathrm{Re}(M){{\tau }_{\varepsilon }} > 0 $, the periodic solution of system (2.1) is stable.

    Proof. Although system (4.1) is two-dimensional, we find that the first equation is independent of $ \theta $, which represents the argument component in polar coordinates, so we only need to consider the first equation of system (4.1). If $ \frac{\mathrm{Re}\left(M \right) \tau _{\varepsilon}}{\mathrm{Re}\left(H \right) \tau _c} < 0 $, $ \dot{r} = \text{Re}\left(M \right) \tau _{\varepsilon}r+\text{Re}\left(H \right) r^3 $ has a nontrivial fixed point $ r^* = \sqrt{-\frac{\mathrm{Re}\left(M \right) \tau _{\varepsilon}}{\mathrm{Re}\left(H \right) \tau _c}} $. At this equilibrium $ r^* $, the eigenvalue is $ \lambda = -2\mathrm{Re}(M)\tau_\epsilon $. According to Lyapunov indirect method, if $ \mathrm{Re}(M)\tau_\epsilon < 0 $, the equilibrium $ r^* $ is unstable. On the contrary, if $ \mathrm{Re}(M)\tau_\epsilon > 0 $, the equilibrium $ r^* $ is stable. Therefore, the periodic solution of system (2.1) is unstable if $ \mathrm{Re}(M){{\tau }_{\varepsilon }} < 0 $ and is stable if $ \mathrm{Re}(M){{\tau }_{\varepsilon }} > 0 $.

    In this section, we first analyze the official data for rationalization and select the parameter values that are realistic. Second, the accuracy of the theoretical analysis is verified by simulating under these parameters. According to the results, we select the most suitable booster vaccination time and provide some reasonable suggestions to control the epidemic.

    In this section, we use statistical methods to analyze the values of parameters, and on this basis, we select two groups of parameters most consistent with reality.

    1) Cure rate: $ u $

    We obtain the average annual cure rates in different countries from the World Health Organization website (https://www.who.int/). By eliminating the missing values and outliers, we get the cure rates for 65 countries and plot the histogram in Figure 2.

    Figure 2.  Histogram of cure rates $u$ in 65 countries.

    We can clearly find that the cure rates of these countries are primarily concentrated in two categories: one is countries with cure rates greater than 0.8, which are mainly concentrated at approximately 0.95. We call it the first echelon. The other group is countries with cure rates between 0.6 and 0.8. We call it the second echelon, typically represented by the United States, with a cure rate of 0.66. Clearly, the higher the cure rate is in a country, the easier it is to control the outbreak. However, since we want to control the epidemic as effectively as possible, we focus on countries in the second echelon, choosing the cure rate $ u = 0.7\in(0.6, 0.8) $.

    2) Vaccination rates: $ \alpha $

    By visiting the official website of Our World in Data (https://ourworldindata.org/covid-vaccinations), we get different vaccination amounts from different countries around the world. Following data collation, the global intercontinental distribution of vaccination was obtained, as shown in Figure 3. According to the law of large numbers, the more experiments there are, the closer the frequency is to the probability, that is, the more accurately and reliably we can approximate the probability with frequency. Based on this, after obtaining the vaccination rates of each continent, we use them as the weights of the vaccination rates of each continent so that the coverage rate is more reliable.

    Figure 3.  Intercontinental distribution of vaccination worldwide (Unit: ten thousands).

    Based on the intercontinental distribution of global vaccination, we can easily see that although 4.4 billion doses have been vaccinated worldwide, the distribution is very uneven. At the intercontinental level, Asia accounts for 65% of global vaccinations with more than 2.8 billion doses, followed by Europe with 15%, North America with 12%, and South America with 6%, while Africa, with nearly 20% of global population, accounts for less than 2% of global vaccinations. Therefore, we need to find a vaccination rate that reflects the global vaccination situation. Our approach is a weighted average of cure rates across continents. Europe and North America are the continents with the highest vaccination rates per 100 people, and Europe has recently surpassed North America with nearly 90 doses per 100 people. South America is at the same level as Asia, with an average of 55 doses per 100 people. Oceania ranks fifth with 32 doses per 100 people, and Africa has only five doses per 100 people, which is much lower than others. Based on the above information, we weight the vaccination rates of each continent and finally obtain the average annual vaccination rate $ \alpha = 0.63 $.

    3) Natural mortality rate: $ d $

    Since 2019, the outbreak of the COVID-19 epidemic has caused thousands of deaths and a new increase in mortality. However, for the natural mortality rate of the population, we no longer have access to data from official statistics due to the interference of deaths from infection. Here, our approach is to count the official population mortality rates (https://databank.worldbank.org/) before the outbreak and obtain an estimate of the natural mortality rate after the COVID-19 outbreak through fitting. With the 95% confidence interval, we select the following fit results based on the merits of the best fit, which are shown in Table 2. We find that the error sum of squares SSE is small, indicating that the error is small and the fit is great.

    Table 2.  Results of fitting.
    Country Function SSE
    Australia $ y_1=0.006549 + 4.366\mathrm{e}^{-5}\cos(0.3491t) + 8.852\mathrm{e}^{-5}\sin(0.3491t) $ $ 1.724\mathrm{e}^{-8} $
    Canada $ y_2= 0.007141 + 7.21\mathrm{e}^{-5}\cos(0.6981t) + 2.463\mathrm{e}^{-6}\sin(0.6981t) $ $ 1.09\mathrm{e}^{-7} $
    Switzerland $ y_3= 0.008143 -0.0002755\cos(0.2715t) + 5.128\mathrm{e}^{-5}\sin(0.2715t) $ $ 4.291\mathrm{e}^{-7} $
    China $ y_4= 0.006841 -0.0002921\cos(0.2853t) +0.0002776\sin(0.2853t) $ $ 1.352\mathrm{e}^{-7} $

     | Show Table
    DownLoad: CSV

    Based on the fitting results, we predict the natural mortality rates of the population after the outbreak, and the results for the next five years are shown in Figure 4.

    Figure 4.  Dates and trends of natural mortality rate $ c $ of Australia, Canada, Switzerland and China.

    We find that the natural population mortality rate is distributed between 0.0065 and 0.0085 for the four representative countries in 2019–2025, so we choose $ d = 0.00714 $ as the natural population mortality rate.

    4) COVID-19 mortality rate: $ c $

    We also obtain the average annual COVID-19 mortality data of different countries from the official website of the World Health Organization (https://www.who.int/). After data processing, we select the rates from 160 countries and draw the corresponding histogram with the beta distribution density curve, as shown in Figure 5. We find that the frequency of the distribution of COVID-19 mortality rates peaks at approximately $ c = 0.01 $ and $ c = 0.02 $. To simulate the real situation because of the different medical capacities in different regions and the different characteristics of mutated virus at different time, we finally choose $ c = 0.01 $ and $ c = 0.02 $ as the COVID-19 mortality rates.

    Figure 5.  Probability distribution of COVID-19 mortality rate $ c $ in 160 countries.

    5) Booster vaccination rates: $ \rho $

    For the booster vaccination every year on average, we use the analytical estimation method. Since the booster vaccination is for people who already have antibodies, that is, the booster vaccination rate is relative to $ R $ and a large portion of $ R $ comes from vaccination, the initial vaccination has already screened out those who could not be vaccinated due to individual reasons. Therefore, $ R $ has a stronger desire to be vaccinated than $ S $, so we believe that the rate of booster vaccination is higher than that of initial vaccination. With the above considerations, we choose the average annual booster vaccination rate $ \rho = 0.7 $.

    6) Infection rate: $ \beta $, $ \gamma $

    Since the global COVID-19 outbreak, countries have paid great attention to the infection and have adopted policies to stop the spread of the epidemic. One of the most important parts is wearing masks in public, which has led to a significant reduction in the infection rate. However, since the epidemic differs in its ability to spread in different countries, the infection rate also varies. Here, we take $ \beta $ to be 0.0009 and 0.00015 respectively. In addition, the immune barrier of $ V $ has not been fully formed although these individuals are vaccinated. Therefore, there is also a specific transmission rate for them, but this rate is significantly smaller than that of susceptible people. Therefore, we take $ \gamma = 0.0002 $ and $ \gamma = 0.00008 $.

    Based on the above discussion and estimating the reality of the situation, we take the following two groups of parameters:

    (Ⅰ): $ \Lambda = 10, $ $ \mu = 0.7, $ $ \alpha = 0.63, $ $ \beta = 0.0009, $ $ \gamma = 0.0002, $ $ d = 0.00714, $ $ \rho = 0.7, $ $ u = 0.7, $ $ c = 0.02, $ $ m = 0.9; $

    (Ⅱ): $ \Lambda = 80, $ $ \mu = 0.8, $ $ \alpha = 0.63, $ $ \beta = 0.00015, $ $ \gamma = 0.00008, $ $ d = 0.00714, $ $ \rho = 0.7, $ $ u = 0.7, $ $ c = 0.01, $ $ m = 0.9. $

    For the group of parameters (Ⅰ): $ \Lambda = 10, \; \mu = 0.7, \; \alpha = 0.63, \; \beta = 0.0009, \; \gamma = 0.0002, \; d = 0.00714 $, $ \rho = 0.7, \; u = 0.7, \; c = 0.02, \; m = 0.9 $, we calculate the disease-free equilibrium $ E_1 = [431, 591, 0, 378] $. We find that the assumptions (H1) and (H2) are not valid, so equilibria $ E_2 $ and $ E_3 $ do not exist. Under the group of parameters, (H3) holds, so $ E_1 = [431, 591, 0, 378] $ is locally asymptotically stable when $ \tau_1 = \tau_2 = 0 $. When $ \tau_1 > 0, \tau_2 = 0 $, we can obtain (H4) holds and $ h(z) $ only has one positive root, $ \omega_{11, 1} = 0.2073 $, $ \sin(\omega_{11, 1}\tau_{11, 1}^{(0)}) = 0.3039 $, $ \cos(\omega_{11, 1}\tau_{11, 1}^{(0)}) = -0.9527 $, $ \tau_{11, 0} = \tau_{11, 1}^{(0)} = 13.66 $. We choose $ \tau_1 = \tau_1^* = 4.5 $ and substitute it into Eqs (9.1)$ - $(9.4), we have $ \omega_{21, 1} = 0.61 $, $ \sin(\omega_{21, 1}\tau_{21, 1}^{(0)}) = 0.5821 $, $ \cos(\omega_{11, 1}\tau_{11, 1}^{(0)}) = -0.8131 $, $ \tau_{21, 0} = \tau_{21, 1}^{(0)} = 1.57 $. If ${\bf(H7)}$ and ${\bf(H8)}$ hold, the equilibrium $ E_1 $ is locally asymptotically stable when $ \tau_2\in\left[\left.0, \tau_{21, 0}\right)\right. $.

    When $ \tau_1 = \tau_2 = 0 $, the equilibrium $ E_1 $ is locally asymptotically stable according to Theorem 3.1. This means that the antibody expires instantaneously, and people maintain the effectiveness of the antibodies through uninterrupted booster shots. First, we choose $ [S(t), V(t), I(t), R(t)]^\mathrm{T} = [500, 600, 100, 500]^\mathrm{T} $ for $ t\in[-\tau_2, 0] $ as the initial function, which represents the initial state of the epidemic and then draw a time history graph in Figure 6(b). We can see that in this case, the epidemic basically reaches stability in 5 years and is eliminated after 10 years. This result is very optimistic. Next, we select three sets of initial functions 1) $ [S(t), V(t), I(t), R(t)]^\mathrm{T} = [500, 600, 100, 500]^\mathrm{T} $, 2) $ [S(t), V(t), I(t), R(t)]^\mathrm{T} = [5000, 6000, 1000, 5000]^\mathrm{T} $, 3) $ [S(t), V(t), I(t), R(t)]^\mathrm{T} = [10, 000, 20, 000, 10, 000, 10, 000]^\mathrm{T} $ for $ t\in[-\tau_2, 0] $, which represent different initial states of the epidemic. We find that when the initial functions are far from equilibrium, the fluctuations are large during the initial years, and the time taken to reach stability is long, as shown in Figure 6(a). However, the epidemic eventually reaches stability regardless of the initial functions, so we conjecture that the equilibrium is not only locally asymptotically stable, but may also be globally stable. This finding suggests that under this set of parameters, even if the epidemic is severe in a given region, the epidemic will eventually reach a manageable level and will be eliminated through uninterrupted booster vaccination. However, as far as reality is concerned, the vaccines we have developed cannot be immediate and uninterrupted, so this situation does not exist.

    Figure 6.  When $ \tau_1 = \tau_2 = 0 $, equilibrium $ E_1 $ of system (2.1) is locally asymptotically stable.

    When $ \tau_1 = 4.5\in(0, \tau_{11, 0}) = (0, 13.66) $, $ \tau_2 = 1\in(0, \tau_{21, 0}) = (0, 1.57) $, the equilibrium $ E_1 $ is locally asymptotically stable according to Theorem 3.3. $ \tau_1 = 4.5 $ means the vaccine will fail 4.5 years after individuals acquires antibodies, and $ \tau_2 = 1 $ means that people begin the booster vaccination after 1 year to cope with the decline in antibodies.

    First, we still choose the initial function $ [S(t), V(t), I(t), R(t)]^\mathrm{T} = [500,600,100,500]^\mathrm{T} $ for $ t\in[-\tau_2, 0] $ and show the variation in the number of people in different compartments over time in Figure 7(b). This figure suggests that receiving the booster vaccinations after one year will lead to a stable epidemic in 10 years and the society will persist in this state consistently. In addition, the number of infected people would eventually approach 0, meaning that the epidemic can eventually be eliminated. That is the most desirable outcome. After that, we choose 1) $ [S(t), V(t), I(t), R(t)]^\mathrm{T} = [500,600,100,500]^\mathrm{T} $, 2) $ [S(t), V(t), I(t), R(t)]^\mathrm{T} = [5000,6000,1000,5000]^\mathrm{T} $, 3) $ [S(t), V(t), I(t), R(t)]^\mathrm{T} = [10,000,20,000,10,000,10,000]^\mathrm{T} $ for $ t\in[-\tau_2, 0] $ as three sets of initial functions that are same as $ \tau_1 = \tau_2 = 0 $, and the results are shown in Figure 7(a). Although the amplitude of fluctuation and the time taken to reach stability increase gradually as the initial function moves away from equilibrium, it will eventually stabilize. Thus, if we carry out booster vaccination after one year, the outbreak will be controlled. However, the time needed to contain the outbreak depends on the current situation. The more severe the epidemic is, the more difficult it will be to control. According to the above results, we guess that the equilibrium $ E_1 $ is not only locally asymptotically stable but also globally stable under this group of parameters.

    Figure 7.  When $ \tau_1 = 4.5 $, $ \tau_2 = 1 $, equilibrium $ E_1 $ of system (2.1) is locally asymptotically stable.

    For the group of parameters (Ⅱ): $ \Lambda = 90, \; \mu = 0.8, \; \alpha = 0.63, \; \beta = 0.00015, \; \gamma = 0.00008, \; d = 0.00714, $ $ \rho = 0.7, \; u = 0.7, \; c = 0.02, \; m = 0.9, $ we find that (H1) holds, so equilibrium $ E_2 $ makes sense and is [2817, 3683, 875, 2605]. (H2) and (H3) do not hold, equilibrium $ E_1 $ is unstable and equilibrium $ E_3 $ is meaningless. Substituting this group of parameters into Eq (3.7), (H4) is satisfied, so equilibrium $ E_2 $ is locally asymptotically stable when $ \tau_1 = \tau_2 = 0 $. When $ \tau_1 > 0, \tau_2 = 0 $, we can get that $ h(z) $ has three roots and $ z_3 = 0.0004 < z_2 = 0.0206 < z_1 = 0.1458 $, $ \omega_{12, 1} = 0.38 $ $ \sin(\omega_{12, 1}\tau_{12, 1}^{(0)}) = 0.4150 $, $ \cos(\omega_{12, 1}\tau_{12, 1}^{(0)}) = -0.9098 $, $ \tau_{12, 0} = \tau_{12, 1}^{(0)} = 7.11 $. Choosing $ \tau_1 = 4 $, we obtain that (H9) and (H10) hold and $ \omega_{22, 1} = 0.66 $, $ \sin(\omega_{22, 1}\tau_{22, 1}^{(0)}) = 0.5738 $, $ \cos(\omega_{22, 1}\tau_{22, 1}^{(0)}) = -0.8190 $, $ \tau_{22, 0} = \tau_{22, 1}^{(0)} = 0.92 $. Substitute the parameters (Ⅱ) into Eqs. (B9) and (B14), we have $ \mathrm{Re}(M) > 0, \mathrm{Re}(H) < 0 $. We can deduce $ \tau_\varepsilon > 0 $, $ \mathrm{Re}(M_k)\tau_\varepsilon > 0 $, the periodic solution is stable based on Theorem 4.1.

    When $ \tau_1 = \tau_2 = 0 $, we choose 1) $ [S(t), V(t), I(t), R(t)]^\mathrm{T} $ = $ [3000,4000,1000,3000]^\mathrm{T} $, 2) $ [S(t), V(t), I(t), R(t)]^\mathrm{T} $ = $ [5000,8000,2000,5000]^\mathrm{T} $, 3) $ [S(t), V(t), I(t), R(t)]^\mathrm{T} $ = $ [10,000,20,000,5000,10,000]^\mathrm{T} $ for $ t\in[-\tau_2, 0] $ as the initial functions separately. We draw a plot of the epidemic over time for the three sets of initial functions and a separate time history graph with initial function 1) in Figure 8.

    Figure 8.  When $ \tau_1 = \tau_2 = 0 $, equilibrium $ E_2 $ of system (2.1) is locally asymptotically stable.

    It is easy to observe from Figure 8(b) that the epidemic reaches stability rapidly within 8 years. However, it cannot be eradicated and coexists with humans for a long time. In Figure 8(a), the initial function 3) takes longer to reach the steady state, and the stabilization effect is not as good as that of 1) and 2), but it eventually reaches stability. Therefore, we conjecture that the equilibrium $ E_2 $ is globally stable. Thus, we can get the following results: it is effective to control the epidemic through uninterrupted booster vaccination, but the time needed to control the epidemic mainly depends on the current state. The more severe the epidemic, the more difficult it will be. This finding shows that booster vaccination is of great significance, but considering that the vaccines developed at present are not instantaneous failures and vaccines cannot be administered in a short period of time, $ \tau_1 = \tau_2 = 0 $ is basically impossible.

    When $ \tau_1 = 4\in(0, \tau_{12, 0}) = (0, 7.11), \tau_2 = 0.8\in(0, \tau_{22, 0}) = (0, 0.92) $, the antibody expires after 4 years and we re-inject booster shots after 0.8 years. We also choose the three sets of initial functions: 1) $ [S(t), V(t), I(t), R(t)]^\mathrm{T} = [500,600,100,500]^\mathrm{T} $, 2) $ [S(t), V(t), I(t), R(t)]^\mathrm{T} = [5000,6000,1000,5000]^\mathrm{T} $, 3) $ [S(t), V(t), I(t), R(t)]^\mathrm{T} = [10,000,20,000,10,000,10,000]^\mathrm{T} $ for $ t\in[-\tau_2, 0] $ and make figures respectively similar to the above situation. The results are shown in Figure 9.

    Figure 9.  When $ \tau_1 = 4, \tau_2 = 0.8 $, equilibrium $ E_2 $ of system (2.1) is locally asymptotically stable.

    With antibody failure after 4 years and booster vaccination after 0.8 years, the system quickly becomes stable when the initial function is chosen near equilibrium according to (b) in Figure 9. From the figure, we also find that although the epidemic cannot be eliminated, it can always stay in a stable state. In Figure 9(a), as the initial function moves away from equilibrium, the epidemic worsens. It fluctuates dramatically in a short period of time and will take decades to stabilize. But no matter what initial function is chosen, the epidemic will eventually be controlled. Therefore, we hypothesize that this stability is global. Even if the current epidemic is not optimistic in some countries, the trend will be better with the booster vaccination in time.

    When $ \tau_1 = 4, \; \tau_2 = 1.11 > \tau_{22, 0} = 0.92 $, which means that the antibody will expire after 3.9 years and booster vaccination will be carried out after 1.11 years, the equilibrium $ E_2 $ of system (2.1) is unstable. According to Theorem 4.1, system (2.1) has a forward Hopf bifurcation near $ \tau_{22, 1}^{(0)} = 0.92 $ and the periodic solution is stable. We choose $ [S(t), V(t), I(t), R(t)]^\mathrm{T} = [3000,4000,1000,3000]^\mathrm{T} $ for $ t\in[-\tau_2, 0] $ as the initial function, and the epidemic situation is shown in Figure 10.

    Figure 10.  When $ \tau_1 = 4, \tau_2 = 1.1 $, equilibrium $ E_2 $ of system (2.1) is unstable.

    In this case, the epidemic will continue to fluctuate, which is not promising for controlling the epidemic.

    Based on the above simulations, we know that the only way to avoid out-of-control epidemic is to carry out booster vaccinations within the critical time. However, when is the best time to vaccinate? Next, we will analyze the optimal time of booster vaccination. First, we take the first set of parameters as the study object given $ \tau_2 = 0.1, 0.4, 0.7, 1.0, 1.3 $ respectively, and then we draw the images of the infected people ($ I $) with time as shown in Figure 11.

    Figure 11.  Variation trends of $ I $ with different $ \tau_2 $.

    In Figure 11(a), we find that the time to eventual stabilization is essentially the same regardless of the booster vaccination time we select, but when we focus on the case near the peak in Figure 11(b), we find that the smaller the booster vaccination time is, the smaller the number of infections will be. This finding reveals that the earlier the booster vaccination is performed, the fewer people will be infected with COVID-19. However, considering human and material resources and financial resources, it is not an easy task to complete community-wide vaccination, so we need to implement booster vaccinationss as soon as possible according to the local reality and within the scope of ensuring the normal operation of the community.

    From a medical point of view, it is generally believed that the immune effect weakens after six months of vaccination, so 6 months after the completion of the initial vaccine is the critical period for booster vaccination [51,52]. For the optimal strategy of vaccine distribution, researchers solve the optimal scheduling of the mixed vaccination strategy by using different optimization methods [53,54,55]. From the point of view of completing one round of vaccination, we need to allow sufficient time for booster vaccination [56]. Taking these considerations into account, we can obtain the optimal booster vaccination time.

    Based on the above simulations and analysis, we have the following results.

    1) When the time of antibody expiration is determined, the epidemic will become increasingly not controllable as the time of booster vaccination increases, and when the time of booster vaccination exceeds the critical value, the system (2.1) will be unstable, and the epidemic will be out of control. However, the time of booster vaccination is not as small as possible. In reality, we need to consider the comprehensive capacity of the country to ensure that booster vaccination time is minimized within the capacity.

    2) After simulating the changes in epidemic development under three different sets of initial functions, we find that the epidemic will eventually be controlled, so we assume that the stability of equilibrium $ E_k\; (k = 1, 2, 3) $ is not only locally asymptotically stable, but may also be globally stable. That is, even if the current epidemic is not optimistic in a certain region, the possible global stability ensures that the epidemic can eventually reach effective control after timely booster vaccination.

    3) When the equilibrium $ E_1 $ satisfies the condition of locally asymptotical stability, the epidemic can eventually be eliminated. If the selected parameters meet the existence and stability of equilibrium $ E_2 $ or $ E_3 $, it can reach a stable state. Although the epidemic cannot be eliminated, to some extent, it can be considered to be effectively controlled.

    4) Our model is flexible. The current state of the epidemic varies from country to country, and we can simulate the epidemic of different countries by changing the parameters. For example, when a mutant strain is found in a country and has been confirmed to spread widely around the country, we can simulate this situation by increasing the infection rate $ \alpha $ and determine the population size by natural population growth $ \Lambda $. If the current vaccine is found to be less resistant to the mutant strain, we can simultaneously increase $ \gamma $, $ \mu $ and determine population size by natural population growth $ \Lambda $. After determining the parameters, we can calculate the critical booster vaccination time. When the booster vaccination time is less than the critical value, the system can reach stability, that is, the current outbreak can be controlled. Based on the above analysis, we can develop a policy applicable to this region that allows to control the epidemic as much as possible within the capacity of the region.

    Based on the above conclusions, we give the following recommendations.

    1) The time of the booster vaccination should be kept within the threshold so that the epidemic can be controlled eventually. However, the overall strength of the area, such as economic development, medical level, etc., also needs to be assessed. For example, in the second group of parameters, $ \tau_{22, 0} $ = 0.92. If this is a well-developed country, we would suggest that the booster shot is needed after 0.5–0.92 years. As the epidemic continues to evolve, we are concerned about more than just one booster shot. Our model can apply to more than one booster shot, gradually developing a presence similar to that of the annual influenza vaccination. So we need booster vaccinations every 0.5–0.92 years to ensure the continued effectiveness of the antibodies.

    2) Even if the development of the epidemic in a given country is not promising, it should not be taken lightly. Although providing booster vaccination at an appropriate time cannot guarantee that the epidemic will reach control within a short period of time, it has a positive effect on the development of the epidemic. Because the possible global stability ensures that the epidemic will eventually be effectively controlled.

    3) It is necessary to detect changes in the epidemic continuously, such as the emergence of new mutant strains and changes in the time of antibody expiration. We can regain a critical booster vaccination time through calculation and then change the booster vaccination strategy in time to cope with the variable COVID-19.

    In this paper, we have constructed an $ SVIR $ model with two time delays to study the suitable time for booster vaccination. We have considered the existence and stability of equilibria in our model. Then, we have also analyzed the existence and dynamic properties of Hopf bifurcation. In Section 5, we chose two groups of parameters by analyzing the data from official websites. Then, we have carried out numerical simulations to verify the analytical results and provide some reasonable suggestions to control the development of the epidemic.

    According to the results of numerical simulations, we found that the epidemic would become increasingly difficult to control as the booster vaccination time $ \tau_2 $ increased. Considering the factors of all aspects, we suggest that it is suitable to vaccine with the booster shots approximately 6–11 months later for the group of parameters (Ⅱ). We also found that the equilibrium of system (2.1) might be globally stable when $ \tau_2 $ is within the threshold value. That means, the booster vaccination is effective for outbreak control, so it is crucial to provide booster vaccination in time even if the current epidemic is serious. In addition, considering that the selection of parameters might change due to virus mutation or other factors, we need to continuously monitor the status of the epidemic in real time and adjust the strategy to ensure that the epidemic can be effectively controlled.

    This study was funded by Fundamental Research Funds for the Central Universities of China (Grant No.2572022DJ06).

    All the authors declare no conflicts of interest in this paper.

    For equilibrium $ E_1 $, we let $ \lambda = \mathrm{i}\omega\; (\omega > 0) $ be a pure imaginary root of Eq (3.3) and substitute it into Eq (3.3). Then separating the real part and imaginary part, we have:

    $ {ω2+b2+μωsin(ωτ1)+μc1cos(ωτ1)=ρc2cos(ωτ2)ρωsin(ωτ2),b1ω+μωcos(ωτ1)μc1sin(ωτ1)=ρc2sin(ωτ2)ρωcos(ωτ2), $ (A1)

    where $ b_1 = m+\alpha+2d, b_2 = (m+d)(\alpha+d), c_1 = m+\alpha+d, c_2 = \alpha+d. $ Eq (A1) leads to:

    $ {sin(ωτ2)=ω[ω2+b2+μωsin(ωτ1)+μc1cos(ωτ1)]ρ(c22+ω2)+c2[b1ω+μωcos(ωτ1)μc1sin(ωτ1)]ρ(c22+ω2),cos(ωτ2)=c2[ω2+b2+μωsin(ωτ1)+μc1cos(ωτ1)]ρ(c22+ω2)+ω[b1ω+μωcos(ωτ1)μc1sin(ωτ1)]ρ(c22+ω2). $ (A2)

    Adding the square of the two equations in Eq (A2), we get:

    $ F1(ω)=ω4+(2b2+b21+μ2ρ2)ω2+b22μ2c21ρ2c22+2sin(ωτ1)[μω3+b2μωb1μωc1]+2cos(ωτ1)[μc1ω2+μb2c1+b1μω2], $ (A3)

    where $ b_1, b_2, c_1, c_2 $ are given in Eq (A1). We suppose:

    $ {\bf (H7)} \; \; \mu^2c_1^2-\rho^2c_2^2+2\mu b_2c_1 < 0. $

    Under (H7), we have $ F_1(0) = \mu^2c_1^2-\rho^2c_2^2+2\mu b_2c_1 < 0 $, $ F_1(\infty) > 0 $. Therefore, there is at least one positive root of $ F_1(\omega) = 0 $. We assume there are $ l $ positive roots of $ F_1(\omega) = 0 $, hence, there is the expression of $ \tau_{21, l}^{(j)}, j = 0, 1, 2, \cdots $.

    $ τ(j)21,l={1ω21,l[arccos(P21,l)+2jπ],           Q21,l0,1ω21,l[2πarccos(P21,l)+2jπ],     Q21,l<0, $ (A4)

    where $ Q_{21, l} = {\left. {\sin(\omega\tau_2)} \right|_{\omega = {\omega _{21, l}}, \tau = \tau _{21, l}^{(j)}}}, {P_{21, l}} = \cos(\omega\tau _{2}) {\left. \right|_{\omega = {\omega _{21, l}}, \tau = \tau _{21, l}^{(j)}}}, $ $ \sin(\omega\tau_2) $ and $ \cos(\omega\tau _{2}) $ are given in Eq (A2).

    Let $ \tau_{21, 0} = \min\; \tau_{21, l}^{(j)} = \tau_{21, l}^{(0)}, j = 0, 1, 2, \cdots $. When $ \tau_2 = \tau_{21, 0} $, Equation (3.5) has a pair of purely imaginary roots $ \pm\mathrm{i}\omega_{21, 0} $. Assume:

    $ {\bf(H8)}\; \; \left.\mathrm{Re}\left(\frac{\mathrm{d}\lambda }{\mathrm{d}\tau_2 }\right)^{-1} \right|_{\tau = \tau _{21, 0}}\neq0. $

    Under ${\bf(H8)}$, the equilibrium $ E_1 $ is locally asymptotically stable when $ \tau_2\in\left[\left.0, \tau_{21, 0}\right)\right. $.

    For equilibrium $ E_k\; (k = 2, 3) $, similar to the above analysis, we let $ \lambda = \mathrm{i}\omega(\omega > 0) $ be a pure imaginary root of characteristic equation (3.4) and substitute it into Eq (3.4). Then, separating the real part and imaginary part, we have:

    $ {ω4ϱ2,kω2+ϱ4,k+μcos(ωτ1)(σ1,kω2+σ3,k)+μsin(ωτ1)(ω3+σ2,kω)=ρ(ς1,kω2+ς3,k)cos(ωτ2)ρ(ω3+ς2,kω)sin(ωτ2),ϱ1,kω3+ϱ3,kω+μcos(ωτ1)(ω3+σ2,kω)μsin(ωτ1)(σ1,kω2+σ3,k)=ρ(ς1,kω2+ς3,k)sin(ωτ2)ρ(ω3+ς2,kω)cos(ωτ2). $ (A5)

    Then, we can obtain:

    $ {sin(ωτ2)=(ω4ϱ2,kω2+ϱ4,k)(ω3+ς2,kω)+(ϱ1,kω3+ϱ3,kω)(ς1,kω2+ς3,k)ρ[(ς1,kω2+ς3,k)2+(ω3+ς2,kω)2]+μcos(ωτ1)[(σ1,kω2+σ3,k)(ω3+ς2,kω)+(ω3+σ2,kω)(ς1,kω2+ς3,k)]ρ[(ς1,kω2+ς3,k)2+(ω3+ς2,kω)2]+μsin(ωτ1)[(ω3+σ2,kω)(ω3+ς2,kω)(σ1,kω2+σ3,k)(ς1,kω2+ς3,k)]ρ[(ς1,kω2+ς3,k)2+(ω3+ς2,kω)2].cos(ωτ2)=(ω4ϱ2,kω2+ϱ4,k)(ς1,kω2+ς3,k)+(ϱ1,kω3+ϱ3,kω)(ω3+ς2,kω)ρ[(ς1,kω2+ς3,k)2+(ω3+ς2,kω)2]μcos(ωτ1)[(σ1,kω2+σ3,k)(ς1,kω2+ς3,k)+(ω3+σ2,kω)(ω3+ς2,kω)]ρ[(ς1,kω2+ς3,k)2+(ω3+ς2,kω)2]+μsin(ωτ1)[(ω3+σ2,kω)(ς1,kω2+ς3,k)+(σ1,kω2+σ3,k)(ω3+ς2,kω)]ρ[(ς1,kω2+ς3,k)2+(ω3+ς2,kω)2]. $ (A6)

    Adding the square of the two equations in (A6), we have:

    $ Fk(ω)=(ω4ϱ2,kω2+ϱ4,k)2+(ϱ1,kω3+ϱ3,kω)2+μ2(σ1,kω2+σ3,k)2+μ2(ω3+σ2,kω)2+2μcos(ωτ1)[(σ1,kω2+σ3,k)(ω4ϱ2,kω2+ϱ4,k)+(ω3+σ2,kω)(ϱ1,kω3+ϱ3,kω)]+2μsin(ωτ1)[(ω3+σ2,kω)(ω4ϱ2,kω2+ϱ4,k)(σ1,kω2+σ3,k)(ϱ1,kω3+ϱ3,kω)]ρ2(ς1,kω2+ς3,k)2ρ2(ω3+ς2,kω)2=0,k=2,3. $ (A7)

    Then we give the following assumption:

    $ {\bf (H9)} \; \; \varrho_{4, k}^2+\mu^2\sigma_{3, k}^2+2\mu \sigma_{3, k}\varrho_{4, k}-\rho^2\varsigma_{3, k}^2 < 0. $

    Under (H9), we can deduce $ F_k(0) < 0 $ and $ F_k(\infty) > 0 $. Thus, $ F_k(\omega) = 0 $ has at least one positive root. We assume that there are $ l $ positive roots of $ F_k(\omega) = 0 $ and denote them as $ \omega_{2k, l} $.

    $ τ(j)2k,l={1ω2k,l[arccos(P2k,l)+2jπ],           Q2k,l0,1ω2k,l[2πarccos(P2k,l)+2jπ],     Q2k,l<0, $ (A8)

    where $ Q_{2k, l} = {\sin(\omega\tau_2)}|_{\omega = {\omega _{2k, l}}, \tau = \tau _{2k, l}^{(j)}}, {P_{2k, l}} = \cos(\omega\tau _{2}) |_{\omega = {\omega _{2k, l}}, \tau = \tau _{2k, l}^{(j)}}, $ $ \sin(\omega\tau_2) $ and $ \cos(\omega\tau _{2}) $ are given in Eq (A6).

    Let $ \tau_{2k, 0} = \min\; \tau_{2k, l}^{(j)}, j = 0, 1, 2, \cdots, k = 2, 3 $. When $ \tau_2 = \tau_{2k, 0} $, Equation (3.5) has a pair of purely imaginary roots $ \pm\mathrm{i}\omega_{2k, 0} $. Assume

    $ {\bf(H10)}\ \left.\mathrm{Re}\left(\frac{\mathrm{d}\lambda}{\mathrm{d}\tau_2}\right)^{ - 1} \right|_{\tau = \tau _{2k, 0}}\neq0. $

    Under ${\bf(H10)}$, the equilibrium $ E_k\; (k = 2, 3) $ is locally asymptotically stable if $ \tau_2\in\left[\left.0, \tau_{2k, 0}\right)\right. $.

    The system (2.1) can be written as follows:

    $ X(t)=AX(t)+BX(tτ1)+CX(tτ2)+F[X(t),X(tτ1),X(tτ2)], $ (B1)

    where

    $ X(t)=(S,V,I,R)T,X(tτ1)=(S(tτ1),V(tτ1),I(tτ1),R(tτ1))T,X(tτ2)=(S(tτ2),V(tτ2),I(tτ2),R(tτ2))T,A=(αβId0βS0αmγIdγV0βIγIβS+rVduc00mud),B=(000μ00000000000μ),C=(0000000ρ0000000ρ),F=(FSFVFIFR)=(βSIγVIβSI+γVI0). $

    Let $ h = (h_{ 1}, h_{ 2}, h_{ 3}, h_{ 4})^T $ be the eigenvector of the linear operator corresponding to the eigenvalue $ \mathrm{i}\omega^* $, and let $ h^* = (h_{ 1}^*, h_{ 2}^*, h_{ 3}^*, h_{ 4}^*)^T $ be the normalized eigenvector of the adjoint operator of the linear operator corresponding to the eigenvalues $ -\mathrm{i}\omega^* $ satisfying the inner product $ < h^*, h > = 1 $. By a simple calculation, we can obtain:

    $ h1=1,h2=uh3+(iω+d+μeiωτ1+ρeiωτ2)h4m,h3=mμeiωτ1βImγI(iω+d+μeiωτ1+ρeiωτ2)(iω+α+βI+d)[(iωβSγV+d+u+c)m+uγI]μeiωτ1+γIβS(iω+d+μeiωτ1+ρeiωτ2),h4=(iω+α+βI+d)+βSh3μeiωτ1,˜h2=(iω+α+βI+d)˜h1βIα,hj=d˜hj,d=1h1¯˜h1+h2¯˜h2+h3¯˜h3+h4¯˜h4,˜h1=βI[um(iω+m+d+γI)+γV]α(iωβSγV+d+u+c+umγI)αβS+[um(iω+m+d+γI)+γV](iω+α+βI+d),˜h3=1,˜h4=(iω+m+d+γI)˜h2γI˜h3m. $ (B2)

    The solution of Eq (B1) $ X(t) $ can be written as:

    $ X(t)=X(T0,T1,T2,)=n=1εnXn(T0,T1,T2,), $ (B3)

    where $ X({{T}_{0}}, {{T}_{1}}, {{T}_{2}}, \cdots) =[S(T0,T1,T2,),V(T0,T1,T2,),I(T0,T1,T2,),R(T0,T1,T2,)]T$,$Xn(T0,T1,T2,) = {{[{{S}_{n}}({{T}_{0}}, {{T}_{1}}, {{T}_{2}}, \cdots), {{V}_{n}}({{T}_{0}}, {{T}_{1}}, {{T}_{2}}, \cdots), {{I}_{n}}({{T}_{0}}, {{T}_{1}}, {{T}_{2}}, \cdots), {{R}_{n}}({{T}_{0}}, {{T}_{1}}, {{T}_{2}}, \cdots)]}^\mathrm{T}} $, $ T_i = \epsilon^i t $, $ i = 0, 1, 2, \cdots $ and $ T_i $ is the scaling transform in the time direction.

    $ X'(t) $ can be written as:

    $ X(t)=dX(t)dt=εdX1dt+ε2dX2dt+ε3dX3dt+=ε(X1T0+εX1T1+ε2X1T2)+ε2(X2T0+εX2T1)+ε3X3T0+=εD0X1+ε2D1X1+ε3D2X1+ε2D0X2+ε3D1X2+ε3D0X3+, $ (B4)

    where $ D_i = \frac{\partial}{\partial T_i} (i = 1, 2, 3, \cdots) $ is differential operator.

    Since we are more concerned about the influence of booster vaccination time, we take $ \tau_2 $ as the bifurcation parameter. We let $ \tau_2 = \tau_c+\varepsilon\tau_\varepsilon $, where $ \tau_c $ is the critical time delay given in Eqs (A4) or (A8), $ \tau_\varepsilon $ is the disturbance parameter and $ \varepsilon $ is the dimensionless scale parameter. Using Taylor expansion of $ X(t-\tau_1^*) $ and $ X(t-\tau_2) $ respectively, we have:

    $ X(tτ1)=εX1,τ1+ε2(X2,τ1D1X1,τ1)+ε3(X3,τ1D1X2,τ1D2X1,τ1)+,X(tτ2)=εX1,τc+ε2X2,τc+ε3X3,τcε2τεD0X1,τcε3τεD0X2,τcε2τcD1X1,τcε3τεD1X1,τcε3τcD2X1,τcε3τcD1X2,τc+, $ (B5)

    where $ X_{j, \tau_1^*} = X_j(T_0-\tau_1^*, T_1, T_2, \cdots) $, $ X_{j, \tau_c} = X_j(T_0-\tau_c, T_1, T_2, \cdots) $, $ j = 1, 2, 3\cdots $. Then we substitute Eqs (B3)–(B5) into Eq (B1). For the $ \varepsilon $-order terms, we have:

    $ {D0S1+αS1+βSI1+βIS1+dS1μR1,τ1=0,D0V1αS1+mV1+γVI1+γIV1+dV1ρR1,τc=0,D0I1βSI1βIS1+dI1+uI1+cI1γVI1γIV1=0,D0R1mV1uI1+dR1+μR1,τ1+ρR1,τc=0. $ (B6)

    Since $ \pm \mathrm{i}{{\omega }^{*}} $ are the eigenvalues of the linear part of Eq (B1), the solution of Eq (B6) can be expressed in the following form:

    $ X1(T1,T2,T3,)=G(T1,T2,T3,)eiωT0h+ˉG(T1,T2,T3,)eiωT0ˉh, $ (B7)

    where $ h $ is given in Eq (B2).

    For the $ \varepsilon^2 $-order terms, we obtain:

    $ {D0S2+αS2+βSI2+βIS2+dS2μR2,τ1=D1S1βS1I1μD1R1,τ1,D0V2αS2+mV2+γVI2+γIV2+dV2ρR2,τc=D1V1γV1I1ρτεD0R1,τcρτcD0R1,τc,D0I2βSI2βIS2+dI2+uI2+cI2γVI2γIV2=D1I1+γV1I1+βS1I1,D0R2mV2uI2+dR2+μR2,τ1+ρR2,τc=D1R1+μD1R1,τ1+ρτεD0R1,τc+ρτcD0R1,τc. $ (B8)

    Then, we substitute Eq (B7) into the right side of Eq (B8) and mark the coefficient before $ \text{e}^{\text{i}\omega^* T_0} $ as vector $ m_1 $. In accordance with the solvability condition $ < h^*, m_1 > = 0 $, we can obtain the expression of $ \frac{\partial G}{\partial T_1} $:

    $ GT1=MτεG, $ (B9)

    where $ M = \frac{\rho \text{i}\omega^* \left(h_{4}\bar{h}_{4}^{*}-h_{4}\bar{h}_{2}^{*} \right)}{1-\mu \text{e}^{-\mathrm{i}\omega^*\tau _1^*}\left(h_{4}\bar{h}_{4}^{*}-h_{4}\bar{h}_{1}^{*} \right) -\rho \tau _c\text{e}^{-\mathrm{i}\omega^*\tau _c}\left(h_{4}\bar{h}_{4}^{*}-h_{4}\bar{h}_{2}^{*} \right)} $.

    We assume that the solution of Eq (B8) will take the following form:

    $ S2=g1e2iωτcT0G2+ˉg1e2iωτcT0ˉG2+l1GˉG,A2=g2e2iωτcT0G2+ˉg2e2iωτcT0ˉG2+l2GˉG,I2=g3e2iωτcT0G2+ˉg3e2iωτcT0ˉG2+l3GˉG,R2=g4e2iωτcT0G2+ˉg4e2iωτcT0ˉG2+l4GˉG. $ (B10)

    Substituting them into Eq (B8), we can solve the expression of $ g_{1} $, $ g_{2} $, $ g_{3} $, $ g_{4} $, $ l_{1} $, $ l_{2} $, $ l_{3} $, $ l_{4} $ from the following equations:

    $ [(2iω+ϵ1)0βSμe2iωτ1α(2iω+ϵ2)γVρe2iωτcβIγI(2iω+ϵ3)00mu(2iω+ϵ4)][g1g2g3g4]=[βh1h3γh2h3γh2h3+βh1h30], $ (B11)
    $ [ϵ10βSμαϵ2γVρβIγIϵ300mu(d+μ+ρ)][l1l2l3l4]=[βh1ˉh3βˉh1h3γˉh2h3γh2ˉh3γˉh2h3+βˉh1h3+γh2ˉh3+βh1ˉh30], $ (B12)

    where $ \epsilon_1 = \alpha +\beta I^*+d, \epsilon_2 = m+d+\gamma I^*, \epsilon_3 = -\beta S^*-\gamma V^*+d+u+c, \epsilon_4 = d+\mu \text{e}^{-\mathrm{i}\omega^*\tau _1^*}+\rho \text{e}^{-\mathrm{i}\omega^*\tau _2}. $ For $ \varepsilon^3 $-term, we have:

    $ {D0S3+αS3+βSI3+βIS3+dS3μR3,τ1=D2S1βS1I2μD1R2,τ1D1S2βS2I1μD2R1,τ1,D0V3αS3+mV3+γVI3+γIV3+dV3ρR3,τc=D2V1γV2I1ρτεD1R1,τcρτcD1R2,τcD1V2γV1I2ρτεD0R2,τcρτcD2R1,τc,D0I3βSI3βIS3+dI3+uI3+cI3γVI3γIV3=D2I1+γV1I2+βS2I1D1I2+γV2I1+βS1I2,D0R3mV3uI3+dR3+μR3,τ1+ρR3,τc=D2R1+μD1R2,τ1+ρτεD0R2,τc+ρτcD1R2,τcD1R2+μD2R1,τ1+ρτεD1R1,τc+ρτcD2R1,τc. $ (B13)

    We substitute Eqs (B7), (B9) and (B10) into the right expression of Eq (B13) and note the coefficient of $ \text{e}^{\text{i}\omega^*T_0} $ as vector $ m_2 $. According to solvability condition $ < h^*, m_2 > = 0 $, we have the expression of $ \frac{\partial G}{\partial T_2} $. Since $ \tau_\varepsilon^2 $ has less impact on the normal form, we can ignore the $ \tau_\varepsilon^2G $ term. Thus, we can obtain:

    $ GT2=HG2ˉG, $ (B14)

    where $ H = \frac{\beta \left(h_{1}l_{3}+\bar{h}_{1}g_{3}+h_{3}l_{1}+\bar{h}_{3}g_{1} \right) \left(-\bar{h}_{1}^{*}+\bar{h}_{3}^{*} \right)}{1+\mu \left(h_{4}\bar{h}_{1}^{*}-h_{4}\bar{h}_{4}^{*} \right) \text{e}^{-\mathrm{i}\omega^*\tau _1^*}+\rho \tau _c\text{e}^{-\mathrm{i}\omega^* \tau _c}\left(h_4\bar{h}_{1}^{*}-h_{4}\bar{h}_{4}^{*} \right)} +\frac{\gamma \left(h_{2}l_{3}+\bar{h}_{2}g_{3}+h_{3}l_{2}+\bar{h}_{3}g_{2} \right) \left(-\bar{h}_{2}^{*}+\bar{h}_{3}^{*} \right)}{1+\mu \left(h_{4}\bar{h}_{1}^{*}-h_{4}\bar{h}_{4}^{*} \right) \text{e}^{-\mathrm{i}\omega^*\tau _1^*}+\rho \tau _c\text{e}^{-i\omega^* \tau _c}\left(h_4\bar{h}_{1}^{*}-h_{4}\bar{h}_{4}^{*} \right)}. $

    Then, we let $ G\rightarrow G/\varepsilon $. Therefore, we obtain the normal form of Hopf bifurcation for system (2.1):

    $ ˙G=MτεG+HG2ˉG, $ (B15)

    where $ M, H $ are given in Eqs (B9) and (B14).


    Abbreviation MRSA: Methicillin resistant ; MSSA: Methicillin susceptible ; MRCoNS: Methicillin reistant Coagulase Negative Staphilococci; PVL: Panton valentine leucocidin; PCR: Polymerase chain reaction; CFU: Colonies forming unit; PBP2a: Penicillin binding protein 2a;

    Conflicts of interests



    All authors declare no conflicts of interest in this paper.

    [1] Gordon RJ, Lowy FD (2008) Pathogenesis of Methicilln-Resistant Staphylococcus aureus Infection. Clin Infect Dis. 46: S350–S359. doi: 10.1086/533591
    [2] Meyer F, Girardot R, Piemont Y, et al. (2009) Analysis of the specificity of Panton-Valentine leucocidin and gamma-hemolysin F component binding. Infect Immun 77: 266–273. doi: 10.1128/IAI.00402-08
    [3] Vandenesch F, Lina G, Gillet Y, et al. (2009) The end of the controversy: Panton Valentine is the culprit. Med Sci 25: 984–986.
    [4] McDonald RR, Antonishyn NA, Hansen T, et al. (2005) Development of a triplex real-time PCR assay for detection of Panton-Valentine leukocidin toxin genes in clinical isolates of methicillin-resistant Staphylococcus aureus. J Clin Microbiol 43: 6147–6149. doi: 10.1128/JCM.43.12.6147-6149.2005
    [5] Shallcross LJ, Fragaszy E, Johnson AM, et al. (2013) The role of Panton-Valentine leucocidin toxin in staphylococcal disease: a systematic review and meta-analysis. Lancet Infect Dis 13: 43–54. doi: 10.1016/S1473-3099(12)70238-4
    [6] Dohin B, Gillet Y, Kohler R, et al. (2007) Pediatric bone and joint infections caused by Panton- Valentine leukocidin-positive Staphylococcus aureus. Pediatr Infect Dis J 26: 1042–1048. doi: 10.1097/INF.0b013e318133a85e
    [7] Gillet Y, Issartel B, Vanhems P, et al. (2002) Association between Staphylococcus aureus strains carrying gene for Panton-Valentine leukocidin and highly lethal necrotising pneumonia in young immunocompetent patients. Lancet 359: 753–759. doi: 10.1016/S0140-6736(02)07877-7
    [8] Lina G, Piémont Y, Godail-Gamot F, et al. (1999) Involvement of Panton-Valentine leukocidin-producing Staphylococcus aureus in primary skin infections and pneumonia. Clin Infect Dis 29: 1128–1132. doi: 10.1086/313461
    [9] Calfee DP, Salgado CD, Milstone AM, et al. (2014) Strategies to Prevent Methicillin-Resistant Staphylococcus aureus Transmission and Infection in Acute Care Hospitals. Infec Cont Hosp Epidemiol 35: 772–796. doi: 10.1086/676534
    [10] Al-Talib H, Yean CY, Al-Khateeb A, et al. (2014) Rapid detection of methicillin-resistant Staphylococcus aureus by a newly developed dry reagent-based polymerase chain reaction assay. J Microbiol Immunol Infection 47: 484–490. doi: 10.1016/j.jmii.2013.06.004
    [11] Johnsson D, Molling P, Stralin K, et al. (2004) Detection of Panton-Valentine leukocidin gene in Staphylococcus aureus by Light Cycler PCR: clinical and epidemiological aspects. Clin Microbiol Infect 10: 884–889. doi: 10.1111/j.1469-0691.2004.00976.x
    [12] Llarrull LI, Fisher JF, Mobashery S (2009) Molecular basis and phenotype of methicillin resistance in Staphylococcus aureus and insights into new β-Lactams that meet the challenge. Antimicrob Agents Chemother 53: 4051–4063. doi: 10.1128/AAC.00084-09
    [13] Teng CS, Lo WT, Wang SR, et al. (2009) The role of antimicrobial therapy for treatment of uncomplicated skin and soft tissue infections from community-acquired methicillin-resistant Staphylococcus aureus in children. J Microbiol Immunol Infect 42: 324–328.
    [14] Boyle-Vavra S, Daum RS (2007) Community-acquired methicillin-resistant Staphylococcus aureus: the role of Panton-Valentine leukocidin. Lab Invest 87: 3–9. doi: 10.1038/labinvest.3700501
    [15] Becker K, Heilmann C, Peters G (2014) Coagulase negative staphylococci. Clin Microbiol Rev 27:870–926. doi: 10.1128/CMR.00109-13
    [16] Stegger M, Andersen PS, Kearns A, et al. (2012) Rapid detection, differentiation and typing of methicillin-resistant Staphylococcus aureus harboring either mecA or the new mecA homologue mecA(LGA251). Clin Microbiol Infect 18: 395–400. doi: 10.1111/j.1469-0691.2011.03715.x
    [17] Söderquist B, Neander M, Dienus O, et al. (2012) Real-time multiplex PCR for direct detection of methicillin-resistant Staphylococcus aureus (MRSA) in clinical samples enriched by broth culture. APMIS 120: 427–432. doi: 10.1111/j.1600-0463.2011.02849.x
    [18] Okolie CE, Wooldridge KG, Turner DP, et al. (2015) Development of a new pentaplex real-time PCR assay for the identification of poly-microbial specimens containing Staphylococcus aureus and other staphylococci, with simultaneous detection so staphylococcal virulence and methicillin resistance markers. Mol Cell Probes 29: 144–150. doi: 10.1016/j.mcp.2015.03.002
  • microbiol-05-02-138-s001.pdf
  • This article has been cited by:

    1. Anuj Kumar, Yasuhiro Takeuchi, Prashant K Srivastava, Stability switches, periodic oscillations and global stability in an infectious disease model with multiple time delays, 2023, 20, 1551-0018, 11000, 10.3934/mbe.2023487
    2. Deepika Solanki, Sumit Kaur Bhatia, Harendra Pal Singh, Praveen Kumar, Dynamic analysis of a communicable disease fractional order model incorporating vaccination and multiple time delays, 2025, 116, 11100168, 147, 10.1016/j.aej.2024.12.055
    3. Mohammed Amine Menouer, Salih Djilali, Anwar Zeb, Ilyas Khan, Abdoalrahman S.A. Omer, The role of imperfect vaccination and treatment in containing an epidemic for a delayed SVITR epidemic model, 2025, 33, 2769-0911, 10.1080/27690911.2025.2461095
  • Reader Comments
  • © 2019 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)
通讯作者: 陈斌, bchen63@163.com
  • 1. 

    沈阳化工大学材料科学与工程学院 沈阳 110142

  1. 本站搜索
  2. 百度学术搜索
  3. 万方数据库搜索
  4. CNKI搜索

Metrics

Article views(10476) PDF downloads(2099) Cited by(36)

Figures and Tables

Figures(3)  /  Tables(1)

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog