Towards Quantitative Microbiome Community Profiling Using Internal Standards

Abstract

An inherent issue in high-throughput rRNA gene tag sequencing microbiome surveys is that they provide compositional data in relative abundances. This often leads to spurious correlations, making the interpretation of relationships to biogeochemical rates challenging. To overcome this issue, we quantitatively estimated the abundance of microorganisms by spiking in known amounts of internal DNA standards. Using a 3-year sample set of diverse microbial communities from the Western Antarctica Peninsula, we demonstrated that the internal standard method yielded community profiles and taxon cooccurrence patterns substantially different from those derived using relative abundances. We found that the method provided results consistent with the traditional CHEMTAX analysis of pigments and total bacterial counts by flow cytometry. Using the internal standard method, we also showed that chloroplast 16S rRNA gene data in microbial surveys can be used to estimate abundances of certain eukaryotic phototrophs such as cryptophytes and diatoms. In Phaeocystis, scatter in the 16S/18S rRNA gene ratio may be explained by physiological adaptation to environmental conditions. We conclude that the internal standard method, when applied to rRNA gene microbial community profiling, is quantitative and that its application will substantially improve our understanding of microbial ecosystems.IMPORTANCE High-throughput-sequencing-based marine microbiome profiling is rapidly expanding and changing how we study the oceans. Although powerful, the technique is not fully quantitative; it provides taxon counts only in relative abundances. In order to address this issue, we present a method to quantitatively estimate microbial abundances per unit volume of seawater filtered by spiking known amounts of internal DNA standards into each sample. We validated this method by comparing the calculated abundances to other independent estimates, including chemical markers (pigments) and total bacterial cell counts by flow cytometry. The internal standard approach allows us to quantitatively estimate and compare marine microbial community profiles, with important implications for linking environmental microbiomes to quantitative processes such as metabolic and biogeochemical rates.

Keywords: amplicon sequencing; community profiling; internal standard; marine microbiome.

FIG 1

(A) Sampling location at the WAP across the gradients of coast to open ocean (y axis) and open water to ice edge (x axis). (B) Relative abundances of 18S rRNA gene from a test run with a 5:1-diluted internal standard using representative samples from the open ocean and the coast. Coastal 2A and 2B are duplicate samples from the same location. The internal standard 18S rRNA gene reads (in red) are a small portion of the total reads and are proportional to the dilution factor.

FIG 2

A 3-year data set of WAP community OTU counts (97% similarity) in RMP (upper panels) and QMP in rRNA gene copies milliliter⁻¹ (lower panels) for eukaryotes (A) and prokaryotes (B). The bar plots present QMP and RMP, with the top 20 OTUs associated with their taxonomic identifications (down to the finest identified level in the SILVA 128 database) and all other OTUs combined into one bin labeled “others” in gray. Stations are ranked in ascending total 16S or 18S rRNA gene counts. Environmental variables are plotted at the bottom: grid station or the approximate distance to shore (kilometers) in red, chlorophyll a concentration (milligrams meter⁻³) in green, and primary production (milligrams meter⁻³ day⁻¹) or bacterial production (leucine incorporation [picomoles liter⁻¹ hour⁻¹]) in blue. QMP captures significant variations across samples and is correlated with environmental variables.

FIG 3

Abundances of cryptophytes (A and B), diatoms (C and D), and Phaeocystis (E and F) estimated by 18S rRNA gene QMP (left panels) and RMP (right panels) compared to abundances estimated by CHEMTAX-HPLC pigment analyses. For significant correlations (P < 0.05), Pearson’s r² values for original data and log₁₀-transformed data (in parentheses) are shown at the bottom of the plot.

FIG 4

Comparison of total bacterial abundances estimated by FCM and 16S rRNA gene QMP. (A) The size of each data point represents the rrn effect, calculated as the averaged rrn of the top 20 classified OTUs in each sample. The color coding represents the particle association effect, calculated as the cumulative cell percentage (after rrn correction) of the top 10 most abundant particle-associated OTUs, identified as the >3-µm fraction bacteria reported by Delmont et al. (41). An exhaustive survey of the particle-associated OTUs is not feasible considering that significant portions of the prokaryotic OTUs are of unknown physiology or are even unclassified in the current rRNA database (SILVA 128). (B) FCM versus rrn and particle association effect-corrected 16S rRNA QMP in cells milliliter⁻¹. The four points inside the gray circle are likely outliers. After excluding these four points, data points are plotted on both linear and log scales.

FIG 5

Taxon cooccurrence matrices based on QMP (rRNA gene copies per milliliter) versus RMP (rRNA gene copy percentage). The top 24 most abundant prokaryotic (16S rRNA gene) OTUs and eukaryotic (18S rRNA gene) OTUs were used to construct the pairwise correlation matrices based on QMP (upper triangle) and RMP (lower triangle). Spearman’s ρ values for significant correlations (P < 0.05) are presented as squares in the heat map. The bottom histograms display the distribution of negative (blue) and positive (red) ρ values for matrices based on QMP and RMP.

FIG 6

Comparison of phytoplankton abundances based on normalized genomic 18S and chloroplast 16S rRNA genes for cryptophytes (A), Fragilariopsis (diatom) (B), Corethron (diatom) (C), and Proboscia (diatom) (D) and for Phaeocystis in relationship to the north-south geographic gradient (color-coded by Palmer LTER line number as an index for open water to ice edge gradient) (E) and Phaeocystis subclade ratios calculated as ln(OTU138/OTU2) (F). OTU2 and OTU138 are the top two Phaeocystis OTUs, comprising 87.27% and 12.70% of the total Phaeocystis 18S rRNA gene counts, respectively. The slope and intercept values are presented in Table S3 in the supplemental material.