Project Data

Location:

Berkeley, CA

Posted:

November 13, 2012

Contact this candidate

Resume:

The ISME Journal (****) *, **** ****

www.nature.com/ismej

SHORT COMMUNICATION

Selection of primers for optimal taxonomic

classification of environmental 16S rRNA

gene sequences

David AW Soergel1, Neelendu Dey2,4,5, Rob Knight3 and Steven E Brenner1

Department of Plant and Microbial Biology, University of California, Berkeley, CA, USA; 2Division of

Gastroenterology, Department of Medicine, University of California, San Francisco, CA, USA and

Howard Hughes Medical Institute and Department of Chemistry and Biochemistry, University of Colorado,

Boulder, CO, USA

Microbial community profiling using 16S rRNA gene sequences requires accurate taxonomy assignments.

Universal primers target conserved sequences and amplify sequences from many taxa, but they provide

variable coverage of different environments, and regions of the rRNA gene differ in taxonomic

informativeness especially when high-throughput short-read sequencing technologies (for example,

454 and Illumina) are used. We introduce a new evaluation procedure that provides an improved

measure of expected taxonomic precision when classifying environmental sequence reads from a

given primer. Applying this measure to thousands of combinations of primers and read lengths,

simulating single-ended and paired-end sequencing, reveals that these choices greatly affect

taxonomic informativeness. The most informative sequence region may differ by environment, partly

due to variable coverage of different environments in reference databases. Using our Rtax method of

classifying paired-end reads, we found that paired-end sequencing provides substantial benefit in

some environments including human gut, but not in others. Optimal primer choice for short reads

totaling 96 nt provides 82 100% of the confident genus classifications available from longer reads.

The ISME Journal (2012) 6, 1440 1444; doi:10.1038/ismej.2011.208; published online 12 January 2012

Subject Category: microbial population and community ecology

Keywords: 16S ribosomal RNA; taxonomy; phylogeny; classification; bacteria; sequencing

Variation in 16S ribosomal gene sequences has been Previous work on taxonomic classification of

used since the mid-1980 s to characterize microbial environmental 16S rRNA gene sequences has

diversity (Stahl et al., 1984). Interest in sequence- focused on whether reference sequences matching

based surveys of environmental microbes has a given query share taxonomic annotations

exploded in recent years with the availability of (Jonasson et al., 2002; Desantis et al., 2006; Sogin

sequencing technologies that produce ever-larger et al., 2006; Wang et al., 2007). Validations of

data sets at ever-decreasing cost; in particular, taxonomic classifiers have typically compared a

the Illumina platform is attractive because of limited range of primers, read lengths and environ-

throughput, despite its short reads (Sogin et al., ments (Sundquist et al., 2007; Huse et al., 2008;

2006; Lazarevic et al., 2009; Claesson et al., 2010; Liu et al., 2008; Wu et al., 2008; Hamp et al., 2009).

Caporaso et al., 2011; Degnan and Ochman, 2012). Reference databases contain many sequences de-

Here, we examine the reliability of assignment of rived from some environments and few associated

novel sequences to known taxa under thousands with others, however (Supplementary Figure S1),

of simulated scenarios, varying primer choice, read leading to substantial variation in classification quality.

length and environment. In addition, the use of leave-one-out cross-validation

at the sequence level (Sundquist et al., 2007; Wang

et al., 2007; Liu et al., 2008; Wu et al., 2008) where

Correspondence: DAW Soergel. Current address: Department of

a single sequence with a known annotation is held

Computer Science, University of Massachusetts, 140 Governors

Drive, Amherst, MA 01003-9264, USA. out from a reference database and classified using

E-mail: *******@**.*****.***

the remainder is problematic: reference sequences

Current address: Center for Genome Sciences and Systems

matching held-out query sequences are likely to

Biology, Washington University School of Medicine, Saint Louis,

originate from the same sample, because natural

MO 63108, USA.

Current address: Division of Gastroenterology, Department of environments contain microdiverse clusters of

Internal Medicine, Washington University School of Medicine,

closely related strains (Acinas et al., 2004).

Saint Louis, MO 63108, USA.

We addressed these issues by simulating

Received 28 March 2011; revised 21 November 2011; accepted 3

truncated reads from eight large environmental data

December 2011; published online 12 January 2012

Primer choice for 16S taxonomic classification

