Metagenomic analysis reveals wastewater treatment plants as hotspots of antibiotic resistance genes and mobile genetic elements

Metagenomic analysis reveals wastewater treatment plants as hotspots of antibiotic resistance genes and mobile genetic elements

Accepted Manuscript Metagenomic analysis reveals wastewater treatment plants as hotspots of antibiotic resistance genes and mobile genetic elements Ji...

4MB Sizes 0 Downloads 16 Views

Accepted Manuscript Metagenomic analysis reveals wastewater treatment plants as hotspots of antibiotic resistance genes and mobile genetic elements Jianhua Guo, Jie Li, Hui Chen, Philip Bond, Zhiguo Yuan PII:

S0043-1354(17)30565-1

DOI:

10.1016/j.watres.2017.07.002

Reference:

WR 13044

To appear in:

Water Research

Received Date: 17 April 2017 Revised Date:

25 June 2017

Accepted Date: 1 July 2017

Please cite this article as: Guo, J., Li, J., Chen, H., Bond, P., Yuan, Z., Metagenomic analysis reveals wastewater treatment plants as hotspots of antibiotic resistance genes and mobile genetic elements, Water Research (2017), doi: 10.1016/j.watres.2017.07.002. This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

ACCEPTED MANUSCRIPT

Animals develop ARB in their guts after taking antibiotics

Person develop ARB in his gut after taking antibiotics Spr gen eads A era R l c o B in mm uni ty

RI PT

on ain ls m re ima an n B c om a R A t fr a me

population

SC

Slaughter house

Hospital

Wastewater produced from various sectors

M AN U

Biosolids reuse

Water reuse

WWTP-ARG Hotspot

Drinking water plant

AC C

EP

TE D

River to sea

Land farming

ACCEPTED MANUSCRIPT

1

Metagenomic analysis reveals wastewater treatment plants as hotspots of antibiotic resistance

2

genes and mobile genetic elements

3

Jianhua Guo1,*, Jie Li2, Hui Chen1, Philip Bond1, Zhiguo Yuan1

4

1

5

QLD 4072, Australia

6

2

7

*Corresponding author: Dr. Jianhua Guo, Phone: +61 (0) 7 3346 3215; Fax: +61 (0) 7 3365 4726;

8

E-mail: [email protected]

RI PT

Advanced Water Management Centre (AWMC), The University of Queensland, St Lucia, Brisbane,

Beijing Genomics Institute (BGI)-Shenzhen, Shenzhen 518083, China

SC

9

Abstract: The intensive use of antibiotics results in their continuous release into the environment and

11

the subsequent widespread occurrence of antibiotic resistant bacteria (ARB) and antibiotic resistance

12

genes (ARGs), and mobile genetic elements (MGEs). This study used Illumina high-throughput

13

sequencing to investigate the broad-spectrum profiles of both ARGs and MGEs in activated sludge

14

and anaerobically digested sludge from a full-scale wastewater treatment plant. A pipeline for

15

identifying antibiotic resistance determinants was developed that consisted of four categories: gene

16

transfer potential, ARG potential, ARGs pathway and ARGs phylogenetic origin. The metagenomic

17

analysis showed that the activated sludge and the digested sludge exhibited different microbial

18

communities and changes in the types and occurrence of ARGs and MGEs. In total, 42 ARGs

19

subtypes were identified in the activated sludge, while 51 ARG subtypes were detected in the digested

20

sludge. Additionally, MGEs including plasmids, transposons, integrons (intI1) and insertion

21

sequences (e.g. ISSsp4, ISMsa21 and ISMba16) were abundant in the two sludge samples. The co-

22

occurrence pattern between ARGs and microbial taxa revealed by network analysis indicated that

23

some environmental bacteria (e.g. Clostridium and Nitrosomonas) might be potential hosts of multiple

24

ARGs. The findings increase our understanding of WWTPs as hotspots of ARGs and MGEs, and

25

contribute towards preventing their release into the downstream environment.

AC C

EP

TE D

M AN U

10

26

1

ACCEPTED MANUSCRIPT

Keywords: Antibiotics; Antibiotic resistant bacteria (ARB); Antibiotic resistance genes (ARGs);

28

Mobile genetic elements (MGEs); Metagenomics; Biological wastewater treatment;

AC C

EP

TE D

M AN U

SC

RI PT

27

2

ACCEPTED MANUSCRIPT

1. Introduction

30

The intensive use of antibiotics for human and veterinary medicine, animal farming, and agri-

31

industrial production results in the continuous release of antibiotics into the environment and

32

subsequently the widespread occurrence of antibiotic resistance in both natural and engineered

33

ecosystems (Alonso et al., 2001; Larson, 2007; Aminov, 2009; Taylor et al., 2011; Udikovic-Kolic et

34

al., 2014). Antibiotic compounds, antibiotic resistant bacteria (ARB) and antibiotic resistance genes

35

(ARGs) are widely found in soil (Torres-Cortés et al., 2011; Joy et al., 2013; Zhu et al., 2013; Chen et

36

al., 2016), surface water (Li et al., 2009a; Pruden et al., 2012; Li et al., 2014), groundwater (Koike et

37

al., 2007; Lapworth et al., 2012), and even deep ocean sediments (Zhang et al., 2009; Chen et al.,

38

2013; Shah et al., 2014). Many ARGs are also found to be transcribed in natural environments

39

(Versluis et al., 2015). Additionally, increasing levels of ARB and ARGs are being detected in clinical

40

settings around the world (Zhang and Li, 2011; Schmieder and Edwards, 2012; Marti et al., 2014).

41

The widespread occurrence and control of ARB and ARGs is a major public health issue and an

42

emerging challenge to deal with worldwide (Pruden, 2014; Berendonk et al., 2015).

M AN U

SC

RI PT

29

Although antibiotic resistance is recognized as a major threat to human public health worldwide, it

44

remains unclear for the diversity, distribution and fate of ARGs in urban water systems (Munck et al.,

45

2015). Some studies have applied amplification-based methods to investigate the occurrence and

46

abundance of various ARGs in natural or engineered ecosystems (Liu et al., 2012; Novo et al., 2013;

47

Calero-Caceres et al., 2014; Mao et al., 2015; Zhang et al., 2016). These studies have been useful to

48

show quantitative polymerase chain reaction (PCR) is an efficient tool to study the widespread

49

occurrence and abundance of ARGs. The PCR detection depends on available primers that are based

50

on known resistance gene and is not suitable for the discovery of novel ARGs (Schmieder and

51

Edwards, 2012). As a powerful technology, metagenomic sequencing can overcome the drawbacks of

52

amplification-based methods (Schmieder and Edwards, 2012; Yang et al., 2014a) and can be used for

53

broad-spectrum screening of ARGs. Recently, there have been a few studies based on metagenomic

54

sequencing to detect ARGs in wastewater treatment plants (WWTPs). These have been conducted to

AC C

EP

TE D

43

3

ACCEPTED MANUSCRIPT

investigate the occurrence and abundance of ARGs in WWTPs (Li et al., 2013; Li et al., 2015; Paiva

56

et al., 2017), to detect differences of ARG abundances in different treatment facilities (Yang et al.,

57

2014a; Bengtsson-Palme et al., 2016), as well as to determine the removal efficiency of ARGs in

58

different wastewater or sludge treatment processes (Christgen et al., 2015). However, few studies have

59

simultaneously investigated the profiles of ARGs and mobile genetic elements (MGEs) (Wang et al.,

60

2013). This is important as a main mechanism for acquired resistance is through horizontal gene

61

transfer of ARGs from the environment or from other bacteria, and the gene transfer is mediated by

62

MGEs such as plasmids, transposons, bacteriophages, insertion sequences and integrons (Stalder et

63

al., 2012; Rizzo et al., 2013). Consequently, it would be advantageous to obtain a more

64

comprehensive picture regarding the fate and transformation of ARGs by simultaneous metagenomic

65

analysis of both ARGs and MGEs in WWTPs. In addition, the most commonly used methods are

66

based on short reads for ARGs annotation, which are not accurate to study antibiotic resistance in the

67

environment (Parnanen et al., 2016). In contrast, if the assembled contigs with long enough sequences

68

are used to the ARG identification, the mobility potential could be assessed.

M AN U

SC

RI PT

55

In this study, we aimed to comprehensively determine the occurrence, abundance and diversity of

70

ARGs and MGEs in a full-scale WWTP treating domestic wastewater. In order to obtain the broad-

71

spectrum profiles of ARGs and MGEs in this WWTP, we applied metagenomic sequencing approach

72

using the Illumina HiSeq 2000 platform and developed an integrated pipeline for identifying antibiotic

73

resistance determinants that consisted of four categories: gene transfer potential, ARG potential,

74

ARGs pathway and ARGs phylogenetic origin. About 9 Gb of high quality sequence data were

75

obtained from two complex sludge samples (aerobic activated sludge (AAS) and anaerobically

76

digested sludge (ADS)), and then the assembled contigs were used for annotation of ARGs and

77

MGEs. The correlation between the occurrence of ARGs and ARB, and the effect of anaerobic

78

digestion for the removal of ARGs and MGEs, were evaluated.

AC C

EP

TE D

69

79 80

2. Materials and Methods

4

ACCEPTED MANUSCRIPT

2.1 Sludge sampling

82

The aerobic activated sludge (AAS, 30 mL) and the anaerobically digested sludge (ADS, 30 mL) were

83

collected from an aeration tank and from a mesophilic anaerobic digester at a full-scale municipal

84

WWTP (Beijing, China) respectively. There was not any reported flu happened in the catchment

85

before sampling (March 2013). This WWTP mainly treats domestic wastewater (without hospital

86

wastewater service), with a mean influent flow of 600,000 m3/day. Bar screens, aerated grit chambers

