A Dissolved Oxygen Threshold for Shifts in Bacterial Community Structure in a Seasonally Hypoxic Estuary

Abstract

Pelagic ecosystems can become depleted of dissolved oxygen as a result of both natural processes and anthropogenic effects. As dissolved oxygen concentration decreases, energy shifts from macrofauna to microorganisms, which persist in these hypoxic zones. Oxygen-limited regions are rapidly expanding globally; however, patterns of microbial communities associated with dissolved oxygen gradients are not yet well understood. To assess the effects of decreasing dissolved oxygen on bacteria, we examined shifts in bacterial community structure over space and time in Hood Canal, Washington, USA−a glacial fjord-like water body that experiences seasonal low dissolved oxygen levels known to be detrimental to fish and other marine organisms. We found a strong negative association between bacterial richness and dissolved oxygen. Bacterial community composition across all samples was also strongly associated with the dissolved oxygen gradient, and significant changes in bacterial community composition occurred at a dissolved oxygen concentration between 5.18 and 7.12 mg O2 L-1. This threshold value of dissolved oxygen is higher than classic definitions of hypoxia (<2.0 mg O2 L-1), suggesting that changes in bacterial communities may precede the detrimental effects on ecologically and economically important macrofauna. Furthermore, bacterial taxa responsible for driving whole community changes across the oxygen gradient are commonly detected in other oxygen-stressed ecosystems, suggesting that the patterns we uncovered in Hood Canal may be relevant in other low oxygen ecosystems.

About the Authors:

Rachel L. Spietz

* E-mail: [email protected]

Current address: School of Oceanography, University of Washington, Seattle, Washington, United States of America

Affiliation: School of Aquatic and Fishery Sciences, University of Washington, Seattle, Washington, United States of America

Cheryl M. Williams

Affiliation: School of Oceanography, University of Washington, Seattle, Washington, United States of America

Gabrielle Rocap

Affiliation: School of Oceanography, University of Washington, Seattle, Washington, United States of America

M. Claire Horner-Devine

Affiliation: School of Aquatic and Fishery Sciences, University of Washington, Seattle, Washington, United States of America

Introduction

Globally, over a quarter-million square kilometers of marine ecosystems are threatened by low dissolved oxygen (DO) levels, or hypoxia, which can result in the exclusion or death of resident macroorganisms creating so-called ‘dead-zones’ [1]. The number of hypoxic coastal areas has grown at an exponential rate of 5.54% year-1 [2], largely due to increased eutrophication [3,4]. Seasonal hypoxia is detrimental to coastal ecosystems and is known to have negative effects on fish [5], crustacea [6], benthic invertebrate communities [7], and trophic interactions and energy flow [8]. However, Bacteria and Archaea remain active in hypoxic waters and perform important functions in mineralization of organic matter and other biogeochemical cycling [1,9].

Extensive research has gone into understanding the thresholds of hypoxia for management strategies aimed at avoiding catastrophic collapses of ecosystems. In a review of over 800 published experiments, Vaquer-Sunyer and Duarte [2] report that published thresholds of hypoxia span a broad range from 0.29 mg O2 L-1 to 4 mg O2 L-1, but on average, sublethal effects on marine benthic macroorganisms occur when DO concentration drops below 2.61 ± 0.17 mg O2 L-1 and lethal effects occur at 2.05 ± 0.09 mg O2 L-1, with many taxa of Crustacea and Fishes demonstrating negative effects of hypoxia near 4 mg O2 L-1. It is worth noting, however, that thresholds for hypoxia likely differ between coastal open ocean hypoxic zones, as the organisms in each may be adapted to different durations of hypoxia [10]. Lacking from our current knowledge is whether a threshold of hypoxia for bacterial communities exists, and, if so, at what DO concentration do bacterial communities experience taxonomic and functional shifts.

Understanding the effect of hypoxia on lower trophic levels including plankton assemblages is crucial for elucidating whole ecosystem effects and interactions among trophic levels. Many microbial studies have focused on truly oxygen deficient zones where anaerobic metabolisms dominate including the Eastern Tropical South Pacific [11] and North Pacific [12], Arabian Sea [13,14], Cariaco Basin [15,16], Baltic Sea [17], Black Sea [18–20], and Saanich Inlet [21]. A recent meta-analysis examining global patterns of microbial community composition in oxygen-minimum zones found that particular bacterial and archaeal taxa, specifically Proteobacteria, Bacteroidetes, marine group A, Actinobacteria, and Planctomycetes, are common in suboxic waters and that these taxa tend to co-occur in a non-random pattern along the oxycline [22]. Only some of the members of these suboxic-associated phyla are known to use electron acceptors alternative to oxygen, which might suggest that mechanisms other than oxygen-depletion drive changes in the bacterial community. Though these previous studies have demonstrated that particular bacterial taxa persist in oxygen-minimum zones, the threshold concentration of DO at which bacterial communities change remains unclear. The bacterial community that exists at dysoxic (0.66–2.96 mg O2 L-1) and suboxic (0.03–0.66 mg O2 L-1) levels has the potential to become increasingly important in global biogeochemical cycles as oxygen-minimum zones in both the open and coastal ocean continue to expand with climate change [23], thus emphasizing the importance of defining a threshold DO concentration for Bacteria.

We examined spatial and temporal variation in bacterial communities in Hood Canal, WA, USA between April and October of 2007. Hood Canal is a long and narrow glacial fjord located 80 miles west of Seattle, Washington, USA (Fig 1A), that seasonally experiences periods of low and even completely depleted DO concentrations as a result of naturally occurring physical and hydrographical conditions [24]. This system has relatively small-scale, enclosed circulation patterns and rapidly changing abiotic conditions. Throughout the year, Hood Canal remains highly stratified, which limits the movement of oxygenated waters [25]. Furthermore, human activity over the past few decades has exacerbated nutrient loading particularly in the southern reaches of Hood Canal [26], in turn increasing the frequency and duration of recent hypoxia events [27]. Fish and macroinvertebrates in Hood Canal are directly affected by hypoxia through mortality or distributional shifts [28,29]. Particularly significant fish kill events near the southern reaches of the canal have occurred in the Falls of 2002, 2003, 2004, and 2006 as well as in the Spring of 2006 when the oxygen-deprived deep-water mass shoaled to the surface as a result of a combination of southerly winds and fresh-water intrusion from the Puget Sound at the North end of the canal [30]. Additionally, mesozooplankton assemblages in Hood Canal demonstrated altered composition and behavior when DO dropped below ~1.5 mg O2 L-1 [31]. These shifts have implications for predator-prey dynamics and trophic energy transfer and suggest that hypoxia may affect the entire biologic community in Hood Canal. We examined the bacterial communities associated with the seasonal variability in abiotic factors in Hood Canal, including depleted DO concentration. We hypothesized that bacterial communities in Hood Canal demonstrate non-random changes in composition and taxa richness related to variation in abiotic factors, and that these community patterns are associated particularly with changing levels of DO.

