Long-term overuse of antibiotics has driven the propagation and spreading of antibiotic resistance genes (ARGs) such as efflux pumps in the environment, which can be transferred to clinically relevant pathogens. This study explored the abundance and diversity of ARGs and mobile genetic elements within bacterial communities from sediments of three Western Balkans glacial lakes: Plav Lake (high impact of human population), Black Lake (medium impact of human population) and Donje Bare Lake (remote lake, minimal impact of human population) via shotgun metagenomics. Assembled metagenomic sequences revealed that Resistance-Nodulation-Division (RND) efflux pumps genes were most abundant in metagenome from the Plav Lake. The Integron Finder bioinformatics tool detected 38 clusters of attC sites lacking integron-integrases (CALIN) elements: 20 from Plav Lake, four from Black Lake and 14 from Donje Bare Lake. A complete integron sequence was recovered only from the assembled metagenome from Plav Lake. Plasmid contents within the metagenomes were similar, with proportions of contigs being plasmid-related: 1.73% for Plav Lake, 1.59% for Black Lake and 1.64% for Donje Bare Lake. The investigation showed that RNDs and mobile genetic elements content correlated with human population impact.

The discovery and introduction of antibiotics into clinical practice enabled the treatment of infectious diseases and facilitated advances in modern medicine such as complex surgery, transplantation and chemotherapy (Aminov 2010). Besides their role in clinical use, antibiotics have been utilised in a vast array of applications, from agriculture and growth-promotion in animals to preserving building materials from contamination (Singer et al. 2016). However, the rapid spread of antibiotic resistance in communities and clinical settings threatens the successful treatment of infectious diseases in modern society (Hawkey 2008). Since the antibiotic resistance crisis emerged (Ventola 2015), it has caused increased mortalities and substantial costs (Thorpe et al. 2018).

Genetic determinants of antibiotic resistance are commonly associated with mobile genetic elements – plasmids and transposons (Gaze et al. 2011). This association enables rapid spread among bacterial species, emphasising the genetic context of antibiotic resistance as a major factor for risk assessment when considering transfer from an environmental source to clinically relevant pathogens (Bengtsson-Palme et al. 2018). The selective pressure of antibiotics derived from anthropogenic activity results in maintaining a wide range of known antibiotic resistance genetic determinants in the environment as well as in the recruitment of novel resistance determinants from environmental bacteria (Larsson et al. 2007). Integrons are a type of genetic assembly platforms often associated with mobile genetic elements (Mazel 2006). They often harbour cassettes encoding antibiotic resistance determinants and have the ability to capture new gene cassettes by the activity of integrases, a family of tyrosine recombinases (Mazel 2006). Prevalence of integrons has been reported to be elevated in anthropogenically impacted sites, indicating a potential correlation between pollution and antibiotic resistance dissemination (Rosewarne et al. 2010).

Lakes are predicted to have the potential to store and accumulate antibiotic resistance genes (ARGs) to a greater extent than rivers (Czekalski et al. 2015), due to increased residence time of contaminants because of water retention pollutants from wastewaters slowly circulating around the lakes (Reuther 2009).

Western Balkans glacial lakes are exposed to anthropogenic activity mostly due to the vicinity of human settlements, agricultural production or tourism. Glacial lakes in Balkan Peninsula could be considered as endangered due to climate changes and their sensitivity to atmospheric conditions. However, studies of the human population impact on bacterial diversity in glacial lakes are rare (Hotaling et al. 2017). Recently, we correlated bacterial composition and diversity in three Western Balkans glacial lakes with a human population in the vicinity by using 16S rRNA gene-based metagenomic analysis (Malešević et al. 2019). Having this in mind, we recognised that bacterial communities of sediments from Western Balkan glacial lakes are unexplored and should be deeply analysed in order to fully assess the human impact by determining the presence of ARGs and mobile genetic elements associated with these metagenome sequences. Bioinformatic analysis of metagenomes was used to directly access the genetic potential of entire microbial communities of glacial lake sediments without cultivation. Metagenomic data sets were examined in order to detect the presence of antibiotic resistance determinants by sequence similarity with available databases.

Sampling sites and sediment collection

Sediment samples analysed in this study were collected in June 2016 from three glacial lakes, Plav Lake, Black Lake and Donje Bare Lake, exposed to different levels of human population impact. Plav Lake (42.595833 N 19.925000 E) is located in Plav municipality of Montenegro and is very exposed to human activity. Sampling was carried out in the vicinity of the urban area of the town of Plav. Black Lake (43.143333 N 19.087500 E) consists of two lakes, Big Lake and Small Lake, and sampling was carried out from Big Lake counterpart. Black Lake is the premium tourist attraction of the Durmitor mountain area in Montenegro and it is easily accessible as it is 3 km distant from the centre of the town of Žabljak. Donje Bare Lake (43.318333 N 18.630833 E) is a lake of Bosnia and Herzegovina and it is located 20 km from the nearest inhabited village. Sampling was described in detail in our previous work (Malešević et al. 2019). Briefly, ten samples of sediment within an area of 10 m2 were taken from each lake. Samples from the same lake were pooled together and after that metagenomic DNA was isolated and sequenced.

DNA extraction and metagenome sequencing

Total DNA was extracted from 1 g of each sediment sample by using a Qiagen PowerSoil DNA isolation kit (Carlsbad, CA, USA) in accordance with the manufacturer's instructions. The amount and purity of isolated DNA were determined by using a NanoVue Spectrophotometer (GE Healthcare, USA) based on the absorbancy of 260 nm (A260). The extracted total DNA was stored at –80 °C prior to analysis. Approximately 5 μg of a DNA aliquot was proceeded for shotgun metagenomic sequencing using an Illumina HiSeq 2000 sequencing platform (100 bp paired-end reads) at the commercial service in Macrogen (Seoul, Republic of Korea).

Metagenomic sequence bioinformatic analysis

The quality of each sequencing library was assessed using FastQC (Andrews 2010). Due to varying depths across mixed organisms, it was difficult to assemble the metagenome using a single strategy. IDBA-UD with multi k-mer mode outperformed the assembly using De Bruijn Graph methods (Peng et al. 2012). In addition, unassembled reads were collected and assembled by Celera Assembler with the Best Overlap Graph (CABOG) (Myers et al. 2000).