87

and primary sedimentation are included as preliminary treatment. An Anaerobic-Anoxic-Oxic (A2O)

88

process was adopted for simultaneous nitrification, denitrification and biological phosphorous

89

removal. The hydraulic retention time and the sludge retention time were 8–10 h and 12–15 days,

90

respectively. The average concentrations of chemical oxygen demand (COD), NH4+-N, and PO43--P in

91

the influent were 305 mg/L, 46.4 mg/L, and 4.6 mg/L, respectively, while the average concentrations

92

of those in the effluent were 45 mg/L, 2.4 mg/L and 0.6 mg/L, respectively. The sludge treatment

93

consisted of thickening tanks, anaerobic mesophilic digestion and dewatering. 800 m3/d sludge was

94

treated through the oval sludge digestion system at 35 °C. During the sampling period, the anaerobic

95

digester demonstrated stable performance in terms of volatile solids destruction (an average of 53%)

96

and the biogas yield (1000 m3/d).

97

2.2 DNA extraction

98

Both AAS and ADS samples were mixed with 100% ethanol at a ratio of 1:1 (volume/volume)

99

immediately after being collected from this full-scale WWTP, then were transferred to the lab on ice.

100

Genomic DNA extraction was conducted within 24 hours after sampling, using the FastDNA SPIN

101

Kit for Soil (QBIOgene Inc., Carlsbad, CA, USA) according to the manufacturer’s instructions. The

102

quality of the DNA was examined using gel electrophoresis (1% agarose), and the concentration and

103

purity were measured using a Qubit Fluorometer (Thermo, USA).

104

2.3 Library preparation, sequencing and read assembly

105

The extracted DNA was transported on dry-ice to the Beijing Genomic Institute at Shenzhen, China.

106

Illumina shotgun DNA library construction and sequencing were conducted as follows. After

AC C

EP

TE D

M AN U

SC

RI PT

81

5

ACCEPTED MANUSCRIPT

fragmentation, a paired-end fragment library with length of ~170 base pairs (bp) was constructed.

108

Adaptor-appended fragments were sequenced using the Illumina HiSeq 2000 platform and 90 bp

109

length reads were generated for each end. Quality filtering was firstly conducted for the two

110

metagenomic datasets before the downstream analysis. This was performed by removing raw reads

111

that: contained more than 3 ambiguous nucleotides, were shorter than 35 bp, had more than 15 bp

112

overlap with adapter sequences, included more than 36 nucleotides with quality value lower than 20,

113

or were potential duplicated reads due to amplification artifacts.

RI PT

107

After quality filtering, clean reads were assembled into contigs using SOAPdenovo assembler (v

115

1.05, set as -p 8 -F -M 3 -D 1 -L 90 -u) (Li et al., 2010). A range of values (23−57) for the parameter

116

K (k-mer size) was investigated to obtain the maximum contig number and the maximum value of

117

N50 (Miller et al., 2010). The optimal value of ‘33’ was selected for both samples after the K-mer

118

analysis (Figure S1, Supporting Information). Only contigs longer than 500 bp were used for further

119

analysis. SOAP2 (Li et al., 2008) was then used to align all reads to the contigs for each sample (set as

120

-s 40 –v 5 –l 32 –r 2 –m 130 –x 210 –M 4 –p 1) (Li et al., 2009b). The successfully aligned reads were

121

assigned as ‘assembled reads’. The bioinformatic framework for quantifying the community structure,

122

the abundance of ARGs, and gene transfer potential is schematically shown in Figure 1. In the present

123

study, DNA assembly was firstly conducted to merge short reads into contigs, then run through the

124

following analyses (the detailed bioinformatic pipeline is described as below). The sequences

125

described in this study were deposited in NCBI Sequence Read Archive with accession numbers

126

(SRP064834).

127

2.4 Gene prediction from contigs

128

ORFs were predicted from contigs using MetaGeneMark (v2.10, default setting) (Zhu et al., 2010).

129

The predicted ORFs longer than 100 nt were translated into protein sequences based on the NCBI

130

translation table 11 (Junjie et al., 2010). Reads were mapped back to the non-redundant ORF set for

131

each sample and the coverage for each ORF was calculated as the number of mapped reads.

132

2.5 Taxonomic classification

AC C

EP

TE D

M AN U

SC

114

6

ACCEPTED MANUSCRIPT

Taxonomic classification of contigs was conducted firstly by searching the reads against the NCBI NT

134

database using SOAP2 (V2.21, set as -s 40 -v 5 -l 35 -r 2 -m 130 -x 204 -M 4 -p 1) (Li et al., 2009b),

135

and then further by MEGAN (V5.3) (Huson and Weber, 2013) using the default settings. Sequences

136

were assigned to NCBI taxonomies with MEGAN by using the lowest common ancestor algorithm

137

and the default parameters.

138

2.6 Functional annotation

139

Protein sequences of the predicted genes were searched against the Non-supervised Orthologous

140

Groups (eggNOG, V3.0) (Powell et al., 2014) and the Kyoto Encyclopedia of Genes and Genomes

141

(KEGG, V59) (Tatusov et al., 2000; Kanehisa et al., 2006) databases using BLASTP (V2.2.26) with

142

the E-value cutoff of 10-5. The abundance of a certain COG or KEGG entry in each sample was

143

calculated by the total found genes weighted by their coverage. The eggNOG database provides

144

multiple sequence alignments and broad functional annotation (Powell et al., 2014). The annotated

145

sequences were sorted into 25 functional categories. For the function class of defense mechanisms, we

146

further investigated specific variations of microbial functions at Level 2.

147

2.7 ARGs and MGE analysis

148

The bioinformatic framework for quantifying ARGs and MGE is illustrated as Figure 1. Metagenomic

149

data relevant to antibiotic resistance determinants were classified into four categories: gene transfer

150

potential, ARG potential, ARGs pathway and ARGs origin. Compared to the framework proposed

151

with three categories by (Port et al., 2014), the fourth category, ARGs origin, relates to identifying

152

potential sources of ARGs through linking with community composition profiling was also

153

incorporated in this study. In addition, we identified these four categories based on gene predicted

154

from contigs, rather than directly based on reads, which would provide more reliable annotation

155

results.

AC C

EP

TE D

M AN U

SC

RI PT

133

156

Regarding ARGs annotation, a local database of ARGs was created by downloading all sequences

157

from Antibiotic Resistance Database-ANNOTation (ARG-ANNOT), which is a new platform to

158

detect existing and putative ARGs in bacterial genomes (Gupta et al., 2014). All quality-filtered

7

ACCEPTED MANUSCRIPT

metagenomic data were compared against the ARG database using BLAST (BLASTx) with a cutoff E

160

value of <10-5. ARG types and subtypes in total annotated ARG sequences in both samples were

161

visualized using Circos (Krzywinski et al., 2009). To predict the phylogenetic origins of the ARGs,

162

the predicted genes were simultaneously searched against NCBI nt databases and the ARG-ANNOT

163

database. Cytoscape, as an open source software platform, was used for visualizing complex networks

164

between ARGs and potential ARB. Here, microorganisms that carry ARG sequence were counted as

165

the potential ARB. The BLAST hit outputs were further filtered to annotate different species carrying

166

ARGs using strict criteria with amino acid identity ≥50% and alignment length ≥50 amino acids.

SC

RI PT

159

All sequence reads were aligned to plasmid sequences available in the NCBI RefSeq database and

168

those in the Genbank nucleotide database annotated as plasmids (Kristiansson et al., 2011). The

169

plasmid coverage was obtained by counting the number of reads aligned at each base. A read was

170

annotated as a plasmid in terms of its best BLAST hit with a sequence similarity of above 90% over

171

an alignment of at least 25 amino acids (Wang et al., 2013; Yang et al., 2014a). In addition, local

172

databases of insertion sequences and integron integrase genes were created by downloading ISs

173

sequences from the ISfinder (Siguier et al., 2006) and integrase genes from the INTEGRALL (Moura

174

et al., 2009), respectively. A read was annotated as an insertion sequence or integrase gene according

175

to its best BLASTn hit with a nucleotide sequence identity of above 90% over an alignment of at least

176

50 bp (Wang et al., 2013). Similarly, a cutoff E value of <10-5 was used in the BLAST searches

177

against these two local databases.

178

AC C

EP

TE D

M AN U

167

179

3. Results

180

3.1 Overall taxonomic structure and community diversity analysis

181

Overall, Illumina sequencing generated more than 4.4 Gb of reads for each sample. After quality

182

filtering, each sample contained more than 4.3 Gb high quality sequence reads, which were then

183

assembled using SOAPdenovo assembly algorithm. About 84571 and 82591 contigs longer than 500

184

bp were assembled from the metagenomic reads of Samples AAS and ADS. The minimum contig

8

ACCEPTED MANUSCRIPT

185

length of 500 bp was chosen in order to ensure a reasonable length for open reading frames (ORF)

186

prediction and obtain a dataset of manageable size. The complexity of the full-scale WWTP

187

challenged the assembly. The percentage of reads that could be assembled into contigs over 500 bp

188

was 26.5% for the AAS and 35.1% for the ADS, respectively. Taxonomic annotation of the metagenomic datasets was conducted to reveal microbial

190

compositions of the AAS and ADS samples. Bacteria were the dominant domain in both AAS and

191

ADS, accounting for 99.5% and 76.9% of the DNA sequences respectively (Figure S2). It should be

192

noted that the abundance of Archaea was higher in ADS (22.0%), compared to that in AAS (0.1%).

193

Sequences from Eukaryota and Viruses accounted for only 0.12% and 0.34% in AAS, and 0% and

