Deploying dengue-suppressing Wolbachia: robust models predict slow but effective spatial spread in Aedes aegypti

A novel strategy for controlling the spread of arboviral diseases such as dengue, Zika and chikungunya is to transform mosquito populations with virus-suppressing Wolbachia. In general, Wolbachia transinfected into mosquitoes induce fitness costs through lower viability or fecundity. These maternally inherited bacteria also produce a frequency-dependent advantage for infected females by inducing cytoplasmic incompatibility (CI), which kills the embryos produced by uninfected females mated to infected males. These competing effects, a frequency-dependent advantage and frequency-independent costs, produce bistable Wolbachia frequency dynamics. Above a threshold frequency, denoted p̂, CI drives fitness-decreasing Wolbachia transinfections through local populations; but below p̂, infection frequencies tend to decline to zero. If p̂ is not too high, CI also drives spatial spread once infections become established over sufficiently large areas. We illustrate how simple models provide testable predictions concerning the spatial and temporal dynamics of Wolbachia introductions, focusing on rate of spatial spread, the shape of spreading waves, and the conditions for initiating spread from local introductions. First, we consider the robustness of diffusion-based predictions to incorporating two important features of wMel-Aedes aegypti biology that may be inconsistent with the diffusion approximations, namely fast local dynamics induced by complete CI (i.e., all embryos produced from incompatible crosses die) and long-tailed, non-Gaussian dispersal. With complete CI, our numerical analyses show that long-tailed dispersal changes wave-width predictions only slightly; but it can significantly reduce wave speed relative to the diffusion prediction; it also allows smaller local introductions to initiate spatial spread. Second, we use approximations for p̂ and dispersal distances to predict the outcome of 2013 releases of wMel-infected Aedes aegypti in Cairns, Australia, Third, we describe new data from Aedes aegypti populations near Cairns, Australia that demonstrate long-distance dispersal and provide an approximate lower bound on p̂ for wMel in northeastern Australia. Finally, we apply our analyses to produce operational guidelines for efficient transformation of vector populations over large areas. We demonstrate that even very slow spatial spread, on the order of 10-20 m/month (as predicted), can produce area-wide population transformation within a few years following initial releases covering about 20-30% of the target area.

pling drift or demographic stochasticity, the variance in position increases by ~1 N in one dimension, where population density is N. In contrast, if f < pf '[0] everywhere, then the wave is pulled at a rate determined by the individuals at the advancing front; individuals behind the wave almost all descend from this small fraction (Brunet et al. 2006). Now, the spread of the wave is sensitive to the form of longrange dispersal, and to the population density. Provided that the dispersal distribution is bounded by an exponential, the wave will settle to a steady shape, with a speed that depends only on the growth rate from low density; random fluctuations will slow it by an amount that depends logarithmically on density (Brunet et al. 2006). If the rate of long-range dispersal is faster than exponential, then the wave will accelerate, though again, random fluctuation will limit this acceleration, and cause the wave to fragment.

