Comparison of Methods to Identify Pathogens and Associated Virulence Functional Genes in Biosolids from Two Different Wastewater Treatment Facilities in Canada

The use of treated municipal wastewater residues (biosolids) as fertilizers is an attractive, inexpensive option for growers and farmers. Various regulatory bodies typically employ indicator organisms (fecal coliforms, E. coli and Salmonella) to assess the adequacy and efficiency of the wastewater treatment process in reducing pathogen loads in the final product. Molecular detection approaches can offer some advantages over culture-based methods as they can simultaneously detect a wider microbial species range, including non-cultivable microorganisms. However, they cannot directly assess the viability of the pathogens. Here, we used bacterial enumeration methods together with molecular methods including qPCR, 16S rRNA and cpn60 gene amplicon sequencing and shotgun metagenomic sequencing to compare pre- and post-treatment biosolids from two Canadian wastewater treatment plants (WWTPs). Our results show that an anaerobic digestion WWTP was unsuccessful at reducing the live indicator organism load (coliforms, generic E. coli and Salmonella) below acceptable regulatory criteria, while biosolids from a dewatering/pelletization WWTP met these criteria. DNA from other pathogens was detected by the molecular methods, but these species were considered less abundant. Clostridium DNA increased significantly following anaerobic digestion treatments. In addition to pathogen DNA, genes related to virulence and antibiotic resistance were identified in treated biosolids. Shotgun metagenomics revealed the widest range of pathogen DNA and, among the approaches used here, was the only approach that could access functional gene information in treated biosolids. Overall, our results highlight the potential usefulness of amplicon sequencing and shotgun metagenomics as complementary screening methods that could be used in parallel with culture-based methods, although more detailed comparisons across a wider range of sites would be needed.


Introduction
Large wastewater treatment facilities start with clarification and end with disinfection of the liquid portion before discharging it into a nearby watercourse. The remaining non-liquid portion, sewage sludge, can undergo different biological as well as physical-chemical treatment processes by means of anaerobic or aerobic digestion, dewatering or pelletization [1]. Municipal biosolids, as defined by the Canadian Council of Ministers of the Environment (CCME), are organic-based products which may be solid, semi-solid or liquid and which are produced from the treatment of municipal sludge. Municipal biosolids are municipal sludge which has been treated to meet to jurisdictional standards, requirements or guidelines including the reduction of pathogens. It is estimated that 0.4 to 8 million tons of municipal biosolids are produced annually in Canada, USA and Europe [2][3][4]. A substantial amount of these biosolids are formulated into fertilizer for land application as a means of waste management [4,5].
The recycling of organic wastes for land application as fertilizers and supplements (e.g., soil amendments) can result in benefits through the suppression of plant diseases [6], return and cycling of nutrients to the soil [7], and improvement of the physical properties of the soil (e.g. moisture absorbance) by increasing the overall organic matter content [8]. In contrast, there may also be risks associated with adding biosolids to soil, since these materials can be a potential source of pathogens, endotoxins and chemicals from industrial and household sources, which could lead to adverse environmental and human health effects [9][10][11]. As such, the benefits must be carefully balanced against the potential safety hazards associated with these materials. Consideration of the sources of waste-derived materials and the level of processing and treatment used during their manufacture are essential in determining the risks, since concerns over plant, animal, and human pathogens can be effectively alleviated with adequate treatment. Although very little is known of public health issues directly linked to pathogens in biosolids [10,12,13], direct contact or contamination of food crops represent two plausible routes whereby pathogens, if present in significant amounts, could affect human health. Pathogens of concern that may be present in sewage include: bacteria (e.g. Salmonella spp, Escherichia coli pathogenic strains, Campylobacter jejuni), viruses (e.g. Adenovirus, Rotavirus, Hepatitis A), protozoa (e.g. Cryptosporidium sp., Entamoeba histolytica, Giardia lamblia), and helminths (e.g. Ascaris lumbricoides, Ascaris suum, Trichuris trichiura) [14][15][16].
Pathogen inactivation is a key goal in biosolids production. Previous studies have shown that pathogens that have survived sewage treatment processes end up in biosolid-amended soils [9,11]. Additionally, protozoan parasites including Cryptosporidium sp. and Giardia sp. were also reported to survive wastewater treatment processes [17]. Various regulatory bodies both domestically (at the provincial, territorial and federal level) as well as internationally, typically employ indicator organisms (fecal coliforms, E. coli and Salmonella) to assess the adequacy and efficiency of the treatment process in reducing pathogen loads in the final product. For example, according to the federal Canadian Food Inspection Agency (CFIA) Salmonella must be absent (non-detectable) and fecal coliform levels must not exceed 1000 MPN/g of dry weight in biosolids that are sold or imported into Canada [18]. These regulated levels vary between Canadian provinces and other countries and sometimes depend on the intended use of biosolids as fertilizers (e.g., food vs. non-food crops). These microbial indicators do not represent a comprehensive list of pathogens found in biosolids, but are used as indicators of treatment efficiency regarding pathogen inactivation. Since culture-based methods are used to enumerate these bacteria, they fail to provide information on non-indicator pathogens as well as viable but non-culturable (VBNC) organisms. Although qPCR-based quantification could circumvent some of these limitations, the large list of potential pathogens would render it a highly laborious and costly process. As such, a more holistic approach is needed to better characterize the pathogen population/load in treated biosolids intended for field application. In the present study, our goals were two-fold: 1) observe the effectiveness of different wastewater treatment processes through changes in the microbial taxonomical and functional community composition in the end product by a genomic approach and 2) compare traditional pathogen detection methods to modern molecular detection methods. To achieve the latter, bacterial enumeration (most probable number, MPN) methods were compared to indirect molecular detection techniques including qPCR, 16S rRNA and cpn60 gene amplicon pyrosequencing and shotgun metagenomic sequencing in their ability to detect pathogens and virulence genes in biosolids obtained from two different WWTP. The WWTP used different treatments, namely anaerobic digestion and dewatering-pelletization, and samples were taken before and after treatment at various time points over the course of one year.