194

0.01% in ADS respectively.

M AN U

SC

RI PT

189

For better understanding the differences of community structure in these two samples, the

196

taxonomic affiliation at phylum level was determined (Figure 2). In the AAS sample, the microbial

197

community was predominantly comprised of four phyla: Proteobacteria (69.9%) and Nitrospirae

198

(15.4%) were major constituents, followed by Bacteroidetes (8.6%) and Actinobacteria (1.5%) (Figure

199

2). In comparison, the most abundant organism populations in the ADS sample were Cloacimonetes,

200

Euryarchaeota, Proteobacteria, Firmicutes, which accounted for 24.6%, 21.9%, 19.6% and 6.6% of the

201

community, respectively. The dominance of Nitrospirae in the AAS likely accounts for the good

202

nitrification performance in this WWTP (the mean ammonium removal efficiency was > 90%), while

203

the proliferation of Euryarchaeota in the ADS supports good methane production in the digester.

EP

AC C

204

TE D

195

205

3.2 Functional profile of microbial communities

206

For functional annotation, reads were assembled into contigs, and subsequently genes were predicted

207

from the contigs and functionally annotated through the KEGG and eggNOG databases. 141,1438 and

208

152,329 genes with average length of more than 600 bp were predicted from the AAS and ADS

209

samples respectively. Most of the functions annotated based on KEGG and eggNOG databases were

210

shared by both samples (Figures S3 and S4).

9

ACCEPTED MANUSCRIPT

A total of 127,836 genes of the AAS sample and 126,356 genes of the ADS sample were annotated

212

using the eggNOG database. Functional gene annotation based on COG classes revealed that the

213

relative distribution of the 25 basic metabolic categories of the anaerobic and aerobic sludge

214

metagenomes was to a great extent similar (Figure S3). In addition to general function prediction,

215

these two samples harbored more genes involved in amino acid transport and metabolism, energy

216

production and conversion, inorganic ion transport and metabolism, transcription, and carbohydrate

217

transport and metabolism (Figure S4). In particular, we further compared functional profiles related to

218

the top-23 significant OGs of defense mechanisms (Figure 3).

SC

RI PT

211

These OGs can be divided into three major functions: (i) ARGs like antimicrobial peptide

220

transporters, acriflavin resistance protein and beta-lactamase (COG0577, COG1619, COG1680,

221

dproNOG00127 and dproNOG01088); (ii) the restriction enzymes (COG1002, COG0286, COG4096,

222

COG0610, COG0732, arCOG00878 and arCOG02632); and (iii) multidrug efflux pump and

223

transporter (COG0841, COG0842, COG1132, COG1132 and arCOG02841) (Figure 3). Comparison

224

of the defense mechanism genes detected in the two samples revealed some differences in

225

abundances. For example, the number of ABC-type multidrug transporters was relatively higher in the

226

ADS in comparison to the AAS. The abundance of types of predicted proteins involved in antibiotic

227

response suggests that the WWTP harbors a diverse range of ARGs.

228

3.3 Occurrence, abundance and diversity of ARGs

229

In order to explore the profile of ARGs present in the municipal WWTP, we conducted BLAST

230

analysis of the AAS and the ADS metagenomes against the ARG-ANNOT database. This analysis

231

indicated that a total of 11,429 reads (0.049%) of the sample AAS and 10,110 reads (0.051%) of the

232

sample ADS were assigned to 10 and 9 types of the known ARGs respectively (Figure 4). The most

233

prominent ARG types found in both samples were genes conferring resistance to fluoroquinolones,

234

tetracycline,

235

aminoglycosides, glycopeptides, phenicols, and trimethoprim (Figure 4). After the anaerobic sludge

236

digestion treatment, the dominant ARG types and subtypes exhibited some changes in abundance. The

AC C

EP

TE D

M AN U

219

beta-lactam,

macrolide-lincosamide-streptogramin

(MLS),

sulfonamides,

10

ACCEPTED MANUSCRIPT

dominant ARG types in the AAS were genes conferring resistance to beta-lactamases (38.1%),

238

fluoroquinolones (22.7%), glycopeptides (15.7%) and sulfonamides (9.2%). In contrast, the dominant

239

types in the ADS were glycopeptides (33.0%), MLS (27.9%), fluoroquinolones (14.2%) and

240

sulfonamides (6.8%). In addition, the abundances of the ARGs relating to confer resistance to

241

tetracyclines and rifampicin were greater in the AAS (6.5% and 1.1%) in comparison to those in the

242

ADS (5.7% and undetectable). In contrast, the abundance of the ARGs relating to aminoglycosides

243

and phenicols were greater in the ADS sample (3.3% and 0.8%) than the AAS (0.6% and 0.4%).

RI PT

237

We further analyzed the sub-types of ARGs in the sludge samples (Figure 5 and Table S1). The

245

dominant subtypes in the AAS were oqxBgb (fluoroquinolones resistance gene, 18.9%), pep-EC

246

(beta-lactamases resistance gene, 11.5%) and penA (beta-lactamases resistance gene, 11.3%), while

247

the dominant subtypes in the ADS were oqxBgb (fluoroquinolones resistance gene, 17.2%), vanR-F

248

(glycopeptides resistance gene, 11.6%) and vanR-M (glycopeptides resistance gene, 10.0%).

249

Additionally, the multidrug resistance genes, tetracycline resistance genes (tet) and sulfonamide

250

resistances genes (sul) had high abundance in both samples, each accounting for around 4% of the

251

reads involved in antibiotic resistance (Table S1). Also, we detected that the ARG subtype number

252

had slightly increased from 42 in the AAS to 51 in the ADS. Moreover, there were 18 shared ARG

253

subtypes detected in both samples (Figure S5 and Table S1).

254

3.4 Correlation between ARGs and bacteria

255

To reveal which microorganisms were carrying the ARG genes, the network analysis approach was

256

used to investigate the co-occurrence patterns between ARG subtypes and microbial taxa (Figure 6).

257

Out of the potential 42 and 41 hosts in AAS and ADS, there were 9-shared microorganisms, although

258

the abundances of these microorganisms varied. Interestingly, we detected potential multiple ARGs in

259

a few hosts. For example, seven ARG subtypes were assigned to Dechloromonas in the AAS sample

260

(pep-EC, oqxA, oqxBgb, vanR-B, vanRc4, vanR-Pt, and dfrA3), while three different ARG subtypes

261

(vanR-Pt, vgaB and tetO) were affiliated with Clostridium in the ADS sample. In addition, we

262

observed potential multidrug resistance in function microorganisms of WWTPs, including

AC C

EP

TE D

M AN U

SC

244

11

ACCEPTED MANUSCRIPT

Nitrosomonas (a typical ammonia-oxidizing bacteria), which potentially carried 5 ARG subtypes

264

(penA, oqxBgb, vanHAc2, vanR-F, and dfrK) in the AAS, and Candidatus ‘Accumulibacter’ (a

265

typical organism for phosphorus removal), which carried 4 ARGs (peb_EC, SFO-1, vanR-B, and

266

vanR-C). The presence of these environmental ARB (Ashbolt et al., 2013) might be associated with

267

two reasons: (1) horizontal gene transfer between pathogenic bacteria and indigenous environmental

268

bacteria might be existed; (2) selection pressures (e.g. antibiotics or heavy metal in wastewater) result

269

in the development of environmental ARB.

270

3.5 Occurrence, abundance and diversity of MGEs

271

In this study, we also determined the presence of MGEs, including plasmids, conjugative transposons,

272

ISs and integrons in the sludge metagenomes (Table 1), which play important roles in the mobility and

273

acquisition of ARG among various microorganisms via the horizontal gene transfer pathway. In total,

274

there were 295,370 reads and 1,443,924 reads belonging to Class II transposons in the sample AAS

275

and ADS, respectively. Conjugative transposons, like conjugative plasmids contain an origin of

276

transfer and the genes required to make the conjugation apparatus (van Hoek et al., 2011). These

277

elements probably contribute as much as plasmids to the spread of ARG in some genera of pathogens.

278

A total of 233,999 reads of the AAS metagenome and 1,064,347 reads of the ADS sequences were

279

matched to long interspersed elements (LINEs) (Table S2). Additionally, a total of 42,528 reads of the

280

AAS and 269,748 reads of the ADS belonged to short interspersed elements (SINEs). If an ARG is on

281

a conjugative or mobilizable plasmid, then this ARG has the potential to transfer among bacteria

282

(Bennett, 2008; van Hoek et al., 2011), so that resistance can spread throughout a given ecosystem. In

283

the metagenomic data, a total of 93,004 reads (0.38%) of the AAS and 31,842 reads (0.13%) of the

284

ADS were matched to plasmids (Table S3).

AC C

EP

TE D

M AN U

SC

RI PT

263

285

Insertion sequences (IS) are segments of bacterial DNA that can transfer from one position on a

286

chromosome to another position on the same chromosome or onto a different chromosome. Our

287

search detected a total of 91,043 reads (0.38%) from the AAS and 45,780 reads (0.19%) from the

288

ADS that were matched to 255 and 142 types of known IS respectively. (Table S4). It was seen that 12

ACCEPTED MANUSCRIPT

the two samples shared 77 common IS types and all the types detected belonged to 21 IS families

290

(Table S4). Among the ISs in the aerobic activated sludge, ISSsp4 (4259 reads, 0.018%) had the

291

highest abundance, followed by IS994 (4121 reads, 0.017%). ISSsp4, a member of the IS21 family, is

292

usually hosted by Sphingomonas sp. (Shintani et al., 2007) , and IS994, a member of the IS3 family, is

293

originally detected from the fish pathogen Renibacterium salmoninarum (Rhodes et al., 2000). Most