Effect of dispersal distribution on wave speed: deterministic models
We illustrate these points using the simple cubic model: where p + q = 1. Note that the model only makes sense for sufficiently small s : it can be taken as a weak-selection approximation to a variety of more detailed models (Barton and Turelli, 2011). However, it is accurate even for quite large s, and so we take it as a surrogate for a much wider class of models.
Because we will be considering the case p / < 0, it is more natural to use the parameter in this section.
Without loss of generality, we assume p / < 1 2 , or = 1 -− 2 p / > 0. If p / > 0, < 1, there is an unstable equilibrium, whereas if p / < 0, > 1, the allele can increase from low frequency. (p / = 0 represents selection on a strictly recessive allele). If < 3 or p / > -−1, then f [p] ≥ pf '[0] for some p, and so a pushed wave solution exists. However, a pulled solution also exists, with speed c = s/∕4 . The transition Because we will be considering the case p / < 0, it is more natural to use the parameter in this section.
Without loss of generality, we assume p / < 1 2 , or = 1 -− 2 p / > 0. If p / > 0, < 1, there is an unstable equilibrium, whereas if p / < 0, > 1, the allele can increase from low frequency. (p / = 0 represents selection on a strictly recessive allele). If < 3 or p / > -−1, then f [p] ≥ pf '[0] for some p, and so a pushed wave solution exists. However, a pulled solution also exists, with speed c = s/∕4 . The transition between a pulled and pushed wave occurs when the speeds of the two solutions intersect, at = 2, or p / = -− 1 2 . For larger , the deterministic solution is "pulled", with speed c = s( -− 1) = 2 s-−p /  (Fisher 1937). Figure 1 shows how the wave speed increases with , measured relative to the prediction for a pushed wave, s/∕4 (horizontal line). For Gaussian dispersal, the speed is very close to the diffusion approximation for <2, and close to the prediction for a pulled wave, c = s( -− 1) , for > 2 (see below). A reflected exponential ("Laplace") distribution, which has substantially fatter tails, is also close to the diffusion approximation for a pulled wave (blue dots to right), but is slightly slower than predicted for a pushed wave for <1. This is because genes that move a long way will be at a selective disadvantage, and so will be lost: the wave speed depends on the shape of the bulk of the dispersal kernel, which for a given variance, is narrower than a Gaussian (Fig. 2). The exponential square root distribution has even fatter tails, and shows a correspondingly slower speed for <1. For ≥2, this fat-tailed distribution gives a higher speed, but still settles to a travelling wave with constant speed. This is at first puzzling, since distributions with fatter than exponential tails are predicted to give an accelerating wave. However, numerical calculations were made with a dispersal distribution truncated at ± 20 demes, or ±10 standard deviations (black dots). If this truncation is increased to 50 or 100 demes (purple and red dots), this makes little difference to the wave speed for < 1 (speeds are indistinguishable for <1), but substantially increases the speed for > 1, as predicted. However, the speed is only increased by ~ 20% even when individuals can disperse out to 100 demes, or ± 50 . Moreover, even for the exponential square root distribution, only an extremely small fraction of the distribution is excluded by truncation at 25 . Arguably, it would be unrealistic to calculate for yet longer tails, since real ranges are finite, and since over such large distances, stochastic effects will dominate even in very large populations, as we discuss below. Figure 1. Wave speed, c, relative to that expected for a pushed wave, s/∕4 , is plotted against = 1 -− 2 p / . The horizontal line at 1 is the expectation for a pushed wave, whilst the downward curve at the right is the ratio expected for a pulled wave, whose speed is c= s( -− 1) . Dots show numerical calculations for three forms of dispersal: Gaussian (green), Laplace (blue), and exponential square root (black); the dispersal distribution is truncated at ±20 demes, and the standard deviation adjusted to equal =2. For the exponential square root, results are also shown for truncation at ±50 demes (purple) and ±100 demes (red). The selection coefficient is adjusted so that the maximum rate of change is always 0.025; this is equivalent to constant directional selection of s = 20 %. With these parameters, the deme spacing has negligible effect.  Table 1. The mass and the fraction of the total variance that is excluded by truncating distributions at 10 ( Fig. 1, green, blue, black), 25 (Fig. 1, purple) and 50 (Fig. 1, red).

