Genotyping
All founder strains were genotyped for the mitochondrial cytochrome c
oxidase I (COI) gene, by extracting DNA from 10 cells per strain and
applying a PCR amplification protocol and sequence analysis (56). At
cycle 30 in the evolutionary experiment, mixes of 50 cells from each
line were processed in the same way and analysed for (multiple) COI
marker signals. This method characterises the line for the most frequent
COI genotype, and it has a resolution threshold of c. 5%, i.e., it can
detect 2-3 cells of a minority genotype in the sample (56). The
sequences used are deposited in GenBank (50, 57).
Range expansion modelOur model is designed to capture the specificities of the evolutionary
experiment and the characteristics of the strains in the founder
population. Thus, we model the population dynamics of Paramecium
strains, assuming logistic growth following the Verhulst equation
expanded to include both intra- and inter-strain competition:
\begin{equation}
\frac{\text{dN}_{\text{i\ }}}{\text{dt}}=(r_{0,i}-\sum_{j}{(\alpha_{\text{ij}}N_{j}}))\ N_{i}\nonumber \\
\end{equation}where \(N_{i}\) is the population size of strain i, \(r_{0,i}\) is its
intrinsic rate of increase and \(\alpha_{\text{ij}}\) as the competition
coefficients. The model is parametrized with the posteriors extracted
from growth curve fitting (r0, \(\overline{N}\)), as
described above. We make the simplifying assumption that intra- and
interspecific competition is of equal strength. We further assume a
quasi-extinction threshold of 0.7 (the bottleneck occurring during
dispersal, we have tested the effect of different quasi-extinction
thresholds ranging from 0.0001 to 0.9), which implies that strains
experience an extinction if they exhibit densities below this value.
We model the community dynamics of the strains for 7 days, followed by a
3h dispersal phase in a 2-patch metapopulation. During this 3h dispersal
phase all strains can disperse from their patch of origin to their
destination patch according to the dispersal rates estimated from
dispersal assay; the model is parametrized with posteriors extracted
from the statistical analysis described above. After the dispersal phase
we follow the patch of origin (residents in the range core treatment),
the destination patch (dispersers in the front treatment) or the
combined patches (dispersers and residents mixed in the control
treatment). We repeat this procedure for a total of 10 iterations. As in
the experiment, we control for densities between rounds of iteration by
selecting the equivalent of 10 mL samples.
This approach allows us to predict, based only on measurements of growth
parameters and dispersal rates of the founder strains, which strains
should predominate in each of the three treatments at the end of the
experiment. It is important to keep in mind that the underlying model is
deterministic. However, since we parametrize the model with draws from
posteriors, our approach takes into account the uncertainty associated
with the data and yields a distribution of likely outcomes, given these
uncertainties. Note that our model depicts a scenario of selection from
standing genetic variation; it does not include mutational change.
Statistical analysisStatistical analyses were performed in R (ver. 4.1.2) and JMP 14 (SAS
Institute Inc. 2018). We analysed dispersal (proportion of dispersers),
using generalised linear models (GLM) with binomial error distribution.
We considered evolutionary treatment (core, front, control),
experimental cycle and line (nested within treatment) as explanatory
variables. We analysed variation in intrinsic population growth rate
(r0) and equilibrium density (\(\overline{N}\); averages
per line and year) using GLMs, with evolutionary treatment, year and
line as explanatory variables. To illustrate how selection acts on
standing genetic variation in our model, we associated the winning
probability of each of the 20 founder strains with their respective
median values of dispersal, r0, and \(\overline{N}\)from the distributions used by the model. For each treatment, we then
performed multiple regressions, with winning probability as response
variable and the three traits as explanatory variables. To investigate
associations between dispersal, r0 and \(\overline{N}\),
we constructed a data matrix based on trait means per year and line (3
years x 15 lines, n = 45), after centering and scaling trait
distributions. One range-front line was lost in year 3, leading to n =
44. We used a Bayesian approach with the “rstan” package version
2.21.2 (58) to estimate pairwise trait correlations (SI Appendix, Fig.
S3). We also performed a principal component analysis (PCA), considering
the joint variation in all three traits.
AcknowledgementsGZ was supported by a grant from the Agence Nationale de Recherche (n°
ANR-20-CE02-0023-01) to OK. This is publication ISEM-XXXX-XXX of the
Institut des Sciences de l’Evolution. Florent Deshors and Thomas
Teissedre assisted with movement and growth assays. Alison Duncan, Flore
Zelé and Alexey Potekhin gave helpful comments for the experimental set
up and the interpretation of the results.