The Spread of Antibiotic Resistance Is Driven by Plasmids Among the Fastest Evolving and of Broadest Host Range

Abstract

Microorganisms endure novel challenges for which other microorganisms in other biomes may have already evolved solutions. This is the case of nosocomial bacteria under antibiotic therapy because antibiotics are of ancient natural origin and resistances to them have previously emerged in environmental bacteria. In such cases, the rate of adaptation crucially depends on the acquisition of genes by horizontal transfer of plasmids from distantly related bacteria in different biomes. We hypothesized that such processes should be driven by plasmids among the most mobile and evolvable. We confirmed these predictions by showing that plasmid species encoding antibiotic resistance are very mobile, have broad host ranges, while showing higher rates of homologous recombination and faster turnover of gene repertoires than the other plasmids. These characteristics remain outstanding when we remove resistance plasmids from our dataset, suggesting that antibiotic resistance genes are preferentially acquired and carried by plasmid species that are intrinsically very mobile and plastic. Evolvability and mobility facilitate the transfer of antibiotic resistance, and presumably of other phenotypes, across distant taxonomic groups and biomes. Hence, plasmid species, and possibly those of other mobile genetic elements, have differentiated and predictable roles in the spread of novel traits.

Introduction

Bacteria adapt quickly to environmental challenges, including stress, antibiotic treatment, or prophage predation. This is facilitated by their ability to acquire genes by horizontal transfer using mobile genetic elements (MGEs) (Brockhurst et al. 2019; Arnold et al. 2022). These elements can be mobilized by conjugation or phage transduction and remain in cells as plasmids (extrachromosomal elements) or integrated in the chromosome (Cabezon et al. 2015). Many potentially adaptive accessory traits are carried by conjugative elements or other MGEs mobilized by them or by phages and their satellites (Rankin et al. 2011). One might expect that such adaptive genes would tend to get fixed in the novel genome. Yet, functions under balancing selection, whose fitness value varies with time, space, or frequency in the population (Hedrick 2007), may be maintained in MGEs because they are constantly gained and lost. Functions may also remain in plasmids because gene expression and evolutionary rates may be modulated by plasmid copy number (San Millan et al. 2016; Ilhan et al. 2019). In all these cases, there is either an intermittent or a constant pressure to maintain a trait in a plasmid that circulates in a community. But certain functions may be predominantly found in MGEs because selection for the trait is recent and genes encoding such a trait were not initially available in the community. In such cases, genes are predominantly in MGEs because they are just arriving in the local community and MGEs are their vehicles.

An example of the importance of MGE-mediated horizontal gene transfer results from the challenge imposed by the introduction of antibiotic therapy to control human-related bacteria in the last few decades (Nicolaou and Rigol 2018). Under such an unprecedented challenge, adaptation required either functional innovation, which is often a slow process, or the acquisition of antibiotic resistance genes from bacteria inhabiting natural environments where antibiotics have been present for a long time (Larsson and Flach 2022). Indeed, horizontal gene transfer, especially of plasmids, is driving the rapid spread of antimicrobial resistance (AMR) genes in the last decades (Castaneda-Barba et al. 2024). Importantly, AMR is driven by plasmids that lacked these genes 70 years ago (Hughes and Datta 1983). Conjugative plasmids were shown to disseminate resistance to penicillins, aminoglycosides, sulfonamides and last-resort antibiotics such as colistin or carbapenems (Bean et al. 2009; Rozwandowicz et al. 2018; Wang et al. 2018; León-Sampedro et al. 2021). They are also frequently associated with the spread of pandemic multidrug-resistant high-risk clones such as Escherichia coli sequence type 131 and Klebsiella pneumoniae sequence type 258 (Deleo et al. 2014; Kondratyeva et al. 2020). Integrative conjugative elements and phage-plasmids also contribute to the spread of AMR genes (Botelho and Schulenburg 2021; Pfeifer et al. 2022). This raises the question of which MGE characteristics favor the spread of novel traits.

