An 18S rRNA Workflow for Characterizing Protists in Sewage, with a Focus on Zoonotic Trichomonads

Abstract

Microbial eukaryotes (protists) are important components of terrestrial and aquatic environments, as well as animal and human microbiomes. Their relationships with metazoa range from mutualistic to parasitic and zoonotic (i.e., transmissible between humans and animals). Despite their ecological importance, our knowledge of protists in urban environments lags behind that of bacteria, largely due to a lack of experimentally validated high-throughput protocols that produce accurate estimates of protist diversity while minimizing non-protist DNA representation. We optimized protocols for detecting zoonotic protists in raw sewage samples, with a focus on trichomonad taxa. First, we investigated the utility of two commonly used variable regions of the 18S rRNA marker gene, V4 and V9, by amplifying and Sanger sequencing 23 different eukaryotic species, including 16 protist species such as Cryptosporidium parvum, Giardia intestinalis, Toxoplasma gondii, and species of trichomonad. Next, we optimized wet-lab methods for sample processing and Illumina sequencing of both regions from raw sewage collected from a private apartment building in New York City. Our results show that both regions are effective at identifying several zoonotic protists that may be present in sewage. A combination of small extractions (1 mL volumes) performed on the same day as sample collection, and the incorporation of a vertebrate blocking primer, is ideal to detect protist taxa of interest and combat the effects of metazoan DNA. We expect that the robust, standardized methods presented in our workflow will be applicable to investigations of protists in other environmental samples, and will help facilitate large-scale investigations of protistan diversity.

Keywords: 18S rRNA amplicon sequencing; Environmental sequencing; Protist; Sewage; Trichomonad; Zoonoses.

Fig. 1

Workflow employed for testing and optimization of protocols to detect microbial eukaryotes in raw sewage. Part 1 was used to evaluate two variable regions of the 18S rRNA gene, V4 and V9, and determine best practices for zoonotic protists. Part 2 was used to optimize the methods in Part 1 and develop experimental methods for Illumina sequencing of raw sewage samples. Methods that were tested but determined suboptimal are indicated in boxes with dotted lines and/or crossed out text

Fig. 2

Phylogenetic tree of trichomonad Sanger sequences from a V4 region and b V9 region. Numbers at nodes correspond to Bayesian posterior probabilities. GenBank accession numbers are shown at the left for each taxon. Colored taxa are identical to each other and brackets indicate OTU membership and corresponding ID number from the supplementary material after clustering at 98% identity

Fig. 3

Alpha diversity analysis for sewage 18S communities using the phylogenetic diversity metric and Shannon Index. All analyses were calculated from replicate OTU tables sampled at a depth of 600,000 sequences. a Alpha diversity of the different extraction volumes for the V4 and V9 regions. Values shown represent the measurements from all 1 and 10 mL (fresh samples only) extraction volumes by region. Error bars represent ± one standard deviation; double asterisks indicates p value <0.001. b Alpha diversity of sewage samples based on extraction volume and blocking primer use for the V9 region. Asterisk indicates p value <0.01

Fig. 4

Relative sequence abundances of eukaryotic taxa across experimental conditions for each region. Each bar represents the average of two replicates per condition per time point. All data shown were from fresh extractions. a V4 region results from the 1 and 10 mL extractions. b V9 region results from the 1 and 10 mL extractions (with (B) and without (N) the blocking primer)

Fig. 5

LEfSe analysis of 18S communities to determine biomarker taxa across sample groupings for each region. All data shown were from fresh extractions. LDA scores for taxa identified as differentially abundant between a 10 and 1 mL for the V4 region and b 1 mL, 1 mL blocking, 10 mL, and 10 mL blocking samples for the V9 region. Taxa are ranked according to their effect size (LDA score) and associated with the group with the highest median