Optimal barrier zones for stopping the invasion of Aedes aegypti mosquitoes via transgenic or sterile insect techniques

Abstract

Biological invasions have dramatically altered the natural world by threatening native species and their communities. Moreover, when the invading species is a vector for human disease, there are further substantive public health and economic impacts. The development of transgenic technologies is being explored in relation to new approaches for the biological control of insect pests. We investigate the use of two control strategies, classical sterile insect techniques and transgenic late-acting bisex lethality (Release of Insects carrying a Dominant Lethal), for controlling invasion of the mosquito Aedes aegypti using a spatial stage-structured mathematical model. In particular, we explore the use of a barrier zone of sterile/transgenic insects to prevent or impede the invasion of mosquitoes. We show that the level of control required is not only highly sensitive to the rate at which the sterile/transgenic males are released in the barrier zone but also to the spatial range of release. Our models characterise how the distribution of sterile/transgenic mosquitoes in the barrier zone can be controlled so as to minimise the number of mass-produced insects required for the arrest of species invasion. We predict that, given unknown rates of mosquito dispersal, management strategies should concentrate on larger release areas rather than more intense release rates for optimal control.

Appendices

Appendix 1: Numerical methods

The delayed PDE equation (4a) was solved by constructing a finite difference method. The kinetics are considered explicitly with a fully implicit treatment of diffusive transport (Morton and Mayers 1994). This numerical algorithm has been validated against the simulations of Seirin-Lee et al. (2010), and it has also been confirmed that refinements in the time steps and grid size do not influence the results presented. In simulations, we used $10^{-3}$ or $10^{-4}$ for the time step and $10^{-3}$ for the spatial step.

Appendix 2: The influence of mosquito diffusion on optimal release ratio, optimal release distance, release effort minima and wave speeds

The region $\Omega$ is chosen sufficiently large to ensure that the boundary conditions do not influence the prospects for wild-type invasion through the barrier zone. Thus, there are two parameters in the model which are changed by a rescaling of length: the diffusion coefficient and the barrier zone size. For example, we can obtain a rescaled diffusion rate, $D_{\mathrm {esti}}=D_{0}/L_{0}^{2}$, for some given diffusion rate, $D_{0}$, and length scale, $L_{0}$. For a different diffusion rate, $D=kD_{0}$, with k constant, we can rescale the spatial length from $L_{0}$ to $L=\sqrt {k}L_{0}$ which leads to the same rescaled diffusion rate, $D_{\mathrm {esti}}$, such that $D_{\mathrm {esti}}=D_{0}/L_{0}^{2}=D/L^{2}$.

Now, we redefine $[\mathrm {EF}]_{\min }$ directly for an arbitrary diffusion rate. From $L=\sqrt {k}L_{0}$ and $D=kD_{0}$, we have

$$ L(\equiv L(D))=\sqrt{k}L_{0}=\sqrt{\frac{D}{D_{0}}}L_{0}.$$

(9)

This shows that an optimal release region size can be described by a minimising value of the release region size, $\bar \gamma$, at the minimum release effort for the diffusion rate, $D_{0}$, such that

$$ \gamma_{\mathrm{s}}^{\mathrm{opt}}=\sqrt{\frac{D}{D_{0}}}\bar\gamma.$$

(10)

Then, the release effort function for the arbitrary diffusion coefficient D is given by

$$ [\mathrm{EF}]_{\min}(D)= \gamma_{\mathrm{s}}^{\mathrm{opt}}\theta^{\mathrm{opt}} N^*=\sqrt{\frac{D}{D_{0}}}\bar\gamma \bar\theta N^*,$$

(11)

where $\bar \theta$ is the minimising value of the release rate ratio at the minimum release effort, $[\mathrm {EF}]_{\min }(D_{0})$, and is not affected by the spatial length scaling so that $\theta ^{\mathrm {opt}}=\bar \theta$. From the numerical results of Fig. 3, we know the values of $D_{0}=1~\text {km}^{2}/\text {day}$ and the detailed values of $(\bar \gamma , \bar \theta )$, so that we obtain Fig. 4 directly from noting

