Bacterial Community Shift in Treated Periodontitis Patients Revealed by Ion Torrent 16S rRNA Gene Amplicon Sequencing

Periodontitis, one of the most common diseases in the world, is caused by a mixture of pathogenic bacteria and inflammatory host responses and often treated by antimicrobials as an adjunct to scaling and root planing (SRP). Our study aims to elucidate explorative and descriptive temporal shifts in bacterial communities between patients treated by SRP alone versus SRP plus antibiotics. This is the first metagenomic study using an Ion Torrent Personal Genome Machine (PGM). Eight subgingival plaque samples from four patients with chronic periodontitis, taken before and two months after intervention were analyzed. Amplicons from the V6 hypervariable region of the 16S rRNA gene were generated and sequenced each on a 314 chip. Sequencing reads were clustered into operational taxonomic units (OTUs, 3% distance), described by community metrics, and taxonomically classified. Reads ranging from 599,933 to 650,416 per sample were clustered into 1,648 to 2,659 non-singleton OTUs, respectively. Increased diversity (Shannon and Simpson) in all samples after therapy was observed regardless of the treatment type whereas richness (ACE) showed no correlation. Taxonomic analysis revealed different microbial shifts between both therapy approaches at all taxonomic levels. Most remarkably, the genera Porphyromonas, Tannerella, Treponema, and Filifactor all harboring periodontal pathogenic species were removed almost only in the group treated with SPR and antibiotics. For the species T. forsythia and P. gingivalis results were corroborated by real-time PCR analysis. In the future, hypothesis free metagenomic analysis could be the key in understanding polymicrobial diseases and be used for therapy monitoring. Therefore, as read length continues to increase and cost to decrease, rapid benchtop sequencers like the PGM might finally be used in routine diagnostic.


Introduction
Chronic periodontitis has a complex etiology and is an endemic inflammatory disease caused by a mixed bacterial infection with pathogens in intraoral biofilms. It is characterized by severe destruction of teeth supporting periodontal tissues, i.e., the connective tissue attachment apparatus and alveolar bone. Advanced forms of periodontitis may lead to edentulism in adulthood despite therapy and cause high costs for prosthetic rehabilitation. Several bacteria, mainly Gram-negative, are shown to be strongly associated with periodontal infections, e.g. Porphyromonas gingivalis, Tannerella forsythia, and Treponema denticola [1]. Of the about 700 species detected in the oral cavity based on culture and molecular analysis, more than half have not yet been cultivated [2,3]. Approximately 400 of these species have been detected in subgingival pockets, although single patients harbor only a fraction of these species. Recently, next-generation sequencing techniques were used to investigate oral samples from healthy individuals enabling for a deeper coverage of the microbial community [4,5].
Mechanical debridement, i.e. the removal of biofilms from tooth surfaces by scaling and root planing (SRP) is the established gold standard in the treatment of periodontal disease [1] but adjunctive systemic antibiotics are frequently used as supportive therapeutic agents. Two systematic reviews about the impact of various adjunctive antibiotic agents in periodontal therapy conclude that antibiotics may offer a clinical benefit [6,7]. However, it is insufficiently studied which periodontitis patients may profit from antibiotic treatment, what the preferred treatment interval might be, and how antibiotics affect the whole subgingival microflora. Therefore, our study aims to elucidate explorative and descriptive temporal shifts in bacterial communities as analyzed by nextgeneration sequencing between patients treated by SRP alone versus SRP with adjunctive antibiotics. To explore richness and diversity of subgingival plaque samples the hypervariable region V6 of the 16S rRNA gene was targeted. By using the fast and per experiment rather affordable Ion Torrent's Personal Genome Machine (PGM) a platform was applied that eventually might also be used in routine diagnostics.