We hypothesize that two characteristics are crucial for gene transfer between distantly related bacteria living in different biomes: the ability of MGEs to acquire, exchange, and deliver genes (genetic plasticity, here considered as evolvability by large genetic changes) and their ability to transfer to distant hosts (mobility and especially host range). The first may result from allelic exchanges between homologous genes in MGEs (homologous recombination). Alternatively, genetic exchanges may be mediated by recombinases, including transposable elements and integrons (Amos et al. 2018; Wang et al. 2024). A recent study reported that 77% of the plasmids carrying AMR genes have insertion sequences (Che et al. 2021), and many studies described intragenomic transfer of resistance genes in composite transposons and integrons (Deng et al. 2015; Wang and Dagan 2024). For instance, the colistin resistance gene mcr-1 is present in many different plasmid types largely because of its association with the ISApl1 transposon (Wang et al. 2018). Recombinases facilitate the translocation of AMR genes in the environmental bacteria encoding antibiotic resistance, their subsequent transfer between plasmids, and their potential future integration in chromosomes (Wang et al. 2022). Of note, plasmids have by far the highest densities of transposable elements among bacterial replicons. Mobility is also essential for genetic transfer and host range is thus expected to contribute to explain the preferential presence of AMR genes in some plasmids. Broad host-range plasmids are more likely to achieve transfer from distantly related environmental bacteria. These variables may combine to produce particularly effective vectors of AMR. For example, the broad host-range plasmid pNDM-1_Dok01 partly responsible for the spread of β-lactamase resistance in E. coli ST38 carries the New Delhi metallo-β-lactamase (NDM-1) in a composite transposon very similar to the ones of distantly related environmental bacteria (Sekizuka et al. 2011).

Comparisons between antibiotic resistance plasmids often reveal a quick pace of change for these elements (Lipworth et al. 2024), whereas other plasmids are remarkably stable (Anda et al. 2023). This led us to hypothesize that plasmids carrying antibiotic resistance in current nosocomial bacteria are among those that are most plastic (in terms of rates of allelic exchanges and pace of change of gene repertoires), mobile (able to transfer or be transferred horizontally), and with broader host ranges. These characteristics would contribute to explain why transfer of a novel resistance occurs rapidly after the introduction of the antibiotic on the market (Hutchings et al. 2019). Our hypothesis might have broader implications: When an unprecedented challenge leads to strong selection for novel traits that are absent in a community, their acquisition from other communities will be driven by MGEs particularly plastic and with broad host ranges. To test our hypothesis, we focus specifically on antibiotic resistance because of its clinical relevance and data availability. We focus on plasmids because they are the key spreaders of AMR and because they are easy to identify and delimit in completely assembled genomes. More precisely, we identified the plasmids carrying and lacking antibiotic resistance genes, clustered them in plasmid taxonomic units (PTUs), and analyzed them in terms of mobility, host range, homologous recombination, and gene turnover. Finally, we integrated these results to rank the relevance of the different variables and to understand if these plasmids are more plastic because they carry resistance genes or because they are part of particularly plastic plasmid species.

Materials and Methods

Plasmid Dataset

