Modeling population dynamics based on experimental trials with genetically modified (RIDL) mosquitoes
Mario A. Natiello Hernán G. Solari






















lyxlist[1]
	
		#1
		 
		 
		 
	
















Modelling population dynamics based on experimental trials with genetically
modified (RIDL) mosquitoes
    Mario A. Natiello^† and Hernán G. Solari^
    December 30, 2023
======================================================================================================



^†Centre for Mathematical Sciences, Lund University, Box
118, S 221 00 LUND, Sweden. mario.natiello@math.lth.se (corresponding
author)

^Departamento de Física, FCEN-UBA and IFIBA-CONICET,
Argentina. solari@df.uba.ar

Recently, the RIDL-SIT technology has been field-tested for control
of Aedes aegypti. The technique consists of releasing genetically
modified mosquitoes carrying a “lethal gene”. In 2016 the World
Health Organization (WHO) and the Pan-American Health Organization
(PAHO) recommended to their constituent countries to test the new
technologies proposed to control Aedes aegypti populations.
However, issues concerning effectiveness and ecological impact have
not been thoroughly studied so far. In order to study these issues,
we develop an ecological model. It presents interdependent dynamics
of mosquito populations and food in a homogeneous setting. Mosquito
populations are described in a stochastic compartmental setup, in
terms of reaction norms depending on the available food in the environment.
The development of the model allows us to indicate some critical biological
knowledge that is missing and could (should) be produced. Hybridisation
levels, release numbers during and after intervention and population
recovery time after the intervention as a function of intervention
duration and target are calculated under different hypotheses with
regard to the fitness of hybrids and compared with two field studies
of actual interventions. This minimal model should serve as a basis
for detailed models when the necessary information to construct them
is produced. For the time being, the model shows that nature will
not clean non-lethal introgressed genes.

Keywords: environmental impact, RIDL-SIT, compartmental model, epidemic
risk, genetic modification, Aedes aegypti



§ INTRODUCTION


The use of the Sterile Insect Technique (SIT) has been proposed to
control populations of the mosquito Aedes aegypti (Ae ae.)
because of public health reasons <cit.>. Ae ae. is
the main vector of several viruses such as the Flavivirus producing
Yellow Fever, Dengue and Zika, as well as the Alphavirus producing
Chikungunya. The idea is not new and follows the success obtained
with the fruit fly (Drosophila sp.) in agricultural settings.

The early attempts to use the method with Ae ae. sterilised
by radiation ended in failure <cit.>. The evolution of SIT
by the introduction of genetically modified insects, GMI, incorporating
an autocidal gene <cit.> (termed RIDL-SIT, Release of Insects
carrying a Dominant Lethal) was later considered to be a more promising
technique for mosquito control than irradiated mosquitoes <cit.>.
The RIDL-SIT technique for Ae ae. was developed recently and
communicated in a series of papers <cit.>
(all references except the first one refer to actual field-tests).

The World Health Organization (WHO) recently promoted the realisation
of evaluations of this technique <cit.> in the wake of the
2015 Zika pandemic in the Americas. More recently, the Food and Drug
Administration (FDA-USA) authorised the realisation of such evaluations
in the US <cit.>. In correspondence with technological and
political developments, the Pan-American Health Organization (PAHO)
endorsed WHO's recommendations for its member states. It is apparent
that there exists a growing interest in performing tests of the RIDL-SIT
technique in the field, despite the limited amount of available
information and the uncertainty of impact assessments. Nevertheless,
warnings were raised by scientists <cit.> that foresaw the
ecological and social dimension of the problem, in correspondence
with the considerations of social scientists regarding technological
fixes <cit.>.

In April 20-21, 2017 a scientific meeting was held at the Ministry
of Public Health of Argentina with the support of PAHO, addressing
the reasonability of conducting such experiments in the country. The
participants were the National health authorities, a PAHO representative,
a group of scientists proposing integrated control technologies (including
SIT), a group of scientists not involved in technological development
and two representatives of the foremost RIDL-SIT company. Despite
the presence of a scientist representing the team developing the RIDL-SIT
technique, some questions concerning impact assessment remained unanswered,
being nevertheless relevant for the scientific community and the society
as a whole, for example:


  * Which lasting modifications of the environment will be produced by
these tests? 

  * Will the population of vectors previous to the environmental intervention
be re-established? If yes, how fast?

The goal of the present work is to address these and related issues
focusing on the technique promoted by the WHO <cit.> concerning
the transgenic mosquito OX513A <cit.>
(in the preliminary work <cit.> it is called LA513A) that
was preliminarily evaluated in <cit.>. We will use the mathematical
methods of stochastic population dynamics, which should in general
be used before any intervention. We investigate the time-evolution
of a mosquito population as a result of the proposed interventions,
focusing on: 


  * The hybridisation of the released mosquitoes with the local populations
as a result of genetic diffusion

  * The expected ratio between released males and local males necessary
to achieve different levels of control

  * The recovery time of the new, hybrid, populations when the intervention
ceases 

  * The reasons for the difficulties found in establishing such control
systems 

To achieve these goals we implement a stochastic model, able to deal
with small population numbers, based upon the available information.
Because of the scarcity of biological information concerning the released
mosquitoes, we have kept the model as general as possible, avoiding
the need to replace missing information with conjectures. We indicate
the most noticeable missing information in our critical review of
the data (see Section <ref>). It is also
worth mentioning that construction of conclusions concerning some
data available in the reports is occasionally inconsistent between
different reports and/or existing biological knowledge.

This manuscript is structured as follows: In Section <ref>
we describe the biological background of the technique, including
a critical review of the available data, as well as the mathematical
methods supporting the simulation. Section <ref> displays
the outcome of the simulations and their relation to relevant biological
questions. In Section <ref> results of this work are compared
with actual field observations, while Section <ref> contains
a general discussion and concluding remarks.



§ BIOLOGICAL BACKGROUND AND MATHEMATICAL METHODS




 §.§ The RIDL-SIT technique


The Sterile Insect Technique (SIT) was introduced in the past
century for insect pest control. The goal is to release sterile males
in sufficient number to reduce (control) the population by reducing
the number of offspring. Sterility was originally achieved by radiation
and it has been successful as part of the measures to control the
Mediterranean fruit fly <cit.> (see below), the screwworm
fly <cit.> and other insects regarded as pests (for a history
see <cit.>). Hence, the idea of controlling populations of
mosquito vectors of viral diseases by sterilisation emerged clearly.
By 1962 Aedes aegypti was considered a mosquito that could
be controlled with the technique. Field experiments were performed
by <cit.> but despite overwhelming releases of sterilised
males in relation to the local mosquito population, no control was
made evident. In consideration of this failure the technique had to
be revised. <cit.> discussed some conditions that a sterile
insect technique should satisfy to be successful; we quote:

Apart from the obvious requisite of mating competitiveness, two criteria
must be satisfied by any sterilisation method: (1) sterility must
not be suppressed by complementation with any genotype present in
field populations. Any variation in the expression of the sterility
factor(s) following interaction with the large and heterogeneous genome
of the field population will quickly lead to selection for non-sensitive
females and fertile matings. (2) Absolute assurance of effectiveness
must result from the sterilisation method in the factory setting.

These authors also recommended that 

a first step for the future deployment of transgenic technology is
to demonstrate the safety of release material in the absence of confounding
concerns about the spread of alleles and drive elements. This is possible
only if the released material is genetically sterile.

Several authors <cit.> have indicated
that density-dependent-effects can/will reduce the efficiency of SIT
programs to the point of making it non-viable.

To overcome in part these population (ecological) effects, the RIDL-SIT
technique was introduced. Transgenic mosquitoes were designed to have
a late acting dominant lethal gene, LG, <cit.> following
a technique developed for the fruit fly <cit.>. The original
strains of modified mosquitoes were constructed using the Rockefeller
strain of Ae ae., a strain that has been kept in captivity
for many years and presents adaptation of its life-cycle to the laboratory
<cit.>, as well as manifestly larger sizes
and fertility than wild strains.

This work focuses on the technique developed upon the transgenic mosquito
OX513A specified in the previous Section. The proposed strategy is
to release solely male mosquitoes outnumbering the wild males and
thus strongly influencing the fecundation of wild females. The offspring
is expected to express the lethal gene and die before reaching fertile
age. Thus, the introduced mosquitoes would, after one generation,
disappear from the natural environment along with their offspring.
The expected outcome in a closed environment (no adult immigration
or emigration) is to control the original wild mosquito population.

The genetic engineering of the mosquitoes proceeds by inoculating
eggs with molecular constructs to “induce tetracycline-repressible
dominant lethality in both males and females” <cit.>, producing
four transgenic lines of modified mosquitoes which appear to be obtained
by modifications at a single site in the genome. Three of these lines
gave the desired mortality. One of these three presented increased
mortality near pupation and was selected for mass production because
of this desired attribute. According to <cit.>, a penetrance
of 93-97% was achieved, meaning that at the level of phenotype a
3-7% of the individuals carrying the gene do not manifest the early
mortality. We found no studies concerning whether the individuals
carrying the gene and arriving to adulthood transmit the trait of
not-manifesting-the-modified-gene. For the case of the Mediterranean
fruit fly, <cit.> followed one generation (n=4590) concluding
that survival to adulthood despite the lethal gene is not an inheritable
trait. Such conclusion is inaccurate, a proper statement is that the
probability of survival being inheritable lies in the interval [0, 5.9·10^-4]
with 95% confidence. The difference matters when this probability
is multiplied by millions of released individuals. Nevertheless, the
authors acknowledge that “the possibility of biochemical resistance
to the lethal effector molecule remains a potential drawback to RIDL
that is not a significant issue for radiation-based sterilization”.

As part of the engineering process, a fluorescent marker was incorporated
to the same modified gene. This marker facilitates identification
of larvae of genetically modified mosquitoes.

The requisite of mating competitiveness <cit.> of the modified
mosquitoes, is monitored by the quantity C, defined as the ratio
between offspring (lower case) to adult males (upper case) of modified
mosquitoes, m/M, normalised to the same ratio for the wild mosquitoes,
w/W: 

    C=(m/M)(w/W)^-1.

Low estimates values of C will have to be compensated by performing
larger releases. The estimation of C performed in cage and field
experiments produced values of 0.059 (Cayman Islands <cit.>)
and 0.031 (Brazil, <cit.>). We leave aside a previous
report for the Cayman Islands of C=0.56 <cit.> obtained
in the context of trying to establish that mating of wild females
with modified males was actually possible.

An additional source of uncertainty is the lifespan of modified adults
in the wild. Mark, release and capture experiments were performed
twice. For Ae ae., the ability to produce reliable data in
such experiments has been questioned. A large dispersion of results
is obtained when trying to estimate dispersal distances, i.e., dispersal
measurements appear to be greatly influenced by the measurement procedure
<cit.>. The reported results appear not to escape
to this rule: survivals of 2 days are reported in <cit.>
and of 5 days in <cit.>. In the first experiment the unmodified
laboratory strain was released as well showing a 2.1 day survival
in average with no significant differences with respect to the modified
strain. The release of a control strain was not repeated in the second
experiment, probably because of the urban setting in which it was
performed. In contrast, the survival of modified mosquitoes tested
against a strain collected in Chiapas (Mexico) termed Latin Wild Type,
LWT, used as the basis of the FDA report <cit.>
contains only experiments in the laboratory, but none in the field.
The experiment confirms the existence of significant differences in
phenotype between LWT and the GMI in terms of fertility and hatch
rate measured at the laboratory, the environment in which the GMI
was developed. We have not found data testing survival rates of modified
mosquitoes released in the field against the local wild population
under similar conditions among the reported literature (<cit.>).



 §.§ Mosquito, genetics and environment


