Escherichia coli K12 (MG1655) was evolved for 14 days to three different antibiotics: Amikacin sulfate (AMK) (Sigma), Piperacillin sulfate (PIP) (Sigma), and Tetracycline hydrochloride (TET) (Sigma), covering three major classes of antibiotics, including both bactericidal and bacteriostatic drugs. The antibiotics were dissolved in water (10 mg/l) and the stock solutions were stored at −20°C. Four replicate lineages were evolved in parallel for each drug. 96-well plates (Almeco), containing 1 ml Mueller-Hinton broth II (MHBII) (Sigma) per well and a 2-fold antibiotic gradient in 10 dilutions, were prepared at the start of the experiment and stored at −20°C. The minimal inhibitory concentration (MIC) of the wild type, as defined by the European Committee on Antibiotic Susceptibility Testing (EUCAST), was located in the second well, allowing growth of the wild type in the first well under sub-inhibitory conditions (exact plate setup and drug concentrations are given in Supplementary Table 1). Plates were defrosted at the day of usage, pre-heated to 37°C, inoculated with 50 μl of freshly growing cells and incubated at ~900 r.p.m. and 37°C for 22 h. One hundred fifty microliters of each well were transferred into a 96-well microtiter plate and the optical density was measured at a wavelength of 600 nm (OD₆₀₀) by an ELx808 Absorbance Reader (BioTek). Based on the OD measurement a cut-off value, that was the minimal growth that clearly set itself apart from the background growth, was chosen to define distinct growth for each drug (Supplementary Figure 1). An OD₆₀₀ > 0.1 corresponding to ~8.0 × 10⁷ CFU/ml defined distinct growth for AMK and TET and an OD₆₀₀ > 0.3 equivalent to about 2.4 × 10⁸ CFU/ml defined growth for PIP due to a background growth level of around OD₆₀₀ = 0.18. Fifty microliters of the well with the highest drug concentration that showed distinct growth in the deep-well plate were used to inoculate a fresh gradient (exact OD₆₀₀ values and corresponding drug concentrations of the well chosen for each transfer are given in Supplementary Table 2). Remaining cells in these wells were mixed to a final glycerol concentration of 20% and stored at −80°C. On each plate 16 wells served as negative control resulting in a total of 448 wells during the course of the experiment, of which 1% showed growth. For each lineage seven colonies were isolated for genomic and phenotypic characterization from the population that had been maintained for two passages at or above the clinical breakpoint as defined by EUCAST for the specific antibiotic. The clinical breakpoint is the drug concentration that is used as a cut-off value to classify pathogens as susceptible or resistant toward a specific drug (Turnidge and Paterson, 2007).

E. coli K12 (MG1655) was not only evolved in drug gradients but also to a daily relative increase of drug concentration. Three different approaches were used (exact drug concentrations for each day for the different increment approaches and drugs are given in Supplementary Table 3). The lineages in the “Increment 100” setting were exposed to a 100% increase in drug concentration. Under this regime the drug concentration was consequently doubled every day, applying a constantly strong selection pressure to the lineages. The clinical breakpoint was supposed to be reached after 7 days of the ALE experiment. “Increment 50” lineages were also exposed to a rather high environmental change rate by growth in a 50% higher drug concentration every day, reaching the clinical breakpoint on the 9th day of the experiment. The drug concentration was raised by 25% for the “Increment 25” lineages, allowing a mild selection and twice as much time to adapt to the clinical breakpoint concentration compared to the “Increment 100” lineages. Eight lineages were evolved in parallel in each setting to AMK, PIP, and TET. The experiment was designed that all experimental setups reached the MIC as defined by EUCAST at the 4th day. As control for media adaptations eight wild type lineages were evolved to the media without antibiotics. The MHBII antibiotic mixture was prepared in falcon tubes for each day, drug and experimental setup in the beginning of the experiment and stored at −20°C. The drug containing media was defrosted on the day of usage and pre-heated to 37°C. The lineages were grown in 1 ml MHBII and antibiotic in a 96-deep-well dish for 22 h at 37°C and ~900 r.p.m. 50 μl of cells were transferred every 22 h. The remaining cells were mixed to a final concentration of 20% glycerol and stored at −80°C. The adaptive evolution was stopped after 14 days when the “Increment 25” lineages had passed the clinical breakpoint. During the experiment ~3% of 1,152 negative controls showed growth. Cells were streaked on LB plates and identified by visual investigation as E. coli. All colonies looked identical, suggesting that there was no contamination. Once the adaptive evolution experiment was ended, lineages that were adapted to the clinical breakpoint were streaked on LB agar and seven isolated colonies were used for further analysis. If lineages died out before they had reached the clinical breakpoint, the last possible time point was chosen.

