Implementation of Control Strategies for Sterile Insect Techniques

TL;DR

This paper develops and analyzes control strategies for eliminating wild mosquito populations using sterile insect techniques, focusing on open-loop, closed-loop, and mixed impulsive release methods with proven convergence conditions.

Contribution

It introduces a comprehensive mathematical framework for designing and analyzing both open-loop and closed-loop control strategies for sterile insect releases, including a novel mixed approach.

Findings

Global convergence to mosquito-free state under specified conditions.

Closed-loop control adapts to real-time population estimates.

Mixed strategy optimizes elimination time and release efficiency.

Abstract

In this paper, we propose a sex-structured entomological model that serves as a basis for design of control strategies relying on releases of sterile male mosquitoes (Aedes spp) and aiming at elimination of the wild vector population in some target locality. We consider different types of releases (constant and periodic impulsive), providing necessary conditions to reach elimination. However, the main part of the paper is focused on the study of the periodic impulsive control in different situations. When the size of wild mosquito population cannot be assessed in real time, we propose the so-called open-loop control strategy that relies on periodic impulsive releases of sterile males with constant release size. Under this control mode, global convergence towards the mosquito-free equilibrium is proved on the grounds of sufficient condition that relates the size and frequency of releases. If periodic assessments (either synchronized with releases or more sparse) of the wild population size are available in real time, we propose the so-called closed-loop control strategy, which is adjustable in accordance with reliable estimations of the wild population sizes. Under this control mode, global convergence to the mosquito-free equilibrium is proved on the grounds of another sufficient condition that relates not only the size and frequency of periodic releases but also the frequency of sparse measurements taken on wild populations. Finally, we propose a mixed control strategy that combines open-loop and closed-loop strategies. This control mode renders the best result, in terms of overall time needed to reach elimination and the number of releases to be effectively carried out during the whole release campaign, while requiring for a reasonable amount of released sterile insects.

Tables

Table 1. Table 1: Parameters of the sex-structured entomological model ( 2 )

Table 2. Table 2: Aedes spp parameters values

Parameter	Value	Description
$\rho$	4.55	Number of eggs a female can deposit
$r$	0.5	$r:\left(1-r\right)$ expresses the primary sex ratio in offspring
$\sigma$	0.05	Regulates the larvae development into adults under
		density dependence and larval competition
$K$	140	Carrying capacity
${\mu }_M$	0.04	Mean mortality rate of wild adult male mosquitoes
${\mu }_F$	0.03	Mean mortality rate of wild adult female mosquitoes
${\mu }_S$	0.04	Mean mortality rate of sterile adult male mosquitoes
$\gamma$	1	Fitness of sterile adult male mosquitoes

Table 3. Table 3: Cumulative number of released sterile males for each open-loop periodic SIT control treatment.

Period (days)	Cumulative Number of released sterile males	Nb of Weeks to reach elimination
$\tau =7$	$\mathrm{𝟗𝟐𝟒},\mathrm{𝟔𝟐𝟕}$	84
$\tau =14$	$942,869$	84

Table 4. Table 4: Cumulative number of released sterile males and number of releases for each closed-loop periodic SIT control treatment when k = 1 5 ​ 𝒩 F 𝑘 1 5 subscript 𝒩 𝐹 k=\dfrac{1}{5}{\cal N}_{F} . See Figure 4 , page 4 .

	Cumulative Nb of		Nb of weeks needed
	released sterile males		to reach elimination		Nb of nonzero releases
	$1$	$4$	$1$	$4$	$4$
$\tau =7$	$2,251,052$	$4,363,430$	$64$	$54$	$34$
$\tau =14$	$2,390,676$	$2,896,835$	64	56	17

Table 5. Table 5: Cumulative number of released sterile males and number of releases for each closed-loop periodic SIT control treatment when k = 0.99 𝒩 F 𝑘 0.99 subscript 𝒩 𝐹 k=\dfrac{0.99}{{\cal N}_{F}} . See Figure 5 , page 5

	Cumulative Nb of		Nb of weeks needed
	released sterile males		to reach elimination		Nb of nonzero releases
	$1$	$4$	$1$	$4$	$4$
$\tau =7$	$794,807$	$1,221,593$	240	58	37
$\tau =14$	$909,344$	$1,043,107$	130	62	20

Table 6. Table 6: Cumulative number of released sterile males and number of releases for each mixed open/closed-loop periodic SIT control treatment when k = 0.2 𝒩 F 𝑘 0.2 subscript 𝒩 𝐹 k=\dfrac{0.2}{{\cal N}_{F}} . See Figure 7 , page 7 .

	Cumulative Nb of		Nb of weeks needed
	released sterile males		to reach elimination		Nb of nonzero releases
	$1$	$4$	$1$	$4$	$4$
$\tau =7$	$450,668$	$534,849$	72	65	53
$\tau =14$	$465,187$	$499,497$	72	66	25

Table 7. Table 7: Cumulative number of released sterile males and number of releases for each mixed open/closed-loop periodic SIT control treatment when k = 0.99 𝒩 F 𝑘 0.99 subscript 𝒩 𝐹 k=\dfrac{0.99}{{\cal N}_{F}} . See Figure 8 , page 8 .

	Cumulative Nb of		Nb of weeks needed
	released sterile males		to reach elimination		Nb of nonzero releases
	$1$	$4$	$1$	$4$	$4$
$\tau =7$	$457,489$	$450,077$	246	69	53
$\tau =14$	$427,701$	$449,059$	136	74	28

Equations