The populations of the tree hole mosquito Aedes aegypti are
mainly limited by the development during the larval stage. The two
main environmental factors that have been so far identified are temperature
<cit.> and food. It is expected that genetics will
play a role but comparative studies of quantitative traits are not
abundant. Studies comparing developmental rates of Ae ae. collected
at different locations in Argentina and those of the Rockefeller strain
(a laboratory strain), indicate <cit.> that there are
significant differences, under the same temperature and laboratory
feeding conditions, between local strains and the Rockefeller strain,
as well as statistically significant differences between local (country
wide) strains, despite being all of them identified in the same branch
of the phylogenetic tree as hybrids between Aedes aegypti aegypti
and Aedes aegypti formosus <cit.>. For example, the
Rockefeller strain spends less time as larvae. Yet, the most substantial
difference corresponds to the total number of eggs and the daily fertility
which is five times larger for the Rockefeller strain than for the
local strains. A recent study performed in Trinidad and Tobago <cit.>
reports a fertility of local females raised in the laboratory, with
local genetics, which is in line with the experiments reported in
<cit.>.

Any model (be it mathematically explicit or not) aimed at exploring
the possible outcome of an environmental intervention, must include
the reaction norms <cit.> of the local, the
released and the resulting hybrid mosquitoes. The traits that have
been identified so far are: hatching inhibition as a response to low
levels of food (the Gillett effect) <cit.>,
duration of the larvae stage <cit.>,
mortality of preimaginal forms <cit.> and
fertility <cit.>, directly associated to the body size of
females. Furthermore, it has been observed that temperature and food
are not independent factors <cit.>.

Needless to say, the feeding behaviour of the larvae will alter the
food resources in the container. Hence, food and phenotype must be
considered as dynamical variables.



 §.§ Simple population dynamic model to evaluate effectiveness and gene
diffusion




  §.§.§ Comments on previous models


In Table <ref> we display a few mathematical
models that have bearing on the issues discussed in this work. The
table classifies the models according to key features of our modelling
strategy. Pop denotes the number of compartments. A large number
of  compartments allows for finer graining, since compartments differentiate
groups of individuals with different characteristics, different dynamical
responses and/or different biology. Var indicates whether modelling
is performed with continuous or discrete variables. Populations are
counted with non-negative integers and are hence intrinsically discrete;
continuous variables may constitute a fair approximation to large
populations, but they are in conflict with the possibility of extinction
of one or more compartments. Dyn indicates the type of dynamics.
All models but the last in the table are deterministic (either differential
equations or difference equations on discrete time). Such strategy
forces in practice the use of continuous variables. SkeeterBuster,
the last entry of the Table has an intermediate approach, being the
time-lapse of the preimaginal stages stochastic. Eco describes
whether the model takes into account the interaction of the populations
with the environment. In fact, model parameters have a biological
meaning, at least metaphorically, but in general they express the
action of the environment on the populations, rather than a true interaction
where the dynamics of the population modifies the environment as well.
Sp indicates whether the spatial distribution is taken into account.
Homogeneous environments can be modelled without spatial resolution,
contrasting with inhomogeneous environments. Int: indicates whether
introgression is considered. GS: denotes the genetic strategy in
the model (LG=lethal gene, Ge=Generic).


Concerning the central issues of this work, these models were not
conceived to address the possibility of hybridisation and its eventual
propagation in time, or the interruption of the interventions, and
are hence insufficient to discuss the questions of interest for this
study. The model in <cit.> includes a discussion
about fitness of the RIDL released males but the context is how to
compensate deficiencies in fitness with a higher release rate. The
model proposed in <cit.> is the only one addressing the LG
strategy. It has a low level of detail since it deals only with one
population. Its goal was to address questions regarding the effectiveness
of the proposed technique. FDA (Food and Drug Administration, US)
<cit.> refers to this model in their assessment of the RIDL-SIT
technique. Finally, most models treat the environment as “fixed”.
Even if biological parameters in the models are e.g., temperature
dependent, the environment is imposed on the populations without feedback,
being in this way outside the dynamical evolution. The possibility
of (local) environmental modification originated in the time-evolution
of the populations is absent.

A separate paragraph deserves the approach of SkeeterBuster
(last entry of the table). The model produces detailed simulations
based upon Focks' model CIMSiM <cit.>. The original
model was revised and improved to incorporate the discussion of some
genetic aspects, mainly the fixation index, F_ST, in relation
to spatial dispersion. The computation is performed on one neutral
locus (two alleles). One locus accounts only for (00,01,11) classes,
meaning that wild and hybrids with no LG are considered identical.
No hybridisation can be tracked in this form. The model is highly
detailed (it tracks individual females) but its biological content
targets other questions as compared with the present work. Additionally,
the information required to run the model is difficult to obtain,
excpet for the area of Iquitos (Perú), which is the main target for
SkeeterBuster. Iquitos has been the locus of a long term research
project in dengue and Aedes aegypti, but even in such conditions,
parametrisation of the model is difficult: We read in the conclusions
of <cit.>, “CIMSiM and DENSiM are most sensitive to
parameters for which little empirical information exists, and although
this information can be extremely valuable for setting research priorities,
it also places constraints on the use and interpretation of model
outputs.” When compared to early versions of aedesBA (a compartmental
model that uses an intermediate level of information) no definite
advantage was observed for any of the models <cit.>. The
difference between modelling philosophies has been briefly discussed
in <cit.>.

Mathematical models do not produce biological data, they articulate
information to produce understanding and insight about what to expect
given the present level of information. However, knowledge at times
takes the form of knowledge of our ignorance: which information do
we need that has not yet been produced? Measurements and experiments
are made to match a preconceived form of articulation. Models call
for new experiments to produce substantially superior forms of articulation
than the irrationality in habits and intuition of (raw) empiricism.

