Robustness and predictability of evolution in bottlenecked populations

Deterministic and stochastic evolutionary processes drive adaptation in natural populations. The strength of each component process is determined by the population size: deterministic components prevail in very large populations, while stochastic components are the driving mechanisms in small ones. Many natural populations, however, experience intermittent periods of growth, moving through states in which either stochastic or deterministic processes prevail. This growth is often countered by population bottlenecks, which abound in both natural and laboratory populations. Here we investigate how population bottlenecks shape the process of adaptation. We demonstrate that adaptive trajectories in populations experiencing regular bottlenecks can be naturally scaled in time units of generations; with this scaling the time courses of adaptation, fitness variance, and genetic diversity all become relatively insensitive to the timing of population bottlenecks, provided the bottleneck size exceeds a few thousand individuals. We also include analyses at the genotype level to investigate the impact of population bottlenecks on the predictability and distribution of evolutionary pathways. Irrespective of the timing of population bottlenecks, we find that predictability increases with population size. We also find that predictability of the adaptive pathways increases in increasingly rugged fitness landscapes. Overall, our work reveals that both the adaptation rate and the predictability of evolutionary trajectories are relatively robust to population bottlenecks.

Figure 1

Global optimum, antipode of global optimum and global minimum fitness values, versus sequence length, $L$. The analytical prediction [Eq. (10), evaluated numerically] for the global optimum fitness, $E\left[f_{\max }\right]$, is shown for $K=L-1$ (magenta solid line), along with the analytical prediction of $E\left[f_{\max }\right]=2/3$ for $K=0$ (magenta dotted line). Analogous analytical predictions for the global minimum fitness are shown for comparison (black lines). Simulation results are shown for comparison for $K=$ 1, 2, 3 and 4 (magenta and black symbols as indicated). Simulation results for the expected fitness of the antipode of the global optimum are shown in blue, along with the asympotic expectation (0.5 for large $L$ and $K=L-1$, dashed line). Simulation results show the mean across 100,000 randomly generated fitness landscapes in each case. Error bars for simulation results are similar to symbol heights and omitted for clarity.

Figure 2

Fitness trajectories, that is, mean population fitness versus time for different bottleneck ratios (upper panels) and fitness versus bottleneck size at different times (lower panels). Time is measured in units of bottlenecks (left panels) and doublings (right panels). The parameter values are mutation rate $U=1\times {10}^{-4}$, sequence size $L=8$, epistasis parameter $K=2$, and $N_f$ is set at $N_f=2^{15}=32768$. The bottleneck sizes are indicated in the legends. In the bottom panels the curves correspond to fixed numbers of bottlenecks (or doublings) as indicated in the legends.

Figure 3

Mean population fitness, fitness variance and change in fitness $\mathrm{\Delta }f$ as a function of time. Time is expressed in units of doublings. The parameter values are $N_f=2^{15}=32768$, mutation rate $U=1\times {10}^{-4}$, sequence size $L=8$ and epistasis and $K=2$. The bottleneck sizes are $N_0=32$ (blue dashed-lines), $N_0=4096$ (orange dashed-lines) and $N_0=16384$ (green dashed-lines). $\mathrm{\Delta }f$ is simply the mean population fitness at time $t+1$ minus the mean population fitness at time $t$, for $t$ in the units indicated.

Figure 4

Hill diversity numbers $D_0$, $D_1$, and $D_2$ versus time. Time is expressed in units of bottlenecks (doublings) on the left (right) panels. The measures of diversity are presented at the end of the growth phase (solid lines) and just after the bottleneck protocol (dashed lines). The parameter values are $N_f=2^{15}=32768$, mutation rate $U=1\times {10}^{-4}$, sequence size $L=8$ and epistasis parameter $K=2$. The bottleneck sizes are $N_0=32$ (blue lines), $N_0=4096$ (orange lines) and $N_0=16384$ (green lines) as indicated in the legends. The symbol $a$ in the legend means just after bottlenecks, whereas $b$ means just before bottlenecks.

Figure 5

Predictability and mean path divergence. In the upper panels both quantities are shown as a function of the population size at the end of the growth phase $N_f$. In these panels the population size after the bottleneck is set at $N_0=32$. In the lower panels both quantities are shown as a function of the population size after the bottleneck $N_0$. In these panels the population size $N_f$ is set at $N_f=4096$, whereas the mutation rate is set at $U=5\times {10}^{-2}$ and sequence size at $L=8$. The epistasis parameter $K$ is displayed in the legends.

Figure 6

Predictability with respect to the ending points, $P_{2,\mathrm{ending}}$. In all plots, the population size at the end of the growth phase is $N_f=32768$. In the left panel, the mutation rate is $U=0.005$ whereas the sequence size $L$ and epistasis parameter $K$ are both varied such that the correlation is $\gamma =0.5$. In the right panel, the sequence size is $L=8$ and the epistasis parameter is set at $K=3$. The Hamming distance from starting points to the global optimum is five, $d_{\mathrm{GO}}=5$. The simulation data refers to an everage over 10 distinct fitness landscapes, and 10 random starting points for each landscape. The dashed-lines correspond to the predictability with respect to ending points for random adaptation walks (RAW) and for S-weighted walks (SWW).

Figure 7

Multidimensional scaling plot of the evolutionary pathways. In the upper panels $N_f$ is set at $N_f=4096$, whereas $N_0=64$ (left upper panel) and $N_0=2048$ (right upper panel). In the lower panels, $N_0$ is set at $N_0=32$, whereas $N_f=128$ (left lower panel) and $N_f=4096$ (right lower panel). The data correspond to a fixed fitness landscape with epistasis parameter $K=2$. Those evolutionary pathways that achieve a frequency higher than 0.05 are highlighted in the plot (dark circles, with numbers indicating path frequency).