294

of members belonging to the family of IS21 or IS3 have active horizontal transmission (Mahillon and

295

Chandler, 1998). In contrast, the ADS was dominated by ISMsa21 (9398 reads, 0.039%), ISMba16

296

(9237 reads, 0.039%) and ISPa38 (2625 reads, 0.011%). Both ISMsa21and ISMba16 belong to the

297

IS200/IS605 family, which was originally detected in Salmonella typhimurium, but related elements

298

have been found in a variety of bacterial and archaeal species (Mahillon and Chandler, 1998; Filee et

299

al., 2007). In this study, these two ISs were carried by the genus of Methanosarcina, which had high

300

abundance in the anaerobically digested sludge. Compared to the IS families of IS21 or IS3, IS200 has

301

a significantly lower horizontal transfer frequency (Mahillon and Chandler, 1998).

M AN U

SC

RI PT

289

Alignment against the INTEGRALL database indicated that a total of 5026 reads (0.021%) of the

303

AAS and 9119 reads (0.038%) of the ADS were assigned to integrase genes (Table S5). The most

304

abundant integrase gene was intI1, accounting for 91.1% and 100% of alignment hits of the AAS and

305

the ADS respectively. Previous studies have also reported the widespread occurrence of integrase

306

genes in WWTPs (Wang et al., 2013), including class 1 integrons carrying various ARGs. Our results

307

also detected a number of reads from both sludges annotated as gene cassettes (Table S6). As multiple

308

cassettes can reside together in an integron, this could be an important mechanism for providing

309

resistance to multiple antibiotics in the sludges. For example in the AAS, 85 reads were annotated as

310

catB8-like/aacA4. This cassette enables resistance to chloramphenicol, amikacin, dibekacin,

311

isepamicin, netilmicin, sisomicin and tobramycin (Moura et al., 2012). Additionally, in the ADS, 611

312

reads were aligned to aadA2, which confers resistance to streptomycin and spectinomycin.

AC C

EP

TE D

302

313 314

4. Discussion

13

ACCEPTED MANUSCRIPT

4.1 WWTPs might be hotspots of ARB and ARGs in the urban water cycle

316

This study adopted high-throughput sequencing to provide comprehensive insight in microbial

317

community structures and the diversity and abundance of ARGs and MGEs occurring in a full-scale

318

WWTP. There is a distinct difference in microbial communities between aerobic activated sludge and

319

the anaerobically digested sludge: Proteobacteria, Nitrospirae, Bacteroidetes and Actinobacteria

320

dominated in the activated sludge, while the most dominant phyla in anaerobic digestion sludge were

321

Cloacimonetes, Euryarchaeota, Proteobacteria, and Firmicutes. Our metagenomic data also suggested

322

the prevalence of a variety of ARGs in the WWTP, where a total of 42 ARGs subtypes belonging to

323

10 ARGs types were identified in the AAS, while 51 subtypes belonging to 9 ARGs types were

324

detected in the ADS. Among them, the prominent ARG types found in both samples were able to

325

resistance to fluoroquinolones, tetracycline, beta-lactam, MLS, sulfonamides, aminoglycosides,

326

glycopeptides, phenicols, and trimethoprim. Fluoroquinolone resistance genes (oqxBgb) had the

327

highest abundance in both aerobic activated sludge (around 18.9% of the total ARGs reads), and the

328

digested sludge (accounting for around 17.2% of the total ARGs reads). The abundances of ARGs

329

detected in this study (0.049% in AAS and 0.051% in ADS) were similar to those detected in other

330

studies of sewage (0.06%) (Yang et al., 2014a) and antibiotic contaminated river sediments (0.02%–

331

1.71%) (Kristiansson et al., 2011). However, this proportion of ARGs was relatively higher than that

332

detected in several previous studies of activated sludge (0.003 – 0.007%) (Zhang and Li, 2011; Yang

333

et al., 2014a), the effluent of WWTPs (0.008 – 0.012%) (Zhang and Li, 2011; Yang et al., 2014a),

334

drinking water (0.0004 – 0.0071%) (Shi et al., 2013), marine water (0.0017%) (Port et al., 2012), and

335

sludge in a Tannery WWTP (0.0081 – 0.0101%) (Wang et al., 2013). The difference of percentage

336

between the studies may be due to different bioinformatic methods or parameters used to estimate the

337

proportion of ARGs. On one hand, the enriched occurrence of these ARGs in this WWTP might be

338

associated with the antibiotics discharged from the upstream household wastewater. This has been

339

detected in the sewage and influent of this WWTP (Sui et al., 2010). The level of antibiotics in the

340

influent would likely affect the microbial community by selecting for ARB. Other studies have also

AC C

EP

TE D

M AN U

SC

RI PT

315

14

ACCEPTED MANUSCRIPT

341

indicated that WWTPs might be a hotspot for the occurrence of ARGs (Zhang and Li, 2011; Yang et

342

al., 2014a; Xu et al., 2015). On the other hand, heavy metals and other emerging contaminations (e.g.

343

triclosan) could also pose selective pressures to stimulate the horizontal transfer of ARGs (Carey et al.,

344

2016). WWTPs might be hotspots for the occurrence of ARGs in urban water systems (Rizzo et al., 2013),

346

and perhaps contribute to their broader distribution. Firstly, WWTPs receive various wastewaters

347

from households, slaughterhouses and hospitals, which contain various antibiotics, and potentially the

348

wastewater microbiota has large proportions of ARB and ARGs. Recent metagenomic analyses are

349

detecting a diverse range of ARGs in the influents of WWTPs (Zhang and Li, 2011; Yang et al.,

350

2014a; Xu et al., 2015). Secondly, niches existing in wastewater or sludge treatment facilities (e.g.

351

antibiotic residues at sub-inhibitory concentrations and high biomass concentrations) might be

352

favorable to stimulate the horizontal transfer of ARGs (Stalder et al., 2012; Rizzo et al., 2013). In this

353

study, our metagenomic analysis identified very diverse and high levels of MGEs, which included

354

integrons, transposons, IS and plasmids, that could facilitate activities of horizontal gene transfer.

355

Finally, the effluent or biosolid discharged by WWTPs will contain ARB and ARGs, which will

356

contaminate and significantly affect the abundance of ARGs in downstream environments. This would

357

include the receiving water bodies and land using the effluent as an irrigation source (Lupo et al.,

358

2012; Pruden et al., 2012; Yang et al., 2014a). Previous reports demonstrate that effluents discharged

359

from WWTPs and landfill disposal of wasted sludge contribute to the widespread occurrence of ARGs

360

in natural environments (LaPara et al., 2011; Munir et al., 2011; Lupo et al., 2012; Pruden et al., 2012;

361

Threedeach et al., 2012). Consequently, effective measures for the removal of ARGs during

362

wastewater treatment need to be developed.

363

4.2 An integrated pipeline for identifying antibiotic resistance determinants

364

A new framework to annotate antibiotic resistance determinants was developed, which include four

365

categories: gene transfer potential, ARG potential, ARGs pathway and ARGs origin (Figure 1). In

366

comparison of the widely used framework in previous studies, this pipeline is able to identify potential

AC C

EP

TE D

M AN U

SC

RI PT

345

15

ACCEPTED MANUSCRIPT

sources of ARGs through linking with community composition profiling based on the assembled

368

contigs. This proposed method provides more reliable annotation results because the assembled

369

contigs were used, rather than short reads.

370

So far, most metagenomic studies are using the short Illumina sequence reads in taxonomic

371

assignment and ARG annotation, which would generate unreliable results. As indicated by our

372

previous results (Guo et al., 2017), merging short reads into contigs increased annotation quality in

373

taxonomic, because of the increased length of the sequences. Yang et al. (2014b) also claimed that

374

using longer sequences not only improved the taxonomic analysis quality, but also increased

375

functional annotation quality. More importantly, our study enables simultaneously investigation of the

376

full profiles of ARGs and MGEs in wastewater treatment systems, which is advantageous in obtaining

377

a comprehensive picture regarding the fate and transformation of ARGs in engineered ecosystems.

378

The ARG identified have the potential resistance to antibiotics, and thus represent potential risks to

379

public health. It should be noted that the lack of functional demonstrations for ARGs in environmental

380

metagenomes are considerable limitations for the characterization of the environmental resistome and

381

the assessment of its clinical relevance, which warrants further studies.

TE D

M AN U

SC

RI PT

367

A few studies have used metagenomic sequencing to detect ARGs in engineered ecosystems, but

383

the correlation of ARGs and their possible hosts through metagenomics analysis is rarely determined

384

(Wang et al., 2013). By determining the taxonomic origin of ARGs, this study found quite diverse

385

possible hosts, which could potentially carry ARGs in both samples. In addition to pathogens

386

(Pseudomonas aeruginosa, Pasteurella multocida and Escherichai coli), a few of microorgnisms

387

belong to indigenous environmental bacteria (e.g. Nitrosomonas and Candidatus ‘Accumulibacter’),

388

implying the occurrence of horizontal gene transfer between pathogen and environmental bacteria

389

within the components of the WWTP. To date, information regarding non-pathogenic environmental

390

species, which can be important carriers of ARGs, is scarce. It is required to assess the potential risk

391

of these non-pathogenic environmental species on the spread of antibiotic resistance. Multidrug-

392

resistant bacteria has become a major public health concern in many countries as there are fewer

AC C

EP

382

16

ACCEPTED MANUSCRIPT

effective antimicrobial agents to treat infections caused by these bacteria (Ibrahim et al., 2012; van der

394

Donk et al., 2012). The potential for these systems to adding to the problem of multiple drug

395

resistance may be a concern for the operation of WWTP. In this study potential multidrug resistance

