
In this paper, we introduce a class of stochastic harvesting population system with Fractional Brownian Motion (FBM), which is still unclear when the stochastic noise has the character of memorability. Stochastic optimal control problems with FBM can not be studied using classical methods, because FBM is neither a Markov pocess nor a semi-martingale. When the external environment impact on the system of FBM, the necessary and sufficient conditions for the optimization are offered through the stochastic maximum principle, Hamilton function and Itˆo formula in our work. To illustrate our study, we provide an example to demonstrate the obtained theoretical results, which is the expansion of certainty population system.
Citation: Chaofeng Zhao, Zhibo Zhai, Qinghui Du. Optimal control of stochastic system with Fractional Brownian Motion[J]. Mathematical Biosciences and Engineering, 2021, 18(5): 5625-5634. doi: 10.3934/mbe.2021284
[1] | Bo Zheng, Lihong Chen, Qiwen Sun . Analyzing the control of dengue by releasing Wolbachia-infected male mosquitoes through a delay differential equation model. Mathematical Biosciences and Engineering, 2019, 16(5): 5531-5550. doi: 10.3934/mbe.2019275 |
[2] | Yunfeng Liu, Guowei Sun, Lin Wang, Zhiming Guo . Establishing Wolbachia in the wild mosquito population: The effects of wind and critical patch size. Mathematical Biosciences and Engineering, 2019, 16(5): 4399-4414. doi: 10.3934/mbe.2019219 |
[3] | Martin Strugarek, Nicolas Vauchelet, Jorge P. Zubelli . Quantifying the survival uncertainty of Wolbachia-infected mosquitoes in a spatial model. Mathematical Biosciences and Engineering, 2018, 15(4): 961-991. doi: 10.3934/mbe.2018043 |
[4] | Mugen Huang, Zifeng Wang, Zixin Nie . A stage structured model for mosquito suppression with immigration. Mathematical Biosciences and Engineering, 2024, 21(11): 7454-7479. doi: 10.3934/mbe.2024328 |
[5] | Bo Zheng, Wenliang Guo, Linchao Hu, Mugen Huang, Jianshe Yu . Complex wolbachia infection dynamics in mosquitoes with imperfect maternal transmission. Mathematical Biosciences and Engineering, 2018, 15(2): 523-541. doi: 10.3934/mbe.2018024 |
[6] | Rajivganthi Chinnathambi, Fathalla A. Rihan . Analysis and control of Aedes Aegypti mosquitoes using sterile-insect techniques with Wolbachia. Mathematical Biosciences and Engineering, 2022, 19(11): 11154-11171. doi: 10.3934/mbe.2022520 |
[7] | Diego Vicencio, Olga Vasilieva, Pedro Gajardo . Monotonicity properties arising in a simple model of Wolbachia invasion for wild mosquito populations. Mathematical Biosciences and Engineering, 2023, 20(1): 1148-1175. doi: 10.3934/mbe.2023053 |
[8] | Luis Almeida, Michel Duprez, Yannick Privat, Nicolas Vauchelet . Mosquito population control strategies for fighting against arboviruses. Mathematical Biosciences and Engineering, 2019, 16(6): 6274-6297. doi: 10.3934/mbe.2019313 |
[9] | Yun Li, Hongyong Zhao, Kai Wang . Dynamics of an impulsive reaction-diffusion mosquitoes model with multiple control measures. Mathematical Biosciences and Engineering, 2023, 20(1): 775-806. doi: 10.3934/mbe.2023036 |
[10] | Chen Liang, Hai-Feng Huo, Hong Xiang . Modelling mosquito population suppression based on competition system with strong and weak Allee effect. Mathematical Biosciences and Engineering, 2024, 21(4): 5227-5249. doi: 10.3934/mbe.2024231 |
In this paper, we introduce a class of stochastic harvesting population system with Fractional Brownian Motion (FBM), which is still unclear when the stochastic noise has the character of memorability. Stochastic optimal control problems with FBM can not be studied using classical methods, because FBM is neither a Markov pocess nor a semi-martingale. When the external environment impact on the system of FBM, the necessary and sufficient conditions for the optimization are offered through the stochastic maximum principle, Hamilton function and Itˆo formula in our work. To illustrate our study, we provide an example to demonstrate the obtained theoretical results, which is the expansion of certainty population system.
In Memory of Geoffrey J. Butler and Herbert I. Freedman
The invasive Asian tiger mosquito Aedes albopictus, originally indigenous to Southeast Asia, has invaded most continents, including Africa, Europe, and the Americas, prompted by international travel and trade, especially the used tire trade [1,2]. As a competent vector of more than 25 viruses such as dengue, Zika, and chikungunya, Aedes albopictus is a species of great medical concern in the world. It is the sole transmission vector of dengue in southern China. In 2014, an unprecedented outbreak of dengue fever hit Guangzhou, the capital city of Guangdong province, with 37,354 laboratory confirmed cases of infection [3].
With no effective therapies or licensed vaccines available, the dominating dengue control strategy has been vector elimination, including community-based source reduction and insecticide spraying. However, the invasiveness of Aedes albopictus and the creation of ubiquitous larval sources make source reduction a challenging task, while heavy insecticide applications induce serious environmental pollution and the insecticide resistance [1,2,4]. Fortunately, novel disease control methods using the endosymbiotic bacterium Wolbachia have been developed since the pioneering study of Xi et al. [5] in 2005. The success of the methods can be attributed mostly to the following facts: Biological safety : Wolbachia are naturally presented in up to 60% of insects, but are not usually found in the Aedes aegypti mosquito that transmits human viruses [6]. Virus control : The infection of some Wolbachia strains such as wMel has shown to block the transmission of dengue and other viruses in both Aedes aegypti and Aedes albopictus [6,7,8,9]. Maternal transmission : The bacterium is transmitted from infected females to the next generation. Cytoplasmic incompatibility (CI): Ⅰf a Wolbachia strain is infected by a male mosquito, but not by a female, then their crossing is incompatible that induces zygotic death and female sterility [5,6,8,9].
These tantalizing properties have engineered two Wolbachia-driven approaches for the elimination of mosquito borne diseases: population replacement and population suppression of wild Aedes mosquitoes. In the first approach, both male and female mosquitoes infected by a Wolbachia strain are released in natural areas. If the released numbers exceed a threshold level, then the reproduction advantage of the infected females driven by the CI mechanism and maternal transmission can facilitate the spread and fixation of Wolbachia in the wild mosquito population [10]. In the second approach, only infected males are released in natural areas. The incompatible crossing between the released males and the wild females induces female sterility and a suppression of the next generation [8,11]. These developments have stimulated extensive research activities in mathematical ecology to study the Wolbachia interfered population dynamics of mosquitoes [10,11,12,13,14]. In recent years, we have estimated the threshold release level for Wolbachia fixation and quantified how this critical level is affected by various factors such as the spatial movement of mosquitoes [15,16], the randomness of climatic conditions [17,18,19], and the leakage of maternal transmission [20,21,22]. More recently, we have also assessed the sensitivity of system parameters on the efficacy of Wolbachia driven Aedes mosquito suppression [23,24,25].
Let μ≥0 denote the mating competitiveness of the male mosquitoes infected by a Wolbachia strain comparing to wild males in the competition for wild female mating. Let ξ∈[0,1] denote the CI intensity – the zygotic death rate from the incompatible crossing of infected males and wild females. In the modeling of Wolbachia interfered mosquito dynamics so far, it has been almost always assumed that released males are equally competitive as wild males with μ=1, and the CI is complete with ξ=1. These assumptions are indeed strongly supported by the experiments in laboratory or even in the semi-field cages [8,11]. However, the field trials of Aedes albopictus population suppression in Guangzhou since 2015 have revealed a significant reduction in the competitiveness with μ ranging from 0.50 to 0.75 in wild areas [8,11,26]. In fact, wild Aedes albopictus mosquitoes in Guangzhou are naturally infected with two Wolbachia strains wAlbA and wAlbB. It is the infection of the third type of Wolbachia strain, wPip, established by embryonic micro-injections, that blocks virus transmission and induces CI in the crossing of triple-infected males and double-infected females [8]. Although the CI intensity is found to be strong with ξ≥0.95, it does indicate a possibility that CI may not be complete.
In this work, we initiate a preliminary assessment on how the reduced mating competitiveness and incomplete CI impair the efficacy of the Wolbachia-driven population suppression of Aedes mosquitoes in natural areas. This is of primary importance for designing more effective mosquito release policies in the planned large scale mosquito suppression campaign. We consider a population of wild adult mosquitoes in a total number A(t) at time t in the studying area, evenly divided in sex, that are interfered by R(t) released male adults infected by a Wolbachia strain not presented in the natural population. If the released males are equally competitive as wild males in mating, then the chance for a female to mate a released male, or the incompatible crossing probability, equals R(t) over the total number A(t)/2+R(t) of all males. In general, as the contribution of released males is scaled by the mating competitiveness μ, the incompatible crossing probability becomes μR(t)/[A(t)/2+μR(t)]. The offsprings produced by a female consist of two parts: those from compatible crossings, and those from incompatible crossings reduced by CI. Assume that the waiting time from the mating to the eclosion of next generation takes τ days, and each female produces b adult offsprings per day on average in compatible crossings. The delay τ changes with the climate and nutrient conditions and varies from 16 days to 66 days for Aedes albopictus in Guangzhou [27,28,29]. By taking account of the compatible and incompatible crossings and the CI intensity into consideration, we obtain the expected number of adults on day t produced by a single female on day t−τ as
b⋅A(t−τ)/2A(t−τ)/2+μR(t−τ)+b(1−ξ)⋅μR(t−τ)A(t−τ)/2+μR(t−τ). | (1.1) |
By multiplying theses terms with the total number of females, A(t−τ)/2, on day t−τ, we obtain the eclosion rate on day t. For the mortality terms, we follow the idea of Herb Freedman in the description of the classical logistic model [30]: For small A(t) the population decays linearly with a minimum mortality rate m>0; for large A(t) the mosquitoes ``compete each other for the limited resources" in the breeding sites in the larval stage and the decay is dominated by a second order term. In summary, we set
dA(t)dt=b2A(t−τ)+2(1−ξ)μR(t−τ)A(t−τ)+2μR(t−τ)A(t−τ)−m(1+A(t)K)A(t). | (1.2) |
The constant K plays the role that characterizes the density restriction as in the logistic model, but is not the carrying capacity due to the complexity of the eclosion terms and the factor m. If we set r(t)=2R(t)/A(t), the ratio of the released male numbers over the wild male numbers, then (1.2) is converted to
dA(t)dt=b21+μ(1−ξ)r(t−τ)1+μr(t−τ)A(t−τ)−m(1+A(t)K)A(t). | (1.3) |
In accordance with their biological meanings, we maintain the following conditions for the system parameters:
b,τ,K>0;m∈(0,1);μ≥0;ξ∈[0,1]. | (1.4) |
The solution of (1.3) subject to the initial condition A(t)=ϕ(t)∈C([t0−τ,t0],(0,∞)) for t∈[t0−τ,t0], t0≥0, will be denoted by A(t)=A(t,t0,ϕ).
We assess the impact of the mating competitiveness μ and the CI intensity ξ on the suppression efficacy of wild Aedes mosquitoes, by analyzing the global dynamics of (1.3) and performing numerical simulations with the experimental data. Our analysis identifies a threshold CI intensity ξ0∈(0,1) that increases in the reproduction of wild mosquitoes. When ξ≤ξ0, the population suppression is improbable no matter how many infected males are released. For ξ>ξ0, we find a threshold release ratio r∗>0: Ⅰf r(t) is kept at a low level with ˉr=sup[0,∞)r(t)<r∗, then A(t) is bounded below by a positive constant depending on ˉr; if r(t)≥r∗ for large t, then limt→∞A(t)=0. Our theorems are inclined to support a more important role of ξ than μ in the population suppression. It is further supported by our numerical examples which show that a slight decrease of ξ from 1 to 0.92 is more devastating than halving μ from 1 to 0.5. The CI intensity ξ0 not only warns an absolute failure of the population suppression when ξ<ξ0, it also predicts a great challenge of mosquito control when ξ is not considerably higher than ξ0. Finally, we estimate the optimal starting date for infected male release to target a more than 95% reduction of the wild population on October 1, normally in the middle of the peak season for mosquito growth and dengue fever transmission in Guangzhou. We find that the optimal date is almost independent of μ but is sensitive to ξ. If CI is complete, then starting about two months ahead can be an optimal option for less financial and labor costs. A slight reduction in the CI intensity may require not only more male releases, but also a considerably earlier starting date.
We estimate the thresholds for population suppression of wild Aedes mosquitoes by studying the global dynamics of (1.3). We begin with the following simple and fundamental result.
Lemma 2.1. For each ϕ∈C([t0−τ,t0],(0,∞)), the solution A(t)=A(t,t0,ϕ) of (1.3) is positive and bounded in t>t0.
Proof. Suppose for contradiction that A(t)>0 does not hold for all t>t0. Then we may find a t1>t0 such that A(t)>0 in [t0−τ,t1) and A(t1)=0. As A vanishes the first time at t1, we have A′(t1)≤0, which contradicts
A′(t1)=b21+μ(1−ξ)r(t1−τ)1+μr(t1−τ)A(t1−τ)>0 |
obtained from (1.3). If A(t) is unbounded, then there exists a sequence {tn}, with tn→∞ as n→∞, such that A′(tn)≥0, A(tn−τ)≤A(tn), A(tn)→∞ as n→∞. Substituting t=tn into (1.3) and applying these conditions yield
m(1+A(tn)K)≤b21+μ(1−ξ)r(tn−τ)1+μr(tn−τ)≤supx≥0b21+μ(1−ξ)x1+μx, |
which is apparently incompatible with the unboundedness of A.
As b counts the average number of adult offspring produced by a single female per day in compatible crossings, while m∈(0,1) is the average daily mortality rate of adults, it holds in normal environmental conditions that b>2m, that is,
ξ0=1−2mb>0. | (2.1) |
As in [25], we will also maintain this condition in the following studies. Interestingly, our next result shows that ξ0 defines a lower bound on the CI intensity below which population suppression is, independent of releasing efforts, absolutely improbable.
Theorem 2.1. If the CI intensity ξ<ξ0, then it holds uniformly that
A_=lim inft→∞A(t)≥bK(ξ0−ξ)2m>0, | (2.2) |
independent of the initial function ϕ∈C([t0−τ,t0],(0,∞)) and the releasing amount R(t).
Proof. We first prove A_>0. Suppose for contradiction that it is not true. Then A_=0. For an arbitrary initial function ϕ∈C([t0−τ,t0],(0,∞)), let ϕ0 denote its absolute minimum value on [t0−τ,t0]. For each n=1,2,⋯, define tn as the least time at which A reaches ϕ0/(n+1). By the assumption A_=0, tn is well-defined, and
A(t)>A(tn)=ϕ0n+1fort∈[t0−τ,tn),A′(tn)≤0,limn→∞tn=∞. | (2.3) |
Let t=tn in (1.3). We then have
m(1+A(tn)K)A(tn)≥b21+μ(1−ξ)r(tn−τ)1+μr(tn−τ)A(tn−τ). |
As A(tn−τ)>A(tn), the inequality remains valid after replacing A(tn−τ) by A(tn). Reorganizing terms slightly leads to
1+A(tn)K≥b2m1+μ(1−ξ)r(tn−τ)1+μr(tn−τ). |
Because [1+μ(1−ξ)x]/(1+μx) decreases in x≥0, by using (2.1) we find further that
1+A(tn)K>b(1−ξ)2m⇒A(tn)>K[b(1−ξ)2m−1]=bK(ξ0−ξ)2m. |
This is obviously inconsistent with (2.3) when n is sufficiently large and verifies A_>0.
By using the fluctuation lemma, see Lemma A.1 in [31], we can find a sequence {sn} such that sn→∞, A(sn)→A_ and A′(sn)→0 as n→∞. Let t=sn in (1.3). We obtain
m(1+A(sn)K)A(sn)≥b(1−ξ)2A(sn−τ)−A′(sn) |
whose infimum limit with n→∞ yields
m(1+A_K)A_≥b(1−ξ)2A_. |
from which (2.2) follows at once.
In the special case of ξ=ξ0, no meaningful conclusion can be made from Theorem 2.1 since the low bound in (2.2) equals zero. The following theorem shows that A(t) is still bounded below by a positive constant which depends on the upper bound of the releasing ratio r(t).
Theorem 2.2. Assume that the CI intensity ξ=ξ0 and ˉr=sup[t0−τ,∞)r(t)<∞. Then
A_=lim inft→∞A(t)≥bKξ02m(1+μˉr)>0 | (2.4) |
for any initial function ϕ∈C([t0−τ,t0],(0,∞)).
Proof. Similar to the proof of Theorem 2.1, we first prove A_>0. In fact, if A_=0, then there is an infinite series {sn} such that A(t)>A(sn) for t∈[t0−τ,sn), A′(sn)≤0 and A(sn)→0 as n→∞. From (1.3), we have
m(1+A(sn)K)A(sn)≥b21+μ(1−ξ0)r(sn−τ)1+μr(sn−τ)A(sn−τ). |
As A(sn−τ)>A(sn), by using (2.1) we obtain
1+A(sn)K>b2m1+2mμbr(sn−τ)1+μr(sn−τ)=b2m+μr(sn−τ)1+μr(sn−τ). |
Hence
A(sn)K>b2m−11+μr(sn−τ)=bξ02m[1+μr(sn−τ)]≥bξ02m(1+μˉr) |
which clearly contradicts the assumption A(sn)→0 and confirms A_>0.
The remaining proof uses the same idea as in the proof of Theorem 2.1 based on the fluctuation lemma with a slightly more complicated calculation. Let {sn} be an infinite sequence with A(sn)→A_ and A′(sn)→0 as n→∞. As [1+μ(1−ξ)x]/(1+μx) decreases in x≥0, substituting t=sn in (1.3) leads to
m(1+A(sn)K)A(sn)≥b21+μ(1−ξ0)ˉr1+μˉrA(sn−τ)−A′(sn). |
By taking the infimum limit, we derive
m(1+A_K)A_≥b21+μ(1−ξ0)ˉr1+μˉrA_ |
and so
1+A_K≥b2m1+μ(1−ξ0)ˉr1+μˉr=b2m+μˉr1+μˉr. |
It follows that
A_≥K(b2m+μˉr1+μˉr−1)=Kb2m−11+μˉr=bKξ02m(1+μˉr). |
This completes the proof.
Suppose that one of the following three extreme conditions holds:
(1) ξ=0 — Wolbachia infection does not modify the reproduction of wild Aedes mosquitoes at all.
(2) μ=0 — Wild female mosquitoes refuse to mate with released males completely.
(3) R(t)≡0 — No males carrying a novel Wolbachia strain are released.
Then System (1.3) reduces to
dA(t)dt=b2A(t−τ)−m(1+A(t)K)A(t), | (2.5) |
which describes the dynamics of wild mosquito population without Wolbachia intervention. Theorem 2.1 applies to (2.5) with ξ=0. Let A∗0 denote the constant in (2.2) in this case:
A∗0=bKξ02m=(b2m−1)K. | (2.6) |
Then Theorem 2.1 gives A_≥A∗0. Furthermore, by applying the proof of Theorem 2.1 in the last paragraph to the upper limit ¯A=lim supt→∞A(t), we can also prove ¯A≤A∗0. It follows that ¯A=A_=A∗0, and so A(t)≡A∗0 is globally asymptotically stable as shown by Yu [25]. It indicates that A∗0 defines the carrying capacity of the wild adult mosquito population.
In view of Theorem 2.1 and the discussion following its proof, it is natural to assume
μ>0,ξ>ξ0,r(t)>0fort≥t0−τ, | (2.7) |
in our search for the conditions ensuring the population suppression. When (2.7) holds,
r∗=b−2mμ[2m−b(1−ξ)]=ξ0μ(ξ−ξ0)>0. | (2.8) |
We show that r∗ defines a lower bound for the releasing efforts, in the sense that if the supremum of r(t) in t≥t0−τ is below r∗, then A(t) is bounded below by a positive constant. Corresponding to a constant release ratio r(t)≡r, (1.3) has an equilibrium point
A∗r=(b2m⋅1+μ(1−ξ)r1+μr−1)K, | (2.9) |
besides A≡0. By comparing (2.6) and (2.9) we see that A∗r reduces to A∗0 when r=0, which confirms the consistency in the definitions (2.6) and (2.9). It is easy to check that A∗r decreases in r≥0, which is positive for r∈[0,r∗), vanishes uniquely at r=r∗, and becomes negative for r>r∗.
Theorem 2.3. Let (2.7) hold. Suppose that ¯r=sup[t0−τ,∞)r(t)<r∗. Then for each ϕ∈C([t0−τ,t0],(0,∞)), the solution A(t)=A(t,t0,ϕ) of (1.3) satisfies
0<A∗¯r≤lim inft→∞A(t)≤lim supt→∞A(t)≤A∗r_, | (2.10) |
where r_=inf[t0−τ,∞)r(t).
Proof. From the assumption ¯r<r∗ and the discussion on A∗r before the statement of this theorem, we see that A∗¯r>0. For an arbitrary positive number a1<A∗¯r and a1<ϕ(t) on [t0−τ,t0], we claim A(t)>a1 for all t≥t0. Indeed, if this is not true, then we may let t1>t0 be the least time at which A=a1. Hence A(t)>a1 in [t0,t1), A(t1)=a1, and A′(t1)≤0. Substituting t=t1 into (1.3) leads to
m(1+a1K)a1≥b21+μ(1−ξ)r(t1−τ)1+μr(t1−τ)A(t1−τ)>b21+μ(1−ξ)¯r1+μ¯ra1. |
Consequently, we derive
1+a1K>b2m⋅1+μ(1−ξ)¯r1+μ¯r⇒a1>(b2m⋅1+μ(1−ξ)¯r1+μ¯r−1)K=A∗¯r. |
This contradicts our assumption a1<A∗¯r and establishes the claim. As a result, it also holds that the lower limit A_≥a1.
We use the fluctuation lemma again to complete the proof for the first part of (2.10). Let {sn} be an increasing and divergent sequence such that A(sn)→A_ and A′(sn)→0 as n→∞. Replacing t by sn in (1.3) gives
m(1+A(sn)K)A(sn)≥b21+μ(1−ξ)¯r1+μ¯rA(sn−τ)−A′(sn). |
By taking the infimum limit, we obtain
m(1+A_K)A_≥b21+μ(1−ξ)¯r1+μ¯rA_, |
from which a simplification of terms gives A_≥A∗¯r.
As ¯A≥A_>0, we can use the fluctuation lemma directly to prove the second part of (2.10). Let {tn} be a divergent sequence such that A(tn)→¯A and A′(tn)→0 as n→∞. Let t=tn in (1.3). We have
m(1+A(tn)K)A(tn)≤b21+μ(1−ξ)r_1+μr_A(tn−τ)−A′(tn). |
By taking the supremum limit, we obtain
m(1+¯AK)¯A≤b21+μ(1−ξ)r_1+μr_¯A. |
With a slight simplification, we find ¯A≤A∗r_ by the definition (2.9) of A∗r.
When r(t)≡r for a constant r∈(0,r∗), an over-simplified assumption on the releasing efforts that the number of released males is proportional to the wild male number, we have ¯r=r_=r and Theorem 2.3 implies limt→∞A(t)=A∗r>0. It indicates that A∗¯r and A∗r_ probably provide the sharpest estimates for the lower and upper limits of A(t).
We now prove that r∗ sets a threshold level of infected male mosquito release for population suppression: As long as r(t)≥r∗, the wild population will be eliminated ultimately. It indicates further that ξ0 defines the threshold CI intensity over which a complete mosquito elimination is ascertained provided additionally that the release ratio r(t)≥r∗.
Theorem 2.4. Let (2.7) hold. Suppose there is T>0 such that r(t)≥r∗ for t≥T. Then for any ϕ∈C([t0−τ,t0],(0,∞)), limt→∞A(t,t0,ϕ)=0.
Proof. It suffices to show that the upper limit ¯A=0. By the fluctuation lemma, there exists a divergent sequence {tn} such that A(tn)→¯A and A′(tn)→0 as n→∞. When n is sufficiently large such that tn≥T+τ, (1.3) gives
m(1+A(tn)K)A(tn)=b21+μ(1−ξ)r(tn−τ)1+μr(tn−τ)A(tn−τ)−A′(tn)≤b21+μ(1−ξ)r∗1+μr∗A(tn−τ)−A′(tn). |
Taking limit gives
m(1+¯AK)¯A≤b21+μ(1−ξ)r∗1+μr∗lim supn→∞A(tn−τ)≤b21+μ(1−ξ)r∗1+μr∗¯A. |
The inequality holds when ¯A=0. If ¯A>0, then
m(1+¯AK)≤b21+μ(1−ξ)r∗1+μr∗⇒¯A≤(b2m1+μ(1−ξ)r∗1+μr∗−1)K=A∗r∗ |
by the definition of A∗r in (2.9). However, as we have noticed that A∗r∗=0 right before the statement of Theorem 2.3, the last inequality leads to ¯A≤0 and a contradiction. Therefore, it must hold that ¯A=0.
In accordance with the basic assumptions (1.4), (2.1), and (2.7), the following condition will be maintained in the discussion below:
b,τ,K>0;m∈(0,1);μ>0;ξ∈(0,1];ξ0=1−2mb>0. | (3.1) |
In Theorems 2.1 – 2.4, we have identified ξ0 as the threshold of the CI intensity, and
r∗=r∗(μ,ξ)=b−2mμ[2m−b(1−ξ)]=ξ0μ(ξ−ξ0)>0 | (3.2) |
as the threshold for r(t)=2R(t)/A(t) – the ratio between the abundance R(t) of released males infected by a novel Wolbachia strain absent in the natural population and the abundance A(t)/2 of wild male adults. As b/m can be interpreted as the net production rate of wild females, we see from (3.1) that the threshold CI intensity ξ0 is an increasing function of the natural reproduction rate of wild mosquitoes. If the infection does not induce complete CI and the CI intensity ξ<ξ0, then the population suppression is improbable no matter how many infected males are released in the wild area. If ξ>ξ0 and for some T>0, r(t)≥r∗ for t≥T, then the wild mosquito population will be eliminated ultimately. If ¯r=sup[t0−τ,∞)r(t)<r∗, then A(t) is bounded below by the constant A∗¯r>0 with
A∗¯r=(b2m⋅1+μ(1−ξ)¯r1+μ¯r−1)K |
as defined in (2.9).
Our theorems are inclined to support a more important role of the CI intensity than the mating competitiveness μ in the population suppression. In contrary to the existence of the threshold CI intensity, our model does not generate a threshold level of μ. Since μ and the release ratio r(t) appear together in two product forms of μr(t−τ), they relate reciprocally and a decrease in μ can be compensated by a propositional increase in r(t−τ). As shown in (3.2), the threshold release ratio r∗ is reciprocal to each of μ and ξ−ξ0. The CI intensity ξ0 not only warns an absolute failure of the population suppression when ξ<ξ0, it also predicts a great challenge of mosquito control when ξ is not considerably higher than ξ0.
For Aedes albopictus in Guangzhou, we have estimated b∈[0.9043,6.4594], τ∈[16,66], and m∈[0.0198,0.1368] in [25], by combining the laboratory [32,33,34,35] and field [27,29,34] data. To make our observations above more specific and transparent, we recall from our discussion on ξ and μ in the introduction and fix
b=2,μ∈[0.5,1],ξ∈[0.91,1],m=0.1,τ=19,K=20,000. | (3.3) |
The constant K does not alter the CI intensity threshold ξ0 or the release ratio threshold r∗ and shows only a minimal impact on the dynamical behavior of A(t). As it scales with the size of the studying area, we take K=20,000 as a representative example. Clearly, the parameters in (3.3) satisfy (3.1) and determine ξ0=0.9.
The dependence of r∗ on ξ∈[0.91,1] and μ∈[0.5,1] is shown in Figure 1A, and its dependence on ξ for fixed μ=0.5,0.75,1, is shown in Figure 1B. It decreases in both ξ and μ with the maximum value r∗(0.5,0.91)=180 and the least value r∗=r∗(1,1)=9. It indicates that for a complete elimination of wild Aedes mosquitoes, at least a ratio 9:1 of the released males over wild males needs to be maintained. Although it appears higher than the 5:1 ratio estimated in [25], the two different ratios do not contradict each other because the latter is estimated at the initial time with r(0)=5 in a constant releasing policy with R(t)≡R(0). As shown in Figure 1B, the threshold ratio r∗ is extremely sensitive to the CI intensity ξ when it is close to the threshold CI intensity level ξ0=0.9. Even with an equal mating competitiveness μ=1, a slight change from ξ=0.92 to ξ=0.91 doubles r∗ from r∗(1,0.92)=45 to r∗(1,0.91)=90. On the other hand, when the CI is close to be complete, the change of r∗ is relatively flat; e.g., 9≤r∗≤18 for ξ∈[0.95,1] and μ=1.
We use numerical examples to demonstrate further that a reduction in the CI intensity ξ causes substantially higher damages than a reduction in the mating competitiveness μ in the population suppression. We compare the dynamical behavior of A(t) with varying ξ or μ but fix the other parameter at 1. In the simulations, we use the same values for b, m, τ, and K as in (3.3). The environmental carrying capacity A∗0 determined by (2.6) equals 9K=180,000. For specificity, we also fix
ϕ(t)=20,000,t∈[−19,0],r(t)≡14. | (3.4) |
The initial population size in the area is 1/9 of the carrying capacity, while the releasing ratio is larger than the threshold value r∗(1,1)=9.
In Figure 2A, the mating competitiveness is fixed at μ=1 and the CI intensity decreases from 1 to 0.95 and 0.92. In agree with Theorem 2.4, A(t)→0 as t→∞ when ξ=1. When ξ decreases to 0.95 or 0.92, a complete population suppression is impossible and Theorem 2.3 implies that A(t)→A∗14≈2667 and A(t)→A∗14≈8267, respectively, as t→∞. Although the wild population is not completely eliminated in the latter two cases, it is suppressed by (180,000−2667)/180,000≈98.52% and (180,000−8267)/180,000≈95.41% comparing to the steady-state (or the carrying capacity) of the wild population without Wolbachia intervention. In Figure 2B, ξ=1 and μ decreases from 1 to 0.75 and 0.5. Interestingly, even when μ is decreased by 25% from 1 to 0.75, it still holds that A(t)→0 as t→∞. For μ=0.5, we have A(t)→A∗14=5000 as t→∞. As 8267>5000, we see that a slight decrease of the CI intensity from 1 to 0.92 is more devastating in the population suppression than halving the mating competitiveness from 1 to 0.5.
Aedes albopictus in Guangzhou overwinters as diapause eggs with few surviving adults in open areas from December to February. With the elevation of temperature and precipitation in the spring, most diapause eggs start hatching at the turn of February to March. The adult population grows rapidly from the middle of March and reaches the first peak of the year in late May or early June [28,29]. The population declines in the hot summer and bounces back to reach the second peak in September or October [26,4]. The high-incidence season of dengue fever overlaps the second peak of mosquito abundance at the turn of September to October [36]. Breaking down the second peak of Aedes mosquitoes provides a temporary and efficient control of dengue transmission. Suppose we target at a more than 95% reduction of the adult mosquitoes on October 1 by a Wolbachia driven approach. When shall the infected male release be started? We note that a perfect answer to this question may not be reached easily as it involves many factors such as weather conditions, financial resource, labors, and the community support. In the remaining discussion, we use a numerical simulation based on our model to give a partial answer.
In the simulation, we use again the values for b, m, τ, and K specified in (3.3), which determine the carrying capacity A∗0=180,000. Set t0=0 on March 1 before the burst of adult mosquito population and let ϕ(t)=10,000 on [−19,0]. With the starting date on the first day of each month from March to September, we estimate the minimum constant releasing ratio and the corresponding number of released males such that the number of the wild adults on October 1 (t≈210) over the carrying capacity 180,000 is less than 5%. For a given pair of μ>0 and ξ∈[0,1], we denote by r=rm(μ,ξ) the minimum constant release ratio, and N=N(μ,ξ) the total number of released males. For μ=0.75 and ξ=1, the temporal profiles of the wild adults interfered by infected male mosquitoes starting on different dates are shown in Figure 3A. We find rm=41.8 when the releasing starts on September 1, which reduces sharply to rm=16.6 when it starts one month early from August 1, and reduces moderately to rm=11.8 when it starts another month early on July 1. Very interestingly, we find that the ratio rm shows no significant reduction further when the starting date is shifted back to the first days of June, May, April, and March. These numbers seem to suggest that starting on August, about two months ahead of the target date, can be an optimal option for less financial and labor costs.
We examine further the variations of the minimum release ratio rm(μ,ξ) and the total number N(μ,ξ) on the starting dates for (μ,ξ)=(1,1),(0.75,1),(1,0.95). The values of rm are depicted in Figure 3B, C, and the values of N are depicted in Figure 3D; see also Table 1 where these numbers are listed. The variations of rm(1,1), N(1,1), and N(0.75,1) on the starting dates all follow the same pattern as we discussed above for (0.75,1): A large reduction is seen when the date is moved from September 1 back to August 1, and only a moderate reduction is shown if the dates are moved back further. In addition, a 25% reduction in the mating competitiveness from μ=1 to μ=0.75 requires only a moderate 30%−35% increase in the releasing amounts. In agree with our discussion in the previous section, a mild 5% reduction in the CI tensity from ξ=1 to ξ=0.95 causes more damages in the suppression than a 25% reduction in μ. The devastating impact is more obvious when the starting date is close to the target date October 1. For (μ,ξ)=(1,0.95), our simulation reveals that it is nearly impossible to meet the suppression goal within one month: For a starting date on September 1, we find that even for the releasing ratio as high as r=10000,16000,20000, the adult abundances on October 1 are still 7.39%, 7.33%, 7.20% over the carrying capacity. As shown in Table 1, rm(1,095)=39 when the release starts on August 1, and reduces sharply to 17.2 if it starts on July 1, but shows no significant reduction in earlier starting dates. These numbers suggest an optimal starting date in July for male mosquito releasing in this case, other than August in the other two cases. In summary, a loss in the CI intensity causes considerably more damages in the suppression of wild mosquito populations than the same magnitude of loss in the mating competitiveness. To reduce the adult population by more than 95% during the dengue fever peak season, an optimal starting date for infected male mosquito release is about three months ahead when the CI intensity is reduced by 5%. In a sharp contrast, when CI is complete, a two month implementation of the control measures is sufficient even if the mating competitiveness is reduced up to 25%.
Start Date | rm(1,1), N(1,1) | rm(0.75,1), N(0.75,1) | rm(1,0.95), N(1,0.95) |
Sept 1 | 32.2, 1.8409×107 | 41.8, 2.4195×107 | |
Aug 1 | 12.6, 1.2538×107 | 16.6, 1.6635×107 | 39, 3.8895×107 |
July 1 | 8.8, 1.1896×107 | 11.8, 1.6112×107 | 17.2, 2.3286×107 |
June 1 | 7.4, 1.1886×107 | 9.8, 1.5848×107 | 12.8, 2.0535×107 |
May 1 | 6.8, 1.1395×107 | 9, 1.5197×107 | 11, 1.8596×107 |
Apr 1 | 6.4, 0.9951×107 | 8.6, 1.3249×107 | 10, 1.5750×107 |
Mar 1 | 6, 0.6489×107 | 8, 0.8651×107 | 9.2, 0.9996×107 |
This work was supported by National Natural Science Foundation of China (11471085, 11631005), Program for Changjiang Scholars and Innovative Research Team in University (IRT_16R16), Natural Science Foundation of Guangdong Province (2017A030310597), and Guangdong University of Finance & Economics Big Data and Educational Statistics Application Laboratory(2017WSYS001).
All authors declare no conflicts of interest in this paper.
[1] |
A. Cadenilas, A Stochastic maximum principle for systems with jumps, with applications to finance, Syst. Control Lett., 47 (2002), 433–444. doi: 10.1016/S0167-6911(02)00231-1
![]() |
[2] |
Q. Zhang, J. Tian, X. Li, The asymptotic stability of numerical analysis for stochastic age-dependent cooperative Lotka-Volterra system, Math. Biosci. Eng., 18 (2021), 1425–1449. doi: 10.3934/mbe.2021074
![]() |
[3] | Q. Zhang, W. Liu, Z. Nie, Existence, uniqueness and exponential stability for stochastic age-dependent population, Appl. Math. Comput., 154 (2004), 183–201. |
[4] | Y. Zhao, S. Yuan, Q. Zhang, Numerical solution of a fuzzy stochastic single species age structure model in a polluted environmen, Appl. Math. Comput., 260 (2015), 385–396. |
[5] |
Y. Du, Q. Zhang, Anke Meyer Bases. The positive numerical solution for stochastic age-dependent capital system based on explicit-implicit algorithm, Appl. Numer. Math., 165 (2021), 198–215. doi: 10.1016/j.apnum.2021.02.015
![]() |
[6] | W. Li, M. Ye, Q. Zhang, Numerical approximation of a stochastic age-structured population model in a polluted environment with Markovian switching, Numer. Math. Part. D. E., 36 (2020), 22488. |
[7] | W. Liu, Q. Zhang, Convergence of numerical solutions to stochastic age-structured system of three species, Appl. Numer. Math., 218 (2011), 3973–3980. |
[8] | Q. Zhang, C. Han, Convergence of numerical solutions to stochastic age-structured population system, Appl. Math. Comput., 118 (2005), 134–146. |
[9] | S. Zhu, Convergence of the semi-implicit euler method for stochastic age-dependent population equations with poisson jumps, Int. J. Biomath., 34 (2009), 2034–2043. |
[10] | Y. Pei, H. Yang, Q. Zhang, Asymptotic mean-square boundedness of the numerical solutions of stochastic age-dependent population equations with Poisson jumps, Appl. Math. Comput., 320 (2018), 524–534. |
[11] | G. W. Weber, E. Savku, I. Baltas, Stochastic optimal control and games in a world of regime switches, paradigm shifts, jumps and delay, 17th Europt Workshop on Advances in Continuous Optimization, (2019). |
[12] | G. W. Weber, E. Savku, Y. Y. Okur, Optimal control of stochastic systems with regime switches, jumps and delay in finance, economics and nature, Conference: Seminar at Department of Mathematics, (2016). |
[13] |
W. Ma, Q. Zhang, C. Han, Numerical analysis for stochastic age-dependent population equations with fractional Brownian motion, Commun. Nonlinear Sci., 17 (2012), 1884–1893. doi: 10.1016/j.cnsns.2011.08.025
![]() |
[14] |
P. E. Kloeden, A. Neuenkirch, R. Pavani, Multilevel monte carlo for stochastic differential equations with additive fractional noise, Ann. Oper. Res., 189 (2011), 255–276. doi: 10.1007/s10479-009-0663-8
![]() |
[15] |
M. Giles, Multilevel monte carlo path simulation, Oper. Res., 56 (2008), 607–617. doi: 10.1287/opre.1070.0496
![]() |
[16] | T. E. Duncan, B. Maslowski, Pasik-Duncan B. Semilinear stochastic equations in a Hilbert space with a fractional Brownian motion, Siam J. Math. Anal., 40 (2008), 2286–2315. |
[17] |
W. Zhou, X. Zhou, J. Yang, Stability analysis and application for delayed neural networks driven by fractional Brownian noise, IEEE T. Neur. Net. Lear., 29 (2018), 1491–1502. doi: 10.1109/TNNLS.2017.2674692
![]() |
[18] | Z. Luo, Optimal harvesting control problem for an age-dependent competing system of n species, Appl. Math. Comput., 183 (2006), 119–127. |
[19] |
C. Zhao, M. Wang, Z. Ping, Optimal control of harvesting for age-dependent predator-prey system, Math. Comput. Model., 42 (2005), 573-584. doi: 10.1016/j.mcm.2004.07.019
![]() |
[20] | J. Chen, Z. He, Optimal control for a class of nonlinear age-distributed population systems, Appl. Math. Comput., 214 (2009), 574–580. |
[21] |
Z. He, D. Ni, S. Wang, Optimal harvesting of a hierarchical age-structured population system, Int. J. Biomath., 12 (2019), 1950091. doi: 10.1142/S1793524519500918
![]() |
[22] | R. Dhayal, M. Malik, S. Abbas, Optimal controls for second-order stochastic differential equations driven by mixed fractional Brownian motion with impulses, Math. Method Appl. Sci., 43 (2020), 4107–4124. |
[23] |
S. Adly, A. Hantoute, M. Théra, Nonsmooth Lyapunov pairs for infinite-dimensional first-order differential inclusions, Nonlinear Anal., 75 (2012), 985–1008. doi: 10.1016/j.na.2010.11.009
![]() |
[24] | X. R. Mao, Stochastic differential equations and applications, Horwood, UK, (2007). |
[25] | D. Bainov, P. Simeonov, Integral inequalities and applications, Kluwer Academic Publishers, (1992). |
[26] | Q. Lü, M. L. Deng, W. Q. Zhu, Stochastic averaging of quasi integrable and resonant Hamiltonian systems excited by fractional Gaussian noise with Hurst index 1/2<H<1, Acta Mech. Solida Sin., 1 (2017), 11–19. |
1. | Yantao Shi, Jianshe Yu, Wolbachia infection enhancing and decaying domains in mosquito population based on discrete models, 2020, 14, 1751-3758, 679, 10.1080/17513758.2020.1805035 | |
2. | A.S. Benedito, C.P. Ferreira, M. Adimy, Malay Banerjee, Modeling the dynamics of Wolbachia-infected and uninfected Aedes aegypti populations by delay differential equations, 2020, 15, 0973-5348, 76, 10.1051/mmnp/2020041 | |
3. | Mugen Huang, Linchao Hu, Modeling the suppression dynamics of Aedes mosquitoes with mating inhomogeneity, 2020, 14, 1751-3758, 656, 10.1080/17513758.2020.1799083 | |
4. | Genghong Lin, Yuanxian Hui, Stability analysis in a mosquito population suppression model, 2020, 14, 1751-3758, 578, 10.1080/17513758.2020.1792565 | |
5. | Yijie Li, Zhiming Guo, Yanyuan Xing, ModelingWolbachiaDiffusion in Mosquito Populations by Discrete Competition Model, 2020, 2020, 1026-0226, 1, 10.1155/2020/8987490 | |
6. | Jianshe Yu, Existence and stability of a unique and exact two periodic orbits for an interactive wild and sterile mosquito model, 2020, 269, 00220396, 10395, 10.1016/j.jde.2020.07.019 | |
7. | Bo Zheng, Jia Li, Jianshe Yu, Existence and stability of periodic solutions in a mosquito population suppression model with time delay, 2022, 315, 00220396, 159, 10.1016/j.jde.2022.01.036 | |
8. | Yun Li, Hongyong Zhao, Kai Wang, Dynamics of an impulsive reaction-diffusion mosquitoes model with multiple control measures, 2022, 20, 1551-0018, 775, 10.3934/mbe.2023036 | |
9. | Jianshe Yu, Jia Li, A delay suppression model with sterile mosquitoes release period equal to wild larvae maturation period, 2022, 84, 0303-6812, 10.1007/s00285-022-01718-2 | |
10. | Mu-gen Huang, Jian-she Yu, Global Asymptotic Stability in a Delay Differential Equation Model for Mosquito Population Suppression, 2022, 38, 0168-9673, 882, 10.1007/s10255-022-1021-8 | |
11. | Yantao Shi, Bo Zheng, Modeling Wolbachia infection frequency in mosquito populations via a continuous periodic switching model, 2023, 12, 2191-950X, 10.1515/anona-2022-0297 | |
12. | Yantao Shi, Bo Zheng, Discrete dynamical models on Wolbachia infection frequency in mosquito populations with biased release ratios, 2022, 16, 1751-3758, 320, 10.1080/17513758.2021.1977400 | |
13. | Mugen Huang, Zifeng Wang, Zixin Nie, A stage structured model for mosquito suppression with immigration, 2024, 21, 1551-0018, 7454, 10.3934/mbe.2024328 | |
14. | Hui Wan, Yin Wu, Guihong Fan, Dan Li, Wolbachia invasion dynamics of a random mosquito population model with imperfect maternal transmission and incomplete CI, 2024, 88, 0303-6812, 10.1007/s00285-024-02094-9 | |
15. | Sourav Kumar Sasmal, Yasuhiro Takeuchi, Yukihiko Nakata, A simple model to control the wild mosquito with sterile release, 2024, 531, 0022247X, 127828, 10.1016/j.jmaa.2023.127828 | |
16. | Mugen Huang, Jianshe Yu, Modeling the Impact of Migration on Mosquito Population Suppression, 2023, 22, 1575-5460, 10.1007/s12346-023-00834-8 | |
17. | Xiaoke Ma, Ying Su, Wolbachia Invasion in Mosquitoes with Incomplete CI, Imperfect Maternal Transmission and Maturation Delay, 2024, 86, 0092-8240, 10.1007/s11538-024-01363-4 | |
18. | Hongpeng Guo, Accelerating mosquito population replacement: a two-dimensional model with periodic release of Wolbachia-infected males, 2025, 2025, 2731-4235, 10.1186/s13662-025-03907-x | |
19. | Xinyu Wang, Liping Wang, Ruizhe Shang, Peng Wu, Dynamic analysis of a mathematical model for wild mosquito population control: combining incompatible and sterile insect techniques, 2025, 140, 2190-5444, 10.1140/epjp/s13360-025-06125-2 | |
20. | Mu-gen Huang, Jian-she Yu, Global asymptotic stability in a delay stage structured model for mosquito population suppression, 2025, 40, 1005-1031, 122, 10.1007/s11766-025-4396-5 |
Start Date | rm(1,1), N(1,1) | rm(0.75,1), N(0.75,1) | rm(1,0.95), N(1,0.95) |
Sept 1 | 32.2, 1.8409×107 | 41.8, 2.4195×107 | |
Aug 1 | 12.6, 1.2538×107 | 16.6, 1.6635×107 | 39, 3.8895×107 |
July 1 | 8.8, 1.1896×107 | 11.8, 1.6112×107 | 17.2, 2.3286×107 |
June 1 | 7.4, 1.1886×107 | 9.8, 1.5848×107 | 12.8, 2.0535×107 |
May 1 | 6.8, 1.1395×107 | 9, 1.5197×107 | 11, 1.8596×107 |
Apr 1 | 6.4, 0.9951×107 | 8.6, 1.3249×107 | 10, 1.5750×107 |
Mar 1 | 6, 0.6489×107 | 8, 0.8651×107 | 9.2, 0.9996×107 |
Start Date | rm(1,1), N(1,1) | rm(0.75,1), N(0.75,1) | rm(1,0.95), N(1,0.95) |
Sept 1 | 32.2, 1.8409×107 | 41.8, 2.4195×107 | |
Aug 1 | 12.6, 1.2538×107 | 16.6, 1.6635×107 | 39, 3.8895×107 |
July 1 | 8.8, 1.1896×107 | 11.8, 1.6112×107 | 17.2, 2.3286×107 |
June 1 | 7.4, 1.1886×107 | 9.8, 1.5848×107 | 12.8, 2.0535×107 |
May 1 | 6.8, 1.1395×107 | 9, 1.5197×107 | 11, 1.8596×107 |
Apr 1 | 6.4, 0.9951×107 | 8.6, 1.3249×107 | 10, 1.5750×107 |
Mar 1 | 6, 0.6489×107 | 8, 0.8651×107 | 9.2, 0.9996×107 |