[Figure omitted. See PDF.]

Fig 1. Geography and biogeochemical conditions in Hood Canal, WA, USA in 2007.

(A) Map of sampling stations within Hood Canal, WA, USA. (B) Contour plots showing range of dissolved oxygen (DO) and chlorophyll a concentrations at Hama Hama and Sister’s Point stations in Hood Canal, WA. Data for high-resolution depth profiles were collected by Oceanic Remote Chemical Analyzer buoys maintained by the Northwest Association of Networked Ocean Observing Systems (http://www.nanoos.org). Note changes in scale on both the x- and y-axes.

https://doi.org/10.1371/journal.pone.0135731.g001

Methods

Sample collection

High resolution depth profiles of DO and chlorophyll a for two mid Hood Canal, WA stations, Hama Hama and Sister’s Point, were obtained via Oceanic Remote Chemical Analyzer (ORCA) buoys maintained by the Northwest Association of Networked Ocean Observing Systems (http://www.nanoos.org) [32]. Water samples were collected during April, June, and October of 2007 into Niskin bottles mounted on a rosette, equipped with a CTD (conductivity, temperature, and depth) sampler (Sea-Bird Electronics, Bellevue, WA, USA) that measured in-situ temperature and salinity. April and October samples were collected from two stations near the middle of Hood Canal: Hama Hama and Sister’s Point. During the June sample collection, two additional stations (Bangor and Lynch Cove) were added in order to survey a full north to south transect of the Canal. At each sampling station and time, samples were collected both at five meters depth (herein referred to as surface samples) as well as ten meters above the bottom substrate (herein referred to as deep samples). These so-called deep samples varied in depth from 140−145 m at Hama Hama, 125 m at Bangor, 40−50 m at Sister’s Point, and 12 m at Lynch Cove, as the depth of the water column decreases near the southern reaches of the canal. The timing of the collection points were intended to coincide with the late spring algal bloom (April), subsequent die-off (June), and a smaller autumn algal bloom (October), which affect DO concentration and nutrient availability (S1 Table). A field permit was not required for the sampling in this study, as there are no protections in place for the waters or organisms that would be impacted by our sampling and sampling did not affect any endangered or protected species.

Immediately upon retrieval of the Niskin bottles, DO samples were collected using surgical tubing into calibrated glass bottles with ground glass stoppers and were analyzed using a modified Winkler titration [33]. The remaining water was used for additional measurements of biotic and abiotic parameters. For inorganic nutrients, 50 mL of water was passed through a 0.45-μm syringe filter and immediately frozen until being analyzed at the University of Washington’s Oceanography Marine Chemistry Lab for NO3- (detection limit 0.15 μM), NO2- (detection limit 0.02 μM), NH4+ (detection limit 0.12 μM), PO43- (detection limit 0.03 μM), and SiOH4 (detection limit 0.59 μM) following the Protocols for the Joint Global Ocean Flux Survey (JGOFS) Core Measurements (1994) [34]. To measure bacterial abundance, 5 mL of water was preserved in formaldehyde with a final concentration of 1% and stored at -80°C until 0.5–1.0 mL aliquots were filtered onto 0.22-μm polycarbonate filters and stained with 4’,6-diamidino-2-phenylindole (DAPI). The DAPI filters were viewed under a Nikon Eclipse 80i with UV light from an X-cite Series EXFO to count the number of cells present in fifteen randomly selected fields using NIS-Elements BR 3.0 software then averaged to quantify the total bacterial abundance. Principal components analyses and all other statistical analyses were conducted in R using the vegan package (http://cran.r-project.org/web/packages/vegan/index.html), unless noted otherwise. Metadata are provided in S1 Table.

DNA extraction and processing

For bacterial community analyses, 300 mL of water was filtered onto 0.22-μm Supor filters and preserved with 3 mL of solution of 0.5 M NaCl, 10 mM Tris, and 100 mM EDTA, then immediately frozen at -80°C. A Qiagen DNeasy Tissure Kit was used to extract DNA from the filters following the manufacturer’s instructions for extraction from gram-negative Bacteria (Qiagen Valencia, CA). The DNA extracts were sent to MBL’s Keck Facility for PCR amplification and sequencing on the GS-FLX Titanium 454 platform (see Sogin, et al., 2006 for previously described methods [35]). Briefly, the V6 region of the 16S rRNA gene was amplified via polymerase chain reaction (PCR) using a pool of bacterial specific primers, excluding amplification of Archaea (see Huse, et al., 2010 for description of V6 primer pool [36]). Sequencing was performed in conjunction with the International Census of Marine Microbes (ICoMM) project (http://icomm.mbl.edu), which uses a standardized 454 pyrosequencing pipeline to sequence marine microbial samples from across the globe [37]. All sequences were deposited in GenBank Sequence Read Archives (www.ncbi.nlm.nih.gov) under accession number SRP001226.

As part of the ICoMM project, the sequence data were initially analyzed using a standardized pipeline [35,36,38]. Sequences were pre-clustered using a single-linkage algorithm to smooth sequencing errors and reduce noise. Sequences were then clustered using average-linkage into operational taxonomic units (OTUs) based on 97% sequence identity. Prior to community analyses, we removed sequences that were highly similar to chloroplast sequences.

Statistical analyses

Abiotic factors, including NO3-, NO2-, NH4+, PO43-, SiO42-, salinity, temperature, and DO concentration, were log-transformed prior to analyses to meet assumptions of normality and homogenization of variance. Euclidean distance between samples was calculated for environmental variables [39]. Principal components analysis (PCA) was performed to assess dominant trends in environmental variables between samples across depths, sites, and seasons. A Monte Carlo randomization test was used to test the significance of each principal component at an alpha level of 0.05.

To minimize the effect of sequencing effort on richness estimates, all samples were randomly subsampled down to 3663 sequences, which was the number of sequences of the sample with the fewest number of sequences (Table 1), the richness measures were calculated, then the process was repeated and averaged over 1000 iterations. Observed taxa richness, or the number of OTUs identified at the pre-defined sequence identity, as well as the Chao richness estimator were calculated using the subsampling method and averaged over 1000 iterations (Table 1). We used multiple regression modeling with a requirement of α < 0.05 for a variable to enter the model to determine which abiotic factors best described the changes in taxa richness, evenness (Pielou’s J), and proportion of functional groups across all samples.

[Figure omitted. See PDF.]

Table 1. Number of sequences obtained in each sample, as well as the observed and Chao estimates of taxa richness and Pielou’s J index of taxa evenness for all Hood Canal samples based on an operational taxonomic unit cutoff of 97% sequence identity.

https://doi.org/10.1371/journal.pone.0135731.t001

We used multivariate statistical analyses to identify important drivers of the relationship between bacterial community composition and environmental factors. We selected a Sorensen abundance-based dissimilarity measure because it incorporates probability of detection of taxa shared between two samples to account for under-sampling, which is common when studying microbial communities [40]. Similarly to the alpha diversity calculations, the bacterial community matrix was randomly subsampled to the number of sequences of the smallest sample then the Sorensen abundance-based similarity was calculated between samples to account for sequencing effort bias [41,42] and averaged over 1000 iterations.

Analysis of similarity (ANOSIM) tests examined the effects of factors including depth, site, season, or DO (high or low) on bacterial community composition of each sample and significance was tested using 999 random permutations of the community composition dataset. A non-parametric iterative BIO-ENV analysis was employed to select the set of abiotic variables most correlated with changes in bacterial community composition [43].

Threshold Indicator Taxa Analysis (TITAN) detects changes in taxa distributions along an environmental gradient that exists over space and/or time [44]. The analysis assesses the synchrony among taxa change points as evidence for community thresholds. Prior to TITAN, the sixteen Hood Canal samples were randomly subsampled to the smallest number of sequences present in a single sample (3663 sequences). In order to accurately calculate an indicator value (IndVal) for each taxon, the taxon must be present in three or more samples. Therefore, taxa with narrow distributions (present in less than three samples) were eliminated from each subsampled dataset. Finally, all bacterial abundance data were log-transformed to down-weight highly abundant types. TITAN was run in R using a minimum split of three samples, 100 random permutations to calculate IndVal p-values, and 100 random permutations for the bootstrapping calculation. Subsampling and TITAN analysis was repeated ten times to confirm that stochastic changes in sampling do not affect the results of TITAN.

To validate the DO threshold identified by TITAN, we used a non-parametric ANOSIM test grouping samples either by depth or DO level split at either 4 mg O2 L-1 or 6 mg O2 L-1. The R-values from ANOSIM, which describes how well samples group by factor based on community composition, were compared to determine the best factor for grouping.

Finally, hierarchical clustering using average linkage was used to cluster samples based on dissimilarity of community composition quantified by the Sorensen abundance-based estimator to identify the DO concentration at which the community composition changed significantly. All three analyses, TITAN, ANOSIM, and clustering were performed in R using scripts from Baker and King [44]; vegan, cluster, and pvclust packages and scripts written for this study, which are available upon request.

A similarity percentage (SIMPER) analysis was used to identify bacterial taxa contributing most to the similarity within groups of samples defined by DO and thus identify taxa strongly associated with either high or low DO [45]. SIMPER analysis identifies taxa that contribute the most to the similarity in community composition between samples of the same group, or similarly, taxa that contribute the most to dissimilarity in community composition across groups of samples. Thus, the taxa identified by SIMPER are significant indicators of a particular factor, in this case low DO. Dissolved oxygen is usually higher in surface waters than in deep waters, so to separate a depth effect on taxa distribution from an oxygen effect, the SIMPER analysis was run a second time using sampling depth as a factor rather than oxygen level. To test the significance of an oxygen effect over a depth effect, a generalized linear model was run using an intercept of 0 and a slope of 1, which would imply a null hypothesis that there is no difference between a depth or DO effect. Studentized residuals were calculated for each taxon, and a Bonferroni correction was applied to calculate a 95% confidence interval because of the high number of taxa or data points. Taxa that demonstrated a significant effect of DO level rather than sampling depth were those that fell above the 95% confidence interval, and thus had a studentized residual greater than 4.2169. For taxa identified to be significant indicators of high or low DO, their taxonomy was further resolved by using the V6-16S rRNA gene sequences to perform both a blastn search against the NCBI nr/nt database in (http://blast.ncbi.nlm.nih.gov/Blast) and also by aligning the V6-16S rRNA gene sequences against the arb-silva database using the online SINA v1.2.11 aligner (http://www.arb-silva.de/aligner/) [46]. All alignments were performed in June and July of 2014.

Results

Environmental context

Using high-resolution measurements of DO and chlorophyll a (Fig 1B), it was apparent that patterns of hypoxia in 2007 were not as extreme as in past years when low DO had caused fish kills. Yet, our samples captured a wide range of abiotic conditions and likely represent typical seasonal variability in Hood Canal. Water samples were collected from Hama Hama and Sister’s Point (Fig 1A) in April, directly following the major spring phytoplankton bloom as indicated by the heightened level of chlorophyll a in the surface water layer at both stations, marked by high DO concentrations (12.32−13.18 mg O2 L-1) and low nutrient concentrations in the surface waters (0.03−0.73 μM NO3-) and low DO concentrations (3.59−3.94 mg O2 L-1) and high nutrient concentrations in deep waters (Fig 1B). June samples were collected at the end of the phytoplankton bloom, and thus the surface waters were characterized by slightly lower DO (8.66−10.25 mg O2 L-1) and nutrient concentrations than in April. At the same time, the deep waters contained low DO (2.95−7.12 mg O2 L-1) and high nutrient concentrations. October samples were collected after a lesser fall phytoplankton bloom. Dissolved oxygen and nutrient concentrations in October were in the mid to high range, relative to previous months, in both water depths, which were less stratified than earlier in the year (S2 Fig and S1 Table). Bacterial abundances across all stations and time points ranged from 2.16 × 105 cells mL-1 to 2.70 × 106 cells mL-1, with an average abundance of 8.46 × 105 cells mL-1 across all samples.

Using the environmental variables measured at the time of water sample collection, a principal components analysis (PCA) was conducted to visualize abiotic factors driving the differences among samples (Fig 2). Nitrate, DO, salinity, and phosphate were most correlated with PCA axis one, which explained 57% of the total measured abiotic variation among samples. Samples appear to fall into two groups along axis one, mostly described by depth. An additional 19% of the total variation in abiotic parameters was explained by PCA axis two, which was most correlated with ammonium concentration.

[Figure omitted. See PDF.]

Fig 2. Principal components analysis of environmental characteristics of the sixteen samples ordinated based on Euclidean distance calculated from environmental factors.

Points are coded by water depth as either deep (filled symbols) or surface (open symbols) and by season as April (square), June (circle), or October (triangle). Nitrate, phosphate, salinity, and dissolved oxygen best explain the separation of samples along axis one, and ammonium alone explains variation along axis two based on the variable loadings on each principal component.

https://doi.org/10.1371/journal.pone.0135731.g002

Bacterial alpha diversity in Hood Canal, WA

Observed bacterial taxa richness, or the number of OTUs (defined by 97% sequence similarity clusters) recovered, in Hood Canal ranged from 155 OTUs to 590 OTUs per sample and was significantly higher in the deep-water samples than in the surface water samples by an average factor of two (F1,14 = 22.55; p < 0.01) (Table 1, Table 2, Fig 3). Neither station (horizontal space) nor season had a significant effect on bacterial richness (results not shown). By applying a forward-selection multiple regression model using all abiotic variables measured in this study, DO concentration alone best explained the variation in observed bacterial richness across samples (adjusted R2 = 0.738, p < 0.01, Fig 3). Taxonomic richness peaked at a DO concentration of 4 mg O2 L-1 (Fig 3).

[Figure omitted. See PDF.]

Fig 3. Bacterial alpha-diversity varies by DO.

A significant relationship exists between DO concentration and bacterial taxa richness. Points are coded by water depth as either deep (filled symbols) or surface (open symbols) and by season as April (square), June (circle), or October (triangle). Dashed line represents the best fit from the linear regression model of bacterial richness as a function of dissolved oxygen (adjusted R2 = 0.7379, p < 0.01).

https://doi.org/10.1371/journal.pone.0135731.g003

[Figure omitted. See PDF.]

Table 2. Continuous abiotic variables affect the diversity and composition of bacterial communities in Hood Canal, WA, USA (α = 0.05 for all tests).

The strength of richness and evenness, measured by Pielou’s J, models was assessed by the adjusted R2 value from multiple linear regressions, and the strength of the composition models was assessed by the Spearman’s rho correlation coefficient from BIO-ENV analyses.

https://doi.org/10.1371/journal.pone.0135731.t002

Across all samples, over 57% of the sequences clustered into only 42 unique OTUs (defined by 97% sequence identity), while less than 2% of the sequences clustered into over 2300 unique OTUs. The abundant OTUs were typically widespread across all the samples, whereas each of the less abundant OTUs was detected in one or two samples. Uneven distribution of individuals into taxa is not unusual in bacterial communities [35], which we further explored by calculating an index of taxon evenness, Pielou’s J (Table 1) that provided an additional window into bacterial community structure. Evenness of taxa was higher in deep waters than in the surface water (F1,14 = 9.517, p < 0.01), while neither season nor sampling site had an effect on the evenness of the bacterial community (Table 2). Using multiple regression modeling, taxa evenness was best explained by PO43- (adjusted R2 = 0.524, p < 0.01, Table 2).

Although an uneven distribution of individuals into taxa is common in bacterial communities [35], changes in the relative proportion of the dominant classes provide insight into factors impacting bacterial community structure (Fig 4). Sampling depth explained a significant amount of the variation in relative proportion of the fifteen most abundant classes across samples (G-test; G = 37.585, p < 0.01), whereas season or horizontal space had no significant effect (alpha < 0.05). Of the most abundant bacterial classes, the majority of sequences recovered from all samples belonged to three classes: Flavobacteria, Alphaproteobacteria, and Gammaproteobacteria. In surface waters these classes comprised 88.5% of sequences, whereas in deep waters they comprised only 68.7% of sequences, which demonstrates higher taxa evenness in deep-water samples. The proportions of Alphaproteobacteria and Gammaproteobacteria were not significantly different between surface and deep waters; however, the proportion of Flavobacterial sequences was significantly higher in the surface water samples compared to deep-water (F1,14 = 21.463; p < 0.01). The majority of Flavobacteria detected in the Hood Canal surface samples are members of the Flavobacteriaceae family, which are known to be highly associated with algal particles and are important players in organic matter degradation as part of the microbial loop [47,48]. The high relative abundance of Flavobacteriaceae in the surface water samples of Hood Canal reflect the considerable algal bloom dynamics in the estuary and a recent pulse of organic matter particles into the euphotic zone. The five most abundant taxa in each sample included common, ubiquitous marine bacterial taxa such as those belonging to SAR11, SAR86, Flavobacteriaceae, Thalassobacter, SAR406, SAR324, and Rhodobacteriacaea. However, as DO concentration decreased across samples, the relative abundances of the dominant types changed as well (Fig 4).

[Figure omitted. See PDF.]

Fig 4. Relative proportions of the most abundant classes change across a dissolved oxygen gradient and are significantly different when grouped by sampling depth.

The samples are ordered by dissolved oxygen concentration (solid line). Abbreviations: Hama Hama (HH), Sister’s Point (SP), Bangor (BA), Lynch Cove (LC), April (APR), June (JUN), October (OCT), Surface (S), Deep (D). Classes that comprised less than 1% of the sample’s total community were condensed into the “other” category.

https://doi.org/10.1371/journal.pone.0135731.g004

Patterns of functionally relevant taxa

In order to evaluate potential functional changes in the bacterial communities, we examined changes in the relative abundances of several groups of bacteria in our dataset known to be important in specific biogeochemical processes (Table 3). Two groups of ammonia oxidizing bacteria (AOB), Nitrosomonas and Nitrosococcus, were detected in Hood Canal. The relative abundance of AOB was higher in deep water than surface water (F1,14 = 18.98, p < 0.01) and also demonstrated a significant negative relationship with DO concentration (adjusted R2 = 0.5765, p < 0.01). Given that we used bacterial-specific 16S rRNA primers, we did not detect ammonia-oxidizing Archaea, which are also known to inhabit Hood Canal [9]. The relative abundance of two groups of nitrite-oxidizing bacteria (NOB) in Hood Canal, Nitrospira and Nitrospina, was also significantly higher in the deep water than in the surface water (F1,14 = 24.9, p < 0.01). Additionally, the relative abundance of NOB was not only strongly associated with a decrease in DO but also with decreased NH4+ concentration and temperature (adjusted R2 = 0.9315, p < 0.01). There were no significant patterns of increased Cyanobacteria abundance by water depth, sampling season, or sample site; however, the relative abundance of this photosynthetic phylum was strongly associated with higher temperatures (adjusted R2 = 0.8812, p < 0.01). The relative abundance of methylotrophic bacteria in Hood Canal, including Methylobacteriaceae, Methylococcaceae, Methylocystaceae, and Methylophilaceae, also demonstrated no significant patterns across water depths, season, or sampling site but their abundances did increase with decreasing DO and higher temperature (adjusted R2 = 0.5984, p < 0.01).

[Figure omitted. See PDF.]

Table 3. Abiotic factors affect the average proportion of various bacterial functional groups.

A forward selection stepwise multiple linear regression modeling approach was employed with a requirement of p < 0.05 to enter the model.

https://doi.org/10.1371/journal.pone.0135731.t003

A bacterial community threshold for DO

Change in bacterial community composition (BCC) was examined across all possible combinations of abiotic variables measured within this study, and DO concentration alone was most correlated with dissimilarities in the composition of bacterial communities among samples (Spearman’s rho = 0.7320, Table 2). Further, we examined the relevance of a threshold of DO concentration for bacterial communities in Hood Canal by applying three separate statistical tests: 1) a parametric assessment of the community data using hierarchical clustering, 2) a non-parametric assessment of the community data using an ANOSIM, and 3) an analysis of community thresholds determined by patterns of individual taxa using TITAN. Using a parametric clustering approach the BCC samples were plotted on a dendrogram based on the dissimilarity in BCC between pairs of samples. The clustering dendrogram reveals two distinct clusters of samples based solely on community composition data without any a priori knowledge of environmental factors (Fig 5A). Most of the surface water samples fell in the high DO cluster, except for the surface sample collected in October, which clustered with the lower DO samples as expected given the DO level at that sample was below 5.18 mg O2 L-1, the level that defined the low DO cluster. Similarly, all deep samples clustered with each other except for the deep sample collected in June from Bangor when DO concentrations were higher than 7.12 mg L-1. The two clusters suggest that DO is a strong factor in shaping BCC.

[Figure omitted. See PDF.]

Fig 5. Both (A) parametric and (B) non-parametric analyses suggest that the bacterial community composition shifts between 5.18 and 7.12 mg O2 L-1.

(A) Dendrogram of parametric hierarchical clustering approach using Sorensen-abundance based community composition dissimilarity. Numbers following Sample ID show the DO concentration in mg L-1. (B) Non-metric multidimensional scaling ordination shows bacterial community composition of samples are more similar to those from the same water depth, season, or DO concentration. Position of samples is calculated using the Sorensen abundance-based similarity index of bacterial community composition and the absolute dissimilarities between samples is condensed onto two-dimensions with a resulting stress value of 0.0679. Points are coded by water depth as either deep (filled symbols) or surface (open symbols) and by season as April (square), June (circle), or October (triangle). The size of the sample point is proportional to DO content with the highest DO represented by the largest sample point.

https://doi.org/10.1371/journal.pone.0135731.g005

In the second analysis, the samples were again plotted based on the pairwise dissimilarities in BCC using non-parametric multidimensional scaling (NDMS). The NMDS further supports that two groups of samples exist and that these two groups are delineated by DO, rather than depth alone (Fig 5B). An analysis of similarity (ANOSIM) confirmed that depth and DO are significant descriptors of the two groups of samples, but the association of DO with community composition is stronger than the association with depth (ANOSIM: RDO = 0.827, pDO < 0.01; Rdepth = 0.6144, pdepth < 0.01).

Lastly, we examined taxon-specific responses to the DO gradient to identify a threshold concentration of DO at which both positively and negatively-associated taxa change, thus giving a community-wide change point. TITAN identified a threshold value of 6.15 mg O2 L-1 for both positively and negatively associated taxa (Fig 6). This implies that a shift in BCC occurs at 6.15 mg O2 L-1 (95% CI [3.62, 8.78]) from a community dominated by taxa that are positively-associated with DO to a community dominated by negatively-associated taxa. Though the confidence interval surrounding the community change-point predicted by TITAN is broad, the other two independent statistical tests (hierarchical clustering and NMDS) also suggest a community change point existing between 5.18 mg O2 L-1 and 7.12 mg O2 L-1.

[Figure omitted. See PDF.]

Fig 6. Bacterial community threshold analysis for DO in Hood Canal.

Taxa-Z plot from TITAN showing the cumulative change points of taxa abundances across the DO gradient. The peak in change-point is considered an approximate community threshold for DO concentration. Z- are taxa that have a negative association with DO and Z+ are taxa that have a positive association with DO.

https://doi.org/10.1371/journal.pone.0135731.g006

Dominant taxa in Hood Canal, WA and indicator taxa of low DO

The relative proportions of the fifteen most abundant bacterial classes changed significantly as DO dropped below 6 mg O2 L-1 (Fig 4, G-test; G = 44.564, p < 0.01). Shifts in the most abundant bacterial classes prompted further investigation into the individual taxa that may be representative of changing DO conditions. By comparing the abundances of taxa across our samples using pyrosequencing data, we identified taxa that may be indicators of low DO in Hood Canal. We identified 23 taxa that had significantly higher DO SIMPER values than depth values at an alpha level of 0.05 (Fig 7A and 7B), suggesting that these taxa drive dissimilarity between high and low DO communities and can thus be considered indicators of different DO conditions. Only two taxa identified from SIMPER analysis were indicators of high DO waters in Hood Canal. One of these taxa, Alphaproteobacteria_03_61, was most similar to members of the family Rhodobacteriacaea that are commonly detected in marine surface waters, particularly during phytoplankton blooms (S2 Table). The remaining 21 taxa were indicators of low DO waters in Hood Canal. The majority of the low DO indicator taxa were common marine bacterial groups including Rhodobacteracaea, SAR11, SAR406, Nitrospina, Methylobacter, and SAR86, which are often found across a wide range of marine habitats (S2 Table). These indicator taxa were abundant OTUs throughout the entire dataset but have significant shifts with decreasing DO, suggesting that these OTUs may be specialized members of common, ubiquitous groups that have undergone niche differentiation.

[Figure omitted. See PDF.]

Fig 7. DO rather than depth affects taxa distributions and twenty-five taxa were identified as indicators of different DO conditions.

(A) Plot of each OTU’s contribution to the dissimilarity between high and low DO groupings of samples versus surface and deep grouping of samples. Red line shows one-to-one line representing equal contribution to both DO effect and depth effect. Taxa falling to the right of the red one-to-one line contribute more to the dissimilarity between high and low DO groupings than depth groupings suggesting a stronger effect of DO than depth. (B) Plot showing the studentized residuals calculated for each OTU to identify significant deviation from the null hypothesis that there is no difference between the OTU’s contribution to dissimilarity between DO and depth groupings. Each red horizontal line shows the cutoff for significance. Points above the upper red line represent OTUs that contribute significant dissimilarity between DO groupings thus are indicators of high or low DO, whereas points below the lower red line are significant indicators of a depth effect.

https://doi.org/10.1371/journal.pone.0135731.g007

Discussion

On a global scale, DO content in coastal waters has decreased in recent decades and is arguably becoming one of the most important indicators of aquatic ecosystem health [1]. We investigated changes in bacterial community diversity and composition along a DO gradient in Hood Canal, WA, and, for multiple indices of bacterial community structure, we found a role for DO that was stronger than that of the other abiotic factors assessed here. Further, we identified a threshold concentration of DO at which the bacterial community composition shifted significantly, suggesting that the structure of bacterial communities is intimately linked to DO.

A strong negative relationship between bacterial richness and DO was detected across all samples (Fig 3), which highlights a pattern consistent with other systems including the seasonally anoxic Canadian fjord, Saanich Inlet [21] and in the open ocean water column in the Eastern Tropical North Pacific [12]. However, when DO concentration falls below hypoxic levels in Saanich Inlet and the Eastern Tropical North Pacific, bacterial richness has been shown to decrease sharply. Expanding our sampling gradient to extremely low DO conditions periodically present in Hood Canal could elucidate such a potential unimodal relationship between DO and richness (Fig 3). Indirect effects of oxygen-depletion on the bacterial community at low DO could explain the increasing richness at low DO. For instance, energy is transferred from macrofauna to microbes as oxygen is depleted as oxygen-dependent organisms either avoid hypoxic areas or die [1,8,49]. This transfer of energy to the bacterial community as DO decreases could also allow for higher taxa richness exemplifying a common theme in ecology where increasing bioavailable energy allows for coexistence of multiple taxa resulting in greater biodiversity [50–52].

Dissolved oxygen played a strong role in shaping bacterial community composition, as well. Our data suggest a significant shift in bacterial community composition at a threshold DO concentration between 5.18 and 7.12 mg O2 L-1, which is not only higher than the published thresholds of 2–4 mg O2 L-1 for higher trophic-level organisms [2], but also much higher than the minimum concentration for aerobic metabolisms (< 96 ng O2 L-1) [53,54] as well as the maximum DO concentration for denitrification (160 μg O2 L-1) [55]. Therefore, the direct effects of oxygen on bacteria may not be the cause of the compositional shift we detected. Alternatively, indirect effects through trophic interactions including altered bacterial grazing or viral infection may favor a different bacterial community, as has been shown under other climate change scenarios such as increased temperature and UV radiation [56]. Additionally, abiotic factors not measured in this study that may co-vary with DO may be causing a compositional shift. For example, the composition and quality of dissolved organic matter can strongly influence bacterial community composition [57]. Nonetheless, evidence for significant changes to the bacterial community at threshold DO concentrations above the conventional definition of hypoxia suggests that impacts of decreasing DO have implications for ecosystem processes and health before sublethal and lethal effects on macrofauna are detected.

We identified 21 indicator taxa of low DO that were present across the majority of samples spanning a full range of DO, but their relative abundances change significantly at the threshold DO concentration. These low-DO indicator taxa were identified as members of common and ubiquitous marine bacteria including Roseobacter, SAR11, SAR202, SAR406, SAR86, OM190, and Nitrospina, which are often detected in coastal regions, estuaries, upwelling zones, and other pelagic oxygen minimum zones (S2 Table). Other OTUs from these common marine groups were also detected throughout our dataset, but were not identified as indicators of low DO. This may suggest that particular low-DO adapted ecotypes of common marine bacteria exist, highlighting the taxonomic and functional redundancy of bacterial communities [58]. Ecotypes of SAR11 are known to be distributed across chemical and physical gradients in marine surface waters suggesting that niche separation has led to the diversification within this group [59,60]. It is possible that similar mechanisms are responsible for shaping the distribution of indicator taxa we detected. Further research into the indicator taxa we detected here will elucidate the adaptations of these potential ecotypes that are responsible for providing a competitive advantage at low DO.

Our initial survey across a gradient of spatial and temporal conditions in Hood Canal has provided a snapshot into the complex relationship between bacterial community structure and DO. An obvious limitation of this study is that a few discrete samples do not capture fine-scale patterns or more extreme hypoxia conditions; yet with even limited sampling we were able to detect a significant shift in community composition at a threshold level of DO. Expanded sampling in both space and time will aid in refining this threshold concentration and the mechanisms for compositional shifts. The long residence time of water in Hood Canal [25] combined with the reduced grazing pressure on bacteria at depth or lower DO as a result of a shift in the distribution of higher trophic levels where DO is depleted [31] could also allow for the detection of more bacterial taxa, many of which may be in a dormant or inactive state. Dormancy can contribute to high taxonomic diversity and also act as a reservoir of functional diversity that can be resuscitated under favorable conditions [61].

As bacterial communities shift in diversity and composition, it is likely that the bacterially mediated nutrient and energy cycles also change, which has implications for changes in global cycles as oxygen minimum zones continue to expand. The hydrography of Hood Canal generates a wide gradient of abiotic factors including DO that change on an annual scale, which establishes a unique natural laboratory in which these questions can be explored. We are only beginning to understand how microbial communities are responding to decreasing DO, but further questions into the effects on individual functional groups of Bacteria will be crucial in predicting the global consequences of expanding oxygen minimum zones.

Supporting Information

[Figure omitted. See PDF.]

S1 Fig. Rarefaction curves of (A) All 16S sequences, (B) Alphaproteobacterial sequences only, (C) Gammaproteobacterial sequences only, and (D) Flavobacterial sequences only.

Each solid line represents the sequence accumulation curve for each sample in the dataset. The dotted line represents a 1:1 line, which is the maximum taxa accumulation per sequence. The dashed, vertical line in each plot represents the number of sequences to which each set was subsampled. Note that plot A is on a different scale from plots B-D.

https://doi.org/10.1371/journal.pone.0135731.s001

(TIF)

S2 Fig. Depth profiles of physical and chemical data recorded during sampling.

(A) April sampling at Hama Hama and Sister’s Point, (B) during June sampling at all stations, and (C) during October sampling at Hama Hama and Sister’s Point.

https://doi.org/10.1371/journal.pone.0135731.s002

(PDF)

S1 Table. Abiotic data collected to accompany each bacterial pyrosequencing sample and used in PCA and regression analyses.

https://doi.org/10.1371/journal.pone.0135731.s003

(PDF)

S2 Table. Indicator taxa of high or low DO conditions in Hood Canal, WA, as identified by SIMPER analysis, and their associated top BLAST matches on NCBI.

“Contrib.depth” and “Contrib.DO” describe the relative contribution of a particular taxon to the dissimilarity between sample groups defined by water depth or DO content, respectively. “Hi” or “Lo DO ave abund” describes the average abundance of a particular taxon in the samples grouped as either high or low DO, respectively.

https://doi.org/10.1371/journal.pone.0135731.s004

(XLSX)

Acknowledgments

The authors wish to thank Allan Devol, Jan Newton, and Wendi Ruef for supplying oxygen and chlorophyll a data from ORCA buoys in Hood Canal, which are operated by the University of Washington and maintained as part of the Northwest Association of Networked Ocean Observing Systems (NANOOS), part of the U.S. Integrated Ocean Observing System (IOOS). We thank Michael Dyen for assistance in sample processing and data collection and KathiJo Jankowski for helpful and insightful comments on the manuscript.

Author Contributions

Conceived and designed the experiments: CW GR MCHD. Performed the experiments: CW GR MCHD. Analyzed the data: RLS GR MCHD. Contributed reagents/materials/analysis tools: RLS CW GR MCHD. Wrote the paper: RLS GR MCHD.

Citation: Spietz RL, Williams CM, Rocap G, Horner-Devine MC (2015) A Dissolved Oxygen Threshold for Shifts in Bacterial Community Structure in a Seasonally Hypoxic Estuary. PLoS ONE 10(8): e0135731. https://doi.org/10.1371/journal.pone.0135731

References

1. Diaz RJ, Rosenberg R (2008) Spreading dead zones and consequences for marine ecosystems. Science 321: 926–929. pmid:18703733

2. Vaquer-Sunyer R, Duarte CM (2008) Thresholds of hypoxia for marine biodiversity. Proc Natl Acad Sci U S A 105: 15452–15457. pmid:18824689

3. Howarth R, Chan F, Conley DJ, Garnier J, Doney SC, Marino R, et al. (2011) Coupled biogeochemical cycles: Eutrophication and hypoxia in temperate estuaries and coastal marine ecosystems. Front Ecol Environ 9: 18–26.

4. Rabalais NN, Cai W-J, Carstensen J, Conley DJ, Fry B, Hu X, et al. (2014) Eutrophication-driven deoxygenation in the coastal ocean. Oceanography 27: 172–183.

5. Eby LA, Crowder LB (2002) Hypoxia-based habitat compression in the Neuse River Estuary: context-dependent shifts in behavioral avoidance thresholds. Can J Fish Aquat Sci 59: 952–965.

6. Cheng W, Liu C-H, Hsu J-P, Chen J-C (2002) Effect of hypoxia on the immune response of giant freshwater prawn Macrobrachium rosenbergii and its susceptibility to pathogen Enterococcus. Fish Shellfish Immunol 13: 351–365. pmid:12458742

7. Jorgensen BB (1980) Seasonal oxygen depletion in the bottom waters of a Danish fjord and its effect on the benthic community. Oikos 34: 68–76.

8. Baird D, Christian RR, Peterson CH, Johnson GA (2004) Consequences of hypoxia on estuarine ecosystem function: Energy diversion from consumers to microbes. Ecol Appl 14: 805–822.

9. Urakawa H, Martens-Habbena W, Stahl DA (2010) High abundance of ammonia-oxidizing Archaea in coastal waters, determined using a modified DNA extraction method. Appl Environ Microbiol 76: 2129–2135. pmid:20118363

10. Hofmann AF, Peltzer ET, Walz PM, Brewer PG (2011) Hypoxia by degrees: Establishing definitions for a changing ocean. Deep Res Part I Oceanogr Res Pap 58: 1212–1226.

11. Stevens H, Ulloa O (2008) Bacterial diversity in the oxygen minimum zone of the Eastern tropical South Pacific. Environ Microbiol 10: 1244–1259. pmid:18294206

12. Beman JM, Carolan MT (2013) Deoxygenation alters bacterial diversity and community composition in the ocean’s largest oxygen minimum zone. Nat Commun 4: 2705. pmid:24162368

13. Fuchs BM, Woebken D, Zubkov MV, Burkill P, Amann R (2005) Molecular identification of picoplankton populations in contrasting waters of the Arabian Sea. Aquat Microb Ecol 39: 145–157.

14. Zubkov MV, Tarran GA, Burkill PH (2006) Bacterioplankton of low and high DNA content in the suboxic waters of the Arabian Sea and the Gulf of Oman: Abundance and amino acid uptake. Aquat Microb Ecol 43: 23–32.

15. Madrid VM, Taylor GT, Scranton MI, Chistoserdov AY (2001) Phylogenetic diversity of bacterial and archaeal vommunities in the anoxic zone of the cariaco basin. Appl Environ Microbiol 67: 1663–1674. pmid:11282619

16. Lin X, Wakeham SG, Putnam IF, Astor YM, Scranton MI, Chistoserdov A, et al. (2006) Comparison of vertical distributions of prokaryotic assemblages in the anoxic Cariaco Basin and Black Sea by use of fluorescence in situ hybridization. Appl Environ Microbiol 72: 2679–2690. pmid:16597973

17. Labrenz M, Jost G, Jürgens K (2007) Distribution of abundant prokaryotic organisms in the water column of the central Baltic Sea with an oxic–anoxic interface. Aquat Microb Ecol 46: 177–190.

18. Vetriani C, Tran H V, Kerkhof LJ (2003) Fingerprinting microbial assemblages from the oxic/anoxic chemocline of the Black Sea. Appl Environ Microbiol 69: 6481–6488. pmid:14602603

19. Wakeham SG, Amann R, Freeman KH, Hopmans EC, Jorgensen BB, Putnam I, et al. (2007) Microbial ecology of the stratified water column of the Black Sea as revealed by a comprehensive biomarker study. Org Geochem 38: 2070–2097.

20. Fuchsman CA, Kirkpatrick JB, Brazelton WJ, Murray JW, Staley JT (2011) Metabolic strategies of free-living and aggregate-associated bacterial communities inferred from biologic and chemical profiles in the Black Sea suboxic zone. FEMS Microbiol Ecol 78: 586–603. pmid:22066565

21. Zaikova E, Walsh DA, Stilwell CP, Mohn WW, Tortell PD, Hallam SJ (2010) Microbial community dynamics in a seasonally anoxic fjord: Saanich Inlet, British Columbia. Environ Microbiol 12: 172–191. pmid:19788414

22. Wright JJ, Konwar KM, Hallam SJ (2012) Microbial ecology of expanding oxygen minimum zones. Nat Rev Microbiol 10: 381–394. pmid:22580367

23. Keeling RE, Körtzinger A, Gruber N (2010) Ocean deoxygenation in a warming world. Ann Rev Mar Sci 2: 199–229. pmid:21141663

24. Moore SK, Mantua NJ, Newton JA., Kawase M, Warner MJ, Kellogg JP (2008) A descriptive analysis of temporal and spatial patterns of variability in Puget Sound oceanographic properties. Estuar Coast Shelf Sci 80: 545–554.

25. Babson AL, Kawase M, MacCready P (2006) Seasonal and interannual variability in the circulation of Puget Sound, Washington: A box model study. Atmosphere-Ocean 44: 29–45.

26. Steinberg PD, Brett MT, Bechtold JS, Richey JE, Porensky LM, Smith SN (2010) The influence of watershed characteristics on nitrogen export to and marine fate in Hood Canal, Washington, USA. Biogeochemistry 106: 415–433.

27. Newton J, Devol A, Ruef W (2009) Measuring nitrate fluxes to assess estuarine eutrophication. Decis Anal 2: 5–6.

28. Froehlich HE, Essington TE, Beaudreau AH, Levin PS (2013) Movement patterns and distributional shifts of Dungeness Crab (Metacarcinus magister) and English Sole (Parophrys vetulus) durings easonal hypoxia. Estuaries and Coasts 37: 449–460.

29. Palsson WA, Pacunski RE, Parra TR, Beam J (2008) The effects of hypoxia on marine fish populations in Southern Hood Canal, Washington. Am Fish Soc Symp 64.

30. Newton J, Bassin C, Devol AH, Richey JE, Kawase M, Warner MJ (2011) Hood Canal Dissolved Oxygen Program Integrated Assessment and Modeling Report 2011.

31. Keister JE, Tuttle LB (2013) Effects of bottom-layer hypoxia on spatial distributions and community structure of mesozooplankton in a sub-estuary of Puget Sound, Washington, U.S.A. Limnol Oceanogr 58: 667–680.

32. Devol AH, Ruef W, Emerson S, Newton J (2007) In situ and remote monitoring of water quality in Puget Sound: the ORCA time-series. Environ Res Eng Manag 1: 19–23.

33. Carpenter JH (1965) The Chesapeake Bay Institute Technique for the Winkler Dissolved Oxygen Method. Limnol Oceanogr 10: 141–143.

34. Protocols for the Joint Global Ocean Flux Survey (JGOFS) Core Measurements (1994). UNESCO 29.

35. Sogin ML, Morrison HG, Huber JA, Welch DM, Huse SM, Neil PR (2006) Microbial diversity in the deep sea and the underexplored “rare biosphere.” Proc Natl Acad Sci 103: 12115–12120. pmid:16880384

36. Huse SM, Welch DM, Morrison HG, Sogin ML (2010) Ironing out the wrinkles in the rare biosphere through improved OTU clustering. Environ Microbiol 12: 1889–1898. pmid:20236171

37. Amaral-Zettler L, Artigas LF, Baross J, Bharathi L, Boetius A, Chandramohan D, et al. (2010) A global census of marine microbes. In: McIntyre , editor. Life in the World’s Oceans: Diversity, Distribution and Abundance. Oxford: Blackwell Publishing Ltd. pp. 223–245.

38. Huse SM, Huber JA, Morrison HG, Sogin ML, Welch DM (2007) Accuracy and quality of massively parallel DNA pyrosequencing. Genome Biol 8: R143. pmid:17659080

39. Legendre P, Legendre L (1998) Numerical Ecology. 2nd Englis. Amsterdam, The Netherlands: Elsevier Science B.V.

40. Chao A, Chazdon RL, Colwell RK, Shen T-J (2006) Abundance-based similarity indices and their estimation when there are unseen species in samples. Biometrics 62: 361–371. pmid:16918900

41. Martiny JBH, Eisen JA, Penn K, Allison SD, Horner-Devine MC (2011) Drivers of bacterial β-diversity depend on spatial scale. Proc Natl Acad Sci 108: 7850–7854. pmid:21518859

42. Horner-Devine MC, Lage M, Hughes JB, Bohannan BJM (2004) A taxa-area relationship for bacteria. Nature 432: 750–753. pmid:15592412

43. Clarke KR, Ainsworth M (1993) A method of linking multivariate community structure to environmental variables. Mar Ecol Prog Ser 92: 205–219.

44. Baker ME, King RS (2010) A new method for detecting and interpreting biodiversity and ecological community thresholds. Methods Ecol Evol 1: 25–37.

45. Clarke K (1993) Non-parametric multivariate analyses of changes in community structure. Aust J Ecol 18: 117–143.

46. Pruesse E, Peplies J, Glöckner FO (2012) SINA: accurate high-throughput multiple sequence alignment of ribosomal RNA genes. Bioinformatics 28: 1823–1829. pmid:22556368

47. Cottrell MT, Kirchman DL (2000) Natural assemblages of marine proteobacteria and members of the Cytophaga-Flavobacter cluster consuming low- and high-molecular-weight dissolved organic matter. Appl Environ Microbiol 66: 1692–1697. pmid:10742262

48. Kirchman D (2002) The ecology of Cytophaga–Flavobacteria in aquatic environments. FEMS Microbiol Ecol 39: 91–100. pmid:19709188

49. Pearson PN, Palmer MR (2000) Atmospheric carbon dioxide concentrations over the past 60 million years: 695–699.

50. Mittelbach GG, Steiner CF, Scheiner SM, Gross KL, Reynolds HL, Waide RB, et al. (2001) What is the observed relationship between species richness and productivity? Ecology 82: 2381–2396.

51. Waide RB, Willig MR, Steiner CF, Mittelbach G, Gough L, Dodson SI, et al. (1999) The relationship between productivity and species richness. Annu Rev Ecol Syst 30: 257–300.

52. Horner-Devine MC, Leibold MA, Smith VH, Bohannan BJM (2003) Bacterial diversity patterns along a gradient of primary productivity. Ecol Lett 6: 613–622.

53. Stolper DA, Revsbech NP, Canfield DE (2010) Aerobic growth at nanomolar oxygen concentrations. Proc Natl Acad Sci U S A 107: 18755–18760. pmid:20974919

54. Brewer PG, Hofman AF, Peltzer ET, Ussler W (2014) Evaluating microbial chemical choices: the ocean chemistry basis for the competition between use of O2 or NO3 as an electron acceptor. Deep Sea Res Part I Oceanogr Res Pap 87.

55. Codispoti LA, Brandes JA, Christensen JP, Devol AH, Naqvi SWA, Paerl HW, et al. (2001) The oceanic fixed nitrogen and nitrous oxide budgets: Moving targets as we enter the anthropocene?. Sci Mar 65: 85–105.

56. Bouvy M, Bettarel Y, Bouvier C, Domaizon I, Jacquet S, Le Floc’h E, et al. (2011) Trophic interactions between viruses, bacteria and nanoflagellates under various nutrient conditions and simulated climate change. Environ Microbiol 13: 1842–1857. pmid:21605305

57. Judd KE, Crump BC, Kling GW (2006) Variation in dissolved organic matter controls bacterial production and community composition. Ecology 87: 2068–2079. pmid:16937646

58. Allison SD, Martiny JBH (2008) Resistance, resilience, and redundancy in microbial communities. PNAS 105.

59. Carlson CA, Morris R, Parsons R, Treusch AH, Giovannoni SJ, Vergin K (2009) Seasonal dynamics of SAR11 populations in the euphotic and mesopelagic zones of the northwestern Sargasso Sea. ISME J 3: 283–295. pmid:19052630

60. Vergin KL, Beszteri B, Monier A, Thrash JC, Temperton B, Treusch AH, et al. (2013) High-resolution SAR11 ecotype dynamics at the Bermuda Atlantic time-series study site by phylogenetic placement of pyrosequences. ISME J 7: 1322–1332. pmid:23466704

61. Lennon JT, Jones SE (2011) Microbial seed banks: the ecological and evolutionary implications of dormancy. Nat Rev Microbiol 9: 119–130. pmid:21233850