The dataset used in this study consists of 32,798 complete genomes from 8,345 species retrieved from NCBI RefSeq database of high-quality complete nonredundant prokaryotic genomes (https://ftp.ncbi.nlm.nih.gov/genomes/refseq/, last accessed in May 2023 (O'Leary et al. 2016)). These genomes contain 38,057 plasmids. To avoid the misidentification of secondary chromosomes as plasmids, 996 plasmids larger than 500 kb were excluded (Harrison et al. 2010; Smillie et al. 2010). The complete list of 37,061 plasmids can be found in supplementary table S1, Supplementary Material online.

Coding Sequence Prediction

We have previously observed that some automatic annotations mis-annotate some genes in MGEs, notably the relaxase, which is essential to study plasmid mobility. To ensure that sequence annotation is homogeneous and appropriate for our purpose, all plasmid sequences were re-annotated using prodigal 2.6.3 (Hyatt et al. 2010) with the “-p meta” option to obtain the putative protein coding genes.

Plasmid Mobility Prediction

Conjugative systems were detected using MacSyFinder v2.07rc (Néron et al. 2023) with the “Plasmids” models of the CONJScan module. Briefly, MacSyFinder uses HMM protein profiles and a set of rules (which are defined in the models) about their presence and proximity. MacSyFinder was used with default parameters, and CONJScan models were ran all at once with the parameters “- all -ordered -circular”. The plasmids were then classed in relation to the conjugation genes they contain (as in Coluzzi et al. 2022). The plasmids encoding a complete mating pair formation system and a relaxase were classed as conjugative (pCONJ). The remaining plasmids that encode a relaxase were classed as mobilizable plasmids (pMOB). Plasmids that are none of the above but have an oriT and may thus be mobilized by a relaxase and a MPF encoded in trans were classed as pOriT. These origins of transfer were identified using blastn (v2.09.0+) (Camacho et al. 2009) with parameters -task blastn-short -evalue 0.01 filtered at identity > 80% and coverage > 80% and a dataset of 91 sequences of different origins of transfer (Ares-Arroyo et al. 2023).

Plasmids were classed as phage-plasmids (P-P) as in Pfeifer et al. (2021). The remaining plasmids were considered as putatively non-transmissible by conjugation (pNT).

Detection of Antibiotic Resistances, Integrons, and Insertion Sequences

Plasmid protein sequences were screened for the genes conferring antibiotic resistance to the 31 AMR drug classes detected by AMRFinderPlus v3.11.14 with default parameters (identity > 90%, e-value < 1e-20, coverage > 50%) (Feldgarden et al. 2021) and the database version 2023-04-17.1. Plasmids encoding at least one antibiotic resistance gene were considered AMR+ plasmids. Through the text, AMR classes and AMR genes correspond to “AMR classes” and “AMR gene symbol” as defined by AMRFinderPlus.

Insertion sequences were detected using ISEScan v1.7.2.3 (Xie and Tang 2017) with default parameters and the “--removeShortIS” option.

Integrons were detected using IntegronFinder v2 (Neron et al. 2022) with the “–local-max” and “—circ” options.

PTU and Host-Range Predictions

The PTUs and host range were assigned to each plasmid with COPLA v1.0 and default parameters (Redondo-Salvo et al. 2021). The host range indicates the ability of a PTU to colonize different bacterial species. Each PTU has a host-range grade that reflects the diversity of bacterial species in which the plasmids within that PTU have been found. Host-range grades vary between I and VI, the higher the grade, the broader the host range, where grades correspond to PTUs:

• Grade I: Plasmids of the PTU are found in strains of the same species.

• Grade II: Plasmids of the PTU are found in species of the same genus.

• Grade III: Plasmids of the PTU are found in genera that belong to the same family.

• Grade IV: Plasmids of the PTU are found in different taxonomic families of the same order.

• Grade V: Plasmids of the PTU are found in different orders of the same class.

• Grade VI: Plasmids of the PTU are found in different classes of the same phylum or even beyond, in different phyla.

Grades V and VI were merged to have enough counts for statistical tests.

PTU Pangenome and Core-Genome Construction and Phylogeny Inference

The pangenome and core genome of each PTU were computed using PanACoTA v.1.3.1 (Perrin and Rocha 2021). Briefly, the pangenome was constructed by clustering all protein sequences in the set of genomes using MMseqs2 (protein identity > 80% (Steinegger and Söding 2017)). We retrieved the core genome from each pangenome. It consists of the gene families present in exactly one copy in at least 80% of all plasmids.

The multiple sequence alignments of the families of the core genes were computed using the “align” PanACoTA v1.3.1 module. Briefly, the protein sequences of the core genes were aligned with MAFFT v7.467 (-auto parameters) (Katoh and Standley 2014) then back-translated to nucleotide alignments (i.e. each amino acid was replaced by the original codon) and concatenated. The phylogenetic inference was done from the resulting multiple alignments using IQ-TREE 2.0.6 (Nguyen et al. 2015) with the ultra-fast bootstrap option (-bb 1000 bootstraps) and with the best fitting model estimated using ModelFinder Plus (-MFP) (Kalyaanamoorthy et al. 2017). The best model for each tree was determined according to the Bayesian information criterion. To obtain robust phylogenetic reconstructions, we followed (Filipski et al. 2014) and removed from the phylogenetic analyses the plasmids encoding less than 50% of the core genes of a PTU.

The average tip-to-root distance for each phylogenetic tree was computed using the “get_distance()” function of ETEToolkit v3.1.2 (Huerta-Cepas et al. 2016).

Rarefaction Curve and Heaps’ Law Fit

Rarefaction curves for the pangenome of each PTU were computed using the package vegan for R, v2.5.6 (https://CRAN.R-project.org/package=vegan) on the presence/absence matrixes of gene families obtained by PanACoTA. As rarefaction curves require a large sample size to be accurate (Newton 2007), we computed them only for PTUs with more than 30 plasmids. The Heaps’ law was fitted to each rarefaction curves (Tettelin et al. 2008).

where F is the number of pangenome families and κ and β are empirical parameters that were inferred using the “optimize.curve_fit” method from scipy v.11.4 for Python (Virtanen et al. 2020).

Assessment of Congruence Between Phylogenetic Trees

The Robinson–Foulds symmetric distance between two trees corresponds to the number of splits or fusions that are required to transform the topology of one tree into another. The normalized Robinson–Foulds symmetric difference metric (nRF) corresponds to the Robinson–Foulds symmetric distance (Robinson and Foulds 1981) divided by the value of the theoretical maximum of the Robinson–Foulds distance. The nRF varies between 0 and 1, where pairs of trees with nRF = 0 have identical topologies and those with nRF = 1 have no common splits. For each PTU with a phylogenetic tree based on more than one core gene multiple alignment, we computed the mean nRF between the tree constructed from all concatenated core genes of the PTU and the tree constructed from each individual gene within the core genome (concatenate-based comparisons). This allows to measure the degree of incongruence between each individual core gene and the set of all core genes. High values of nRF are thus an indication of recombination. Of note, such tests are usually made with trees that do not include the genes being tested, e.g. robust phylogenies produced by external data. Yet, such trees are unavailable for plasmids. Given the recombination rates we observe, they are probably impossible to identify in many PTUs thus also precluding the use of recombination detection methods based on the bacterial clonal frame (Didelot and Falush 2007). The use of the focal gene in the comparison (as part of the genes that were used to build the core genes phylogeny) renders the method conservative, i.e. it will tend to produce lower rates of incongruence. To assess this potential problem, we also computed the nRF between the tree of every gene of the core genome of the PTU (all-pairs comparisons, self-comparisons are omitted, i.e. tree of the gene A is not compared to itself). In this case, there is no assumption of a core gene single phylogeny. These scores were computed using ete-compare method from the ETEToolkit v3.1.2 (Huerta-Cepas et al. 2016). To ensure only the robust parts of the trees were considered, all branches with a bootstrap support inferior to 75 were discarded, with the “min_support_src” and “min_support_ref” options.

wGRR Calculation

We estimated the weighted gene repertoire relatedness (wGRR) values between each pair of plasmids in our dataset using the wGRR software (https://gitlab.pasteur.fr/jgugliel/wgrr) with default parameters as described in Cury et al. (2018). Concretely, we identified plasmid proteins that are significantly similar (identity ≥ 35%, coverage ≥ 50%) among all pairs of plasmid sequences with MMseqs2 search v. 9-d36de with the parameter “sensitivity 7.5 (-s 7.5)” and “max-seqs 100000” (to retrieve all possible hits). We then identified the best bidirectional hits (BBH) between pairs of plasmids which correspond to couples of proteins that are reciprocal best matches. Finally, we use the BBH to calculate the wGRR.

where A_i and B_i are the ith BBH pair of P total pairs, id(A_i, B_i) is the identity between the pair, and min(#A, #B) is the number of genes of the plasmid having the fewest. The wGRR varies between 0 and 1. A wGRR close to 1 means that both plasmids are very similar (all genes in a plasmid have very similar best bidirectional hits in the other plasmid) and a wGRR of 0 means that plasmids have no homologs.

For each PTU, we computed the median and mean wGRR values between every possible pair of plasmids within the PTU. A high mean wGRR value means that most of the plasmids in the PTU are closely related and a small mean wGRR value means that most plasmids in the PTU are more distantly related.

Statistics

Analysis of variance, χ² tests, and ordinal and nominal logistic regressions were performed using the “Fit model” module of JMP v17. All others statistical tests were performed within R (v. 4.2.3).

Results

Plasmid Groups Are Heterogeneously Enriched in ARGs

We retrieved all plasmids sequenced within complete genomes of the nonredundant database RefSeq of the NCBI to study their evolution. We excluded the few very large elements that could be secondary chromosomes or chromids (>500 kb) (supplementary fig. S1, Supplementary Material online, supplementary table S1, Supplementary Material online), re-annotated the remaining 37,061 plasmids, and screened them for antibiotic resistance genes. We used conservative methods to focus on genes with high similarity to known antibiotic resistance genes in nosocomial pathogens (identified with AMRFinderPlus). About a quarter of them (9,384) encode at least one AMR gene (Fig. 1a and b). These plasmids are 60% larger than the others (Fig. 1c). We analyzed in detail the 392 unique antibiotic resistance genes (determined as “gene symbol” by AMRFinderPlus) which were present in a total of 46,472 times in these plasmids (supplementary table S1, Supplementary Material online). These genes conferred resistance to 31 drug classes (determined as “AMR class” by AMRFinderPlus). The most common antibiotic resistance classes found in the dataset target beta-lactams (6,032 plasmids), aminoglycosides (5,250 plasmids) and sulfonamides (3,725 plasmids) (Fig. 1d). These genes often co-occur with those of other AMR classes (Fig. 1d). We found 895 different combinations of AMR classes in plasmids, some of which are more often found in combination than alone (Fig. 1d). Accordingly, plasmids with antibiotic resistance genes have a median of three different AMR classes.

To assess if different combinations of antibiotic resistance classes could often be found in closely related plasmids, we identified the PTUs of the 37,061 plasmids. PTUs are cohesive groups of plasmids that share high average nucleotide identity and approximately similar size. Almost half (17,314) of the plasmids and two thirds (6,321) of those containing resistance genes could be assigned to one of the previously described 355 PTUs (supplementary table S1, Supplementary Material online (Redondo-Salvo et al. 2021)). Two hundred twenty-three (63%) PTUs lacked plasmids with antibiotic resistance genes. Among the others, the plasmids encoding resistance genes had widely variable frequencies, some of them only containing such plasmids (Fig. 2, supplementary fig. S2, Supplementary Material online, supplementary table S2, Supplementary Material online). Most of these PTUs contained several combinations of AMR classes. The PTU-HI2 (including pSTM6 and relatives) is an extreme case: 98% of its plasmids have resistance genes and these plasmids harbor 197 different combinations of AMR classes (supplementary fig. S3, Supplementary Material online). The high diversity observed among closely related plasmids (within the same PTU) suggests that resistance genes were acquired in multiple events.

Resistance Plasmids and Their PTUs Are Mobile and Have Broad Host Ranges

To test the hypothesis that plasmids encoding antibiotic resistance tend to have broad host ranges, we compared the PTUs containing such plasmids with the others. For these analyses, we used a previously published classification of host ranges of PTUs (Redondo-Salvo et al. 2020). In that study, host ranges were classed in six groups from I (narrowest) and VI (broadest). To have enough sample size for statistical analyses, we put together the classes V and VI. We found that the frequency of PTUs encoding AMR was higher for the broadest host ranges (nominal logistic regression P < 0.0001, Fig. 3a, supplementary table S3, Supplementary Material online). As we feared that the high number of replicons in certain PTUs could inflate the PTU host range, we made a complementary analysis where we added the number of plasmids in each PTU in the statistical model. The association between PTU containing resistance plasmids and broad host range remained significant (nominal logistic regression P < 0.0001, supplementary table S3, Supplementary Material online), and is thus not an artifact caused by the different number of plasmids per PTU.

To study the mobility of plasmids in light of the carriage of AMR, we classified the former as: phage-plasmids (P-P), conjugative plasmids (pCONJ), mobilizable plasmids (encoding a relaxase) (pMOB), plasmids with an origin of transfer for conjugation (pOriT), and putatively non-transmissible plasmids (pNT, see Materials and Methods). We then identified the most frequent mobility type for each PTU and assigned it to the group. In 97% of the cases, the major category accounted for more than 50% of the plasmids of the group and often much more (supplementary fig. S4, Supplementary Material online). For instance, in 232 out of the 355 PTU (65%), more than 90% of plasmids had the same classification in terms of mobility. Our analysis shows that PTUs with a majority of conjugative plasmids are those with the highest frequency of plasmids carrying antibiotic resistance genes, followed by diverse types of mobilizable plasmids and phage-plasmids. The plasmids putatively non-transmissible had the lowest frequency of antibiotic resistance (Fig. 3b, supplementary table S4, Supplementary Material online, P < 0.0001). Taken together, these results suggest that plasmids encoding antibiotic resistance tend to be mobile and have broader host ranges than the others.

Antibiotic Resistance Plasmids Recombine More

Plasmids may evolve rapidly by exchanging genes and alleles with other plasmids or chromosomes. Homologous recombination is usually regarded as a mechanism of allelic variation, which can speed up adaptation. But recombination involving multiple genes may also result in gene acquisition or loss. The existence of numerous plasmids in some PTUs allows to identify homologous recombination in their core genes. We found 270 PTUs having between 1 and 381 gene families present in more than 80% of the plasmids (defined as core genes of the PTU, supplementary fig. S5, Supplementary Material online). As expected, the number of core genes in a PTU correlates with the mean size of the plasmids in the PTU (Spearman correlation: 0.7033, P < 0.0001). Importantly, there is no significant difference in the number of core genes between PTUs with or without AMR genes (Wilcoxon test, P = 0.482). While this may seem surprising given that resistance plasmids are larger (see above), many of the smallest plasmids could not be placed in PTUs and the number of core genes is shaped by the size of the element and its rate of variation of gene repertoires (see below). Hence, whatever differences one may find between these two groups of PTUs, they are not caused by differences in the number of their core genes.

We analyzed the rates of recombination using two variants of a phylogenetic method. In method I, we compared the evolutionary history of each core gene in the PTU with that of the concatenate of the core genes of the PTU (henceforth called the consensus phylogeny). This was done using the normalized Robinson–Foulds (nRF) distances which measure the similarity between two phylogenetic trees (Fig. 4). This analysis involved the analysis of 11,294 gene trees from the 250 PTUs with more than one core gene. Seventeen PTUs had congruent trees, i.e. the evolutionary histories of individual core genes were similar from the consensus, suggesting little or no recombination (mean nRF = 0, supplementary table S2, Supplementary Material online). Most of the others showed evidence of extensive recombination. The average incongruence was high (median nRF of 0.505), i.e. the evolutionary histories of the individual core genes were quite different from the consensus. If there are high recombination rates in plasmids, then the phylogeny of the concatenate of the core genes may not be reliable. Hence, we developed a complementary approach (method II, Fig. 4) where we compared the average nRF between every pair of gene families in the PTU. This analysis does not require the use of a consensus trees and results in similar conclusions (median nRF of 0.59, supplementary table S5, Supplementary Material online), confirming that recombination is frequent in plasmid core genes.

To test the idea that resistance plasmids recombine at higher rates than the others, we compared the average incongruence (obtained using method I) between the PTUs with or without AMR genes. In agreement with our hypothesis, incongruences were higher in PTUs containing resistance plasmids (Fig. 4, nominal logistic regression, P < 0.0001, supplementary table S6, Supplementary Material online). A complementary statistical model including the number of core genes and the size of the plasmids as covariates confirmed the significance of the positive association between recombination (measured by nRF) and the presence of antibiotic resistance in PTUs (P < 0.0001, same test, supplementary table S6, Supplementary Material online). This result is not caused by different depths in the phylogenetic trees of the two groups of PTUs, since they are not significantly different in this respect (supplementary fig. S6, Supplementary Material online). A similar type of association between nRF and AMR was found with the method II (P < 0.0001, same test, supplementary table S5, Supplementary Material online). Overall, these results show that PTUs with antibiotic resistance plasmids engage in allelic exchanges at higher rates than the others.

Gene Turnover Is Higher in Plasmid Species Encoding AMR

If plasmids with antibiotic resistance genes are more plastic, then their PTUs might have a higher diversity of gene families, i.e. larger pangenomes. To test this hypothesis, we computed the rarefaction curves of the pangenomes for PTUs with more than 30 plasmids and fitted them with the Heaps’ law (see Materials and Methods). We then analyzed the β parameter that varies between 0 (close, invariant genome) and 1 (open, highly variable genome) (Tettelin et al. 2008) (Fig. 5a). The pangenomes of PTUs are open (median β is 0.526), and those with antibiotic resistance genes are more open than the others (β higher by 17.5%, Fig. 5a, P = 0.0007, nominal logistic regression, supplementary table S7, Supplementary Material online). As the number of plasmids in the PTU and the size of the plasmids could have an impact on the saturation curves, we made an additional model including both parameters. This complementary analysis confirmed a significantly higher average value of β in the PTUs containing AMR genes (P = 0.004, same test, supplementary table S7, Supplementary Material online). These results indicate that plasmids with antibiotic resistance genes have a more open pangenome.

If the pangenomes of PTUs with antibiotic resistance plasmids are more open, then these plasmids should be more diverse in terms of accessory genes. To test this hypothesis, we compared the divergence in gene repertoires between all pairs of plasmids of the same PTU with the weighted gene repertoire relatedness index (wGRR, see Materials and Methods). The wGRR varies between 0 and 1. A wGRR close to 1 means that both plasmids are very similar (all genes in a plasmid have a very similar BBH in the other) and a wGRR of 0 means that plasmids lack homologs. As expected, PTUs with antibiotic resistance genes had a significantly smaller wGRR than the others (Fig. 5b). When we fitted a nominal logistic model between the AMR status of a PTU and the mean wGRR, smaller values of wGRR were associated with PTUs with antibiotic resistance genes (nominal logistic regression, P = 0.0006), even controlling for the number of core genes and the number of plasmids in the PTU (supplementary table S8, Supplementary Material online). Hence, despite the high similarity of core genes within PTUs, their plasmids have highly variable gene repertoires. PTUs encoding AMR are even more variable than the others.

Acquisition and exchange of AMR genes in plasmids is often driven by integrons and/or transposons (see Introduction). Indeed, we found that PTUs encoding AMR harbor higher densities of transposable elements and are the only ones with integrons (Fig. 5c and supplementary fig. S7, Supplementary Material online). These results suggest that, while having similar numbers of core genes, the PTUs with antibiotic resistance genes are enriched in genetic elements that increase the rates of gain and loss (turnover) of accessory genes.

PTUs with antibiotic resistance genes have higher recombination rates in core genes, higher gene repertoire diversity, and higher gene turnover. Since these factors are not independent, we assessed their joint effect. All three variables are significantly and positively associated with PTUs with antibiotic resistance genes (P < 0.0001, nominal logistic regression, supplementary table S9, Supplementary Material online), suggesting that all these factors favor the acquisition of antibiotic resistance genes. We then compared the Wald χ^² and the odds ratios of these factors to determine which had the strongest effect. The openness of the pangenome had the strongest effect (beta_χ² = 11.58, odds ratio: 7,995.991), followed by the gene repertoire diversity (wGRR_χ² = 8.79, odds ratio: 0.003466), and the core-genome recombination rate (nRF_χ² = 5.38, odds ratio: 70.631). This suggests that while all these factors contribute to the acquisition of AMR genes by plasmids, the capacity to quickly acquire new genes (resulting in a more open pangenome and higher wGRR) is particularly important.

Evolvability and Mobility Are Intrinsic to the AMR PTUs

The previous results showed that PTUs with antibiotic resistance genes are more plastic than the others. This can be interpreted in two ways. (i) Plasmids that are intrinsically more plastic and mobile are more likely to acquire resistance genes and spread them. (ii) Resistance itself favors rapid genetic diversification and mobility. The second hypothesis is consistent with the association of resistance genes with recombinases (IS, integrons), although the way resistance could favor plasmid mobility is less obvious. To test the first hypothesis, we retained the classification of PTUs into AMR+ and AMR− based on whether they originally contained resistance plasmids. However, for this analysis, we removed the resistance plasmids from the AMR+ PTUs. As a result, we still have two groups—AMR+ and AMR−—but the AMR+ PTUs no longer contain any resistance plasmids (Fig. 6). AMR+ PTUs still contained more conjugative elements and had broader host ranges than the others (supplementary table S10, Supplementary Material online, P < 0.0001). They also had higher recombination rates, higher mean nRF values (P = 0.0170, nominal logistic regression, supplementary table S10, Supplementary Material online, Fig. 6a), lower wGRR, and more open pangenomes (P_beta = 0.0027, P_wGRR = 0.0006, same test, supplementary table S10, Supplementary Material online, Fig. 6b and c). The density of transposable elements was higher in these PTUs (P = 0.01025, Fig. 6d). Hence, the rapid evolution of antibiotic resistance plasmids seems intrinsic to the plasmid families and is not just a consequence of the presence of resistance genes in some plasmids of the PTU. This suggests a causality: More plastic and mobile plasmid families are more likely to acquire antibiotic resistance genes and spread them among human-related bacteria.

Discussion

Here, we investigated how bacteria acquire new traits via horizontal gene transfer, by focusing on the characteristics of MGEs that enable them to tackle unprecedented challenges. Antibiotic resistance in nosocomials differs from traits like virulence, competition, or phage immunity by its extreme novelty. This is why we focused on known antibiotic resistance genes of relevance in the clinic, as they are recent in these bacteria. Environmental bacteria with other antibiotic resistance genes may endure different dynamics because we cannot ensure such trait is recent in their habitat. It is not surprising that conjugative plasmids have the highest frequency of such AMR genes in our dataset, since they are mobile, tend to have broad host range, can transfer large amounts of genetic information per event, and tolerate large integrations of genetic material (Forster et al. 2022). Plasmids that are not autonomously conjugative but can be mobilized by conjugation also account for a large fraction of transfer of resistance genes, and some of them have among the broadest known host ranges (Meyer 2009). Future work will be needed to include the substantial number of plasmids, mostly small, that could not be classed in PTUs. Overall, our results are in line with the hypothesis that broad host range conjugation has a key role in driving the initial stages of transfer of novel traits between distantly related bacteria.

By focusing on homogeneous groups of closely related plasmids (PTUs), we could analyze these plasmids in a phylogenetic context to quantify their rates of change in gene repertoires and of sequence recombination. Unfortunately, linking all the PTUs by a phylogeny is impossible due to the lack of homologous genes. Hence, we treated the PTUs as effectively independent units. These analyses confirmed that families of plasmids encoding antibiotic resistance have larger pangenomes and faster-changing gene repertoires than the others. We cannot exclude the possibility that the bias in genome databases toward closely related antibiotic resistant bacteria may result in larger numbers of plasmids in PTUs encoding AMR. This should make our analysis conservative, since more closely related plasmids are more likely to reveal higher average values of gene repertoire relatedness (wGRR) and we found the opposite trend. It could also be argued that antibiotic resistance plasmids inflate estimates of the PTUs rate of diversification of gene repertoires because we know these genes were acquired recently. Yet, the analysis of PTUs after removing antibiotic resistance plasmids corroborated our hypothesis that plasmid evolvability favors AMR acquisition.

The acquisition of AMR is likely facilitated by the abundance of transposable elements and integrons in plasmids, as previously described for resistance plasmids (Che et al. 2021) and here further extended to their PTUs (even when removing AMR plasmids from the analyses). The reasons for the higher frequency of transposable elements in some plasmids are poorly understood. It may be related with the availability of more neutral integration sites in the larger conjugative plasmids or with preferential transposition during plasmid conjugation (Ton-Hoang et al. 2010). It is important to note that transposable elements may also facilitate the transfer of resistance from plasmids to chromosomes and this might stabilize resistance in the genome at a later stage (Kadibalban et al. 2024). It is thus possible that the same MGEs that facilitated the acquisition of antibiotic resistance by highly mobile and plastic plasmids will subsequently have a key role in the transfer of these genes to more stable regions of the bacterial genome.

The relevance of homologous recombination in the acquisition of AMR by plasmids has rarely been studied, even though it was previously observed that homologous recombination was frequent in core genes of IncP-1 plasmids (Norberg et al. 2011). Actually, the extent of recombination, measured by phylogenetic incongruence, was very different between PTUs, and was particularly high in PTUs encoding AMR. It's important to note that homologous recombination may result in the acquisition of intervening genes when recombination tracts are long (Doroghazi and Buckley 2010). Hence, homologous recombination may be part of the mechanisms allowing the plasmids to acquire resistance genes. The abundance of transposable elements in AMR plasmids may also increase the chances of DNA breaks that stimulate repair by homologous recombination (Mahillon and Chandler 1998; Li and Heyer 2008), and thereby stimulate allelic exchanges. Plasmid copy number, when high, may also result in higher frequency of genetic exchanges (Nicoloff et al. 2024). Finally, the mobility of antibiotic resistance plasmids may increase the likelihood of their occasional co-occurrence with other plasmids, which may result in more frequent and varied genetic exchanges (Wang et al. 2024). In the future, it would be interesting to compare these results with those of plasmids carrying other types of traits like virulence or defense to understand if given time antibiotic resistance genes will tend to move toward plasmids with other characteristics such as low cost and high stability.

Phages can disperse freely in the environment and might be better at transferring novel genes between biomes. But temperate phages tend to have narrow host ranges caused by their dependency on specific receptor binding proteins (de Jonge et al. 2019). These phages are also less tolerant to gene acquisition than conjugative elements because increase in genome size rapidly results in defaults in DNA packaging in the viral particle (Hatfull and Hendrix 2011). Our analyses highlighted that variability in gene repertoire changes seems to be the most important trait for the acquisition of AMR genes. This may explain why integrative phages tend to lack such genes (Enault et al. 2017) since they are very constrained in size. In contrast, phage-plasmids play a role in the spread of antibiotic resistance possibly because they may accept novel genes more easily (they are larger), they have more transposable elements and integrons, and they exchange genes with plasmids at larger rates (Pfeifer et al. 2022; Pfeifer and Rocha 2024). Interestingly, P1-like phage-plasmids, which are particularly large and of broad host range (Murooka and Harada 1979; Pfeifer and Rocha 2024), are among the ones having most resistances (supplementary fig. S3, Supplementary Material online). This suggests that host range and evolvability also facilitate the transfer of novel traits by phages. Of note, while conjugation seems to have a key role in the initial transfer of traits across biomes and distantly related bacteria (Che et al. 2021), once the antibiotic resistance genes are acquired by nosocomial bacteria, phage-mediated transduction and natural transformation may further spread them between closely related bacteria (Winter et al. 2021).

PTUs encoding antibiotic resistance are associated with conjugation, broad host range, and several forms of genetic diversification, even when we removed antibiotic resistance plasmids from the dataset. This suggests that these are intrinsic characteristics of the plasmid species that acquire AMR, and not (only) the result of antibiotic resistance acquisition. That some plasmid species are more prone to the transfer of adaptive genes from distant bacteria in different biomes, implies that plasmids have differentiated roles in the processes of horizontal gene transfer. This may be of importance to study the transfer of traits tackling other unprecedented challenges, such as response to pollution or climate change.

Supplementary Material

Supplementary material is available at Molecular Biology and Evolution online.

Acknowledgments

We thank Manuel Ares-Arroyo and Marie Touchon for discussions and comments on the manuscript. We thank Eugen Pfeifer for providing the phage-plasmid classification. This work used the computational and storage services (TARS cluster) provided by the IT department at Institut Pasteur, Paris. Work was supported by the Laboratoire d’Excellence IBEID Integrative Biology of Emerging Infectious Diseases [ANR-10-LABX-62-IBEID], Equipe FRM (Fondation pour la Recherche Médicale) [EQU201903007835], the TRANSFOCONFLICT project [ANR20-CE12-0004-02] and the INCEPTION project [PIA/ANR-16-CONV-0005].

Data Availability

All major data, or identifiers to public repositories with the data, are incorporated into the article and its online Supplementary material. Other data, e.g. intermediary results files, underlying this article will be shared on reasonable request to the corresponding author.

References