16S rRNA Amplicon Sequencing for Epidemiological Surveys of Bacteria in Wildlife

Raw sequencing results.The sequencing of 1,569 PCR products in two MiSeq runs generated a total of 23,698,561 raw paired-end sequence reads (251 bp) of the V4 region of the 16S rRNA gene. Because we made PCR replicates for each rodent sample, and because we included several controls in each sequencing run, we had more PCR products (N = 1,569) than rodent samples (N = 711) (see summary in Table S1 in the supplemental material and complete information by sample and run in Table S2). Overall, 99% of the PCRs generated more than 3,000 raw reads (mean, 11,908 reads; standard deviation, 6,062 reads). The raw sequence files have been deposited in FASTQ format in the Dryad Digital Repository (http://dx.doi.org/10.5061/dryad.m3p7d) (41).

Using mothur v1.34 (42) and the MiSeq standard operating procedure (http://www.mothur.org/wiki/MiSeq_SOP), we removed 20.1% of the paired-end reads because they were misassembled, 1.5% of sequences because they were misaligned, 2.6% because they were chimeric, and 0.2% because they were nonbacterial. The remaining reads were grouped into operational taxonomic units (OTUs) with a divergence threshold of 3%. Bioinformatics analysis identified 13,296 OTUs, corresponding to totals of 7,960,533 sequences in run 1 and 6,687,060 sequences in run 2.

Taxonomic assignment of sequences.We used the Bayesian classifier (bootstrap cutoff = 80%) implemented in mothur with Silva SSU Ref database v119 (42) as a reference for the taxonomic assignment of OTUs. The 50 most abundant OTUs accounted for 89% (minimum, 15,284 sequences; maximum, 2,206,731 sequences) of the total sequence data set (see Table S3 in the supplemental material). The accuracy of taxonomic assignment (to the genus level) was assessed with positive controls for PCR, corresponding to DNA extracts from laboratory isolates of Bartonella taylorii, Borrelia burgdorferi, and Mycoplasma mycoides (PC_{Bartonella_t}, PC_{Borrelia_b}, and PC_{Mycoplasma_m}, respectively), which were correctly assigned to a single OTU corresponding to the appropriate reference sequences (Table 2). Note that the sequences of PC_{Mycoplasma_m} were assigned to Entomoplasmataceae rather than Mycoplasmataceae because of a frequent taxonomic error reflected in most databases, including Silva (43). This problem might also affect other taxa. We therefore recommend systematically carrying out a blast analysis against the sequences of taxa of interest in GenBank to confirm the taxonomic assignment obtained with the 16S databases. Finally, we assumed that the small number of sequences per sample might limit the completeness of bacterial detection (36). For this reason, we discarded seven rodent samples (2 M. erythroleucus and 5 M. domesticus) yielding fewer than 500 sequences for at least one of the two PCR replicates (1% of the samples).

Filtering for reagent contaminants.16S rRNA amplicon sequencing data may be affected by the contamination of reagents (23). We therefore filtered the data, using negative controls for extraction (NC_ext), corresponding to extraction without the addition of a tissue sample, and negative controls for PCR (NC_PCR), corresponding to PCR mixtures to which no DNA was added. We observed between 2,843 and 8,967 sequences in the NC_ext and between 5,100 and 9,145 sequences in the NC_PCR. On the basis of their high number of reads in negative controls, we identified 13 contaminant genera, including Pseudomonas, Acinetobacter, Herbaspirillum, Streptococcus, Pelomonas, Brevibacterium, Brachybacterium, Dietzia, Brevundimonas, Delftia, Comamonas, Corynebacterium, and Geodermatophilus, some of them having been previously identified in other studies (23). These contaminants accounted for 29% of the sequences in the data set (Fig. 2). They also differed between MiSeq runs: Pseudomonas, Pelomonas, and Herbaspirillum predominated in run 1, whereas Brevibacterium, Brachybacterium, and Dietzia predominated in run 2 (see Table S4 and Fig. S1 in the supplemental material). This difference probably reflects the use of two different PCR kits manufactured several months apart (Qiagen technical service, personal communication). The majority of the other contaminants, such as Streptococcus, most likely originated from the DNA extraction kits used, as they were detected in abundance in the negative controls for extraction (NC_ext).

• Open in new tab
• Download powerpoint

FIG 2

Taxonomic assignment of the V4 16S rRNA sequences in wild rodents and in negative controls for extraction and PCR. The histograms show the percentages of sequences for the most abundant bacterial genera in the two MiSeq runs combined. Notice the presence in the controls of several bacterial genera, which was likely due to the inherent contamination of laboratory reagents by bacterial DNA (termed “contaminant genera”). These contaminant genera are also present (to a lesser extent) in the rodent samples. The insertions represent the proportion of sequences from rodent samples which were incorrectly assigned to the controls. See Fig. S1 for separate histograms for the two MiSeq runs.

Genera identified as contaminants were then simply removed from the sample data set. Note, however, that the exclusion of these results does not rule out the possibility that our samples represented true rodent infections (at least for some of them, such as those by species of the genus Streptococcus, which contains both saprophytic and pathogenic species). However, as mentioned by Razzauti et al. (8), distinguishing between those two possibilities seems difficult, if not impossible. Faced with this lack of certainty, it was most prudent to simply remove these taxa from the sample data set. These results highlight the importance of carrying out systematic negative-control procedures to filter the taxa concerned in order to prevent inappropriate data interpretation, particularly for the Streptococcus genus, which contains a number of important pathogenic species. The use of DNA-free reagents would improve the quality of sequencing data and likely increase the depth of sequencing of the samples.

After filtering for reagent contaminants was performed as described above, 12 OTUs, belonging to 7 genera for which at least one species or one strain is known to be pathogenic in mammals (and that are therefore referenced here as “pathogenic genera”), accounted for 66% of the sequences identified in wild-rodent samples for the two MiSeq runs combined (Fig. 2). These genera are Bartonella, Borrelia, Ehrlichia, Mycoplasma, Orientia, Rickettsia, and Streptobacillus. Six different OTUs were obtained for Mycoplasma (Mycoplasma_OTU_1 to Mycoplasma_OTU_6) and one OTU each for the other genera (Table 2). Finally, the precise significance of the remaining 34% of sequences, which potentially corresponded to commensal bacteria (Bacteroidales, Bacteroides, Enterobacteriaceae, Helicobacter, Lactobacillus), unknown pathogens, undetected contaminants, or undetected sequencing errors, was undetermined.

Filtering for false-positive results.mothur analysis produced a table of abundance, giving the number of sequences for each OTU in each PCR product (data have been deposited in the Dryad Digital Repository [http://dx.doi.org/10.5061/dryad.m3p7d]) (41). The multiple biases present during the experimental steps and data processing steps listed in Table 1 made it impossible to infer prevalence and cooccurrence directly from the table of sequence presence/absence in the PCR products. We suggest filtering the data with estimates of the different biases calculated from the multiple controls introduced during the process. This strategy involves calculating sequence number thresholds from our bias estimates. Two different thresholds were set for each of the 12 OTUs and two MiSeq runs. We then discarded positive results associated with sequence counts below the threshold (Fig. 1).

Threshold T_CC: filtering for cross-contamination.One source of false positives is cross-contamination between samples processed in parallel (Table 1). Negative controls for dissection (NC_mus), consisting of the spleens of healthy laboratory mice manipulated during sessions of wild-rodent dissection, and negative controls for extraction (NC_ext) and PCR (NC_PCR) were used, together with positive controls for PCR (PC_{Bartonella_t}, PC_{Borrelia_b}, and PC_{Mycoplasma_m}), to estimate levels of cross-contamination. For each sequencing run, we calculated the maximal number of sequences for the 12 pathogenic OTUs in the negative and positive controls. These numbers ranged from 0 to 115 sequences, depending on the OTU and the run considered (Table 2), and we used them to establish values for OTU-specific thresholds of cross-contamination (T_CC) for each run. For example, in sequencing run 2, the highest number of sequences in a control for Mycoplasma_OTU_2 was 115 (in an NC_ext). Therefore, we established the threshold value at 115 sequences for this OTU in sequencing run 2. Thus, PCR products with fewer than 115 sequences for the Mycoplasma_OTU_2 in sequencing run 2 were considered to represent false-positive results for this OTU. The use of these T_CC values led to 0% to 69% of the positive results being discarded, corresponding to only 0% to 0.14% of the sequences, depending on the OTU considered (Fig. 3; see also Table S5 in the supplemental material). A PCR product may be positive for several bacteria in cases of coinfection. In such cases, the use of a T_CC makes it possible to discard the positive result for one bacterium while retaining positive results for other bacteria.

• Open in new tab
• Download powerpoint

FIG 3

Numbers of rodents yielding positive results (positive rodents), and of sequences from positive rodents, removed for each OTU at each step in data filtering. These findings demonstrate that the positive rodents filtered out corresponded to only a very small number of sequences. (Left panel) The histogram shows the number of positive rodents discarded because of likely cross-contamination, false index pairing, and failure to replicate in both PCRs, as well as the positive results retained at the end of data filtering (in green). (Right panel) The histogram shows the number of sequences corresponding to the same class of positive rodents. Note that several positive results may be recorded for the same rodent in cases of coinfection.

Threshold T_FA: filtering out incorrectly assigned sequences.Another source of false positives is the incorrect assignment of sequences to a PCR product (Table 1). This phenomenon may be due either to cross-contamination between indices during the experiment or to the generation of mixed clusters during the sequencing (27).

First, the cross-contamination of indices may happen during the preparation of indexed primer microplates. This cross-contamination was estimated using negative-control index pairs (NC_index) corresponding to particular index pairs not used to identify the samples. NC_index returned very low (1 to 12) read numbers, suggesting that there was little or no cross-contamination between indices in our experiment.

Second, the occurrence of mixed clusters during the sequencing of multiplexed samples was previously reported by Kircher et al. (27). Mixed clusters on an Illumina flow cell surface are described by Kircher et al. (27) as the predominant source of error of sequence assignment to a PCR product. The impact of this phenomenon on our experiment was estimated using “alien” positive controls (PC_alien) (sequenced in parallel with the rodent samples): PC_{Mycoplasma_m}, corresponding to the DNA of Mycoplasma mycoides, which cannot infect rodents; and PC_{Borrelia_b}, containing the DNA of Borrelia burgdorferi, which is not present in Africa. Neither of these bacterial species can survive in abiotic environments, so the presence of their sequences in African rodent PCR products indicates a sequence assignment error due to false index pairing (27). Using PC_{Mycoplasma_m}, we obtained an estimate of the global false-index-pairing rate of 0.14% (i.e., 398 of 280,151 sequences of the Mycoplasma mycoides OTU were assigned to samples other than PC_{Mycoplasma_m}). Using PC_{Borrelia_b}, we obtained an estimate of 0.22% (534 of 238,772 sequences of the Borrelia burgdorferi OTU were assigned to samples other than PC_{Borrelia_b}). These values are very close to the estimate of 0.3% obtained by Kircher et al. (27). Close examination of the distribution of misassigned sequences within the PCR 96-well microplates showed that all PCR products with misassigned sequences had one index in common with either PC_{Mycoplasma_m} or PC_{Borrelia_b} (see Fig. S2 in the supplemental material).

We then estimated the impact of false index pairing for each PCR product by calculating the maximal number of sequences of alien bacteria assigned to PCR products other than the corresponding PC. These numbers ranged from 28 to 43, depending on the positive control for run 1 (Table 2)—run 2 was discarded because of the low values of the numbers of sequences, which were obtained likely due to the fact that DNAs of PC were diluted 100-fold in run 2 (see Table S1 in the supplemental material). We then estimated a false-assignment rate for each PCR product (R_fa) by dividing the numbers given above by the total number of sequences from alien bacteria in sequencing run 1. R_fa was estimated for PC_{Mycoplasma_m} and PC_{Borrelia_b} separately. R_fa reached 0.010% and 0.018% for PC_{Mycoplasma_m} and PC_{Borrelia_b}, respectively. We adopted a conservative approach by fixing the R_fa value at 0.020%. This number signifies that each PCR product may receive a maximum of 0.020% of the total number of sequences of an OTU present in a run due to false index pairing. Moreover, the number of sequences for a specific OTU misassigned to a PCR product should increase with the total number of sequences of the OTU in the MiSeq run. We therefore defined the second threshold (T_FA) as the total number of sequences in the run for an OTU multiplied by R_fa. T_FA values varied with the abundance of each OTU in the sequencing run (Table 2). Because the abundance of each OTU varied from one sequencing run to the next, T_FA also varied according to the sequencing run. The use of the T_FA led to 0% to 87% of positive results being discarded. This corresponded to 0% to 0.71% of the sequences, depending on the OTU (Fig. 3; see also Table S5 in the supplemental material).

Validation performed with PCR replicates.Random contamination may occur during the preparation of PCR 96-well microplates. These contaminants may affect some of the wells, but not those for the negative controls, leading to the generation of false-positive results. We thus adopted a conservative approach in which we considered rodents to be positive for a given OTU only if both PCR replicates were considered to represent positive results after the filtering steps described above were performed. The relevance of this strategy was supported by the strong correlation between the numbers of sequences for the two PCR replicates for each rodent (R² > 0.90) (Fig. 4; see also Fig. S3 in the supplemental material). At this stage, 673 positive results for 419 rodents were validated for both replicates (note that a rodent sample may be positive for several bacterial species and may thus be counted several times), whereas only 52 positive results were discarded because the result for the other replicate was negative. At this final validation step, 0% to 60% of the positive results for a given OTU were discarded, corresponding to only 0% to 7.17% of the sequences (Fig. 3; see also Table S5 and Table S6 in the supplemental material). Note that the number of replicates may be increased, as previously described for the strategy of Gómez-Díaz et al (44).

• Open in new tab
• Download powerpoint

FIG 4

Plots of the number of sequences [log (x + 1) scale] from bacterial OTUs in both PCR replicates (PCR1 and PCR2) of the 348 wild rodents analyzed in the first MiSeq run. Note that each rodent was tested with two replicate PCRs. Green points correspond to rodents with two positive results after filtering; red points correspond to rodents with one positive result and one negative result; and blue points correspond to rodents with two negative results. The light blue area and lines correspond to threshold values used for the data filtering: samples below the lines were filtered out. See Fig. S3 for plots corresponding to the second MiSeq run.

Postfiltering results.Finally, the proportions of rodents positive for a given OTU filtered out by the complete filtering approach ranged from 6% to 86%, depending on the OTU, corresponding to only 1% of the total sequences (Fig. 3). Indeed, our filtering strategy mostly excluded rodents with a small number of sequences for the OTU concerned. These rodents were considered to represent false-positive results.

Refining bacterial taxonomic identification.We refined the taxonomic identification of the 12 bacterial OTUs through phylogenetic and blast analyses. We were able to identify the bacteria present down to the genus level, and, in some cases, we could even identify the most likely species (Table 3; see also Fig. S4 in the supplemental material). For instance, the sequences of the six Mycoplasma OTUs were consistent with three different species—M. haemomuris for OTU_1 and OTU_3, M. coccoides for OTU_4, OTU_5, and OTU_6, and M. sp. nov. (45) for OTU_2—with high percentages of sequence identity (≥93%) and bootstrap values (≥80%). All three of these species belong to the Hemoplasma group, whose members are known to infect mice, rats, and other mammals (46, 47) and are thought to cause anemia in humans (48, 49). The Borrelia sequences grouped with three different species of the relapsing fever group (B. crocidurae, B. duttonii, and B. recurrentis) with a high percentage of identity (100%) and a bootstrap value of 71%. In West Africa, B. crocidurae causes severe borreliosis, a rodent-borne disease transmitted by ticks and lice (50). The Ehrlichia sequences were 100% identical to and clustered with the sequences of the recently described “Candidatus Ehrlichia khabarensis” species isolated from voles and shrews in the Far East of the Russian Federation (51). The Rickettsia sequences were 100% identical to the sequence of R. typhi, a species of the typhus group responsible for murine typhus (52), but species of this clade were differentiated from many other Rickettsia species only weakly (bootstrap support of 61%). The most likely species corresponding to the sequences of the Streptobacillus OTU was S. moniliformis, with a high percentage of identity (100%) and a bootstrap value of 100%. This bacterium is common in rats and mice and causes a form of rat-bite fever, Haverhill fever (53). The Orientia sequences corresponded to O. chuto, with a high percentage of identity (100%) and a bootstrap value of 77%. This species was recently isolated from a patient infected in Dubai (54). Finally, accurate species determination was not possible for Bartonella, as the 16S rRNA gene does not resolve the species of this genus well (55). Indeed, the sequences from the Bartonella OTU detected in our rodents corresponded to at least seven different species (B. elizabethae, B. japonica, B. pachyuromydis, B. queenslandis, B. rattaustraliani, B. tribocorum, and B. vinsonii) and a putative new species recently identified in Senegalese rodents (56).

These findings demonstrate the considerable potential of 16S rRNA amplicon sequencing for the rapid identification of zoonotic agents in wildlife, provided that the postsequencing data are cleaned beforehand. Borrelia (50) and Bartonella (56) were the only ones of the seven pathogenic bacterial genera detected here in Senegalese rodents to have been reported as present in rodents from West Africa before. The other bacterial genera identified here have previously been reported to be present in rodents only in other parts of Africa or on other continents. Streptobacillus moniliformis has recently been detected in rodents from South Africa (57), and there have been a few reports of human streptobacillosis in Kenya (58) and Nigeria (59). R. typhi was recently detected in rats found in Zaire, in Central Africa (60), and human seropositivity for this bacterium has been reported in coastal regions of West Africa (61). With the exception of one study in Egypt published some time ago (62), Mycoplasma spp. have never before been reported in African rodents. Several species of Ehrlichia (from the E. canis group: E. chaffeensis, E. ruminantium, E. muris, and E. ewingii) have been characterized in West Africa, but only in ticks from cattle (63), and a report of possible cases of human ehrlichioses in this region was previously published (64). Finally, the present study reports the first identification of Orientia in African rodents (9). There have already been a few reports of suspected human infection with this bacterium in Zaire, Cameroon, Kenya, and Tanzania (65).

Estimating prevalence and coinfection.After performing data filtering, we were able to estimate the prevalence in rodent populations and to assess coinfection in individual rodents for the 12 bacterial OTUs. Rates of bacterial prevalence differed considerably between rodent species (Table 3). Bartonella spp. were highly prevalent in the two multimammate rat species M. natalensis (79%) and M. erythroleucus (27%); Orientia spp. were prevalent in the house mouse species M. musculus (22%); and Ehrlichia spp. occurred frequently in only one of the two rats of the multimammate species M. erythroleucus (18%). In contrast, the prevalence of Streptobacillus and Rickettsia was low in all rodent species (<5%). Coinfection was common, as 184 rodents (26%) were found to be coinfected with bacteria from two (19%), three (5%), four (2%), or five (0.1%) different bacterial pathogens.

Interestingly, several Mycoplasma OTUs appeared to be specific to a rodent genus or species (Table 3 and Fig. 5A). OTU_2, putatively identified as a recently described lineage isolated from brown rat, Rattus norvegicus (45), was specifically associated with R. rattus in this study. Of the OTUs related to M. coccoides, OTU_4 was found exclusively in R. rattus, whereas OTU_5 and OTU_6 seemed to be specific to the two multimammate rats (M. erythroleucus and M. natalensis). Comparative phylogenies of Mycoplasma OTUs and rodents showed that R. rattus, which is phylogenetically more distantly related to the other three rodents, contained a Mycoplasma community different from that in the rodents in the Mus-Mastomys clade (Fig. 5A). Pathogen prevalences also differed considerably between sites, as shown for the six Mycoplasma OTUs (Fig. 5B). This suggests that the infection risks for animals and humans vary greatly according to environmental characteristics and/or biotic features potentially related to recent changes in the distribution of rodent species in Senegal (66, 67).

• Open in new tab
• Download powerpoint

FIG 5

Prevalence of Mycoplasma lineages in Senegalese rodents, by site, and phylogenetic associations between Mycoplasma lineages and rodent species. (A) Comparison of phylogenetic trees based on the 16S rRNA V4 sequences of Mycoplasma and on the mitochondrial cytochrome b gene and the two nuclear gene fragments (IRBP exon 1 and GHR) for rodents (the tree was drawn based on data from reference 92). Lines link the Mycoplasma lineages detected in the various rodent species (for a minimum site prevalence exceeding 10%). The numbers next to the branches are bootstrap values (shown only if >70%). (B) Plots of OTU prevalences, with 95% confidence intervals calculated by Sterne’s exact method (93) according to rodent species and site (see reference 69 for more information about site codes and their geographic locations). The gray bars on the x axis indicate sites from which the rodent species concerned is absent.

Perspectives. (i) Recommendation for future experiments.Our experiments demonstrated the need to include many different kinds of controls, at different steps, in order to avoid data misinterpretation. In particular, alien positive controls are important for establishing threshold values for OTU positivity. These alien positive controls should include taxa that are distant enough from those of potential pathogens to avoid potential confusion between sequences of alien controls and sequences that result from actual infection of rodent samples. Ideally, one should choose alien positive controls from bacterial genera which are not able to infect the study’s host species. In our study, the use of Mycoplasma and Borrelia species as alien positive controls was not ideal because both genera are potential rodent pathogens. Thankfully, the species used as alien controls could be easily distinguished from the species found in rodents on the basis of the phylogenetic analyses of the V4 sequences. However, on the basis of our experience, we recommend using bacterial genera phylogenetically distant from pathogenic genera as alien controls when possible.

The inclusion of negative controls of PCR and DNA extraction in studies based on massive sequencing of 16S rRNA amplicons had long been overlooked until a report published by Salter in 2014 (23) demonstrated the high pollution of laboratory reagents by bacterial DNA. Most studies published prior to that report reported no extraction controls in their protocols. Here, we performed one negative control for extraction per DNA extraction microplate; with each run consisting of four DNA extraction microplates, and each control having been analyzed in two replicate experiments, we had a total of 8 negative controls for extraction per run which were analyzed twice. On the basis of our experience, we recommend performing at least this number of extraction controls per run. Further increases in the number of extraction controls per microplate would further improve the efficiency of data filtering and so the quality of the data produced.

The protocol of PCR amplification is also of importance for ensuring data quality. In our study, we built separate amplicon libraries for each sample separately and used very long (5-min) PCR elongation times in order to mitigate the formation of chimeric reads (18). High numbers of PCR cycles are also known to increase chimera formation, yet, as mentioned by Schnell et al. (68), this parameter is mainly critical only when bulk amplification of pools of tagged/indexed amplicons is performed (e.g., when using an Illumina TrueSeq library preparation kit). As we used separate amplicon libraries for each sample, we believe that the relatively high number of PCR cycles that we used (40 cycles) had minimal impact on chimera formation and that this protocol ensures the absence of chimeric sequences between samples. We had chosen to maximize the number of cycles to enhance our ability to detect pathogenic bacteria, which are sometimes of low quantity in animal samples. Fine-tuning the balance between these parameters deserves further study.

Moreover, we targeted the spleen to detect bacterial infections in our study based on the fact that this organ is known to filter microbial cells in mammals. However, we lack the data to be certain that the spleen is the best organ for use in giving a global picture of bacterial infection in rodents (and, more broadly, in vertebrates). We are currently conducting new experiments to address this issue.

Finally, in our experiments, about a third of OTU sequences were attributable neither to contamination nor to (known) pathogenic genera. We currently have no precise idea of the significance of the presence of these OTUs in the rodent spleens. Some of these OTUs could be linked to further undetected biases during data generation; in spite of all the precautions we have implemented here, other biases may still elude detection. Such biases could explain the very high numbers of rare OTUs (11,947 OTUs corresponding to <100 reads), which together represent more than 88% of the total number of OTUs but less than 1% of the total number of sequences (with the two runs combined).

Additionally, the presence of an OTU in a rodent spleen does not necessarily imply that the OTU is pathogenic. We know little about the microbiome of healthy organs of vertebrates, and yet the sharp increase in the number of microbiome studies over the last few years has led to the discovery that communities of microbiota appear to be specific to each part of the vertebrate’s body, including internal tissues and blood (69). The OTUs detected in rodent’s spleen could thus simply be part of the healthy microbiome of the organ. These issues deserve better documentation. Our results thus pave the way for future research on unknown bacterial pathogens and on the microbiome of healthy organs in vertebrates.

(ii) Improving HTS for epidemiological surveillance.The screening strategy described here has the considerable advantage of being nonspecific, making it possible to detect unanticipated or novel bacteria. Razzauti et al. (8) recently showed that the sensitivity of 16S rRNA amplicon sequencing on the MiSeq platform was equivalent to that of whole-RNA sequencing (RNA-seq) on the HiSeq platform for detecting bacteria in rodent samples. However, little is known about the comparative sensitivities of HTS approaches relative to the sensitivity of quantitative PCR (qPCR) performed with specific primers, the current gold standard for bacterial detection within biological samples. Additional studies are required to address this issue. Moreover, as 16S rRNA amplicon sequencing is based on a short sequence, it does not yield results sufficiently high in resolution to distinguish between species in some bacterial genera, such as Bartonella, or to distinguish between pathogenic and nonpathogenic strains within the same bacterial species. To get this information, we thus need to follow up the 16S rRNA amplicon sequencing with complementary laboratory work. Whole-genome shotgun or RNA-seq techniques provide longer sequences, through the production of longer reads or the assembly of contigs, and they might therefore increase the accuracy of species detection (70). However, these techniques would be harder to adapt for the extensive multiplexing of samples (8). Other methods could be used to assign sequences to bacterial species or strains for samples found to be positive for a bacterial genus following the 16S rRNA screening. For example, positive PCR assays could be carried out with bacterial genus-specific primers, followed by amplicon sequencing, as commonly used in multilocus sequence analysis (MLSA) strategies (71), or high-throughput microfluidic qPCR assays based on bacterial species-specific primers could be used (72). High-throughput amplicon sequencing approaches could be fine-tuned to amplify several genes for species-level assignment, such as the gltA gene used by Gutiérrez et al. (73) for the Bartonella genus, in parallel with the 16S rRNA-V4 region.

This strategy could also easily be adapted for other microbes, such as protists, fungi, and even viruses, provided that universal primers are available for their detection (see references 74 and 75 for protists and fungi and reference 76 for degenerate virus family-level primers for viruses). Finally, our filtering method could also be translated to any other postsequencing data set of indexed or tagged amplicons in the framework of environmental studies (e.g., metabarcoding for diet analysis and biodiversity monitoring [77], the detection of rare somatic mutations [78], or the genotyping of highly polymorphic genes [e.g., major histocompatibility complex {MHC} or HLA typing] [79, 80]).

(iii) Monitoring the risk of zoonotic diseases.Highly successful synanthropic wildlife species, such as the rodents studied here, will probably play an increasingly important role in the transmission of zoonotic diseases (81). Many rodent-borne pathogens cause only mild or undifferentiated disease in healthy people, and these illnesses are often misdiagnosed and underreported (53, 82–85). The information about pathogen circulation and transmission risks in West Africa provided by this study is important in terms of human health policy. We show that rodents carry species of seven major pathogenic bacterial genera: Borrelia, Bartonella, Mycoplasma, Ehrlichia, Rickettsia, Streptobacillus, and Orientia. The last five of these genera have never before been reported in West African rodents. The data generated with our HTS approach could also be used to assess zoonotic risks and to formulate appropriate public health strategies involving the focusing of continued pathogen surveillance and disease-monitoring programs on specific geographic areas or rodent species likely to be involved in zoonotic pathogen circulation, for example.