Isolated colonies were used to inoculate a 96-well microtiter plate containing 150 μl MHBII. About 10⁵ cells of an overnight culture were transferred with a 96-pin replicator to 10 dilutions of a 2-fold drug gradient spanning from 0.5 to 256 mg/l AMK and 0.25 to 128 mg/l of PIP or TET, respectively. For each lineage one isolated colony was tested in two technical replicates against all three antibiotics. Eight inoculated wells containing MHBII served as positive control while eight wells only filled with MHBII served as negative control. The plates were incubated at 37°C and 900 r.p.m. for 18 h and subsequently the OD₆₀₀ was measured by an ELx808 Absorbance Reader (BioTek). The data was further analyzed using R (Team R Core, 2014). The average OD₆₀₀ values of the negative controls were subtracted from all remaining OD₆₀₀ values. Percent inhibition was calculated by subtraction of the OD₆₀₀ values divided through the average of the OD₆₀₀ values of the positive controls from 1.

A dose-response curve was fitted to the values using a logistic model from the drc package (Ritz and Streibig, 2005), with x for the molar drug concentration and default values for the other variables, where b describes the steepness of the curve, c and d the lower and upper asymptotes and e the effective dose (Munck et al., 2014; Ritz et al., 2015):

Dose-response curves were plotted with the package ggplot (Wickham, 2009). The drug concentration causing 85% growth inhibition (IC₈₅) was calculated with the inverted function, normalized to the wild type, and plotted grouped by drug using ggplot (Wickham, 2009). A non-parametric distribution of parallel lineages was observed for some experimental setups and drugs wherefore the non-parametric Kruskal-Wallis one-way analysis of variance was applied in R to detect significant (P < 0.05) differences between the four experimental setups adapted to each drug.

A 96-well microtiter plate, containing 200 μl MHBII per well, was inoculated with cells in exponential growth phase using a 96-pin replicator. All seven isolated colonies per lineages were included in the growth measurement. OD₆₀₀ was measured in a ELx808 Absorbance Reader (BioTek) every 5 min for 10 h at 37°C and 650 r.p.m. The data was analyzed with R (Team R Core, 2014). The growth rate was calculated based on the steepest part of the growth curve during exponential growth. The doubling time was normalized to the wild type and plotted grouped by drug with the package ggplot (Wickham, 2009). To test whether the observed differences in growth rate were significant (P < 0.05) the non-parametric Kruskal-Wallis one-way analysis of variance was applied in R.

One isolated colony of each of the 92 lineages, that was also used to measure the IC₈₅ and growth rate, was grown in LB and DNA was extracted with the A & A Genomic Mini kit (A&A Biotechnology). The DNA was sent to Macrogen, who prepared genomic libraries with the TruSeq DNA Nano (550 bp) kit from Illumina and sequenced them by Illumina MiSeq 300 bp paired ends.

The reads were aligned to E. coli K12 U00096 reference genome with CLC Genomics Workbench (Escherichia coli Str. K-12 Substr, 2014). On average each base was covered at least 37.5 times (Supplementary Table 4). SNP and INDEL sites were determined with CLC Genomics Workbench. Only loci with a phred score of 30 at the variable position as well as at the three neighboring bases occurring at least with a frequency of 80% were included in the analysis. Single nucleotide polymorphisms (SNPs) that were detected in all lineages including the media adapted wild type were considered mutations that had occurred before the start of the experiment and were therefore excluded. If two lineages shared two identical SNPs they were considered cross-contaminated and one of them was excluded from further analysis.