DAW Soergel et al

1441

sets of near-full-length 16S rRNA gene sequences 64 nt, 80 nt, 96 nt, 112 nt, 128 nt, 260 nt, 400 nt,

extracted from GreenGenes (Supplementary Table 800 nt and full-length), with the constraint that read

S1), using pairs of 44 universal primers commonly length could not exceed amplicon length for each

found in the literature (Supplementary Tables S2 primer pair, produced 6617 single-end and 3061

and S3). These were selected from an initial set of 94 paired-end datasets per environment.

primers by the criterion that each primer had to Reference databases were constructed by holding

match at least 40% of the sequences in at least one of out each entire study in turn from GreenGenes,

the chosen environmental samples. Single-end clustering the remainder at 99% using UCLUST

reads were tested from each primer with all viable (Edgar, 2010), and selecting one representative

amplification partners (794 combinations), and sequence per cluster (see Supplementary Methods

paired-end reads were tested using all 374 viable for details). Each query fragment was then matched

pairings of the 22 forward and 22 reverse primers. against remaining representative sequences using

Simulations using 11 read lengths (32 nt, 48 nt, USEARCH (Edgar, 2010), configured to penalize

Figure 1 Classification performance, at three levels of estimated accuracy (Supplementary Methods), of 6617 possible choices of

amplification primer, sequencing primer and read length for single-ended reads from different environments (left portion of each panel)

and 3061 possible choices of primer pair and read length for paired-end reads (right portion). Combinations of primers and read lengths

are sorted on the x axis according to a measure of overall classification performance (Supplementary Methods). Stacked bars show the

proportion of non-chimeric, non-unique sequences from each sample not the proportion of the total sample that can be classified to

each taxonomic level for each combination. See Supplementary Figure S1 and Supplementary Table S1 for the excluded proportion of

novel (and thus a priori unclassifiable) sequences in each sample. The top of each colored section indicates how much of the sample can

be classified to the given level or better. Primer miss (black) indicates sequences that did not match a given primer and so would not be

amplified. Classifications more specific than the genus level are exceedingly rare and so are not visible here. Horizontal lines indicate the

maximum proportion of each sample classifiable to the genus level using 96 nt or less of sequence (i.e., with an optimal choice of primer

or primer pair; see also Supplementary Tables S4 and S5), showing that short reads from the best primers frequently but not always

provide taxonomic information nearly matching that obtained from longer read lengths. Full-size versions of these panels are available

in the supplementary data.

The ISME Journal

Primer choice for 16S taxonomic classification

DAW Soergel et al

1442

Table 1 Genus classification rates for optimal choices of primers, grouped by total read length

Hypersaline Mat

Grassland Soil

Steer Rumen

Dust & Skin

Termite Gut

Human Gut

Ocean

Coral

Total

nucleotides Read length Forward primer Reverse primer

32 single 32 (E786F end) E826R 0 5 2 5 3 25 42 83

no improvement

single 64 (end) E533Ra **-**-**-**-**-** 76 83

pair 3 2 E9 69F E1492 R 17 24 74 90 78

80 single 80 E341F (E1406R E533Ra end) **-**-**-**-**-** 66 45

pair 4 8 E3 41F E926Ra **-**-**-**-**-** 73 64

single 96 E517F U515F (end) **-**-**-**-**-** 65 32

pair 4 8 E3 41F E1064 R **-**-**-**-**-** 80 42

96 single 96 E3 4 1 F (end) **-**-**-**-**-** 80 45

single 96 E343F U341F (end) **-**-**-**-**-** 80 46

pair 4 8 E517F U515F E92 6 R a **-**-**-**-**-** 70 45

pair 4 8 E9 69F E1492 R 16 16 76 92 65

single 112 E517F U515F (end) **-**-**-**-**-** 66 42

112

single 112 (end) E926Ra **-**-**-**-**-** 80 40

pair 6 4 E3 41F E1406 R **-**-**-**-**-** 79 46

pair 6 4 E517F U515F E1 4 0 6 R **-**-**-**-**-** 72 47

percentage of sample classified to genus level

single 128 E341F (E1406R E533Ra end) **-**-**-**-**-** 73 71

pair 6 4 E517F U515F E1 4 0 7 R **-**-**-**-**-** 72 47

pair 64 E343F E1406R U1406R **-**-**-**-**-** 80 46

pair 6 4 E3 43F E926Ra **-**-**-**-**-** 86 42

128

single 128 E5 1 7 F (end) **-**-**-**-**-** 62 43

single 128 (end) E1406 R **-**-**-**-**-** 85 44

pair 6 4 U519F E1 4 0 6 R **-**-**-**-**-** 71 13

single 128 (end) E357R **-**-**-**-**-** 71 42

pair 6 4 E3 41F E343 F E1 4 9 2 R 39 32 78 92 82

pair 6 4 E517F U515F E1 4 9 2 R 36 34 78 94 78

pair 80 E517F U515F E1406R E1407R U1406R **-**-**-**-**-** 80 45

pair 80 E341F E1406R E1407R **-**-**-**-**-** 86 46

160 pair 8 0 E517F U515F E9 2 6 R a **-**-**-**-**-** 76 43

pair 8 0 E517F U515F E1 4 9 2 R 39 36 79 95 81

pair 8 0 E3 41F E1492 R 39 36 77 91 85

pair 9 6 E5 17F E1407 R **-**-**-**-**-** 80 45

pair 9 6 E5 17F E1406 R **-**-**-**-**-** 79 45

192

pair 9 6 E5 17F E926Ra **-**-**-**-**-** 76 42

pair 9 6 E3 41F E926Ra **-**-**-**-**-** 88 41

pair 112 E341F E1406R E1407R U1406R **-**-**-**-**-** 87 45

pair 112 E5 1 7 F E1406 R **-**-**-**-**-** 80 43

224

pair 112 E5 1 7 F E926Ra **-**-**-**-**-** 74 42

classification

pair 112 U519F E1 4 0 7 R **-**-**-**-**-** 81 13

estimated

accuracy:

pair 128 E517F U515F E1406R U1406R **-**-**-**-**-** 79 44

pair 128 E341F E1406R E1407R **-**-**-**-**-** 85 47

256

pair 128 E517F U515F E9 2 6 R a **-**-**-**-**-** 76 43

pair 128 U519F E1 4 0 6 R **-**-**-**-**-** 79 13

bold: >= 95%

normal: 80% - 95%

italic:

single 260 (end) E1406 R **-**-**-**-**-** 86 44

260

single 260 (end) E926Ra **-**-**-**-**-** 84 41

no improvement

400

pair 260 E517F U515F E1 4 0 6 R **-**-**-**-**-** 78 43

pair 260 E341F E1406R E1407R **-**-**-**-**-** 86 44

520

pair 260 E3 4 1 F E926Ra **-**-**-**-**-** 83 43

pair 260 E3 4 1 F E1492 R 36 33 80 90 89

pair 400 E3 4 1 F E1406 R **-**-**-**-**-** 87 43

800

pair 400 E517F U515F E9 2 6 R a **-**-**-**-**-** 70 42

1600 no improvement

Maximum **-**-**-**-**-** 89 83

Thousands of combinations that produce suboptimal results are not shown (see Supplementary Methods). Cells are colored on a gradient from

worst (red) to best (green) per column. Estimated classification accuracy (Supplementary Methods) is indicated by bold or italic font. Primers in

parentheses are used in single-ended experiments for amplification but not sequencing. Primers appearing together perform equivalently; that is,

for a given row, any choice among the given sequencing and amplification primers will produce the same result. End indicates an end primer

such as E8F, E1406R, U1406R, E1407R, E1492R or E1506R. Primer E1492R could not be tested in three datasets because sequences were not long

enough; the corresponding cells remain blank.

The ISME Journal

Primer choice for 16S taxonomic classification

DAW Soergel et al

1443

indels and mismatches equally. Clusters were then extent to which it overlaps any of the classical

selected that matched within 0.5% identity of the V-regions .

best hit (hits o80% identity were disregarded). For No one combination of primers and read length

paired-end query sequences, our Rtax procedure works best in all environments, but near-optimal

(Supplementary Methods) selected those reference performance in six out of the eight environments is

clusters that matched both reads simultaneously available using paired-end 80 nt reads from primers

with an average percent identity within 0.5% such as E517F, U515F or E341F paired with E1406R

identity of the maximum. Taxonomic classifications or closely related primers (Table 1). However,

were made at each level by retaining annotations practical considerations such as ability to amplify

agreeing among 450% of the clusters (including low-biomass samples will sometimes influence which

those with no annotation in the denominator); primers are used. For instance, short amplicons

these generally extended at best to the genus may be preferred because these are less subject

level, because the reference database provides few to length heterogeneity biases and chimera forma-

species-level annotations. tion. Similarly, short single-ended sequences are

Sequences from novel taxa (or sequences that less subject to errors due to chimeras, simply

appear novel due to sequencing error or chimerism) because they are less likely to span a breakpoint.

clearly cannot be correctly classified; however, such Classification performance for experimental choices

sequences may constitute a substantial proportion of matching such constraints can be found in the

a given sample (Supplementary Figure S1 and supplementary data.

Supplementary Table S1). The version of Green- The choice of reference database and taxonomy

Genes that we used excluded taxa (defined by 97% can have a dramatic impact on the resulting

identity) that were unique to a single sample, as one classification accuracy. In this study, we used the

of the several strategies to remove chimeras. These current GreenGenes taxonomy, which has been

unique sequences were therefore excluded from our filtered to remove chimeras and where the taxo-

query sets. Thus the classification rates we report nomic annotations are comprehensive and consis-

represent the proportion of non-chimeric, non- tent with the phylogenetic tree (McDonald et al.,

unique sequences that can be classified to each 2011). Experiments using a previous version of the

rank. If an environmental sample is not similarly GreenGenes taxonomy lacking these features

filtered before classification, then the classifiable yielded far poorer accuracy (data not shown). In

proportion (that is, taken with respect to the total addition, bolstering areas of low coverage in refer-

sample) will be correspondingly lower. ence databases will substantially improve classifier

Classification rate and accuracy vary widely performance. For instance, taxa in the hypersaline

among environments and sequence regions, for mat, coral and grassland soil samples were under-

several reasons: (1) the reference database provides represented in the reference database (Supple-

different levels of coverage of each environment, mentary Figure S1), and presumably as a conse-

(2) no primer is truly universal and different quence classifications of sequences from those

primers (and pairs) hit different proportions of samples were less likely to prove correct (Figure 1).

sequences in each environment and (3) the targeted Additional data sets from poorly sampled environ-

regions are variably informative. Figure 1 shows ments will also help to distinguish chimeric from

proportions of sequences from each environment legitimate but novel sequences.

classified to each rank, for all 9678 single-ended and In combination, these results indicate that taxo-

paired-end primer and read length combinations. nomic classifications of short reads especially

Horizontal panels compare unfiltered results to genus-level classifications should be treated with

classifications passing 80% and 95% estimated skepticism, unless the specific combination of

accuracy filters (see Supplementary Methods), primer, read length, environmental source, reference

showing that most classifications can be made with database and assignment method has been thor-

high accuracy when optimal primers are chosen. oughly validated. At the same time, optimal choices

Remarkably, only 96 nt of sequence (taken as a single of these parameters allow high classification rates

read or as a pair of 48 nt reads) can provide 82 100% and high accuracy. Thus, large-scale projects such as

of the 80% accurate genus classifications available the Earth Microbiome Project (Gilbert et al., 2010),

from any read length (Supplementary Table S4). which aims to collect and analyze samples from tens

Paired-end sequencing can provide substantial of thousands of microbial habitats around the globe,

gains in classification rate for some but not all may reasonably proceed with standardized primer

environments and read lengths. Paired-end classi- choices and short reads.

fications are typically more accurate than those

made from single reads, and so are more likely

to pass the 95% estimated accuracy filter (Supple-

Acknowledgements

mentary Table S5). Another surprise is that hyper-

variable regions need not be specifically targeted, This research used ShaRCS, UC Shared Research Comput-

as there is no obvious relationship between ing Services Cluster, which is technically supported by

taxonomic informativeness of a region and the multiple UC IT divisions and managed by the University

The ISME Journal

Primer choice for 16S taxonomic classification

DAW Soergel et al

1444

of California, Office of the President. This work was taxonomy using SSU rRNA hypervariable tag sequen-

supported in part by the National Institutes of Health cing. PLoS Genet 4: e1000255.

(HG4872; HG4866; T32 institutional training grant Jonasson J, Olofsson M, Monstein HJ. (2002). Classifica-

DK00700733) and the Howard Hughes Medical Institute. tion, identification and subtyping of bacteria based on

pyrosequencing and signature matching of 16S rDNA

fragments. APMIS 110: 263 272.

Lazarevic V, Whiteson K, Huse S, Hernandez D, Farinelli L,

Osteras M et al. (2009). Metagenomic study of the oral

microbiota by Illumina high-throughput sequencing.

References J Microbiol Methods 79: 266 271.

Liu CH, Lee SM, Vanlare JM, Kasper DL, Mazmanian SK.

Acinas SG, Klepac-Ceraj V, Hunt DE, Pharino C, Ceraj I,

(2008). Regulation of surface architecture by symbiotic

Distel DL et al. (2004). Fine-scale phylogenetic

bacteria mediates host colonization. Proc Natl Acad

architecture of a complex bacterial community.

Sci USA 105: 3951 3956.

Nature 430: 551 554.

McDonald D, Price MN, Goodrich J, Nawrocki EP,

Caporaso JG, Lauber CL, Walters WA, Berg-Lyons D,

DeSantis TZ, Probst A et al. (2011). An improved

Lozupone CA, Turnbaugh PJ et al. (2011). Global

Greengenes taxonomy with explicit ranks for ecologi-

patterns of 16S rRNA diversity at a depth of millions

cal and evolutionary analyses of bacteria and archaea.

of sequences per sample. Proc Natl Acad Sci USA

ISME J; e-pub ahead of print 1 December 2011;

108(Suppl 1): 4516 4522.

doi:10.1038/ismej.2011.139.

Claesson MJ, Wang Q, O Sullivan O, Greene-Diniz R, Cole

Sogin ML, Morrison HG, Huber JA, Welch DM, Huse SM,

JR, Ross RP et al. (2010). Comparison of two next-

Neal PR et al. (2006). Microbial diversity in the deep

generation sequencing technologies for resolving

sea and the underexplored 00 rare biosphere00 . Proc Natl

highly complex microbiota composition using tandem

Acad Sci USA 103: 121**-*****.

variable 16S rRNA gene regions. Nucleic Acids Res 38:

Stahl DA, Lane DJ, Olsen GJ, Pace NR. (1984). Analysis of

e200.

hydrothermal vent-associated symbionts by ribosomal

Desantis TZ, Hugenholtz P, Larsen N, Rojas M, Brodie EL,

RNA sequences. Science 224: 409 411.

Keller K et al. (2006). Greengenes, a chimera-checked

Sundquist A, Bigdeli S, Jalili R, Druzin ML, Waller S,

16S rRNA gene database and workbench compatible

Pullen KM et al. (2007). Bacterial flora-typing with

with ARB. Appl Environ Microbiol 72: 5069 5072.

targeted, chip-based pyrosequencing. BMC Microbiol

Degnan PH, Ochman H. (2012). Illumina-based analysis

7: 108.

of microbial community diversity. The ISME J 6:

Wang Q, Garrity GM, Tiedje JM, Cole JR. (2007).

183 194.

Naive Bayesian classifier for rapid assignment of

Edgar RC. (2010). Search and clustering orders of

rRNA sequences into the new bacterial taxonomy.

magnitude faster than BLAST. Bioinformatics 26:

Appl Environ Microbiol 73: 5261 5267.

2460 2461.

Wu D, Hartman A, Ward N, Eisen JA. (2008). An automated

Gilbert JA, Meyer F, Antonopoulos D, Balaji P, Brown CT,

phylogenetic tree-based small subunit rRNA taxonomy

Brown CT et al. (2010). Meeting report: the terabase

and alignment pipeline (STAP). PLoS One 3: e2566.

metagenomics workshop and the vision of an Earth

Microbiome Project. Stand Genomic Sci 3: 243 248.

Hamp TJ, Jones WJ, Fodor AA. (2009). Effects of experi-

This work is licensed under the Creative

mental choices and analysis noise on surveys of

CommonsAttribution-NonCommercial-Share

the 00 rare biosphere00 . Appl Environ Microbiol 75:

Alike 3.0 Unported License. To view a copy of this

3263 3270.

license, visit http://creativecommons.org/licenses/

Huse SM, Dethlefsen L, Huber JA, Welch DM, Relman DA,

by-nc-sa/3.0/

Sogin ML. (2008). Exploring microbial diversity and

Supplementary Information accompanies the paper on The ISME Journal website (http://www.nature.com/ismej)

The ISME Journal

1444

Contact this candidate