Modeling confinement and reversibility of threshold-dependent gene drive systems in spatially-explicit Aedes aegypti populations
BMC Biology volume 18, Article number: 50 (2020)
The discovery of CRISPR-based gene editing and its application to homing-based gene drive systems has been greeted with excitement, for its potential to control mosquito-borne diseases on a wide scale, and concern, for the invasiveness and potential irreversibility of a release. Gene drive systems that display threshold-dependent behavior could potentially be used during the trial phase of this technology, or when localized control is otherwise desired, as simple models predict them to spread into partially isolated populations in a confineable manner, and to be reversible through releases of wild-type organisms. Here, we model hypothetical releases of two recently engineered threshold-dependent gene drive systems—reciprocal chromosomal translocations and a form of toxin-antidote-based underdominance known as UDMEL—to explore their ability to be confined and remediated.
We simulate releases of Aedes aegypti, the mosquito vector of dengue, Zika, and other arboviruses, in Yorkeys Knob, a suburb of Cairns, Australia, where previous biological control interventions have been undertaken on this species. We monitor spread to the neighboring suburb of Trinity Park to assess confinement. Results suggest that translocations could be introduced on a suburban scale, and remediated through releases of non-disease-transmitting male mosquitoes with release sizes on the scale of what has been previously implemented. UDMEL requires fewer releases to introduce, but more releases to remediate, including of females capable of disease transmission. Both systems are expected to be confineable to the release site; however, spillover of translocations into neighboring populations is less likely.
Our analysis supports the use of translocations as a threshold-dependent drive system capable of spreading disease-refractory genes into Ae. aegypti populations in a confineable and reversible manner. It also highlights increased release requirements when incorporating life history and population structure into models. As the technology nears implementation, further ecological work will be essential to enhance model predictions in preparation for field trials.
The discovery of CRISPR and its application as a gene editing tool has enabled gene drive systems to be engineered with much greater ease [1, 2]. Recent attention has focused on homing-based drive systems and their potential to control mosquito-borne diseases on a wide scale, either by spreading disease-refractory genes  or by spreading genes that confer a fitness load or sex bias and thereby suppress mosquito populations [4, 5]. The increased ease of gene editing has also advanced the entire field of gene drive, including systems appropriate during the trial phase of the technology . Such systems would ideally be capable of enacting local population control by (a) effectively spreading into populations to achieve the desired epidemiological effect, and (b) being recallable from the environment in the event of unwanted consequences, public disfavor, or the end of a trial period.
As gene drive technology has progressed, a number of systems have been proposed with the potential to enact localized population control without spreading on a wide scale [6, 7]. Sterile male releases provide one option , a recent version of which is based on the same molecular components as CRISPR gene drive systems [8, 9]. At the interface between homing-based and localized suppression systems, an autosomal X chromosome-shredding system has been proposed that induces a transient male sex bias and hence population suppression before being selected out of the population . Population modification drive systems that display transient drive activity before being eliminated by virtue of a fitness cost could also spread disease-refractory genes into populations in a localized way. Examples of this variety of drive system include split-gene drive , daisy drive , and killer-rescue systems . Each system has its own strengths and weaknesses and could be suited to a different situation. In this paper, we theoretically explore the potential for two recently engineered threshold-dependent gene drive systems to achieve localized and reversible population modification in structured populations—reciprocal chromosomal translocations  and a toxin-antidote-based system known as UDMEL .
Threshold-dependent gene drive systems must exceed a critical threshold frequency in a population in order to spread. Based on this dynamic, simple population models, in which two randomly mating populations exchange small numbers of migrants with each other, predict that these systems can be released at high frequencies in one population and spread to near-fixation there, but never take off in the neighboring population because they do not exceed the required threshold there [16, 17]. These systems can also be eliminated through dilution to sub-threshold levels with wild-type organisms at the release site, making them excellent candidates for the trial phase of a population modification gene drive strategy, or when localized population modification is desired . Elimination of a non-driving transgene can in fact be more difficult, as the dynamics of threshold-dependent systems actively drive them out of populations at sub-threshold levels. However, whether these dynamics hold in real ecosystems depends crucially on the dispersal patterns and population structure of the species being considered. First steps towards modeling the spatial dynamics of these systems have been taken by Champer et al. , who model spatially structured releases of various threshold-dependent systems without life history, and Huang et al. , who model engineered underdominance  on a grid-based landscape incorporating life history for Aedes aegypti, the mosquito vector of dengue, Zika, and other arboviruses.
Here, we present a detailed ecological analysis of the expected population dynamics of two recently engineered threshold-dependent drive systems, translocations and UDMEL, in Ae. aegypt in a well-characterized landscape—Yorkeys Knob, a suburb ~ 17 km northwest of Cairns, Australia (Fig. 1c)—suitable for confineable and reversible releases. Yorkeys Knob and the nearby town of Gordonvale were field sites for releases of Wolbachia-infected mosquitoes in 2011 , and the prevalence of Wolbachia infection over time provided information on the number of adult Ae. aegypti mosquitoes per household and other mosquito demographic parameters for that location , as well as an opportunity to validate our modeling framework. Wolbachia is an intracellular bacterium and biocontrol agent that biases inheritance in its favor if infected females are present, and blocks transmission of dengue and other arboviruses . Yorkeys Knob is a partially isolated suburb, separated by a 1–2-km-wide uninhabited, vegetated area from the nearest suburb, Trinity Park, a control site for the Wolbachia trial. This allowed us to simulate trials of transgenic mosquitoes in a well-characterized population, while also theoretically exploring their potential of spread to a neighboring community.
A number of other ecological details are relevant to the spread of threshold-dependent gene drive systems that have not been considered in previous modeling studies. Perhaps of greatest importance, the frequency of the introduced transgene in the mating pool is markedly different from the frequency of introduced adults. It is typical to release only male mosquitoes as part of an intervention, as only females are involved in human disease transmission. Life cycle and mating structure therefore become relevant, as immature life stages are not available for mating, and female adults are thought to mostly mate only once soon after emergence . This means that many of the released adult males will not find a mating partner, and hence, larger releases will be required to exceed threshold frequencies than predicted in simple population frequency models.
The nature of mosquito dispersal behavior is also relevant to the spatial dispersal of transgenes. Our species of interest, Ae. aegypti, is understood to display leptokurtic dispersal behavior in a suburban setting, in which mosquitoes tend to remain in the same household for the majority of their lifespan, while a few mosquitoes disperse over larger distances . With these landscape, dispersal, and life cycle considerations accounted for, we theoretically explore the ability to drive two threshold-dependent systems, translocations and UDMEL, into populations of Ae. aegypti in one community, Yorkeys Knob, without them spreading in significant numbers to a neighboring community, Trinity Park, and to be remediated from Yorkeys Knob at the end of the simulated trial period.
We use the Mosquito Gene Drive Explorer (MGDrivE) modeling framework  to model the spread of translocations and UDMEL through spatially structured mosquito populations (Fig. 1). This is a genetic and spatial extension of the lumped age-class model of mosquito ecology  modified and applied by Deredec et al.  and Marshall et al.  to the spread of homing gene drive systems. The framework incorporates the egg, larval, pupal, and adult life stages, with egg genotypes being determined by maternal and paternal genotypes and the allelic inheritance pattern of the gene drive system. Spatial dynamics are accommodated through a metapopulation structure in which lumped age-class models run in parallel and migrants are exchanged between populations according to a zero-inflated exponential dispersal kernel with parameters defined in Additional file 1: Table S1 [22, 30,31,32,33,34,35,36,37,38]. Further details of the framework are described in the “Methods” section.
Applying the MGDrivE modeling framework to our research questions, we incorporate the inheritance patterns of reciprocal chromosomal translocations and UDMEL into the inheritance module of the model (Fig. 1a, b, Additional file 2: Fig. S1), the life cycle parameters of Aedes aegypti (Additional file 1: Table S1) into the life history module, and the distribution of households in Yorkeys Knob (923 households) and Trinity Park (1301 households) along with their expected mosquito population sizes and movement rates between them into the landscape module (Fig. 1c). The suburb of Trinity Park served as a control site for field releases of Wolbachia-infected mosquitoes, to quantify the extent to which the Wolbachia infection could spread from one community to another, and plays a similar role for our simulated releases of threshold-dependent gene drive systems.
The inheritance patterns that result from chromosomal translocations are depicted in Fig. 1a. Chromosomal translocations result from the mutual exchange between terminal segments of two non-homologous chromosomes. When translocation heterozygotes mate, several crosses result in unbalanced genotypes and hence unviable offspring, resulting in a heterozygote reproductive disadvantage. This results in bistable, threshold-dependent population dynamics, confirmed in laboratory drive experiments . The inheritance patterns produced by the UDMEL system are depicted in Fig. 1b. This system consists of two unlinked constructs, each possessing a maternally expressed toxin active during oogenesis and a zygotically active antidote expressed by the opposite construct. Offspring are more likely to have the antidote(s) to the maternal toxin(s) when transgenes are present at moderate-to-high population frequency. This produces threshold-dependent dynamics since, above a critical frequency, the selective advantage of the antidotes outweighs the selective disadvantage of the toxins, and below the critical frequency, the selective disadvantage of the toxins dominates. Mathematical models predict a threshold frequency of ~ 24% in the absence of an additional fitness cost—a result consistent with laboratory drive experiments for a construct having a modest additional fitness cost .
Using data from field trials of Wolbachia-infected Ae. aegypti mosquitoes in Yorkeys Knob and Gordonvale, Australia , we validated our modeling framework prior to application to other threshold-dependent drive systems. Wolbachia biases the offspring ratio in favor of those carrying Wolbachia through a mechanism known as cytoplasmic incompatibility, in which offspring of matings between infected males and uninfected females result in the death of some or all progeny, while matings involving infected females tend to produce infected offspring . For the wMel strain of Wolbachia that was used in the Australian field trials, incompatible crosses produce no viable offspring, and Wolbachia is inherited by all offspring of infected females. There is also a fitness cost associated with Wolbachia infection, the value of which has been estimated between 0 and 20% [21, 40, 41].
We used Wolbachia surveillance data from Fig. 1 of Hoffmann et al.  to validate our model framework and parameter values. Based on the description of the field trials in Hoffmann et al. , we simulated weekly releases of 20 Wolbachia-infected mosquitoes (10 female and 10 male) per household at a coverage of 30% over 10 weeks in the spatially explicit landscapes of Yorkeys Knob and Gordonvale with the exception that in Gordonvale, the fifth release was postponed by a week due to a tropical cyclone. Model predictions were calculated for a variety of literature-based values of adult mortality rate and Wolbachia-associated fitness cost [21, 33, 40,41,42] and were compared to observed Wolbachia prevalence over time. We found that model predictions most closely matched field data for a baseline adult mortality rate of 0.090 per day  and that predictions matched field data quite well for both 5% and 10% fitness costs, with a 10% fitness cost being closer to that estimated elsewhere [21, 40, 41] (Additional file 3: Fig. S2). Model predictions in Fig. 2 use these parameter values, along with others listed in Additional file 1: Table S1, and their agreement with the observed field data provides good validation of our modeling framework.
Population replacement and remediation for translocations
The use of translocations for transforming pest populations was initially suggested by Serebrovskii  and later Curtis  for the introduction of disease-refractory genes into mosquito populations. A number of models have been proposed to describe their spread through randomly mating populations [14, 16, 45, 46]; however, with one recent exception addressing spatial structure , these have largely ignored insect life history and mating structure. Such models suggest that the translocation need only exceed a population frequency of 50%, in the absence of a fitness cost associated with the translocation, to spread to fixation in a population, which could conceivably be achieved through a single seeding release round. Here, we find that incorporating life history and population structure into mosquito population dynamic models significantly increases release requirements.
In Figs. 3 and 4, based on the precedent set by the 2011 Wolbcahia field trial , we consider weekly releases of 20 adult Ae. aegypti males homozygous for the translocation per household for given durations and coverage levels, where coverage level is the proportion of households that receive the releases. Releases are simulated in the community of Yorkeys Knob, in which prior releases of Wolbachia-infected mosquitoes suggested a local population of ~ 15 adult Ae. aegypti per household , and for mosquito movement rates inferred from previous studies [25, 37, 38] (Additional file 1: Table S1). For a coverage level of 100%, and in the absence of a fitness cost, four weekly releases of 20 Ae. aegypti males (~ 3:1 released to local males) are required for the translocation to spread to fixation in the community (Fig. 3), as opposed to the single release expected when ignoring life history and population structure . As coverage is reduced to 50%, the required number of releases increases to 7, and for a coverage level of 25%, as seen for the World Mosquito Program in Yorkeys Knob, the required number of releases increases to 16 (Figs. 3 and 4). Although large, these releases are achievable, considering the much larger releases conducted for sterile insect programs .
To simulate remediation of a translocation, we consider weekly releases of 20 adult Ae. aegypti wild-type males in the community of Yorkeys Knob, whereby the translocation has already reached fixation in that community. In the absence of a fitness cost associated with the translocation, translocations are symmetrical in their threshold dynamics, and so, for a coverage level of 100%, four weekly releases are required for the translocation to be completely remediated from the community, and for a coverage of 25%, 16 weekly releases are required for the translocation to be completely remediated (Figs. 3 and 4). Encouraging features of these results are that (i) remediation can be achieved through releases of non-biting, non-disease-transmitting males; (ii) release sizes are achievable; and (iii) despite the spatial household structure, both replacement and remediation are complete within the community. The time to replacement is highly dependent on the coverage level and number of releases, but is reasonably quick given sufficient releases. At a coverage of 50%, 20 weekly releases led to the translocation spreading to a frequency > 99% within half a year of the final release (or within 300 days of the first release). For equivalent wild-type releases, this is the same as the time to > 99% elimination.
Population replacement and remediation for UDMEL
UDMEL was the first synthetic gene drive system to be engineered that displays threshold-dependent dynamics . The system consists of two unlinked constructs, each possessing a maternally expressed toxin active during oogenesis and a zygotically active antidote expressed by the opposite construct (Fig. 1b). At low population frequencies, the maternal toxin confers a significant selective disadvantage, leading to elimination, while at high population frequencies, the zygotic antidote confers a selective benefit in the context of a prevalent toxin, leading to fixation. The dynamics of this system in randomly mating populations have been characterized by Akbari et al. , suggesting that the system need only exceed a population frequency of ~ 24%, in the absence of a fitness cost, to spread to fixation, while the wild-type must exceed a population frequency of ~ 76% to eliminate the construct. Both replacement and remediation should therefore be achievable with 1–2 releases of transgenic and wild-type organisms, respectively; however, as for translocations, we find that incorporating life history and population structure into our models increases release requirements in both cases.
In Figs. 3 and 4, we consider weekly releases of 20 adult Ae. aegypti males homozygous at both loci for the UDMEL system in the community of Yorkeys Knob. The lower threshold for UDMEL as compared to translocations means that replacement is much easier to achieve for UDMEL. For a coverage level of 50% or higher, and in the absence of a fitness cost, a single release of 20 Ae. aegypti males leads to the UDMEL system spreading to fixation throughout the community (Fig. 3). As coverage is reduced to 25%, the required number of releases to achieve fixation increases to two (Figs. 3 and 4). As for translocations, the time to replacement is highly dependent on the coverage level and number of releases. From Fig. 3, it is apparent that UDMEL reaches total allele fixation slowly, although the number of individuals having at least one copy of the transgene increases quickly. At a coverage of 50%, 10 weekly releases lead to wild-type individuals falling to a frequency < 2% within 2.6 years of the final release (or within 3 years of the first release).
Remediation, however, is more difficult to achieve for UDMEL compared to translocations due to the higher threshold that wild-type organisms must exceed to eliminate UDMEL. Additionally, wild-type females must be included in the releases to propagate the wild-type allele because, assuming continued functioning of UDMEL components, the maternal toxins of females having UDMEL at both loci kill all offspring that do not inherit UDMEL at both loci. To simulate remediation, we first consider weekly releases of 10 adult Ae. aegypti wild-type females and 10 adult males in the community of Yorkeys Knob. In the absence of a fitness cost associated with the UDMEL construct, and for a coverage level of 75%, nine weekly releases are required for a reduction in UDMEL allele frequency over the first year (Fig. 3); however, a closer inspection of the simulation results reveals that complete remediation of UDMEL from the community is not possible even with 20 releases, as the UDMEL allele frequency bounces back. Comparison of these results to those for a panmictic population with a population size equal to that of Yorkeys Knob reveals that complete remediation of UDMEL can occur at coverage levels as low as 25% (for 16 or more weekly releases) (Additional file 4: Fig. S3). Inspection of the spatially explicit simulation results suggests that the rebound in UDMEL allele frequency in the structured population is due to UDMEL remaining at super-threshold levels after the wild-type releases in a small number of households, and slowly recolonizing the landscape following that. Complete remediation of UDMEL is possible, however, for 15 or more releases at a coverage level of 100% (Fig. 4). These results make a strong case for translocations as preferred systems to introduce transgenes in a local and reversible way as (i) remediation of UDMEL requires releases of biting, vector-competent females and (ii) release requirements for these biting, vector-competent females are burdensomely high due to the high threshold that must be surpassed consistently throughout a spatially structured population.
Confinement of translocations and UDMEL to release site
Confinement of translocations and UDMEL to partially isolated populations has previously been modeled by Marshall and Hay  and Akbari et al. . In both cases, two randomly mating populations were modeled that exchange migrants at given rates. Population structure was otherwise ignored, as was mosquito life history. Results from these analyses suggest that translocations would spread and remain confined to populations for which the migration rate is less than ~ 5.8% per mosquito per generation , and that UDMEL would remain confined to populations for which the migration rate is less than ~ 1.6% per mosquito per generation . These migration rates are relatively low; however, this may be beneficial for the types of landscapes we are considering here, whereby the system may spread between neighboring households, but not from one suburb to another. Recently, Champer et al.  showed that translocations would remain confined to and persist in a population connected to another by a “migration corridor” under a range of parameter values.
For our landscape of interest—the suburbs of Yorkeys Knob and Trinity Park—it is very unlikely that Ae. aegypti mosquitoes will travel from one suburb to another by their own flight. Extrapolating the exponential dispersal kernel used in our simulations, fitted to data from mark-release-recapture experiments collated by Guerra et al.,  suggests these events to be negligible, before accounting for the fact that the intervening vegetated area may serve as a barrier to Ae. aegypti flight . Furthermore, rare migrant mosquitoes are unlikely to cause the threshold frequency for either drive system to be exceeded, thus making spatial spread due to such movements unlikely. In considering confineability to the release suburb, we therefore model “batch migration,” in which several mosquitoes are carried, perhaps by a vehicle, from one community to another at once. Batch migration events could be thought of as several adult mosquitoes being carried at once, or perhaps more likely, as a larval breeding site, such as a tyre, being carried from one household to another, with several adults emerging from the tyre following transport. We model batch migration events as occurring between randomly chosen households, and vary the number of daily migration events and the effective number of adults carried per event. For computational simplicity, we focus on migration events from Yorkeys Knob, in which either system has already reached fixation, to households in Trinity Park, which is initially fixed for wild-type mosquitoes.
In Fig. 4e, f, we see that both the number and size of daily batch migration events affect the chance of either system establishing itself in the neighboring suburb, Trinity Park. For translocations, ~ 16 daily migration events of batches of 5 adults are required for spread in Trinity Park. For batches of 10 adults, ~ 9 daily migration events are required, and for batches of 20 adults, ~ 5 daily migration events are required. For UDMEL, ~ 3 daily migration events of batches of 5 adults are required for spread in Trinity Park, and for batches of 10 adults, ~ 2 daily migration events are required.
These results continue to make a strong case for translocations as preferred systems to introduce transgenes in a local and reversible way as (i) many more batch migration events are required to lead to spread for translocations as opposed to UDMEL and (ii) the rate of migration events required for translocations to spread is higher than what would be expected between these communities. Specifically, Wolbachia releases in Yorkeys Knob in 2011 provide evidence for occasional batch migrations to the nearby suburb of Holloways Beach; however, the spatio-temporal pattern of Wolbachia spread, as inferred from monitored trap data, suggests only ~ 1–2 batch migration events over the course of a month, consisting of less than 5 adult females per event .
A theoretical study by Khamis et al.  on toxin-antidote-based underdominance gene drive systems, similar to UDMEL but for which the toxins are zygotic rather than maternal , found that the gene drive threshold frequency is highly sensitive to (i) the increase in adult mortality rate in organisms having the transgene, (ii) the duration of the larval life stage, and (iii) the parameters determining the character or strength of larval density dependence. In Fig. 5 and Additional file 5: Fig. S4, we explore the sensitivity of our model outcomes of replacement, remediation, and confinement for translocations and UDMEL as we vary (i) the duration of the larval life stage, (ii) the baseline adult mortality rate, (iii) the fitness cost associated with the gene drive system, and (iv) the mean adult dispersal distance. For translocations, we model a 10% fitness cost as a 10% reduction in mean adult lifespan for organisms homozygous for the translocation and a 5% reduction for organisms heterozygous for the translocation. For UDMEL, since its inheritance bias is induced through the action of maternal toxins, we model a 10% fitness cost as a 10% reduction in female fecundity for organisms homozygous for UDMEL at both loci, with 2.5% additive fitness costs contributed by each transgenic allele.
For translocations, the associated fitness cost had the greatest impact on the release scheme required for the system to be fixed or remediated from the population, given the life parameters considered (Fig. 5). A 10% fitness cost led to ~ 10 weekly releases at a coverage of 50% being required for the translocation to reach fixation (an increase of 3 releases), while a 20% fitness cost led to ~ 13 weekly releases being required (an increase of 6 releases). Small changes in the duration of the larval life stage had minor impacts on the release requirements, with an increase in larval lifespan of 2 days leading to one more weekly release being required for the translocation to reach fixation, and vice versa. A 2% change in the baseline adult mortality rate and 50% change in the mean migration distance had negligible impact on release requirements. Remediation, on the other hand, requires fewer wild-type releases when there is a fitness cost associated with the translocation. A 10% fitness cost led to 5 weekly releases at a coverage of 50% being sufficient to eliminate the translocation (a decrease of 2 releases), and a 20% fitness cost led to 4 weekly releases at a coverage of 50% being sufficient for elimination (a decrease of 3 releases). Small changes in the duration of the larval life stage had minor impacts on the wild-type release requirements for elimination, with an increase in larval lifespan of 2 days or a 2% decrease in the adult mortality rate leading to one fewer release being required.
The sensitivity of our predictions regarding confinement to the release site is of particular interest, as invasion of a neighboring community may be more likely under some parameter values than others. Fortunately, a fitness cost associated with the translocation leads to a higher threshold and hence more batch migration events required for invasion of a neighboring community. A 10% fitness cost led to ~ 2–3 additional daily migration events of 10 adults required for spread to Trinity Park, and a 20% fitness cost led to ~ 6–7 additional daily migration events required (Fig. 5). Also noteworthy, a 2% increase in the adult mortality rate led to ~ 2 fewer daily migration events required for spread to Trinity Park—i.e., ~ 7 migration events for batches of 10 adults and ~ 14 migration events for batches of 5 adults. While still above inferred batch migration rates, this highlights that there could exist parameter sets beyond those explored for which invasion is feasible.
UDMEL displays similar parameter sensitivities regarding fixation and batch migration outcomes as for translocations, with the exception that these outcomes are less sensitive to fitness costs (Additional file 5: Fig. S4), likely due to the fact that fitness is accommodated through a reduction in female fecundity rather than an increase in adult mortality. A 20% fitness cost led to ~ 1 additional weekly release being required for the system to spread to fixation, whether at a coverage of 25% or 50%. Similarly, for an invasion of Trinity Park, a 10% fitness cost required ~ 1 additional daily migration event of 5 adults, and a 20% fitness cost required ~ 2 additional daily migration events. Of note, a 2% increase in the adult mortality rate or a 2-day increase in the duration of the larval stage led to ~ 1 fewer daily migration event required for spread to Trinity Park, making this now very achievable—i.e., ~ 2–3 migration events for batches of 5 adults and ~ 1–2 migration events for batches of 10 adults.
Finally, we conducted an analysis of the sensitivity of our results to population structure, exploring the impact of (i) removing all population structure by treating Yorkeys Knob and Trinity Park as randomly mixing populations and (ii) incorporating heterogeneity in mosquito household population size. Results of the comparison to panmictic populations are depicted in Additional file 6: Fig. S5. Previously, we had seen that introducing population structure greatly increases the release requirements to eliminate UDMEL from a community (Additional file 4: Fig. S3). The trend of higher release requirements in structured populations is also seen for translocations, although to a lesser extent, with one additional release required for either replacement or remediation at a coverage of 75%, two additional releases required at a coverage of 50%, and 7–8 additional releases required at a coverage of 25%. Invasion of a neighboring population, on the other hand, requires ~ 1–2 fewer daily migration events in structured populations for batches of 10 adults having either the translocation or UDMEL.
Incorporating heterogeneity in household mosquito population size, we retain a mean of 15 adults per household, as inferred from Wolbachia field trial data in Yorkeys Knob , and distribute population sizes across households according to a zero-inflated, truncated exponential distribution with 55% of households having no mosquitoes and none having more than 45 adults. This distributional form, including zero inflation, is based on results of a large field survey conducted across a set of households in Kamphaeng Phet province, Thailand . Including this source of heterogeneity substantially increases release requirements for replacement and remediation with translocations, with 3–4 additional releases required at a coverage of 75%, 5–6 additional releases required at a coverage of 50%, and 20 releases being insufficient at a coverage of 25% (Additional file 7: Fig. S6). Release requirements are marginally increased for UDMEL, with ~ 1–2 additional releases required at coverages of 25–100%. These results are likely due to the threshold frequency being more difficult to exceed in households with large numbers of mosquitoes, and this being less of an issue for UDMEL due to its lower threshold frequency. Fortunately, population size heterogeneity makes confinement more promising for both systems, increasing the required number of daily migration events for batches of 10 adults by 3–4 for translocations and by 1–2 for UDMEL.
The idea of using threshold-dependent gene drive systems to replace local populations of disease vectors with varieties that are unable to transmit diseases has been discussed for over half a century now, since Curtis  famously proposed the use of translocations to control diseases transmitted by mosquitoes. While Curtis had been primarily concerned with introducing and spreading genes into a population, as the technology nears implementation, confining them is also becoming a significant concern. As CRISPR-based homing gene drive technology edges closer to field application, concerns are increasingly being raised regarding the invasiveness of these systems [50, 51], and systems such as split drive , daisy drive , and threshold-dependent underdominant systems  are gaining interest, at least during the trial phase of population replacement technology . In this paper, we model the introduction of two drive systems, chromosomal translocations and UDMEL, that have been engineered in the laboratory and shown to display threshold-dependent spread [14, 15]. While previous papers have described the population dynamics of these two systems in randomly mating populations ignoring life history [14,15,16, 44], with one recent paper including spatial structure , we present the first analysis of these systems in a spatially structured population including mosquito life history and reflecting a well-characterized landscape where field trials could conceivably be conducted .
Our results provide strong support for the use of translocations to implement confineable and reversible population replacement in structured Ae. aegypti populations. Regarding reversibility, translocations are preferable to UDMEL as (i) they can be remediated through releases of non-disease-transmitting male Ae. aegypti and (ii) required releases sizes are achievable (~ 10 weekly releases at a coverage level of 50%). UDMEL requires less effort to introduce into a population, but is much more difficult to remove once it has been introduced, requiring a very large number of both males and disease-transmitting females to be released. This highlights the benefit of a ~ 50% threshold for reversible population replacement: the symmetry allows both replacement and remediation to be achieved with similar effort. Extreme underdominance is another example of system with a 50% threshold [20, 52]. Regarding confineability, translocations again outperform UDMEL as ~ 16 daily migration events of batches of 5 Ae. aegypti adults are required for translocations to spread to the neighboring suburb of Trinity Park, while UDMEL can spread to Trinity Park given only ~ 3 daily migration events (or 2–3 daily migration events for alternative model parameterizations). The true batch migration rate between suburbs is expected to be smaller than either of these ; however, the rate required for translocations to spread is highly unlikely to be reached, while the rate for UDMEL is conceivable.
As with any modeling study, there are limitations inherent in our analysis. Several of the parameters we assumed to be constant here would indeed be dynamic in a real intervention scenario. At the genetic level, lab experiments suggest non-outbred individuals homozygous for the translocation had a fitness cost that largely disappeared once offspring were produced that were the product of at least one wild-type individual . Models fitted to data from UDMEL drive experiments also suggested dynamic fitness costs that depended on the frequency of transgenic organisms in the population . Another recent modeling study highlights the possibility of toxin and antidote mutational breakdown for underdominance constructs; however, this is expected to be impactful over a larger timescale than considered here (hundreds of generations) . At the ecological level, our model of Ae. aegypti life history , based on the lumped age-class model of Hancock and Godfray , assumes the existence of a constant equilibrium population size and other constant ecological parameters, such as adult death rate and larval development times. These parameters have indeed been shown to vary in space and time, and in response to local mosquito density [54, 55], which our sensitivity analyses suggest could have significant impacts on release thresholds and gene drive outcomes  (Figs. 5 & S4). At the landscape level, we have assumed a relatively homogenous distribution of mosquitoes per household, and movement rates between households that depend only on distance and household distribution. Extensive landscape heterogeneities have been shown to slow and alter the spread of Wolbachia [41, 56] and would likely impact the spread of translocations and UDMEL as well. Future work that helps to characterize the environmental drivers of mosquito population dynamics will inform iterative model development to address this.
In conclusion, our analysis supports the use of translocations as a threshold-dependent drive system capable of spreading disease-refractory genes into structured Ae. aegypti populations in a confineable and reversible manner. If such a system were engineered in Ae. aegypti, it would be an excellent candidate for the introduction of disease-refractory genes during the trial phase of population replacement technology, or whenever a localized release was otherwise desired. As the technology nears implementation, further ecological work characterizing the density dependencies, seasonality, and spatial heterogeneities of Ae. aegypti populations will be essential to enhance model predictions in preparation for field trials.
We used the MGDrivE framework  (https://marshalllab.github.io/MGDrivE/) to simulate releases of adult Ae. aegypti males homozygous for one of two threshold-dependent gene drive systems—reciprocal chromosomal translocations or UDMEL—in the community of Yorkeys Knob in Queensland, Australia. To simulate remediation, we modeled releases of wild-type adult Ae. aegypti into populations in Yorkeys Knob already fixed for the gene drive system. Houses receiving releases were randomly chosen from a uniform distribution. These were conserved within simulation runs, but varied between runs. To determine confineability, we simulated batch migration events from Yorkeys Knob (fixed for the gene drive system) to the neighboring community of Trinity Park (initially wild-type). The MGDrivE framework models the egg, larval, pupal, and adult (male and female) mosquito life stages implementing a daily time step, overlapping generations and a mating structure in which adult males mate throughout their lifetime, while adult females mate once upon emergence, retaining the genetic material of the adult male with whom they mate for the duration of their adult lifespan. Density-independent mortality rates for the juvenile life stages are assumed to be identical and are chosen for consistency with the population growth rate in the absence of density-dependent mortality. Additional density-dependent mortality occurs at the larval stage, the form of which is taken from Deredec et al. . Full details of the modeling framework are available in the S1 Text of Sánchez et al.  and in the software documentation available at https://marshalllab.github.io/MGDrivE/docs/reference/. Parameters describing Ae. aegypti life history and the gene drive systems and landscape of interest are listed in Additional file 1: Table S1.
The inheritance patterns for reciprocal chromosomal translocations (depicted in Fig. 1a) and UDMEL (depicted in Additional file 2: Fig. S1) are modeled within the inheritance module of the MGDrivE framework , and their impacts on female fecundity and adult lifespan are implemented in the life history module. The distribution of households in Yorkeys Knob, Trinity Park, and Gordonvale was taken from OpenStreetMap (https://www.openstreetmap.org/) (Fig. 1c). We implement the stochastic version of the MGDrivE framework to capture the randomness associated with events that occur in small populations, such as households, which serve as nodes in the landscape modeled here. In the stochastic implementation of the model, the number of eggs produced per day by females follows a Poisson distribution, the number of eggs having each genotype follows a multinomial distribution, all survival/death events follow a Bernoulli distribution, and female mate choice follows a multinomial distribution with probabilities given by the relative frequency of each adult male genotype in the population.
Availability of data and materials
The distribution of households in Yorkeys Knob, Trinity Park, and Gordonvale was obtained from OpenStreetMap (https://www.openstreetmap.org/). The software package, MGDrivE, is available on GitHub (https://marshalllab.github.io/MGDrivE/) . All information required to reproduce the simulations is available in the GitHub deposition, the manuscript, and its additional files.
Doudna JA, Charpentier E. Genome editing. The new frontier of genome engineering with CRISPR-Cas9. Science. 2014;346(6213):1258096.
Champer J, Buchman A, Akbari OS. Cheating evolution: engineering gene drives to manipulate the fate of wild populations. Nat Rev Genet. 2016;17:146–59.
Gantz VM, Jasinskiene N, Tatarenkova O, Fazekas A, Macias VM, Bier E, et al. Highly efficient Cas9-mediated gene drive for population modification of the malaria vector mosquito Anopheles stephensi. Proc Natl Acad Sci U S A. 2015;112:E6736–43.
Hammond A, Galizi R, Kyrou K, Simoni A, Siniscalchi C, Katsanos D, et al. A CRISPR-Cas9 gene drive system targeting female reproduction in the malaria mosquito vector Anopheles gambiae. Nat Biotechnol. 2016;34:78–83.
Kyrou K, Hammond AM, Galizi R, Kranjc N, Burt A, Beaghton AK, et al. A CRISPR-Cas9 gene drive targeting doublesex causes complete population suppression in caged Anopheles gambiae mosquitoes. Nat Biotechnol. 2018;36:1062–6.
Marshall JM, Akbari OS. Can CRISPR-based gene drive be confined in the wild? A question for molecular and population biology. ACS Chem Biol. 2018;13:424–30.
Leftwich PT, Edgington MP, Harvey-Samuel T, Carabajal Paladino LZ, Norman VC, Alphey L. Recent advances in threshold-dependent gene drives for mosquitoes. Biochem Soc Trans. 2018;46:1203–12.
Alphey L. Re-engineering the sterile insect technique. Insect Biochem Mol Biol. 2002;32:1243–7.
Kandul NP, Liu J, Sanchez CHM, Wu SL, Marshall JM, Akbari OS. Transforming insect population control with precision guided sterile males with demonstration in flies. Nat Commun. 2019;10:84.
Galizi R, Doyle LA, Menichelli M, Bernardini F, Deredec A, Burt A, et al. A synthetic sex ratio distortion system for the control of the human malaria mosquito. Nat Commun. 2014;5:3977.
Akbari OS, Bellen HJ, Bier E, Bullock SL, Burt A, Church GM, et al. Safeguarding gene drive experiments in the laboratory. Science. 2015;349:927–9.
Noble C, Min J, Olejarz J, Buchthal J, Chavez A, Smidler AL, et al. Daisy-chain gene drives for the alteration of local populations. Proc Natl Acad Sci U S A. 2019;116:8275–82.
Gould F, Huang Y, Legros M, Lloyd AL. A killer-rescue system for self-limiting gene drive of anti-pathogen constructs. Proc Biol Sci. 2008;275:2823–9.
Buchman AB, Ivy T, Marshall JM, Akbari OS, Hay BA. Engineered reciprocal chromosome translocations drive high threshold, reversible population replacement in drosophila. ACS Synth Biol. 2018;7:1359–70.
Akbari OS, Matzen KD, Marshall JM, Huang H, Ward CM, Hay BA. A synthetic gene drive system for local, reversible modification and suppression of insect populations. Curr Biol. 2013;23:671–7.
Marshall JM, Hay BA. Confinement of gene drive systems to local populations: a comparative analysis. J Theor Biol. 2012;294:153–71.
Altrock PM, Traulsen A, Reeves RG, Reed FA. Using underdominance to bi-stably transform local populations. J Theor Biol. 2010;267:62–75.
Champer J, Zhao J, Champer S, Liu J, Messer PW. Population dynamics of underdominance gene drive systems in continuous space. bioRxiv. 2019. https://doi.org/10.1101/449355.
Huang Y, Lloyd AL, Legros M, Gould F. Gene-drive into insect populations with age and spatial structure: a theoretical assessment. Evol Appl. 2011;4:415–28.
Davis S, Bax N, Grewe P. Engineered underdominance allows efficient and economical introgression of traits into pest populations. J Theor Biol. 2001;212:83–98.
Hoffmann AA, Montgomery BL, Popovici J, Iturbe-Ormaetxe I, Johnson PH, Muzzi F, et al. Successful establishment of Wolbachia in Aedes populations to suppress dengue transmission. Nature. 2011;476:454–7.
Ritchie SA, Montgomery BL, Hoffmann AA. Novel estimates of Aedes aegypti (Diptera: Culicidae) population size and adult survival based on Wolbachia releases. J Med Entomol. 2013;50:624–31.
Walker T, Johnson PH, Moreira LA, Iturbe-Ormaetxe I, Frentiu FD, McMeniman CJ, et al. The wMel Wolbachia strain blocks dengue and invades caged Aedes aegypti populations. Nature. 2011;476:450–3.
Hartberg WK. Observations on the mating behaviour of Aedes aegypti in nature. Bull World Health Organ. 1971;45:847–50.
Schmidt TL, Barton NH, Rašić G, Turley AP, Montgomery BL, Iturbe-Ormaetxe I, et al. Local introduction and heterogeneous spatial spread of dengue-suppressing Wolbachia through an urban population of Aedes aegypti. PLoS Biol. 2017;15:e2001894.
Sánchez CHM, Wu SL, Bennett JB, Marshall JM. MGDrivE: a modular simulation framework for the spread of gene drives through spatially-explicit mosquito populations. Methods Ecol Evol. 2019; https://doi.org/10.1111/2041-210X.13318.
Hancock PA, Godfray HCJ. Application of the lumped age-class technique to studying the dynamics of malaria-mosquito-human interactions. Malar J. 2007;6:98.
Deredec A, Godfray HCJ, Burt A. Requirements for effective malaria control with homing endonuclease genes. Proc Natl Acad Sci U S A. 2011;108:E874–80.
Marshall JM, Buchman A, Sánchez CHM, Akbari OS. Overcoming evolved resistance to population-suppressing homing-based gene drives. Sci Rep. 2017;7:3776.
Otero M, Solari HG, Schweigmann N. A stochastic population dynamics model for Aedes aegypti: formulation and application to a city with temperate climate. Bull Math Biol. 2006;68(8):1945–974.
Eisen L, et al. The impact of temperature on the bionomics of Aedes (Stegomyia) aegypti with special reference to the cool geographic range margins. J Med Entomol. 2014;51:496–516.
Simoy MI, Simoy MV, Canziani GA. The effect of temperature on the population dynamics of Aedes aegypti. Ecological Modelling. 2015;314:100–10.
Focks DA, Haile DG, Daniels E, Mount GA. Dynamic life table model for Aedes aegypti (Diptera: Culicidae): analysis of the literature and model development. J Med Entomol. 1993;30(6):1003–17.
Horsfall WR (1955) Mosquitoes: Their Bionomics and Relation to Disease. Hafner, London.
Fay RW. The biology and bionomics of Aedes aegypti in the laboratory. Mosq News. 1964;24:300–308.
Schmidt TL, Filipović I, Hoffmann AA, Rašić G. Fine-scale landscape genomics helps explain the slow spatial spread of Wolbachia through the Aedes aegypti population in Cairns, Australia. Heredity. 2018;120(5):386–95.
Harrington LC, Edman JD, Kittayapong P, Coleman RC, Clark GG, Sithiprasasna R, et al. Dispersal of the dengue vector Aedes aegypti within and between rural communities. Am J Trop Med Hyg. 2005;72(2):209–20.
Guerra CA, Reiner RC Jr, Perkins TA, Lindsay SW, Midega JT, Brady OJ, et al. A global assembly of adult female mosquito mark-release-recapture data to inform the control of mosquito-borne pathogens. Parasites & Vectors. 2014;7(1):276.
Turelli M, Hoffmann AA. Microbe-induced cytoplasmic incompatibility as a mechanism for introducing transgenes into arthropod populations. Insect Mol Biol. 1999;8:243–55.
Ross PA, Endersby NM, Hoffmann AA. Costs of three Wolbachia infections on the survival of Aedes aegypti larvae under starvation conditions. PLoS Negl Trop Dis. 2016;10:e0004320.
Hancock PA, Ritchie SA, Koenraadt CJM, Scott TW, Hoffmann AA, Godfray HCJ. Predicting the spatial dynamics of Wolbachia infections in Aedes aegypti arbovirus vector populations in heterogeneous landscapes. J Appl Ecol. 2019; https://doi.org/10.1111/1365-2664.13423.
Mordecai EA, Cohen JM, Evans MV, Gudapati P, Johnson LR, Lippi CA, et al. Detecting the impact of temperature on transmission of Zika, dengue, and chikungunya using mechanistic models. PLoS Negl Trop Dis. 2017;11:e0005568.
Serebrovskii AS. On the possibility of a new method for the control of insect pests. Zool Zhurnal. 1940;19:618–30.
Curtis CF. Possible use of translocations to fix desirable genes in insect pest populations. Nature. 1968;218:368–9.
Curtis CF, Robinson AS. Computer simulation of the use of double translocations for pest control. Genetics. 1971;69:97–113.
Wright S. On the probability of fixation of reciprocal translocations. Am Nat. 1941;75:513–22.
Carvalho DO, McKemey AR, Garziera L, Lacroix R, Donnelly CA, Alphey L, et al. Suppression of a field population of Aedes aegypti in Brazil by sustained release of transgenic male mosquitoes. PLoS Negl Trop Dis. 2015;9:e0003864.
Khamis D, El Mouden C, Kura K, Bonsall MB. Ecological effects on underdominance threshold drives for vector control. J Theor Biol. 2018;456:1–15.
Koenraadt CJ, Aldstadt J, Kijchalao U, Sithiprasasna R, Getis A, Jones JW, et al. Spatial and temporal patterns in pupal and adult production of the dengue vector Aedes aegypti in Kamphaeng Phet, Thailand. Am J Trop Med Hyg. 2008;79:230–8.
Marshall JM. The effect of gene drive on containment of transgenic mosquitoes. J Theor Biol. 2009;258:250–65.
Noble C, Adlam B, Church GM, Esvelt KM, Nowak MA. Current CRISPR gene drive systems are likely to be highly invasive in wild populations. Elife. 2018;7. https://doi.org/10.7554/eLife.33423.
Maselko M, Heinsch SC, Chacón JM, Harcombe WR, Smanski MJ. Engineering species-like barriers to sexual reproduction. Nat Commun. 2017;8:883.
Edgington MP, Alphey LS. Modeling the mutation and reversal of engineered underdominance gene drives. J Theor Biol. 2019;479:14–21.
Hancock PA, White VL, Callahan AG, Godfray CHJ, Hoffmann AA, Ritchie SA. Density-dependent population dynamics in Aedes aegypti slow the spread of w Mel Wolbachia. J Appl Ecol. 2016;53:785–93.
Ross PA, Endersby NM, Yeap HL, Hoffmann AA. Larval competition extends developmental time and decreases adult size of wMelPop Wolbachia-infected Aedes aegypti. Am J Trop Med Hyg. 2014;91:198–205.
Hancock PA, Godfray HCJ. Modelling the spread of Wolbachia in spatially heterogeneous environments. J R Soc Interface. 2012;9:3045–54.
The authors would like to thank Dr. Gregory Lanzaro, Dr. Yoosook Lee, Dr. Tomás León, and Ms. Partow Imani for discussions on Aedes aegypti life history and dispersal behavior.
This work was supported by a DARPA Safe Genes Program Grant (HR0011-17-2-0047), awarded to OSA and JMM, and funds from the Innovative Genomics Institute, awarded to JMM.
Ethics approval and consent to participate
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Life history, population size and movement parameters for Aedes aegypti in Cairns, Australia.
Complete inheritance pattern of UDMEL. UDMEL is composed of two unlinked constructs (here referred to as A and B), each consisting of a maternally-expressed toxin and a zygotically-expressed antidote for the toxin on the opposite construct (see Fig. 1b). The cross here represents matings between all nine possible parental genotypes (“+” represents the wild-type allele, and “A” and “B” represent alleles corresponding to the two UDMEL constructs). Offspring lacking the antidotes to the maternal toxins produced by their mother are unviable (shaded). At high population frequencies, the selective advantage on the constructs, by virtue of the antidotes, outweighs the fitness load due to the toxins, and hence results in frequency-dependent spread.
Comparison of model predictions to Wolbachia field trial data from Yorkeys Knob and Gordonvale, Australia for a variety of parameterizations. Field observations of Wolbachia population frequency are depicted in light blue, with 95% binomial confidence intervals based on the frequency and sample size reported for the 2011 field trials in Yorkeys Knob and Gordonvale. Model predictions are depicted for an analogous release scheme consisting of 20 Wolbachia-infected mosquitoes (10 female and 10 male) per household at a coverage of 30% over 10 weeks with the exception that, in Gordonvale, the fifth release was postponed by a week due to a tropical cyclone. Parameter values are listed in Additional file 1: Table S1. Model predictions more closely match observed field data for a baseline adult mortality rate of 0.090 per day (Top), as compared to one of 0.050 per day (Bottom). Observations match predictions well for both 5% and 10% fitness costs associated with Wolbachia. Agreement between observations and predictions is strong, providing good validation for the modeling framework.
Remediation results for UDMEL in spatially-explicit and panmictic populations. Time-series results are shown for a given number of weekly releases of 20 adult wild-type Ae. aegypti per household (10 female and 10 male) with the intent of remediation in the community of Yorkeys Knob (Fig. 1c), and at given coverage levels, where coverage is the proportion of households that receive the releases. (Top) Remediation in the spatially-explicit population is extremely difficult. At a coverage level of 75%, nine weekly releases are required for a reduction in UDMEL allele frequency over the first year; however, complete remediation is not possible even with 20 releases. Complete remediation is possible at a coverage level of 100% for 15 or more weekly releases. (Bottom) Remediation in the panmictic population is much less demanding. At a coverage level of 25%, it can be achieved with 16 or more releases, and at a coverage level of 50%, it can be achieved with nine or more releases.
Sensitivity of model outcomes for replacement and confinement of UDMEL. Changes are depicted in the proportion of 50 stochastic simulations that result in fixation for replacement and confinement of UDMEL. Proportions are compared to those in the second row of Fig. 4 as we vary: i) the mean dispersal distance of adult mosquitoes (+/- 50%), ii) the duration of the larval life stage (+/- 2 days), iii) the baseline adult mortality rate (+/- 2%), and iv) the fitness cost associated with being homozygous for the translocation (+10% or +20%). UDMEL displays similar parameter sensitivities regarding fixation and batch migration outcomes as for translocations (Fig. 5), with the exception that these outcomes are less sensitive to fitness costs, likely due to the fact that fitness is accommodated through a reduction in female fecundity rather than an increase in adult mortality.
Sensitivity of model outcomes for translocations and UDMEL comparing spatially-explicit and panmictic populations. Changes are depicted in the proportion of 50 stochastic simulations that result in fixation for replacement, remediation and confinement of translocations and UDMEL. Proportions are compared to those in Fig. 4 as we simulate a model where Yorkeys Knob and Trinity Park are panmictic populations of equivalent size to their spatially-explicit versions. Introducing population structure greatly increases the release requirements to remediate UDMEL (as seen in Additional file 4: Fig. S3), and substantially increases the release requirements for replacement or remediation of translocations. Invasion of a neighboring population, on the other hand, requires moderately fewer daily migration events in structured populations for both translocation and UDMEL.
Sensitivity of model outcomes for translocations and UDMEL comparing spatially-explicit populations with and without heterogeneity in household mosquito population size. Changes are depicted in the proportion of 50 stochastic simulations that result in fixation for replacement, remediation and confinement of translocations and UDMEL. Proportions are compared to those in Fig. 4 as we simulate a model where household mosquito population size is distributed according to a zero-inflated, truncated exponential distribution with a mean of 15 adults, 55% of households having no mosquitoes, and none having more than 45 adults. Introducing household population size heterogeneity substantially increases release requirements for replacement and remediation with translocations, and marginally increases release requirements for replacement with UDMEL. Fortunately, population size heterogeneity makes confinement moderately more promising for both systems.
About this article
Cite this article
Sánchez C., H.M., Bennett, J.B., Wu, S.L. et al. Modeling confinement and reversibility of threshold-dependent gene drive systems in spatially-explicit Aedes aegypti populations. BMC Biol 18, 50 (2020). https://doi.org/10.1186/s12915-020-0759-9