$$ \text{[EF]}_{\min}(D)=\bar\gamma_{\mathrm{s}} \bar\theta N^*\sqrt{D},$$

(12)

with D measured in units of square kilometre per day.

Invasion speed and the diffusion rate

Analogous observations apply for the invasive wave speed of pest insects in the absence of control. In particular, noting independence with respect to the overall domain size, the only length scale in the current context within the model occurs in the diffusion coefficient. This is of degree two with respect to powers of the length dimension whereas the speed of invasion is of degree one; thus, the invasion speed must scale linearly with the square root of the diffusion coefficient. Any other relation will conflict with the need for dimensional consistency.

Appendix 3: The release effort function

The sterile/transgenic males are released locally in space and continuously in time with constant rate, $\kappa$; this is described by the release function

$$\begin{array}{rll}\kappa(x)&=& \theta N^{*} \chi(x),\\\chi(x) &=& \left\{\begin{array}{cl}1 & x \in A=[\bar x-\gamma_{\mathrm{s}}/2, \bar x+\gamma_{\mathrm{s}}/2] \\0 & x\in \Omega \backslash A\end{array}\right. ,\end{array}$$

(13)

where $\gamma _{\mathrm {s}} \ll |\Omega |$ is the release region size and $\bar {x}$ is the release region centre.

Here, we show that the release effort function, defined by $[\mathrm {EF}]=\gamma _{\mathrm {s}}\theta N^*$, is also proportional to the number of control males present in the environment after the initial transient dynamics, in particular,

$$ [\mathrm{EF}]=\mu\lim_{t \to \infty}\int_{\Omega}S \,\mathrm{d}x.$$

(14)

Integrating Eq. 4b,

$$\begin{array}{rll} \frac{\partial{S(x, t)}}{\partial{t}} &=&D_{S}\frac{\partial^{2}{S(x, t)}}{\partial{x}^{2}} \\ && + \kappa(x)-\mu S(x, t), x\in\Omega, t\in(0, +\infty], \end{array}$$

(15)

over the spatial variable, with use of Eq. 13, we have

$$ \int_{\Omega}\frac{\partial S}{\partial t}\,\mathrm{d}x=\int_{\Omega}D_{S} \frac{\partial^{2} S}{\partial x^{2}}\,\mathrm{d}x +\int_{\Omega}[\kappa(x)-\mu S]\,\mathrm{d}x,$$

(16)

and thus

$$ \frac{\text{d}} {\text{d}{t}}\int_{\Omega}S\,\mathrm{d}x=\theta N^* \gamma_{\mathrm{s}}-\mu \int_{\Omega}S\,\mathrm{d}x.$$

(17)

This is a linear ordinary differential equation in $\int _{\Omega }S\,\mathrm {d}x$, which can be readily be solved and, for general initial conditions, this yields

$$\lim\limits_{t \to \infty}\int_{\Omega}S\,\mathrm{d}x=\frac{\theta N^* \gamma_{\mathrm{s}}}{\mu}=\frac{1}{\mu}[\mathrm{EF}].$$

Hence, we obtain the required result, Eq. 14.

Appendix 4: Equilibrium relativity on varying r and β

In Fig. 6, we vary the birth rate and density-dependent coefficient for four different release rates (rows) and for the two control strategies (columns) and calculate the equilibrium relative to the no-strategy equilibrium in the non-spatial model. As one might expect, lower values of density dependence and birth rate lead to a lower required release rate for population control. Also, as one can see, the transition across parameter space is relatively smooth, with no counterintuitive behaviour for intermediate values. Similar qualitative results may be obtained for the spatial model, but we omit these for the sake of brevity.

Fig. 6

The equilibrium relative to the no-strategy equilibrium in the non-spatial model when varying the birth rate, r, and density-dependent coefficient, $\beta$, for each of four different release rates (rows) and two control strategies, SIT vs. RIDL (columns)