Patients, Sample Collection and DNA Isolation
In total, samples from four patients were taken before (pre) and two months after (post) intervention. Patients were part of a doubleblind, parallel group, randomized, placebo-controlled multi-center efficacy study (ISRCTN64254080). The study was approved by the Medical Ethical Committee of the University Münster (ref: 2006-474-f-A) and a signed written informed consent was obtained from each patient. The four selected patients fulfilled the study's inclusion criteria and were part of one stratum (non-smokers, generalized severe chronic periodontitis, i.e. more than 38% of sites with pocket probing depths of 6 mm or more). The intervention consisted either of mechanical debridement plus 500 mg amoxicillin and 400 mg metronidazole three times daily for 7 days (experimental, Ex; two patients) or mechanical debridement alone (control, Co; two patients). For microbiological analysis, subgingival plaque samples from each patient were pooled. Specimens were taken from four different teeth with at least one site having a pocket probing depth (PPD) of $6 mm. Sample teeth were randomly selected before intervention to assure equal distribution throughout the mouth. At the four sample sites supragingival plaque was gently removed, teeth were air-dried and isolated with cotton rolls. One sterile paper point (ISO45, Roeko Dental, Langenau, Germany) was inserted for 10 seconds in each site and all paper points were removed and pooled in one transport tube. Samples were stored at 220uC until further use. Bacterial genomic DNA was isolated and purified with the QiaAmp Mini DNA-Isolation Kit (Qiagen, Hilden, Germany). The protocol followed the manufacturer's instructions with minor modifications, i.e., a pretreatment step with proteinase K (20 mg/ ml) at 56uC was carried out overnight. After isolation, the purified DNA was eluted in 200 ml of elution buffer. Quality and purity of the isolated genomic DNA was confirmed by agarose gelelectrophoresis and spectrophotometry on the NanoDrop 2000 device (Fisher Scientific, Schwerte, Germany). DNA concentration was estimated with the Qubit 2.0 instrument applying the Qubit dsDNA HS Assay (Life Technologies, Invitrogen division, Darmstadt, Germany).

16S Primers and Amplicon Library Generation
PCR amplification of the 16S rRNA hypervariable region V6 was performed with a pool of six degenerated forward and six degenerated reverse primers targeting bacteria and archaea as described by Huber at al. [8]. The 59-ends of the forward primers were fused with the A-Adaptor plus key sequence, whereas the reverse primers were fused with the truncated Pi-adapter sequence (trP1), respectively. The primers were diluted in molecular grade laboratory water and pooled equimolar. For amplicon library preparation 4 ng of each genomic DNA, 1 U Platinum Taq DNA polymerase High Fidelity, 5 mM dNTPs, 2 mM MgCl 2 (all Life Technologies), and 10 pmol primer-mix were used per 25 ml amplification reaction. The PCR conditions were as follows: 94uC for 3 min, followed by 36 cycles of 94uC for 15 s, 58uC for 15 s, 68uC for 10 s, and a final elongation step of 68uC for 30 s. Amplicon product purification was performed via gel-electrophoresis on a 1.5% Tris Borat EDTA agarose gel [9] stained with SybrGold (Life Technologies) using the peqGold S gel extraction kit (Peqlab Biotechnologie GmbH, Erlangen, Germany). Control of purification and determination of exact fragment sizes was carried out with the Caliper GX system using the HT DNA High Sensitivity LabChip Kit (Caliper Life Sciences GmbH, Mainz, Germany). Amplicon library concentration was estimated with the Qubit 2.0 instrument using the Qubit dsDNA HS assay (Life Technologies). From the concentration and the average size of each amplicon library, the amount of DNA fragments per micro liter was calculated and libraries were diluted to 2.8610 8 DNA molecules per micro liter prior to clonal amplification.

Emulsion PCR and Sequencing
The emulsion PCR was carried out applying the Ion XPress Template kit V2.0 (Life Technologies) as described in the appropriate user Guide

