1. Introduction
The global incidence of dengue is placing more than 50% of the world’s population at risk, and affects more than 100 countries around the world [
1]. As a mosquito-borne disease, dengue viruses are transmitted through the bites of infected
Aedes female mosquitoes, including
Aedes aegypti and
Aedes albopictus, which act as the main vectors of dengue. There are four serotypes of dengue viruses. As the repeated proof of the antibody-dependent enhancement phenomenon between different serotypes, the application of the dengue vaccine has been halted [
2,
3]. Furthermore, no effective drugs are available for dengue viruses, especially for severe cases, including dengue hemorrhagic fever and dengue shock syndrome. Hence, the most straightforward way to control dengue is to eradicate mosquitoes by spraying insecticide or removing mosquito breeding sites, but neither method works effectively due to the appearance of insecticide resistance and difficulty in eradicating breeding sites [
4].
Luckily, the World Mosquito Program’s
Wolbachia release method is helping communities worldwide to control and prevent the spread of mosquito-borne diseases. The endosymbiotic bacteria
Wolbachia occurs naturally in 60% of insect species, which are safe for humans and the environment. As the successful establishment of
Wolbachia in
Aedes mosquitoes through microinjection [
5,
6],
Wolbachia-driven mosquito control techniques have become one of the effective biological control methods. The
Wolbachia bacteria is maternally inherited, which can induce a special mechanism called “cytoplasmic incompatibility”(CI), resulting in early developmental arrest of the embryo when infected male mosquitoes mate with uninfected female mosquitoes [
5,
6,
7,
8,
9]. By using CI, and the maternal transmission of
Wolbachia, one of the World Mosquito Program’s
Wolbachia release methods aims to replace the wild vector mosquitoes with
Wolbachia infected ones, whose capability of transmitting diseases has been greatly reduced owing to the
Wolbachia infection [
5,
6,
7,
8,
9]. Taking the release in Australia as an example, back to 2011,
Wolbachia-infected mosquitoes were released to replace the natural
Aedes populations in Cairns, reaching a near-fixation of
Wolbachia in a few months [
9,
10]. Similar field trials have also been carried out in Brazil, Vietnam, Indonesia, etc. Post surveillance has shown that in areas where
Wolbachia is self-sustaining at a high level, there have been no dengue outbreaks.
Prior to area-wide releases, cage experiments are essential for designing effective and affordable release strategies. Motivated by the experimental observation in [
7] that supplemental male release in every generation can accelerate the
Wolbachia invasion, we develop a general discrete model which is based on the classic models in [
11,
12,
13,
14,
15,
16], which focused on
Wolbachia spread dynamics in cage mosquito populations without supplemental releases. All these discrete models only include a one-time release of both
Wolbachia-infected females and males, which generated an unstable equilibrium, below which the
Wolbachia infections will disappear, and above which,
Wolbachia invasion is successful.
The generalized model developed in
Section 2 covers the models in [
11,
12,
13,
14,
15,
16] as special cases. Particularly, we focus on a periodic and impulsive release strategy, where
Wolbachia-infected males are released with the release ratio
during the first
N generations, and the release ratio from generation
to generation
T is
. This leads to a discrete switching model. Very recently, this release strategy was investigated thoroughly, in [
17,
18,
19], by using ordinary differential equation models, leaving the study of discrete models, which are well suited for cage experiment absence. In
Section 3, sufficient conditions to guarantee the existence of at least one
T-periodic solution are obtained for the discrete switching model. Although our numerical example shows that the
T-periodic solution is also unique when it exists, we failed to prove it theoretically for general
N and
T, and methods used in [
17,
18,
19] are not applicable for discrete models. However, we only succeed to prove that this numerical observation holds for the special case
and
. The whole paper is enveloped with a brief discussion about the insight of our main results in
Section 4.
2. Model Development
Motivated on experimental observations in [
5,
6,
7,
8,
9], and our previous works [
20,
21,
22,
23], we assume that (i) the maternal leakage rate is
: among offspring from infected females, where
of them are infected, and
of them are uninfected; (ii) CI intensity is
: among eggs produced from wild pregnant female mosquitoes mated with infected males,
of them will not hatch, and the hatch rate of other three compatible crosses is 1; (iii) equal sex determination: the birth ratio of male-to-female is 1:1, irrelevant of their mothers infection status. Thus, we assume that at the
n-th generation, the number of infected female/male mosquitoes is
, and the number of uninfected female/male mosquitoes is
. The infection frequency is then defined by
. At the
-th generation, the dynamic of the mosquitoe population is interfered upon by the supplemental release of infected males, whose abundance is assumed to be
. Here,
can be interpreted as the release ratio of released males to the total number of female or male mosquitoes.
Although the supplemental release of the infected males causes temporary biased male-to-female ratios at the
n-th generation, we assume that the sex ratio of the newly borne population at the
-th generation is 1:1, due to the assumption of equal sex determination. With the release of these infected males, the production rate, as well as the infection frequency, is reduced. Since infected females are not supplementally released, the proportions of infected and uninfected females at the
n-th generation are still
and
, respectively. With
infected males added at the
n-th generation, the proportion of infected males goes from
to:
and the proportion of uninfected males goes from
to:
At the
-th generation, the infection status of infected offspring is independent on the male partner’s infection status, the proportion of infected individuals is
, where
represents the fitness cost of infected females relative to uninfected ones, that is, if the fitness of
Wolbachia-free mosquitoes is 1, then the fitness of infected mosquitoes is reduced to
. We should mention here, that in our model,
counts the fitness costs caused by
Wolbachia infection, such as fecundity reduction, egg viability, hatch rate decrease, and mean longevity [
5,
6,
7,
8,
9]. Meanwhile, the maternal transmission leakage rate
reduces the proportion of infected mosquitoes from
to
. And the contribution of infected mothers to uninfected proportion is
.
On the uninfected offspring, we need to consider CI. Under the random mating assumption [
11,
13,
14,
16], the mating probability for a female mating with an uninfected male is
, and the hatch rate is 1, which is reduced to
if the father is infected, with the mating probability then being
. Hence, the contribution of uninfected mothers to the proportion of uninfected mosquitoes at the
-th generation is:
To sum up, we have our general model as:
with
.
The generalized model (
1) covers the models in [
11,
12,
13,
14,
15,
16] as special cases. For example, the first discrete model studying
Wolbachia infection frequency with laboratory data and parameters was developed in 1959 by Caspari and Watson [
11], who were motivated by the CI mechanism in mosquitoes. The model in [
11] is a special case of model (
1) with
and
, which produces an unstable equilibrium point
such that only when the initial infection frequency surpasses
can
Wolbachia invasion be guaranteed. After that, the discrete model for
Wolbachia spread dynamics has been put on the shelf for more than 15 years until 1975. Fine [
12,
13] observed the maternal leakage from
Wolbachia-infected mothers and introduced the maternal leakage rate
into the model in [
11] as the third crucial parameter. Again, the model in [
12,
13] is a special case of model (
1) with
. For
, since the 1990’s, systematic studies have been carried out by Hoffmann and Turelli on
Wolbachia spread dynamics in
Drosophila simulans [
14,
15,
16]. The dynamics can be completely determined by the position of the unstable equilibrium frequency, which can be explicitly expressed in terms of
,
and
. Once prior information on
,
and
is known, the unstable equilibrium frequency determines the fate of the
Wolbachia invasion. After that, discrete models in [
11,
12,
13,
14,
15,
16] have been seldom revisited due to their clear dynamics, until 2019, when we introduced the notion of supplementary releases into the discrete model [
24] for the first time, as in model (
1), with
.
Although model (
1) is simple in its form, it can generate much richer dynamics than the models in previous works. In the current work, we begin with the easiest case with
and
, i.e., the CI is complete and the maternal transmission is perfect. In this case, model (
1) is reduced to:
with
. Particularly, when
, model (
2) is further reduced to:
By defining:
From Theorems 3.1–3.3 in [
24] with
and
, we can obtain the following lemma.
Lemma 1 ([
24]).
For , model (3) has a zero equilibrium point and a Wolbachia fixation equilibrium point . Furthermore,- (1)
When , model (3) has an unstable equilibrium point , which serves as a threshold value on the initial frequency : if , then as , whereas for , as . - (2)
When , we have for any , i.e., the equilibrium is globally asymptotically stable.
Here, we aim to study the case that a periodic and impulsive release strategy, in each release period
T,
Wolbachia-infected male mosquitoes are released with
during the first
N generations for
, and
for
. Then the model switches between:
for
, and:
for
. We assume that
holds without loss of generality throughout the paper.
It’s easy to see that when , we get . Therefore, for any , is strictly monotonously increasing to 1, and we can reach our first result as follows.
Theorem 1. When , the trivial equilibrium 0 of model (4) and (5) is unstable, and Wolbachia fixation equilibrium point 1 is globally asymptotically stable. 3. At Least One Nontrivial -Periodic Solution to (4) and (5)
Let
denote the solution of model (
4) and (
5) with initial value
. Then
satisfies (
4) with the initial value
and
, which offers:
by induction. Initiated from
,
satisfies (
5) with
Hence,
Thus, the dynamics of model (
4) and (
5) depends on the switching frquency between Equations (
4) and (
5). If models (
4) and (
5) have a
T-periodic solution initiated from
u, then
. Therefore, finding the periodic solutions of models (
4) and (
5) is equivalent to finding the fixed points of
h. From (
6), (
7), and the facts that
and
, we get:
which means that 0 and 1 are two fixed points of
h, or equivalently, models (
4) and (
5) has 0 and 1 as two trivial
T-periodic solutions.
Next, we seek sufficient conditions to guarantee that models (
4) and (
5) have at least one nontrivial
T-periodic solution. We begin with the case
.
In this case, Lemma 1 implies that both (
4) and (
5) produce bistable dynamics with the existence of unstable equilibrium points, denoted by
and
with
for
, respectively. Since
, we reach
. Furthermore, Lemma 1 implies that,
The continuity of
in
u guarantees that
h has at least one fixed point in
. Therefore, we have:
Lemma 2. When , model (4) and (5) has at least one nontrivial T-periodic solution initiated from . For the case
, similarly from Lemma 1, we have:
However, the sign of
is indefinite for
. This leads us to calculate
. Since
, we reach:
where:
Since:
for
with
. Therefore, from (
12),
Similarly,
Hence, combining (
11), (
13), and (
14), we have:
which implies that when
lie in:
we have
. Therefore, there exists
such that:
where
is sufficiently small. Recalling that
, there exists at least one
such that:
The above analysis leads to another sufficient condition for the existence of nontrivial
T-periodic solutions to the models (
4) and (
5) as follows.
Lemma 3. If and , then models (4) and (5) have at least one nontrivial T-periodic solution with u lying in , where δ is defined in (16). Combining the results shown in Lemmas 2 and 3, sufficient conditions to guarantee the existence of nontrivial
T-periodic solutions to the models (
4) and (
5) are stated as follows.
Theorem 2. Models (4) and (5) have at least one nontrivial T-periodic solution, if one of the following two conditions holds,where D is defined in (15). The following example displays that when the conditions in (
17) hold, the
T-periodic solution not only exists, but also may be unique.
Example 1. Let . Then . When we take and , models (4) and (5) read:Letting and . Then (i) in (17) holds. By randomly selecting 400 points for , we solve (18) in Matlab and plot against u in Figure 1A. The initial value for the unique nontrivial periodic solution is dotted in red. Similarly, when taking and , we have:and hence (ii) in (17) holds. Figure 1B again displays that the nontrivial fixed point of lying in is unique. Although the uniqueness of the nontrivial periodic solutions in (
18) is verified as the uniqueness of the nontrivial fixed points in
, we failed to theoretically prove this observation for the general
N and
T. In the next section, we prove the uniqueness of the nontrivial fixed points of
, or equivalently, the uniqueness of the nontrivial periodic solutions to models (
4) and (
5) for the special case with
and
. Methods to prove the conclusion for the general
N and
T are still pending.
A Unique Nontrivial T-Periodic Solution for and
We claim that when
and
, condition (
17) is also sufficient to guarantee that the nontrivial periodic solution to the models (
4) and (
5) is unique, i.e.,
Theorem 3. For and , model (4) and (5) admits a unique nontrivial periodic solution initiated from , when (17) holds. To make the proof of Theorem 3 compact, we provide the following two lemmas as the preliminary analysis. The first one is devoted to the calculation at u satisfying .
Lemma 4. For , we have:where: Proof. Let
be the solution of models (
4) and (
5) initiated from
. Therefore from (
4) and (
6), we have:
and from (
5) and (
7), we get:
Then, taking the derivative on both sides of (
21) with respect to
u, we get:
Similarly, taking the derivative on both sides of (
22) with respect to
, we get:
Combining (
23) and (
24), for
, we reach:
where:
and:
which leads to (
19), proving the lemma. □
The uniqueness of the nontrivial periodic solutions to models (
4) and (
5) stated in Theorem 3 will be proved by the uniqueness of the nontrivial fixed points of the Poincaré map
. To this end, the comparison between
and 1 for
will be visited often. Tedious but simple calculations offer:
where:
with:
Similarly, we have:
Hence, the relation between and 1 for is totally determined by the sign of the quartic polynomial . The next lemma shows that when , is strictly decreasing in u.
Lemma 5. When , we have for all .
Proof. We firstly prove that:
In fact, since:
we get:
and:
where:
and:
Therefore,
and (
27) holds.
We next prove that
for all
. From:
where
are defined in (
26). Hence,
whose sign is indefinite. However,
Furthermore, the
x-coordinate of the vertex of
, denoted by
, is:
It is easy to see that
since
. Regarding the sign of
, we have:
Therefore we have
when
. The facts
,
and
imply that the function
defined in (
29) belongs to one of the following cases.
- (i)
If and , then for .
- (ii)
If and , then there exists a unique such that for and for .
- (iii)
If , then for .
The function is strictly increasing from to for for case (i). For case (ii), function decreases from to , and then increases from to for . For case (iii), function is strictly decreasing from to for . For all these three cases, we have for . The proof is completed. □
Proof of Theorem 3 We assume that models (
4) and (
5) have another nontrivial
T-periodic solution initiated from
by contradiction. See
Figure 2 for illustration.
Then two cases occur:
On case (i), we have
,
, which lead to
and
, this contradicts the conclusion in Lemma 5. In case (ii), we get
,
, which goes well with the monotonicity of function
. To achieve a contradiction, we rewrite:
to consider the monotonicity of these functions in
.
Then, from (
22) and (
32), we have:
which implies that
is strictly increasing in
. Hence, the function
is also strictly increasing in
.
Finding
small enough, such that
, from the monotonicity of
respect of
, we observe that the curve
totally lies below
except at
, which admits three zeros,
with:
See
Figure 2B for illustration. This implies that
and
, which contradicts the monotonicity of
in
u, proving the Theorem. □
4. Discussion
The endosymbiotic bacterium Wolbachia was first identified in 1924, which naturally infects 60% of insect species. In mosquitoes, Wolbachia often induces CI which results in early embryonic death when Wolbachia-infected males mate with uninfected females. However, Wolbachia-infected females produce viable and infected embryos owing to the mechanism of maternal transmission, irrelevant of the paternal infection status. The CI mechanism and maternal transmission of Wolbachia bring a reproductive advantage of infected mosquitoes over uninfected ones. Regarding Wolbachia release, the infected males can sterilize uninfected females through CI, and the release of infected females can increase Wolbachia infection frequency in mosquito populations through maternal transmission. To date, two release strategies have been applied in field trials: population suppression by only releasing infected males to sterilize wild females, and population replacement by releasing both infected females and males so that wild mosquitoes are replaced by infected ones.
The field trials of both population suppression and population replacement have been implemented in several countries, including Australia, Brazil, and China, which shows that
Wolbachia release is becoming a safe and sustainable approach to controlling mosquito populations and mosquito-borne disease transmission. Particularly, since March 2015, by combining the incompatible and sterile insect techniques (IIT-SIT),
Wolbachia-infected male mosquitoes have been released on two isolated islands in south of Guangzhou city. The implementation of IIT-SIT in 2016 and 2017 eradicated more than 95% of wild-type
Aedes albopictus populations [
25]. As collaborates, we took charge in mathematical modeling in this project, and our discrete mathematical model driven by semi-field cage experiments
“accurately described and predicted target population dynamics in the semi-field cage experiments, and supported the notion that a 5:1 over-flooding ratio of HC (the name of the Wolbachia strain) to wild-type males is sufficient for effective population suppression and/or elimination”.
The notion of
supplemental releases is motivated by the implementation of incompatible and sterile insect techniques combined to eliminate mosquitoes in Guangzhou, Guangdong China since 2015 [
25], where supplementary releases are frequently carried out to guarantee the population suppression. Very recently, we summarized works on
Wolbachia spread dynamics since 1959 in [
26] and proposed some open questions on the periodic and impulsive release strategies. We took our first try in [
27] to answer one of them. This is a continual study aiming to count the exact numbers of nontrivial periodic solutions to the discrete switching model, which is a challenging task to achieve due to the lack of effective methods.
In [
17,
18,
19], the authors developed mosquito population suppression models by switching to ordinary differential equation (ODE) models in each release waiting period
T. The global dynamics of these models were elegantly described by using the corresponding Pioncaré maps
, where
u is the initial point. By solving the initial value problems in the first
T-periodic interval
, the explicit expression of
was obtained, and then the existence, and uniqueness of the periodic solutions to the ODE models was equivalent to the existence and uniqueness of the fixed points to the map
. Particularly, at
,
can also be explicitly expressed in
u, which is a key step to showing the uniqueness of the periodic solutions to the ODE models.
Following this procedure, we obtained sufficient conditions to guarantee that the models (
4) and (
5) have at least one nontrivial periodic solution in Theorem 2. However, this method of proving the uniqueness does not work for the discrete models developed in this paper except for the case
and
under which the comparison between
and 1 can be translated to the comparison between a quartic polynomial
and 0. By developing the qualitative property of
, sufficient conditions to guarantee the existence and uniqueness of the discrete switching model were obtained in Theorem 3. Unfortunately, for the general
N and
T, as shown in (
7), the Poincaré map
is a composite function whose explicit expression is absent. Furthermore, when calculating
for the general
N and
T, we find that it is an impossible mission to get a similar expression as shown in (
19). This drives us to seek new methods to treat the general cases in our future research.