396

was also observed in these environmental bacteria mentioned above. It should be noted that the

397

taxonomic origin of ARGs was conducted through simultaneously searching the predicted genes

398

against both the taxonomic database and the ARG database. The obtained results could indicate the

399

possible host information of AGRs if the predicted genes mapped to both the ARGs and microbial

400

taxa sequences. The correlation between ARGs and bacteria revealed by this approach need to be

401

further validated using culture-based approaches or function metagenomics (Mullany, 2014).

402

4.3 Role of anaerobic sludge digestion on removal of ARGs and MGEs

403

It is hypothesized that the niches existing in anaerobic digesters (e.g. pathogen removal and DNA

404

hydrolysis under starvation conditions) might reduce ARGs abundance and decrease the horizontal

405

transfer of ARGs that is facilitated by MGEs. In this study, we also compared the types and

406

abundance of ARGs in sludge samples prior to and following anaerobic digestion. It seems that the

407

sludge digestion process had a mixed effect. Following the digestion there was a loss of some ARGs,

408

but also a number of unique ARGs emerged. These lost ARGs (in total 29) included resistance genes

409

for beta-lactamases, trimethoprim, glycopeptides, rifampicin, aminoglycosides and phenicol (Table

410

S1). In contrast, following the anaerobic digestion, 33 ARGs were exclusively detected in the ADS,

411

which might imply the sludge treatment induced the occurrence of ARGs. These post digestion ARGs

412

included resistance genes for MLS, aminoglycosides, beta-lactamases, trimethoprim, phenicols,

413

glycopeptides and sulfonamides (Table S1). Among the 18 persistent ARG subtypes detected in both

414

samples, only 8 ARGs had reduced proportions after the anaerobic digestion process. This included

415

the fluoroquinolone resistance protein (oqxBgb), the glycopeptide resistance protein (vanR-C),

416

tetracycline resistance genes (tetC and tet-41), beta-lactamase resistance protein (OXA-53), and

417

sulfonamide resistance genes (sul1 and sul2). In addition to the subtypes exclusively emerging after

418

the anaerobic digestion, described above, around 10 persistent ARG subtypes were enriched in the

AC C

EP

TE D

M AN U

SC

RI PT

393

17

ACCEPTED MANUSCRIPT

419

sludge after the digestion (Table S1). The different percentage and compositions of ARGs might be

420

attributed to the changes of the community structure caused by sludge digestion and consequently the

421

proportional representation of taxa that harbor these genes and mobile elements. To date the effects of the anaerobic digestion process on the removal or promotion of MGEs have

423

received limited attention. In this study, metagenomic analysis identified very high levels of several

424

types of MGEs for horizontal gene transfer that included integrons, conjugative transposons, plasmids

425

and insertion sequences. Our results showed that the anaerobic digestion process has potential to

426