The modelling situation presents the dilemma: should we refrain from
modelling because of the lack of information to run our most detailed
codes? Such a restriction has as consequence that models which are
too limited to address the task (e.g., <cit.>), become the
only decision-making contribution provided by the modelling community.
The best is the enemy of the good, is a popular aphorism and Italian
proverb (il meglio è l'inimico del bene) quoted by Voltaire.
In what follows we present a viable option under the current circumstances.



  §.§.§ Lack of information and modelling decisions.


In the light of the limited or uncertain information concerning the
RIDL-SIT technique, some modelling decisions were made, namely:


  * The lack of information regarding the reaction norms to environmental
conditions such as the availability of food have been decisive in
the modelling decision of not considering high-level models such as
<cit.> which incorporate phenotype variation in terms of
environmental conditions. High-level models are sensitive to differences
found in the phenotype for eggs collected at different places <cit.>.

  * The effect of the lethal gene seems to be apparent not before L4
stage and it is strongest at pupae level (see Table 1 in <cit.>).
Hence, the relation between mortality to emergence related to unmodified
pupae is adjusted so that for individuals carrying the lethal gene
the probability of death as pupae is 19/20=0.95 and the
probability of adult emergence is 1/20=0.05, since these
are the only competing events in an exponential race <cit.>
at pupa stage.

  * Another substantial missing information concerns the contributions
of dead larvae and pupae as a secondary source of food. Such contribution
was incorporated in earlier models <cit.> based upon results
in <cit.> and amount to a 12% recovery (the product of
40% nutritive value and 30% conversion of food in weight) according
to these sources. Recycling of dead larvae was assumed to influence
food availability during the developmental stages (see Section <ref>
for a description of the food cycle).

  * Heritability of the survival trait despite the presence of the lethal
gene has been set to zero.

  * Mating effectiveness as compared with the wild strain has been rounded
up to C=0.06.

  * We have made no difference in the survival of adult males of different
strains (other than the food deficit effect for individuals raised
in the wild). None of the reported experiments checks lifespan of
modified males vs. lifespan of wild males in the same wild environment
under similar conditions. In our simulations, that target for a pre-established
reduction level of the wild population, the mortality of released
males will affect the release rate since a comparatively high proportion
of males must be maintained throughout the process. Eventually, larger
release rates may enter in conflict with the factory production capabilities.
The assumption made here underestimates the strain put on the factory.

  * Since there is no information regarding reaction norms of hybrids
between local strains of Ae ae. and the Rockefeller strain,
we have initially assigned to all wild and hybrid mosquitoes the same
reaction to food deficit stress using the set of reaction norms available
for mosquitoes collected at Buenos Aires (Argentina) <cit.>
corresponding, in the present model, to: pupation rate, fertility,
mortality at the preimaginal stadiums, as well as the Gillett effect
as quantitatively described in <cit.> and the mortality of
adults. However, released adult mosquitoes (overwhelmingly males)
are not raised under food deficit stress and their mortality rate
is not corrected by food deficit but fixed as indicated in the previous
point.

  * In Subsection <ref> we explore to some extent the dependence
of hybridisation with respect to hybrid fitness in two ways: (a) We
test two values for the mating efficiency of heterozygous males carrying
one LG, O1, setting it equal to that of the wild strain or, alternatively,
to that of homozygous males (carrying two copies of the LG). (b) For
adult females, we tested doubling the mortality rate for the O1
hybrids. In all cases, there is no point in considering hybrids not
carrying a LG worse fitted than O1 mosquitoes. For hybrids not
carrying a LG (named Rock in Table <ref>)
we considered a linear interpolation between the two extreme cases
just described, modifying at the same time the mating efficiency of
males and the mortality rate of females.

A model is not just a mere transcription of experimental findings
but a logical construction as well. In our schematic construction
we have considered the two food regimes described in <cit.>
and the dependence of phenomenological traits with the logarithm of
the food concentration. A region where the life form is impossible
was considered for food concentrations below 2^-8 of the optimal
concentration values obtained in the laboratory (which correspond
to 1 in these units). The region of high mortality was considered
in the range [2^-8,2^-6) while the favourable region spans
the interval [2^-6,1]. In the most unfavourable region the only
life form possible are dormant eggs inhibited from hatching by the
Gillett effect (see Section <ref>), all
other preimaginal forms have probability one of dying. In the transition
region, the Gillett effect diminishes from total inhibition to no
inhibition and the mortality moves from probability one to the values
measured in laboratory experiments <cit.>. The fertility
is considered to be proportional to the mass excess to the minimal
form of the adult found at the 2^-8 food concentration where
it takes a zero value and grows to reach the measured value for Ae
ae. in the region of optimal laboratory conditions <cit.>.
Laboratory results indicate that the increase in size as measured
by the wing length is linear with the logarithm of the food concentration
<cit.> all along the considered region. Using the raw data
of the experiment reported in <cit.> we determined
that the event rate as larvae (considering both death and pupation)
increased linearly (with the logarithm of food concentration) in the
favourable region and was approximately constant in the transition
region.



  §.§.§ Stochastic Population Dynamics


This approach to population studies assumes that the dynamics is adequately
described by classifying the population into compartments X_j,
j=1,⋯,N (described by nonnegative integer numbers), each compartment
counting the number of individuals of a given type as a function of
time. For example, mosquito subpopulations in the present study are
represented by the class D with elements {egg, larvae,
pupae, adult female, adult male}. Each element in D will
be further subdivided with respect to its genetic type according to
G, with elements {wild, hybrid, modified heterozygous,
modified homozygous}. Thus, the index j can be thought
of as listing the elements of the external product D× G. Changes
in time occur stochastically and are identified by events,
labelled by the index α, where α=1,⋯,E. The event
α produces an instantaneous change in population compartment
j given by the integer δ_j^α. Note that events
where an individual undergoes a transition from one compartment to
another (e.g., pupation where a larva becomes pupa) have exactly two
nonzero indices δ_j^α: For j=l (larva) δ_l^α=-1,
while for j=p (pupa) δ_p^α=+1. Also, death events
have δ=-1 for the corresponding compartment. Birth events
(oviposition) may have δ>1(see Subsection <ref>).
The time-evolution of the population, given an initial condition (here
labelled as t=0) can be summarised as:

    X_j(t)=X_j(0)+∑_α=1^Eδ_j^αn_α(t)

where n_α(t) is the number of events of type α that
occurred up to time t (a stochastic variable). This approach is
called a Markov Jump Process (Markov since the only necessary
information to compute the next change is the knowledge of the present
state; Jump since the changes that take place are described
by the integer δ_j^α indicating a change in the number
of individuals in compartment j). In the mathematical literature
attention to these processes starts with Kolmogorov’s foundational
work <cit.> and its further elaboration by Feller
<cit.>. A substantial effort to relate the stochastic
description to deterministic equations was performed by Kurtz <cit.>.

Traditionally, attention was focused on the dynamics in population
space, rather than event space. In the case of probability rates that
are linear in the subpopulations both approaches are quite interchangeable<cit.>.
However, event space is more appealing for a general understanding.
For example, when the dynamics is affected by introducing a treatment
or strategy regarded as an event (think e.g., of vaccination in an
epidemic disease) it will be important for policy reasons (logistics,
costs, etc.) to count the number of occurrences of that event, rather
than attempting to sort it out by reversing the information contained
in the population set X={ X_j(t), j=1,⋯,N}.

The Markov assumption leads to defining probabilities per unit time
W_α (in the sequel called rates or probability
rates) for each event. In this work it will suffice to consider that
these probability rates depend only on the population X at each
given time (other dependencies will be handled along the way). In
a small time-interval h the probability of occurrence/non occurrence
of event α is the pair (W_α(X)h+o(h),1-W_α(X)h+o(h)),
where o(h) is a quantity that goes to zero with h faster than
linearly.

The process is further modelled implementing the Feller-Kendall algorithm
<cit.> that can be summarised in the following steps
<cit.>:


  * Starting at t=0 or immediately after an event has occurred, the
waiting time to the next event is exponentially distributed, with
rate R=∑_αW_α(X) equal to the sum of all event
probability rates, where X refers to the population values after
the most recent event has occurred (or to X(0)). The time τ
to the next event is simulated by picking an exponential random deviate
with rate R.

  * Given that an event has occurred at time τ, the probability
of it being event α is

    P(α/event at τ)=W_α(X)/∑_αW_α(X)=W_α(X)/R.

 Choosing an arbitrary order for the events, the event that took place
at τ is simulated by picking a uniformly distributed random
number Y in [0,R] and checking to which event in the order it
corresponds. For example, if ∑_α=1^KW_α(X)<Y<∑_α=1^K+1W_α(X),
we assign the occurrence to event K+1.

  * Populations are updated according to δ_j^α for the
assigned event and a new cycle 1-3 is computed until the final simulation
time is reached.




  §.§.§ Schematic description of the model


The setup in consideration is that within a limited natural population
of Ae ae., male adult mosquitoes are introduced (with an error
of 1 female every 4300 releases <cit.>). These mosquitoes
differ in two ways from the wild population. Firstly, they have two
copies of a dominant LG (not present in the wild population) that
induces massive death of the offspring before adult stage. The released
mosquitoes have been previously raised in the laboratory, from the
Rockefeller strain of Ae ae., a laboratory-adapted strain.
These mosquitoes are genetically different from the original wild
population. If the totality of the offspring does not die before reproduction,
there will be some degree of population mixing and genetically modified
individuals (relative to the original wild population) could arise.

Hence, to describe this problem a model has to be designed covering
three goals: 


  * Describe the wild population with sufficient accuracy taking care
of environmental constraints. 

  * Take care of the evolution of the LG in the offspring along generations.
A released homozygous mosquito mating with a wild female (where the
corresponding gene is not lethal) will yield heterozygous offspring
and any eventual subsequent offspring will propagate this special
gene according to the mendelian rules.

  * Take care of the mixing between the natural population and the released,
genetically different, mosquitoes (mixing of other genes apart from
the lethal one).

Within this approach, the specification of the model requires (a)
to identify the populations (compartments) participating in the dynamics
and (b) to specify the events that define the dynamics, their reach,
biological content and environmental conditions. 

As mentioned in the previous Subsection, the subpopulations involved
are chosen according to the developmental stages (egg, larva, pupa,
male adult and female adult). In the absence of released mosquitoes
the events governing the population dynamics are egg hatching and
mortality, larval pupation and mortality, pupa emergence and mortality,
fecundation (mating) and oviposition and finally adult mortality.
These events take care of modelling goal i.

All compartments are further classified according to their genetical
content. There is a compartment for the wild lineage and three for
the mixed lineage where 0, 1 or 2 copies of the LG may be
present. This population space allows for dealing with modelling goal
ii.

Finally, each subpopulation has attached a number in [0,1] representing
the average percentage of foreign genetic material. Released mosquitoes
contribute with 1, wild mosquitoes with 0 and the mixing propagates
according to the law of independent assortment of alleles <cit.>.
This takes care of modelling goal iii.

A detailed scheme of compartments and relating events is shown in
Figure <ref>.






  §.§.§ Food dynamics


Most rates, especially larval pupation and mortality rates, depend
on the availability of food. Following the model in <cit.>
we assume that the totality of oviposition sites can optimally host
l_opt larvae (grown to their potential size), that may feed unrestrictedly.
At the same time, food in the sites is produced by (environmental)
bacterial activity, it degrades at a certain rate u (also depending
on the environment, in this work u=1/days) and it is consumed
at a pace depending on the amount of larvae. Dead larvae and pupae
are to some extent recycled as available food, as mentioned in Section
<ref>. Let P_f/C_f denote
the ratio of produced to consumed food. Initially, we set its value
to the deterministic equilibrium value, an environment-dependent estimate.
We denote the leftover food (relative to consumed food) as L_f.
In conditions of food scarcity, L_f=0. The dead larvae or pupae
recycled as food is denoted as X_l. This last parameter is initially
zero and it is updated at every larval or pupal mortality event (see
below). The ratio P_f/C_f indicates the environmental
scarcity conditions: if lesser than one, there is no food available
for the larvae to reach the potential weight at pupation (this ratio
enters in several steps of the modelling procedure). The food-cycle
for each simulation time-step τ, where the initial total larva
population is L and the ideal estimated capacity of the oviposition
sites is l_opt, proceeds as follows: 

    L_f    =     max(P_f/C_f-1,0)
    P_f/C_f    =     L_f exp(-uτ)+l_opt+X_l/L
    
    X_l    =     X_lexp(-uτ)

At each time-step, the values of L and X_l are updated and
all above quantities recalculated. The influence of the food in mosquito
dynamics is described by larval event rate LE(x), mortality probability
as larvae and pupae ML(x) , hatching inhibition GL(x) and fertility
rate FT(x)

    LE(x)     =    max(0.025 , 0.25269+0.031974 x)
    
    ML(x)     =    
    1     x<-8
    max(0.033 , 0.033-0.4835 (x+6.))     x≥-8
    
    
    GL(x)     =    
    1     x≥-6
    
    (1.+0.5 (x+6.))     -8≤ x<-6
    
    0     x<-8
    
    
    FT(x)     =    max(0. , 0.127742 (x+8.))

where x=log_2(Cf/Pf). The fertility rate changes the
fertility of the emerging females, and then the average fertility
of females (see below). The values for the reaction norms correspond
to the values reported in <cit.>. In addition, the average
lifetime as adult is corrected as m_a=m_a0/1-ML(x)
obtaining in this form a reasonable correspondence between times measured
at the laboratory (see for example <cit.>) and in field experiments
(for example <cit.>). Released male mosquitoes, that had
no food deficit while raised, were set to have average lifetime m_a0
(see values below).



  §.§.§ Specific details


The theoretical basis of the model when it comes to the population
dynamics of Ae ae. refers to <cit.>.
The problem is quasi-linear <cit.> in the populations,
meaning that (a) all events that decrease a subpopulation produce
a decrease in one unit in only one population compartment and (b)
all non constant rates can be written in the form W_α=m_α^iX_i,
for the relevant population X_i. However, the coefficients m_α^i
need not be all constant. 

The events hatching and pupation are assumed to have
the same coefficient for all four compartments of the participating
subpopulation. All these events pick one individual from the original
subpopulation and passes it to the corresponding subpopulation of
the next stage. Pupa emergence is less likely for the compartments
carrying 1 or 2 copies of the LG (see below for the explicit
values). Emergence events are randomly assigned to yield 50% males
and females. Each mortality event simply decreases the corresponding
subpopulation in one unit. Except in the case of pupa mortality, where
the effect of the LG also distinguishes two cases, the mortality coefficient
is taken to be the same along each stage. The pupa mortality is much
higher for the compartments having 1 or 2 copies of the LG (see
Appendix I for details of the estimates). As indicated in Section
<ref>, adult mortality was assumed to
be equal for all strains. This modelling decision may overestimate
the size of non wild subpopulations, an effect that may be mitigated
by larger releases. 

The release was acted weekly during the treatment period, following
the implementations of this technique discussed above, that only report
weekly releases. The release goals of <cit.> were a relation
of 10:1 in adult males (later increased to 25:1) and 50%
fluorescence in the larvae from ovitrap-collected eggs. The initial
size of the release is therefore set to be 11 times larger than
the adult male population in our computations. Later release events
are increased or reduced adaptively every week, in order to adjust
the weekly proportion of eggs with LG against a target value. There
is however a maximum release size corresponding to the maximal capacity
of the production plant. In this work it is set to 100000 individuals
per release. Release is intended to concern only modified males, but
a few females are released as well, since the accuracy of the separation
technique is limited (one female every 4300 releases <cit.>).

The model simulates a closed environment, i.e., no emigration/immigration
of adult individuals to and from adjacent areas is taken into account.
Clearly, such effects could only delay or reduce any effects of the
release since immigrants, at least initially, would be wild individuals
while emigrants could be any among released, wild and mixed individuals. 

In this implementation, fecundation (mating) is not treated as a separate
event. Instead, we assume that given that the adult male population
is nonzero, all females are fecundated immediately after emergence.
Thus, each emerged female has two associated indices, corresponding
to the adult female compartment (inherited from the pupa) and the
fecundating male compartment. The latter is chosen randomly among
the four male subpopulations, weighing the size of the subpopulation
with two copies of the LG with a factor 0.06 as described in the
previous Sections. We assume hence that this reduced effectiveness
is coupled to a high proportion of non-wild genes and to the homozygous
property.

The index pair given by fecundation is used to propagate to the offspring
both the LG and the percentage of introduced genes from the released
population (according to the independent assortment of alleles <cit.>,
i.e., with the average of parental percentages). The propagation takes
place in the oviposition event that adds individuals to the egg population
(distributed among the four egg compartments according to mendelian
rules) without decreasing any other subpopulation. The amount of eggs
laid by each female at each oviposition depends on temperature and
on the size of the female, which in turn depends on the environmental
food availability. In average, wild mosquitoes collected in Argentina
and reared in an environment with optimal food availability (P_f/C_f=1)
at T=26 C, lay 2 eggs per day <cit.>. Hence, the
number of eggs per oviposition is δ_egg=2/m_ovi(26) FT(1)FT(P_f/C_f)≡ c_layFT(P_f/C_f).
Once a female has emerged, its size is fixed as well as its egg-laying
capacity, both by the value of P_f/C_f at emergence
(which varies among events according to food dynamics). Thus, at each
emergence event, δ_egg is updated as follows:

    δ_egg:=δ_eggF+c_layFT(P_f/C_f)/F+1,

where F is the total number of adult females existing previous
to the actual emergence. This is a sort of “moving average”.

Dead larvae and pupae are partially recycled by the environment within
the hatching sites. For each death event we add to X_l the amount
Δ l=xlFT(P_f/C_f)/FT(1). The quotient
estimates the reduced size of each dead larva or pupae relative to
ideal food conditions. The values xl=0, 0.12 (as suggested by
results in <cit.> and <cit.>, see Section <ref>)
and 0.5 were tested.

Most rates are known to be temperature dependent <cit.>.
In this work we fixed temperature at the value T=26 C, corresponding
to a year-round stable tropical environment. The coefficients for
the different rates (in units of (days)^-1) are displayed in
Table <ref>, along with the references where the model
structure and rates have been discussed:






  §.§.§ Hybridisation degree


The inheritance of genotype is acted at fecundation. On one hand,
the different subpopulations are classified in four compartments,
depending on they being “purely wild” or mixed (e.g., a wild female
mating a released male), with 0, 1 or 2 copies of the LG (here
labelled W, M_0, M_1, M_2). When the treatment starts,
the first generation offspring of the released males with wild females
will invariably land on compartment M_1 but later generations
may mix further, since pre-adult mortality of the mixed offspring
is not total. In any case, the offspring falls always in one of the
four compartments.

Regarding the rest of the genetic material, the average degree of
mixing R of each compartment is computed in several steps. This
value changes for each individual at fecundation/oviposition. Further,
when new individuals enter a compartment, the compartmental average
for R changes accordingly and this modification propagates along
the life cycle. Wild population has initial mixing value 0, while
released mosquitoes have value 1, corresponding to 100% Rockefeller-strain
genetic material. Following the law of independent assortment of alleles,
offspring has a mixing equal to the average of the parental values,
R_o=1/2(R_m+R_f). Laid eggs (eq. <ref>),
however, distribute among compartments according to mendelian laws.
Finally, the average mixing of the affected compartments is consequently
moved as, 

    R_n:=R_nX(n)+R_oE/X(n)+E,

where E is the amount of eggs (with inherited mixing R_0)
belonging to compartment n and X_nis the preexistent number
of eggs in that compartment (with previous average mixing R_n).



§ RESULTS


The code used for the simulations is deposited on GitHub. See <https://github.com/MNatiello/gene-diffusion-mosquito>

Model simulations were performed using the following scheme: 


  * The initial wild populations and environmental food conditions were
set to the deterministic equilibrium value (see Appendix II for a
description).

  * For a transient time t_T=156 weeks (about 3 years), the wild
population was allowed to proceed without interference. The population
data was saved to start the next simulation.

  * After the previous transient, intervention was started, for a duration
t_I=52 or 260 weeks for the results of Subsection <ref>
and t_I=116 weeks for those in Subsection <ref>.

  * After finishing the intervention, the system was simulated for another
156 weeks.

  * The exhibited results are the average of 10 runs of the previous
scheme.

  * The short treatment was targeted to a weekly fluorescence ratio of
0.5 and the longer ones to ratios of 0.5 and 0.65 between
eggs hatching to fluorescent vs total number of hatched eggs (per
week).

  * The recycle coefficient xl at equilibrium was set to 0.12 for
the short treatment while for the long treatment we illustrate the
dynamics with xl=0, 0.12 and 0.5.

  * There is no allowance for differentiated hybrid adult fitness (other
than in Subsection <ref>).




 §.§ Summary of results


We illustrate the simulations results with a few graphs indicating
the time-evolution of the displayed properties. The abscissa of all
graphs has units of days, counted from the beginning of the simulation.
The transient, treatment and final portions are evident.

In Figure <ref> we display the time-evolution of the
number of adult females carrying no LG for the short treatment (1
year), where the ratio of eggs with LG to total number of eggs each
week is targeted to 50%. The recycling coefficient of dead larvae
and pupae is set to 0.12. The drop in population size during the
intervention is evident. Initially, the drop is more intense since
the adjusting mechanism to target fluorescent eggs is slow (see below
for a discussion). For comparison, we plot one single simulation run
superimposed to the average of 10 simulation runs.


Figure <ref> displays the time-evolution of the number
of adult females carrying no LG for the long treatment (5 years),
where the ratio of eggs with LG to total number of eggs each week
is targeted to 50% or 65%. Same considerations as in the previous
picture hold. The three curves correspond to different values of the
recycling coefficient of dead larvae and pupae (see caption).



The target goal we attempt to simulate follows <cit.> and
aims to obtain a 50 or 65% proportion of laid eggs every week carrying one or two
copies of the LG. The idea of monitoring eggs in order to regulate
the treatment may be simple, but the control procedure is quite involved.
On one hand, eggs have to be hatched in order to observe the larval
fluorescence, meaning that the information about the system in actual
implementations is obtained with certain delay (of the order of days).
On the other hand, the control is acted by modifying the release,
i.e., mainly adult males. In the model, it takes two (target 50%)
or three (target 65%) weeks before the release of adults propagates
to a modified proportion of eggs of about half the desired target.

The simulations display a reduction of adult females (and all other
subpopulations except adult males) during the treatment period. This
reduction could be said to be environment dependent, since it varies
with the ability of the oviposition/breeding sites to recycle dead
larvae and pupae into more food. There is an initial drop passing
the target since the control procedure is indirect (eggs are monitored
but the only modification to the system acts through the introduction
of adult males). The other target in <cit.> is satisfied
in excess, as seen in Figure <ref>, left. The adult male
ratio is permanently above 10:1 during the treatment, while in
the most demanding treatments it can be as high as 40:1. 

After finishing the treatment, the system returns to a situation that
is comparable with the initial equilibrium condition, the main difference
being displayed in Figure <ref>, right. This figure displays
the time-evolution of the proportion of alien genes in the population,
a measure for genetic diffusion. Only the recycling coefficient 0.12
is plotted, since the differences with the other situations are minor.
Three of the four female compartments die out, remaining only the
hybrid compartment with no LG. The new equilibrium condition, although
being quantitatively similar to the one before treatment, bears now
different mosquitoes, the difference being larger for the longer and
more demanding treatments.


 The fact that the system returns to its previous equilibrium value
can be understood in different levels. From the simulation viewpoint,
the system was assumed to be in a stable equilibrium before treatment.
Despite the more or less drastic influence of the treatment, it eventually
returns to this equilibrium when left alone. However, simulations
with large weekly releases (about eight times larger than the maintenance
release for the 0.5 target) under a sufficiently long time (more
than one year) shift the final equilibrium to a smaller permanent
population with 100% Rockefeller hybridisation homozygous in the
LG. In Appendix D we comment on the effect of constant weekly release
of different sizes.

From the modelling viewpoint, we pointed out before that there exists
no knowledge about the adaptability of the modified hybrid mosquito
populations to the actual environment. This possibility is consequently
not considered in the model. Similarly, it is not known whether mosquitoes
with 50% of the genetic content of the Rockefeller strain differ
from the wild strain in their behaviour and their efficiency as vectors
of viral diseases.



 §.§ Dependence of the degree of hybridisation with hybrid fitness.


The target value was adjusted to obtain a control (reduction) of the
viable-egg population close to 70%, resulting in a value of 0.65
for our targeting strategy. In Table <ref>
we show the achieved population control, the hybridisation level,
R, and its range (given as the extreme values of R over 10 runs
of the simulation) for different control targets and recycling level
of food, while using different hypotheses regarding the fitness of
the hybrid adult population. O1 are hybrids with one LG and Rock
is the hybrid with no LG. At the end of the runs all the population
was of Rock type. The intervention were chosen to be 116 weeks
long in accordance with the intervention at Juazeiro <cit.>. 






 §.§ Report of simulated scenarios: Questions and answers.


We summarise our findings with a set of questions and answers. The
first two questions are those present in the Introduction and originally
triggering this study.
00.00.0000

  A. Which lasting modifications of the environment will be produced
by these tests? 

The main modification is the hybridisation of the released mosquitoes
with the local populations as a result of genetic diffusion. The degree
of hybridisation depends both on the duration of the treatment and
its intensity (target), as well as on the (environmental) differences
in food dynamics.

  B. Will the population of vectors previous to the environmental
intervention be re-established?

The recovery time of populations can be counted in months. There is
no difference in population sizes before and after the treatment for
the targeted programmes.

  C. What is the expected ratio between released males and local
males necessary to achieve different levels of control? 

The ratio ranges from 12:1 to 40:1 to sustain a target of 0.5
to 0.65 of the eggs having the modified genetics. Ratios larger
than 100:1 are needed to detect sizeable effects at the beginning
of the intervention. This can be compared with the reports in <cit.>
where the test started setting a ratio of 10:1 in a comparatively
large area but it soon became clear that a smaller (or much smaller)
area and a higher ratio (25:1) were necessary in order to observe
some effects. Recall that the simulations assume a closed environment.
Field trials are likely to require higher ratios because of the possible
immigration of wild adults into the intervention region (see point
F.) and the emigration of released mosquitoes.

  D. How are the results affected by the limitations of a hypothetical
mosquito factory? 

In the simulations we have arbitrarily chosen to limit the production
to 100000 modified mosquitoes per week. Some explorations have
managed to pass this limit. The present simulations always remain
below the limit. In practice, a smaller area requires fewer mosquitoes.
Therefore, the production limit is not an issue at this level. Repeatedly
hitting the limit would slow down or reduce the possibility of attaining
the targets of the treatment.

  E. How does recycling of dead larvae and pupae affect the dynamics?


Dead individuals could be reprocessed in different ways depending
on e.g., temperature, bacterial contents, etc. This issue is not intrinsic
to the RIDL-SIT technique apart from the fact that this or any other
technique must act on a given natural environment. The model allows
for sensing this issue to some extent. The effects appear to be less
important regarding the genetic diffusion but sensibly larger regarding
the temporary reduction in subpopulation sizes during the treatments.

  F. Is it possible to eliminate the mosquito?

The model does not consider immigration (or emigration) to (from)
the intervention region. Under these conditions, the targeted treatment
will not eliminate the mosquito but rather replace it with some hybrid
strain. Simulations releasing weekly all the production of the factory
suggest that a 100% Rockefeller strain with two copies of the
LG may replace the wild population (remaining stable) if the treatment
is long enough (see Appendix D). Analytical results confirm this observation.
Population levels are expected to be somewhat smaller and hence mosquitoes
larger than in the untreated situation.

  G. How does the technique compare to other infrastructural measures
when it comes to costs, permanent benefits, risks and mosquito populations?

The model does not contemplate social and economic costs (in part
because of lack of information). The simulations suggest that the
technique will become a permanent patch, this is, there will
be no end to the release programme and as such it represents a structural
cost that should be compared with other long-term policies (such as
running water, a sewage system, proper handling of water reservoirs,
etc.) that set the possibility of a healthier environment.




§ COMPARISON WITH FIELD OBSERVATIONS


Since the initial submission of the present work, two relevant articles
reporting on the evolution of Aedes aegypti populations during
and after interventions have been published. It is relevant then to
compare quantitatively the field observations with our calculations.

The first work, <cit.>, follows the release of genetically
modified mosquitoes in Jacobina and Juazeiro, two locations in Brazil.
The reported reduction in egg population is of 70%, and we will
compare their results with our runs for a 65% target (achieving
68% reduction). Releases lasted for 10 months (Juazeiro) and
27 months (Jacobina). The authors report that reduction of the
population lasted for 17 weeks and 32 weeks after the end of
the releases. Actually, the data for Jacobina presents the first important
peak 14 weeks after the end of the releases. As the authors discuss,
these differences are expected because the times at which the intervention
ceased corresponds to different weather conditions for each case,
namely a favourable season (rainy, Juazeiro) and an unfavourable one
(dry, Jacobina). The present model for the 65% reduction target
estimates 13-18 weeks for a full recovery of the population
under optimal weather conditions, the shorter times achieved when
O1 adults are considered unfitted (the results present little sensibility
to the fitness of Rock hybrid type). In our model, a decrease of
90% in the number of females carrying one or two LG's is observed
after 70 days (~2.5 months) (for O1 unfitted
adults but fitted Rock hybrids) while field data reports that “transgenic
larvae were detected up to” 2 months (Jacobina) and 5 months
(Juazeiro) after the release interruption. The recovery of the population
during the rainy (favourable) season observed in the field is somewhat
faster than in the model. The eradication of the population with LG
is produced by the competition with the wild type, a faster recovery
of the latter implies a faster eradication of larvae carrying the
LG, thus our overestimation of the time after the intervention when
they are detectable is linked to the same causes as the recovery time.
When O1 adults are considered fitted, the O1 eggs are no longer
present after 154 days, a number harder to consider in agreement
with field data.

The second work, <cit.>, measured hybridisation in Jacobina.
We will compare their results with a simulation for a 27 months
intervention (116 weeks) at 65% reduction level (see Table
<ref>). They report “The degree of introgression
is not trivial. Depending on sample and criterion used to define unambiguous
introgression, from about 10% to 60% of all individuals have
some OX513A genome”. Their procedure is to compute an h-value,
namely the fraction of modified genes out of a set of 21770 well
validated, biallelic, mendelian genes presenting polymorphism. h-values
larger than 0.02 or 0.04 are taken to be “unambiguous introgressed
individuals”. The h-value may be compared with the average degree
of hybridisation R in our simulations.

The R value corresponds to an upper limit since the migration flux
favours a lower level of hybridisation measurable in the field and
the simulated intervention considers that the weather was favourable
to mosquito development all the time, but this is not true for the
case studies. Of the different cases of fitness of the hybrid population
studied, we can consider that the introgression simulated with O1
mosquitoes being as fit as wild mosquitoes produces an R value
which is too high compared to field results. When mosquitoes carrying
the LG are considered unfitted, it is required that the Rock mosquitoes
(hybrids with no LG) have a better fitness than the O1 in order
to match the observations. From our simulations and the reported field
studies it is not possible to conclude that the Rock hybrids are
poor competitors in front of the Wild type.

In summary, the results are consistent with hybrids carrying the LG
being substantially unfit but first generations of Rock hybrids
not being substantially different, in terms of fitness, to the Wild
type. At this point it is worth to consider Figure 2 in <cit.>
where an increase in population is observed during the weeks 18
(2013) and 16 (2014) during the intervention time, and an increase
in population after week 16 (2014), when the intervention had already
ended. This is, in the field data, the pre-intervention population
is lower than the post-intervention population. There are a number
of hypotheses that can account for such changes, and they are worth
to  be investigated, such as: Changing weather conditions, improved
genetics, changes in standard control measures applied and changes
in people's behaviour (relaxing environmental sanity).

Considering the limited information available about the local and
genetically modified populations, that the trials were performed in
a heterogeneous and open area in contrast with the model, and also
that weather data was not introduced in our simulations (and is not
available), the performance of our minimal-model must be considered
fairly good since it produces values in the same range as the corresponding
measurements for all the variables considered. 



§ FINAL DISCUSSION AND CONCLUSIONS


We have introduced a minimal model that allows study of hybridisation
produced by interventions with a particular form of genetically modified
mosquitoes. Several dynamically relevant factors had to be omitted
because of the lack of empirical information. In this respect, the
minimal model can serve as the basis of a more detailed model when,
and if, the information becomes available.

Introgression of genes (as an independent action) in a population
subject to natural selection results always in a population equally
or better adapted than the pre-existing population. This is a basic
theorem (supported on the mathematics behind genetic propagation)
in as much as fitness can be associated with the genetics of the population.
Only the quality and quantity of the improvement remains to be ascertained
in each situation. The same theorem indicates that the result of lowering
genetic variability will in general decrease fitness, the limiting
case being a null decrease. Hence, if controlling the population with
a SIT technique results in a decrease of variability, the total outcome
of the RIDL-SIT experiment is completely uncertain.

When results from the model are compared to recent interventions we
observe that they are quantitatively compatible when hybrids carrying
the LG are considered to have a mating efficiency equal to mosquitoes
carrying two copies of the LG and an increased mortality rate, while
hybrids carrying no LG are closer in mating and survival to the wild
population. This situation implies that if F1 males (all of them
pertaining to the O1 class of the model) are as unfitted as released
males and F1 females (and all members of the O1 class) die twice
as faster that wild females (an exaggeration), genetic diffusion still
reaches important levels. The fitness relative to the Wild population
cannot be asserted by the present methods. The intervention is not
self-cleaning and will leave behind genetic contamination with the
unavoidable consequences stated by the theorem just mentioned.

The present model is a “proof of concept” rather than a “production
model”. One of the goals of this work has been to show how to model
the problem of genetic diffusion and to highlight its relevance. Mosquito
models can be highly elaborated and as such they are sensible to local
differences in populations or environment. A compartmental model intending
to forecast the outcome of an intervention can be conceived by merging
the modelling strategies developed in this work into a detailed model.
The information needed to asses particular interventions includes
the reaction norms for developmental treats as a function of genetics
(wild (local) mosquitoes, F1 hybrids, F1b hybrids), concentration
of food and temperature. The mating efficiency of F1 and F1b
males is needed as well (along with more standard information such
as several years of weather records: at least rain and average daily
temperature).

The processes of development of models and of introduction of new
technologies present some parallels that are worth discussing. In
the first place, both processes respond to the different attitudes
towards knowledge and learning, or more generally to the epistemic
frame(s) in action.

From pre-modelling (intuition-based) conclusions we enter the realm
of theoretical conclusions based on mathematical models. A first stage
in model development can be called the chimeric level where
the modelled mosquito only roughly resembles the real one (take this
expression as self-criticism <cit.>), but at least
the time-scales and life-cycle characteristics are present. Models
mature into ecologically oriented models <cit.>
fundamentally including phenotypic plasticity <cit.> and
focusing on one strain at a time (subpopulations of mosquitoes in
different cities may be strikingly different <cit.>
and even within the same city <cit.>, different genetics
is associated to different breeding sites) and may further progress
towards evolutionary models <cit.> accounting for the evolution
of phenotypic plasticity. In this higher level we still have everything
to learn.

In parallel, technologies start as laboratory tests, essentially ignorant
of the hazards and singularities of life outside the laboratory. At
this level, technologies are short-sighted; they cannot look ahead
the time and space scales of the controlled laboratory tests. In order
to mature into reliable tools, they too have to develop, including
insight about the mutual influences between technologies and environment.
This triggers the need of a research programme and an action policy
which eventually may outcome a strategy to handle the challenges of
the environment and of (co)evolution.

The perspectives of the RIDL-SIT technique in Aedes aegypti
considered before this work are clearly at the first stage. We have
pointed out in the previous sections issues about missing information
on both the life-cycle of the modified mosquitoes and on their interaction
with the environment and with the wild strain. Moreover, already at
the level of bibliographical background research we note that important
requirements mentioned e.g., in <cit.> are not fulfilled.

The information so far produced corresponds to the development of
the technology and extrapolations of laboratory results based upon
intuitions. We have termed this epistemic frame the technological
frame <cit.>. This issue has been discussed as:
“typological versus population thinking in biology“ as well <cit.>.
Typological thinking usually implies the fundamental decision of not
exploring the differences between a wild environment and the laboratory.

In this work we implemented an environment-aware model of the life-cycle
of a strain of Aedes aegypti compatible with field observations
<cit.> also introducing some genetics since
the issue at stake with the RIDL-SIT technique is to take advantage
of genetic differences between mosquito strains. We have demonstrated
that, with the help of mathematical models, it is possible to detect
some of the information needed to make less risky decisions than just
“trial and error” (which has been the method of technology in
all civilisations but it is not the method of science), and that there
is more to monitor than the success of the method considered only
in the restricted terms of the developers of the technology. The model
suggests that not even in a minimal area protected from immigration
of wild specimens the technique is able to eliminate the mosquito
population. Even worse, the output of the treatment is to produce
a new, hybrid, mosquito strain, whose characteristics are unknown.
The RIDL-SIT strategy can be viewed as a technological patch for an
unsustainable control programme, a conclusion that has been reached
previously from social considerations <cit.>. Indeed, the
results of this work indicate that there is a long path that needs
to be travelled before moving into field tests.

The limited nature of the model manifests the need of understanding
evolutionary issues in order to proceed further: How might the mosquito
change under selective pressure? How will the differences between
hybrid and wild mosquitoes manifest themselves? The results of hybridisation
cannot be forecast in the present context of lack of information,
yet this difficulty should not be rephrased as “there are no consequences
expected from hybridisation”. In all cases, the larger genetic pool
available will result in a better matching of the local forms of the
mosquito to the environment. This situation opens a large range of
possibilities we only can speculate about. They go from a return to
the pre-release equivalent epidemic risk to situations of increased
risk. Let us consider two possible mechanisms for the latter. If the
larger fertility provided by the Rockefeller background of the released
mosquitoes results in larger fertility of the new hybrid population
established in place of the original wild population, without altering
the available biomass of food, a larger number of adults of smaller
size than the original wild population is to be expected. A larger
number of vectors facilitates the propagation of mosquito transmitted
diseases, while the favourable or unfavourable incidence of the size
of the vector is a matter of discussion <cit.>. In the same
form, a co-evolution of bacteria in which Ae ae. feeds prompted
by the increase in mortality in the breeding sites is expected to
result in substantially larger populations because of the sensitivity
detected to “recycling” dead pupae and larvae. Thus, an “arms
race” scenario might arise (increased population retaliated by increased
releases). The epidemiological situation after the intervention could
then be worse than before the intervention. Some detailed models like
SkeeterBuster predict overshooting of mosquito populations going above
the equilibrium value after interventions, others like aedesBA do
not <cit.>, it is a matter under discussion that depends
on modelling decisions under uncertainty. We must conclude that the
possibility of a renewed risk for epidemics after the intervention
has ceased should be taken into account. The observed increases in
population during and after the Juazeiro tests should be addressed,
they cannot be ignored.

This work, as much as previous experience in modelling Aedes
aegypti populations, indicates that local differences in environment
and mosquito's phenotype are reflected in the population dynamics
in sensible forms. This observation warns us about the risks of attributing
to the species what is in fact the result of the genetics and environment.
Mosquito strains are different, they interact in (so far) unexpected
ways and are subject to a large, and mostly unknown today, number
of influences. The technological frame (or typological thinking) excels
at the time of proposing interventions but ecological or population
thinking must be used in evaluating the proposals. In the same direction,
it would be a mistake to extrapolate the situation studied in this
work to other conceivable cases of RIDL-SIT; the evaluation of every
case must be performed based on actual measurements of biological
traits, and such experimental values should be produced. Possible
progress is then conceived as a dialogue between both biological views.
The present work indicates that mathematical modelling can be an important
form to explore proposals and possible outcomes.

The standard analogy between controlling fruit flies and controlling
mosquitoes such as Ae ae. is not a very close one from an ecological
point of view. A very obvious difference is the landscape, going from
agricultural to urban, yet a deeper difference consists in that Ae
ae. is an established domestic insect limited only by the environment,
while the fruit fly, when production has not been abandoned, is not
limited by the available environment, but rather the fly is just attempting
to colonise it, being in short number relative to the carrying capacity.
Current knowledge in SIT indicates that “releasing sterile insects
routinely in certain areas may be more expedient to prevent establishment
of major pests than eliminating them after they become established”
<cit.>.



§ ACKNOWLEDGEMENTS


HGS acknowledges support from the University of Buenos Aires under
grants 20020130100778BA and 20020130100361BA, as well as support from
the MOSTICAW (STIC-AMSUD) collaboration. MAN acknowledges support
from SveFUM and from Kungliga Fysiografiska Sällskapet i Lund.




80


urlstyle

[Alphey et al.(2010)Alphey, Benedict, Bellini, Clark, Dame, Service,
  and Dobson]alph10
Luke Alphey, Mark Benedict, Romeo Bellini, Gary G Clark, David A Dame, Mike W
  Service, and Stephen L Dobson.
Sterile-insect methods for control of mosquito-borne diseases: an
  analysis.
Vector-Borne and Zoonotic Diseases, 100 (3):0
  295–311, 2010.

[Morlan et al.(1962)Morlan, McCray Jr, and Kilpatrick]morl62
HB Morlan, EM McCray Jr, and JW Kilpatrick.
Field tests with sexually sterile males for control of Aedes
  aegypti.
Mosquito News, 220 (3):0 295–300, 1962.

[Gong et al.(2005)Gong, Epton, Fu, Scaife, Hiscox, Condon, Condon,
  Morrison, Kelly, Dafa'alla, et al.]gong05
Peng Gong, Matthew J Epton, Guoliang Fu, Sarah Scaife, Alexandra Hiscox,
  Kirsty C Condon, George C Condon, Neil I Morrison, David W Kelly, Tarig
  Dafa'alla, et al.
A dominant lethal genetic system for autocidal control of the
  mediterranean fruitfly.
Nature biotechnology, 230 (4):0 453, 2005.

[Benedict and Robinson(2003)]bene03
Mark Q Benedict and Alan S Robinson.
The first releases of transgenic mosquitoes: an argument for the
  sterile insect technique.
Trends in parasitology, 190 (8):0 349–355,
  2003.

[Phuc et al.(2007)Phuc, Andreasen, Burton, Vass, Epton, Pape, Fu,
  Condon, Scaife, Donnelly, et al.]phuc07
Hoang Kim Phuc, Morten H Andreasen, Rosemary S Burton, Céline Vass,
  Matthew J Epton, Gavin Pape, Guoliang Fu, Kirsty C Condon, Sarah Scaife,
  Christl A Donnelly, et al.
Late-acting dominant lethal genetic systems and mosquito control.
BMC Biology, 50 (1):0 11, 2007.

[Harris et al.(2011)Harris, Nimmo, McKemey, Kelly, Scaife, Donnelly,
  Beech, Petrie, and Alphey]harr11
Angela F Harris, Derric Nimmo, Andrew R McKemey, Nick Kelly, Sarah Scaife,
  Christl A Donnelly, Camilla Beech, William D Petrie, and Luke Alphey.
Field performance of engineered male mosquitoes.
Nature biotechnology, 290 (11):0 1034–1037,
  2011.

[Harris et al.(2012)Harris, McKemey, Nimmo, Curtis, Black, Morgan,
  Oviedo, Lacroix, Naish, Morrison, et al.]harr12
Angela F Harris, Andrew R McKemey, Derric Nimmo, Zoe Curtis, Isaac Black,
  Siân A Morgan, Marco Neira Oviedo, Renaud Lacroix, Neil Naish, Neil I
  Morrison, et al.
Successful suppression of a field mosquito population by sustained
  release of engineered male mosquitoes.
Nature biotechnology, 300 (9):0 828–830,
  2012.

[Lacroix et al.(2012)Lacroix, McKemey, Raduan, Wee, Ming, Ney, AA,
  Salman, Subramaniam, Nordin, et al.]lacr12
Renaud Lacroix, Andrew R McKemey, Norzahira Raduan, Lim Kwee Wee, Wong Hong
  Ming, Teoh Guat Ney, Siti Rahidah AA, Sawaluddin Salman, Selvi Subramaniam,
  Oreenaiza Nordin, et al.
Open field release of genetically engineered sterile male Aedes
  aegypti in Malaysia.
PloS one, 70 (8):0 e42771, 2012.

[Carvalho et al.(2015)Carvalho, McKemey, Garziera, Lacroix, Donnelly,
  Alphey, Malavasi, and Capurro]carv15
Danilo O Carvalho, Andrew R McKemey, Luiza Garziera, Renaud Lacroix, Christl A
  Donnelly, Luke Alphey, Aldo Malavasi, and Margareth L Capurro.
Suppression of a field population of Aedes aegypti in Brazil
  by sustained release of transgenic male mosquitoes.
PLoS Negl Trop Dis, 90 (7):0 e0003864, 2015.

[Winskill et al.(2015)Winskill, Carvalho, Capurro, Alphey, Donnelly,
  and McKemey]wins15
Peter Winskill, Danilo O Carvalho, Margareth L Capurro, Luke Alphey, Christl A
  Donnelly, and Andrew R McKemey.
Dispersal of engineered male Aedes aegypti mosquitoes.
PLoS neglected tropical diseases, 90 (11):0
  e0004156, 2015.

[World Health Organization(2016)]who16b
World Health Organization.
Promising new tools to fight Aedes mosquitoes.
Bulletin of the World Health Organ 2016, 94:0
  562–563, 2016.
http://dx.doi.org/10.2471/BLT.16.020816.

[Center for Veterinary Medicine(2016)]fda16
Center for Veterinary Medicine.
Environmental assessment for investigational use of Aedes
  aegypti ox513a.
Technical report, United States Food and Drug Administration,
  Department of Health and Human Services, 2016.

[The Entomological Society of America (ESA) and Sociedade
  Entomológica do Brasil (SEB)(2016)]alag16
The Entomological Society of America (ESA) and Sociedade Entomológica do
  Brasil (SEB).
Aedes aegypti and the zika virus: Symptoms of a wider problem.
In Summit on the Aedes aegypti Crisis in the Americas,
  2016.

[Reis-Castro and Hendrickx(2013)]reis13
Luisa Reis-Castro and Kim Hendrickx.
Winged promises: exploring the discourse on transgenic mosquitoes in
  Brazil.
Technology in Society, 350 (2):0 118–128,
  2013.

[Enkerlin et al.(2015)Enkerlin, Gutiérrez-Ruelas, Cortes, Roldan,
  Midgarden, Lira, López, Hendrichs, Liedo, and Arriaga]enke15
Walther Enkerlin, José Manuel Gutiérrez-Ruelas, Antonio Villaseñor
  Cortes, Edgar Cotoc Roldan, David Midgarden, Estuardo Lira, Jose Luis Zavala
  López, Jorge Hendrichs, Pablo Liedo, and Francisco Javier Trujillo
  Arriaga.
Area freedom in Mexico from Mediterranean fruit fly (diptera:
  Tephritidae): a review of over 30 years of a successful containment program
  using an integrated area-wide sit approach.
Florida entomologist, 980 (2):0 665–681,
  2015.

[Hendrichs et al.(2005)Hendrichs, Vreysen, Enkerlin, and Cayol]hend05
J Hendrichs, MJB Vreysen, WR Enkerlin, and JP Cayol.
Strategic options in using sterile insects for area-wide integrated
  pest management.
In Sterile Insect Technique, pages 563–600. Springer, 2005.

[Klassen and Curtis(2005)]klas05
W Klassen and CF Curtis.
History of the sterile insect technique.
In Sterile Insect Technique, pages 3–36. Springer, 2005.

[Weidhaas et al.(1974)Weidhaas, Breeland, Lofgren, Dame, and
  Kaiser]weid74
DE Weidhaas, SG Breeland, CS Lofgren, DA Dame, and R Kaiser.
Release of chemosterilized males for the control of anopheles
  albimanus in El Salvador.
The American journal of tropical medicine and hygiene,
  230 (2):0 298–308, 1974.

[Rogers and Randolph(1984)]roge84
David J Rogers and Sarah E Randolph.
From a case study to a theoretical basis for tsetse control.
International Journal of Tropical Insect Science, 50
  (5):0 419–423, 1984.

[Dobson et al.(2002)Dobson, Fox, and Jiggins]dobs02
Stephen L Dobson, Charles W Fox, and Francis M Jiggins.
The effect of wolbachia-induced cytoplasmic incompatibility on host
  population size in natural and manipulated systems.
Proceedings of the Royal Society of London B: Biological
  Sciences, 2690 (1490):0 437–445, 2002.

[Tejerina et al.(2009)Tejerina, Almeida, and Almirón]Tejerina2009
Edmundo Fabricio Tejerina, Francisco Felipe Ludueña Almeida, and
  Walter Ricardo Almirón.
Bionomics of Aedes aegypti subpopulations (Diptera:
  Culicidae) from Misiones province, northeastern Argentina.
Acta Tropica, 1090 (1):0 45–49, Jan 2009.
10.1016/j.actatropica.2008.09.014.
URL <http://dx.doi.org/10.1016/j.actatropica.2008.09.014>.

[Grech et al.(2010)Grech, Ludueña-Almeida, and Almirón]Grech2010
Marta Gladys Grech, Francisco Ludueña-Almeida, and Walter Ricardo Almirón.
Bionomics of Aedes aegypti subpopulations (Diptera:
  Culicidae) from Argentina.
Journal of Vector Ecology, 350 (2):0 277–285,
  Dec 2010.
10.1111/j.1948-7134.2010.00083.x.
URL <http://dx.doi.org/10.1111/j.1948-7134.2010.00083.x>.

[Otero et al.(2008)Otero, Schweigmann, and Solari]oter08
M Otero, N Schweigmann, and H G Solari.
A stochastic spatial dynamical model for Aedes aegypti.
Bulletin of Mathematical Biology, 70:0 1297–1325,
  2008.

[Bergero et al.(2013)Bergero, Ruggerio, Lombardo, Schweigmann, and
  Solari]berg13
P Bergero, C Ruggerio, R Lombardo, N Schweigmann, and H Solari.
Dispersal of Aedes aegypti: field study in temperate areas
  and statistical approach.
Journal of Vector Borne Diseases, 50:0 163–170, 2013.
Journal of Vector Borne Diseases.

[Sharpe and DeMichele(1977)]shar77
P J H Sharpe and D W DeMichele.
Reaction kinetics of poikilotherm development.
Journal of Theoretical Biology, 64:0 649–670, 1977.

[Rueda et al.(1990)Rueda, Patel, Axtell, and Stinner]rued90
L. M. Rueda, K. J. Patel, R. C. Axtell, and R. E. Stinner.
Temperature-dependent development and survival rates of Culex
  quinquefasciatus and Aedes aegypti (Diptera: Culicidae).
Journal of Medical Entomology, 27:0 892–898, 1990.

[Gloria-Soria et al.(2016)Gloria-Soria, Ayala, Bheecarry,
  Calderon-Arguedas, Chadee, Chiappero, Coetzee, Elahee, Fernandez-Salas,
  Kamal, Kamgang, Khater, Kramer, Kramer, Lopez-Solis, Lutomiah, Martins Jr,
  Micieli, Paupy, Ponlawat, Rahola, Rasheed, Richardson, Saleh, Sanchez-Casas,
  Seixas, Sousa, Tabachnick, Adriana Troyo, and Powell]glor16
Andrea Gloria-Soria, Diego Ayala, Ambicadutt Bheecarry, Olger
  Calderon-Arguedas, Dave D Chadee, Marina Chiappero, Maureen Coetzee,
  Khouaildi Bin Elahee, Ildefonso Fernandez-Salas, Hany A Kamal, Basile
  Kamgang, Emad I. M. Khater, Laura D. Kramer, Vicki Kramer, Alma Lopez-Solis,
  Joel Lutomiah, Martins Jr, Maria Victoria Micieli, Christophe Paupy,
  Alongkot Ponlawat, Nil Rahola, Syed Basit Rasheed, Joshua B. Richardson,
  Amag A. Saleh, Rosa Maria Sanchez-Casas, Gonçalo Seixas, Carla A. Sousa,
  Walter J. Tabachnick, Adriana Adriana Troyo, and Jeffrey R. Powell.
Global genetic diversity of Aedes aegypti.
Molecular ecology, 250 (21):0 5377–5395,
  2016.

[Chadee and Martinez(2016)]chad16
Dave D Chadee and Raymond Martinez.
Aedes aegypti (l.) in Latin American and Caribbean
  region: With growing evidence for vector adaptation to climate change?
Acta tropica, 156:0 137–143, 2016.

[Hamilton(2011)]hamilton2011population
Matthew Hamilton.
Population genetics.
John Wiley & Sons, 2011.

[Gillett(1955a)]gill55
J D Gillett.
Variation in the hatching-response of Aedes eggs (Diptera:
  Culicidae).
Bulletin of the World Health Organization, 46:0
  241–255, 1955a.

[Gillett(1955b)]gill55b
J D Gillett.
The inherited basis of variation in the hatching response of
  Aedes eggs (Diptera: Culicidae).
Bulletin of the World Health Organization, 46:0
  255–265, 1955b.

[Gillett(1959)]gill59
J D Gillett.
Control of hatching in prediapause eggs of Aedes mosquitoes.
Nature, 104:0 1621–1623, 1959.

[Gillett et al.(1977)Gillett, Roman, and Phillips]gillet77
J D Gillett, E A Roman, and V Phillips.
Erratic hatching in aedes eggs: A new interpretation.
Proc. R. Soc. Lond. B, 196:0 223–232, 1977.

[Livdahl and Edgerly(1987)]livd87
Todd P Livdahl and J S Edgerly.
Egg hatching inhibition: Field evidence for population regulation in
  a treehole mosquito.
Ecological Entomology, 12:0 395–399, 1987.

[Edgerly and Marvier(1992)]edge92
Janice S Edgerly and Michelle A Marvier.
To hatch or not to hatch? egg hatch response to larval density and to
  larval contact in a treehole mosquito.
Ecological Entomology, 17:0 28–32, 1992.

[Edgerly et al.(1993)Edgerly, Willey, and Livdhal]edge93
J S Edgerly, M S Willey, and T P Livdhal.
The community ecology of Aedes egg hatching: Implications
  for mosquito invasion.
Ecological Entomology, 18:0 123–128, 1993.

[Southwood et al.(1972)Southwood, Murdie, Yasuno, Tonn, and
  Reader]sout72
T. R. E. Southwood, G. Murdie, M. Yasuno, R. J. Tonn, and P. M. Reader.
Studies on the life budget of Aedes aegypti in Wat
  Samphaya Bangkok Thailand.
Bulletin of the World Health Organization, 46:0
  211–226, 1972.

[Focks et al.(1993a)Focks, Haile, Daniels, and
  Moun]fock93a
D A Focks, D C Haile, E Daniels, and G A Moun.
Dynamics life table model for Aedes aegypti: Analysis of the
  literature and model development.
Journal of Medical Entomology, 30:0 1003–1018,
  1993a.

[Maciá(2006)]Macia2006
Arnaldo Maciá.
Differences in performance of Aedes aegypti larvae raised at
  different densities in tires and ovitraps under field conditions in
  Argentina.
Journal of Vector Ecology, 310 (2):0 371–377,
  Dec 2006.

[Maciá(2009)]maci09
A Maciá.
Effects of larval crowding on development time, survival and weight
  at metamorphosis in Aedes aegypti (Diptera: Culicidae).
Revista Sociedad Entomológica Argentina, 68 (1-2):0
  107–114, 2009.

[Romeo Aznar et al.(2015)Romeo Aznar, De Majo, Fischer,
  Natiello, and Solari]rome15
Victoria Romeo Aznar, María Sol De Majo, Sylvia Fischer, Mario A
  Natiello, and Hernán G Solari.
A model for the development of  Aedes (Stegomyia) aegypti
  (and other insects) as a function of the available food.
Journal of Theoretical Biology, 365:0 311–324, 2015.
URL <http://dx.doi.org/10.1016/j.jtbi.2014.10.016>.

[Arrivillaga and Barrera(2004)]arri04
J Arrivillaga and R Barrera.
Food as a limiting factor for Aedes aegypti in water-storage
  containers.
Journal of Vector Ecology, 29:0 11–20, 2004.

[Padmanabha et al.(2011)Padmanabha, Bolker, Lord, Rubio, and
  Lounibos]Padmanabha2011
H. Padmanabha, B. Bolker, C. C. Lord, C. Rubio, and L. P. Lounibos.
Food availability alters the effects of larval temperature on aedes
  aegypti growth.
Journal of Medical Entomology, 480 (5):0
  974–984, Sep 2011.

[Esteva and Yang(2005)]este05
Lourdes Esteva and Hyun Mo Yang.
Mathematical model to assess the control of aedes aegypti mosquitoes
  by the sterile insect technique.
Mathematical biosciences, 1980 (2):0 132–147,
  2005.

[Yakob et al.(2008)Yakob, Alphey, and Bonsall]yakob2008aedes
Laith Yakob, Luke Alphey, and Michael B Bonsall.
Aedes aegypti control: the concomitant role of competition, space and
  transgenic technologies.
Journal of applied Ecology, 450 (4):0
  1258–1265, 2008.

[Li(2009)]li2009simple
Jia Li.
Simple stage-structured models for wild and transgenic mosquito
  populations.
Journal of Difference Equations and Applications, 150
  (4):0 327–347, 2009.

[White et al.(2010)White, Rohani, and Sait]white2010modelling
Steven M White, Pejman Rohani, and Steven M Sait.
Modelling pulsed releases for sterile insect techniques: fitness
  costs of sterile and transgenic males and the effects on mosquito dynamics.
Journal of Applied Ecology, 470 (6):0
  1329–1339, 2010.

[Lee et al.(2013a)Lee, Baker, Gaffney, and
  White]lee2013modelling
S Seirin Lee, RE Baker, EA Gaffney, and SM White.
Modelling aedes aegypti mosquito control via transgenic and sterile
  insect techniques: Endemics and emerging outbreaks.
Journal of theoretical biology, 331:0 78–90,
  2013a.

[Lee et al.(2013b)Lee, Baker, Gaffney, and
  White]lee2013optimal
S Seirin Lee, Ruth E Baker, Eamonn A Gaffney, and Steven M White.
Optimal barrier zones for stopping the invasion of aedes aegypti
  mosquitoes via transgenic or sterile insect techniques.
Theoretical ecology, 60 (4):0 427–442,
  2013b.

[Magori et al.(2009)Magori, Legros, Puente, Focks, Scott, Lloyd, and
  Gould]mago09
Krisztian Magori, Mathieu Legros, Molly E. Puente, Dana A. Focks, Thomas W.
  Scott, Alun L. Lloyd, and Fred Gould.
Skeeter buster: A stochastic, spatially explicit modeling tool for
  studying Aedes aegypti population replacement and population
  suppression strategies.
PLoS Neglected Tropical Diseases, 30 (9):0
  e508, 09 2009.
10.1371/journal.pntd.0000508.
URL <http://dx.doi.org/10.1371%2Fjournal.pntd.0000508>.

[Legros et al.(2012)Legros, Xu, Okamoto, Scott, Morrison, Lloyd, and
  Gould]legros2012assessing
Mathieu Legros, Chonggang Xu, Kenichi Okamoto, Thomas W Scott, Amy C Morrison,
  Alun L Lloyd, and Fred Gould.
Assessing the feasibility of controlling aedes aegypti with
  transgenic methods: a model-based evaluation.
PloS one, 70 (12):0 e52235, 2012.

[Robert et al.(2014)Robert, Legros, Facchinelli, Valerio, Ramsey,
  Scott, Gould, and Lloyd]robert2014mathematical
Michael A Robert, Mathieu Legros, Luca Facchinelli, Laura Valerio, Janine M
  Ramsey, Thomas W Scott, Fred Gould, and Alun L Lloyd.
Mathematical models as aids for design and development of
  experiments: the case of transgenic mosquitoes.
Journal of medical entomology, 490 (6):0
  1177–1188, 2014.

[Focks et al.(1993b)Focks, Haile, Daniels, and
  Mount]fock93b
D A Focks, D C Haile, E Daniels, and G A Mount.
Dynamic life table model for Aedes aegypti: Simulations
  results.
Journal of Medical Entomology, 30:0 1019–1029,
  1993b.

[Ellis et al.(2011)Ellis, Garcia, Focks, Morrison, and
  Scott]Ellis2011
Alicia M Ellis, Andres J Garcia, Dana A Focks, Amy C Morrison, and Thomas W
  Scott.
Parameterization and sensitivity analysis of a complex simulation
  model for mosquito population dynamics, dengue transmission, and their
  control.
American Journal of Tropical Medicine and Hygiene, 850
  (2):0 257–264, Aug 2011.
10.4269/ajtmh.2011.10-0516.
URL <http://dx.doi.org/10.4269/ajtmh.2011.10-0516>.

[Legros et al.(2016)Legros, Otero, Romeo Aznar, Solari, Gould, and
  Lloyd]legr16
Mathieu Legros, Marcelo Otero, Victoria Romeo Aznar, Hernan Solari, Fred Gould,
  and Alun L Lloyd.
Comparison of two detailed models of aedes aegypti population
  dynamics.
Ecosphere, 70 (10), 2016.

[Romeo Aznar et al.(2013)Romeo Aznar, Otero, de Majo, Fischer, and
  Solari]rome13
V Romeo Aznar, M J Otero, M S de Majo, S Fischer, and H G Solari.
Modelling the complex hatching and development of Aedes
  aegypti in temperated climates.
Ecological Modelling, 253:0 44–55, 2013.
URL <http://dx.doi.org/10.1016/j.ecolmodel.2012.12.004>.

[Durrett(2001)]durr01
R Durrett.
Essentials of Stochastic Processes.
Springer Verlag, New York, 2001.

[Gilpin and McClelland(1979)]gilp79
Michael E Gilpin and GAH McClelland.
Systems analysis of the yellow fever mosquito Aedes aegypti.
Fortschritte der Zoologie, 25:0 355:388, 1979.

[Romeo Aznar(2015)]rome15b
Victoria Romeo Aznar.
Biología teórica, modelo y experimentos aplicados al
  entendimiento de la dinámica poblacional del mosquito Aedes
  aegypti.
PhD thesis, Facultad de Ciencias Exactas y Naturales, Universidad de
  Buenos Aires, 2015.
URL
  <http://digital.bl.fcen.uba.ar/gsdl-282/cgi-bin/library.cgi?a=d   c=tesis   d=Tesis_5982_RomeoAznar>.

[Kolmogoroff(1931)]Kolmogoroff1931
A Kolmogoroff.
Über die analytischen methoden in der wahrscheinlichkeitsrechnung.
Mathematische Annalen, 104:0 415–458, 1931.
ISSN 0025-5831.
URL <http://dx.doi.org/10.1007/BF01457949>.
10.1007/BF01457949.

[Feller(1940)]Feller1940
Willy Feller.
On the integro-differential equations of purely discontinuous
  Markoff processes.
Transactions of the American Mathematical Society, 480
  (3):0 488–515, 1940.
ISSN 00029947.
URL <http://www.jstor.org/stable/1990095>.

[Ethier and Kurtz(1986)]ethi86
S N Ethier and T G Kurtz.
Markov Processes.
John Wiley and Sons, New York, 1986.

[Kurtz(1970)]kurt70
T G Kurtz.
Solutions of ordinary differential equations as limits of pure jump
  Markov processes.
J. Appl. Prob., 7:0 49–58, 1970.

[Kurtz(1978)]kurt78
T. G. Kurtz.
Strong approximation theorems for density dependent Markov chains.
Stochastic Processes and Applications, 6:0 223, 1978.

[Solari and Natiello(2014)]sola14
H G Solari and M A Natiello.
Linear processes in stochastic population dynamics: Theory and
  application to insect development.
The Scientific World Journal - Journal of Probability and
  Statistics., ID 873624:0 1–15, 2014.
URL <http://downloads.hindawi.com/journals/tswj/aip/873624.pdf>.

[Kendall(1949)]kend49
David G Kendall.
Stochastic processes and population growth.
Journal of the Royal Statistical Society. Series B
  (Methodological), 11:0 230–282, 1949.

[Kendall(1950)]kend50
David G. Kendall.
An artificial realization of a simple “birth-and-death” process.
Journal of the Royal Statistical Society. Series B
  (Methodological), 12:0 116–119, 1950.

[Otero et al.(2006)Otero, Solari, and Schweigmann]oter06
M Otero, H G Solari, and N Schweigmann.
A stochastic population dynamic model for Aedes aegypti:
  Formulation and application to a city with temperate climate.
Bulletin of Mathematical Biology, 68:0 1945–1974,
  2006.

[Trpis(1972)]trpi72
M. Trpis.
Dry season survival of Aedes aegypti eggs in various breeding
  sites in the Dar es Salaam area, Tanzania.
Bull. W.H.O., 47:0 433–437, 1972.

[Garziera et al.(2017)Garziera, Pedrosa, de Souza, Gómez, Moreira,
  Virginio, Capurro, and Carvalho]garz17
Luiza Garziera, Michelle Cristine Pedrosa, Fabricio Almeida de Souza, Maylen
  Gómez, Márcia Bento Moreira, Jair Fernandes Virginio, Margareth Lara
  Capurro, and Danilo Oliveira Carvalho.
Effect of interruption of over-flooding releases of transgenic
  mosquitoes over wild population of a edes aegypti: two case studies in
  brazil.
Entomologia Experimentalis et Applicata, 1640
  (3):0 327–339, 2017.

[Evans et al.(2019)Evans, Kotsakiozi, Costa-da Silva, Ioshino,
  Garziera, Pedrosa, Malavasi, Virginio, Capurro, and Powell]evan19
Benjamin R Evans, Panayiota Kotsakiozi, Andre Luis Costa-da Silva,
  Rafaella Sayuri Ioshino, Luiza Garziera, Michele C Pedrosa, Aldo Malavasi,
  Jair F Virginio, Margareth L Capurro, and Jeffrey R Powell.
Transgenic aedes aegypti mosquitoes transfer genes into a natural
  population.
Scientific reports, 90 (1):0 1–6, 2019.

[Reed et al.(2010)Reed, Waples, Schindler, Hard, and Kinnison]reed10
Thomas E Reed, Robin S Waples, Daniel E Schindler, Jeffrey J Hard, and
  Michael T Kinnison.
Phenotypic plasticity and population viability: the importance of
  environmental predictability.
Proc Biol Sci, 2770 (1699):0 3391–3400, Nov
  2010.
10.1098/rspb.2010.0771.
URL <http://dx.doi.org/10.1098/rspb.2010.0771>.

[Paupy et al.(2004)Paupy, Chantha, Huber, Lecoz, Reynes, Rodhain, and
  Failloux]paup04
Christophe Paupy, Ngan Chantha, Karine Huber, Nicolas Lecoz, Jean-marc Reynes,
  Francois Rodhain, and Anna-Bella Failloux.
Influence of breeding sites features on genetic differentiation of
  Aedes aegypti populations analyzed on a local scale in Phnom Penh
  municipality of Cambodia.
The American journal of tropical medicine and hygiene,
  710 (1):0 73–81, 2004.

[Scheiner(1993)]sche93
Samuel M Scheiner.
Genetics and evolution of phenotypic plasticity.
Annual review of ecology and systematics, 240
  (1):0 35–68, 1993.

[Solari(2016)]sola16b
Hernán G Solari.
Recognizing the epistemic frame in opportunity of the zika pandemic.
In Congrès mondial pour la pensée complexe (Les défis d'un
  monde globalisé), 2016.
URL
  <https://www.reseau-canope.fr/fileadmin/user_upload/Projets/pensee_complexe/solari_recognizing_epistemic_frame_opportunity_zika_pandemic.pdf>.

[Solari(2017)]sola17
Hernán G Solari.
Epistemic considerations on research about flavivirus induced fevers.
2017.
URL
  <http://users.df.uba.ar/solari/papers/sola17-ddenuefever.pdf>.

[Powell(2018)]powe18
Jeffrey Powell.
Genetic variation in insect vectors: Death of typology?
Insects, 90 (4):0 139, 2018.

[Juliano et al.(2014)Juliano, Ribeiro, Maciel-de Freitas, Castro,
  Codeço, Lourenço-de Oliveira, and Lounibos]juli14
Steven A Juliano, Gabriel Sylvestre Ribeiro, Rafael Maciel-de Freitas,
  Márcia G Castro, Claudia Codeço, Ricardo Lourenço-de
  Oliveira, and L Philip Lounibos.
She’sa femme fatale: low-density larval development produces good
  disease vectors.
Memorias do Instituto Oswaldo Cruz, 1090 (8):0
  1070–1077, 2014.

[Verhulst(1838)]verh38
PF Verhulst.
Notice sur la loi que la population suit dans son accroissement.
Correspondence Mathematic et Physic, 10:0 113–121,
  1838.

[Solari and Natiello(2003)]sola03
H G Solari and M A Natiello.
Stochastic population dynamics: the Poisson approximation.
Physical Review E, 67:0 031918, 2003.








§ LETHAL-GENE ENHANCED MORTALITY AT PUPAE LEVEL 


The assumptions behind the modified death and emergence rates for
pupae are as follows. For the model, pupae can undergo only two events,
namely emergence m_p→ a and death m_p. We assume that
the time to the next event is not altered with respect to the wild
population, so in both cases this time responds to an exponential
distribution with parameter R_pupa=m_p→ a+m_p=m_pl→ a+m_pl
(where l stands for “lethal”). Further, we assume that 5%
of pupae carrying the LG emerge while the other 95% dies. Hence,

    m_pl→ a/m_pl→ a+m_pl=0.05

 and therefore

    m_pl→ a+m_pl/m_pl→ a=20  or   m_pl=19m_pl→ a

Hence,

    m_pl→ a    =    1/20(m_p→ a+m_p),  and
    
    m_pl    =    19/20(m_p→ a+m_p).





§ DETERMINISTIC EQUILIBRIUM POPULATION


Deterministic population equations as those are frequently seen in
epidemiological modelling, have an older history <cit.> than
the stochastic approach. However, deterministic models have a more
limited value, since they average away relevant issues of the dynamics.
The most conspicuous limitation is that these models cannot properly
describe few individuals, as would be the case in problems where there
is a potential extinction risk (such as ours). In certain situations
<cit.> the deterministic dynamics can be seen as a
sort of large-population limit of the stochastic approach. However,
since the wild population is not spontaneously facing extinction,
the deterministic equilibrium values are a reasonable starting point
for our simulations.

To produce deterministic equations out of the stochastic setup amounts
to combine the rates W_α and associated modifications δ_j^α
populationwise.

The general expression for the deterministic dynamics of each population
species in the absence of released individuals reads, 

    dX_j/dt=∑_αW_α(X_1,⋯,X_n)δ_j^α=∑_αm_α^iX_iδ_j^α

where the sum over α is relevant only for those events modifying
population X_j and X_i is the specific subpopulation associated
to the event with coefficient m_α^i. At this point recall
that fecundation is not treated as an event, i.e., it is assumed to
occur in a time-scale that is negligible in front of the life-span
of adult females. Letting E, L, P, F denote the size of the
wild populations of eggs, larvae, pupae and females, we obtain, 


    dE/dt    =     m_ovi c_layFT(P_f/C_f) F-(m_e→ l+m_e)E
    dL/dt    =     m_e→ lE-LE(P_f/C_f) L
    dP/dt    =     LE(P_f/C_f) (1-ML(P_f/C_f)) L-0.5787 P
    dF/dt    =    1/20.5787 (1-ML(P_f/C_f))P-0.04/1-ML(P_f/C_f)F

The food equilibrium condition influences all equations and also the
larvae population L. Being the deterministic system quasilinear
<cit.> in the populations, the equilibrium is obtained when
the determinant of the coefficient matrix is zero, determining in
such a form the food level of the equilibrium. In terms of exponential
races, the equilibrium is independent of the average time required
for the race and depends only on the probability of occurrence of
one or another outcome. The population number is determined by the
production of food by the environment. The equilibrium values were
obtained solving the equilibrium food dynamics to obtain the equilibrium
P_f/C_f value and the equilibrium larvae population
L. Subsequently, all other populations were obtained by solving
dE/dt=0, dP/dt=0, dF/dt=0 as a function
of L and P_f/C_f . 



§ REACTION NORMS DETERMINED FOR THIS WORK


The data available from the experiment reported in <cit.>
was used to produce the estimates for the dependence of body size
(surrogated by wing length) and the rate for events (pupation or death)
as larvae. The results are displayed in Figure <ref>.
The data displays clearly two different responses to food depletion:
the starvation region where the probability of mortality increases
but the average time for the next event in the larvae compartment
does not change significantly; and the scarcity region where mortality
is not increased but the average time between events (the reciprocal
of the rate) changes monotonically. The left panel of Figure <ref>
shows the wing length as a function of the logarithm of the fraction
of optimal food density. Wing length, body size and fertility are
monotonically related. The right panel shows the event rate for larvae
as a function of the same environmental variable.






§ CONSTANT VS TARGETED RELEASE


The goal of control strategies is to reduce mosquito populations to
a level in which epidemic outbreaks are not possible. Releasing genetically
modified mosquitoes in excess of those needed to reach the target
represents a misuse of costly resources. A compromise must be taken
between the area under control and the degree of control (meaning
the number of mosquitoes per human being for example).

In contrast, in a trial phase, one may want to test how to achieve
maximal effect. To depict this situation we simulated five “large”
constant weekly releases (larger or much larger than the maintenance
release targeted in Figures <ref>, <ref>
and <ref>), acting under two years, for the case of 12%
food recycling. Only for the larger values the wild and hybrid (with
no LG) populations are outnumbered to extinction in the present modelling
conditions (closed environment with no immigration/emigration of individuals,
no heterogeneities), while a smaller population of 100% Rockefeller
with 2 LG individuals persists for over three years after the intervention.
Results are displayed in Figure <ref>. However, this
situation is unrealistic in field conditions where immigration/emigration
cannot be controlled.