Sequence Analysis
Raw sequencing reads were consecutively checked for different quality criteria. Low quality reads were removed as follows: i) reads with no match against the trP1 adapter allowing for up to three errors; ii) reads not matching the degenerated PCR primers with at most two errors; iii) sequences containing ambiguously called bases (N); and iv) reads with an average quality score below the 5 th percentile of all reads (95 percent confidence interval). Correctly detected primer sequences were trimmed off during this procedure. Primer matching and trimming were performed by extensions to the Conveyor workflow engine [10] in combination with Cutadapt [11], requiring a minimum primer overlap of 10 bases and penalizing initial gaps in the semi-global alignment. Derived high quality sequences were screened for artificial chimeric formations using the UCHIME algorithm [12]. In order to de-noise sequencing data, reads were further applied to the Single Linkage Preclustering (SLP) method [13] implemented in mothur [14]. In this process, sequences having a pairwise distance less than 2% were clustered using a modified single linkage algorithm and merged while retaining cluster sizes as counts. Following analyses were performed using only high quality reads de-noised by SLP.
Subsequently, operational taxonomic unit (OTU) analysis was done on a clustering basis for each sample individually. Reads were clustered at 97% identity threshold using the complete linkage clustering algorithm implemented in ESPRIT [15]. To meet the known problem of species overestimation [16], reads from singleton OTUs were excluded based on results of another upstream clustering, as recommended [17]. Rarefaction curves and the species richness estimators ACE [18] were calculated by ESPRIT. The abundance accounting Shannon and Simpson diversity indices [19] were computed using the statistical software suite R v 2.9.10 [20] and the vegan R-package [21]. Finally, assignment of taxonomic annotation to high quality reads was achieved using the GAST [22] pipeline with default parameters in combination with the SILVA rRNA database (version 102) [23]. In brief, GAST compares each sequence to a reference database, selects the closest reference(s), and assigns the taxonomy based on the lowest common ancestor for a two-third majority of the references. A schematic overview of the experimental work-flow and applied bioinformatic procedure is given in Figure S1.

