Parameters | Biological meaning | Values | Citations |
μ | Death rate | 0.02 | [2,10] |
ν | Vaccination probabilities | 0.8 | [10] |
κ | Coefficient for boosting | 0.03 | estimation |
η | Loss of immunity rate | 0.16 | [2] |
γ | Recovery rate | 0.5 | estimation |
The collision of a vehicle is often a process of strong nonlinearity and large deformation, and the finite element method is time-consuming. For this reason, a method of collision inversion research is proposed. Through the definition of the forward and inverse problems, the forward problem is inversely solved from the perspective of the inverse problem, and the collision process can be predicted quickly while the accuracy is ensured. In this paper, the idea of inversion is introduced into the collision test of the front bumper of the vehicle. First, a finite element model is established based on the geometric model of the bumper. The accuracy of the finite element model is verified by comparing the results of the drop weight test with the simulation results of the finite element model. Then, using the built simulation model, spring mass model and drop weight test, the inversion research of vehicle collision and the front bumper of a vehicle is carried out. The inversion research of the bumper first inverts the collision course in the inversion algorithm derived from the vehicle collision, and then compares the collision course derived from the inversion algorithm formula with the simulation test results. The research results show that the change trends of the time-velocity curve and the time-deformation curve obtained by the collision inversion algorithm are basically consistent with the simulation test results, indicating that the collision inversion algorithm can realize the rapid prediction of the collision course and improve derivation efficiency significantly, and it provides a new idea. Finally, under the constant energy condition of the drop weight test E = mgh, through the relationship between energy and deformation, it is concluded that the depth deformation of low-speed collision of the front bumper is greater than that of high-speed collision.
Citation: Miao Luo, Yousong Chen, Dawei Gao, Lijun Wang. Inversion study of vehicle frontal collision and front bumper collision[J]. Electronic Research Archive, 2023, 31(2): 776-792. doi: 10.3934/era.2023039
[1] | Ceyu Lei, Xiaoling Han, Weiming Wang . Bifurcation analysis and chaos control of a discrete-time prey-predator model with fear factor. Mathematical Biosciences and Engineering, 2022, 19(7): 6659-6679. doi: 10.3934/mbe.2022313 |
[2] | Xiaoling Han, Xiongxiong Du . Dynamics study of nonlinear discrete predator-prey system with Michaelis-Menten type harvesting. Mathematical Biosciences and Engineering, 2023, 20(9): 16939-16961. doi: 10.3934/mbe.2023755 |
[3] | Jinxing Zhao, Yuanfu Shao . Bifurcations of a prey-predator system with fear, refuge and additional food. Mathematical Biosciences and Engineering, 2023, 20(2): 3700-3720. doi: 10.3934/mbe.2023173 |
[4] | Mianjian Ruan, Xianyi Li, Bo Sun . More complex dynamics in a discrete prey-predator model with the Allee effect in prey. Mathematical Biosciences and Engineering, 2023, 20(11): 19584-19616. doi: 10.3934/mbe.2023868 |
[5] | Yuanfu Shao . Bifurcations of a delayed predator-prey system with fear, refuge for prey and additional food for predator. Mathematical Biosciences and Engineering, 2023, 20(4): 7429-7452. doi: 10.3934/mbe.2023322 |
[6] | Xiaoyuan Chang, Junjie Wei . Stability and Hopf bifurcation in a diffusivepredator-prey system incorporating a prey refuge. Mathematical Biosciences and Engineering, 2013, 10(4): 979-996. doi: 10.3934/mbe.2013.10.979 |
[7] | Christian Cortés García, Jasmidt Vera Cuenca . Impact of alternative food on predator diet in a Leslie-Gower model with prey refuge and Holling Ⅱ functional response. Mathematical Biosciences and Engineering, 2023, 20(8): 13681-13703. doi: 10.3934/mbe.2023610 |
[8] | Christian Cortés García . Bifurcations in a discontinuous Leslie-Gower model with harvesting and alternative food for predators and constant prey refuge at low density. Mathematical Biosciences and Engineering, 2022, 19(12): 14029-14055. doi: 10.3934/mbe.2022653 |
[9] | Rajalakshmi Manoharan, Reenu Rani, Ali Moussaoui . Predator-prey dynamics with refuge, alternate food, and harvesting strategies in a patchy habitat. Mathematical Biosciences and Engineering, 2025, 22(4): 810-845. doi: 10.3934/mbe.2025029 |
[10] | Tingting Ma, Xinzhu Meng . Global analysis and Hopf-bifurcation in a cross-diffusion prey-predator system with fear effect and predator cannibalism. Mathematical Biosciences and Engineering, 2022, 19(6): 6040-6071. doi: 10.3934/mbe.2022282 |
The collision of a vehicle is often a process of strong nonlinearity and large deformation, and the finite element method is time-consuming. For this reason, a method of collision inversion research is proposed. Through the definition of the forward and inverse problems, the forward problem is inversely solved from the perspective of the inverse problem, and the collision process can be predicted quickly while the accuracy is ensured. In this paper, the idea of inversion is introduced into the collision test of the front bumper of the vehicle. First, a finite element model is established based on the geometric model of the bumper. The accuracy of the finite element model is verified by comparing the results of the drop weight test with the simulation results of the finite element model. Then, using the built simulation model, spring mass model and drop weight test, the inversion research of vehicle collision and the front bumper of a vehicle is carried out. The inversion research of the bumper first inverts the collision course in the inversion algorithm derived from the vehicle collision, and then compares the collision course derived from the inversion algorithm formula with the simulation test results. The research results show that the change trends of the time-velocity curve and the time-deformation curve obtained by the collision inversion algorithm are basically consistent with the simulation test results, indicating that the collision inversion algorithm can realize the rapid prediction of the collision course and improve derivation efficiency significantly, and it provides a new idea. Finally, under the constant energy condition of the drop weight test E = mgh, through the relationship between energy and deformation, it is concluded that the depth deformation of low-speed collision of the front bumper is greater than that of high-speed collision.
Pertussis, commonly recognized as whooping cough, is an extremely contagious respiratory infection caused by Bordetella pertussis bacteria [1]. The infectious disease is primarily spread through direct contact, such as respiratory droplets. It primarily affects the respiratory system and is characterized by severe coughing. Pertussis predominantly affects infants and young children but can also impact adolescents and adults, often leading to severe complications, particularly in vulnerable populations. The early signs of pertussis resemble common cold symptoms and are easily misdiagnosed. Therefore, it is difficult to quickly identify pertussis in the early stages. The treatment of pertussis mainly relies on antibiotics and symptomatic treatment. Currently, there is no specific drug for pertussis, resulting in no effective elimination of the pathogen. According to data from the Chinese Center for Disease Control and Prevention, there were nearly 500,000 cases of pertussis infections in China in 2024.
It is worth noting that over the past six decades, injectable vaccines have become an effective strategy to prevent diseases. Before the pertussis vaccine was widely available, pertussis caused a significant number of childhood deaths around the world every year. As the pertussis vaccine has been promoted globally, the mortality rate has decreased dramatically. However, in recent decades, there has been a resurgence in pertussis cases [2,3,4,5,6,7,8]. In 2014, 24.1 million cases of pertussis were reported worldwide, resulting in 160,700 deaths among children under 5 years old, particularly in developing countries [9]. This trend can be attributed to two main factors. On the one hand, vaccination rates in developing countries are generally low, resulting in insufficient vaccine coverage to form effective herd immunity. Herd immunity is the idea that when enough of the population is vaccinated, the entire community becomes more resistant to the disease, thereby protecting individuals who cannot be vaccinated, such as newborns or those with immune deficiencies. However, if vaccine coverage is insufficient, the virus can still spread among people, leading to outbreaks [10]. On the other hand, the immunity provided by the vaccine is not permanent, and the effectiveness of the vaccine diminishes over time, especially for adult and adolescent populations, where the protective effect is relatively weak. This means that even if children are fully vaccinated, they may lose their immunity to pertussis over time, increasing the risk of reinfection. It is also important to note that individuals who recover from an infection may gain immunity, but that does not offer complete protection [11]. The impact of immunity on disease transmission has gained great attention, while there are relatively few works related to immune waning. Based on the above background, this paper focuses on a pertussis model with long-term immune waning and natural immune boosting.
Mathematical models are essential for characterizing transmission dynamics of infectious diseases. Over the past few decades, mathematical models of pertussis have provided important insights into the occurrence of new pertussis cases. Tian and Wang [12] proposed a pertussis model incorporating recessive infection, thereby laying a solid foundation for understanding the impact of latent infection on pertussis recurrence. Given the multifaceted nature of pertussis infection, a model accounting for multiple infection pathways was given by [13]. Research on the effects of vaccination on pertussis dynamics has also garnered significant attention. Rozhnova and Nunes [2] focused on the long-term trends of pertussis prior to vaccination, aiming to minimize the number of free parameters. Safan et al.[14] established a pertussis model based on data from newborns, which concluded that symptomatic infections can be completely eradicated through vaccination. Conversely, Águas proposed that the resurgence of pertussis may result from reduced transmission, independent of vaccination [15]. Most existing studies [1,12,13,14,15] have concentrated on deterministic models, often neglecting the influence of environmental noise and its role in pertussis dynamics. Meanwhile, the dynamic behavior of most epidemics (such as pertussis) is greatly influenced by random factors [16,17]. Due to stochastic models having advantages in adapting to complex future conditions [18,19,20,21,22,23,24,25], we apply nonlinear perturbations [21] to a deterministic pertussis model. During the modeling process, we mainly considered the influence of environmental factors such as temperature and humidity changes on the transmission of pertussis. The transmission of pertussis has seasonal variations [26,27], mainly reflected in terms of temperature and humidity. In order to more accurately capture the seasonal transmission characteristics of pertussis, i.e., peak of infection in summer and winter, while infections in spring and autumn are relatively rare, we introduced a Markov switching process [25] in the model, aiming to simulate the influence of seasonal fluctuations of environmental factors (such as temperature and humidity) on the transmission of pertussis.
As far as we are aware, up to now, a pertussis model taking both immune and random factors into consideration has not received any attention. In order to fill this gap, this paper proposes a stochastic pertussis model with Markov switching. We apply Lyapunov functions and the Has'minskii theorem to analyze the pertussis model. The primary contributions and innovations of this paper are described below:
● The transmission rate is more sensitive to the external environment. A stochastic pertussis model with Markov switching of transmission rate is discussed in this paper. This model is an extension of Lavine et al. [10].
● In response to the recurrence and state-switching phenomena of pertussis, it is essential to implement appropriate measures. To this end, sufficient conditions are presented to ensure the extinction of a stochastic pertussis model with Markov switching. Additionally, considering persistence of the pertussis is equally important. We derive persistence conditions for pertussis by proving ergodicity.
The remainder of this paper is structured as follows: Section 2 introduces the mathematical model and necessary lemmas. In Section 3, we establish sufficient conditions for disease extinction. The existence of a unique ergodic stationary distribution (ESD) is demonstrates in Section 4. Section 5 shows numerical simulations to validate our theoretical results and explore the effects of random perturbations and natural immune boosting on disease transmission. Section 6 offers a discussion of our conclusions.
Lavine et al.[10] proposed a deterministic model that is outlined as follows:
{dS(t)dt=μ(1−ν)−(μ+βI(t))S(t)+ηW(t),dI(t)dt=βS(t)I(t)−(μ+γ)I(t),dW(t)dt=ηR(t)−(μ+η+κβI(t))W(t),dR(t)dt=κβI(t)W(t)+γI(t)+μν−(μ+η)R(t), | (2.1) |
where S(t), I(t), W(t), and R(t) are the susceptible compartment, the infectious compartment, the waning compartment, and the recovered compartment at time t, respectively. It is worth noting that R(t) consists of individuals who have recently recovered or been boosted, possessing a strong immune response. W(t) includes individuals who maintain immunity against infection but can further strengthen their immune response upon re-exposure κβI(t)W(t). If individuals in W(t) are not re-exposed, they will eventually lose their immunity entirely and revert to the susceptible group ηW(t). Besides that, μ represents the death rate, ν denotes the vaccination probabilities specific to different age groups, β characterizes the transmission rate, η represents the rate at which immunity wanes when there is no immune boosting, γ represents the rate of recovery, and κ represents the coefficient for boosting. Denoting N(t)=S(t)+I(t)+W(t)+R(t) as the total polulation, and on the basis of (2.1), we can obtain dN(t)=(μ−μN(t))dt. By integrating, we can get
N(t)=1+(N(0)−1)e−μt, | (2.2) |
it will come out N(t)≤1, while N(0)≤1.
Denote the possible region of the determinist system (2.1)
Γ∗=(S(t),I(t),W(t),R(t)∈R4+|0<S(t)+I(t)+W(t)+R(t)≤1). | (2.3) |
The system (2.1) possesses a disease-free equilibrium point E0=(S0,I0,W0,R0)=(μ(1−ν)+ηW0μ,0,ημν(μ+η)2,μνμ+η). Combining the theories of Ma et al.[18] and Sahu et al.[19], the basic reproduction number R0 of the system (2.1) is obtained as follows:
R0=βS0μ+γ=β(μ+γ)((1−ν)+η2ν(μ+η)2). |
By applying the global stability theory of equilibrium, system (2.1) has the following threshold dynamics:
● E0 is globally asymptotically stable [10], which implies the number of infectious will decrease to zero.
● If R0>1, the endemic equilibrium E∗=(S∗,I∗,W∗,R∗) is globally asymptotically stable and satisfies the following equation:
{S∗=μ+γβ,W∗=βμ(ν−1)+(μ+βI∗)(μ+ν)ηβ,R∗=(μ+η+κβI∗)W∗η,κβI∗W∗+γI∗+μν−(μ+η)R∗=0. |
The proof of asymptotic stablility of E∗ is similar to that in Theorem 3.2 of [20], so we omit it here.
Environmental noise is ubiquitous and can significantly impact the transmission of pertussis. The deterministic model (2.1) has inherent limitations in accurately predicting future dynamics. Linear perturbation, a simplifying assumption commonly employed, is commonly used to describe disease spread [29,30,31]. Based on [21,23,25,32], we use nonlinear perturbation for the system (2.1). This study investigates the natural mortality rate of the system (2.1) under nonlinear perturbations, such as μS(t)→μ+(σii+σijS(t))S(t)dBi(t), which is given by
{dS(t)=[μ(1−ν)−(μ+βI(t))S(t)+ηW(t)]dt+(σ11S(t)+σ12)S(t)dB1(t),dI(t)=[βS(t)I(t)−(μ+γ)I(t)]dt+(σ21I(t)+σ22)I(t)dB2(t),dW(t)=[ηR(t)−(μ+η+κβI(t))W(t)]dt+(σ31W(t)+σ32)W(t)dB3(t),dR(t)=[κβI(t)W(t)+γI(t)+μν−(μ+η)R(t)]dt+(σ41R(t)+σ42)R(t)dB4(t), | (2.4) |
where Bi(t)(i=1,2,…,4) refers to four standard Brownian motions, each independent. σ2i1>0, σ2i2>0 describe the intensities of nonlinear and linear perturbations on Bi(t), respectively, i=1,2,…,4. According to Mao [18], the processes Bi(t)(i=1,2,…,4) are defined on a complete probability space {Ω,F, {Ft}t>0, P}, where the filtration {Ft}t≥0 is right continuous and increasing, provided that F0 contains all P−null sets.
In the actual pertussis transmission model, changes in the external environment, especially alterations in temperature and humidity, can lead to variations in parameters. In the modeling process, fixing parameters cannot represent such state changes. Similarly, introducing random noise σ cannot accurately depict such changes. The mortality rate μ, the immunity wanes η, and vaccination probabilities ν are less affected by environmental changes, while the contact rate β is more sensitive to environmental changes [33]. In this paper, we introduce Markov switching to describe the state transitions caused by environmental changes and simplify the model by only considering the influence of the state change on β, described in detail as follows:
{dS(t)=[μ(1−ν)−(μ+β(r(t))I(t))S(t)+ηW(t)]dt+(σ11(r(t))S(t)+σ12(r(t)))S(t)dB1(t),dI(t)=[β(r(t))S(t)I(t)−(μ+γ)I(t)]dt+(σ21(r(t))I(t)+σ22(r(t)))I(t)dB2(t),dW(t)=[ηR(t)−(μ+η+κβ(r(t))I(t))W(t)]dt+(σ31(r(t))W(t)+σ32(r(t)))W(t)dB3(t),dR(t)=[κβ(r(t))I(t)W(t)+γI(t)+μν−(μ+η)R(t)]dt+(σ41(r(t))R(t)+σ42(r(t)))R(t)dB4(t), | (2.5) |
where (r(t))t≥0 is a right-continuous Markov chain taking values in a finite space N={1,2,…,N}. The following is transition probability:
P{r(t+Δ)=j|r(t)=i}={qijΔ+o(Δ),if i≠j,1+qiiΔ+o(Δ),if i=j, |
where the time increment Δ>0, and qij>0 is the transition coefficient from i to j if i≠j. For any i∈N, if j=i, qii=−∑i≠jqij. o(Δ) signifies limΔ→0o(Δ)/Δ=0. In addition, (r(t))t≥0 is irreducible and has a unique stationary distribution π={π1,π2,…,πN}∈R1×n, which is determined by πΓ=0 subject to ∑Nk=1πk=1, πk>0.
We need to give several mathematical notations that will be used in this paper. Let an N-dimensional vector →A=(A(1),A(2),…,A(N))T, assuming that ˇA=maxk∈S{A(k)}, ˆA=mink∈S{A(k)}.
To investigate the dynamical behavior of the stochastic epidemic model (2.5), it is essential to establish the existence of a global positive solution. Since S(t), I(t), W(t), and R(t) represent numbers of different kinds of individuals, these values should remain non-negative. And based on the background of pertussis transmission, we assume μ≥0,ν≥0,β≥0,γ≥0,κ≥0,η≥0. In this section, we will show that system (2.5) possesses a unique global positive solution.
Theorem 2.1. For any initial condition (S(0),I(0),W(0),R(0),r(0))∈R4+×N, system (2.5) admits a unique solution (S(t),I(t),W(t),R(t),r(t)) on t≥0, and the solution will remain in R4+×N with probability one, a.s..
Proof. We omit the proof, which is similar to Zu et al. [24] and Zhou et al. [25].
Remark 2.1. The existence and uniqueness of the global positive solution are essential for ensuring the reliability of predictions in the stochastic pertussis model. Without establishing these fundamental properties, the model's predictive power would be compromised, as multiple or non-positive solutions could arise, leading to ambiguous or incorrect interpretations of the disease dynamics.
The extinction of diseases is an important topic discussed in the field of biomathematics. In this section, we seek to establish the sufficient conditions for the extinction of random systems (2.5). To this end, a value ˜R0 associated with R0 is defined below
˜R0=∑Ni=1πiβ(i)μ+γ+∑Ni=1πiσ22(i)2. |
Theorem 3.1. For any initial value (S(0),I(0),W(0),R(0))∈R4+×N, let (S(t),I(t),W(t),R(t)) be the solution of stochastic system (2.5). If ˜R0<1 holds, then the disease of the stochastic system (2.5) will go extinct.
Proof. Utilizing the Itˆo's formula yields
d(lnI(t))=1I(t)[β(r(t))S(t)I(t)−(μ+γ)I(t)]dt−(σ21(r(t))I(t)2+σ22(r(t)I(t))22I(t)2dt+1I(t)(σ21(r(t))I(t)2+σ22(r(t))I(t))dB2(t)=[β(r(t))S(t)−(μ+γ+(σ21(r(t))I(t)+σ22(r(t))22)]dt+(σ21(r(t))I(t)+σ22(r(t)))dB2(t). | (3.1) |
Integrate both sides of Eq (3.1) from 0 to t, we can obtain
lnI(t)t≤lnI(0)t+1t∫t0[β(r(u))S(u)−(μ+γ+(σ21(r(u))I(t)+σ22(r(u))22)]du+Φ(t)t≤lnI(0)t+1t∫t0[β(r(u))N(u)−(μ+γ+σ22(r(u))22)]du+Φ(t)ta.s.≤lnI(0)t+1t∫t0[β(r(u))−(μ+γ+σ22(r(u))22)]du+Φ(t)ta.s. | (3.2) |
where Φ(t)=∫t0(σ21(r(u))I+σ22(r(u)))dB2(u) is a continuous real-valued local martingale and its quadratic variation is
⟨Φ,Φ⟩t=∫t0(σ21(r(u))I+σ22(r(u)))2du. |
Applying the exponential martingale inequality of [18, Theorem 7.4], one gets
P{sup0≤t≤T[Φ(t)−ϵ2∫t0(σ21(r(u))I+σ22(r(u)))2du]>2ϵlnk}≤1k2, |
wherein 0<ϵ<1, k is a random positive number. Utilizing the Borel-Cantelli's lemma [18, Lemma 2.4] for every ω∈Ω, there exists a integer k0 such that all k≥k0,
Φ(t)≤ϵ2∫t0(σ21(r(u))I+σ22(r(u)))2du+2ϵlnk. | (3.3) |
Combining (3.2) and (3.3) for all t∈(k−1,k], and letting ϵ→0, one has
lnI(t)t≤lnI(0)t+1t∫t0[β(r(u))−(μ+γ+σ22(r(u))22)]du+ϵ2t∫t0(σ21(r(u))I+σ22(r(u)))2du+2lnkϵt≤lnI(0)t+1t∫t0[β(r(u))−(μ+γ+1t∫t0σ22(r(u))22)]du+2lnkϵ(k−1). |
If k→+∞, then t→+∞, one obtains that lnI(0)t→0, lnkk−1→0. By combining the dominated convergence theorem and the ergodic properties of Markov chains, one can gain
lim supt→∞1t∫t0σ22(r(u))22du=N∑i=1πiσ22(i)2,lim supt→∞1tβ(r(u))du=N∑i=1πiβ(i). | (3.4) |
Combining (3.4) and ˜R0<1, we have
lim supt→∞lnI(t)t≤lim supt→∞1t∫t0[β−(μ+γ+σ22(r(u))22)]du=N∑i=1πiβ(i)−(μ+γ+N∑i=1πiσ22(i)2)=(μ+γ+N∑i=1πiσ22(i)2)[∑Ni=1πiβ(i)(μ+γ+∑Ni=1πiσ22(i)2)−1]=(μ+γ+N∑i=1πiσ22(i)2)(˜R0−1)<0, | (3.5) |
which shows limt→∞I(t)=0.
Remark 3.1. In Theorem 3.1, we utilize the property of S(t)≤N(t) almost surely during stochastic perturbations. When the external environment remains stable and the transmission coefficient does not undergo sudden changes N=1, R0 will be a threshold result common to deterministic system (2.1) and stochastic system (2.5). Based on ˜R0, the measures of isolating close contacts and wearing masks can effectively reduce transmission rate and curb the spread of pertussis.
Whether a disease will persist is a hot topic of current research [23,24,25,36,37,38]. Based on the Has'minskii theorem [22], an endemic equilibrium does not exist in the stochastic epidemic model, but an ESD can still exist. The ESD reflects the persistence of disease transmission within the community [16,39].
Denote ¯μ1=μ+∑Nk=1πkσ212(k)2+2√μ∑Nk=1πkσ11(k)σ12(k)+23√μ2∑Nk=1πkσ211(k), μ2=μ+γ+∑Nk=1πkσ222(k)2, ¯μ3=μ+η+∑Nk=1πkσ232(k)2, μ4=μ+η+∑Nk=1πkσ242(k)2, Re0=ˆβSeμ2.
Theorem 4.1. Suppose that Re0>1. Then the solution (S(t),I(t),W(t),R(t),r(t)) of system (2.5) has a distinct invariant distribution regardless of the initial state (S(0),I(0),W(0),R(0),r(0))∈R4+×N.
Proof. The diffusion matrix related to system (2.5) is illustrated by
A=((σ11(k)S2+σ12(k)S)20000(σ21(k)E2+σ22(k)E)20000(σ31(k)I2+σ32(k)I)20000(σ41(k)R2+σ42(k)R⋁)2). |
Opt for ˜M=min(S,I,W,R)∈˜Dϵ⊂R4+{(σ11(k)S2+σ12(k)S)2,(σ21(k)E2+σ22(k)E)2,(σ31(k)I2+σ32(k)I)2,(σ41(k)A2+σ42(k)A)2}>0, ensuring
4∑i,j=1aij(S,I,W,R)ζiζj=(σ11(k)S2+σ12(k)S)2ζ21+(σ21(k)I2+σ22(k)I)2ζ22+(σ31(k)W2+σ32(k)W)2ζ23+(σ41(k)R2+σ42(k)R)2ζ24≥˜M‖ζ‖2, |
wherein ζ:=(ζ1,ζ2,ζ3,ζ4). Inspired by [25], we will employ the generalized θ-stochastic criterion method. To verify condition (ⅱ), we have the following four steps.
Step 1. For all τ∈(0,1), here are some continuous τ-stochastic Lyapunov functions to consider:
μ1=μ+N∑k=1πkσ212(k)2+2√μ∑Nk=1πkσ11(k)σ12(k)1−τ+23√μ2∑Nk=1πkσ211(k)(1−τ)2,μ3=μ+η+N∑k=1πkσ232(k)2+τ2σ231(k)6,Re0(τ)=ˆβSeμ2. |
We observe that μi (i = 1, …, 4) are all functions of the given variable τ∈(0,1) that increase monotonically. Furthermore, it is found that
infτ∈(0,1)μi=limτ→0+μi=μ∗i,i=1,3.supτ∈(0,1)Re0(τ)=limτ→0+Re0(τ)=Re0. |
It should be pointed out that limτ→0+Se=¯Se.
Step 2. To counteract the effects of nonlinear perturbations and Markov switching, we propose the following Lyapunov functions:
V1=−lnS+2∑j=1ζj(S+ξj)ττ+ω1(k),V2=−lnI+ω2(k),V3=−lnW+ζ3(W+τ)ττ+ω3(k),V4=−lnR+ω4(k),∀k∈N, |
where ζ1,ζ2,ζ3, and ξ1,ξ2 are determined by (4.2), (4.10), and (4.5). ωi(k),i=1,2,3,4 are determined subsequently.
Utilizing the Itˆo′s formula to V1, combining
a3⩾(a−12)(a2+1)anda4⩾14(3a2−1)(a2+1). |
and (a+b)2≤2(a2+b2), we get
LV1=−1S[μ(1−ν)−(μ+β(k)I)S+ηW]+2∑j=1ζj(S+ξj)τ−1[μ(1−ν)−(μ+βI)S+ηW]+12(σ11(k)S+σ12(k))2−2∑j=1ζj(1−τ)2(S+ξj)2−τ(σ11(k)S2+σ12(k)S)2+∑l∈Nqklω1(l)≤−μS+μνS+μ+ˇβI−ηWS+σ212(k)2+σ11(k)σ12(k)S+σ211(k)2S2+2∑j=1μζjξ1−τj−2∑j=1ζjξτ−2j(1−τ)2(1+Sξj)2−τ(σ11(k)S2+σ12(k)S)2+∑l∈Nqklω1(l)≤−μS+μνS+μ+ˇβI−ηWS+σ212(k)2+σ11(k)σ12(k)S+σ211(k)2S2+2∑j=1μζjξ1−τj−ζ1ξτ−21(1−τ)σ11(k)σ12(k)S3(1+Sξ1)2−ζ2ξτ−22(1−τ)σ211(k)S42(1+Sξ2)2+∑l∈Nqklω1(l)≤−μS+μνS+μ+ˇβI−ηWS+σ212(k)2+σ11(k)σ12(k)S+σ211(k)2S2+2∑j=1μζjξ1−τj−ζ1ξτ+11(1−τ)σ11(k)σ12(k)(Sξ1)32[1+(Sξ1)2]−ζ2ξτ+22(1−τ)σ211(k)(Sξ2)44[1+(Sξ2)2]+∑l∈Nqklω1(l)≤−μS+μνS+μ+ˇβI−ηWS+σ212(k)2+σ11(k)σ12(k)S+σ211(k)2S2+2∑j=1μζjξ1−τj−ζ1ξτ+11(1−τ)σ11(k)σ12(k)(Sξ1−12)2−ζ2ξτ+22(1−τ)σ211(k)[3(Sξ2)2−1]16+∑l∈Nqklω1(l) | (4.1) |
Choosing
ζ1=2(1−τ)ξτ1,ζ2=83(1−τ)ξτ2, | (4.2) |
we gain
LV1≤−μS+μνS+μ+ˇβI−ηWS+σ212(k)2+2μ(1−τ)ξ1+ξ1σ11(k)σ12(k)2+8μ3(1−τ)ξ2+ξ22σ211(k)6+∑l∈Nqklω1(l). | (4.3) |
Let Z1(k):=σ212(k)2+2μ(1−τ)ξ1+ξ1σ11(k)σ12(k)2+8μ3(1−τ)ξ2+ξ22σ211(k)6. Motivated by [38], due to Γ being irreducible, for →Z1=(Z1(1),Z1(2),…,Z1(N))τ, we can get a vector →ω1=(ω1(1),ω1(2),…,ω1(N))τ so that the following Poisson system Γ→ω1=∑Nl=1πlZ1(l)−→Z1 is true, which indicates
Z1(k)+∑l∈Mqklω1(l)=N∑k=1πkZ1(k),∀k∈N. | (4.4) |
Combining (4.3) and (4.4), it is clear that
LV1≤−μS+μνS+μ+ˇβI−ηWS+N∑k=1πkσ212(k)2+2μ(1−τ)ξ1+ξ1∑Nk=1πkσ11(k)σ12(k)2+8μ3(1−τ)ξ2+ξ22∑Nk=1πkσ211(k)6. |
In order to ensure 2μ(1−τ)ξ1+ξ1∑Nk=1πkσ11(k)σ12(k)2+8μ3(1−τ)ξ2+ξ22∑Nk=1πkσ211(k)6 reaches a minimum value, we opt for
ξ1=2√μ(1−τ)∑Nk=1πkσ11(k)σ12(k),ξ2=23√μ(1−τ)∑Nk=1πkσ211(k). | (4.5) |
We obtain
LV1≤−μS+μνS+ˇβI−ηWS+μ1. | (4.6) |
Employing the Itˆo's formula to V2 and referring to the similar methods described earlier, we get
LV2=−[β(k)S−(μ+γ)]+12(σ21(k)I+σ22(k))2+∑l∈Nqklω2(l)≤−ˆβS+μ+γ+ˇσ2212I2+ˇσ21ˇσ22I+σ222(k)2+∑l∈Nqklω2(l):=−ˆβS+ˇσ21ˇσ22I+Z2(k)+ˇσ2212I2+∑l∈Nqklω2(l). | (4.7) |
where Z2(k)=μ+γ+σ222(k)2, a vector →ω2=(ω2(1),ω2(2),…,ω2(N))τ is chosen to meet the Poisson system Γ→ω2=∑Nl=1πlZ2(l)−→Z2, where→Z2=(→Z2(1),→Z2(2),…,→Z2(N))τ. Therefore
LV2≤−ˆβS+ˇσ21ˇσ22I+ˇσ2212I2+μ2. | (4.8) |
In light of the similar steps mentioned in (4.7), we can derive
LV3≤−ηRW+μ+η+κˇβI+(σ31(k)W+σ32(k))22+ζ3ηRτ1−τ−ζ3(1−τ)ττ−2σ231(k)W42(1+(Wτ))2−τ+∑l∈Nqklω3(l)≤−ηRW+μ+η+κˇβI+ˇσ31ˇσ32W+σ232(k)2+σ231(k)2W2+ζ3ηRτ1−τ−ζ3(1−τ)ττ+2σ231(k)(Wτ)44(1+(Wτ)2)+∑l∈Nqklω3(l). | (4.9) |
We select
ζ3=83(1−τ)ττ, | (4.10) |
and then derive
LV3≤−ηRW+μ+η+κˇβI+ˇσ31ˇσ32W+σ232(k)2+σ231(k)2W2+8ηR3τ(1−τ)−τ2σ231(k)[3(Wτ)2−1]6+∑l∈Nqklω3(l)≤−ηRW+κˇβI+ˇσ31ˇσ32W+8ηR3τ(1−τ)+μ+η+σ232(k)2+τ2σ231(k)6+∑l∈Nqklω3(l). | (4.11) |
Let Z3:=μ+η+σ232(k)2+τ2σ231(k)6, a vector →ω3=(ω3(1),ω3(2),⋯,ω3(N))τ can be found that satisfies the Poisson system Γ→ω3=∑Nl=1πlZ3(l)−→Z3, wherein →Z3=(Z3(1),Z3(2),⋯,Z3(N))τ. Hence
LV3≤−ηRW+κˇβI+ˇσ31ˇσ32W+8ηR3τ(1−τ)+μ+η+N∑k=1πkσ232(k)2+τ2σ231(k)6=−ηRW+κˇβI+ˇσ31ˇσ32W+8ηR3τ(1−τ)+μ3. | (4.12) |
Applying the Itˆo's formula to V4, we can get
LV4≤−κˆβIWR−γIR−μνR+μ+η+(σ41(k)R+σ42(k))22++∑l∈Nqklω4(l)≤−κˆβIWR−γIR−μνR+μ+η+ˇσ41ˇσ42R+σ242(k)2+ˇσ2412R2+∑l∈Nqklω4(l). | (4.13) |
Define Z4=μ+η+σ242(k)2. Using the same method mentioned in (4.4), one can obtain
LV4≤−κˆβIWR−γIR−μνR+ˇσ41ˇσ42R+ˇσ2412R2+μ+η+N∑k=1πkσ242(k)2=−κˆβIWR−γIR−μνR+ˇσ41ˇσ42R+ˇσ2412R2+μ4. | (4.14) |
Step 3. Define two Lyapunov functions Y1,Y2:R5+×N→R that play important roles
Y1=α1SeV1+1SeV2+α2WeV3+α3ReV4,Y2=Y1+ϕW+ϕR, | (4.15) |
wherein the constants (α1,α2,α3) and ϕ, which are positive, are derived from (4.25) and (4.28), respectively. Suppose Se,Ie,We,Re satisfy the following equations:
{−ηRe+μ3We=0,−κˆβIeWe−γIe−μν+μ4Re=0,−μ+μν−ηWe+μ1Se+ˇβIeSe=0. | (4.16) |
Let Ie=1, meanwhile, we can get
Se=μ+ηWe−μνμ1+ˇβ>0,We=η(μν+γ)μ3μ4−ηκˆβ>0,Re=μ3(μν+γ)μ3μ4−ηκˆβ>0. |
Take a standard transformation of variables (S,I,W,R)
ˆS:=SSe,ˆI:=IIe,ˆW:=WWe,ˆR:=RRe. | (4.17) |
Considering that
−ˆβ+μ2Se=−μ2Se(ˆβSeμ2−1)=−μ2Se(Re0−1). | (4.18) |
From the inequality lna≤a−1(∀a>0), using the formulas (4.6), (4.16), and (4.17), we get
L(SeV1)≤Se(−μS+μνS+ˇβI−ηWS+μ1)≤Se(−μSeˆS+μνSeˆS+ˇβ(1+lnI)−ηWeˆWSeˆS+μ1)≤Se[−μSe(1+ln1ˆS)+μνSe(1+ln1ˆS)+ˇβ+ˇβlnI−ηWeSe(1+lnˆWˆS)+μ1]=−μ+μlnˆS+μν−μνlnˆS+ˇβSe+ˇβSelnI−ηWe−ηWelnˆW+ηWelnˆS+μ1Se=μlnˆS+ˇβSelnI+ηWelnˆS−μνlnˆS−ηWelnˆW. | (4.19) |
Using the same method as above, together with formula (4.8) and (4.18), we can obtain
L(1SeV2)≤1Se[−ˆβS+ˇσ21ˇσ22I+ˇσ2212I2+μ2]≤1Se[−ˆβSe(1+lnˆS)+ˇσ21ˇσ22I+ˇσ2212I2+μ2]=−ˆβ−ˆβlnˆS+ˇσ21ˇσ22ISe+ˇσ2212SeI2+μ2Se=−μ2Se(Re0−1)−ˆβlnˆS+ˇσ21ˇσ22ISe+ˇσ2212SeI2. | (4.20) |
Based on (4.12) and (4.16), we derive
L(WeV3)≤We[−ηRW+κˇβI+ˇσ31ˇσ32W+8ηR3τ(1−τ)+μ3]≤We[−ηReWe(1+lnˆRˆW)+κˇβI+ˇσ31ˇσ32W+8ηR3τ(1−τ)+μ3]=−ηRe−ηRelnˆR+ηRelnˆW+ˇσ31ˇσ32WWe+κˇβIWe+8ηRWe3τ(1−τ)+μ3We=−ηRelnˆR+ηRelnˆW+ˇσ31ˇσ32WWe+κˇβIWe+8ηRWe3τ(1−τ). | (4.21) |
According to (4.14) and (4.16), we have
L(ReV4)≤Re(−κˆβIWR−γIR−μνR+ˇσ41ˇσ42R+ˇσ2412R2+μ4)≤Re[−κˆβWeRe(1+lnIˆWˆR)−γRe(1+lnIˆR)−μνRe(1+ln1ˆR)+ˇσ41ˇσ42R+ˇσ2412R2+μ4]=−κˆβWe−κˆβWelnˆW−κˆβWelnI+κˆβWelnˆR−γ−γlnI+γlnˆR−μν+μνlnˆR+ˇσ41ˇσ42RRe+ˇσ2412R2Re+μ4Re=−κˆβWelnˆW−κˆβWelnI+κˆβWelnˆR−γlnI+γlnˆR+μνlnˆR+ˇσ41ˇσ42RRe+ˇσ2412R2Re. | (4.22) |
Substituting (4.19)–(4.22) into Y1, we obtain
LY1≤−μ2Se(Re0−1)−βlnˆS+ˇσ21ˇσ22ISe+ˇσ2212SeI2+α1(μlnˆS+ˇβSelnI+ηWelnˆS−μνlnˆS−ηWelnˆW)+α2[−ηRelnˆR+ηRelnˆW+ˇσ31ˇσ32WWe+κˇβIWe+8ηRWe3τ(1−τ)]+α3(−κˆβWelnˆW−κˆβWelnI+κˆβWelnˆR−γlnI+γlnˆR+μνlnˆR+ˇσ41ˇσ42RRe+ˇσ2412R2Re)=−μ2Se(Re0−1)+[α1(μ+ηWe−μν)−β]lnˆS+[α1ˇβSe−α3(γ+κˆβWe)]lnI+[α2ηRe−α1ηWe−α3κˆβWe]lnˆW+[α3(κˆβWe+μν+γ)−α2ηRe]lnˆR+[ˇσ21ˇσ22Se+α1ˇβSe+α2κˇβWe]I+α2ˇσ31ˇσ32WeW+[α3ˇσ41ˇσ42Re+α28ηWe3τ(1−τ)]R+ˇσ2412ReR2+ˇσ2212SeI2. | (4.23) |
Let
{α1(μ+ηWe−μν)−β=0,α1ˇβSe−α3(γ+κˆβWe)=0,α2ηRe−α1ηWe−α3κˆβWe=0,α3(κˆβWe+μν+γ)−α2ηRe=0. | (4.24) |
We derive through calculation
α1=βμ+ηWe−μν>0,α3=α1ˇβSeγ+κˆβWe>0,α2=α3(κˆβWe+μν+γ)ηRe>0. | (4.25) |
For the sake of convenience, we define χ1=ˇσ21ˇσ22Se+α1ˇβSe+α2κˇβWe,χ2=α2ˇσ31ˇσ32We,χ3=α3ˇσ41ˇσ42Re+α28ηWe3τ(1−τ). Combining (4.23) and (4.24), we get
LY1≤−μ2Se(Re0−1)+χ1I+χ2W+χ3R+ˇσ2412ReR2+ˇσ2212SeI2. | (4.26) |
By (4.15) and (4.26), we have
LY2≤−μ2Se(Re0−1)+χ1I+χ2W+χ3R+ϕ[ηR−(μ+η+κβI)W]+ϕ[κβIW+γI+μν−(μ+η)R]+ˇσ2412ReR2+ˇσ2212SeI2=−μ2Se(Re0−1)+(χ1+ϕγ)I+[χ2−ϕ(μ+η)]W+(χ3−ϕμ)R+ϕμν+ˇσ2412ReR2+ˇσ2212SeI2, | (4.27) |
wherein ϕ fulfills the following two equations:
{χ2−ϕ(μ+η)=0χ3−ϕμ=0. | (4.28) |
Substituting (4.27) into (4.28), we obtain
LY2≤−μ2Se(Re0−1)+(χ1+ϕγ)I+ϕμν+ˇσ2412ReR2+ˇσ2212SeI2. | (4.29) |
Step 4. Define the following two functions Y3, Y4:
Y3=(ˆσ11S+ˆσ12)ττ+(ˆσ21I+ˆσ22)ττ+(ˆσ31W+ˆσ32)ττ+(ˆσ41R+ˆσ42)ττ,Y4=−lnS−lnW−lnR. |
L(Y3+Y4)≤ˆσ11ˆστ−112μ+ˆσ41ˆστ−142μν+3μ+2η+ˇσ11ˇσ12S+(ˆσ41ˆστ−142γ+β+κˇβ)I+(ˆσ11ˆστ−112η+ˇσ31ˇσ32)W+(ˆσ31ˆστ−132η+ˇσ41ˇσ42)R+ˆσ21ˆστ−122βSI+ˆσ41ˆστ−142κˇβIW+ˇσ2122+ˇσ2322+ˇσ2422+ˇσ2112S2+ˇσ2312W2+ˇσ2412R2−μS−ηRW−γIR−(1−τ)ˆσ2+τ11S2+τ2−(1−τ)ˆσ2+τ21I2+τ2−(1−τ)ˆσ2+τ31W2+τ2−(1−τ)ˆσ2+τ41R2+τ2≤ˆσ11ˆστ−112μ+ˆσ41ˆστ−142μν+3μ+2η+Q1+ˇσ2122+ˇσ2322+ˇσ2422−μS−ηRW−γIR−(1−τ)ˆσ2+τ11S2+τ4−(1−τ)ˆσ2+τ21I2+τ4−(1−τ)ˆσ2+τ31W2+τ4−(1−τ)ˆσ2+τ41R2+τ4, | (4.30) |
where
Q1=sup(S,I,W,R)∈R4+{ˇσ11ˇσ12S+(ˆσ41ˆστ−142γ+β+κˇβ)I+(ˆσ11ˆστ−112η+ˇσ31ˇσ32)W+(ˆσ31ˆστ−132η+ˇσ41ˇσ42)R+ˆσ21ˆστ−122βSI+ˆσ41ˆστ−142κˇβIW+ˇσ2112S2+ˇσ2312W2+ˇσ2412R2−(1−τ)ˆσ2+τ11S2+τ4−(1−τ)ˆσ2+τ21I2+τ4−(1−τ)ˆσ2+τ31W2+τ4−(1−τ)ˆσ2+τ41R2+τ4}<+∞. |
Constructing the following form
V(S,I,W,R,k)=(M0Y2+Y3+Y4)−¯Y, |
wherein M0>0 is large enough to fulfill the following inequality:
−M0μ2Se(Re0−1)+M0ϕμν+ˆσ11ˆστ−112μ+ˆσ41ˆστ−142μν+3μ+2η+Q1+ˇσ2122+ˇσ2322+ˇσ2422≤−2, | (4.31) |
LV≤M0[−μ2Se(Re0−1)+(χ1+ϕγ)I+ϕμν+ˇσ2412ReR2+ˇσ2212SeI2]+ˆσ11ˆστ−112μ+ˆσ41ˆστ−142μν+3μ+2η+Q1+ˇσ2122+ˇσ2322+ˇσ2422−μS−ηRW−γIR−(1−τ)ˆσ2+τ11S2+τ4−(1−τ)ˆσ2+τ21I2+τ4−(1−τ)ˆσ2+τ31W2+τ4−(1−τ)ˆσ2+τ41R2+τ4≤−2+M0[(χ1+ϕγ)I+ˇσ2412ReR2+ˇσ2212SeI2]−μS−ηRW−γIR−(1−τ)ˆσ2+τ11S2+τ4−(1−τ)ˆσ2+τ21I2+τ4−(1−τ)ˆσ2+τ31W2+τ4−(1−τ)ˆσ2+τ41R2+τ4. | (4.32) |
We define the following compact set:
D={(S,I,W,R)∈R|ϵ≤S≤1ϵ,ϵ≤I≤1ϵ,ϵ3≤W≤1ϵ3,ϵ2≤R≤1ϵ2}, |
where ϵ<1 is a small enough positive constant that satisfies the following inequalities:
−2−(1−τ)min{ˆσ2+τ11,ˆσ2+τ21}8(1ϵ)2+τ+Q2≤−1, | (4.33) |
−2−(1−τ)ˆσ2+τ314(1ϵ)6+3τ+Q2≤−1, | (4.34) |
−2−(1−τ)ˆσ2+τ414(1ϵ)4+2τ+Q2≤−1, | (4.35) |
−2−min{μ,η,γ}ϵ+Q2≤−1, | (4.36) |
−2+M0(χ1+ϕγ+ˇσ241Re+ˇσ2212Se)ϵ≤−1, | (4.37) |
where Q2=sup(I,R)∈R2+{M0[(χ1+ϕγ)I+ˇσ2412ReR2+ˇσ2212SeI2]}. Next, we partition R4+∖Dϵ into eight domains.
Dc1={(S,I,W,R)∈R5+|S>1ϵ},Dc2={(S,I,W,R)∈R5+|I>1ϵ},Dc3={(S,I,W,R)∈R5+|W>1ϵ3},Dc4={(S,I,W,R)∈R5+|R>1ϵ2},Dc5={(S,I,W,R)∈R5+|S<ϵ},Dc6={(S,I,W,R)∈R5+|I<ϵ,R≥ϵ2},Dc7={(S,I,W,R)∈R5+|W<ϵ3R≥ϵ2},Dc8={(S,I,W,R)∈R5+|R<ϵ2,I≥ϵ}. |
It is clear that Dcϵ=∪8i=1Dci. After this, to demonstrate LV≤−1.
Case 1. Suppose (S,I,W,R)∈Dc1, in consideration of (4.32) and (4.33), we gain
LV≤−2−(1−τ)ˆσ2+τ118S2+τ+Q2≤−2−(1−τ)min{ˆσ2+τ11,ˆσ2+τ21}8(1ϵ)2+τ+Q2≤−1. |
Case 2. If (S,I,W,R)∈Dc2, combining (4.32) and (4.33), we have
LV≤−2−(1−τ)ˆσ2+τ218I2+τ+Q2≤−2−(1−τ)min{ˆσ2+τ11,ˆσ2+τ21}8(1ϵ)2+τ+Q2≤−1. |
Case 3. If (S,I,W,R)∈Dc3, in consideration of (4.32) and (4.34), we obtain
LV≤−2−(1−τ)ˆσ2+τ314W2+τ+Q2≤−2−(1−τ)ˆσ2+τ314(1ϵ)6+3τ+Q2≤−1. |
Case 4. If (S,I,W,R)∈Dc4, given (4.32) and (4.35), we get
LV≤−2−(1−τ)ˆσ2+τ414R2+τ+Q2≤−2−(1−τ)ˆσ2+τ414(1ϵ)4+2τ+Q2≤−1. |
Case 5. If (S,I,W,R)∈Dc5, by (4.32) and (4.36), we gain
LV≤−2−μS+Q2≤−2−μϵ+Q2≤−1≤−2−min{μ,η,γ}ϵ+Q2≤−1. |
Case 6. If (S,I,W,R)∈Dc6, combining (4.32) and (4.36), one gets
LV≤−2+M0[(χ1+ϕγ)I+ˇσ2412ReR2+ˇσ2212SeI2]≤−2+M0[(χ1+ϕγ)ϵ+ˇσ2412Reϵ4+ˇσ2212Seϵ2]≤−2+M0(χ1+ϕγ+ˇσ241Re+ˇσ2212Se)ϵ≤−1. |
Case 7. If (S,I,W,R)∈Dc7, in view of (4.32) and (4.36), we acquire
LV≤−2−ηRW+Q2≤−2−ηϵ2ϵ3+Q2≤−1≤−2−min{μ,η,γ}ϵ+Q2≤−1. |
Case 8. If (S,I,W,R)∈Dc8, based on (4.32) and (4.36), we have
LV≤−2−γIR+Q2≤−2−γϵϵ2+Q2≤−1≤−2−min{μ,η,γ}ϵ+Q2≤−1. |
To summarize, we establish that LV≤−1 for all (S,I,W,R)∈Dcϵ. Thus, the second condition of Has'minskii theorem is met. Given that Re0>1, the solution of system (2.5) possesses a unique ESD π(⋅).
Remark 4.1. We can calculate Re0≤R0. Theorem 4.1 implies that as the noise intensity approaches zero, there is not infection and N=1, system (2.5) transforms into system (2.1), and the basic production number R0 serves as a common threshold. We obtain that when the noise intensity is relatively low, pertussis still prevails. Integrating the aforementioned Theorems 4.1 and 3.1, it is sensible to strengthen external intervention, wear masks, and appropriately control the number of students in the classroom to mitigate the spread of pertussis.
In this section, numerical simulations are utilized to validate our conclusions. Using Milstein's higher-order method [40], the discretization equation associated with system (2.5) is given by:
{Sj+1=Sj+[μ(1−ν)−(μ+β(k)I)Sj+ηWj]Δt+(σ11(k)Sj+σ12(k))Sj√Δtξj+12[2σ211(k)(Sj)3+3σ11(k)σ12(k)(Sj)2+σ212(k)(Sj)](ξ2j−1)Δt,Ij+1=Ij+[β(k)SjIj−(μ+γ)Ij]Δt+(σ21(k)Ij+σ22(k))Ij√Δtηj+12[2σ221(k)(Ij)3+3σ21(k)σ22(k)(Ij)2+σ222(k)(Ij)](η2j−1)Δt,Wj+1=Wj+[ηRj−(μ+η+κβ(k)I)Wj]Δt+(σ31(k)Wj+σ32(k))Wj√Δtζj+12[2σ231(k)(Wj)3+3σ31(k)σ32(k)(Wj)2+σ232(k)(Wj)](ζ2j−1)Δt,Rj+1=Rj+[κβ(k)IjWj+γIj+μν−(μ+η)Rj]Δt+(σ41(k)Rj+σ42(k))Rj√Δtϑj+12[2σ241(k)(Rj)3+3σ41(k)σ42(k)(Rj)2+σ242(k)(Wj)](ϑ2j−1)Δt, |
of which k∈N, Δt>0 represents the size of a single iteration step and random variables ξj, ηj, ζj and ϑj are four independent variables that adhere to a Gaussian distribution N(0,1) for j=1,2,…,n, respectively. We consider r(t) as a right-continuous Markov chain defined on the state space N={1,2}, with the corresponding generator
Γ=(−0.40.40.5−0.5). |
By solving πΓ=0 and π1+π2=1, the stationary distribution of Γ is obtained as π=(π1,π2)=(59,49). The initial values are chosen as (S(0),I(0),W(0),R(0))=(0.45,0.35,0.15,0.05) and the primary parameters in system (2.5) are set as follows:
Parameters | Biological meaning | Values | Citations |
μ | Death rate | 0.02 | [2,10] |
ν | Vaccination probabilities | 0.8 | [10] |
κ | Coefficient for boosting | 0.03 | estimation |
η | Loss of immunity rate | 0.16 | [2] |
γ | Recovery rate | 0.5 | estimation |
It is worth noting that, we estimate the parameters of model (2.1) by using the least squares method, and take the year as the unit of measurement. Based on [2] and [10], μ=0.02, ν=0.8, η=0.16 can be chosen, respectively. Then, by analyzing the data of pertussis infections from 2006 to 2019 (selected from the National Health Commission of the People's Republic of China, see Figure 1), κ=0.03 and γ=0.5 have been selected. In Figure 1, the susceptible initial value is set to 1030, the infectious initial value is set to 2100, the waning initial value is set to 3000, and the recovered initial value is set to 3000.
By direct calculation, R0=3.2004>1. We focus on three aspects, as follows:
● The system (2.5) shows both a stationary distribution and ergodic behavior when Re0>1.
● The dynamical characteristics of system (2.5) for Re0<1.
● The effect of stochastic perturbations on the disease dynamics within system (2.5).
Example 5.1. We opt for the transmission rates (β1,β2)=(1,2) and the stochastic noises (σij(1),σij(2))=(2×10−4,6×10−4) for any i=1,2,3,4;j=1,2, respectively. The calculations show that Re0=1.6654>1. Moreover, since R0>1, we can obtain the existence of an endemic equilibrium for deterministic system (2.1). Figure 2 illustrates that the stochastic system (2.5) admits a unique global positive solution, which follows a unique ESD π(⋅), signifying the persistence of the disease in a community. On the other hand, we notice a slight difference between the stochastic solution and the solution without randomness. Meanwhile, Figures 3 and 4 depict the corresponding fitting density function and frequency histogram, as well as the movement of Markov chain, respectively.
Example 5.2. Let the noises (σij(1),σij(2))=(2×10−1,6×10−1) for each i=1,2,3,4; j=1,2. By means of direct derivation, we get Re0=0.4512<1. Therefore, we cannot determine the existence of an ESD of system (2.5). As shown in Figure 5, we can obtain that disease of system (2.5) will become extinct in the long run. Figure 6 shows the corresponding fitted density functions.
Example 5.3. To demonstrate the influence of stochastic perturbations, we have selected three distinct noise levels σ1=(σij(1),σij(2))=(0.002,0.006), σ2=(σij(1),σij(2))=(0.040,0.080), σ3=(σij(1),σij(2))=(0.300,0.500), wherein i=1,2,3,4,j=1,2. Figure 7 shows the trajectories of S(t),I(t),W(t),R(t), respectively. It is evident that the values of S(t),I(t), W(t),R(t) exhibit significant volatility as the noise intensity increases. Moreover, combining Figures 2, 5, and 7, we can observe that larger perturbations will drive the disease to extinction, whereas smaller perturbations may lead to its persistence within the community.
Example 5.4. (Sensitivity analysis) Under the conditions of satisfying 3.1 and 4.1, we can find that the basic reproduction number R0 is a crucial value for the eradication and persistence of pertussis disease. Therefore, Figure 8 shows sensitivity analysis on the reproductive number R0 of the disease. It is evident that transmission rates β and loss of immunity rate η are positively correlated with R0. Meanwhile, there is a negative correlation between recovery rate γ, vaccination probabilities ν, death rate μ, and R0. Physical exercise to enhance immunity and increase vaccination rates are beneficial for controlling the spread of pertussis.
In this paper, we primarily focus on the dynamics of pertussis in nature. By integrating natural immune enhancement and long-term immune waning, a stochastic pertussis system with Markov switching (2.5) is constructed and analyzed. Following the ideas of papers [24,25], the existence and uniqueness of the global positive solution are proven. Next, the exponential extinction of pertussis is shown in Theorem 3.1. Moreover, inspired by [25], θ−stochastic Lyapunov functions are used to establish the sufficient condition Re0>1, which guarantees the existence and uniqueness of an ESD for system (2.5) in Theorem 4.1. Finally, several examples are provided to confirm our theoretical findings (Figures 2–7). Based on these examples, we infer that stochastic noise is negatively correlated with the persistence of the disease. Furthermore, some questions deserve further study, for example, considering a stochastic pertussis model with semi-Markovian switching and impulses. The density function of a stochastic delay pertussis model with semi-Markov switching and Ornstein-Uhlenbeck process [41] may also be analyzed.
The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.
The research was supported by the National Natural Science Foundation of China (Nos.12161068, 12261069, 12301630).
The authors declare there is no conflicts of interest.
[1] |
G. Tang, J. Liu, Z. She, J. Li, Analysis and optimization on front bumper crash structural performance of a blade electric vehicle, Automob. Technol., 47 (2022), 64–67. https://doi.org/10.16638/j.cnki.1671-7988.2022.007.013 doi: 10.16638/j.cnki.1671-7988.2022.007.013
![]() |
[2] |
S. Wang, P. Xu, D. Wang, W. Guo, Q. Che, Inversion prediction method of subway composite energy absorption structure parameters, J. Railway Sci. Eng., 19 (2022), 3087–3095. https://doi.org/10.19713/j.cnki.43-1423/u.T20211358 doi: 10.19713/j.cnki.43-1423/u.T20211358
![]() |
[3] | L. Zhang, X. Zhao, L. Gan, H. Li, Y. Zheng, C. Zhou, Finite element simulation of ship collision accident inversion, Navig. China, 41 (2018), 78–83. |
[4] | M. Duan, Study on ship bridge collision test based on ship equivalent model, J. Chongqing Jiaotong Univ., 2017. |
[5] | Y. Liu, Research on Ice Load Inversion Method and Application of Polar Sailing Ships, M.A thesis, Harbin Engineering University, 2017. |
[6] | L. Han, The Research on Inversion Problems in Sheet Metal Forming Processes Based on Neural Network, M.A thesis, Hunan University, 2006. |
[7] | S. Yang, A Inversion Method of Material Parameters about the CFRP Used in Vehicles, M.A thesis, Hunan University, 2015. |
[8] | J. Yang, H. Xu, J. Zhao, J. Lei, Research on composite control of lateral position of lane line based on inversion and sliding mode, J. Highway Transp. Res. Dev., 32 (2015), 153–158. |
[9] |
G. E. Backus, J. F. Gilbert, Numerical applications of a formalism for geophysical inverse problems, Geophys. J. Int., 13 (1967), 247–276. https://doi.org/10.1111/j.1365-246X.1967.tb02159.x doi: 10.1111/j.1365-246X.1967.tb02159.x
![]() |
[10] |
G. E. Backus, J. F. Gilbert, The resolving power of gross earth data, Geophys. J. Int., 16 (1968), 169–205. https://doi.org/10.1111/j.1365-246X.1968.tb00216.x doi: 10.1111/j.1365-246X.1968.tb00216.x
![]() |
[11] |
G. E. Backus, J. F. Gilbert, Uniqueness in the inversion of inaccurate gross earth data, R. Soc., 266 (1970), 123–192. https://doi.org/10.1098/rsta.1970.0005 doi: 10.1098/rsta.1970.0005
![]() |
[12] | H. Gao, Research on Key Technologies of Vehicle Contact Collision Simulation, Ph.D thesis, Hunan University, 2008. |
[13] | M. Huang, Vehicle Crash Mechanics, New York, CRA Press, (2002), 69–138. https://doi.org/10.1201/9781420041866 |
[14] | S. Qiu, Vehicle Collision Safety Engineering, Beijing Institute of Technology Press, 2016. |
[15] |
X. Wang, W. Lin, G. Liu, S. Ma, R. Tong, Research progress of impact force based on Newton's impact recovery coefficient evaluation, Mech. Sci. Technol. Aerosp. Eng., 39 (2020), 1526–1533. https://doi.org/10.13433/j.cnki.1003-8728.20200179 doi: 10.13433/j.cnki.1003-8728.20200179
![]() |
1. | Mohammed Alsubhi, Rizwan Ahmed, Ibrahim Alraddadi, Faisal Alsharif, Muhammad Imran, Stability and bifurcation analysis of a discrete-time plant-herbivore model with harvesting effect, 2024, 9, 2473-6988, 20014, 10.3934/math.2024976 | |
2. | Anum Zehra, Parvaiz Ahmad Naik, Ali Hasan, Muhammad Farman, Kottakkaran Sooppy Nisar, Faryal Chaudhry, Zhengxin Huang, Physiological and chaos effect on dynamics of neurological disorder with memory effect of fractional operator: A mathematical study, 2024, 250, 01692607, 108190, 10.1016/j.cmpb.2024.108190 | |
3. | Aqeel Ahmad, Muhammad Farman, Parvaiz Ahmad Naik, Khurram Faiz, Abdul Ghaffar, Evren Hincal, Muhammad Umer Saleem, Analytical analysis and bifurcation of pine wilt dynamical transmission with host vector and nonlinear incidence using sustainable fractional approach, 2024, 11, 26668181, 100830, 10.1016/j.padiff.2024.100830 | |
4. | Parvaiz Ahmad Naik, Rizwan Ahmed, Aniqa Faizan, Theoretical and Numerical Bifurcation Analysis of a Discrete Predator–Prey System of Ricker Type with Weak Allee Effect, 2024, 23, 1575-5460, 10.1007/s12346-024-01124-7 | |
5. | Kottakkaran Sooppy Nisar, A constructive numerical approach to solve the Fractional Modified Camassa–Holm equation, 2024, 106, 11100168, 19, 10.1016/j.aej.2024.06.076 | |
6. | Chih-Wen Chang, Zohaib Ali Qureshi, Sania Qureshi, Asif Ali Shaikh, Muhammad Yaqoob Shahani, Real-Data-Based Study on Divorce Dynamics and Elimination Strategies Using Nonlinear Differential Equations, 2024, 12, 2227-7390, 2552, 10.3390/math12162552 | |
7. | Lana Abdelhaq, Sondos M. Syam, Muhammad I. Syam, An efficient numerical method for two-dimensional fractional integro-differential equations with modified Atangana–Baleanu fractional derivative using operational matrix approach, 2024, 11, 26668181, 100824, 10.1016/j.padiff.2024.100824 | |
8. | Parvaiz Ahmad Naik, Yashra Javaid, Rizwan Ahmed, Zohreh Eskandari, Abdul Hamid Ganie, Stability and bifurcation analysis of a population dynamic model with Allee effect via piecewise constant argument method, 2024, 70, 1598-5865, 4189, 10.1007/s12190-024-02119-y | |
9. | Mominul Islam, M. Ali Akbar, A study on the fractional-order COVID-19 SEIQR model and parameter analysis using homotopy perturbation method, 2024, 12, 26668181, 100960, 10.1016/j.padiff.2024.100960 | |
10. | Cahit Köme, Yasin Yazlik, Stability, bifurcation analysis and chaos control in a discrete predator–prey system incorporating prey immigration, 2024, 70, 1598-5865, 5213, 10.1007/s12190-024-02230-0 | |
11. | R. Prem Kumar, G.S. Mahapatra, P.K. Santra, Dynamical analysis of SARS-CoV-2-Dengue co-infection mathematical model with optimum control and sensitivity analyses, 2024, 80, 14681218, 104175, 10.1016/j.nonrwa.2024.104175 | |
12. | Parvaiz Ahmad Naik, Zohreh Eskandari, Mehmet Yavuz, Zhengxin Huang, Bifurcation results and chaos in a two-dimensional predator-prey model incorporating Holling-type response function on the predator, 2024, 0, 1937-1632, 0, 10.3934/dcdss.2024045 | |
13. | Parvaiz Ahmad Naik, Muhammad Farman, Anum Zehra, Kottakkaran Sooppy Nisar, Evren Hincal, Analysis and modeling with fractal-fractional operator for an epidemic model with reference to COVID-19 modeling, 2024, 10, 26668181, 100663, 10.1016/j.padiff.2024.100663 | |
14. | Hafiz Muhammad Shahbaz, Iftikhar Ahmad, Muhammad Asif Zahoor Raja, Hira Ilyas, Kottakkaran Sooppy Nisar, Muhammad Shoaib, A novel design of recurrent neural network to investigate the heat transmission of radiative Casson nanofluid flow consisting of carbon nanotubes (CNTs) across a curved stretchable surface, 2024, 104, 0044-2267, 10.1002/zamm.202400104 | |
15. | Ning Tian, Xiaoqi Liu, Rui Kang, Cheng Peng, Jiaxi Li, Shang Gao, Noise-to-State Stability of Random Coupled Kuramoto Oscillators via Feedback Control, 2024, 12, 2227-7390, 3715, 10.3390/math12233715 | |
16. | Yu Mu, Wing-Cheong Lo, Yuanshun Tan, Zijian Liu, Hybrid control for the prey in a spatial prey-predator model with cooperative hunting and fear effect time lag, 2025, 491, 00963003, 129217, 10.1016/j.amc.2024.129217 | |
17. | Allah Ditta, Parvaiz Ahmad Naik, Rizwan Ahmed, Zhengxin Huang, Exploring periodic behavior and dynamical analysis in a harvested discrete-time commensalism system, 2025, 13, 2195-268X, 10.1007/s40435-024-01551-z | |
18. | Subarna Roy, Subhas Khajanchi, Pankaj Kumar Tiwari, Fear and its carry-over effects in a generalist predator–prey system featuring cooperative hunting, 2025, 1598-5865, 10.1007/s12190-024-02346-3 | |
19. | Wei Li, Chunrui Zhang, Codimension-1 and Codimension-2 Bifurcations Analysis of Discrete Predator–Prey Model with Herd Behavior, 2025, 35, 0218-1274, 10.1142/S0218127425500105 | |
20. | Ibraheem M. Alsulami, Rizwan Ahmed, Faraha Ashraf, Exploring complex dynamics in a Ricker type predator–prey model with prey refuge, 2025, 35, 1054-1500, 10.1063/5.0232030 | |
21. | Parvaiz Ahmad Naik, Muhammad Farman, Muhammad Umer Saleem, Zhengxin Huang, Hijaz Ahmad, Muhammad Sultan, 2025, Chapter 11, 978-981-97-6793-9, 163, 10.1007/978-981-97-6794-6_11 | |
22. | Pinar Baydemir, Huseyin Merdan, Bifurcation analysis, chaos control, and FAST approach for the complex dynamics of a discrete-time predator–prey system with a weak Allee effect, 2025, 196, 09600779, 116317, 10.1016/j.chaos.2025.116317 | |
23. | Sangeeta Kumari, Sidharth Menon, Abhirami K, Dynamical system of quokka population depicting Fennecaphobia by Vulpes vulpes, 2025, 22, 1551-0018, 1342, 10.3934/mbe.2025050 | |
24. | Anuj Kumar Umrao, Prashant K. Srivastava, Exploring predator–prey dynamics: Integrating competitor predators, harvesting delay and fear effect on prey with a modified Beddington–DeAngelis functional response, 2025, 86, 14681218, 104391, 10.1016/j.nonrwa.2025.104391 | |
25. | Rabia Mehdi, Ranchao Wu, Zakia Hammouch, Chaotic dynamics and control in a discrete predation model with Holling II and prey refuge effects, 2025, 35, 1054-1500, 10.1063/5.0260237 |