Population suppression across the study area by an ‘ideal’ drive allele
We first consider a drive allele that is ‘ideal’, in the sense that its fitness effects are completely recessive and it causes no parental effects on its offspring (i.e. heterozygous females are always fully fertile). In line with experimental data, we assume high but not complete homing (95% in both males and females [11]) and that half of non-homed alleles form non-functional resistance alleles ( [17]; henceforth ‘r2 alleles’). We suppose the drive allele has no effect on the sex ratio of offspring, though below we will consider transgenes that induce male-biased sex ratios in addition to their effects on female fertility. Absent from this model is the potential emergence of a resistant allele that restores dsx function (an ‘r1 allele’). We assume such alleles will not emerge in the case of dsx-targeting transgenes due to the highly conserved sequence of the dsx gene [11], though the potential impact of very rare r1 alleles will be addressed in the ‘Discussion’ section.
The non-spatial model predicts such a drive allele will rise to 96.2% from a small initial introduction with the remaining percentage made up of wildtype alleles (1.3%) and r2 alleles (2.5%) which persist in a mutation-selection balance. Such a predominance of non-functional alleles is predicted to result in a genetic load of 0.974, leading to population extinction given our assumed population growth rate when rare of 18.9 (see the ‘Methods’ section).
We simulated the addition of heterozygous adult male mosquitoes each year to sites chosen at random and calculated the population suppression as the reduction in the number of biting female mosquitoes relative to a non-intervention baseline (see the ‘Methods’ section). The number of release sites per year is important in determining how quickly widespread suppression is achieved, while not affecting the eventual suppression level (Additional file 1: Fig. S1). The number of mosquitoes liberated per release, by contrast, makes little difference to either short- or long-term suppression. Significant suppression can be achieved within 4 years if there are at least a few hundred release sites per year. If releases of 5000 males are made at 1% (434) of sites each year, the female population is reduced by ~ 95% after about 4 years (~ 95% suppression; Fig. 1a; 4-year suppression ranged from 92.6% to 95.8% with mean 94.6% from 10 simulations). In the first few years after releases begin, somewhat greater suppression is achieved if the release sites are distributed on a regular grid rather than chosen at random (this resulted in 95.5% (94.5–97.2%) suppression). This difference diminishes with time, however, and is not apparent after 8 years.
Note that females that lack a functional dsx gene are ignored when we calculate suppression because they are unable to take blood meals. In the above example with randomly located release sites, the inclusion of dsx-negative females would lower the population suppression predicted by the spatial model to ~ 92%. We discuss below the differences in suppression predicted by the spatial and non-spatial models.
Costly drive alleles
Somatic expression of Cas9
Next, consider the possibility that there is some somatic expression of the Cas9 nuclease encoded by the drive allele, yet there is no parental deposition of the protein. We assume the somatic expression of Cas9 reduces the fertility of drive heterozygous females because it partly or completely prevents the production of the dsx gene product in cells that require it. The non-spatial model predicts that such a drive allele will cause population extinction provided the reduction in heterozygous fertility is less than ~ 0.44. For costs greater than this, suppression is predicted to decline to a limit of 0 if heterozygous females are completely infertile (black line in Fig. 1b).
Somatic expression reduces the impact of the drive for two reasons. First, there is a direct cost to the transmission of the drive allele from one generation to the next, because drive heterozygous females will have fewer offspring. Second, the cost increases the relative fitness of r2 alleles over the drive allele, because r2/wildtype heterozygous females are assumed to be fully fertile, and so the frequency of r2 alleles increases at the expense of drive alleles. We find the latter cost has a greater effect on the level of predicted suppression; if no r2 alleles are produced, population extinction is predicted by the non-spatial model for costs less than ~ 0.96, a higher threshold (compare the grey and black lines in Fig. 1b). Both factors act to slow the spread of the drive allele, however, so that as rates of somatic expression increase so does the time it takes to reach maximum population suppression.
The spatial model predicts a similar response to somatic Cas9 expression, though long-term suppression is consistently a few percentage points lower than the genetic load predicted by the non-spatial model (Fig. 1b). Stable levels of suppression are attained after about 4 years when somatic expression costs are 0.4 or less, while for higher costs, the time to equilibrium requires a few more years.
Parental deposition of Cas9
Now consider the possibility that Cas9 enters the zygote through deposition in the gametes so that female fertility is reduced when either the mother or the inseminating father carries the drive allele (with the paternal and maternal effects possibly being different). We assume that there are no additional genotype-specific costs, so that drive heterozygous females have wildtype fitness. Analysis of the non-spatial model suggests that deposition costs can prevent the transgene bringing about population extinction, though only if the effects of paternal deposition are relatively strong (Fig. 1c). Thus, in the absence of maternal effects, paternal effects must cause a greater than ~ 90% reduction in fitness to stop extinction, and as maternal effect costs rise to 100%, this threshold drops to ~ 80%. The greater importance of paternal over maternal deposition is because both drive homozygous and heterozygous males are assumed to be fully fertile, so individuals are more likely to have a father than a mother with the drive allele (see also [9]). Parental effects have less impact than somatic expression, because they affect females with r2 alleles as well as females with drive alleles.
As parental deposition costs increase from zero, there is a small increase in population suppression (genetic load) due to reduced female fecundity (Additional file 1, Fig. S2). However, for higher costs, there is an abrupt transition in dynamics with much reduced genetic load and hence the absence of population extinction. For the case of equal paternal and maternal costs, this threshold occurs when costs exceed 0.813 (Additional file 1, Fig. S2). Mathematically, the system displays two equilibria, with low and high genetic loads. In the absence of, or with moderate, deposition costs, the high and low genetic load equilibria are respectively stable and unstable, but this switches at the threshold giving the abrupt transition. Biologically, below the threshold, the selective advantage that drive alleles have over wildtype alleles due to homing is sufficient to outweigh the costs of parental deposition, but above the threshold, the balance tips the other way allowing the wildtype and r2 alleles to invade. Note that high deposition costs can have a large effect in slowing down the approach to equilibrium (Additional file 1: Fig. S3).
Parental deposition has a greater effect in reducing population suppression in the spatial compared to the non-spatial model, though again paternal deposition costs are more important than maternal. Thus, to achieve at least 90% suppression after 8 years of releases (roughly equivalent to population extinction in the non-spatial context), the reduction in fitness due to paternal deposition must be less than 60% (no cost of maternal deposition) to 40% (maternal deposition renders females infertile). We hypothesise that the dynamics of the spatial model are also affected by underlying alternative stable and unstable equilibria whose precise properties are influenced by local environmental factors in a way that lowers the threshold of the transition to the low genetic load equilibrium. As in the non-spatial model, we find deposition costs slow the effect of the transgene in achieving population suppression.
Combined somatic expression and parental costs
We also consider the possibility of both somatic expression and parental Cas9 deposition acting in concert, assuming the costs of the latter to be the same for maternal and paternal depositions (Fig. 1d). The non-spatial model predicts both types of cost act to reduce eventual suppression though, as discussed above, the somatic costs have a somewhat greater effect. The spatial model, by contrast, predicts that the cost of Cas9 deposition is approximately equal to the cost of somatic expression. Again, we hypothesis that the differences are due to subtle interactions between the genetic dynamics and spatially heterogeneous environmental factors.
Paternal male bias
We next suppose a paternally expressed male-biasing component has been added to the transgene, in line with the construct developed by Simoni et al. [12]. We varied the extent of male bias from 0 to 1 (when all progeny of transgenic males are male), for three levels of fitness cost to the heterozygous females (Fig. 2). The non-spatial and spatial models both indicate that paternal male bias will be of little benefit to a strong drive allele. Indeed, high levels of male bias reduce the potential of the transgene to suppress populations because this results in the creation of fewer transgene homozygous females (see also [12]). However, the suppression potential of moderate and weak transgenes can be substantially enhanced by paternal male bias. This is because the presence of transgenic males skews the population sex ratio towards males, thus reducing both the number of adult females and the population growth rate. In the limit of paternal bias equalling one, there is no difference between drive alleles that differ in heterozygous female fitness cost, because heterozygous females are not produced in this scenario.
Spatial variation
To understand why the spatial models consistently predict lower population suppression than the non-spatial models, we next consider the geographic variation in suppression for three constructs that again differ in their fitness costs experienced by heterozygous females (Fig. 3; assuming no paternal male bias). Though these alleles differ considerably in average suppression after 8 years (from 38.5% for the weak allele to 95.8% for the strong), they share a remarkably similar pattern of the strongest suppression in the least seasonal regions. However, the apparently consistent role of seasonality masks a number of qualitative differences in how strong and weak drive alleles influence local population dynamics. To illustrate these differences, we study four exemplar sites, shown as points a–d in Fig. 3. Site ‘a’, in Mali, is highly seasonal and also somewhat isolated from other settlements; site ‘b’, in Niger, is highly seasonal and highly connected; site ‘c’, on the Benin-Burkina Faso border, is somewhat isolated yet due to its proximity to the Pendjali river is assumed to have year-round larval habitats; and finally, site ‘d’, which is in Ghana, is both highly connected and assumed to have year-round larval habitat due to a mild dry season in this region. Figure 4 shows the typical simulated dynamics at each site and for each of the three drive alleles described above.
If the drive allele is strongly suppressing (left column of Fig. 4), colonisation-extinction dynamics occur in much of the study area (and are seen in all example sites except site ‘d’). The colonisation and extinction cycles are fastest in sites with strongly seasonal environments (sites ‘a’ and ‘b’) and most regular in sites within locally dense networks of human populations (site ‘b’). In these densely populated regions, frequent dispersal of mosquitoes between sites ensures their dynamics are correlated and similar levels of suppression are observed in all locations. In more sparsely populated regions, there is greater variation in dynamics across sites: populations are unaffected in some locations because the transgene has not established or has become extinct, whilst sites that have been successfully colonised by the transgene and become extinct may remain empty for several years. Suppression is therefore less variable in regions with a high density of human settlements (Additional file 1: Fig. S4). In regions that both lack a severe dry season and have a high human density, we find strong and continuous suppression is typical rather than colonisation-extinction dynamics (site ‘d’). The high human density facilitates immigration of mosquitoes from adjoining regions, which in this case is sufficient to maintain population persistence.
If the drive allele is moderately or weakly suppressing, populations rarely become extinct even in severely seasonal populations (middle and right columns of Fig. 4). In highly seasonal sites, however, populations may become so small in the dry seasons that allele frequencies are strongly affected by genetic drift (sites ‘a’ and ‘b’). During these times, the drive allele is at risk of becoming either lost or reduced to a low frequency. Loss of the drive allele appears more common in remote than well-connected sites, since a remote population may receive only a small number of mosquitoes carrying the allele in any rainy season (compare the performance of a weak allele in sites ‘a’ and ‘b’). Whether or not the drive allele tends to survive the dry season, however, the consequences of ‘dry season drift’ is to disrupt the suppression effect of the drive allele and so to reduce its efficiency. In populations with mild seasonality, alleles with medium and low strength drive cause stable suppression, as predicted by the non-spatial model (sites ‘c’ and ‘d’).
Robustness of the results
The results presented so far are based on our best estimates of model parameters. To explore the robustness of our results to uncertainty in these parameters, we next investigate the importance of four factors: mosquito dry season ecology, homing rate, dispersal rate, and egg laying rate (Fig. 5).
Dry season ecology
So far, we have assumed all human settlements contain a small amount of permanent larval habitat in addition to habitat associated with rainfall or the proximity of rivers and other water bodies. This ensures that mosquito populations are maintained in nearly all sites in the absence of drive alleles, as observed in the field. However, it is very difficult to find larval sites in many areas during the dry season, so we now suppose that populations are instead maintained by either the dispersal of adult females using high-altitude winds or the aestivation of adult females. Following our previous studies using the same underlying model [14, 18], long-distance dispersal is modelled by assuming a fraction of adult female mosquitoes are redistributed by seasonal prevailing winds. Aestivation is modelled by assuming that a fraction of adult female mosquitoes are dormant each dry season, meaning their mortality risk is reduced but they do not lay eggs. For each process, we set parameters that ensured mosquitoes were present in most human settlements in the rainy season (on average > 99.4% sites).
The different assumptions about dry season ecology had only modest effects on the average level of population suppression after 8 years across the study area (Fig. 5a). The inclusion of aestivation and long-distance dispersal both resulted in somewhat higher suppression than assuming the presence of permanent water bodies in the case of the strong and medium strength drive allele, while long-distance dispersal resulted in somewhat higher suppression in the case of the weak drive allele (Additional file 2: Table S1). The small effect sizes are chiefly because the dry season assumptions are only important in regions where seasonality is strong. Long-distance dispersal makes little difference to the suppression caused by a strong drive allele, yet tends to increase the long-term suppression caused by medium and weak drive alleles in seasonal locations (Additional file 1: Fig. S5). This is because dispersal connects large year-round populations to sites with strong seasonality, and thus reduces the role of dry season stochasticity.
Aestivation tends to reduce the initial impact of medium and weak drive alleles in seasonal locations, because the inactivity of mosquitoes over the dry season slows the spread of the drive allele (Additional file 1: Fig. S6). This effect diminishes as the drive allele spreads throughout the study area, and in most locations, aestivation either increases or has no effect on suppression after 8 years. The precise effect of aestivation on highly seasonal locations depends on how it influences the fluctuations in population size in those locations. In moderately seasonal locations, aestivation increases suppression because it reduces the dry season bottlenecks in population size that promote stochasticity. In the most strongly seasonal regions, however, aestivation results in lower suppression even after 12 years (blue areas in Fig. S6). In these locations, our model predicts the pre-intervention populations are smaller if mosquitoes rely on aestivation rather than permanent water, and the overall extent of stochasticity is thus greater. These results are sensitive to the parameters we use to model dry season ecologies and so should be treated with caution, though they emphasise that the effects of dry season ecology may be quite subtle in the most strongly seasonal locations.
Homing rate
While it is relatively straightforward to measure homing rates in a laboratory (e.g. [10, 11, 19]), it is more difficult to determine how consistent homing frequencies will be throughout the range of genetic backgrounds a drive allele will potentially encounter. The relationship between the homing rate and the impact of a gene drive is also of interest to molecular biologists who may be able to manipulate this parameter.
The non-spatial model predicts that population suppression will increase with the homing rate, a relationship also shown by the spatial model with weak or medium-strength drives (Fig. 5b). A different pattern is seen for strong drive alleles, which over the range we model (0.7–1) are less sensitive to homing rate and show peak suppression at an intermediate value (0.9). The reduction in suppression at very high homing rates occurs because strong drive causes populations to become extinct more rapidly, which can give rise to extinction-colonisation dynamics in locations that would otherwise be permanently suppressed.
Dispersal
The best estimates of inter-village Anopheles dispersal come from three West African mark-release-recapture experiments using locally caught An. Gambiae s.l. mosquitoes, whose measurements imply dispersal rates of 0.005–0.034 inter-village movements per mosquito per day in our parameterisation [18]. We took the lower value of 0.005 to be our default dispersal rate because the experiments were conducted in closely neighbouring villages, yet the possibility that this is an order of magnitude too high or too low cannot be discounted on the basis of so few data. We find that the predicted suppression after 8 years is not sensitive to the assumed dispersal rate provided it is above a threshold of approximately 0.001 movements per mosquito per day (Fig. 5c). Suppression is reduced if dispersal is below this threshold, for example, by approximately 15% if the dispersal rate is 0.0005 movements per mosquito per day.
The dispersal rate plays a larger role in determining the time to reach equilibrium (Additional file 1: Fig. S7). Unsurprisingly, drive alleles take longer to have an impact across the study area for lower dispersal rates, though the size of this effect depends on the strength of the drive allele. For a given dispersal rate, strongly suppressing drive alleles spread faster than weakly suppressing alleles, and they are also less affected by uncertainty in the dispersal rate.
Population growth rate
Finally, we vary the egg laying rate of wildtype female mosquitoes, keeping the relative fecundities of different genotypes the same. This parameter is directly proportional to the intrinsic population growth rate, Rm. Populations with low growth rates are expected to be more affected by suppression drive alleles, because they have less capacity to withstand reductions in fertility (e.g. [7]). This explains the clear negative relationship between growth rate and suppression seen for the medium and low strength drive alleles and for the strong drive allele if the growth rate is small (Fig. 5d). For higher growth rates (> 8), we see a modest positive relationship between growth rate and suppression for the strong allele. In the presence of a strong drive allele, a small growth rate hastens the extinction of local populations, which both increases suppression locally and also promotes colonisation and extinction dynamics. These counterbalancing factors give rise to the shallow trough-shape of the relationship.