$$|z{|}_{-}:=\{\begin{array}{ll}z& \text{ if }z\le 0\\ 0& \text{otherwise}\end{array}\text{ and }|z{|}_{+}:=\{\begin{array}{ll}z& \text{ if }z\ge 0\\ 0& \text{otherwise}\end{array}$$

$$[left=\text{empheqlbrace}]\dot{M}$$

$$\dot{F}$$

$${\mu }_M\ge {\mu }_F.$$

$$\mathcal{D}=\{(M,F):0\le M\le C,0\le F\le C\}$$

$${\mathcal{N}}_F:=\frac{(1-r)\rho }{{\mu }_F},{\mathcal{N}}_M:=\frac{r\rho }{{\mu }_M}.$$

$$F^{\ast }=\frac{{\mathcal{N}}_F}{{\mathcal{N}}_F+{\mathcal{N}}_M}\frac{1}{\beta }\ln {\mathcal{N}}_F,M^{\ast }=\frac{{\mathcal{N}}_M}{{\mathcal{N}}_F+{\mathcal{N}}_M}\frac{1}{\beta }\ln {\mathcal{N}}_F.$$

$$\dot{F}=\left((1-r)\rho e^{-\beta(M+F)}-\mu_{F}\right)F\leq\big{(}(1-r)\rho-\mu_{F}\big{)}F$$

$$J(M,F)=\left(\begin{array}[]{cc}-\beta r\rho Fe^{-\beta(M+F)}-\mu_{M}&r\rho(1-\beta F)e^{-\beta(M+F)}\\ &\\ -\beta(1-r)\rho Fe^{-\beta(M+F)}&(1-r)\rho(1-\beta F)e^{-\beta(M+F)}-\mu_{F}\end{array}\right)$$

$$J(E_{0}^{*})=\left(\begin{array}[]{cc}-\mu_{M}&r\rho\\ 0&(1-r)\rho-\mu_{F}\end{array}\right),$$

$$J\big{(}E^{*}\big{)}=\left(\begin{array}[]{cc}-\dfrac{\beta r\rho}{{\cal N}_{F}}F^{*}-\mu_{M}&\dfrac{r\rho}{{\cal N}_{F}}\big{(}1-\beta F^{*}\big{)}\\ &\\ -\dfrac{\beta(1-r)\rho}{{\cal N}_{F}}F^{*}&-\dfrac{\beta(1-r)\rho}{{\cal N}_{F}}F^{*}\end{array}\right).$$

$$\det J\big{(}E^{*}\big{)}=\dfrac{\beta}{{\cal N}_{F}}(1-r)\rho F^{*}\left(\mu_{M}+\dfrac{r\rho}{{\cal N}_{F}}\right)>0$$

$${\psi }_1(F):=\frac{1}{F},f_1(M,F):=r\rho F{e}^{-\beta (M+F)}-{\mu }_{M}M,g_1(M,F):=(1-r)\rho F{e}^{-\beta (M+F)}-{\mu }_{F}F,$$

$$D_{1}(M,F):=\dfrac{\partial}{\partial M}\Big{(}\psi_{1}(F)f_{1}(M,F)\Big{)}+\dfrac{\partial}{\partial F}\Big{(}\psi(F)g_{1}(M,F)\Big{)}.$$

$$\dfrac{\partial}{\partial M}\Big{(}\psi_{1}(F)f_{1}(M,F)\Big{)}=-\beta r\rho e^{-\beta(M+F)}-\dfrac{\mu_{M}}{F},\qquad\dfrac{\partial}{\partial F}\Big{(}\psi_{1}(F)g_{1}(M,F)\Big{)}=-\beta(1-r)\rho e^{-\beta(M+F)},$$

$$D_1(M,F)=-\beta \rho e^{-\beta (M+F)}-\frac{{\mu }_M}{F}<0$$

$$M(t)=e^{-{\mu }_{M}t}M(0)+r\rho \underset{0}{\overset{t}{\int }}{e}^{-{\mu }_M(t-s)}F(s)e^{-\beta (M(s)+F(s))}ds>0$$

$$\frac{d}{dt}\left(\frac{F}{M}\right)=\frac{F}{M}\left({\mu }_M-{\mu }_F+\rho e^{-\beta (M+F)}\left(1-r-r\frac{F}{M}\right)\right)>({\mu }_M-{\mu }_F)\frac{F}{M}\text{ if }\frac{F}{M}\le \frac{1-r}{r}.$$

$$\forall t\ge T,\frac{F}{M}>\frac{1-r}{r}.$$

$$\dot{F}=\left((1-r)\rho e^{-\beta (M+F)}-{\mu }_F\right)F\ge \left((1-r)\rho e^{-\frac{\beta }{1-r}F}-{\mu }_F\right)F$$

$$\underset{t\to +\infty }{\mathrm{liminf}}F\ge \frac{1-r}{\beta }\ln {\mathcal{N}}_F>0.$$

$$\forall t\ge T^{\prime },F{e}^{-\beta F}\ge \delta >0.$$

$$\dot{M}\ge r\rho \delta e^{-\beta M}-{\mu }_{M}M,$$

$$[left=\text{empheqlbrace}]\dot{M}$$

$$\dot{F}$$

$${\dot{M}}_S$$

$${\mu }_S\ge {\mu }_M.$$

$$[left=\text{empheqlbrace}]\dot{M}$$

$$\dot{F}$$

$${\Lambda }^{crit}:=2\frac{{\mu }_S}{\beta \gamma }\frac{{\varphi }^{crit}({\mathcal{N}}_F)}{1+\frac{{\mathcal{N}}_F}{{\mathcal{N}}_M}},$$

$$1+\varphi \left(1+\sqrt{1+\frac{2}{\varphi }}\right)={\mathcal{N}}_F\exp \left(-\frac{2}{1+\sqrt{1+\frac{2}{\varphi }}}\right).$$

$$r\rho \frac{F^{\ast }}{M^{\ast }+\gamma M_S^{\ast }}{e}^{-\beta (M^{\ast }+F^{\ast })}={\mu }_M,(1-r)\rho \frac{M^{\ast }}{M^{\ast }+\gamma M_S^{\ast }}{e}^{-\beta (M^{\ast }+F^{\ast })}={\mu }_F.$$

Full text

Keywords: Sterile Insect Technique, periodic impulsive control, open-loop and closed-loop control, global stability, exponential convergence, saturated control.

1. 1 Introduction

2. 2 A sex-structured entomological model

3. 3 Elimination with constant releases of sterile insects

4. 3.1 Existence of positive equilibria

5. 3.2 Asymptotic stability of the equilibria

6. 4 Elimination with periodic impulsive releases of sterile insects

7. 5 Elimination by feedback control

8. 5.1 Principle of the method

9. 5.1.1 Step 1 – Setting directly the sterile population level

10. 5.1.2 Step 2 – Shaping an impulsive control compliant with Step 1

11. 5.2 Stabilization result

12. 5.2.1 Synchronized measurements and releases

13. 5.2.2 Sparse measurements

14. 6 Mixed control strategies

15. 7 Numerical illustrations

16. 8 Conclusion

1 Introduction

Since decades, the control of vector-borne diseases has been a major issue in Southern countries. It recently became a major issue in Northern countries too. Indeed, the rapid expansion of air travel networks connecting regions of endemic vector-borne diseases to Northern countries, the rapid invasion and establishment of mosquitoes population, like Aedes albopictus, in Northern hemisphere have amplified the risk of Zika, Dengue, or Chikungunya epidemics222See, for instance, the most recent distribution map of Aedes albopictus provided by ECDC (European Centre for Disease Prevention and Control, https://ecdc.europa.eu/en/publications-data/aedes-albopictus-current-known-distribution-june-2018).

For decades, chemical control was the main tool to control or eradicate mosquitoes. Taken into account resistance development and the impact of insecticides on the biodiversity, other alternatives have been developed, such as biological control tools, like the Sterile Insect Technique (SIT).

Sterile Insect Technique (SIT) is a promising control method that has been first studied by E. Knipling and collaborators and first experimented successfully in the early 1950’s by nearly eradicating screw-worm population in Florida. Since then, SIT has been applied on different pest and disease vectors (see [10] for an overall presentation of SIT and its applications).

The classical SIT relies on massive releases of males sterilized by ionizing radiations. However, another technique, called the Wolbachia technique, is under consideration because Wolbachia [15] is a symbiotic bacterium that infects many Arthropods, including some mosquito species in nature. These bacteria have many particular properties, including one that is very useful for vector control: the cytoplasmic incompatibility (CI) property [3, 25]. CI can be used for two different control strategies:

• •

Incompatible Insect Technique (IIT): males infected with CI-inducing Wolbachia produce altered sperms that cannot successfully fertilize uninfected eggs. This can result in a progressive reduction of the target population. Thus, IIT can be seen as equivalent to classical SIT.

• •

Population Replacement (PR): in this case, males and females, infected with CI-inducing Wolbachia, are released in a susceptible (uninfected) population, such that Wolbachia-infected females will produce more offspring than uninfected females. Because Wolbachia is maternally inherited, this will result in a population replacement by Wolbachia-infected mosquitoes (such replacements or invasions have been observed in natural population, see [23] for the example of Californian Culex pipiens). Recent studies have shown that PR may be very interesting with Aedes aegypti, shortening their lifespan (see for instance [24]), or more interesting, cutting down their competence for dengue virus transmission [20]. However, it is also acknowledged that Wolbachia infection can have fitness costs, so that the introgression of Wolbachia into the field can fail [24].

Based on these biological properties, classical SIT and IIT (see [7, 8, 9, 16, 19, 27] and references therein) or population replacement (see [4, 5, 11, 12, 13, 17, 21, 24, 28] and references therein) have been modeled and studied theoretically in a large number of papers, in order to derive results to explain the success or not of these strategies using discrete, continuous or hybrid modeling approaches, temporal and spatio-temporal models. More recently, the theory of monotone dynamical systems [26] has been applied efficiently to study SIT [1, 27] or population replacement [2, 18] systems.

In this paper, we derive and study a dynamical system to model the release and elimination process for SIT/IIT. We analyze and compare constant continuous/periodic impulsive releases and derive conditions relating the sizes and frequency of the releases that are sufficient to ensure successful elimination. Such conditions enable the design of SIT-control strategies with constant or variable number of sterile males to be released that drive the wild population of mosquitoes towards elimination. Among all the previous strategies, we are also able to derive the (best) strategy that needs to release the least amount of sterile males to reach elimination. This can be of utmost importance for field applications.

The outline of the paper is as follows. In Section 2, we first develop and briefly study a simple entomological model that describes the natural evolution of mosquitoes. Then, in Section 3, we introduce a constant continuous SIT-control and estimate the size of constant releases that ensures global elimination of wild mosquitoes in the target locality. In Section 4, periodic impulsive SIT-control with constant impulse amplitude is considered, and a sufficient condition relating the size and frequency of periodic releases is derived to ensure global convergence towards the mosquito-free equilibrium. This condition enables the design of open-loop (or feedforward) strategies that ensure mosquito elimination in finite time and without assessing the size of wild mosquito population. Alternatively, Section 5 is focused on the design of closed-loop (or feedback) SIT-control strategies, which are achievable when periodic estimations (either synchronized with releases or more sparse) of the wild population size are available in real time. In such situation, the release amplitude is computed on the basis of these measurements. Thorough analysis of the feedback SIT-control implementation mode leads to another sufficient condition to reach mosquito elimination. This condition relates not only the size and frequency of periodic releases but also the frequency of sparse measurements. Finally, in Section 6 we propose a mixed control strategy for periodic impulsive SIT-control. The latter is essentially based on the use of the smallest of the release values proposed by the previous open-loop and closed-loop strategies. It turns out that this control mode renders the best result from multiple perspectives: in terms of overall time needed to reach elimination and of peak-value of the input control, but also in terms of total amount of released sterile insects and of number of releases to be effectively carried out during a whole SIT-control campaign. The paper ends with numerical simulations highlighting the key features and outcomes of periodic impulsive SIT-control strategies (Section 7) followed by discussion and conclusions.

Notations

For any $z\in\mathbb{R}$, define the decomposition $z=|z|_{-}+|z|_{+}$ in negative and positive parts, fulfilling:

[TABLE]

2 A sex-structured entomological model

We consider the following 2-dimensional system to model the dynamics of mosquito populations. It involves two state variables, the number of males $M$ and the number of females $F$.

[TABLE]

All the parameters are positive, and listed in Table 1. The model assumes that all females are equally able to mate. It includes direct and/or indirect competition effect at different stages (larvae, pupae, adults), through the parameter $\beta$. The latter may be seen as the ratio, $\frac{\sigma}{K},$ between $\sigma$, a quantity characterizing the transition between larvae and adults under density dependence and larval competition, and a carrying capacity $K$, typically proportional to the breeding sites capacity. The primary sex ratio in offspring is denoted $r$, and $\rho$ represents the mean number of eggs that a single female can deposit in average per day. Last, $\mu_{M}$ and $\mu_{F}$ represent, respectively, the mean death rate of male and female adult mosquitoes. As a rule, it is observed that the male mortality is larger, and we assume throughout the paper that:

[TABLE]

Existence and uniqueness of the solutions of the Cauchy problem for dynamical system (2) follow from standard theorems, ensuring continuous differentiability of the latter. System (2) is dissipative: there exists a bounded positively invariant set $\mathcal{D}$ with the property that, for any bounded set in $E\subset\mathbb{R}_{+}^{2}$, there exists $t^{*}=t(\mathcal{D},E)$ such that \big{(}M(0),F(0)\big{)}\in E implies \big{(}M(t),F(t)\big{)}\in\mathcal{D} for all $t>t^{*}$. The set $\mathcal{D}$ is called an absorbing set. In our case, it may be taken, e.g., as:

[TABLE]

for some $C>0$.

Remark 1**.**

Population models of the form $\dot{N}=B(N)N-\mu N$ for several birth rate functions, including $B(N)=e^{-\beta N}$, have been studied in [6]. Maturation delay can also be included [6].

Obviously $E^{*}_{0}=(0,0)$ is a trivial equilibrium of system (2), called the mosquito-free equilibrium. Being the state to which one desires to drag the system by adequate releases of sterile insects, it will play a central role in the sequel. Denote for future use

[TABLE]

These positive constants represent basic offspring numbers related to the wild female and male populations, respectively. The first of them governs the number of equilibria, as stated by the following result, whose proof presents no difficulty and is left to the reader.

Theorem 1** (Equilibria of the entomological model).**

**

• •

If ${\cal N}_{F}<1$, then system (2) possesses $E^{*}_{0}$ as unique equilibrium.

• •

If ${\cal N}_{F}>1$, then system (2) also possesses a unique positive equilibrium $E^{*}:=(M^{*},F^{*})$, namely

[TABLE]

Notice that the total population at the nonzero equilibrium is given by $M^{*}+F^{*}=\dfrac{1}{\beta}\ln{\cal N}_{F}$. It depends upon the basic offspring number and the competition parameter $\beta$. As an example, mechanical control through reduction of the breeding sites induces increase of $\beta$ and consequently decrease of the population at equilibrium. Analogously, altering biological parameters modifies the basic offspring number, and therefore the size of the population.

The stability of the equilibria is addressed by the following result.

Theorem 2** (Stability properties of the entomological model).**

**

• •

If ${\cal N}_{F}<1$, then the (unique) equilibrium $E^{*}_{0}$ is Globally Asymptotically Stable (GAS).

• •

If ${\cal N}_{F}>1$, then $E^{*}_{0}$ is unstable, and $E^{*}$ is GAS in $\mathcal{D}\setminus{\{(M,0),M\in\mathbb{R}_{+}\}}$.

Figure 1 shows the convergence of all trajectories to the positive equilibrium in a case where ${\cal N}_{F}>1$ (the pertinent case for the applications we have in mind).

Proof of Theorem 2.

$\bullet$ Assume first ${\cal N}_{F}<1$. Rewriting equation (2b) as follows:

[TABLE]

one deduces that $\dot{F}<-\varepsilon F$ for some positive $\varepsilon$. The state variable $F$ being nonnegative, it then converges to [math]. Using now equation (2a), we deduce that $M$ converges to [math] too, and the GAS of $E^{*}_{0}$ follows.

$\bullet$ Assume now that ${\cal N}_{F}>1$. Let us compute the Jacobian matrix related to entomological system (2), page 2:

[TABLE]

so that

[TABLE]

from which we deduce that $E_{0}^{*}$ is locally asymptotically stable (LAS) if ${\cal N}_{F}<1$.

For the positive equilibrium $E^{*}$, using the fact that $e^{-\beta(M^{*}+F^{*})}=\dfrac{1}{{\cal N}_{F}}$, we have:

[TABLE]

Obviously {\mathrm{\tt trace}}\big{\{}J(E^{*})\big{\}}<0 and

[TABLE]

so that $E^{*}$ is LAS when ${\cal N}_{F}>1$.

Using Dulac criterion [22], we now show that system (2) has no closed orbits wholly contained in the attracting set $\cal D$. Indeed, setting

[TABLE]

let us study the sign of the function

[TABLE]

We have

[TABLE]

and thus

[TABLE]

for all $(M,F)\in\mathcal{D}$ such that $F>0$. Therefore, Dulac criterion [22] applies, demonstrating that system (2) possesses no nonconstant periodic solutions. Thus, using the fact that $E^{*}$ is LAS, by the Poincaré-Bendixson theorem, all trajectories in $\mathcal{D}\setminus\{(M,0)\ :\ M\geq 0\}$ converge towards $E^{*}$.

Convergence towards $E^{*}_{0}$ clearly occurs in absence of females, i.e. when $F(0)=0$. (Notice that for this reason, the point $E^{*}_{0}$ cannot be repulsive.) Consider on the contrary a trajectory such that $F(0)>0$. We will show that convergence to $E^{*}_{0}$ is impossible, so convergence towards $E^{*}$ occurs. First of all, one deduces from (2) and the continuity of $F$ that

[TABLE]

for any $t>0$. The ratio $\dfrac{F}{M}$ is therefore well defined and remains positive along this trajectory. It is moreover continuously differentiable, and

[TABLE]

From (3), it is deduced immediately that there exists for this trajectory a real number $T\geq 0$, such that

[TABLE]

Then it holds for any $t\geq T$ that

[TABLE]

and thus

[TABLE]

As the compact set $\mathcal{D}$ is absorbing, the trajectory is bounded. We deduce from this and (5) the existence of certain $\delta>0$ and $T^{\prime}\geq T$ (whose precise values depend upon the considered trajectory), such that

[TABLE]

Now, we have for any $t\geq T^{\prime}$

[TABLE]

which is strictly positive in a neighborhood of $M=0$. The trajectory under study therefore stays at a positive distance from the point $E^{*}_{0}$, and, being convergent, has to converge to the other equilibrium, namely $E^{*}$. This shows that any trajectory departing with $F(0)>0$ converges towards $E^{*}$, and finally concludes the proof of Theorem 2. ∎

3 Elimination with constant releases of sterile insects

We now extend system (2), in order to incorporate continuous, constant releases driven by an equation for $M_{S}$, the number of sterile males:

[TABLE]

The positive constants $\mu_{S}$ and $\gamma$ represent, respectively, the mortality rate of sterile insects, and their relative reproductive efficiency or fitness (compared to the wild males), which is usually smaller than 1. The nonnegative quantity $\Lambda$ is the number of sterile insects released per time unit. It is taken constant over time in the present section. The other parameters are the same as for model (2), see Table 1.

The mortality of the sterile males is usually larger than that of wild males, so in complement to (3), we also have:

[TABLE]

Assuming $t$ large enough, we may suppose $M_{S}(t)$ at its equilibrium value $M_{S}^{*}:=\dfrac{\Lambda}{\mu_{S}}$, and the previous system then reduces to

[TABLE]

System (8) is dissipative too, with all trajectories converging towards the same set $\mathcal{D}$ introduced in the previous section. It admits the same mosquito-free equilibrium $E^{*}_{0}$.

We are interested here in the issues of existence and stability of positive equilibria. Driven by the application in view, we assume that the mosquito population is viable (that is ${\cal N}_{F}>1$, see Theorem 1), and focus on conditions sufficient for its elimination.

3.1 Existence of positive equilibria

The mosquito-free equilibrium $E^{*}_{0}$ is always an equilibrium of system (8). The following result is concerned with possible supplementary equilibria.

Theorem 3** (Existence of positive equilibria for the SIT entomological model with constant releases).**

*Assume ${\cal N}_{F}>1$. Then *

• •

there exists $\Lambda^{crit}>0$ such that system (6) admits two positive distinct equilibria if $0<\Lambda<\Lambda^{crit}$, one positive equilibrium if $\Lambda=\Lambda^{crit}$, and no positive equilibrium if $\Lambda>\Lambda^{crit}$;

• •

the value of $\Lambda^{crit}$ is uniquely determined by the formula

[TABLE]

where $\phi^{crit}:=\phi^{crit}({\cal N}_{F})$ is the unique positive solution to the equation

[TABLE]

Theorem 3 provides a characterization of the constant release rate above which no positive equilibrium may appear. We prove in the next section (Section 3.2) that in such a situation, convergence towards the mosquito-free equilibrium $E^{*}_{0}$ occurs, that ensures elimination of the wild population.

Proof of Theorem 3.

Clearly, nullity of $M$ at equilibrium is equivalent to nullity of $F$. In order to find possible nonzero equilibria, let $(M^{*},F^{*})$ with $M^{*}>0,F^{*}>0$ be one of them. The populations at equilibrium have to fulfill:

[TABLE]

In particular, we have, for ${\cal N}_{F},{\cal N}_{M}$ defined in (4),

[TABLE]

which imply the relation:

[TABLE]

Injecting this value in the first equation of (11), the number of males $M^{*}$ at equilibrium has to verify the equation

[TABLE]

or again

[TABLE]

The study of equation (12) is done through the following result, whose proof is given in Appendix.

Lemma 1**.**

Let $a,b,c$ be positive constants, with $b>1$. Then the equation

[TABLE]

admits a unique positive root, denoted $\phi^{crit}$. Moreover, the equation

[TABLE]

admits two positive distinct roots if $0<ac<2\phi^{crit}$; one positive root if $ac=2\phi^{crit}$; no positive root otherwise.

Using Lemma 1 with

[TABLE]

one deduces that equation (12) admits exactly one positive root when the root of (10) is equal to $\phi^{crit}=\frac{ac}{2}=\frac{1}{2}\beta\gamma\left(1+\frac{{\cal N}_{F}}{{\cal N}_{M}}\right)\dfrac{\Lambda^{crit}}{\mu_{S}}$, which implies (9) and thus achieves the proof of Theorem 3. ∎

3.2 Asymptotic stability of the equilibria

Assume ${\cal N}_{F}>1$. We first study the asymptotic stability of the mosquito-free equilibrium $E^{*}_{0}$ in the case where it is the unique equilibrium, that is when $\Lambda>\Lambda^{crit}$.

Theorem 4** (Stability of the mosquito-free equilibrium of the SIT entomological model with constant releases).**

If system (6) admits no positive equilibrium (that is, if $\Lambda>\Lambda^{crit}$), then the mosquito-free equilibrium $E^{*}_{0}$ is globally exponentially stable.

Proof of Theorem 4.

The Jacobian matrix of the reduced system (8) is equal to

[TABLE]

Its value at the mosquito-free equilibrium $E^{*}_{0}$ is just ${\mathrm{\tt diag}}\{-\mu_{M};-\mu_{F}\}$, which guarantees local asymptotic stability at this point.

We use again Dulac criterion to show that system (8) has no closed orbits wholly contained in the set $\mathcal{D}$. We set

[TABLE]

and then study the sign of the function

[TABLE]

As

[TABLE]

one has

[TABLE]

for all $(M,F)\in\mathcal{D}$ such that $F>0$. Thus, by the Poincaré-Bendixson theorem, since $E^{*}_{0}$ is the only asymptotically stable equilibrium, all trajectories in $\mathcal{D}$ approach the equilibrium $E^{*}_{0}$. This concludes the proof of Theorem 4. ∎

On the other hand, when $\Lambda<\Lambda^{crit}$ is not large enough and system (8) admits two distinct positive equilibria $E^{*}_{1}<E^{*}_{2}$, one may show by studying the spectrum of the Jacobian matrices that $E^{*}_{0}=(0,0)$ and $E^{*}_{2}=(M_{2}^{*},F_{2}^{*})$ are locally asymptotically stable. It is likely that this case presents bistability and that $E^{*}_{1}=(M_{1}^{*},F_{1}^{*})$ is unstable, with the basin of attraction of $E^{*}_{0}$ containing the interval $[0,E_{1}^{*}):=\{(M,F)\in\mathbb{R}_{+}^{2}\ :\ 0\leq M<M^{*}_{1},\ 0\leq F<F^{*}_{1}\}$, and the basin of attraction of $E_{2}^{*}$ containing the interval $(E_{2}^{*},\mathbf{\infty}):=\{(M,F)\in\mathbb{R}_{+}^{2}\ :\ M>M^{*}_{2},\ F>F^{*}_{2}\}$. This is at least what is suggested by the vector field illustrating this situation presented in Figure 2. It is worth noting that when $\Lambda\to\Lambda^{crit}$ from below, we have $E^{*}_{1}\to E^{*}_{2}$ and the two positive equilibria merge.

4 Elimination with periodic impulsive releases of sterile insects

We now consider periodic impulsive releases $\Lambda(t)$, modeled by the following variant of system (6):

[TABLE]

We choose in this section $\Lambda_{n}$ constant, and drop consequently the subindex $n$. For such release schedule, it is clear that the function $M_{S}$ converges when $t\to+\infty$ towards the periodic solution

[TABLE]

We therefore introduce now the following periodic system:

[TABLE]

Existence and uniqueness of continuously differentiable solutions of system (17) on the interval $[0,+\infty)$ may be shown by standard arguments, as well as the forward invariance of the positive orthant. Notice that the mosquito-free equilibrium $E^{*}_{0}$ previously introduced is still an equilibrium of (17). We are interested here in studying the conditions under which $E^{*}_{0}$ is globally asymptotically stable. For future use, we note that the mean value of $1/M_{S}^{\text{per}}$ corresponding to (16) verifies:

[TABLE]

Theorem 5** (Sufficient condition for elimination by periodic impulses).**

For any given $\tau>0$, assume that $\Lambda$ is chosen such that

[TABLE]

Then every solution of system (17) converges globally exponentially to the mosquito-free equilibrium $E^{*}_{0}$.

The previous result provides a simple sufficient condition for stabilization of the mosquito-free equilibrium, through an adequate choice of the amplitude of the releases, $\Lambda$, for given period $\tau$.

Remark 2**.**

When $r=1-r$ and ${\cal N}_{F}>{\cal N}_{M}$ (which is the case of the application we are interested in), the expression of $\Lambda_{per}^{crit}$ simplifies as follows:

[TABLE]

The function $\tau\mapsto\dfrac{2\cosh\left(\mu_{S}\tau\right)-1}{\mu_{S}\tau^{2}}$ is increasing and tends towards $\mu_{S}$ when $\tau\to 0$. Making $\tau\to 0_{+}$, we derive the following sufficient condition for stabilization:

[TABLE]

to be compared to $\Lambda^{crit}=2\dfrac{\mu_{S}}{\beta\gamma}\frac{\phi^{crit}({\cal N}_{F})}{1+\frac{{\cal N}_{F}}{{\cal N}_{M}}}$ (see Theorem 3).

Proof of Theorem 5.

First rewrite (17) as

[TABLE]

in order to emphasize the factorization of $M$ and $F$.

$\bullet$ 1. Notice that, for any $M,F\geq 0$ and any $t\geq 0$,

[TABLE]

where we write for simplicity

[TABLE]

One then deduces from (20b) that, for any $n\in\mathbb{N}$,

[TABLE]

Therefore, the sequence \big{\{}F(n\tau)\big{\}}_{n\in\mathbb{N}} decreases towards [math], provided that

[TABLE]

that is,

[TABLE]

This is sufficient to ensure that $F$ converges towards [math], and this induces the same behavior for $M$: condition (23) implies that $E^{*}_{0}$ is GAS.

$\bullet$ 2. The same argument may be conducted from (20a) rather than (20b), leading to:

[TABLE]

Global asymptotic stability is thereby guaranteed if

[TABLE]

$\bullet$ 3. Define the positive definite function

[TABLE]

and write its derivative along the trajectories of (17) as

[TABLE]

On the one hand, we have

[TABLE]

On the other hand,

[TABLE]

Coming back to (27), we deduce that

[TABLE]

One may conclude that $E^{*}_{0}$ is GAS provided that

[TABLE]

that is,

[TABLE]

$\bullet$ 4. Finally, putting together the sufficient conditions in (23), (25) and (28) yields the following sufficient condition for global asymptotic stability of $E^{*}_{0}$:

[TABLE]

Expressing the mean value as a function of $\Lambda$ with the help of (18), one establishes that $E^{*}_{0}$ is GAS if

[TABLE]

which is exactly the formula (19). This concludes the proof of Theorem 5. ∎

Remark 3**.**

A rough upper bound estimate for $\Lambda^{crit}_{per}$ can be obtained using the result from the constant continuous release case: if $\Lambda$ is chosen such that $\Lambda>\Lambda^{crit}:=2\dfrac{\mu_{S}}{\beta\gamma}\dfrac{\phi^{crit}({\cal N}_{F})}{1+\frac{{\cal N}_{F}}{{\cal N}_{M}}}$, then $E^{*}_{0}$ is GAS. Thus, using a comparison principle, a sufficient condition to ensure global asymptotic stability of $E^{*}_{0}$ is to choose

[TABLE]

where $\underline{M}_{S}^{per}=\min\limits_{t\in[0,\tau]}M_{S}^{per}(t)=\tau\Lambda\dfrac{e^{-\mu_{S}\tau}}{1-e^{-\mu_{S}\tau}}$. Thus, we derive that, for a given $\tau$, if

[TABLE]

then $E^{*}_{0}$ is GAS. When $\tau\rightarrow 0^{+}$, we recover the result for the constant continuous release (cf. Theorem 3).

5 Elimination by feedback control

We now assume that measurements are available, providing real time estimates of the number of wild males and females $M(t),F(t)$, at least for any $t=n\tau,n\in\mathbb{N}$. One thus has the possibility to choose the number $\tau\Lambda_{n}$ of mosquitoes released at time $n\tau$ in view of this information: this is a closed-loop control option. We study in the sequel this strategy.

5.1 Principle of the method

The principle of the stabilization method that we introduce now is based on two steps. The first one (Section 5.1.1) consists in solving the stabilization problem under the hypothesis that one can directly actuate on $M_{S}$. The second one (Section 5.1.2) consists in showing how to realize, through adequate choice of $\Lambda_{n}$, the prescribed behavior of $M_{S}$ defined in Step 1. The formal statement and proof are provided later, in Section 5.2.

5.1.1 Step 1 – Setting directly the sterile population level

We first suppose to be capable of directly controlling the quantity $M_{S}$. We will rely on the following key property.

Proposition 1**.**

Let $k$ be a real number such that

[TABLE]

Assume that

[TABLE]

Then every solution of (6a)-(6b) converges exponentially to $E_{0}^{*}$.

The idea behind formula (30) is quite natural: it suffice to impose a fixed upper bound $k$ on the ratio $\dfrac{M}{M+\gamma M^{*}_{S}}$ in order to make the ‘apparent’ basic offspring number $k{\cal N}_{F}$ smaller than 1, and consequently to render inviable the wild population. Notice that this condition corresponds exactly to the stability of the system linearized around the origin. It may be excessively demanding for large population sizes, as it ignores the effects of competition modeled by the exponential term. We shall come back to this point in Section 6 and introduce saturation.

Proof of Proposition 1.

From equations (6a) and (6b) we have

[TABLE]

The linear autonomous system

[TABLE]

is monotone (it involves a Metzler matrix) and may thus serve as a comparison system for the evolution of (6a)-(6b). Thus, it is deduced that

[TABLE]

where $(M^{\prime},F^{\prime})$ is the solution of (33) generated by the same initial values as the underlying solution $(M,F)$ of (6a)-(6b).

On the other hand, system (33) is asymptotically stable when (30) holds. In other words, $M^{\prime}(t)$ and $F^{\prime}(t)$ converge to $E_{0}^{*}$ asymptotically. In consequence, $M(t)$ and $F(t)$ also converge to $E_{0}^{*}$ asymptotically when (30) is in force. This achieves the proof of Proposition 1. ∎

5.1.2 Step 2 – Shaping an impulsive control compliant with Step 1

We now want to ensure that condition (31) is fulfilled, through an adequate choice of the impulse amplitude $\Lambda_{n}$. In virtue of (15c)-(15d), the value of $M_{S}$ on the interval \big{(}n\tau,(n+1)\tau\big{]} is given by

[TABLE]

and we would like to choose $\Lambda_{n}$ in such a way that (31) stays in force. However, instead of computing the (nonlinear) evolution of $M(t)$ on the interval \big{(}n\tau,(n+1)\tau\big{)}, we will impose, rather than (31), the stronger condition

[TABLE]

where $M^{\prime}(t)$ refers to the super-solution of $M(t)$ introduced in the proof of Proposition 1. (Notice that the conservatism introduced in this step remains reasonable when the original nonlinear system evolves in region where $\beta(M+F)\ll 1$.) Thereby, we can solve (33) explicitly on \big{(}n\tau,(n+1)\tau\big{]} using the following result.

Lemma 2**.**

The solution of system (33) on \big{(}n\tau,(n+1)\tau\big{]} with initial values \big{(}M^{\prime}(n\tau),F^{\prime}(n\tau)\big{)}=\big{(}M(n\tau),F(n\tau)\big{)} is given by

[TABLE]

The proof of Lemma 2 presents no difficulty and is left to the reader.

All the components of the matrix in (36) are nonnegative provided that $\mu_{F}$, $\mu_{M}$, $\rho$ and $k$ are chosen such that $\mu_{F}-\mu_{M}-(1-r)\rho k\leq 0$. It is worthwhile to recall that $\mu_{F}\leq\mu_{M}$ (see (3), page 3); therefore, the former condition is verified for any positive $\rho$ and $k$.

We now come back to the control synthesis. Using (34) and (36), condition (35) is fulfilled provided that, on any interval \big{(}n\tau,(n+1)\tau\big{]},

[TABLE]

This condition is equivalent to

[TABLE]

for any $s\in[0,\tau]$. In virtue of the relationships (3) and (7), the right-hand side of previou inequality (38) is increasing in $s$. Therefore, condition (38) has to be checked only for $s=\tau$.

5.2 Stabilization result

5.2.1 Synchronized measurements and releases

We now state and prove the stabilization result suggested by the previous considerations.

Theorem 6** (Sufficient condition for stabilization by impulsive feedback control).**

*Assume that nonnegative releases $\Lambda_{n}$, $n\in\mathbb{N}$, are used in accordance with the following constraint: *

[TABLE]

If moreover

[TABLE]

then the series of impulses $\sum\limits_{n=0}^{+\infty}\Lambda_{n}$ converges.

Implementing of the previous control law necessitates the measurement of $M(n\tau),F(n\tau)$ (or their upper estimates), and of $M_{S}(t)$ (or its lower estimate). A possibility to have (39a) fulfilled, is to ignore the population of sterile males already present at time $n\tau$ and to take simply the linear control law

[TABLE]

Notice that this expression corresponds to the value in the right-hand side of (39c).

Proof of Theorem 6.

When \big{(}M(n\tau),F(n\tau)\big{)}=(0,0), an impulsion $\Lambda_{n}$ has no effect on the evolution of $(M,F)$: the origin is an equilibrium point of system (15). We now consider the case \big{(}M(n\tau),F(n\tau)\big{)}\neq(0,0).

$\bullet$ 1. Assume first that (39a) is fulfilled with a strict inequality. By construction, one has:

[TABLE]

where $(M^{\prime},F^{\prime})$ stands for solution of (33) departing from \big{(}M(n\tau),F(n\tau)\big{)} at time $n\tau$.

We will first establish that this implies:

[TABLE]

For this, let $t_{0}$ be any element of \big{[}n\tau,(n+1)\tau\big{)} such that $M(t_{0})\leq M^{\prime}(t_{0})$, $F(t_{0})\leq F^{\prime}(t_{0})$ with at least one equality. Let us show the existence of $t_{1}$ such that $t_{0}<t_{1}<(n+1)\tau$ and

[TABLE]

Indeed, due to (40) and by definition of $t_{0}$, one has

[TABLE]

where we write by convention $M_{S}(t_{0}):=M_{S}(n\tau^{+})$ when $t_{0}=n\tau$. By continuity of the functions $M(t)$ and $M_{S}(t)$ on the open interval \big{(}n\tau,(n+1)\tau\big{)}, there thus exists $t_{1}$ such that $t_{0}<t_{1}<(n+1)\tau$ and

[TABLE]

In such conditions, it can be shown as in Proposition 1 that \big{(}M^{\prime}(t),F^{\prime}(t)\big{)}\geq\big{(}M(t),F(t)\big{)} for any $t\in(t_{0},t_{1})$, and even that \big{(}M^{\prime}(t),F^{\prime}(t)\big{)}>\big{(}M(t),F(t)\big{)}, because the functions defining the right-hand sides of (15a) and (15b) take on strictly smaller values than those defining the right-hand sides of (33). Therefore, for any t_{0}\in\big{\{}n\tau^{+}\big{\}}\cup\big{(}n\tau,(n+1)\tau\big{)}, there exists $t_{1}>t_{0}$ such that (42) holds.

From (42) and the fact that \big{(}M(n\tau),F(n\tau)\big{)}=\big{(}M^{\prime}(n\tau),F^{\prime}(n\tau)\big{)}, one deduces that (42) is true for $t_{1}=(n+1)\tau$, and therefore that (41) is true. Finally, putting together (40) and (41) yields the following key property:

[TABLE]

$\bullet$ 2. Assume now that (39a) is fulfilled (with the original non-strict inequality). Considering values of $\Lambda_{n}$ converging from above towards the quantity in the right-hand side of this inequality and relying on the continuity of the flow with respect to $\Lambda_{n}$, yields instead of (43) the non-strict inequality:

[TABLE]

$\bullet$ 3. Consider now the positive semidefinite function

[TABLE]

In view of (44), we have that for any t\in\big{(}n\tau,(n+1)\tau\big{]} it holds that

[TABLE]

Therefore,

[TABLE]

Due to (39b), there exists $\varepsilon>0$ such that

[TABLE]

and then $\dot{F}\leq-\varepsilon F$. This property ensures that $F(t)$ decreases with time, and converges exponentially towards 0. It is then deduced from (15a) that $M(t)$ also converges exponentially towards 0: overall, $(M(t),F(t))$ converges towards $E_{0}^{*}$.

$\bullet$ 4. Last, choose now $\Lambda_{n}$ fulfilling (39a) and (39c). From the property of exponential stability previously demonstrated, there exist $C,\varepsilon>0$ such that $M(t)<Ce^{-\varepsilon t}$ and $F(t)<Ce^{-\varepsilon t}$ for any $t\geq 0$. We then deduce that

[TABLE]

and one gets by summation

[TABLE]

This shows the convergence of the series and concludes the proof of Theorem 6.

∎

Remark 4**.**

Since $\Lambda_{n}$ is chosen nonnegative, inequality (39a) practically becomes

[TABLE]

It means that the release of sterile males at time $t=n\tau$ is not (really) necessary, if the sterile males population is large enough, i.e.

[TABLE]

Using this result, one may avoid unnecessary releases, thereby reducing the overall cumulative number of released males and the underlying cost of SIT control.

5.2.2 Sparse measurements

The feedback control approach requires to assess the size of mosquito population at every $t\in\tau\mathbb{N}$. Using Mark-Release-Recapture (MRR) technique [14], it is possible to estimate (roughly) the wild population, by completing sparse measurements with a period $p\tau$ for some $p\in\mathbb{N}^{*}:=\mathbb{N}\setminus\{0\}$, and computing the sizes of $(p-1)$ intermediate releases using the last sampled information between two measurements. However, the previous operation is long and costly, and reducing the frequency of realization of this protocol is legitimate.

The following result adapts the control laws given in Theorem 6 to sparse measurements.

Theorem 7** (Stabilization by impulsive control with sparse measurements).**

Let $p\in\mathbb{N}^{*}$ and set $k$ complying with (39b). Assume that, for any $m=0,1,\dots,p-1$, nonnegative releases $\Lambda_{np+m}$, $n\in\mathbb{N}$, are accomplished in accordance with the following constraint:

[TABLE]

Then every solution of system (15) converges exponentially towards $E_{0}^{*}$, with a convergence speed bounded from below by a value independent of the initial condition.

If moreover

[TABLE]

then the series of impulses $\sum\Lambda_{n}$ converges.

Notice that Theorem 7 represents an extension of Theorem 6, recovered in the case $p=1,m=0$: in this case, (46) boils down to (39a).

Proof of Theorem 7.

The demonstration comes from a slight adaptation of the proof of Theorem 6. Indeed, it suffices to verify that, under the conditions in Theorem 7, property (37) holds on the interval $(np\tau,(n+1)p\tau]$, of length $p\tau$. Let $m\in\{0,1,\dots,p-1\}$. One has for any $s\in(0,\tau]$ that

[TABLE]

Inequality (37) is thus true on \big{(}(np+m)\tau,(np+m+1)\tau\big{]} if and only if it is imposed that, for any $m\in\{0,1,\dots,p-1\}$ and any $s\in(0,\tau]$,

[TABLE]

that is,

[TABLE]

In virtue of the relationships (3) and (7), the right-hand side of (5.2.2) is an increasing function of $s$. Therefore, (5.2.2) is more restrictive when taken at $s=\tau$. This yields (46) and shows the first part of the result. The convergence of the series of impulses is demonstrated similarly to Theorem 6. ∎

Remark 5**.**

Again, if $\Lambda_{n}$ is chosen nonnegative and using the fact that the size of releases is constant between two MRR experiments, inequality (46) leads to the following practical choice of $\Lambda_{np+m}$:

[TABLE]

6 Mixed control strategies

The results obtained in the previous sections for open-loop and closed-loop SIT control allow us to compare several SIT release strategies. Here, we consider only periodic impulsive control, which is more realistic than continuous control.

The open-loop approach (developed in Section 4), is based on the determination of a sufficient size of sterile males to be released, in order to eradicate the wild population. This choice is made according to (19). Under this approach, even though the previous formula is ‘tight’, the same amount of sterile insects is used during the whole release campaign.

On the contrary, the closed-loop control approach (exposed in Section 5) is based on estimates of the wild population and thereby it enables fitting the release sizes. As evidenced by (39a), under this approach the released volume is essentially chosen as proportional to the measured population. However, this condition is certainly too demanding for large values of $M,F$ (see the comments preceding Lemma 2). Taking advantage of the apparent complementarity of the two approaches, we propose here mixed impulsive control strategies, combining the two previous modes. They gather the advantages of both approaches, guaranteeing convergence to the mosquito-free equilibrium with releases that remain bounded (like the periodic impulsive control strategies, Section 4) and vanishing with the wild population (like the feedback control strategies, Section 5).

Theorem 8**.**

Let $p\in\mathbb{N}^{*}$. Assume that, for any $n\in\mathbb{N}$, $\Lambda_{n}$ is chosen at least equal to the smallest of the right-hand side of (46) and of a positive constant $\bar{\Lambda}$. Then every solution of system (15) converges globally exponentially to $E_{0}^{*}$, in any of the following situations:

$\bullet$* Case 1.*

[TABLE]

$\bullet$* Case 2.*

[TABLE]

The interest of the previous result is of course to consider the smallest of the two values of $\bar{\Lambda}$: it results in saturated control laws.

The main issue of the proof (presented below) is to overcome the possible occurrence of infinitely many switches between the two modes. The demonstration is based on the use of common Lyapunov functions, that decrease along the trajectories of the system, regardless of the mode in use. Different Lyapunov functions are required for the two cases.

Remark 6**.**

Notice that

[TABLE]

so the condition on $k$ contained in (51) is more restrictive than the one in (50).

Remark 7**.**

The values of $\bar{\Lambda}$ that appear in (50) and (51) are two of the three that appear in (19), corresponding to (23) and (28) in the proof of Theorem 5, page 5. See the proof for more explanations.

Proof of Theorem 8.

We consider here the case where $p=1$, and the case with $p>1$ is treated in a similar way.

$\bullet$ 1. For the Case 1, consider the positive semidefinite function $V(M,F):=F$ introduced earlier by (45), page 45. For the sake of simplicity, we will denote V(t):=V\big{(}M(t),F(t)\big{)}. As shown in the proof of Theorem 5, item 1, it holds that

[TABLE]

for a certain $\varepsilon>0$ (independent of $n$) when $\Lambda_{n}$ is at least equal to $\bar{\Lambda}$ given in (50). On the other hand, it is shown in the proof of Theorem 6, item 3, that $V$ decreases exponentially when $\Lambda_{n}$ is chosen according to (39a) (which is (46) in the case $p=1$). Therefore, regardless of the mode commutations, $V(t)$, and thus $F(t)$, converges exponentially towards zero for every trajectory. As substantiated in the proof of Theorem 6, this is sufficient to deduce the convergence of $M(t)$ towards zero. Thereby, Theorem 8 is proved in the Case 1.

$\bullet$ 2. For the Case 2, let $\mathcal{V}$ be the positive definite function $\mathcal{V}(M,F):=\frac{1}{2}(M^{2}+F^{2})$ introduced by (26), page 26. It was shown in the proof of Theorem 5, item 3, that property (53) also holds for some $\varepsilon>0$ when $\Lambda_{n}$ is chosen according to (39a).

On the other hand, when $\Lambda_{n}$ is taken smaller than the value in (19), due to Theorem 6, one has for all t\in\big{(}n\tau,(n+1)\tau)\big{]}, see (37), that

[TABLE]

Therefore, on the same interval, it holds:

[TABLE]

The reduced discriminant of the previous quadratic form is

[TABLE]

which is negative when $k$ is taken according to (51). In such case, $\dot{\mathcal{V}}$ is negative definite. One concludes that $\mathcal{V}$ decreases exponentially to zero, and this ensures the global exponential stability of the mosquito-free equilibrium $E_{0}^{*}$. The result is thus also proved in the Case 2. This achieves the proof of Theorem 8. ∎

7 Numerical illustrations

The values of the vital characteristics of the mosquitoes which are used in the simulations are summarized in Table 2.

With the above numbers, we have here for the global competition coefficient $\beta=\dfrac{\sigma}{K}=3.57\times 10^{-4}$, and for the basic offspring numbers ${\cal N}_{F}\approx 75.83$ and ${\cal N}_{M}\approx 56.87$. At equilibrium, the mosquito population is thus $E^{*}=(M^{*},F^{*})$ with $M^{*}\approx 6,925$, $F^{*}\approx 5,194$ individuals per hectare.

We now present simulations related to the different impulsive strategies developed in the previous sections.

For open-loop periodic impulsive releases carried out every $7$ (resp. $14$) days, we consider the optimal value given in (19), page 19, to estimate the number of sterile males to release, that is, $7\times 1,573=11,011$ (resp. $14\times 1,604=22,456$) sterile males per hectare and per week (resp. every two weeks).

The simulations run as long as $F(t)$ is greater than a threshold value, here $10^{-1}$, below which we assume that elimination has been reached.

The corresponding simulations are given in Figure 3. In Table 3, we summarize the cumulative number of sterile males as well as the number of releases needed to reach nearly “elimination”. While, as expected, the total number of released sterile males is lower for $\tau=7$, there is no gain in terms of treatment duration. Thus, taking into account the cost of each release and also the risk of failure during the transport, it seems preferable to consider the lower number of releases, and thus to choose $\tau=14$.

The closed-loop approach can be used to reduce the cumulative number of released sterile insects and the number of effective releases. Further on, we consider several sub-cases.

We first consider measurements of the wild population every $\tau$ days or every $p\tau$ days for prescribed $p$ (here typically $p=4$). Also, we take several values for the gain: smaller $k$ provides faster convergence towards $E^{*}_{0}$ – at the price of large peak values of $\Lambda_{n}$ –, while the convergence slows down as $k$ approaches $1/{\cal N}_{F}$ with moderate values of $\Lambda_{n}$. We will consider for practical applications two particular values of $k$, namely

[TABLE]

The size $\Lambda_{n}$ of the $n$-th release is taken equal to the right-hand side of formula (39a) for $p=1$ (of (46) for $p=4$): if, at the moment of the estimate, the size of the sterile male population is sufficiently large, $\Lambda_{n}$ may be null or small.

Simulations presented in Figures 4 (page 4) and 5 (page 5) clearly show that the choice of $k$ and $p$, as well as the period $\tau$ of the releases play an important role in the convergence of the wild population to $E^{*}_{0}$. Tables 4 and 5 provide the total cumulative number of released sterile males, the number of weeks of SIT treatment needed to reach elimination, and the number of effective (that is nonzero) releases. For instance, when $(\tau,p)=(14,4)$ and $k=\dfrac{0.2}{{\cal N}_{F}}$ is relatively small, elimination of wild mosquitoes can be achieved in $56$ weeks, with only $17$ effective releases, as shown in Fig. 6(b), page 6. However, this option requires to release significant number of sterile insects per hectare (close to $2.9\times 10^{6}$ for the whole treatment).

For the larger $k=\dfrac{0.99}{{\cal N}_{F}}$ and with $(\tau,p)=(7,1)$ (see Figure 5(a)), the convergence is slower: more than $240$ weeks of SIT treatment are required to reach nearly elimination. For $p=4$ (see Figure 5(b)), the wild population is close to extinction after $58$ weeks of SIT treatment. However, based on Table 5, it seems that the choice $(\tau,p)=(14,4)$ leads to the best result in terms of timing (62 weeks) and also in terms of cumulative size encompassing 20 effective releases.

The parameter $k$ is of main importance: when $p=4$, while the number of weeks to reach elimination is quite similar for both values of $\tau$, the cumulative number of released sterile males is clearly smaller when $k$ is closer to $1/{\cal N}_{F}$.

We now consider mixed control strategies as exposed in Section 6. In Figures 7 and 8 (pages 7 and 8, respectively) we derive the simulations with the two underlying values of $k$ given in (56).

Except for the case with $(\tau,p)=(7,1)$ and $k=\dfrac{0.99}{{\cal N}_{F}}$ (see Table 6, page 6), where the convergence to $E^{*}_{0}$ is slow, it turns out that the mixed open/closed-loop control strategies derive the best results, not only in terms of releases number but also in terms of overall cumulative number of sterile males to be released.

According to Tables 6 and 7 (pages 6 and 7, respectively), for both values of $k$, the optimal solution would be to release sterile insects every 2 weeks with population assessments carried out by MRR experiments every 4 weeks ($p=4$). In addition, and thanks to (49), Figure 9 displays the release sizes $\Lambda_{n}$ for each mixed strategy. It clearly shows that during the first releases, $\Lambda_{n}=\Lambda^{crit}_{per}$. Further, when the wild population drops below a certain threshold, the feedback control occurs or not, depending on the (estimated) size of the sterile male population. That is why in Tables 6 and 7, we derive the number of effective releases (only for the case $p=4$) which confirms that the best combination is $(\tau,p)=(14,4)$, regardless of the value of $k$.

For the mixed open/closed-loop periodic impulsive SIT control, the choice of $k$ does not matter compared to the closed-loop control only. Our preliminary results thereby indicate that a mixed SIT control option with $(\tau=14,p=4)$ leads to the best strategy in terms of the total number of released sterile males and also in terms of effective releases number.

8 Conclusion

In this work, we studied various strategies to control mosquito population using SIT: open-loop and closed-loop periodic impulsive control strategies, as well as their combination (mixed open/closed-loop strategy). For the open-loop strategy (that is usually considered during field experiments) we found the minimal number of sterile males to be released every $\tau$ days in order to reach elimination of wild mosquitoes. This number is constant and relatively low. The question of determining a stopping time for the release campaigns is not simple, but clearly of primordial importance, as premature ending ruins the preceding efforts.

On the contrary, the feedback SIT-control commences with relatively abundant releases and their amplitude steadily declines with the wild population size until fading away and vanishing when the system converges towards the desired mosquito-free state. This closed-loop control strategy requires to assess the current size of the wild population (using MRR experiments, for instance).

Finally, we proposed a mixed control strategy, combining open-loop and closed-loop strategies. This control input mode renders the best result, and turns out rather meaningful from the experimental standpoint: the control input is launched at the open-loop mode during first weeks (initial phase) and then is shifted to the closed-loop mode (final phase), once the size of wild population exhibits steady decline. With this approach, the gain in terms of release pick-value, number of nonzero releases, and overall cumulative volume is clearly visible. This is due to the fact that initial phase of control action is done at the open-loop mode, i.e. by performing releases of sterile males regardless of the current size of wild population, what induces an essential decline of the wild population before switching to the closed-loop control mode. Even considering the shown simulations in terms of cost, the mixed control seems to be definitively the best choice with a release every two weeks and a population estimate every four weeks.

Knowledge of the cost of each stage of the SIT control (mass rearing, sterilization either by irradiation or using Wolbachia, transportation to the target locality, wild population measurements with MRR techniques, and other necessary supplies) will allow to estimate more precisely and optimize the treatment cost, and thus to make the most appropriate choices from an economical point of view.

As a last remark, we notice that, from a mathematical point of view, the use of closed-loop methods, as well as the fact that the proof of their effectiveness is based on argument of monotonicity, are certainly able to guarantee robustness of the proposed closed-loop algorithms with respect to several uncertainties present in the problem under study. In particular, it is believed that the framework developed here could most certainly be extended to consider the effects of modeling and measurement errors, as well as imprecision and delay in the control-loop.

Acknowledgments

Support from the Colciencias – ECOS-Nord Program (Colombia: Project CI-71089; France: Project C17M01) is kindly acknowledged. DC and OV were supported by the inter-institutional cooperation program MathAmsud (18-MATH-05). This study was part of the Phase 2A ‘SIT feasibility project against Aedes albopictus in Reunion Island’, jointly funded by the French Ministry of Health (Convention 3800/TIS) and the European Regional Development Fund (ERDF) (Convention No. 2012-32122 and Convention No. GURDTI 2017-0583-0001899). YD was (partially) supported by the DST/NRF SARChI Chair in Mathematical Models and Methods in Biosciences and Bioengineering at the University of Pretoria (grant 82770). YD also acknowledges the support of the Visiting Professor Granting Scheme from the Office of the Deputy Vice-Chancellor for Research Office of the University of Pretoria.

Appendix: Proof of Lemma 1, page 1

First, it is easy to check that function $f$, defined in (14), page 14, is first decreasing and then increasing, and, thus, may solely have no root, one root or two zeros.

On the other hand, the number of roots of $f(x)=0$ is clearly non-increasing with respect to $a>0$: it has two roots for ‘small values’ of $a$, no root for ‘large values’ of $a$, and exactly one root for a certain critical value $a^{crit}$ separating the two previous regions. This critical value is characterized by the fact that it possesses a double root $x^{crit}>0$, such that $f(x^{crit})=f^{\prime}(x^{crit})=0$, that is:

[TABLE]

Eliminating the exponential term in the previous formulas yields the second-order polynomial equation in $\dfrac{1}{x^{crit}}$

[TABLE]

Its unique positive root is

[TABLE]

Introducing this expression back in (A-1), leads to

[TABLE]

Thus $\phi^{crit}:=\dfrac{a^{crit}c}{2}$ is solution of (10), page 10, with $\phi^{crit}:=\dfrac{a^{crit}c}{2}$, such that, at the critical point, the parameters $a^{crit},b,c$ are interrelated.

For positive values of $a$ smaller than $a^{crit}$, the equation $f(x)=0$ has two roots, and no root whenever $a>a^{crit}$. This achieves the proof of Lemma 1.