Gaussian
Laplace ExpSqrt Wave speed with infinite variance In Fig. 1, we compared dispersal distributions that were scaled to have the same variance. This approach fails for distributions such as the Cauchy, that fall away so slowly that they have infinite variance. However, we can find the Cauchy distribution that yields the same wave speed as a Gaussian with (say) = 2 (Fig. 3). On an infinite range, a Cauchy distribution leads to a wave with everincreasing speed. However, the Cauchy distribution must necessarily be truncated for numerical calculations, and for any given truncation point, the wave does settle to a steady speed, which can be matched to that of a Gaussian. Figure 3 plots the standard deviation of the truncated Cauchy against = 1 -− 2 p / , for truncation at 10 … 100 (black … red dots). For pushed waves ( < 2), fat tails slow down the wave, because long-distance migrants are lost, and so the truncated Cauchy must have a higher variance to maintain the same speed ( < 2, left of Fig. 3). For a pulled wave, however, a truncated Cauchy with a very small variance can maintain high speed, provided that it is truncated sufficiently far out. In the most extreme case shown here (truncation at 100 , =4; bottom right of Fig. 3), 0nly a fraction 7.9*⋆10 -−8 disperse away from their native deme, contributing variance ~10 -−5 , and yet this maintains a wave with the same speed as Gaussian dispersal with standard deviation =2. (Note that for these extreme cases, the choice of deme spacing does affect the results, though not qualitatively).
to that of a Gaussian. Figure 3 plots the standard deviation of the truncated Cauchy against = 1 -− 2 p , for truncation at 10 … 100 (black … red dots). For pushed waves ( < 2), fat tails slow down the wave, because long-distance migrants are lost, and so the truncated Cauchy must have a higher variance to maintain the same speed ( < 2, left of Fig. 3). For a pulled wave, however, a truncated Cauchy with a very small variance can maintain high speed, provided that it is truncated sufficiently far out. In the most extreme case shown here (truncation at 100 , =4; bottom right of Fig. 3), 0nly a fraction 7.9*⋆10 -−8 disperse away from their native deme, contributing variance ~10 -−5 , and yet this maintains a wave with the same speed as Gaussian dispersal with standard deviation =2. (Note that for these extreme cases, the choice of deme spacing does affect the results, though not qualitatively). Figure 2. The standard deviation, *⋆ , of the truncated Cauchy distribution that gives the same speed as a Gaussian with =2 (horizontal line); this is plotted against . The distribution is truncated at 10 , 25 , 50 , 75 , 100 (black, blue, brown, purple, red).

Random fluctuations
With long-tailed dispersal, deterministic models are misleading, because the numbers of long-range migrants will be small even in very dense populations. Indeed, even with Gaussian dispersal, finite population size has an appreciable effect on the speed of a pulled wave, because the dynamics depend on regions where alleles or individuals are very rare. Here, we simulate a finite population on an infinite range. The population is represented as a list of the numbers of the focal allele in each deme; all demes to the left are assumed to be fixed, and all to the right are assumed to be at zero. In each generation, and in each deme, there is selection followed by random sampling of 2N genes, following the Wright-Fisher model. Dispersal is implemented as follows. For each of the n p demes that are currently polymorphic, the 2N genes are drawn from a parent population according to the dispersal function; the continuous dispersal distribution is rounded to the nearest integer. (For =2, as used here, this causes negligible error). These genes can be from arbitrarily far away. The same procedure is followed for 1000 demes to the left of the current set, and 1000 demes to the right. Finally, the list is trimmed to remove fixed demes on the left and on the right (but tracking the location of the leftmost polymorphic deme, which will tend to move to the right as the wave advances. This scheme is extremely efficient, since the random distances moved can be taken as a single draw of 2Nn p random numbers from the dispersal distribution. The only approximations are that individuals live on a discrete grid, and that only ±1000 demes are followed on either side. Table 2 summarises the effects of finite population size. For a pushed wave ( =0.5), random drift has little effect, unless deme size is very small (2 N = 5). For a pulled wave, there is a substantial slowing, even for demes of 2 N = 1000 individuals. Note that this stochastic effect means that long range dispersal no longer increases the speed of advance, even for the exponential square root. There is no sign that rare propagules move far ahead of the initial wave, producing the kind of stochastic advance that is expected for very fat-tailed distribution. The key conclusion is that random fluctuations have very little effect on the speed of "pushed" waves. The slowing effect is greater when waves are "pulled", and there is long-range dispersal. However, even then the effect is not as dramatic as suggested by the asymptotic theory. Table 2 summarises the effects of finite population size. For a pushed wave ( =0.5), random drift has little effect, unless deme size is very small (2 N = 5). For a pulled wave, there is a substantial slowing, even for demes of 2 N = 1000 individuals. Note that this stochastic effect means that long range dispersal no longer increases the speed of advance, even for the exponential square root. There is no sign that rare propagules move far ahead of the initial wave, producing the kind of stochastic advance that is expected for very fat-tailed distribution. The key conclusion is that random fluctuations have very little effect on the speed of "pushed" waves. The slowing effect is greater when waves are "pulled", and there is long-range dispersal. However, even then the effect is not as dramatic as suggested by the asymptotic theory. Table 2. The wave speed in a finite population. This is estimated as the mean speed over 2000 generations. The last row shows the results for the deterministic model ( Fig. 1) with k max = 20 for the Gaussian and Laplace distributions, and k m = 200 for the exponential square root.