Research articlesEvolution in interacting species alters predator life-history traits, behaviour and morphology in experimental microbial communities
Abstract
Predator–prey interactions heavily influence the dynamics of many ecosystems. An increasing body of evidence suggests that rapid evolution and coevolution can alter these interactions, with important ecological implications, by acting on traits determining fitness, including reproduction, anti-predatory defence and foraging efficiency. However, most studies to date have focused only on evolution in the prey species, and the predator traits in (co)evolving systems remain poorly understood. Here, we investigated changes in predator traits after approximately 600 generations in a predator–prey (ciliate–bacteria) evolutionary experiment. Predators independently evolved on seven different prey species, allowing generalization of the predator's evolutionary response. We used highly resolved automated image analysis to quantify changes in predator life history, morphology and behaviour. Consistent with previous studies, we found that prey evolution impaired growth of the predator, although the effect depended on the prey species. By contrast, predator evolution did not cause a clear increase in predator growth when feeding on ancestral prey. However, predator evolution affected morphology and behaviour, increasing size, speed and directionality of movement, which have all been linked to higher prey search efficiency. These results show that in (co)evolving systems, predator adaptation can occur in traits relevant to foraging efficiency without translating into an increased ability of the predator to grow on the ancestral prey type.
1. Introduction
Predator–prey interactions are ubiquitous across ecosystems. Predation has been widely studied at an ecological level [1–3], and recent research also shows that this interaction can be strongly altered by rapid evolution of anti-predatory defence in the prey [4] as well as by counter-adaptations in the predator [5–7], even though selection may be asymmetric, resulting in slower evolutionary change for the predator [8]. Moreover, owing to population growth–defence trade-offs, rapid evolution of the prey and adaptation to predation can result in frequency-dependent selection of defended and undefended prey types as a function of predator population size [9–11], an example of eco-evolutionary feedback dynamics. Common to this spectrum of evolutionary, coevolutionary and eco-evolutionary dynamics is that these dynamics are all driven by natural selection acting on fitness-relevant traits.
Predation can be described by three main phases, namely prey search, capture and ingestion [12]. These three phases are shaped by key traits in predator–prey systems, including those influencing offence and defence level, and all these traits can be subject to evolutionary change [13]. The offence level is determined by sensory faculties and speed enabling location and capture of prey, and defence level by the capacity for predator avoidance and escape prior to ingestion as well as physico-chemical obstruction of ingestion and digestion [12]. Adaptations in defence and offence, in turn, combined with associated trade-offs, modulate the reproduction (i.e. life-history traits) of both parties [14]. Examples abound of the study of the different phases of predation, and adaptation in both predator and prey life-history traits. For example, the timing and population dynamics of many insectivorous bird species are tightly coupled to the dynamics of their prey insect species [15]. Olive baboon sleeping site choice and behaviour (sharing sleeping sites between multiple baboon groups) in Kenya were recently linked to decreased contact and capture rate by leopards [16]. Coevolution has been hypothesized to occur between northern Pacific rattlesnakes and California ground squirrels whereby venom resistance in squirrels is matched by increased venom effectiveness in rattlesnakes based on field data supportive of local adaptation of the traits [17].
The empirical examples of evolving predator–prey interactions described above cannot be used to experimentally investigate (co)evolution in predator–prey systems owing to the long generation times of the species. By contrast, microbial systems offer a unique opportunity to study predator–prey dynamics, as they include efficient (high prey capture rate) predators and allow for high replication as well as experimental approaches capturing both ecological and evolutionary dynamics. Microbial predator–prey systems show many key characteristics found also in other predator–prey systems, such as offence by speed [18] and defence by avoidance of detection [19], escape [20] or physico-chemical obstruction of ingestion or digestion (for an overview, see [12]). Defence level has also been demonstrated to evolve in controlled set-ups [21,22]. However, to our knowledge, there exist little to no empirical studies examining offence mechanisms subject to rapid evolution in microbial predator–prey systems.
Here, we employed an experimental evolution approach to test the influence of approximately 600 generations of predator–prey interaction on predator traits, using a microbial (ciliate–bacteria) model system. Since predator–prey dynamics are characterized by the intrinsically linked dynamics of both interaction partners, we inspected the influence of both prey and predator evolution on predator traits. To find general patterns in predator traits independently of any specific prey species, as most predators have multiple prey species [23], we used seven different prey species that were all separately evolved with the predator. We expected rapid evolution of anti-predatory defence in the prey to cause impairment of predator growth [7,14]. We expected predator evolution to be weaker, in line with the life–dinner principle [8,24] positing that the prey experiences stronger selection pressure since its survival (life) directly depends on defence, while the predator can afford a certain measure of unsuccessful prey encounters (dinner postponement). Asymmetric selection can result in dynamics other than classic arms race dynamics such as frequency-dependent cycling of traits [5], which have also been observed in microbial predator–prey systems [22]. Nevertheless, instead of escalation where predators alone impose selection pressure, we expected to also observe predator evolution, since coevolution has been demonstrated to occur in bacteria–ciliate systems, in line with the Red Queen hypothesis [7,14,25].
2. Material and Methods
We studied the evolutionary dynamics of one focal predator species (the ciliate Tetrahymena thermophila) and seven of its bacterial prey species in all seven combinations of predator–prey species communities, as well as dynamics in prey species populations only. We ran predator–prey evolutionary experiments over about 600 predator generations, and assessed evolutionary effects on life history, morphology and behaviour using common garden experiments.
(a) Strains and culture conditions
The seven prey species used in this study are listed in table 1. In addition to four taxa previously used as models in predator–prey studies, three strains were chosen based on representing genera associated with ciliate predators in natural habitats or potentially exhibiting different anti-predatory defence mechanisms (table 1). Since each strain represent a single genus, strains are referred to by their genus names in the text.
| straina | rationale for species selection |
|---|---|
| Escherichia coli ATCC 11303 | model prey [26] |
| Janthinobacterium lividum HAMBI 1919 | pre-/post-ingestion defence: toxin release [12] |
| Sphingomonas capsulata HAMBI 103 | model prey [27] |
| Brevundimonas diminuta HAMBI 18 | realistic habitat [28] |
| Pseudomonas fluorescens SBW25 [29] | model prey [30] |
| Comamonas testosteroni HAMBI 403 | pre-ingestion defence: oversize [12] |
| Serratia marcescens ATCC 13880 | model prey [27] |
We used a single strain of the asexually reproducing ciliate T. thermophila 1630/1U (CCAP) [31] as a generalist predator capable of consuming all the prey species. Tetrahymena thermophila is a ciliate species characterized by a facultative sexual reproductive cycle and nuclear dualism, where the cells contain a small diploid non-expressed germline nucleus (micronucleus) and a larger highly polyploid somatic nucleus (macronucleus), derived from the micronucleus after sexual reproduction [32]. Only the macronuclear DNA is expressed and hence determines the phenotypic characteristics of Tetrahymena cells [32]. The micronucleus is only relevant for sexual reproduction. The species can be maintained either under settings of recurrent sexual reproduction, or as asexual lineages only. The Tetrahymena strain used in our experiment had been maintained in serial propagation for many years before the experiments. Sexual reproduction only occurs when induced by starvation [33], and because this was not the case during its long-term maintenance, the strain only underwent asexual reproduction. During asexual reproduction, micronuclei and macronuclei divide independently from each other [32]. It has been noted that, when cultured for a long time asexually, the micronuclei can degrade [34] and have subsequent negative effects on the genotype's fitness during a possible sexual reproduction, or even lead to genotypes losing their ability to sexually reproduce. However, given that micronuclei are never expressed and only play a role in sexual reproduction [32], and also given that we do not induce or study the genotype's ability to reproduce sexually, this possible degradation of the micronucleus does not have consequences on fitness as measured in our setting. We also note that it is a common practice to use Tetrahymena cell lines with non-functional micronuclei, as described in the standard handbook for Tetrahymena cell biology work [34]. In all of these cases, the serial propagation is not problematic as long as one is not inducing sexual reproduction. Hence, any evolution observed at the predator level in this experiment stems from either mutations or selection on existing variation in the macronuclear DNA. Furthermore, as the macronucleus is highly polyploid (n = 45), and chromosomes divide randomly during asexual reproduction [32], cells are relatively buffered to the effects of single maladaptive mutations, and can undergo relatively rapid purging of maladaptive mutations or selection for increased copies of adaptive mutations. This, together with the absence of sexual reproduction, which can be affected by serial propagation [34], makes it highly unlikely that the serial propagation set-up in the experiment would itself strongly influence the evolutionary dynamics of the predator.
Prior to the experiments, all bacterial stocks were kept at –80°C and ciliate stocks were cultured axenically in proteose peptone yeast extract (PPY) medium containing 20 g of proteose peptone and 2.5 g of yeast extract in 1 l of deionized water. During the evolutionary experiment, cultures were kept at 28°C (±0.1°C) with shaking at 50 r.p.m.
(b) Predator–prey evolutionary experiment
The evolutionary experiment was started using a small aliquot (20 µl) of a 48 h bacterial culture started from a single colony and 10 000 ciliate cells (approx. 1700 cells ml–1) from an axenic culture. Each bacterial strain was cultured alone and together with the ciliate predator (three replicates each, with the exception of six replicates for Comamonas) in batch cultures of 20 ml glass vials containing 6 ml of 5% King's B (KB) medium, with 1% weekly transfer to fresh medium.
Every four transfers (28 days), bacterial and predator densities were estimated using optical density (1 ml sample at 600 nm wavelength) as a proxy for bacterial biomass and direct ciliate counts (5 × 0.5 µl droplets using light microscopy) as used in this context and described previously [30,35,36], and samples were freeze-stored with glycerol at –20°C for later analysis. Since predators do not survive freeze-storage in these conditions, at time points 52 and 89 weeks, predator cultures were made axenic by transferring 400 µl into 100 ml of PPY medium containing an antibiotic cocktail (42, 50, 50 and 33 µg ml–1 of kanamycin, rifampicin, streptomycin and tetracycline, respectively) and stored in liquid nitrogen. Axenicity was controlled for by plating on agar plates containing 50% PPY medium, on which all the experimental bacterial strains grow. The liquid nitrogen storage protocol was modified from a previously used protocol [34] and included starving a dense ciliate culture in 10 mM Tris-HCl solution (pH 7) for 2–3 days, centrifugation (1700g, 8 min, 4°C), resuspension of the pellet in 1 ml of leftover supernatant and the addition of 4 ml of sterile 10% dimethyl sulfoxide (DMSO). The resultant solution was transferred to cryotubes in 0.3 ml lots, and frozen in a –20°C freezer at a rate of –1°C min−1 using a Mr Frosty™ Freezing Container (Thermo Scientific) for cell preservation before transferring to liquid nitrogen.
(c) Sample collection and preparation
We isolated the populations for the current experiment at time point 89 weeks (approx. 20 months). With the minimal assumption that populations multiply by 100-fold (dilution rate) until reaching the stationary phase, each weekly transfer interval represents 6.64 generations for both prey and predator [37], constituting a total minimum of approximately 600 generations. Community dynamics are shown in electronic supplementary material, figures S1 and S2 and demonstrate clear differences in population size between different prey species.
Bacteria were restored from freeze-storage by transferring 20 µl into 5 ml of 5% KB medium and culturing for 72 h. Predators were restored from liquid nitrogen by thawing cryotubes in a 42°C water bath for 15 s, followed by the addition of 1 ml of 42°C PPY medium. The cryotube contents were then transferred to a Petri dish containing PPY medium at room temperature. Upon reaching a high density (approx. 48 h), predators were transferred to 100 ml of PPY medium and cultured to a high density (approx. 7 days). To ensure that the antibiotic treatment or the liquid nitrogen storage and revival procedures do not contribute to potential differences between the ancestral predator and evolved predator lines, the axenic ancestral predator was subjected to identical procedures and was revived at the same time as the evolved lines. These culturing steps representing over 10 generations should remove the influence of non-genetic changes in predator traits caused by phenotypic plasticity [38].
(d) Physiological measurements
To test bacterial and ciliate performance and traits, we used a combination of automated video analysis, optical density measurements and flow cytometry. To separate evolutionary responses at the predator and prey level, we tested performance of both evolved and ancestral bacteria with evolved and ancestral ciliates for all evolved lines reciprocally. To do so, we prepared 12 50 ml Falcon® tubes by adding 20 ml of 5% KB medium. Three of these were inoculated with ancestral bacteria and ancestral ciliates, three with ancestral bacteria and evolved ciliates, three with evolved bacteria and ancestral ciliates and the remaining three with evolved bacteria and evolved ciliates. We placed the Falcon® tubes in a 28°C incubator, rotating on a shaker at 120 r.p.m. After inoculation, the samples were left to grow for a period of 12 days, to allow populations to grow to equilibrium density. Over the course of these 12 days, we took a total of 10 samples from each culture for analysing population density dynamics of bacteria and ciliates, and morphological and behavioural metrics for the ciliates. We sampled cultures by gently shaking the culture, to ensure it was well mixed and subsequently pipetting out 200 µl from the mixed culture.
(e) Bacterial density measurements
Bacterial density was determined both through measurement of optical density and through flow cytometry. Flow cytometric analyses were based on established protocols [39,40] that facilitate distinction between living bacterial cells and background signals (e.g. dead cells or abiotic matter). For flow cytometry, we sampled 50 µl of all cultures, diluted the samples 1 : 1000 using filtered Evian water and transferred 180 µl of the diluted samples to a 96-well-plate. We then added 20 µl of SybrGreen to stain the cells and measured bacterial cell counts using a BD Accuri™ C6 flow cytometer. As the inner diameter of the needle from the flow cytometer was 20 µm, and hence smaller than typical ciliate cell sizes, it is highly unlikely that ciliate cells were accidentally measured during flow cytometry. Also, given that bacterial densities were typically between one and five orders of magnitude larger than ciliate densities, even an occasional measurement of ciliate cells would have a negligible effect on bacterial density estimates. The full protocol can be found in the electronic supplementary material. For optical density measurement, we sampled 50 µl of all cultures, diluted 1 : 10 using filtered Evian water, and measured absorbance at 600 nm using a SpectroMax 190 plate reader.
(f) Ciliate density and trait measurements
For measuring ciliate density, we performed video analysis [41] using the BEMOVI R-package [42]. We followed a previously established method [43] where we took a 20 s video (25 frames s−1, 500 frames) of a standardized volume using a Leica M165FC stereomicroscope with circular lighting and mounted Hamamatsu Orca Flash 4.0 camera. We then analysed the videos using BEMOVI [42,44], which returns information on the cell density, morphological traits (longest and shortest cell axis length) and movement metrics (gross speed and net speed of cells, as well as turning angle distribution). The video analysis script, including used parameter values, can be found in the electronic supplementary material.
(g) Data analysis
All statistical analyses were done using the R statistical software (v. 3.5.1) [45]. To obtain the reported F- and p-values for predator traits, we performed ANOVA for the best linear models constructed for the different traits as described below.
(i) Predator trait space
To visualize whether the full set of trait data displayed structure depending on the evolutionary history of the predator and prey species, t-distributed stochastic neighbour embedding (t-SNE) was performed for each prey species separately using the Rtsne package [46] with a perplexity parameter of 3 owing to small sample size.
(ii) Beverton–Holt model fit
For analysing the population growth dynamics of the ciliates, we implemented the Beverton–Holt population growth model [47] (electronic supplementary material, figure S3) using a Bayesian framework in RStan [48], following methods used by the authors in [49,50] . This function has the form of
The priors (lognormal distribution) of the model were chosen in such a way that the mean estimates lay close to the overall observed means, but were broad enough so the model was not constrained too strongly:
-
— equilibrium population density K: ln(K) ∼ normal (9.21, 0.5),
-
— intrinsic rate of increase r0: lnr0 ∼ normal (–2.3, 0.5),
-
— rate of mortality d: lnd ∼ normal (–2.3, 0.5).
Models were run with a warm-up of 2000 iterations and a chain length of 8000 iterations.
(iii) Life-history trait analysis
We analysed the estimates of the life-history traits obtained from the Beverton–Holt model fit (r0, α and K) using linear models and model selection. We first constructed a full model with life-history traits being a function of bacterial evolutionary history (evolved/ancestor), ciliate evolutionary history (evolved/ancestor) and bacterial species (seven species factors) in a full interaction model. Next, we used automated bidirectional model selection using the step function (stats package v. 3.5.1) to find the best model. To avoid bias due to starting point, we fitted the model starting from both the intercept model and the full model, and if model selection resulted in different models, we used sample-size adjusted Akaike Information Criterion (AICc) comparison (MuMIn R-package, v. 1.42.1 [51]) to select the model with the smallest AICc value.
(iv) Morphological and behavioural trait analysis
Morphological and behavioural data were available for every time point during the growth curve, and since we know these traits can be plastically strongly affected by density [52,53], we had to take density into account in the model. We hence separated the analysis into two steps: first, we identified key points in the growth curves (early phase, mid-log phase and equilibrium density phase) and analysed the traits for these particular points. Secondly, we fitted models over all data, but taking bacterial (using flow cytometry data) and ciliate densities into account as covariates in the statistical analysis.
We defined the early phase as the second time point in the time series, equilibrium density phase as the first time point where density was larger than 99% of K, or alternatively the highest density, and the mid-log phase as the point between the early and equilibrium density phase where density was closest to 50% of K. We then created statistical models for the traits (major cell axis size, gross speed of cells and turning angle distribution) as a function of bacterial evolutionary history (evolved/ancestor), ciliate evolutionary history (evolved/ancestor) and bacterial species (seven species factors), including a full interaction for the data at the particular time point. Next, we used automated bidirectional model selection to find the best-fitting model. This was done separately for all three phases (early, mid-log and equilibrium density phases). We again performed model selection starting from both the intercept model and full model, and compared the two models using AICc comparison to identify the best model.
We then created models using all the data, where we fitted major cell axis size, gross speed and turning angle distribution as a function of bacterial evolutionary history (evolved/ancestor), ciliate evolutionary history (evolved/ancestor) and bacterial species (seven species factors), ciliate population density (ln-transformed, continuous) and bacterial population density (ln-transformed, continuous), including a full interaction. For turning angle, we also did a log10 transformation of the turning angle distributions, as fitting the model on untransformed data leads to a strong deviation on the qqplot. Next, we used automated bidirectional model selection using the step function starting from intercept model and full model, and compared the two models using AICc comparison to select the best model.
3. Results
The t-SNE maps (figure 1) showed that the evolutionary history of the predator and prey species frequently resulted in predator divergence in trait space. Importantly, this divergence evolved from a single ancestral predator population, which was subjected to co-culture with different prey species. The full results for all statistical analyses presented below to assess this divergence in detail are available in the electronic supplementary material.
Figure 1. t-SNE map of contribution of predator and prey evolutionary history to predator divergence in trait space. The traits included in the analysis encompass life history (intrinsic growth rate, equilibrium density and competitive ability), morphology (cell size and biovolume) and behaviour (speed and cell turning angle distribution). (Online version in colour.)
Prey evolution drove changes in the life-history traits of the predator, including intrinsic rate of increase (r0), equilibrium density (K) and competitive ability (α), although the presence and strength of the effect depended on the bacterial species (ANOVA, r0: prey evolution F1,78 = 15.32, p < 0.001; prey evolution × prey species F6,78 = 9.03, p < 0.001; K: prey evolution × prey species F6,80 = 13.7, p < 0.001; α: prey evolution F1,78 = 4.79, p = 0.031; prey evolution × prey species F6,78 = 5.40, p < 0.001; electronic supplementary material, tables S1–S3 and S7–S9; figure 2). The intrinsic rate of increase of ciliates (r0) was generally lower in the presence of evolved bacterial prey compared with ancestral prey, with the notable exception of Serratia, where intrinsic rate of increase was higher in the presence of evolved prey (table 2 and figure 2). For three species (Brevundimonas, Janthinobacterium and Pseudomonas), evolved predators had a higher intrinsic rate of increase (r0) on evolved prey compared with ancestral prey (figure 2). Changes in population equilibrium density (K) were highly dependent on species, with four species (Brevundimonas, Comamonas, Janthinobacterium and Serratia) showing higher population equilibrium density in the presence of evolved prey compared with ancestral prey, and the remaining three (Escherichia, Pseudomonas and Sphingomonas) showing decreased population equilibrium density in the presence of evolved prey compared with ancestral prey. Competitive ability (α) typically decreased in the presence of evolved prey compared with ancestral prey, with the exception of Pseudomonas, where competitive ability was higher in the presence of evolved bacteria compared with ancestral bacteria. Notably, for Escherichia, Janthinobacterium and Serratia, the competitive ability (α) of evolved predators was higher in the presence of evolved prey compared with ancestral prey (figure 2).
Figure 2. Reaction norms showing the effect of evolving predator–prey interaction on life-history traits of predator (data points with linear model estimate ±95% confidence intervals; N = 3, except 6 for Comamonas). The life-history traits for predators are parameters of Beverton–Holt continuous-time population models fitted to data, and include intrinsic growth rate (r0), equilibrium density (K) and competitive ability (α). The reaction norms for predators (one strain of the ciliate Tetrahymena thermophila) feeding on ancestral or evolved prey (seven bacterial strains indicated by genus name) are depicted separately for ancestral and evolved predators (colour coding). Predators evolved with a particular prey taxon have always been coupled with ancestral or evolved populations of the same taxon, while the ancestral predator is the same for all prey taxa. (Online version in colour.)
| prey species | predator evolution | r0-ratio | K-ratio | α-ratio |
|---|---|---|---|---|
| Escherichia | ancestor | 0.788 | 0.885 | 0.881 |
| Escherichia | evolved | 0.943 | 0.885 | 1.08 |
| Janthinobacterium | ancestor | 0.912 | 1.06 | 0.849 |
| Janthinobacterium | evolved | 1.09 | 1.06 | 1.04 |
| Sphingomonas | ancestor | 0.381 | 0.517 | 0.730 |
| Sphingomonas | evolved | 0.457 | 0.517 | 0.893 |
| Brevundimonas | ancestor | 0.974 | 1.18 | 0.815 |
| Brevundimonas | evolved | 1.17 | 1.18 | 0.997 |
| Pseudomonas | ancestor | 0.904 | 0.835 | 1.07 |
| Pseudomonas | evolved | 1.08 | 0.835 | 1.31 |
| Comamonas | ancestor | 0.475 | 1.26 | 0.374 |
| Comamonas | evolved | 0.569 | 1.26 | 0.457 |
| Serratia | ancestor | 1.09 | 1.16 | 0.930 |
| Serratia | evolved | 1.31 | 1.16 | 1.14 |
In contrast with life-history traits, which were affected by prey evolution alone, morphological and behavioural traits of the predator were affected by predator evolution (figure 3). However, the effect size of predator evolution was also strongly dependent on predator density (for the movement metrics gross speed and turning angles) or both predator and prey density (for the biovolume metric cell size). Evolved predators were slightly but significantly larger than ancestral predators (ANOVA: predator evolution F1,767 = 7.87, p = 0.005). Although there was a significant effect indicating that this was modulated by the evolutionary history of the prey (ANOVA: prey evolution F1,767 = 4.85, p = 0.033), the associated effect size was much smaller than predator evolution. On average, evolved predators were 39.12 µm larger than ancestral predators, and predators were on average 1.629 µm smaller in the presence of evolved prey compared with ancestral prey. The effect of predator evolution also depended strongly on prey densities (ANOVA: log prey density × predator evolution F1,767 = 6.87, p = 0.009; figure 3). The strongest differences in cell size between ancestral and evolved predators were observed at low prey densities (cell sizes 1.2–1.3 times larger for evolved compared with ancestral ciliates), whereas the effects were negligible at high prey densities (approximately equal size for evolved and ancestral ciliates; electronic supplementary material, tables S4 and S10 and figures S4–S6; figure 3).
Figure 3. Ratios of the predicted trait values of the linear models (cell size, gross cell speed and turning angle) for the evolved predator divided by the ancestral predator at different prey densities (5, 50 and 95% quantiles) and predator densities (5, 50 and 95% quantiles). Ratios represent how ciliate traits differ between evolved and ancestral ciliates, with a value of 1 meaning evolved and ancestral ciliates are identical, values larger than 1 meaning higher trait values for evolved strains, and values smaller than 1 higher trait values for ancestral ciliates. (Online version in colour.)
The gross movement speed of predators depended on the interplay between predator density and predator or prey evolutionary history. Evolved predators had, on average, up to 1.25 times higher speed compared with ancestral predators. However, this effect occurred for evolved predators at high predator densities, whereas at low predator densities, movement speed was approximately similar for ancestral and evolved ciliates (ANOVA: predator density F1,763 = 116.20, p < 0.001; predator evolution F1,763 = 1.90, p = 0.239; predator evolution × predator density F1,763 = 4.36, p = 0.037; figure 3). This effect was partially counteracted by prey evolution by driving speed to a lower rate at increasing predator densities (ANOVA: prey evolution F1,763 = 2.17, p = 0.141; prey evolution × predator density F1,763 = 5.46, p = 0.020). The movement speed of ciliate cells was also dependent on the identity of the prey species, with ciliates moving slower when subjected to three particular prey species (Janthinobacterium, Pseudomonas and Serratia; ANOVA: prey evolution F6,763 = 9.11, p < 0.001; electronic supplementary material, tables S5 and S11 and figure S7). Finally, predator evolution altered cell turning angle distribution across prey species such that evolved predator lines moved in straighter trajectories (ANOVA: predator evolution F1,56 = 10.15, p = 0.001). This effect was again highly dependent on predator population size, with evolved predators turning at approximately 0.92 times the turning rate of ancestral predators at low predator density, but turning equally as much at high predator density (ANOVA: predator density F1,763 = 33.90, p < 0.001; predator evolution × predator density F1,763 = 5.44, p = 0.02; figure 3). The effect of predator population size was also dependent on prey species, such that for three prey species (Janthinobacterium, Pseudomonas and Serratia), evolved predators moved even straighter (less turning) at higher predator densities (ANOVA: predator density × prey species F1,763 = 6.76, p < 0.001; electronic supplementary material, tables S6 and S12 and figures S8–S10).
4. Discussion
We quantified the contribution of predator and prey evolution to predator trait change across seven different prey species in a 20-month (approx. 600 predator generations) co-culture experiment. Prey evolution frequently led to changes in predator life-history traits, decreasing intrinsic growth rate, equilibrium density or competitive ability, while not affecting morphological or behavioural traits in the predator. Interestingly, the strength of the effect and the life-history trait affected depended on the prey species. These results may be influenced by different growth dynamics, defence levels or defence mechanisms of the different prey species (table 1; electronic supplementary material, figures S1 and S2) [12].
For two of the predator life-history traits, intrinsic rate of increase (r0) and competitive ability (α), the trait was impaired, with evolved compared with ancestral prey in all except for two cases (Serratia for r0 and Pseudomonas for α). This could be caused by any mechanism of prey defence evolution decreasing effective prey population size or increasing prey handling time, including cell aggregation of bacterial prey, frequently shown under ciliate predation [54,55]. While a similar result was also observed for population equilibrium density (K) with three prey species (Escherichia, Pseudomonas and Sphingomonas), intriguingly, the remaining four prey species (Brevundimonas, Comamonas, Janthinobacterium and Serratia) showed higher K in the presence of evolved compared with ancestral prey. This counterintuitive result may be caused by resource use evolution, which can occur rapidly in bacterial evolutionary experiments [37] but differ in magnitude between bacterial (i.e. prey) species. In this situation, a sufficient increase in prey population size could sustain a higher predator population size despite anti-predatory defence evolution.
Consistent with the Red Queen hypothesis, evolved predators displayed both behavioural and morphological changes linked to prey foraging efficiency. Increased swimming speed and body size were observed for evolved predators with certain prey species, and predators evolved to swim in straighter trajectories across the different prey species. Increased swimming speed and decreased cell turning (i.e. moving in straighter trajectories) have both been linked to prey search efficiency [18,56,57], and in line with this, ciliates have been shown to display decreased cell turning and increased speed at low food concentrations [58]. The role of increased body size is less clear but may also be related to increased prey search efficiency since swimming speed can be a function of body size [18,56]. All these evolutionary trait changes in the predator are consistent with being adaptations to decreased food availability owing to anti-predatory defence evolution in the prey species.
Interestingly, against our expectation based on the Red Queen hypothesis, we did not find detectable levels of adaptation in predator life-history traits when prey-evolved predators fed on their respective ancestral prey species. This could be indicative of asymmetry of selection [5,22] such that predators experience weaker selection pressure compared with prey owing to the life–dinner principle [8], whereby prey species rely on adaptation (needed to stay alive) more strongly than predators (needed to increase energy uptake). Asymmetric evolutionary change for ciliate predators could also result from smaller population size (in the order of 104 ml–1 for ciliates compared with 108 ml–1 for bacteria), larger genome size (greater than 100 Mb for T. thermophila compared with less than 10 Mb for bacteria) or more complex genomic architecture limiting adaptive mutation supply compared with the bacterial prey [59].
There are two ways asymmetric selection could account for our unexpected result regarding the lack of evolution in ciliate life-history traits. First, the offence-related traits (morphology and behaviour) where predator evolution was observed may simply not have improved sufficiently to be detectable as increased predator growth on ancestral prey using our methods. Although the culture conditions were mostly identical between the serial passage experiment and ciliate physiology measurements (same culture medium, temperature, covering 7-day time span representing serial passage culture cycle), it is also possible that minor differences in experimental conditions (different culture vials, volumes and shaking parameters) or the revival of ciliates from liquid nitrogen storage could have introduced noise in the data, masking ciliate evolution in life-history traits. Second, rapid evolution in the prey species may have changed basic features of the prey population early on in the experiment, such as causing cell aggregation, which is widely documented to evolve rapidly in similar set-ups [21,22,54,55]. An improved ability of the predator to feed on defended prey with altered characteristics may not allow for an improved ability to also feed on ancestral prey. For instance, higher speed and directionality of movement may be useful when feeding on unevenly distributed prey aggregates while not causing a benefit when feeding on prey as homogeneously distributed single cells (food being always closely available). Alternatively, as a more complex explanation, a steepening growth–offence trade-off during coevolution [14] could cause stunted growth in coevolved high-offence-level predators, which may, therefore, only display a net fitness improvement against prey in a recent evolutionary state. Since our sample material represents a snapshot from the endpoint of a long-term (co)evolutionary experiment, further experiments would be needed to assess the dynamics of predator trait change over time to test these hypotheses.
Our findings have implications for interpreting data from (co)evolving predator–prey systems. First, the pronounced impairment of predator growth traits upon prey evolution together with the lack of clear improvements in the ability of evolved predators to feed on ancestral prey types support the asymmetric selection hypothesis. Second, the occurrence of predator evolution in other key traits for predator–prey interaction despite this suggests that tracking ecological changes alone may result in an underestimation of predator evolution [60,61]. A deeper understanding of predator–prey evolutionary dynamics is, therefore, likely to critically depend on the identification and examination of key traits for the interaction, preferably over time and including both interaction partners.
Data accessibility
All code and pre-processed data needed to reproduce the ecological and evolutionary analyses are available from the Dryad Digital Repository: https://doi.org/10.5061/dryad.08kprr4zr [62].
Authors' contributions
J.C. and T.H. designed the coevolutionary experiment. J.C. performed and managed the experiment and wrote the draft manuscript. F.M., E.A.F. and F.A. designed and performed physiological measurements. F.M. and J.C. analysed data. All authors interpreted results and participated in improving the manuscript.
Competing interests
We declare we have no competing interests.
Funding
This work was funded by the Academy of Finland (T.H.; project no. 106993), the University Research Priority Program (URPP) ‘Evolution in Action’ of the University of Zurich and the Swiss National Science Foundation (grant no. PP00P3_179089, to F.A.) and the Jenny and Antti Wihuri Foundation (grant no. 00190040, to J.C.).
Acknowledgements
We thank Veera Partanen for technical help with maintaining the coevolutionary experiment and reviving samples for the physiological measurements. We thank Samuel Hürlemann for help during the laboratory work. We thank two reviewers for constructive comments on a previous version of the manuscript. This is publication ISEM-2020-16 of the Institut des Sciences de l'Evolution–Montpellier.