decrease the abundance of MGEs (in particular for plasmids and insertion sequences, which is

427

consistent with a previous study (Miller et al., 2016). They also reported the anaerobic sludge

428

digestion lowered the occurrence of plasmids, thus reducing the transfer risk of ARGs caused by

429

conjugative plasmids. The abundance of plasmids in activated sludge (93,004 reads matching to

430

known plasmids, 0.38%) was higher than that in digested sludge (31,842 reads, 0.13%). Similarly, a

431

total of 255 types of insertion sequences (91,043 reads, 0.38%) detected in the activated sludge, were

432

reduced to only 142 types (45,780 reads, 0.19%) detected in the digested sludge. We detected

433

insertion sequences belonging to IS families IS21 or IS3 had reduced occurrence following the

434

digestion of the sludge. These MGE are most active for horizontal gene transmission (Mahillon and

435

Chandler, 1998) and such reduction of these could effect the horizontal gene transmission rates of the

436

subsequently disposed sludge.

EP

TE D

M AN U

SC

RI PT

422

Several studies have been conducted to investigate the role of sludge digestion on affecting ARGs

438

abundance and reducing horizontal transfer of ARGs, however, contrasting results have been obtained

439

so far. Ghosh et al. (2009) monitored the fate of three tetracycline resistance genes (tetA, tetO and

440

tetX) by q-PCR in a full-scale anaerobic digester and concluded that the thermophilic anaerobic

441

digester exhibited better performance than that of a mesophilic anaerobic digester for removing these

442

tetracycline resistance genes. This might be attribute to incoming ARB death and narrowed host range

443

for horizontal gene transfer under favourable operating conditions in thermophilic digesters (e.g.

444

higher temperature) (Miller et al., 2016). Additionally, they demonstrated that the thermophilic

AC C

437

18

ACCEPTED MANUSCRIPT

treatment was more efficient in reducing integrons (intI1) compared to the mesophilic digestion.

446

Similarly, Diehl and LaPara (2010) studied the removal of five tetracycline-related ARGs (tetA, tetL,

447

tetO, tetW, and tetX) in lab-scale anaerobic digesters. They also suggested the removal of these ARGs

448

was enhanced with the increasing digestion temperature. In contrast to the above q-PCR studies,

449

Zhang et al. (2015) applied metagenomic sequencing to compare the removal of ARGs in lab-scale

450

thermophilic and mesophilic anaerobic digestion reactors. They showed that several ARGs (aadA,

451

macB, and sul1) were enriched in the thermophilic anaerobic digester, while the abundance of

452

resistance genes of erythromycin esterase type I, sul1, and tetM were increased in the mesophilic

453

anaerobic digester. Also based on metagenomic sequencing, Yang et al. (2014a) indicated that

454

anaerobically digested sludge might carry more ARGs than aerobic activated sludge. The effects of

455

sludge digestion on the abundances of ARGs and MGE is likely associated with operational

456

conditions that would include the feed sludge type, reactor temperature, and hydraulic retention time.

457

It cannot be excluded that the lower abundances of ARGs and MGE are caused by the different

458

microbial community structure in anaerobic digestion sludge, compared to aerobic activated sludge.

459

Obviously more studies are required to better reveal the fate of various ARG subtypes and MGE in the

460

anaerobic digestion process.

EP

461

TE D

M AN U

SC

RI PT

445

5. Conclusions

463

This study adopted high-throughput metagenomic sequencing to provide insights in microbial

464

community structures and the diversity and abundance of ARGs and MGEs occurring in a full-scale

465

WWTP. A new framework to annotate antibiotic resistance determinants was proposed that consisted

466

of four categories: gene transfer potential, ARG potential, ARGs pathway and ARGs origin. In

467

comparison of the widely used framework in previous studies, we are able to identify potential

468

sources of ARGs through linking with community composition profiling based on the predicted genes.

469

Our results demonstrate that the activated sludge and the digested sludge exhibited different microbial

470

communities and changes in the types and occurrence of ARGs and MGEs. A variety of ARGs was

AC C

462

19

ACCEPTED MANUSCRIPT

found both in samples of aerobic activated sludge and in anaerobically digested sludge, including

472

ARGs encoding resistance for fluoroquinolones, tetracycline, beta-lactam, MLS, sulfonamides,

473

aminoglycosides, glycopeptides, phenicols, and trimethoprim. It was seen that MGEs were abundant

474

in the two sludge samples, and these are likely important for the acquisition and mobility of various

475

ARGs among the bacterial species. This study also reported potential multidrug resistance was also

476

detected in environmental bacteria. It is evident that further studies are required to better understand

477

the fate of these genetic elements during wastewater treatment processes in efforts to prevent their

478

release from WWTP effluents and from the reuse of treated biosolids.

SC

RI PT

471

480

Conflict of Interest

481

The authors declare no conflict of interest.

482

M AN U

479

Acknowledgements

484

This work is supported by Australian Research Council Future Fellowship (FT170100196). Jianhua

485

Guo also acknowledges the support by the University of Queensland ECR Biomedical Project. The

486

authors would like to thank Mr. Wenkui Dai at Beijing Genomics Institute (Shenzhen) for suggestion

487

about sequencing.

488

Supplementary Information accompanies this paper on the website.

490

EP

AC C

489

TE D

483

491

References

492 493 494 495 496 497 498 499 500

Alonso, A., Sánchez, P., Martínez, J.L., 2001. Environmental selection of antibiotic resistance genes. Environmental Microbiology 3 (1), 1-9. Aminov, R.I., 2009. The role of antibiotics and antibiotic resistance in nature. Environmental Microbiology 11 (12), 2970-2988. Ashbolt, N.J., Amezquita, A., Backhaus, T., Borriello, P., Brandt, K.K., Collignon, P., Coors, A., Finley, R., Gaze, W.H., Heberer, T., Lawrence, J.R., Larsson, D.G.J., McEwen, S.A., Ryan, J.J., Schonfeld, J., Silley, P., Snape, J.R., Van den Eede, C., Topp, E., 2013. Human Health Risk Assessment (HHRA) for Environmental Development and Transfer of Antibiotic Resistance. Environmental Health Perspectives 121 (9), 993-1001. 20

ACCEPTED MANUSCRIPT

EP

TE D

M AN U

SC

RI PT

Bengtsson-Palme, J., Hammaren, R., Pal, C., Ostman, M., Bjorlenius, B., Flach, C.F., Fick, J., Kristiansson, E., Tysklind, M., Larsson, D.G.J., 2016. Elucidating selection processes for antibiotic resisitance in sewage treatment plants using metagenomics. Science of the Total Environment 572, 697-712. Bennett, P.M., 2008. Plasmid encoded antibiotic resistance: acquisition and transfer of antibiotic resistance genes in bacteria. British Journal of Pharmacology 153 (Suppl 1), S347-S357. Berendonk, T.U., Manaia, C.M., Merlin, C., Fatta-Kassinos, D., Cytryn, E., Walsh, F., Buergmann, H., Sorum, H., Norstrom, M., Pons, M.-N., Kreuzinger, N., Huovinen, P., Stefani, S., Schwartz, T., Kisand, V., Baquero, F., Luis Martinez, J., 2015. Tackling antibiotic resistance: the environmental framework. Nature Reviews Microbiology 13 (5), 310-317. Calero-Caceres, W., Melgarejo, A., Colomer-Lluch, M., Stoll, C., Lucena, F., Jofre, J., Muniesa, M., 2014. Sludge as a potential important source of antibiotic resistance genes in both the bacterial and bacteriophage fractions. Environmental Science & Technology 48 (13), 7602-7611. Carey, D.E., Zitomer, D.H., Hristova, K.R., Kappell, A.D., McNamara, P.J., 2016. Triclocarban influences antibiotic resistance and alters anaerobic digester microbial community structure. Environmental Science & Technology 50 (1), 126-134. Chen, B., Yang, Y., Liang, X., Yu, K., Zhang, T., Li, X., 2013. Metagenomic Profiles of Antibiotic Resistance Genes (ARGs) between Human Impacted Estuary and Deep Ocean Sediments. Environmental Science & Technology 47 (22), 12753-12760. Chen, B., Yuan, K., Chen, X., Yang, Y., Zhang, T., Wang, Y., Luan, T., Zou, S., Li, X., 2016. Metagenomic analysis revealing antibiotic resistance genes (ARGs) and their genetic compartments in the Tibetan environment. Environmental Science & Technology. Christgen, B., Yang, Y., Ahammad, S.Z., Li, B., Rodriquez, D.C., Zhang, T., Graham, D.W., 2015. Metagenomics shows that low-energy anaerobic−aerobic treatment reactors reduce antibiotic resistance gene levels from domestic wastewater. Environmental Science & Technology 49 (4), 2577-2584. Diehl, D.L., LaPara, T.M., 2010. Effect of temperature on the fate of genes encoding tetracycline resistance and the integrase of Class 1 Integrons within anaerobic and aerobic digesters treating municipal wastewater solids. Environmental Science & Technology 44 (23), 91289133. Filee, J., Siguier, P., Chandler, M., 2007. Insertion sequence diversity in Archaea. Microbiology and Molecular Biology Reviews 71 (1), 121-157. Ghosh, S., Ramsden, S.J., LaPara, T.M., 2009. The role of anaerobic digestion in controlling the release of tetracycline resistance genes and class 1 integrons from municipal wastewater treatment plants. Applied Microbiology and Biotechnology 84 (4), 791-796. Guo, J., Ni, B.-J., Han, X., Chen, X., Bond, P., Peng, Y., Yuan, Z., 2017. Unraveling microbial structure and diversity of activated sludge in a full-scale simultaneous nitrogen and phosphorus removal plant using metagenomic sequencing. Enzyme and Microbial Technology 102, 16-25. Gupta, S.K., Padmanabhan, B.R., Diene, S.M., Lopez-Rojas, R., Kempf, M., Landraud, L., Rolain, J.M., 2014. ARG-ANNOT, a New Bioinformatic Tool To Discover Antibiotic Resistance Genes in Bacterial Genomes. Antimicrobial Agents and Chemotherapy 58 (1), 212-220. Huson, D.H., Weber, N., 2013. Microbial community analysis using MEGAN. in: E.F. DeLong (Ed.) Microbial metagenomics, metatranscriptomics, and metaproteomics, pp. 465-485. Ibrahim, M.E., Bilal, N.E., Hamid, M.E., 2012. Increased multi-drug resistant Escherichia coli from hospitals in Khartoum state, Sudan. African Health Sciences 12 (3), 368-375. Joy, S.R., Bartelt-Hunt, S.L., Snow, D.D., Gilley, J.E., Woodbury, B.L., Parker, D.B., Marx, D.B., Li, X., 2013. Fate and transport of antimicrobials and antimicrobial resistance genes in soil and runoff following land application of swine manure slurry. Environmental Science & Technology 47 (21), 12081-12088. Junjie, Q., Ruiqiang, L., Raes, J., Arumugam, M., Burgdorf, K.S., Manichanh, C., Nielsen, T., Pons, N., Levenez, F., Yamada, T., Mende, D.R., Junhua, L., Junming, X., Shaochuan, L.,

AC C

501 502 503 504 505 506 507 508 509 510 511 512 513 514 515 516 517 518 519 520 521 522 523 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552

21

ACCEPTED MANUSCRIPT

EP

TE D

M AN U

SC

RI PT

Dongfang, L., Jianjun, C., Bo, W., Huiqing, L., Huisong, Z., Yinlong, X., Tap, J., Lepage, P., Bertalan, M., Batto, J.M., Hansen, T., Paslier, D.l., Linneberg, A., Nielsen, H.B., Pelletier, E., Renault, P., Sicheritz-Ponten, T., Turner, K., Hongmei, Z., Chang, Y., Shengting, L., Min, J., Yan, Z., Yingrui, L., Xiuqing, Z., Songgang, L., Nan, Q., Huanming, Y., Jian, W., Brunak, S., Dore, J., Guarner, F., Kristiansen, K., Pedersen, O., Parkhill, J., Weissenbach, J., Bork, P., Ehrlich, S.D., Jun, W., Meta, H.I.T.C., 2010. A human gut microbial gene catalogue established by metagenomic sequencing. Nature, UK 464 (7285), 59-65. Kanehisa, M., Goto, S., Hattori, M., Aoki-Kinoshita, K.F., Itoh, M., Kawashima, S., Katayama, T., Araki, M., Hirakawa, M., 2006. From genomics to chemical genomics: new developments in KEGG. Nucleic Acids. Res. 34, D354-D357. Koike, S., Krapac, I.G., Oliver, H.D., Yannarell, A.C., Chee-Sanford, J.C., Aminov, R.I., Mackie, R.I., 2007. Monitoring and source tracking of tetracycline resistance genes in lagoons and groundwater adjacent to swine production facilities over a 3-year period. Applied and Environmental Microbiology 73 (15), 4813-4823. Kristiansson, E., Fick, J., Janzon, A., Grabic, R., Rutgersson, C., Weijdegard, B., Soderstrom, H., Larsson, D.G.J., 2011. Pyrosequencing of antibiotic-contaminated river sediments reveals high levels of resistance and gene transfer elements. Plos One 6 (2). Krzywinski, M., Schein, J., Birol, I., Connors, J., Gascoyne, R., Horsman, D., Jones, S.J., Marra, M.A., 2009. Circos: An information aesthetic for comparative genomics. Genome Res. 19 (9), 1639-1645. LaPara, T.M., Burch, T.R., McNamara, P.J., Tan, D.T., Yan, M., Eichmiller, J.J., 2011. TertiaryTreated Municipal Wastewater is a Significant Point Source of Antibiotic Resistance Genes into Duluth-Superior Harbor. Environmental Science & Technology 45 (22), 9543-9549. Lapworth, D.J., Baran, N., Stuart, M.E., Ward, R.S., 2012. Emerging organic contaminants in groundwater: A review of sources, fate and occurrence. Environmental Pollution 163, 287303. Larson, E., 2007. Community factors in the development of antibiotic resistance Annual Review of Public Health, pp. 435-447. Li, B., Yang, Y., Ma, L., Ju, F., Guo, F., Tiedje, J.M., Zhang, T., 2015. Metagenomic and network analysis reveal wide distribution and co-occurrence of environmental antibiotic resistance genes. ISME Journal 9 (11), 2490-2502. Li, B., Zhang, X., Guo, F., Wu, W., Zhang, T., 2013. Characterization of tetracycline resistant bacterial community in saline activated sludge using batch stress incubation with highthroughput sequencing analysis. Water Research 47 (13), 4207-4216. Li, D., Yang, M., Hu, J., Zhang, J., Liu, R., Gu, X., Zhang, Y., Wang, Z., 2009a. Antibiotic-resistance profile in environmental bacteria isolated from penicillin production wastewater treatment plant and the receiving river. Environmental Microbiology 11 (6), 1506-1517. Li, R., Li, Y., Kristiansen, K., Wang, J., 2008. SOAP: short oligonucleotide alignment program. Bioinformatics 24 (5), 713-714. Li, R., Yu, C., Li, Y., Lam, T.-W., Yiu, S.-M., Kristiansen, K., Wang, J., 2009b. SOAP2: an improved ultrafast tool for short read alignment. Bioinformatics 25 (15), 1966-1967. Li, R., Zhu, H., Ruan, J., Qian, W., Fang, X., Shi, Z., Li, Y., Li, S., Shan, G., Kristiansen, K., Li, S., Yang, H., Wang, J., Wang, J., 2010. De novo assembly of human genomes with massively parallel short read sequencing. Genome Res. 20 (2), 265-272. Li, X., Watanabe, N., Xiao, C., Harter, T., McCowan, B., Liu, Y., Atwill, E.R., 2014. Antibioticresistant E. coli in surface water and groundwater in dairy operations in Northern California. Environmental Monitoring and Assessment 186 (2), 1253-1260. Liu, M., Zhang, Y., Yang, M., Tian, Z., Ren, L., Zhang, S., 2012. Abundance and distribution of tetracycline resistance genes and mobile elements in an oxytetracycline production wastewater treatment system. Environmental Science & Technology 46 (14), 7551-7557.

AC C

553 554 555 556 557 558 559 560 561 562 563 564 565 566 567 568 569 570 571 572 573 574 575 576 577 578 579 580 581 582 583 584 585 586 587 588 589 590 591 592 593 594 595 596 597 598 599 600 601 602

22

ACCEPTED MANUSCRIPT

EP

TE D

M AN U

SC

RI PT

Lupo, A., Coyne, S., Berendonk, T.U., 2012. Origin and evolution of antibiotic resistance: the common mechanisms of emergence and spread in water bodies. Frontiers in Microbiology 3, 18. Mahillon, J., Chandler, M., 1998. Insertion Sequences. Microbiology and Molecular Biology Reviews 62 (3), 725-774. Mao, D., Yu, S., Rysz, M., Luo, Y., Yang, F., Li, F., Hou, J., Mu, Q., Alvarez, P.J.J., 2015. Prevalence and proliferation of antibiotic resistance genes in two municipal wastewater treatment plants. Water Research 85, 458-466. Marti, E., Variatza, E., Luis Balcazar, J., 2014. The role of aquatic ecosystems as reservoirs of antibiotic resistance. Trends in Microbiology 22 (1), 36-41. Miller, J.H., Novak, J.T., Knocke, W.R., Pruden, A., 2016. Survival of antibiotic resistant bacteria and horizontal gene transfer control antibiotic resistance gene content in anaerobic digesters. Frontiers in Microbiology 7, 263. Miller, J.R., Koren, S., Sutton, G., 2010. Assembly algorithms for next-generation sequencing data. Genomics. 95 (6), 315-327. Moura, A., Pereira, C., Henriques, I., Correia, A., 2012. Novel gene cassettes and integrons in antibiotic-resistant bacteria isolated from urban wastewaters. Research in Microbiology 163 (2), 92-100. Moura, A., Soares, M., Pereira, C., Leitão, N., Henriques, I., Correia, A., 2009. INTEGRALL: a database and search engine for integrons, integrases and gene cassettes. Bioinformatics 25 (8), 1096-1098. Mullany, P., 2014. Functional metagenomics for the investigation of antibiotic resistance. Virulence 5 (3), 443-447. Munck, C., Albertsen, M., Telke, A., Ellabaan, M., Nielsen, P.H., Sommer, M.O.A., 2015. Limited dissemination of the wastewater treatment plant core resistome. Nature Communications 6. Munir, M., Wong, K., Xagoraraki, I., 2011. Release of antibiotic resistant bacteria and genes in the effluent and biosolids of five wastewater utilities in Michigan. Water Research 45 (2), 681693. Novo, A., André, S., Viana, P., Nunes, O.C., Manaia, C.M., 2013. Antibiotic resistance, antimicrobial residues and bacterial community composition in urban wastewater. Water Research 47 (5), 1875-1887. Paiva, M.C., Reis, M.P., Costa, P.S., Dias, M.F., Bleicher, L., Scholte, L.L.S., Nardi, R.M.D., Nascimento, A.M.A., 2017. Identification of new bacteria harboring qnrS and aac(6 '-Ib/cr and mutations possibly involved in fluoroquinolone resistance in raw sewage and activated sludge samples from a full-scale WWTP. Water Research 110, 27-37. Parnanen, K., Karkman, A., Tamminen, M., Lyra, C., Hultman, J., Paulin, L., Virta, M., 2016. Evaluating the mobility potential of antibiotic resistance genes in environmental resistomes without metagenomics. Scientific Reports 6. 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. Environmental Health Perspectives 122 (3), 222-228. Port, J.A., Wallace, J.C., Griffith, W.C., Faustman, E.M., 2012. Metagenomic profiling of microbial composition and antibiotic resistance determinants in Puget Sound. Plos One 7 (10). Powell, S., Forslund, K., Szklarczyk, D., Trachana, K., Roth, A., Huerta-Cepas, J., Gabaldón, T., Rattei, T., Creevey, C., Kuhn, M., Jensen, L.J., von Mering, C., Bork, P., 2014. eggNOG v4.0: nested orthology inference across 3686 organisms. Nucleic Acids Research 42 (D1), D231D239. Pruden, A., 2014. Balancing water sustainability and public health goals in the face of growing concerns about antibiotic resistance. Environmental Science & Technology 48 (1), 5-14. Pruden, A., Arabi, M., Storteboom, H.N., 2012. Correlation between upstream human activities and riverine antibiotic resistance genes. Environmental Science & Technology 46 (21), 1154111549.

AC C

603 604 605 606 607 608 609 610 611 612 613 614 615 616 617 618 619 620 621 622 623 624 625 626 627 628 629 630 631 632 633 634 635 636 637 638 639 640 641 642 643 644 645 646 647 648 649 650 651 652 653 654

23

ACCEPTED MANUSCRIPT

EP

TE D

M AN U

SC

RI PT

Rhodes, L.D., Grayson, T.H., Alexander, S.M., Strom, M.S., 2000. Description and characterization of IS994, a putative IS3 family insertion sequence from the salmon pathogen, Renibacterium salmoninarum. Gene 244 (1–2), 97-107. Rizzo, L., Manaia, C., Merlin, C., Schwartz, T., Dagot, C., Ploy, M.C., Michael, I., Fatta-Kassinos, D., 2013. Urban wastewater treatment plants as hotspots for antibiotic resistant bacteria and genes spread into the environment: A review. Science of the Total Environment 447, 345-360. Schmieder, R., Edwards, R., 2012. Insights into antibiotic resistance through metagenomic approaches. Future Microbiology 7 (1), 73-89. Shah, S.Q.A., Cabello, F.C., L'Abée-Lund, T.M., Tomova, A., Godfrey, H.P., Buschmann, A.H., Sørum, H., 2014. Antimicrobial resistance and antimicrobial resistance genes in marine bacteria from salmon aquaculture and non-aquaculture sites. Environmental Microbiology 16 (5), 1310-1320. Shi, P., Jia, S., Zhang, X.-X., Zhang, T., Cheng, S., Li, A., 2013. Metagenomic insights into chlorination effects on microbial antibiotic resistance in drinking water. Water Research 47 (1), 111-120. Shintani, M., Urata, M., Inoue, K., Eto, K., Habe, H., Ornori, T., Yamane, H., Nojiri, H., 2007. The Sphingomonas plasmid pCAR3 is involved in complete mineralization of carbazole. Journal of Bacteriology 189 (5), 2007-2020. Siguier, P., Perochon, J., Lestrade, L., Mahillon, J., Chandler, M., 2006. ISfinder: the reference centre for bacterial insertion sequences. Nucleic Acids Research 34 (suppl 1), D32-D36. Stalder, T., Barraud, O., Casellas, M., Dagot, C., Ploy, M.-C., 2012. Integron involvement in environmental spread of antibiotic resistance. Frontiers in Microbiology 3. Sui, Q., Huang, J., Deng, S., Yu, G., Fan, Q., 2010. Occurrence and removal of pharmaceuticals, caffeine and DEET in wastewater treatment plants of Beijing, China. Water Research 44 (2), 417-426. Tatusov, R.L., Galperin, M.Y., Natale, D.A., Koonin, E.V., 2000. The COG database: a tool for genome-scale analysis of protein functions and evolution. Nucleic Acids. Res. 28 (1), 33-36. Taylor, N.G.H., Verner-Jeffreys, D.W., Baker-Austin, C., 2011. Aquatic systems: maintaining, mixing and mobilising antimicrobial resistance? Trends in Ecology & Evolution 26 (6), 278-284. Threedeach, S., Chiemchaisri, W., Watanabe, T., Chiemchaisri, C., Honda, R., Yamamoto, K., 2012. Antibiotic resistance of Escherichia coli in leachates from municipal solid waste landfills: Comparison between semi-aerobic and anaerobic operations. Bioresource Technology 113, 253-258. Torres-Cortés, G., Millán, V., Ramírez-Saad, H.C., Nisa-Martínez, R., Toro, N., Martínez-Abarca, F., 2011. Characterization of novel antibiotic resistance genes identified by functional metagenomics on soil samples. Environmental Microbiology 13 (4), 1101-1114. Udikovic-Kolic, N., Wichmann, F., Broderick, N.A., Handelsman, J., 2014. Bloom of resident antibiotic-resistant bacteria in soil following manure fertilization. Proceedings of the National Academy of Sciences of the United States of America 111 (42), 15202-15207. van der Donk, C.F.M., van de Bovenkamp, J.H.B., De Brauwer, E.I.G.B., De Mol, P., Feldhoff, K.H., Kalka-Moll, W.M., Nys, S., Thoelen, I., Trienekens, T.A.M., Stobberingh, E.E., 2012. Antimicrobial Resistance and Spread of Multi Drug Resistant Escherichia coli Isolates Collected from Nine Urology Services in the Euregion Meuse-Rhine. PLoS ONE 7 (10), e47707. van Hoek, A.H.A.M., Mevius, D., Guerra, B., Mullany, P., Roberts, A.P., Aarts, H.J.M., 2011. Acquired antibiotic resistance genes: an overview. Frontiers in Microbiology 2. Versluis, D., D'Andrea, M.M., Garcia, J.R., Leimena, M.M., Hugenholtz, F., Zhang, J., Ozturk, B., Nylund, L., Sipkema, D., van Schaik, W., de Vos, W.M., Kleerebezem, M., Smidt, H., van Passel, M.W.J., 2015. Mining microbial metatranscriptomes for expression of antibiotic resistance genes under natural conditions. Scientific Reports 5.

AC C

655 656 657 658 659 660 661 662 663 664 665 666 667 668 669 670 671 672 673 674 675 676 677 678 679 680 681 682 683 684 685 686 687 688 689 690 691 692 693 694 695 696 697 698 699 700 701 702 703 704

24

RI PT

SC

M AN U

TE D EP

735

ACCEPTED MANUSCRIPT

Wang, Z., Zhang, X.-X., Huang, K., Miao, Y., Shi, P., Liu, B., Long, C., Li, A., 2013. Metagenomic profiling of antibiotic resistance genes and mobile genetic elements in a tannery wastewater treatment plant. Plos One 8 (10). Xu, J., Xu, Y., Wang, H., Guo, C., Qiu, H., He, Y., Zhang, Y., Li, X., Meng, W., 2015. Occurrence of antibiotics and antibiotic resistance genes in a sewage treatment plant and its effluentreceiving river. Chemosphere 119 (0), 1379-1385. Yang, Y., Li, B., Zou, S., Fang, H.H.P., Zhang, T., 2014a. Fate of antibiotic resistance genes in sewage treatment plant revealed by metagenomic approach. Water Research 62 (0), 97-106. Yang, Y., Yu, K., Xia, Y., Lau, F.T.K., Tang, D.T.W., Fung, W.C., Fang, H.H.P., Zhang, T., 2014b. Metagenomic analysis of sludge from full-scale anaerobic digesters operated in municipal wastewater treatment plants. Appl. Microbiol. Biot. 98 (12), 5709-18. Zhang, J., Chen, M., Sui, Q., Tong, J., Jiang, C., Lu, X., Zhang, Y., Wei, Y., 2016. Impacts of addition of natural zeolite or a nitrification inhibitor on antibiotic resistance genes during sludge composting. Water Research 91, 339-349. Zhang, T., Li, B., 2011. Occurrence, transformation, and fate of antibiotics in municipal wastewater treatment plants. Critical Reviews in Environmental Science and Technology 41 (11), 951998. Zhang, T., Yang, Y., Pruden, A., 2015. Effect of temperature on removal of antibiotic resistance genes by anaerobic digestion of activated sludge revealed by metagenomic approach. Applied Microbiology and Biotechnology, 1-9. Zhang, X.-X., Zhang, T., Fang, H.H.P., 2009. Antibiotic resistance genes in water environment. Applied Microbiology and Biotechnology 82 (3), 397-414. Zhu, W., Lomsadze, A., Borodovsky, M., 2010. Ab initio gene identification in metagenomic sequences. Nucleic Acids. Res. 38 (12), e132. Zhu, Y.G., Johnson, T.A., Su, J.Q., Qiao, M., Guo, G.X., Stedtfeld, R.D., Hashsham, S.A., Tiedje, J.M., 2013. Diverse and abundant antibiotic resistance genes in Chinese swine farms. Proceedings of the National Academy of Sciences of the United States of America 110 (9), 3435-3440.

AC C

705 706 707 708 709 710 711 712 713 714 715 716 717 718 719 720 721 722 723 724 725 726 727 728 729 730 731 732 733 734

25

ACCEPTED MANUSCRIPT

Figure captions

737

Table 1. Summary of MGE abundances, including plasmids, conjugative transposons, ISs and

738

integrons in two sludge metagenomes (figures in bracket show percentages)

739

Figure 1. Metagenomic data relevant to antibiotic resistance contaminants were classified into four

740

categories: ARGs origin, gene transfer potential, ARG potential, and resistance pathway.

741

Figure 2. Circular representation of microbial communities in AAS and ADS samples at phylum level,

742

visualized by Circos software (Krzywinski et al., 2009). The outmost two circles list names of two

743

samples and microorganism composition. The connecting lines inside the circle link the phylum to the

744

samples and the width of the line indicates the relative abundance of different phyla in the two

745

samples. The outset cycles were coloured according to the software default setting.

746

Figure 3. Relative distribution of the predicted proteins in the functional class of defense mechanisms

747

in both the aerobic activated sludge (AAS) and anaerobically digested sludge (ADS). Metagenomic

748

data were annotated against the eggNOG database at a cutoff of E-value < 10-5.

749

Figure 4. Distributions of ARG types (in relation to antibiotic types) detected in the AAS and ADS

750

samples. The data were visualized using Circos (Krzywinski et al., 2009). MLS denotes macrolide-

751

lincosamide-streptogramin. The outmost two circles list names of two samples and each ARG type.

752

The third circle represents the reads number of ARG types. The length of the bars on the outer ring

753

represented the percentage of ARGs in two samples (left side of the diagram) and correlates the

754

percentages of respective ARG types in the AAS and ADS (right side of the diagram). The outset

755

cycles were coloured according to the software default setting.

756

Figure 5. Abundances of the ARGs subtypes in both AAS and ADS samples (the color intensity in

757

each panel shows the ARGs gene log10 number)

758

Figure 6. Network analysis about the correlation between ARG subtypes and taxonomical origin in

759

the AAS and ADS samples. (The colorful nodes stand for various ARG subtypes, while the grey and

760

white nodes indicate diverse microorganisms. The node size shows the abundance of ARGs and

761

microorganisms. A connection represents the correlation between ARG subtype and potential ARB).

AC C

EP

TE D

M AN U

SC

RI PT

736

26

ACCEPTED MANUSCRIPT Table 1. Summary of MGE abundances, including plasmids, conjugative transposons, ISs and integrons in two sludge metagenomes (figures in bracket show percentages) AAS 93,004 reads (0.38%) 295,370 reads (1.24%) 91,043 reads (0.38%) 5026 reads (0.021%)

ADS 31,842 reads (0.13%) 1,443,924 reads (6.04%) 45,780 reads (0.19%) 9119 reads (0.038%)

AC C

EP

TE D

M AN U

SC

RI PT

MGEs Plasmids Conjugative transposons Insert sequences Integrons

ACCEPTED MANUSCRIPT Cleaned metagenomic sequence reads Illumina HiSeq 2000

De novo assembly Sequences

Samples from WWTPs Gene prediction from contigs

NCBI protein database

NCBI nucleotide database

Local ARG database

Gene transfer potential

ARGs origin

ARGs potential

Plasmids

Community structure

ARGs

Transposons Insertion sequences Integrase

KEGG

Network analysis Possible ARG hosts

SC

Isfinder & INTEGRALL

M AN U

NCBI plasmid database

RI PT

Function annotation

eggNOG

ARGs pathway

Resistance pathway

Figure 1. Metagenomic data relevant to antibiotic resistance contaminants were classified into four categories: ARGs origin, gene transfer potential, ARG potential,

EP AC C

 

TE D

and resistance pathway.

EP

TE D

M AN U

SC

RI PT

ACCEPTED MANUSCRIPT

AC C

Figure 2. Circular representation of microbial communities in AAS and ADS samples at phylum level, visualized by Circos software (Krzywinski et al., 2009). The outmost two circles list names of two samples and microorganism composition. The connecting lines inside the circle link the phylum to the samples and the width of the line indicates the relative abundance of different phyla in the two samples. The outset cycles were coloured according to the software default setting.

ACCEPTED MANUSCRIPT

Protein involved in methylation Type II restriction enzyme, methylase subunits ABC transporter, ATP-binding protein Protein involved in drug transport Type I restriction-modification system protein Hydrophobe/Amphiphile efflux-1 (HAE1) protein Type I restriction-modification system methyltransferase subunit ABC-type antimicrobial peptide transport system, permease component Beta-lactamase class C and other penicillin binding proteins Restriction protein Restriction endonuclease ABC-type antimicrobial peptide transport system, ATPase component ABC-type multidrug transport system, ATPase component ABC-2 type transporter protein Type I site-specific deoxyribonuclease ABC-type multidrug transport system, ATPase and permease components Beta-Lactamase Type III restriction protein Type I site-specific restriction-modification system, R subunit and related helicases Type I restriction protein Acriflavin resistance protein Cation/multidrug efflux pump ABC transporter protein

500

1000

M AN U

0

SC

RI PT

AAS ADS

1500

2000

2500

3000

3500

Number of matched genes

Figure 3. Relative distribution of the predicted proteins in the functional class of defense mechanisms in both the aerobic activated sludge (AAS) and anaerobically digested sludge (ADS). Metagenomic data were annotated

AC C

EP

TE D

against the eggNOG database at a cutoff of E-value < 10-5.

4000

TE D

M AN U

SC

RI PT

ACCEPTED MANUSCRIPT

EP

Figure 4. Distributions of ARG types (in relation to antibiotic types) detected in the

AC C

AAS and ADS samples. The data were visualized using Circos (Krzywinski et al., 2009). MLS denotes macrolide-lincosamide-streptogramin. The outmost two circles list names of two samples and each ARG type. The third circle represents the reads number of ARG types. The bar width of the bars between ARG types and samples correlates to the percentages of respective ARG types in the AAS and ADS samples. The outset cycles were coloured according to the software default setting.

SC

RI PT

ACCEPTED MANUSCRIPT

M AN U

Figure 5. Abundances of the ARGs subtypes in both AAS and ADS samples (the color intensity in each panel shows the ARGs gene log10

AC C

EP

TE D

number)

ACCEPTED MANUSCRIPT

beta-lactam

fluoroquinolones

glycopeptides

MLS

rifampicin

tetracyclines

trimethoprim

Genus

Others

phenicols

(A) AAS

M AN U

SC

RI PT

aminoglycosides

(B) ADS

TE D

Figure 6. Network analysis about the correlation between ARG subtypes and taxonomical origin in the AAS and ADS samples. (The colorful nodes stand for various ARG subtypes, while the grey and white nodes indicate diverse microorganisms. The node size shows the abundance of ARGs and

AC C

EP

microorganisms. A connection represents the correlation between ARG subtype and potential ARB).

ACCEPTED MANUSCRIPT Highlights •

Metagenomic sequencing was used to investigate profiles of ARGs and MGEs. Wastewater treatment plants might be hotspots of ARGs and MGEs.



Some environmental bacteria might be potential hosts of multiple ARGs.

AC C

EP

TE D

M AN U

SC

RI PT