Results and Discussion
Next-generation sequencing technologies extended our understanding of microbial communities and were applied in a vast number of metagenomic studies. The recently introduced Ion Torrent PGM constitutes a new type of next-generation benchtop sequencing platform. Main advantages of this machine are the small size, low acquisition price, and fast turnaround time. Here, we used the PGM to analyze subgingival plaque samples of four patients for shifts in the microbial community in response to treatment by SRP and SRP with adjunctive antibiotics. Demographic data, baseline and post intervention clinical parameters (pocket probing depth, bleeding on probing (BOP), and plaque index (PI)) are listed in Table S1. All clinical parameters improved after intervention in the control and experimental group markedly, i.e. shallow pockets (, = 3 mm) increased whereas initially deep sites decreased. No obvious differences between both groups could be observed. Due to limited sample size and short follow-up statistical inference was not performed. However, it is documented in numerous studies and systematic reviews that the administration of adjunctive antibiotics result in further improvements of clinical parameters in comparison to SRP treatment alone e.g. [6,7,24,25].
For the 16S rRNA gene amplicon sequencing we targeted the V6 hypervariable region, corresponding to positions 984-1047 (Escherichia coli). Although being a rather short region, it has been shown to enable successful discrimination between different phylotypes [26]. By using the Silva database it was demonstrated in an in-silico analysis that more than 97% of V6 tags could be unambiguously mapped at genus level [22]. In addition, V6 was not only the hypervariable 16S rRNA gene region used in the first NGS metagenomic amplicon study [27] but is still most frequently applied in environmental and medical microbiology projects e.g. [4,28]. Furthermore, amplicons of single hypervariable regions missing conserved intermediate stretches and short targets in general are less prone to form chimeras, therefore contributing to higher data quality [29]. In total, 4,892,600 sequence reads, ranging from 599,933 to 650,416 per sample and 314 chip were generated on the Ion Torrent PGM machine ( Table 1). Downstream sequence processing and quality filtering removed about 12.5% of the raw sequencing data. Most sequences were removed from further analysis due to untraceable sequencing or amplification primers (,8%) and significant low quality scores (,4.5%). Chimeric sequences were found only in very small quantities (less than 1% per sample). Sequences of each data set were further de-noised by SLP, by which between 9% and 14% of all reads were compensated for potential sequencing errors resulting in 5,210 to 8,462 unique tag sequences per sample ( Table 1).
Clustering of de-noised high quality reads at species-level, defined by 3% sequence difference, generated between 3,499 and 6,058 OTUs per specimen ( Table 2). Exclusion of singleton OTUs prior to computation of community metrics decreased OTU counts to about 50% while affecting only 0.5% of the sequencing data. Non-singleton OTUs ranged from 1,648 to 2,659 counts per sample. ACE species richness estimations were almost identical to the OTU numbers demonstrating a sufficient sampling effort ( Table 2). This is furthermore supported by rarefaction analysis where all curves almost reached saturated plateau phase ( Figure 1). No correlation between the number of OTUs and sampling time or treatment type could be observed. Our phylotype estimates were about a log order higher than those determined by cloning and traditional Sanger sequencing [2]. When analyzing pooled supragingival plaque samples from 89 healthy adults by 454 pyrosequencing, Keijser et al. reported 10,052 OTUs [4] at 3% difference. The first metagenomic study using Illumina's highthroughput sequencing technology by Lazarevic et al. [5] inferred about 8,000 different phylotypes at 3% OTU distance from pooled oral samples of three healthy adults. In contrast, Griffen et al. [30] found between 100 and 300 species OTUs in single individuals when analyzing subgingival samples of 29 subjects with chronic periodontitis by 454 pyrosequencing. However, in their study OTUs were not defined by hierarchical clustering but based on similarity to sequences of a curated oral 16S rDNA database [31]. This demonstrates that OTU profiles and metrics based thereupon are highly dependent on data characteristics and applied processing methods. However, community comparisons by species richness and diversity estimates are still meaningful if data were generated and analyzed in the same way. In our study, Shannon and Simpson indices showed a consistent pattern, i.e. increased diversity in all samples after therapy regardless of the treatment type, with sample Ex-Pre-1 representing the least (Shannon: 3.68; Simpson: 0.11) and Co-Post-2 the most diverse (Shannon: 4.91; Simpson: 0.01) community ( Table 2). As diversity also considers OTU abundances, our results reflect more evenly distributed abundance profiles after treatment, i.e. prevalent organisms of a Figure 1. Rarefaction analysis for each sample based on clustering of high quality de-noised reads. In total samples from four patients were analyzed two times, i.e., before (continuous lines) and two months after intervention (dashed lines). Two control group patients (Co) were treated by SRP alone whereas two experimental group patients (Ex) received SRP plus adjunctive antibiotics. Lines are colored according to the sample origin. doi:10.1371/journal.pone.0041606.g001 mature biofilm were suppressed by therapy inducing favorable growing conditions for rare and new species. This is also supported by Pielou's evenness values that were indicating more equally distributed OTU abundances among community members in post-treatment samples. It is interesting to note that Kumar at al. when comparing subjects with healthy implants and periimplantitis using 16S pyrosequencing found a significant increase in species diversity in the healthy subgroup [32].
Almost all community members in our samples originated from the bacteria domain. The most prevalent phyla in all samples were: Bacteroidetes, Firmicutes, Fusobacteria, Proteobacteria, Actinobacteria, Spirochaetes, and Synergistetes ( Table 1) comprising together more than 99.8% of all reads. Irrespective of treatment type a shift from Gram-negative Bacteroidetes to Gram-positive Firmicutes was observed after treatment. At genus level, sequencing reads could be assigned to 330 individual genera, of which 38 were present in all samples with a relative abundance of more than one percent in at least one sample (Figure 2). On average, in eight percent of all cases a successful classification at genus level was not possible, resulting in a classification at a higher taxonomic rank (data not shown).
Specimens taken prior to any treatment showed a highly similar community composition by means of taxonomic groups and abundance proportions. Bacteroidetes, Firmicutes, Fusobacteria, and Spirochaetes were the most prevalent phyla ( Table 3), which is in accordance with previous studies [2,4,5]. At genus level the predominant genera were Porphyromonas, Prevotella, Treponema, Fusobacterium, and Tannerella (Figure 2), occupying about 67% of the community sequences. Among these, three genera are of particular interest, as they include the species Porphyromonas gingivalis, Treponema denticola, and Tannerella forsythia, proposed to form the pathogenic 'red complex' consortium in periodontitis [33]. Furthermore, the genus Prevotella also includes several known periodontitis associated bacteria (Prevotella nigrescens, Prevotella intermedia, Prevotella melaninogenica) [33,34]. Other genera related to pathogenic organisms could also be identified, but either only in specific samples (e.g. Aggregatibacter) or at lower abundances (e.g. Parvimonas). For post-treatment samples the three most prevalent genera were Prevotella, Streptococcus, and Fusobacterium constituting about 37 percent of the sequences. Specimens of the experimental group showed markedly decreased abundances of the genera Treponema and Filifactor with frequencies dropping down to complete absence for Filifactor in one patient. In the control group these genera remained relatively stable. The genera Porphyromonas and Tannerella were suppressed in both groups. However, in the experimental group they were virtually absent ( Figure 1). Due to the limited read length of most next-generation sequencing platforms, taxonomic resolution down to species level is impossible. Therefore, we additionally screened all samples for P. gingivalis and T. forsythia using species-specific real-time PCR. A positive treatment effect indicated by a reduction of genomic equivalents of both species could be observed irrespectively of the treatment type (Table 4). However, in the experimental group P. gingivalis and T. forsythia were also suppressed to virtual absence, thereby corroborating our NGS data. Our observed bacterial community shifts are in accordance with Valenza et al. [35], who did a longitudinal capillary sequencing rDNA cloning study and reported that two months after treatment with SRP plus amoxicillin and metronidazole, disease associated genera were suppressed (e.g. Treponema, Porphyromonas, Filifactor, Tannerella). However, these therapy effects were reversible after twelve months. Based on cultivation and on species specific PCR assays the very same effect of systemic adjunctive antibiotics on P. gingivalis and T. forsythia was previously reported [24,25].
In our study several genera had increased abundance values in the control and experimental group after treatment, i.e. Prevotella, Selenomonas, Streptococcus, Actinomyces, and Rothia. Of them, Streptococcus, Actinomyces, and Rothia comprise known early colonizer species [36]. The finding that the antibiotic treatment seemed to have selectively eradicated some genera such as Porphyromonas while other genera such as Prevotella were not affected is difficult to explain. In principle the effects might be due to different recolonization dynamics or antibiotic susceptibility. Since taxonomic resolution down to species-level is impossible we were not able to determine if the sequencing reads counting for the genus Prevotella before and after treatment belonged to the same species. Prevotella loescheii is known to be an early colonizing partner of receptor-bearing streptococci (first colonizers) [36] therefore a shift from periodontitis associated Prevotella species towards the early colonizing species could be conceivable. It is known that the gingival crevicular concentration of metronidazole and amoxicillin is at least ten times the in-vitro minimal inhibitory concentration (MIC 90 ) of P. gingivalis, P. intermedia and P. nigrescens [37]. However, some Prevotella species might have developed active resistance mechanisms, e.g. Bernal et al. reported that about 26% of P. intermediaand P. nigrescens isolates are beta-lactamase positive and thereby at least not susceptible to amoxicillin [38]. In addition to the genera Selenomonas and Streptococcus, Valenza et al. [35] also observed for Prevotella an increase in abundance early after intervention. In contrast to them we were able to relate effects after treatment either to antibiotic or mechanical debridement therapy by inclusion of a control group. Whether the detection of selected periodontal pathogens is useful for treatment decision remains questionable. For example, by testing for A. actinomycetemcomitans, Fusobacterium nucleatum spp., P. gingivalis, P. intermedia, T. denticola, and T. forsythia Cionca at al. showed that it was of no predictive value in identifying subjects who would benefit from adjunctive antibiotic treatment [39]. In another recent study the authors suggested a correlation of these species with pocket probing depth rather than periodontal diagnosis [40]. Therefore, analyzing hypothesis free the complete microbial community by NGS might be the key in understanding a chronic disease like periodontitis, as the 'pathogen' might be a disturbed microbial community rather than single organisms [41]. If the in this study observed differences in community profiles and species diversity of periodontitis patients after treatment implies a diagnostic useful prognostic factor for treatment success warrants further experiments with larger sample size, longer follow-up, and a treatment failure subgroup. As performance characteristics like reading length and error rate of the PGM further improves, more specific, e.g. by sequencing longer stretches of adjacent hypervariable regions [31], and less error prone descriptions of the microbial communities will be possible [42].
In conclusion, this translational metagenomic study gives new insights into differences of the subgingival microflora of periodontitis patients for two different treatment regimes. The study results indicate that changes in community profiles and metrics are potentially more diagnostically predictive for treatment success than the detection of selected periodontal pathogens. Monitoring changes in the subgingival microflora with platforms like the PGM might finally pave the way towards routine personalized medicine. Figure S1 Schematic overview of the applied experimental and bioinformatic steps throughout the analysis. For a more detailed description of each individual step please see the Methods section. (PDF)

Supporting Information
Table S1 Demographic patient data and clinical parameters at baseline and post intervention. (PDF)