Large deletions were identified with a workflow in CLC workbench. The reads of all genomes were assembled de novo and used as reference genomes to map the reads from one of the media adapted wild type strains. The reads from the media adapted strain that did not map to the de novo assembled antibiotic evolved genome were collected and also de novo assembled into contigs. Contigs larger than 1 kb with a coverage > 30 were considered large deletions.

We downloaded both the complete genome of E. coli K12 MG1655, with accession U00096, and the corresponding ORFs from the NCBI Nucleotide Archive. The ORFs were then clustered together using cd-hit (Li and Godzik, 2006), allowing ORFs with at least 90% identity and coverage to be in the same cluster. The reads of the sequenced strains were filtered using the FASTX-Toolkit package and a minimum quality of 30. The quality-filtered reads were blasted against the clustered ORFs by setting the word_size to 10 and the e-value to 100 allowing accurate short sequence mapping. Reads mapping to at least two ORFs of different clusters with a minimum of 30 and a maximum of 70% coverage of each cluster were further considered potential large insertion reads. These reads were blasted against the whole genome to confirm that the two genes were not adjacent on the genome. Reads that covered the reference genome with more than 90% were expected to contain two adjacent genes, wherefore they were excluded from the analysis. It was further explored whether the remaining hits overlapped with insertion elements. Reads that aligned to one ORF annotated as insertion element were included for further analysis. Only insertion reads that were detected in more than 10 reads were considered as large DNA insertions.

Gene duplications were identified with CLC Genomics Workbench using the coverage analysis tool by calling for regions with at least 700 bp of significantly (P < 0.001) increased coverage. The identified regions were mapped to the reference genome U00096 in R and genes that overlapped with at least 95% with the region of high coverage were identified as duplicated genes.

To compare the effect of the ALE methodology on final pheno- and genotypes we evolved E. coli to three different antibiotics using four different ALE strategies. The antibiotics chosen for this experiment represent three major groups of antibiotics. Two of the drugs, amikacin (AMK) and tetracycline (TET), target the ribosome with the former being bactericidal and the latter being bacteriostatic. The third drug, piperacillin (PIP), is a bactericidal drug targeting cell-wall biosynthesis. The four different selection regimes can be divided into two categories: (1) A gradient approach in which the population that tolerates the highest drug concentration is passed to a fresh drug gradient every 22 h; and (2) an increment approach in which the evolving population is passed every 22 h to a new drug concentration increased by fixed increments (Figure 1).

A 20-fold dilution of the gradient evolved populations was chosen based on a model by Wahl et al. (2002), to allow a high variation in the transferred population and to base fixation of mutations on the optimal adaptation rather then on limiting bottlenecks. According to Wahl et al. (2002) a 10-fold dilution of the population would be optimal but previous experiments in the laboratory showed a small increase in the inhibitory drug concentration for AMK and TET when a 10-fold dilution was used. To avoid an inoculum effect, which is a significant increase in the inhibitory concentration caused by a larger amount of organisms in the inoculum (Brook, 1989), we chose a 20-fold dilution (Supplementary Figure 2).