Study sites and sample collection and characterization
Two Canadian biosolid treatment facilities named A and C by the authors were sampled at three time intervals in one year. The owners of the sites gave their consent to carry out this study on these sites. The Plant A treats waste activated sludge by anaerobic digestion and dewatering process with end product of wet pellets. Plant C treats waste activated sludge in a dewatering/pelletization process involving belt-filter press and a final process of pelletization by thermal drier (250-450°C at the entry, and 80-130°C at the exit). Samples were taken prior to and just after treatment in triplicate resulting in a total of 36 samples. The samples were transported on ice and were stored at 4°C (culture methods) or -20°C (DNA extraction) immediately after receiving until further use. The biosolid samples were labeled using the facility letter (A or C), treatment, and sampling date ( Table 1). The moisture content (moisture %) of the samples was determined using an electronic moisture analyzer (IR-35 Moisture Analyzer,

Culture-based methods
Fecal coliforms and E. coli in each sample were evaluated using a most probable number assay (MPN) according to method MFHPB-19 [19]. Briefly, 90 ml of peptone water was added to 10g of sample (dry pellets were crushed before the addition of peptone water) followed by homogenization. Ten-fold serial dilutions of the suspension were made using peptone water. One ml aliquots of each dilution were inoculated into five tubes of lauryl sulfate tryptose (LST) broth. The production of gas in LST broth indicated a presumptive positive test for coliforms.  [20]. The original qualitative method was modified to be used as an MPN procedure. Fifty ml of nutrient broth (NB) was added to 25 g of sample and incubated for 1 h at 35°C followed by addition of 175 ml of NB. Ten-fold serial dilutions of the homogenates were made using NB and the tubes (5 for each dilution) were incubated for 18-24 h at 35°C. A portion of each broth culture was inoculated onto Modified Semi-solid Rappaport Vassiliadis (MSRV) agar, and the plates incubated for 72 hours at 42°C. Presumptive Salmonella positive MSRV cultures were subcultured onto MacConkey agar. The selected colonies from MacConkey agar were initially tested for Salmonella using confirmatory test media including the agar slants of triple sugar iron (TSI) and lysine iron agar (LIA) and urea. The suspect isolates were purified using xylose lysine tergitol-4 (XLT-4) or xylose lysine deoxycholate (XLD) agar, and confirmed using biochemical (API 20E) and serological tests including agglutination (Oxoid Salmonella latex test, Oxoid) and an enzyme linked immune-sorbent assay (ELISA) [21]. Salmonella Typhimurium ATCC 14028 or S. Berta ATCC 8392 was used as a positive control, and E. coli ATCC 11775 or Enterobacter cloacae ATCC 1307 used as a negative control. The MPN of Salmonella per g of dry sample was calculated based on dry weight.

DNA extraction
DNA was extracted using the PowerMax 1 Soil DNA Isolation Kit following the manufacturer's instructions (MoBio, Carlsbad, CA). The final purified DNA was then used for PCR amplification or stored at -20°C.

16S rRNA and cpn60 gene sequencing
In order to taxonomically identify microbes potentially present in biosolids and observe shifts in community composition, amplicon sequencing of two distinct marker genes was carried out. For the 16S rRNA gene, eight libraries were prepared using two different universal bacterial primer sets on four different samples ("Plant A August 2009 before-after" and "Plant C May 2009 before-after"). The following primers were used: V1-V3: forward 5' CGTATCGCCTCCCTCGCGCCATCAGACGAGTGCGTAGTTTGATCCTGGCTCAG-3', reverse 5'-CTATGCGCCTTGCCAGCCCGCTCAGACGCTCGACACATTACCGCGGCTGCTGG-3' and V3-5: forward 5'-CGTATCGCCTCCCTCGCGCCATCAGAGACGCACTCGCCTACGGGAGG CAGCAG-3' reverse 5'-CTATGCGCCTTGCCAGCCCGCTCAGAGCACTGTAGCCGT CAATTCMTTTRAGT-3' where the italic sequence represents the sample specific multiplex identifier, the bold sequence represents the template specific sequences and the remaining sequence is the 454 adapter A (forward) and adapter B (reverse).
For the cpn60 gene, four libraries were prepared using the "Plant A August 2009 beforeafter" and "Plant C May 2009 before-after" samples. The following universal primers were used: for "Plant A August 2009 before" and "Plant C May 2009 before", forward: 5'-CGTATCGCCTCCCTCGCGCCATCAGATCAGACACGGCIGGIGAYGGNACNACNAC3', reverse: The same primers were used for "Plant A August 2009 after" and "Plant C May 2009 after" except that the multiplex identifiers (in italics) were replaced by the following: forward primer: CGTGTCTCTA and reverse primer: CTCGCGTGTC. The PCR conditions were as follows: the mixture in a 50ul final volume contained 50ng DNA, 25pmol of each primer for 16S or cpn60, 1X final Taq polymerase buffer, and 2.5 units of Taq polymerase (New England BioLabs Ltd, Pickering, ON, Canada) and 1 μl of 10 mM deoxynucleoside triphosphates. PCR cycling consisted of 94°C for 5 minutes, denaturation at 94°C for 30 seconds, annealing at 56°C for 30 seconds for 16S V3-1, 50°C for 16S V5-3 and 55°C for cpn60, with an extension at 72°C for 45 seconds after 35 cycles. A final extension at 72°C was added for 7 minutes.

Real-time quantitative PCR
The abundance of E. coli was quantified using two distinct sets of primers targeting the beta-Dglucuronidase gene (uidA): uidA1-F: CAGCAATTGCCCGGCTTTCTTGTA, uidA1-R: GGCATT CAGTCTGGATCGCGAAA (generating a fragment of 83bp) and uidA2-F: GTATCGGTGT GAGCGTCGCAG and uidA2-R: GCGTGGTGATGTGGAGTATTGCC (generating a fragment of 154 bp). For Salmonella quantification, two sets of primers targeting the invasion gene A (invA and sal) were used: invA-F: GATTCTGGTACTAATGGTGATGATC, invA-R: GCCAGGCTAT CGCCAATAAC (generating a fragment of 287 bp) and sal-F: GCGTTCTGAACCTTTGGTAATAA, and sal-R: CGTTCGGGCAATTCGTTA (generating a fragment of 102bp) [22]. Real-time quantitative PCR (qPCR) amplification was performed using a Rotor Gene 3000 instrument (Corbett Research, Mortlake, NSW, Australia) using a QuantiTect SYBR Green PCR master mix (Qiagen) in a 20 μl volume containing 10 pmol of each primer and a final MgCl 2 concentration of 2.5mM for uidA1 and 3.5mM for uidA2, invA and sal. The amplification conditions were as follows: 95°C for 15 minutes, followed by 40 cycles of 95°C for 10 seconds, 55°C for 15 seconds and 72°C for 15 seconds. Fluorescence was measured at the end of each cycle at 72°C and a melting curve analysis (65-95°C) was performed at the end of the amplification procedure.
Standard curves were generated from PCR fragments using the above mentioned primer sets and genomic DNA from E. coli K12 and Salmonella enterica serovar Typhimurium. Amplicons from the different primer sets were cloned using the pGEM-T easy Vector System (Promega U.S. Madison, WI) and transformed into E. coli (JM109). Recombinant plasmids were isolated and linearized with ScaI (New England BioLabs) quantified by PicoGreen (Invitrogen) and used to generate standard curves with serial dilutions. Standard curve efficiencies were all between 0.95 and 1.0 with R 2 value between 0.99 and 1.00.

Sequence data analysis
16S sequence data were primarily analyzed through the RDP pyrosequencing pipeline (http:// pyro.cme.msu.edu/). The sequences were deconvoluted and binned according to their multiplex identifier (only accepting perfect matches), and the multiplex identifier and the forward primer were trimmed using the 'Pipeline Initial Process' tool. This resulted in eight distinct data sets. Using the 'Pipeline Initial Process' tool, all sequences that contained undetermined bases (N) or were shorter than 150 bp were removed. The data sets were submitted to the RDP Classifier tool using a bootstrap cutoff of 80%. The datasets from the two different 16S primers were then compared and since no large variation were observed (data not shown), the datasets were pooled for all downstream analyses. The cpn60 sequences were subjected to the same initial process using the tools available in the RDP pyrosequencing pipeline. The four sequence datasets were then compared to the cpnDB (http://cpndb.cbr.nrc.ca/) using blastn [23]. Blast results were analyzed using MEGAN to place the sequences in the NCBI taxonomy using a lowest common ancestor algorithm [24].
For shotgun metagenomic datasets, replicate sequences that resulted from the attachment of DNA to beads during emulsion PCR and were not derived independently from the environmental data were removed from the dataset using the method of Gomez-Alvarez [25]. Sequences were then submitted to MG-RAST v. 2.0 for automated annotation [26]. The abundance of different taxa given by MG-RAST was further normalized to the individual genome size by dividing the number of hits by the individual genome sizes. This normalization is necessary due to variability in genome size between different organisms as larger genomes generate more reads even though the organism is not more abundant in the sample. Data mining efforts were specifically focused on known pathogens and virulence genes. Results are presented as the raw number of sequences related to a particular function or taxon, or the relative abundance of taxa calculated from the number of reads mapping to this taxon divided by the total number of reads in the dataset.

Statistical analyses
Two proportion z-tests were performed according to Wang [27]. Principal coordinate analysis (PCoA) were performed based on Bray-Curtis distance calculated from genus relative abundance in R [28] using the vegan package [29]. Spearman rank-order correlations were performed in R.

Physical characteristics of biosolid samples
A total of 36 samples were taken at three different dates from two wastewater treatment plants. The two different plants used different methods to treat the biosolids: Plant A used anaerobic digestion, which lowered the water content in the biosolids from 95% to approximately 70%, while Plant C used a dewatering/pelletization treatment which reduced the biosolids water content from 63-97% to approximately 10% (Table 1).

Taxonomical shifts following biosolids treatment
One pre/post sample pair from Plant A (August 4, 2009) and another pair from Plant C (May 12, 2009) were selected for amplicon sequencing (16S and cpn60). For both plants, large changes were observed at the phylum/class level following treatment independent of the primer pair used (Fig 1a). For Plant A, Proteobacteria dominated the microbial community before treatment and was replaced by Firmicutes, Bacteroidetes and Chloroflexi following treatment (Fig 1a). For Plant C, the shifts were less drastic, with decreases in the relative abundance of Proteobacteria and increases in Actinobacteria, Bacteroidetes, and Firmicutes following treatment (Fig 1a). When looking at lower taxonomical levels using Unifrac analyses, a similar pattern emerged with the anaerobic digestion from Plant A causing stronger shifts in microbial communities (dots more distant) and the dewatering/pelletization from Plant C causing less dramatic shifts (Fig 1b). Similar results were obtained from cpn60 amplicon analyses for Plant A (Fig 2) and Plant C (not shown). Shifts in the dominant genera were also observed for the two plants following treatments using 16S rRNA and cpn60 gene sequencing. For Plant A, the community shifted from an Acidovorax and Novosphingobium dominated community ( Table 2) to one dominated by anaerobes and syntrophic bacteria like Syntrophus, Sedimentibacter, Prevotella, Clostridium and Thermovirga (Table 3). For Plant C, the shifts were less evident, with a community that was generally dominated by Paludibacter and Microbacterium both before and after the biosolid treatment (Table 4).
For potential pathogens, we also looked at the species level using the metagenomics datasets from Plant A. Anaerobic digestion was successful at reducing the relative abundance of several species, including E. coli, Legionella and Pseudomonas species (Table 5). For E. coli, this relative decrease was also seen for most samples examined by qPCR (uidA1 shown in Table 6, with uidA2 showing identical trends) and MPN (Table 6). MPN analyses also highlighted a general decrease in Salmonella abundance following treatments (Table 6), while Salmonella qPCR (invA and Sal) was below the qPCR detection limit for all samples. However, anaerobic digestion increased the relative abundance of Campylobacter, Chlamydia, Clostridium, Enterococcus, Listeria and Staphylococcus species DNA in the metagenomic datasets (Table 5).

Functional shifts following biosolids treatment
Within the metagenomic datasets (Plant A, August 4, 2009), we focused our attention on functional genes related to pathogenicity and virulence ("Virulence" Subsystem hierarchy 1 category in MG-RAST). Although all subsystems were still detectable after the treatment, several of the "Virulence" subsystems decreased significantly, especially in the "Resistance to antibiotics and toxic compounds" Subsystem level 2 category (Table 7). In contrast, some subsystems were significantly more abundant following treatment. Some relatively abundant (more than 500 hits) subsystems like "Multidrug Resistance Efflux Pumps" and "Resistance to fluoroquinolones" were relatively more abundant in the biosolids following treatment ( Table 7). The total relative abundance of "Virulence" related reads decreased following treatment, from 7.4% to 5.5% of total classified reads.

Method comparison
The detection of selected pathogens was compared for the different methods used for the samples taken on August 4, 2009 in Plant A (Table 8). E. coli was detected using MPN, qPCR and metagenomics sequencing, but not by amplicon sequencing, while Salmonella was detected by MPN and metagenomic sequencing ( Table 8). Most of the other genera containing potential pathogens were only detected by metagenomic sequencing, with some exceptions like Clostridium that was also detected by amplicon sequencing in pre-and post-treatment samples ( Table 8).
All quantification methods revealed low abundance of coliforms, E. coli and Salmonella, being often below the detection limit, especially in the case of Salmonella (Table 6). The quantification of E. coli by MPN and qPCR methods was compared by Spearman correlation analysis, while the tests were not performed for Salmonella as the qPCR results were below detection limits in all cases. The abundance of E. coli measured by qPCR and MPN were not significantly correlated (uidA1 vs. MPN: r s = 0.165, P = 0.261; uidA2 vs. MPN: r s = 0.218, P = 0.137). The two qPCR quantification methods used (uidA1 and uidA2) were significantly correlated to each other (r s = 0.839, P<0.0001).
The community composition patterns were compared between the datasets from metagenomic and amplicon sequencing (two 16S primer pairs and one cpn60 primer pair). The four different methods gave relatively similar patterns at the phylum/class level, especially when comparing the pre-treatment vs. the post-treatment samples (Fig 2a). When looking at the genus relative abundance using principal coordinate analysis (PCoA) of Bray-Curtis distances, the main difference observed was between pre-and post-treatment samples, with pre-and post-treatment samples being separated on the first axis of the ordination for both amplicon (16S and cpn60) and metagenomic datasets (Fig 2b). The samples analysed using metagenomic sequencing were less differentiated than the other samples (Fig 2b). For the 16S rRNA gene amplicon sequencing, both primer pairs clustered tightly together (Fig 2b). Similarly, when also including samples from Plant C, using different primer pairs resulted in very little difference (Fig 1b). When comparing the most abundant genera, the different sequencing methods generally resulted in similar results, even though the order and relative abundance of the dominant genera changed (Tables 2, 3 and 4). In most cases, the metagenomic analysis resulted in the largest differences (Tables 2 and 3). One interesting difference is that the metagenomic analysis resulted in a more even distribution of the different genera, with relative abundances never exceeding 5.1% as compared to 33% for amplicon sequencing (Tables 2 and 3).

Pathogens and virulence genes in treated biosolids
After treatment, the biosolids showed a drastic reduction in E. coli and Salmonella content when looking at the enumeration of classic pathogens on growth media. Although the starting communities differed slightly, the different sewage treatments did affect the communities differently. From the results of the classic enumeration methods, the anaerobic digestion method of Plant A appeared less efficient in reducing the pathogen load in biosolids as compared to the dewatering/pelletization method of Plant C. Plant A exceeded the Canadian Food Inspection Agency (CFIA) criteria for Salmonella levels which must be non-detectable (at 3 out of 3 sampling dates, post-treatment biosolids had detectable Salmonella) and the level of fecal coliforms, which must not exceed 1000 MPN/g of the total dry weight (at 3 out of 3 sampling dates, treated biosolids exceeded that value) [18]. By contrast, treated biosolids from Plant C met the CFIA indicator standards at all sampling dates. The dewatering/pelletization treatment also resulted in a larger and more consistent decrease in extractable DNA than the anaerobic digestion treatment. However, molecular methods highlighted stronger shifts in community composition following anaerobic digestion as compared to dewatering/pelletization. These larger microbial community composition shifts at Plant A as compared to Plant C might have been caused by the differences in the initial communities between the two plants, but, alternatively, it might have been a direct cause of the anaerobic environment itself, as many of the dominant groups of bacteria after anaerobic digestion were from known anaerobic microorganisms within the Firmicutes and Bacteroidetes phyla.
Molecular methods identified DNA originating from many genera, containing pathogenic species and virulence genes, that increased in their relative abundance following biosolid treatments. For instance, sewage treatment at Plant A increased the relative abundance of Clostridium DNA in biosolids whereas Clostridium DNA was also detected before and after the dewatering/pelletization treatment at Plant C. Even though many species of Clostridium are non-pathogenic, some are the causative agents of diseases in humans, including botulism, tetanus and enterocolitis. As such, their persistence in biosolids could be a cause for concern if these organisms were present in a viable state and at concentrations that could cause disease. Consistent with the work presented here, DNA related to Clostridium was detectable by 16S rRNA gene pyrosequencing after biosolids treatment [30]. Clostridia are obligate anaerobes which produce endospores, therefore their increase following anaerobic digestion at Plant A is not surprising since endospores can resist this type of treatment. The detection of DNA from potential pathogenic organisms other than Salmonella and fecal coliforms in post-treated samples emphasize the need for further research in order to determine the validity and applicability Detection of Pathogens and Virulence Genes in Biosolids of indicator organisms currently used for regulatory purposes. The shotgun metagenomic approach also allowed the detection of numerous virulence-related genes in the processed biosolids of Plant A, including mobile genetic elements and antibiotic resistance genes, as previously reported in treated biosolids [31].

Emerging methods to detect pathogens in biosolids
Culture-based methods, although the simplest and most inexpensive way to detect live pathogens, cannot detect viable but non-culturable (VBNC) bacteria which can potentially reanimate in anaerobically digested biosolids. [32,33]. In contrast, at the molecular level, qPCR can detect non-culturable bacteria and could be used to detect the presence of a specific pathogen's gene in biosolids [34,35]. However, the wide range of possible pathogenic targets again renders such Detection of Pathogens and Virulence Genes in Biosolids a method rather cumbersome. In our study, qPCR was not able to detect Salmonella in many samples, while it was detected by culture-based methods, cpn60 pyrosequencing and metagenomic sequencing. This was surprising since the presence of VBNC bacteria, dead cells and naked DNA should have made the qPCR method detect more of this pathogen [34,36]. Previous studies had already reported a lower sensitivity for qPCR than for culture methods [37], mainly related to the starting sample size. Further studies using spiked samples would be necessary to identify which method was the most sensitive, specific and reproducible.
In this study, we compared amplicon sequencing of cpn60 and two regions of the 16S rRNA gene [38] and shotgun metagenomics [39][40][41][42] as potential alternatives to culture-based approaches for the detection of genetic material from pathogens in biosolids. The region of the 16S rRNA gene that was sequenced did not have a strong influence on the community composition at the genus and phylum levels. The main advantage of shotgun metagenomics, apart from avoiding any possible amplification bias, is that detection is not limited to the targeted  ND: not detected, below qPCR detection limit. doi:10.1371/journal.pone.0153554.t008 Detection of Pathogens and Virulence Genes in Biosolids organisms. Culture methods and qPCR can only detect the DNA of organisms that are targeted [34,35,43], while metagenomics and amplicon sequencing can detect DNA from all organisms present in a sample (including eukaryotes and viruses). For instance viruses, which can be an indicator for biosolid treatment efficiency, can represent up to 10-14% of total sequences in metagenomic datasets [44] but were not detected using other methods. Another striking example is that, even with a general decrease of E. coli, several other potential human pathogen DNAs like those from Clostridium increased significantly as previously reported using 16S rRNA pyrosequencing [30]. However, standard protocols for amplicon and metagenomic sequencing are only semi-quantitative (relative abundance), reducing the utility of the data. Using internal spiked standards like previously done in metatranscriptomics [45,46] could help solve this issue. In the present study, if pathogen monitoring methods were constrained only for E. coli detection using culture-based methods, our data would have falsely indicated that the biosolids treatments were highly efficient in reducing pathogen loads in treated biosolids. Another major advantage of shotgun metagenomic sequencing is the ability to detect the whole complement of functional genes in an environmental sample; however, it is closely linked to the quality of the databases used. As the number of annotated species increases and are deposited into relevant databases, the similarity search against annotated genes will be more successful. It is especially true in the case of biosolids metagenomics, as human pathogens are one of the most sequenced and accurately annotated type of microorganism in current databases. Another critical aspect of current metagenomic studies is the read length used which 1) often precludes the definition of the context surrounding the genes detected, 2) cannot link the reads to an organism or a function with a very high level of certainty. The increased throughput of long-read sequencers could probably solve some of these issues in the near future. Our data suggests that shotgun metagenomic sequencing could be an excellent supplementary method to culture-based methods for pathogen detection, which is the primary method used by the Canadian Food Inspection Agency to detect pathogens. Metagenomics has not only detected all the specific DNA of the organisms detected by other methods, but it had also detected scores of potential pathogens and virulence related genes that the other methods were unable to detect. However, this approach incurs significant costs and requires advanced analytical expertise. Until the cost of library preparation and sequencing decreases further and analyses become routine, our data shows that other powerful approaches like amplicon sequencing will work with the appropriate depth of sequencing. DNA-based methods, as used in this study, can indirectly determine the presence of pathogens through their DNA and virulence-associated factors, but cannot directly determine live pathogen counts. In future studies, it will be important to assess the presence of live pathogenic microbial cells, particularly viable but non-culturable organisms, by using more quantitative molecular approaches such as propidium monoazide (PMA) treatment [47,48] in order to discriminate between live and dead cells.