Metagenomes were screened for the presence of the genetic determinants of antibiotic resistance by the publicly available database – The Comprehensive Antibiotic Resistance Database (CARD, https://card.mcmaster.ca/) (McArthur et al. 2013) performing BLAST searches of metagenomes contigs against the CARD database. Results were visualised as a heat map constructed using Rx64 3.5.1. The results were represented on the relative scale ranging from 0 (blue) as the lowest values, progressing to white, then to 100 (red) as the highest values.

Identification and analyses of integrons, cassette arrays and clusters of attC sites lacking integron-integrases (CALIN) elements was performed using Integron Finder (https://github.com/gem-pasteur/Integron_Finder) (Cury et al. 2016) and Linux command line. Genes encoding antibiotic resistance were detected in CALIN elements by performing BlastX searches against the NCBI Proteins database.

For identification of genes encoding various classes of integrases (intI1, intI2 and intI3 genes), BlastX was used on the Linux command line to search against the database made of known intI sequences.

A blast search of reads against the NCBI Plasmid database that already provides a list of plasmids from the RefSeq database (www.ncbi.nlm.nih.gov/genome/plasmids/) was performed in Linux command line. Sequences with identity >90% and hit length >89 nt were further manually checked by BlastN and contigs with the following thresholds: identity >90% and alignment length with at least 90% of query sequence length, were considered as plasmid sequences.

The protein sequences containing known plasmid replicon domains were extracted from the Pfam database v32.0 (El-Gebali et al. 2019) and then formatted as a BLAST database file by formatdb. Pfam families were derived from Schlüter et al. (2008) and Ma et al. (2012). Searching for putative plasmid sequences of all three metagenomes was carried out using BlastX against our Pfam-derived database of plasmid-associated proteins (E-value <1e-11) (Altschul et al. 1990). The scaffold contigs with one or more hits to plasmid protein genes were assigned as putative plasmid sequences.

Quantitative PCR (qPCR)

Quantification of antibiotic resistance and integrase genes was performed by qPCR method. Among twelve sets of primers used, six pairs detected genes that confer resistance to the following antibiotic classes: β-lactams (blaTEM gene), tetracyclines (tetC gene), aminoglycosides (aacA and aadA genes), chloramphenicol (catA1 gene) and macrolides (ermB gene) (Aminov et al. 2002; Mroczkowska & Barlow 2008; Muziasari et al. 2017). Also, the multiple antibiotic efflux systems (mexA, mexB and mexD genes) were quantified, as well as genes encoding different classes of integrases (intI1, intI2 and intI3 genes) (Muziasari et al. 2017). The qPCR was performed with KAPA SYBR Fast qPCR Kit (KAPA Biosystems, MA, USA) in a 7500 Real Time PCR System thermocycler (Applied Biosystems Thermo Fischer Scientific, MA, USA). The amplification was carried out according to conditions described previously (Aminov et al. 2002; Mroczkowska & Barlow 2008; Muziasari et al. 2017) with a limitation of threshold cycle (CT value) of 30. Normalisation was performed against the 16S rRNA gene using the 2–ΔΔCT method (relative) (Livak & Schmittgen 2001). The results obtained for Black Lake and Donje Bare Lake were compared to those of Plav Lake. In addition, the copy number of targeted genes were expressed relative to 16S rRNA gene copies by the 2–ΔCT method (Schmittgen & Livak 2008). Experiments were carried out in triplicate.

Statistical analysis

To determine whether significant differences exist in the presence of the genetic determinants of antibiotic resistance between three metagenomes from different lakes, samples were analysed by Analysis of Variance (ANOVA) carried out using Rx64 3.5.1 (www.rstudio.com/). Mean values were compared using the least significant difference (LSD) test, and the level of significance was determined at p < 0.05. Principal component analysis (PCA) was performed in order to emphasise the variation of detected ARGs between metagenomes using the package FactomineR of the R software. The statistical significant differences of results obtained by qPCR were determined using Student's t-test. Values at p ≤ 0.05 were considered to be statistically significant.

Nucleotide sequence accession number

Metagenomes sequences data derived from bacterial communities from sediments of Western Balkan glacial lakes were submitted to the GenBank database at NCBI. The sediment metagenome whole genome shotgun (WGS) project from Plav Lake has the project accession number PDVJ00000000. The first version of the project (01) has the accession number PDVJ01000000, and consists of sequences PDVJ01000001–PDVJ01092463 (www.ncbi.nlm.nih.gov/Traces/wgs/?val=PDVJ01#contigs). Sediment metagenome from Black Lake is available under accession number PDVI01000000, and consists of sequences PDVI01000001-PDVI01028961 (www.ncbi.nlm.nih.gov/Traces/wgs/?val=PDVI01#contigs), while metagenome sequences from Donje Bare Lake sediments can be reached under accession number PDVH01000000, and consists of sequences PDVH01000001-PDVH01074614 (www.ncbi.nlm.nih.gov/Traces/wgs/?val=PDVH01#contigs).

Sequence analysis

Metagenomic DNA was sequenced with paired-end 2 × 100 bp reads on an Illumina HiSeq 2000 platform (Macrogen, Seoul, Republic of Korea). A total of 325,366,292 sequence reads were generated from metagenomes. The assembly of reads resulted in 196,038 contigs ranging from 201 to 188,848 bp with an average GC content of 56.65% (Figure S1 in Supplementary Materials). The number of total reads and contigs and GC content (%) for each metagenome sample are presented in Figure S1. The contig dataset was used to determine the functional analysis.

Identification of ARGs in metagenomes using CARD bioinformatics database

The diversity of antibiotic resistance genetic determinants was detected by the CARD bioinformatics database which provides data, models and algorithms related to the molecular basis of antimicrobial resistance (McArthur et al. 2013). Detected ARGs (Figure S2) were classified into several groups: genes encoding efflux pumps of Resistance-Nodulation-Division (RND) superfamily (including: mexF, mexY, mexB, mexD, mexI, smeB, smeE, amrB, ceoB, acrD, mdtB and mdtF) and genes which confer resistance to the following antibiotic classes: β-lactams (blaTEM gene), tetracyclines (tetC gene), aminoglycosides (aacA and aadA genes), chloramphenicol (catA1 gene) and macrolides (ermB gene) (Figure S3). RND efflux pumps genes were the most abundant in Plav Lake metagenome, followed by ermB, blaTEM,aacA and aadA genes. In metagenome from Black Lake, genes for RND efflux pumps and ermB were detected, while in metagenome from Donje Bare Lake all tested resistance genes were found, with the blaTEM gene being the most abundant (Figure S3). The distribution of detected genes in metagenomes is presented on a heat map (Figure 1). The most abundant gene was blaTEM with the highest prevalence in Donje Bare Lake (Figure 1).

Figure 1

Heat map demonstrating the distribution of genetic determinants of antibiotic resistance according to the CARD bioinformatics database. Results were approximated on the relative scale ranging from 0 (blue) as the lowest values, progressing to white, then to 100 (red) as the highest values. Please refer to the online version of this paper to see this figure in colour: http://dx.doi.10.2166/wh.2020.227.

Figure 1

Heat map demonstrating the distribution of genetic determinants of antibiotic resistance according to the CARD bioinformatics database. Results were approximated on the relative scale ranging from 0 (blue) as the lowest values, progressing to white, then to 100 (red) as the highest values. Please refer to the online version of this paper to see this figure in colour: http://dx.doi.10.2166/wh.2020.227.

Close modal

The distribution of genes coding for RND efflux pumps is presented in Figure 2. The most abundant genes encoding RND efflux pumps were mexD, mexF and smeB, all present in metagenome sequence originated from Plav Lake.

Figure 2

Number of hits with RND efflux pumps associated genes according to the CARD database.

Figure 2

Number of hits with RND efflux pumps associated genes according to the CARD database.

Close modal

The PCA plot based on the presence of all identified genetic determinants of antibiotic resistance by CARD database is presented in Figure 3. The dimension 1 (PC1) describes 69.35% of variability between metagenome from Donje Bare Lake (characterised by a positive coordinate on the axis) to metagenomes from Plav Lake and Black Lake (characterised by a negative coordinate on the axis). The dimension 2 (PC2) describes 30.65% of variability between metagenome from Plav Lake to metagenomes Donje Bare Lake and Black Lake (Figure 3(a)). The position of metagenome from Donje Bare Lake was determined by the presence of blaTEM,tetC and catA1 genes, while metagenome from Plav Lake was positioned according to the presence of RND genes (Figure 3(b)).

Figure 3

Individual factor map (a) and variables factor map (b) of principal component analysis (PCA) on different genetic determinants of antibiotic resistance detected in metagenomes by CARD database.

Figure 3

Individual factor map (a) and variables factor map (b) of principal component analysis (PCA) on different genetic determinants of antibiotic resistance detected in metagenomes by CARD database.

Close modal

Identification of integrons and cassette arrays by the integron finder (IF) program

Identification and analysis of integrons and cassette arrays in metagenomes were performed by the Integron Finder (IF) program. IF detects: (1) complete integrons, which include integron integrase gene and at least one attC site; (2) In0 elements, composed of an integron integrase and no attC sites; and (3) the clusters of attC sites lacking integron-integrases (CALIN) elements (Cury et al. 2016). Complete integron was detected only in metagenome from the Plav Lake, located on contig PDVJ01003763.1 and it carried gene for bacterial virulence factor gamma-glutamyl transpeptidase (GGT). Complete integron belongs to the Pseudoxanthomonas suwonensis strain, a class of Gammaproteobacteria, and contained the integrase gene belonging to the phage integrase family. The number of detected In0 elements was three (contigs: PDVJ01015765.1, PDVJ01025354.1, PDVJ01032088.1), none and two (contigs: PDVH01014680.1, PDVH01049549.1) from Plav Lake, Black Lake and Donje Bare Lake, respectively, and all In0 elements contained gene coding integron integrase which belong to the phage integrase family. The total number of detected CALIN elements with two or more attC sites was 38: 20, 4 and 14 from Plav Lake, Black Lake and Donje Bare Lake, respectively. A blast search was performed for all 38 CALIN elements and the results showed that each sample of metagenome sequences had one CALIN element coding antibiotic resistance determinants. CALIN elements coding antibiotic resistance determinants detected by Integron Finder (IF) bioinformatics tool are presented in Figure S4. The CALIN element of Plav Lake metagenome (contig PDVJ01003568.1) carried the gene for N-acetyltransferase enzyme, which is involved in resistance to aminoglycosides. CALIN from Black Lake metagenome (contig PDVI01012771.1) encoded the vicinal oxygen chelate (VOC) family of enzymes which confer resistance to bleomycin and CALIN from Donje Bare Lake metagenome (contig PDVH01000596.1) encoded metallo-β-lactamase enzyme that catalyse the hydrolysis of a broad range of β-lactam drugs (Figure S4).

BlastX analyses of metagenomes for the presence of different integrase gene classes (intI1, intI2 and intI3) are presented in Table 1. All three integrases were most abundant in the metagenome of Plav Lake.

Table 1

Number of hits with different classes of integrase genes

SampleintI1intI2intI3
Plav Lake 186 210 196 
Black Lake 35 17 33 
Donje Bare Lake 135 147 135 
Total 356 374 364 
SampleintI1intI2intI3
Plav Lake 186 210 196 
Black Lake 35 17 33 
Donje Bare Lake 135 147 135 
Total 356 374 364 

Detection of plasmid sequences in metagenome data

Using the BlastN tool, 61 contigs (0.031%) were detected as hits with known plasmid sequences, 35 (0.12%) contigs from Black Lake, 18 (0.019%) contigs from Plav Lake and eight (0.01%) from Donje Bare Lake (Figure S5).

Pfam domains located on putative plasmid sequences

In order to identify plasmid-like sequences within metagenomes, metagenome sequences were searched for the presence of protein sequences for basic plasmid functions (replication, mobilisation, stability, transfer and partitioning), domains commonly used to classify plasmids (Petersen 2011). Pfam families covering known plasmid replicon domains are presented in Table 2. A total of 3292 contigs were identified that encode proteins showing significant homologies (E-value <1e-11) to plasmid related proteins (Table 2) that corresponds to approximately 1.68% of all contigs in three metagenome samples. The percentage of plasmid related contigs within metagenomes was similar and was 1.73% for Plav Lake, 1.59% for Black Lake and 1.64% for Donje Bare Lake. The highest number of contigs (1.08%) was related to Ribon-helix-helix protein (Pfam family PF01402.16) and Firmicutes plasmid replication protein (Pfam family PF05732.6) (0.52%) (Table 2).

Table 2

The number of contigs encoding different plasmid replication-related domains

Pfam familyPfam nameDescriptionPlav LakeBlack LakeDonje Bare Lake
PF01446.12 Rep_1 Replication protein 
PF01719.12 Rep_2 Plasmid replication protein 
PF01051.16 Rep_3 Initiator replication protein 
PF05732.6 RepL Firmicutes plasmid replication protein 536 119 364 
PF07042.6 TrfA TrfA protein 30 21 29 
PF04796.7 RepA_C Plasmid encoded RepA protein 
PF02486.14 Rep_trans Replication initiation factor 
PF01402.16 RHH_1 Ribbon-helix-helix protein, copG family 1008 311 815 
PF01815.11 Rop Rop protein 
PF03428.8 RP-C Replication protein C N-terminal domain 12 
PF10134.4 RPA Replication initiator protein A 
PF06970.6 RepA_N Replication initiator protein A (RepA) N-terminus 
PF06504.6 RepC Replication protein C (RepC) 
PF03090.12 Replicase Replication initiator protein 
Pfam familyPfam nameDescriptionPlav LakeBlack LakeDonje Bare Lake
PF01446.12 Rep_1 Replication protein 
PF01719.12 Rep_2 Plasmid replication protein 
PF01051.16 Rep_3 Initiator replication protein 
PF05732.6 RepL Firmicutes plasmid replication protein 536 119 364 
PF07042.6 TrfA TrfA protein 30 21 29 
PF04796.7 RepA_C Plasmid encoded RepA protein 
PF02486.14 Rep_trans Replication initiation factor 
PF01402.16 RHH_1 Ribbon-helix-helix protein, copG family 1008 311 815 
PF01815.11 Rop Rop protein 
PF03428.8 RP-C Replication protein C N-terminal domain 12 
PF10134.4 RPA Replication initiator protein A 
PF06970.6 RepA_N Replication initiator protein A (RepA) N-terminus 
PF06504.6 RepC Replication protein C (RepC) 
PF03090.12 Replicase Replication initiator protein 

Quantification of antibiotic resistance and integrase genes by qPCR

The quantitative PCR (qPCR) method enabled the detection and quantification of antibiotic resistant and integrase genes in sampled lakes. In the metagenomes of all three lakes the presence of seven genes (blaTEM, catA, ermB, mexA, mexB, intI1 and intI3) was confirmed, but not for the other five selected genes (tetC, aacA, aadA, mexD and intI2). The primer sequences are given in Table 3.

Table 3

List of primers used in this study

Target genePrimer sequence (5′ − 3′)Reference
blaTEM Fw CCGAAGGAGCTAACCGCTT
Rev GTGTCACGCTCGTCGTTTG 
Mroczkowska & Barlow (2008)  
tetC Fw GCGGGATATCGTCCATTCCG
Rev GCGTAGAGGATCCACAGGACG 
Aminov et al. (2002)  
aacA Fw AGAGCCTTGGGAAGATGAAGTTT
Rev TTGATCCATACCATAGACTATCTCATCA 
Muziasari et al. (2017)  
aadA Fw GTTGTGCACGACGACATCATT
Rev GGCTCGAAGATACCTGCAAGAA 
Muziasari et al. (2017)  
catA1 Fw GGGTGAGTTTCACCAGTTTTGATT
Rev CACCTTGTCGCCTTGCGTATA 
Muziasari et al. (2017)  
ermB Fw TAAAGGGCATTTAACGACGAAACT
Rev TTTATACCTCTGTTTGTTAGGGAATTGAA 
Muziasari et al. (2017)  
mexA Fw AGGACAACGCTATGCAACGAA
Rev CCGGAAAGGGCCGAAAT 
Muziasari et al. (2017)  
mexB Fw CTGGAGATCGACGACGAGAAG
Rev GAAATCGTTGACGTAGCTGGAA 
Muziasari et al. (2017)  
mexD Fw TTGCCACTGGCTTTCATGAG
Rev CACTGCGGAGAACTGTCTGTAGA 
Muziasari et al. (2017)  
intI1 Fw GCCTTGATGTTACCCGAGAG
Rev GATCGGTCGAATGCGTGT 
Muziasari et al. (2017)  
intI2 Fw GACGGCTACCCTCTGTTATCTC
Rev GCCACCACTTGTTTGAGGA 
Muziasari et al. (2017)  
intI3 Fw GCCACCACTTGTTTGAGGA
Rev GGATGTCTGTGCCTGCTTG 
Muziasari et al. (2017)  
16S rRNA gene Fw CAGCTCGTGTCGTGAGATGT
Rev CGTAAGGGCCATGATGACTT 
Higgins et al. (2004)  
Target genePrimer sequence (5′ − 3′)Reference
blaTEM Fw CCGAAGGAGCTAACCGCTT
Rev GTGTCACGCTCGTCGTTTG 
Mroczkowska & Barlow (2008)  
tetC Fw GCGGGATATCGTCCATTCCG
Rev GCGTAGAGGATCCACAGGACG 
Aminov et al. (2002)  
aacA Fw AGAGCCTTGGGAAGATGAAGTTT
Rev TTGATCCATACCATAGACTATCTCATCA 
Muziasari et al. (2017)  
aadA Fw GTTGTGCACGACGACATCATT
Rev GGCTCGAAGATACCTGCAAGAA 
Muziasari et al. (2017)  
catA1 Fw GGGTGAGTTTCACCAGTTTTGATT
Rev CACCTTGTCGCCTTGCGTATA 
Muziasari et al. (2017)  
ermB Fw TAAAGGGCATTTAACGACGAAACT
Rev TTTATACCTCTGTTTGTTAGGGAATTGAA 
Muziasari et al. (2017)  
mexA Fw AGGACAACGCTATGCAACGAA
Rev CCGGAAAGGGCCGAAAT 
Muziasari et al. (2017)  
mexB Fw CTGGAGATCGACGACGAGAAG
Rev GAAATCGTTGACGTAGCTGGAA 
Muziasari et al. (2017)  
mexD Fw TTGCCACTGGCTTTCATGAG
Rev CACTGCGGAGAACTGTCTGTAGA 
Muziasari et al. (2017)  
intI1 Fw GCCTTGATGTTACCCGAGAG
Rev GATCGGTCGAATGCGTGT 
Muziasari et al. (2017)  
intI2 Fw GACGGCTACCCTCTGTTATCTC
Rev GCCACCACTTGTTTGAGGA 
Muziasari et al. (2017)  
intI3 Fw GCCACCACTTGTTTGAGGA
Rev GGATGTCTGTGCCTGCTTG 
Muziasari et al. (2017)  
16S rRNA gene Fw CAGCTCGTGTCGTGAGATGT
Rev CGTAAGGGCCATGATGACTT 
Higgins et al. (2004)  

The quantity of individual genes in the metagenomes of Donje Bare Lake and Black Lake was normalised against the quantity in Plav Lake. The ratio of the blaTEM, catA1 and ermB genes between the lakes was similar with an average of 1:9.41 (from 9.16 to 9.65) for Donje Bare Lake and 1:0.35 (from 0.235 to 0.43) for Black Lake (Figure 4). The gene encoding MexB was on average 40% more abundant in Donje Bare Lake and Black Lake than Plav Lake. The quantity of the intI1 gene was 0.12 times higher in Plav Lake compared to other lakes, while the expression of intI3 gene was 1.5 times higher in Donje Bare Lake and 0.5 times lower in Black Lake compared to Plav Lake. The quantities of the mexA and intI1 genes were not significantly different between the analysed lakes.

Figure 4

Quantification of antibiotic resistance and integrase genes by qPCR. Quantity of antibiotic resistance and integrase genes from Black Lake and Donje Bare Lake normalised against values obtained for Plav Lake. Values are the means from results obtained in triplicate. Error bars represent the standard deviation of the mean value. Student's t test was used to compare quantity of the target genes in Plav and Black or Donje Bare Lake (*p< 0.05, **p< 0.01, ***p< 0.001).

Figure 4

Quantification of antibiotic resistance and integrase genes by qPCR. Quantity of antibiotic resistance and integrase genes from Black Lake and Donje Bare Lake normalised against values obtained for Plav Lake. Values are the means from results obtained in triplicate. Error bars represent the standard deviation of the mean value. Student's t test was used to compare quantity of the target genes in Plav and Black or Donje Bare Lake (*p< 0.05, **p< 0.01, ***p< 0.001).

Close modal

The copy numbers of antibiotic resistance and integrase genes compared to the 16S rRNA gene for metagenomes are presented in Table 4. Among seven detected genes the blaTEM gene was the most prevalent within all three lakes. Generally, the numbers of ARGs copies were higher than those obtained for integrase genes.

Table 4

The copy numbers of target genes compared to 16S rRNA gene

Target geneCopy numbers
Plav LakeDonje Bare LakeBlack Lake
blaTEM 1.21 × 10–2 1.13 × 10–1 5.07 × 10–3 
catA1 1.71 × 10–4 1.52 × 10–3 6.31 × 10–5 
ermB 3.55 × 10–3 3.18 × 10–2 8.01 × 10–4 
mexA 3.43 × 10–3 3.62 × 10–3 1.72 × 10–3 
mexB 2.3 × 10–3 3.25 × 10–3 3.16 × 10–3 
intI1 1.64 × 10–4 1.26 × 10–4 1.36 × 10–4 
intI3 1.72 × 10–4 2.62 × 10–4 9.13 × 10–5 
Target geneCopy numbers
Plav LakeDonje Bare LakeBlack Lake
blaTEM 1.21 × 10–2 1.13 × 10–1 5.07 × 10–3 
catA1 1.71 × 10–4 1.52 × 10–3 6.31 × 10–5 
ermB 3.55 × 10–3 3.18 × 10–2 8.01 × 10–4 
mexA 3.43 × 10–3 3.62 × 10–3 1.72 × 10–3 
mexB 2.3 × 10–3 3.25 × 10–3 3.16 × 10–3 
intI1 1.64 × 10–4 1.26 × 10–4 1.36 × 10–4 
intI3 1.72 × 10–4 2.62 × 10–4 9.13 × 10–5 

The PCA plot based on quantified ARGs by qPCR is presented in Figure 5. Dimension 1 (PC1) describes 73.29% of variability between metagenome from Donje Bare Lake (characterised by a positive coordinate on the axis) to metagenomes from Plav Lake and Black Lake (characterised by a negative coordinate on the axis). Dimension 2 (PC2) describes 26.71% of variability between metagenome from Black Lake and Donje Bare Lake to metagenome from Plav Lake (Figure 5(a)). The position of metagenome from Donje Bare Lake was mostly determined by the number of blaTEM,ermB and catA1 genes, and less by the number of mexA and mexB genes (Figure 5(b)). The metagenomes from Black Lake and Plav Lake were positioned according to the number of mexB genes and mexA genes, respectively (Figure 5(b)).

Figure 5

Individual factor map (a) and variables factor map (b) of principal component analysis (PCA) on ARGs quantified by qPCR.

Figure 5

Individual factor map (a) and variables factor map (b) of principal component analysis (PCA) on ARGs quantified by qPCR.

Close modal

This study investigated the abundance and diversity of antibiotic resistance genetic determinants, as well as mobile elements involved in the horizontal gene transfer of resistance to antibiotics, in bacterial communities of three Western Balkans glacial lakes with different input of non-treated sewage and agricultural waste. Low levels of antibiotics and other pharmaceuticals are regularly released into water environments via wastewater, and the concern is that such environmental contamination may serve to create hotspots for antibiotic resistance gene selection and dissemination. ARGs have been regarded as the emerging pollutants and have become an important theme in environmental science. In lakes, pollutants slowly circulate around and their residence time increases due to the longer water retention time relative to rivers (Yang et al. 2018). Because of these characteristics lakes have greater potential to accumulate ARGs. Humans can be exposed to ARGs in lakes through ingestion of water containing ARGs, ingestion of aquatic products or direct contact with the lake such as recreation or swimming (Eckert et al. 2018). In that sense, monitoring of ARGs and mobile elements from bacterial communities in lakes, which can be transferred to pathogenic and human commensal bacterial, is of great importance for public health (Port et al. 2014). A better identification of the natural carriers of ARGs will help to develop strategies to limit resistance spread to pathogens (Stalder et al. 2019).

The risk of infection acquired in natural freshwater environments, like lakes, is higher for those who enjoy swimming and recreation. Recreational water-borne illnesses have been considered as fairly common because water laden with microorganisms or contaminated by human activity gain access to healthy tissues through the skin and body orifices (Ayi 2015). Recreational water-borne illnesses continue to be on the rise and although these infections are often mild, they can be life-threatening or cause major outbreaks (Ayi 2015; Graciaa et al. 2018; Liu et al. 2018). Having this in mind, we recognised that bacterial communities of sediments from Western Balkan glacial lakes are unexplored and should be deeply analysed in order to fully assess the human impact and risk for human health.

Correlation between variations in bacterial composition and diversity with human population impact in sediments of three Western Balkans glacial lakes was previously established (Malešević et al. 2019). Thus, our goal was to determine if resistome and plasmidome within those bacterial communities vary with variation of human population in the vicinity of the lakes, and how their presence can influence public health as components of resistome and plasmidome can be transferred into microbial communities affecting humans (Cantón 2009).

Sediments of Western glacial lakes analysed in this study differed in exposure to the human population, chemical properties, as well as in the composition of bacterial communities (Malešević et al. 2019). However, previous studies suggested that sediment characteristics have no significant effect on the distribution of ARGs (Yang et al. 2017). Another difference between analysed sediments was the size and activity of human population in the vicinity of the lakes (Malešević et al. 2019). The town of Plav is located on the banks of the Plav Lake, and human activities generate domestic and industrial waste, including waste from land farming. Black Lake and Donje Bare Lake are not exposed to effluents generated during land farming, however, Black Lake is a tourist recreational centre. It is considered that Plav Lake is the most exposed to factors that cause pollution of soils, tributary streams, groundwater and lake directly, with a variety of chemical and biological contaminants including organic compounds and pathogens (Malešević et al. 2019). Previous studies indicate that agricultural impacts are the main drivers of the increase in ARGs in lakes and that built-up land use is an important factor driving the distribution of antibiotic resistance in lakes (Czekalski et al. 2015; Yang et al. 2017). However, no data related to the presence of antibiotics and the amount of disposed antibiotics are available for these three glacial lakes. Since Plav is a relatively small town, we could assume that human population creates low selective pressure by antibiotics use in medicine or agriculture.

It is generally considered that ARGs in natural environments that are not under high selective pressure of antibiotics are involved in antibiotic detoxification in the producer organisms (Pang et al. 1994), in the resistance to toxic compounds produced in these environments (Bais et al. 2005) or are involved in cellular processes in the bacteria.

Multi-drug-resistance (MDR) efflux pumps are known to efflux the quinolones (Alonso et al. 1999) which are of synthetic origin, indicating that resistance is not the primary function of MDR pumps. It is well recognised that polyselective efflux pumps like RND pumps that were the most common in the metagenome of Plav Lake, which is the most exposed to human influence, are able to recognise and expel a large panel of molecules including antimicrobials, detergents, heavy metals and organic pollutants (Piddock 2006; Poole 2007; Li & Nikaido 2009). It is tempting to conclude that the detoxifying function of MDR efflux pumps contributed to their enrichment in the Plav Lake which is most exposed to human activities. Indeed, efflux pumps were less frequent in metagenomes of Donje Bare Lake and Black Lake. However, on the other hand, RND efflux pumps which evolved independently of antibiotic consumption are detoxification elements in the environmental conditions and are known to have an important role in resistance to antimicrobials and bacterial pathogenicity in clinical settings (Martinez et al. 2009).

Another hallmark of analysed sediments was the presence of TEM β-lactamases. The TEM β-lactamases represent the clinically significant family of β-lactamases that is commonly encountered in Gram-negative bacteria, and beside TEM-1 and TEM-2, are considered as extended spectrum β-lactamases, of which some are IRT β-lactamases inhibitor-resistant enzymes (Bradford 2001). Among the β-lactamases, the TEM family is of overriding significance with regard to diversity, prevalence and distribution. Previous studies showed high numbers of blaTEM gene copies in surface water and sediments (Lachmayr et al. 2009; Bogaerts et al. 2016; Singh et al. 2016). The high abundance of blaTEM genes in both the Donje Bare Lake and Plav Lake could be explained by their ubiquitous presence as housekeeping genes, which has previously been shown to occur frequently among soil bacteria (Demaneche et al. 2008). Although β-lactamases genes are thought to have originally resided solely on bacterial chromosomes, they are often found on plasmids in specimens collected after the introduction of clinical β-lactam antibiotics. The mobility of these genes is usually connected to transposons or integrons, and they can be also found as components of multiresistance transposons (Lachmayr et al. 2009).

The class 1 integron integrase gene intI1 was detected and quantified in all samples, with the highest abundance of the gene in the bacterial community of Plav Lake sediment. The number of hits for the intI1 gene was significantly higher in the metagenome of Plav Lake sediments compared to others. Also, the trend of increase of the intI1 gene copies that is in correlation with human population influence could be noticed, as determined by quantitative PCR. This could indicate that the human population affects the presence of the intI1 gene in analysed samples. This is consistent with several earlier studies that have shown that integron abundances increase in anthropogenically impacted environments (Rosewarne et al. 2010; Gaze et al. 2011; Gillings et al. 2015). Using the Integron Finder program, the complete integron was detected only in metagenome from Plav Lake and it was placed on contig that shares high similarity to the chromosome of Pseudoxanthomonas suwonensis and encodes bacterial virulence factor GGT. Although previous studies were focused on mobile integrons as important carriers of cassettes encoding antibiotic resistance (Rowe-Magnus & Mazel 2002), today it is well known that chromosomal integrons can carry cassettes encoding virulence factors, secreted proteins and toxin-antitoxin modules (Rowe-Magnus et al. 2003; Boucher et al. 2006). Our findings fit into previous observations that integrons lacking antibiotic resistance determinants can be present in natural populations (Gillings et al. 2009) and that integrons can be found in the chromosome of Gammaproteobacteria (Boucher et al. 2006). We also found few resistance genes in CALIN elements which fits previous suggestions that integrons in natural populations carry a much broader set of adaptive traits than just antibiotic resistance (Cury et al. 2016).

Within our study, one of the goals was to obtain an insight into the plasmidome of glacial lakes sediments bacterial communities. Since plasmids play an important role in the evolution and environmental adaption of microorganisms, they are widely represented in bacterial species and contribute to horizontal gene transfer and thus have a significant ecological and health impact (Thomas & Nielsen 2005). Plasmid-mediated gene transfer plays an important role not only in the mobilisation and dissemination of ARGs, but also in the spread of degradative pathways and virulence determinants of pathogens (Smalla et al. 2015). The development of whole genome sequencing technologies enabled the examination of plasmid DNA content in bacterial communities from environmental samples without cultivation. On the other hand, metagenome datasets are mostly a mixture of chromosomes and plasmids, and as the chromosomal DNA is predominant it is not easy to identify plasmid sequences (Dib et al. 2015). Plasmid sequences in the metagenomes data of three analysed glacial lakes were identified using NCBI Plasmid database and only 61 contigs of the total number of contigs were detected as hits with known plasmid sequences. This result is probably due to the domination of large plasmids in the NCBI plasmid database, while the smaller ones are underrepresented (Jørgensen et al. 2014). As small plasmids could potentially be a reservoir of genetic shuffling (Jørgensen et al. 2014), we have tried to identify plasmid-like sequences and putative plasmids. In this regard, assembled contigs were further searched against the Pfam database (Petersen 2011) and a total of 3292 contigs were detected to encode different plasmid replication-related domains with the highest number of contigs (1608) belonging to the metagenome from Plav Lake. The most common plasmid replication-related domains in observed data sets were RHH_1 (Ribon-helix-helix protein) with a biological role in the regulation of transcription and RepL (Firmicute plasmid replication protein) that is involved in plasmid replication and maintenance. The highest abundance of these two replication domains was detected in Plav Lake, while other replication domains appeared with low frequency in all analysed metagenomes. As plasmids facilitate rapid evolution and adaption of their hosts to the changes of environmental conditions, the highest abundance of plasmid-like sequences was expected in Plav Lake, due to the highest anthropogenic impact to this lake within analysed samples.

The investigation showed that resistome and plasmidome within the analysed three metagenomes differ with variation of human population influence in the vicinity of the lakes and that human population could be one of the main factors that shape the abundance of ARGs and mobile genetic elements within microbial communities from sediments of Western Balkans glacial lakes.

This work was supported by an FEMS research grant for Early Career Scientists (B Filipic), the Ministry of Education, Science and Technological Development of the Republic of Serbia, Serbia (Grant No. 173019) and grant No. CRP/SRB15-02 funded by the International Centre for Genetic Engineering and Biotechnology, Trieste, Italy.

Concept and design of the experiments: BJ, MK, BF. Performance of the experiments: KN, MM, NM. Analysis of the data: BF, DS. Funding of the experiments: MK. All of the authors were involved in drafting and/or revising of the manuscript; approving the final version of the manuscript and agreed to be accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved.

The Supplementary Material for this paper is available online at https://dx.doi.org/10.2166/wh.2020.227.

Altschul
S. F.
Gish
W.
Miller
W.
Myers
E. W.
Lipman
D. J.
1990
Basic local alignment search tool
.
J. Mol. Biol.
215
,
403
410
.
Aminov
R. I.
Chee-Sanford
J. C.
Garrigues
N.
Teferedegne
B.
Krapac
I. J.
White
B. A.
Mackie
R. I.
2002
Development, validation, and application of PCR primers for detection of tetracycline efflux genes of Gram-negative bacteria
.
Appl. Environ. Microbiol.
68
(
4
),
1786
1793
.
Andrews
S.
2010
FastQC: A Quality Control Tool for High Throughput Sequence Data
.
Available at: www.bioinformatics.babraham.ac.uk/projects/fastqc (accessed 9 February 2019)
.
Ayi
B.
2015
Infections acquired via fresh water: from lakes to hot tubs
.
Microbiol. Spectr.
3
(
6
).
doi:10.1128/microbiolspec.IOL5-0019-2015
.
Bais
H. P.
Prithiviraj
B.
Jha
A. K.
Ausubel
F. M.
Vivanco
J. M.
2005
Mediation of pathogen resistance byexudation of antimicrobials from roots
.
Nature
434
,
217
221
.
Bengtsson-Palme
J.
Kristiansson
E.
Larsson
D. G. J.
2018
Environmental factors influencing the development and spread of antibiotic resistance
.
FEMS Microbiol. Rev.
42
(
1
).
doi:10.1093/femsre/fux053
.
Boucher
Y.
Nesbø
C. L.
Joss
M. J.
Robinson
A.
Mabbutt
B. C.
Gillings
M. R.
Doolittle
W. F.
Stokes
H. W.
2006
Recovery and evolutionary analysis of complete integron gene cassette arrays from Vibrio
.
BMC Evol. Biol.
6
(
3
).
doi:10.1186/1471-2148-6-3
.
Cury
J.
Jové
T.
Touchon
M.
Néron
B.
Rocha
E. P.
2016
Identification and analysis of integrons and cassette arrays in bacterial genomes
.
Nucleic Acids Res.
44
(
10
),
4539
4550
.
Demaneche
S.
Sanguin
H.
Pote
J.
Navarro
E.
Bernillon
D.
Mavingui
P.
Wildi
W.
Vogel
T. M.
Simonet
P.
2008
Antibiotic-resistant soil bacteria in transgenic plant fields
.
Proc. Natl. Acad. Sci. USA
105
,
3957
3962
.
Dib
J. R.
Wagenknecht
M.
Farías
M. E.
Meinhardt
F.
2015
Strategies and approaches in plasmidome studies – uncovering plasmid diversity disregarding of linear elements
.
Front. Microbiol.
26
(
6
),
463
.
El-Gebali
S.
Mistry
J.
Bateman
A.
Eddy
S. R.
Luciani
A.
Potter
S. C.
Qureshi
M.
Richardson
L. J.
Salazar
G. A.
Smart
A.
Sonnhammer
E. L. L.
Hirsh
L.
Paladin
L.
Piovesan
D.
Tosatto
S. C. E.
Finn
R. D.
2019
The Pfam protein families database in 2019
.
Nucleic Acids Res.
47
(
D1
),
D427
D432
.
Gaze
W. H.
Zhang
L.
Abdouslam
N. A.
Hawkey
P. M.
Calvo-Bado
L.
Royle
J.
Brown
H.
Davis
S.
Kay
P.
Boxall
A. B.
Wellington
E. M.
2011
Impacts of anthropogenic activity on the ecology of class 1 integrons and integron-associated genes in the environment
.
ISME J.
5
(
8
),
1253
1261
.
Gillings
M. R.
Gaze
W. H.
Pruden
A.
Smalla
K.
Tiedje
J. M.
Zhu
Y. G.
2015
Using the class 1 integron-integrase gene as a proxy for anthropogenic pollution
.
ISME J.
9
(
6
),
1269
1279
.
Graciaa
D. S.
Cope
J. R.
Roberts
V. A.
Cikesh
B. L.
Kahler
A. M.
Vigar
M.
Hilborn
E. D.
Wade
T. J.
Backer
L. C.
Montgomery
S. P.
Secor
W. E.
Hill
V. R.
Beach
M. J.
Fullerton
K. E.
Yoder
J. S.
Hlavsa
M. C.
2018
Outbreaks associated with untreated recreational water – United States, 2000–2014
.
MMWR Morb. Mortal. Wkly. Rep.
67
(
25
),
701
706
.
Hawkey
P. M.
2008
The growing burden of antimicrobial resistance
.
J. Antimicrob. Chemother.
62
(
s1
),
i1
i9
.
Jørgensen
T. S.
Xu
Z.
Hansen
M. A.
Sørensen
S. J.
Hansen
L. H.
2014
Hundreds of circular novel plasmids and DNA elements identified in a rat cecum metamobilome
.
PLoS One
9
(
2
),
e87924
.
Lachmayr
K. L.
Kerkhof
L. J.
Dirienzo
A. G.
Cavanaugh
C. M.
Ford
T. E.
2009
Quantifying nonspecific TEM β-lactamase (blaTEM) genes in a wastewater stream
.
Appl. Environ. Microb.
75
(
1
),
203
211
.
Larsson
D. G. J.
de Pedro
C.
Paxeus
N.
2007
Effluent from drug manufactures contains extremely high levels of pharmaceuticals
.
J. Hazard. Mater.
148
,
751
755
.
Malešević
M.
Mirković
N.
Lozo
J.
Novović
K.
Filipić
B.
Kojić
M.
Jovčić
B.
2019
Bacterial diversity among the sediments of glacial lakes in the western Balkans: exploring the impact of human population
.
Geomicrobiol. J.
36
,
261
270
.
Martinez
J. L.
Sánchez
M. B.
Martínez-Solano
L.
Hernandez
A.
Garmendia
L.
Fajardo
A.
Alvarez-Ortega
C.
2009
Functional role of bacterial multidrug efflux pumps in microbial natural ecosystems
.
FEMS Microbiol. Rev.
33
(
2
),
430
449
.
Mazel
D.
2006
Integrons: agents of bacterial evolution
.
Nat. Rev. Microbiol.
4
(
8
),
608
620
.
McArthur
A. G.
Waglechner
N.
Nizam
F.
Yan
A.
Azad
M. A.
Baylay
A. J.
Bhullar
K.
Canova
M. J.
De Pascale
G.
Ejim
L.
Kalan
L.
King
A. M.
Koteva
K.
Morar
M.
Mulvey
M. R.
O'Brien
J. S.
Pawlowski
A. C.
Piddock
L. J.
Spanogiannopoulos
P.
Sutherland
A. D.
Tang
I.
Taylor
P. L.
Thaker
M.
Wang
W.
Yan
M.
Yu
T.
Wright
G. D.
2013
The comprehensive antibiotic resistance database
.
Antimicrob. Agents Chemotherapy
57
(
7
),
3348
3357
.
Mroczkowska
J. E.
Barlow
M.
2008
Fitness trade-offs in blaTEM evolution
.
Antimicrob. Agents Chemother.
52
(
7
),
2340
2345
.
Muziasari
W. I.
Pitkänen
L. K.
Sørum
H.
Stedtfeld
R. D.
Tiedje
J. M.
Virta
M.
2017
The resistome of farmed fish feces contributes to the enrichment of antibiotic resistance genes in sediments below Baltic sea fish farms
.
Front. Microbiol.
7
,
2137
.
Myers
E. W.
Sutton
G. G.
Delcher
A. L.
Dew
I. M.
Fasulo
D. P.
Flanigan
M. J.
Kravitz
S. A.
Mobarry
C. M.
Reinert
K. H.
Remington
K. A.
Anson
E. L.
Bolanos
R. A.
Chou
H. H.
Jordan
C. M.
Halpern
A. L.
Lonardi
S.
Beasley
E. M.
Brandon
R. C.
Chen
L.
Dunn
P. J.
Lai
Z.
Liang
Y.
Nusskern
D. R.
Zhan
M.
Zhang
Q.
Zheng
X.
Rubin
G. M.
Adams
M. D.
Venter
J. C.
2000
A whole-genome assembly of Drosophila
.
Science
287
(
5461
),
2196
2204
.
Pang
Y.
Brown
B. A.
Steingrube
V. A.
Wallace
R. J.
Jr.
Roberts
M. C.
1994
Tetracycline resistance determinants in Mycobacterium and Streptomyces species
.
Antimicrob. Agents Chemother.
38
(
6
),
1408
1412
.
Petersen
J.
2011
Phylogeny and compatibility: plasmid classification in the genomics era
.
Arch. Microbiol.
193
,
313
321
.
Port
J. A.
Cullen
A. C.
Wallace
J. C.
Smith
M. N.
Faustman
E. M.
2014
Metagenomic frameworks for monitoring antibiotic resistance in aquatic environments
.
Environ. Health Perspect.
122
(
3
),
222
228
.
Reuther
R.
2009
Lake and river sediment monitoring. Environmental monitoring
. In:
Encyclopedia of Life Support Systems (EOLSS)
,
Vol. 2
(
Inyang
H. I.
Daniels
J. L.
, eds.
EOLSS Publishers
,
Paris
,
France
, pp.
120
147
.
Rowe-Magnus
D. A.
Mazel
D.
2002
The role of integrons in antibiotic resistance gene capture
.
Int. J. Med. Microbiol.
292
(
2
),
115
125
.
Rowe-Magnus
D. A.
Guerout
A. M.
Biskri
L.
Bouige
P.
Mazel
D.
2003
Comparative analysis of superintegrons: engineering extensive genetic diversity in the Vibrionaceae
.
Genome Res.
13
(
3
),
428
442
.
Schlüter
A.
Krause
L.
Szczepanowski
R.
Goesmann
A.
Pühler
A.
2008
Genetic diversity and composition of a plasmid metagenome from a wastewater treatment plant
.
J. Biotechnol.
136
(
1–2
),
65
76
.
Schmittgen
T. D.
Livak
K. J.
2008
Analyzing real-time PCR data by the comparative CT method
.
Nat. Protoc.
3
(
6
),
1101
1008
.
Smalla
K.
Jechalke
S.
Top
E. M.
2015
Plasmid detection, characterization, and ecology
.
Microbiol. Spectr.
3
(
1
), PLAS-0038-2014. https://doi.org/10.1128/microbiolspec.PLAS-0038-2014
Stalder
T.
Press
M. O.
Sullivan
S.
Liachko
I.
Top
E. M.
2019
Linking the resistome and plasmidome to the microbiome
.
ISME J.
13
(
10
),
2437
2446
.
Ventola
C. L.
2015
The antibiotic resistance crisis: part 1: causes and threats
.
P T
40
(
4
),
277
283
.

Supplementary data