After 14 days of adaptive evolution in the gradient setup, the populations exposed to AMK tolerated on average 512 mg/l of the drug, corresponding to a 170-fold increase compared to the media adapted wild type (P = 2.89654E⁻²⁷, student's t-test) (Figure 2A). The PIP evolved lineages grew in drug concentrations of 192 mg/l on day 14 of the ALE experiment, equal to a 80-fold increase compared to the media adapted wild type (P = 0.00109527, student's t-test). However, large oscillations in resistance were observed for the PIP evolved lineages during the course of the experiment (Figure 2B). This variation could be explained by an inoculum effect, which is more frequently observed for beta-lactam antibiotics (Eng et al., 1984, 1985; Brook, 1989). Lineages evolved to TET did not reach the same drug tolerance compared to AMK or PIP evolved lineages, but still grew in 15 mg/l TET, exceeding the media adapted wild type inhibitory concentration (IC) by 15 times (P = 2.04072E⁻⁰⁶, student's t-test) (Figure 2C). These values are in accordance with previous findings where the IC₉₀ values (the drug concentration at which growth of 90% of the population is inhibited) of isolated colonies were determined after 14 days of adaptive evolution in a gradient system using the same strain and drugs (Munck et al., 2014). Only PIP evolved lineages appear more resistant in the present study, which is attributed to the inoculum effect caused by the larger passaging volume.

In the increment approach three different rates of environmental change were tested (Figure 1), for which the selective pressure (e.g., antibiotic concentration) was increased by 100, 50, or 25% every day. Similar to the gradient approach the 100% increment setup applies a high selection pressure with the risk of exceeding the adaptive potential of the bacterium leading to extinction of the lineages (Figure 1).

We speculated that the increment-evolved lineages would die out when the drug concentrations they were exposed to exceeded the adaptation level of the gradient evolved lineages as a maximal selection pressure was applied to the gradient evolved lineages. Our observations during the ALE experiment support this hypothesis (Figures 2A–C). Whenever the antibiotic concentration of the drugs in the increment evolution experiments exceeded the maximal concentration of at least one of the gradient evolved lineages, some of the parallel increment lineages became extinct. For instance, the remaining 100% increment lineages adapted to AMK became extinct after they grew in 32 mg/l on day 8 of the experiment but one of the gradient adapted lineages only grew in 16 mg/l. All increment lineages died out when exposed to a drug concentration above the drug concentration to which all of the gradient lineages were adapted to (Figures 2A–C). For example all remaining increment 50% lineages died out when they passed an PIP concentration of ~32 mg/l on day 13 of the adaptive evolution experiment when lineages adapted in the gradient approach only grew in PIP concentrations of 8 mg/l. Since the antibiotic exposure for the 100 and 50% increment lineages exceeded the maximum evolutionary potential exhibited by the relevant gradient evolved lineages all 100 and 50% increment lineages went extinct before the end of the experiment (Figures 2A–C). Notably, only 4 out of the 24 lineages that were exposed to 100% increments reached the clinical breakpoint before extinction. A decrease in population density (OD₆₀₀) often preceded extinction (Supplementary Figure 3). In contrast, all 50 and 25% lineages reached the clinical breakpoint. The ability to adapt appeared to be drug specific. The different adaptation potentials are reflected in the extinction of the increment lineages. For instance, 100 and 50% increment lineages died out later when adapted to AMK compared to TET (Figures 2A–C).

We were interested in observing how the different rates of environmental change affected the final genotypes and phenotypes. By design, lineages would be adapted to different antibiotic concentrations at the end of the experiment, depending on the specific ALE approach. Accordingly, we decided to compare the lineages when they had reached the clinical breakpoint (defined by EUCAST, 2016). Gradient evolved lineages were analyzed at the time point when they had reached or exceeded the clinical breakpoint for 2 consecutive days. Increment evolved lineages were analyzed at the time point when the lineage had reached the clinical breakpoint. As most of the 100% increment lineages failed to reach the clinical breakpoint we excluded these lineages from the analysis.

One colony was obtained from each lineage at the time point where the population had reached the clinical breakpoint and the antibiotic tolerance was determined. The IC₈₅ values were normalized to the average IC₈₅ of the media adapted strains (Figure 3A).

For AMK and TET adapted strains the resistance level of the gradient evolved strains was above the clinical breakpoint (Figure 3A). In contrast, only one of the strains adapted to PIP was above the clinical breakpoint. High fluctuations in resistance level were observed in the PIP adapted lineages (Figure 2B) suggesting that an inoculum effect rather than real adaptation contributed to the population tolerance. The inoculum used in this study corresponded to about 10⁸ CFU/ml and did not indicate inoculum effect in previous experiments (Supplementary Figure 2). Yet, it is exceeding the reported CFU/ml concentration causing inoculum effect for PIP (Bryson and Brogden, 2012).

The resistance levels of the 25% increment evolved strains displayed a normal distribution around the clinical breakpoint for TET and PIP and were above the clinical breakpoint for AMK, whereas the 50% increment strains displayed a slightly lower tolerance (Figure 3A). However, when comparing gradient and increment adapted strains, only TET evolved strains showed a significant (P < 0.05 Kruskal-Wallis one-way analysis of variance) difference in their resistance levels (Figure 3A).

Many resistance-conferring mutations are known to confer a fitness cost, which can often be detected by a reduced growth rate (Linkevicius et al., 2013). Since the selection regime seems to influence the fitness of the resulting lineages (Lindsey et al., 2013), we measured the growth rate of the same isolated colonies that were used for the IC₈₅ determination and an additional six isolated colonies for each lineage, resulting in 28 clones for the gradient approaches and 56 clones for the increment experiments (Figure 3B). Adaptation to AMK generally seemed to be connected with a reduced growth rate compared to the other drugs.

For all three drugs the 25% increment strains grew significantly (Kruskal-Wallis one-way analysis of variance P < 0.05) faster than the gradient adapted strains (Figure 3B). The growth advantage of the increment lineages could be due to a larger number of generations that they underwent compared to the gradient evolved lineages providing better opportunity for fitter mutants to outcompete resistant mutants with larger fitness costs and to accumulate compensatory mutations that can balance the fitness costs of resistance conferring mutations. However, the time that a population was evolved for and the doubling time are not significantly correlated (R = −0.019, Pearson's product-moment correlation coefficient, P = 0.71). Accordingly, it is likely that the shorter doubling time of the increment lineages is due to a lower selection pressure toward drug resistance, resulting in an increased selection for high growth rate. This finding is in line with previous studies reporting that E. coli lineages evolved to rifampicin under sudden drug increase have a significantly reduced growth rate compared to lineages evolved to more gradual drug increases (Lindsey et al., 2013) and that E. coli lineages adapted to 22 different antibiotics under mild selection have an elevated growth rate compared to lineages evolved under strong selection regimes (Oz et al., 2014). Yet, no correlation was found between the resistance level and the growth rate suggesting that a mutation that confers high-level resistance is not necessarily linked to a high fitness cost and vice versa (Supplementary Figure 4).

Antibiotic resistant bacteria often show cross-resistance to similar drugs (Szybalski and Bryson, 1952). Interestingly, increased susceptibility toward other antibiotics can also frequently be observed (Szybalski and Bryson, 1952; Imamovic and Sommer, 2013; Lázár et al., 2014; Munck et al., 2014), a phenomenon commonly referred to as collateral sensitivity. In order to test if the cross-resistance and collateral sensitivity phenotype is influenced by the selection regime we determined the drug resistance profiles for each evolved strain toward each of the three drugs tested (Figure 4). All lineages adapted to AMK showed collateral sensitivity toward PIP. However, the gradient and 25% increment AMK adapted strains differed in their collateral sensitivity toward TET (P < 0.05, Kruskal-Wallis one-way) (Figure 4). It should be noted that the end point of the gradient adapted AMK lineages showed collateral sensitivity toward TET in accordance with previous studies (Munck et al., 2014), suggesting that the number of generations that a lineage was allowed to undergo before testing the collateral sensitivity was important. Isolates from gradient and increment lineages all showed collateral sensitivity against AMK when evolved to TET and cross-resistance between PIP and TET (Figure 4). However, strains adapted to PIP in the gradient approach were slightly less resistant (P < 0.05, Kruskal-Wallis one-way) to TET than the strains adapted in the 25% increment approaches (Figure 4). While, Oz et al. (2014) highlights differences in collateral sensitivity in mildly and strongly selected lineages, our results are in line with findings by Lázár et al. (2014) suggesting that phenotypic similarities dominate over differences.

Overall it can be concluded that the cross-resistance was very similar between gradient and increment approaches. The main phenotypic difference between gradient and increment evolved lineages is a slower growth rate of the gradient evolved lineages.

The strains used for IC₈₅ determination and growth rate measurements were sequenced in order to uncover the underlying genetic changes. We identified a total of 173 mutations across 92 sequenced strains (Supplementary Table 5). Large insertions and deletions made up 26.5% of the total number of mutations (Supplementary Figure 5A). These larger genetic rearrangements are frequently overlooked but can play important roles in the genetic adaptation process. Two of the eight parallel strains adapted in a 25% increment to AMK have three identical SNPs in common, suggesting potential cross-contamination between the lineages. Therefore, only one of the strains was used for the following analysis.

On average we identified about two mutations in each strain across the different experiments (Supplementary Figure 5B). Even though the 25% increment lineages were evolved for more generations until they reached the clinical breakpoint there was no significant difference in the number of mutations between experimental setups (P > 0.5, Students t-test) and the number of mutations in the sequenced isolates did not correlate significantly with the number of generations (R = 0.19, Pearson's correlation, P = 0.097).

Whether a mutation confers resistance, compensates fitness costs of other mutations or hitchhikes with a resistance mutation is difficult to determine without re-introducing specific mutations alone and in combinations into the non-evolved wild type. However, if a gene is mutated in more than one independent strain it is likely that the mutation was selected for (Lieberman et al., 2011; Yang et al., 2011; Sandberg et al., 2014). We filtered our dataset according to this criterion and found that 88.8% of genes mutated in the gradient evolved strains were also mutated in the increment strains (Figure 5). Except mutations in the ATP synthase gamma chain (atpG) and the cytochrome bo(3) ubiquinol oxidase subunit 2 (cyoA), all mutated genes of the gradient strains in the filtered dataset have also been found to be mutated in the 25% increment strains (Figure 5). The clones carrying one of the two mutations have an average doubling time that is four times higher than the wild type and twice as high as the average of all strains adapted to AMK. Therefore, it is likely that these mutations come with a high fitness cost, and accordingly were not fixed in the 25% increment lineages.

Interestingly, the increment-adapted strains carried not only most of the mutations found in the gradient adapted strains, but also many mutations that were solely identified in the increment strains (Figure 5). Such mutated genes are for example nuoB, nuoC, and nuoH, subunits of the NADH-quinone oxidreductases that shuttle electrons from NADH to quinones in the respiratory chain that have been identified to confer resistance toward AMK in previous studies (Kohanski et al., 2007; Schurek et al., 2008; Girgis et al., 2009; Wong et al., 2014). In addition, these mutations were linked to the collateral-sensitivity phenotypes of aminoglycosides toward many other classes of antibiotics (Lázár et al., 2014), suggesting that mutations in these genes are relevant for the collateral sensitivity phenotype toward TET, that was not observed in the gradient lineages when they reached the clinical breakpoint, but only in strains isolated from the end point of the gradient evolved lineages.

Mutations in the genes fusA, sbmA as well as in two different two component systems, cpxRA and arcAB appeared to be the dominating mutations in all strains adapted to AMK (Figure 5). All mutations have been previously linked to AMK resistance (Laviña et al., 1986; Busse et al., 1992; Johanson and Hughes, 1994; Salomón and Farías, 1995; Macvanin and Hughes, 2005; Kohanski et al., 2008, 2010; Pena-Miller et al., 2013; Lázár et al., 2014; Munck et al., 2014). Mutations in the elongation factor G encoding gene fusA have been shown to result in collateral sensitivity toward beta-lactam antibiotics, as observed in this study for PIP (Macvanin and Hughes, 2005).

acrR was found to be the predominantly mutated gene in the PIP evolved lineages regardless of the experimental setup (Figure 5). Mutations in acrR as well as in marR and rob, also identified to be mutated in the increment strains adapted to PIP, lead to the multiple antibiotic resistance (mar) phenotype, which was described to confer resistance toward a variety of drugs including beta-lactam antibiotics and tetracyclines. This finding explains the cross-resistance observed for TET and PIP evolved lineages in this experiment (George and Levy, 1983; Cohen et al., 1989; Ariza et al., 1995; Maneewannakul and Levy, 1996; Oethinger et al., 1998). The lack of mutations in marR and rob in the gradient adapted strains might account for the difference in cross-resistance toward TET compared to the 25% increment adapted strains. However, it can be speculated that these mutations would also occur in a gradient system if the inoculum effect can be avoided, since they were observed previously in an experiment following the gradient approach (Munck et al., 2014). Another frequently observed mutation in PIP adapted strains affects the drug target, the peptidoglycan synthase ftsI (penicillin-binding protein 3) (Figure 5) (Matic et al., 2003; Blázquez et al., 2006). Interestingly, mutations in cpxA were solely found in 25 and 50% increment strains adapted to PIP (Figure 5). Mutations in this gene can confer up to 2-fold increases in resistance to beta-lactam antibiotics (Srinivasan et al., 2012; Bernal, 2014). Since a 2-fold increase in drug resistance is moderately low, it can explain why the mutation was only found in 25 and 50% increment lineages that were exposed to low antibiotic concentrations and why it was not identified in the gradient or 100% increment lineages.

The lineages adapted to TET showed, similar to PIP, mutations in genes belonging to the mar phenotype (Figure 5). In addition to mutations belonging to the mar phenotype, two other mutated genes, dksA and waaP, were identified which were previously only indirectly linked to antibiotic susceptibility (Yethon et al., 1998; Yethon and Whitfield, 2001; Hansen et al., 2008; Tamae et al., 2008; Liu et al., 2010). Interestingly, the only gene duplications observed in this experiment were all found in three different lineages in the 50% increment strains adapted to TET. Two genes, yicS and yibT, with uncharacterized gene products, were duplicated as well as phoU, whose deletion mutant was more susceptible toward antibiotics suggesting a potential role in antibiotic tolerance (Li and Zhang, 2007).

The significant difference in resistance between gradient and increment evolved lineages might be explained by the abundance of mutations in marR. Almost all strains adapted under the gradient approach carry a mutation in marR, whereas this genetic change was only observed in a few lineages evolved in the increment regime (Figure 5). Lineages harboring mutated marR were overall 15% more resistant to TET than all TET adapted lineages on average. However, they also had an increased doubling time by ~27% compared to the average. ΔmarR mutants were previously linked to an impaired fitness (Marcusson et al., 2009), suggesting that mutations in marR are more likely to dominate a population under strong selection.

The clones adapted under the gradient approach seem to have fewer mutations in the filtered data set in comparison to the increment lineages. However, they often carry mutations that were only detected once in the whole experiment (Supplementary Table 5), therefore the gradient adapted lineages display a higher diversity in unique mutations. In order to quantify the similarity and dissimilarity between genotypes in the different experimental setups, we did a pairwise comparison between all strains using the unfiltered data set. The overlap of mutated genes between the pairs was calculated in percent of the total number of mutated genes found in the two strains. We chose to analyze the similarity on the gene level and not on SNP or gene family level, since it was suggested as appropriate measure to detect parallel or convergent evolution (Achaz et al., 2014). The genetic similarity within the gradient evolved replicates was on average around 30–50% (Figure 6). The strains adapted under the 25 and 50% increments were about 45 and 30% similar to each other when adapted to AMK and PIP and only around 20 and 10% alike when evolved to TET (Figure 6). Interestingly, the genetic similarity of strains evolved using different selection regimens was comparable to the similarity within replicates from the same selection regimen (Figure 6). The similarity of the gradient and the 25% increment strains was maximal about 3% below and 12% above the group internal similarity of either the gradient or the increment 25 strains (Figure 6). This result underlines that the genetic similarity between the different selection regimens is similar to the genetic similarity observed between parallel lineages that were evolved under identical conditions.

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.