Microb Ecol DOI 10.1007/s00248-009-9525-9
MICROBIOLOGY OF AQUATIC SYSTEMS
Habitat Heterogeneity and Associated Microbial Community Structure in a Small-Scale Floodplain Hyporheic Flow Path Jennifer L. Lowell & Nathan Gordon & Dale Engstrom & Jack A. Stanford & William E. Holben & James E. Gannon
Received: 16 December 2008 / Accepted: 22 April 2009 # Springer Science + Business Media, LLC 2009
Abstract The Nyack floodplain is located on the Middle Fork of the Flathead River, an unregulated, pristine, fifthorder stream in Montana, USA, bordering Glacier National Park. The hyporheic zone is a nutritionally heterogeneous floodplain component harboring a diverse array of microbial assemblages essential in fluvial biogeochemical cycling, riverine ecosystem productivity, and trophic interactions. Despite these functions, microbial community structure in pristine hyporheic systems is not well characterized. The current study was designed to assess whether physical habitat heterogeneity within the hyporheic zone of the Nyack floodplain was sufficient to drive bacterial β diversity between three different hyporheic flow path locations. Habitat heterogeneity was assessed by measuring soluble reactive phosphorous, nitrate, dissolved organic carbon, dissolved oxygen, and soluble total nitrogen levels seasonally at surface water infiltration, advection, and exfiltration zones. Significant spatial differences were detected in J. L. Lowell (*) : N. Gordon : W. E. Holben : J. E. Gannon Microbial Ecology Program, Division of Biological Sciences, The University of Montana, 32 Campus Dr. 4824, Missoula, MT 59812-1002, USA e-mail: [email protected]
J. L. Lowell : W. E. Holben Montana—Ecology of Infectious Diseases Program, The University of Montana, Missoula, MT 59812-1002, USA D. Engstrom Department of Geosciences, The University of Montana, Missoula, MT 59812-1002, USA J. A. Stanford Flathead Lake Biological Station, The University of Montana, Polson, MT 59860, USA
dissolved oxygen and nitrate levels, and seasonal differences were detected in dissolved oxygen, nitrate, and dissolved organic carbon levels. Denaturing gradient gel electrophoresis (DGGE) and cell counts indicated that bacterial diversity increased with abundance, and DGGE fingerprints covaried with nitrate levels where water infiltrated the hyporheic zone. The ribosomal gene phylogeny revealed that hyporheic habitat heterogeneity was sufficient to drive β diversity between bacterial assemblages. Phylogenetic (P) tests detected sequence disparity between the flow path locations. Small distinct lineages of Firmicutes, Actinomycetes, Planctomycetes, and Acidobacteria defined the infiltration zone and α- and β-proteobacterial lineages delineated the exfiltration and advection zone communities. These data suggest that spatial habitat heterogeneity drives hyporheic microbial community development and that attempts to understand functional differences between bacteria inhabiting nutritionally heterogeneous hyporheic environments might begin by focusing on the biology of these taxa.
Introduction Fluvial systems are dynamic environments shaped by interactions between channel-associated geomorphic processes, the alluvial aquifer, and riparian vegetation, all of which drive lateral and vertical connectivity within a floodplain [39, 48, 49]. The hyporheic zone is a unique heterotrophic component of fluvial ecosystems, comprising the interface between open surface water and groundwater [17, 45, 50, 51]. It is within this spatially fluctuating ecotone that the exchange of nutrients and biota between surface water, groundwater, the alluvial aquifer, and the riparian zone occurs . The resulting mosaic of habitats and nutrient gradients are capable of supporting a diverse
J. L. Lowell et al.
array of microorganisms important for biogeochemical cycling between floodplain components [14–16, 23]. As a result, productivity of the hyporheic zone can be orders of magnitude greater than that in benthic sediments  accounting for 76–96% of floodplain ecosystem respiration , and playing a crucial role in nutrient supply and biotic productivity in surface waters [21, 23]. Despite these contributions from hyporheic biofilms, most studies addressing microbial community structure in streams have targeted benthic biofilms [2, 5, 9, 19, 27, 34, 35]. The few existing studies of hyporheic microbial community structure have focused on contaminated systems , artificial systems [3, 4, 6] and geographically distinct hyporheic sites , leaving the drivers of hyporheic biofilm community structure within natural floodplain ecosystems virtually unknown. In the current study, we measured nutrient and dissolved oxygen levels to assess habitat heterogeneity along a 100-m hyporheic flow path in a pristine northwestern Montana floodplain. These values were compared to bacterial abundance, 16S rRNA gene sequence data, and 16S polymerase chain reaction (PCR)-DGGE patterns. We hypothesized that bacterial assemblages in the hyporheic zone would correlate to spatial habitat heterogeneity creating β diversity (i.e., differences in communities between heterogeneous environments) among bacterial communities. A redundancy analysis and tests for correlation were used to relate differences in bacterial community structure to environmental variability. Additionally, we incorporated phylogenetic information into our community analyses, accounting for the degree of sequence divergence and relative sequence abundance, to make a quantitative assessment of bacterial β diversity [32, 36]. Finally, we discuss the characteristics of the bacterial taxa identified using 16S rRNA gene sequence data and the implications of these taxa for hyporheic zone productivity.
Materials and Methods Study Site The Nyack Floodplain (hereafter, Nyack) located on the Middle Fork of the Flathead River is a fifth-order unregulated stream in Northwestern Montana (48°26′ 53.77″ N, 113°49′28.28″ W). The valley encompassing the Nyack borders Glacier National Park and the Bob Marshall/Great Bear wilderness area and is considered a nearly pristine river catchment. The Nyack is approximately 8 km long and 1 km wide, is situated under a riparian forest, and contains 80% of regional vascular plant flora and at least 70% of regional aquatic invertebrate species, including 80 hyporheic species [22, 46]. Hyporheic flow
paths exist at two distinct spatial scales: at the large floodplain scale (kilometers) and at the scale of individual gravel bars (meters) . Our study site was located on a smaller scale flow path (~100 m) designated as the Movie Road gravel bar (hereafter, MR) within the Nyack (48°26′ 53.77″ N, 113°48′33.41″ W). Infiltration, river-dominated lateral flow of groundwater (advection), and exfiltration zones were previously identified by measuring hydraulic head gradients through piezometer placements in the near bank and river bottom around the gravel bar perimeter. Water and Sediment Sampling Groups of five 1-in. polyvinyl chloride (PVC) sampling wells were installed using a Geoprobe model 5400 at each of the three hydraulic zones (infiltration, advection, and exfiltration) across the gravel bar (n = 3, N = 15). Each group of five wells was placed perpendicular to the flow path in order to obtain five replicate samples for each time point and location (Fig. 1). Bulk sediment was collected from the banks of the MR gravel bar and hand-sieved using methods previously described . The 1.7–2.36-mm size sediment fraction was transported to the laboratory, washed with Milli-Q water, autoclaved at 120°C for 120 min, and baked overnight in a drying oven at 265°C. This procedure was repeated three times on separate days. Sterile sediment (10 g) was packed into acid-washed PVC tubes 9 cm long × 2.1 cm in diameter and slotted with 0.15 cm wide slots . The ends of each tube were plugged with autoclaved, size 16D silicone stoppers (Saint-Gobain, Akron, OH) to retain the sediments. Four tubes of sediment were connected end to end using acid-washed zip ties and installed in each of the 15 1-in. wells. The tubes were suspended within the top meter of the saturated zone using 60-pound-test fishing line attached to the well cap. Following sampler placement in the wells, biofilms were allowed to colonize the sediment for 20 weeks. This time interval was selected based on prior experiments (performed at this site using the same techniques) that compared microbial communities colonizing implanted samplers to those on fresh sediment collected from adjacent pits of equal depth and sampling time interval. DGGE pattern matching analysis showed that, after 16–20 weeks of incubation in situ, sediment community composition in the samplers was highly similar to that on freshly taken native sediment samples . Following the colonization period, one sediment sampler from each well was collected every 16 weeks over a 64-week period, for a total of 60 sediment samples. Sediment samples were placed in sterile “Twirl’Em” bags (Fisher Scientific, Pittsburg, PA) and transported to the laboratory on wet ice. Samples to be used in molecular analyses were immediately frozen at −70°C before DNA extraction.
Hyporheic Flow Path Microbial Communities
nitrogen (STN) were measured during the spring and early fall sampling times only. Water chemistries were measured using the following calorimetric methods according to “The standard methods for the examination of water and wastewater, American Public Health Association guidelines” (http://www.apha.org): SRP, ascorbic acid colorimetry; NO3−, Cd reduction–azole colorimetry; and STN, filtration–persulfate digestion–CD reduction–azole colorimetry. DOC was measured using a Tekmar–Dohrmann carbon analyzer (Teledyne Tekmar, Mason, OH), and DO was measured using an YSI Model 85 handheld DO meter (Yellow Springs Instruments Inc., Yellow Springs, OH). Direct Microscopic Enumeration of Sediment Biofilm Bacterial Communities
Figure 1 The top panel shows the location of the Nyack floodplain within MT. The bottom panel shows the MR gravel bar expanded from the Nyack floodplain map. M1 A-E infiltration wells, M2 D-H advection wells, M3 D-H exfiltration wells
Additionally, wells were purged for several minutes, and two water samples per well were collected. The first were collected in sterile 250-ml plastic bottles for inorganic nutrient analyses and the second in sterile 250-ml glass bottles for dissolved organic carbon (DOC) analyses. Samples were either acidified with phosphoric acid or filtered through 0.2-μm nylon membrane filters (Whatman, Florham Park, NJ) and frozen for nutrient measurements. Sampling dates were July 8, 2005 (summer), November 10, 2005 (late fall), May 8, 2006 (spring), and September 15, 2006 (early fall). Winter samples were not taken, as this area receives very heavy snowfall. Environmental Variables To assess habitat heterogeneity, soluble reactive phosphorous (SRP), nitrate (NO3−), and DOC were measured for all sampling times. Dissolved oxygen (DO) and soluble total
Bacteria in each sediment sample were enumerated (in triplicate) by direct epifluorescent microscopy after staining with 4′,6-diamidino-2-phenylindole (DAPI) . Treatments were carried out in sterile 15-ml centrifuge tubes. Two grams of sediment (wet weight) were suspended in filter-sterilized formation water containing 0.01 M tetrasodium pyrophosphate  and fixed with formalin (3% final concentration). Samples were mixed thoroughly for 1 min, held on ice for 30 min, mixed again for 1 min, and sonicated for 20 min using a Branson 3210R-DTH ultrasonic cleaner (Fredricksburg, VA). After sonication, samples were vortex-mixed for 1 min and centrifuged at 750×g for 10 min using a Sorvall RC 5B Plus (Thermo Scientific, Waltham, MA) to remove large particulates . Subsamples from each supernatant were transferred to amber microcentrifuge tubes (Fisher Scientific, Waltham, MA) and stained for 30 min with DAPI at a final concentration of 10.0 mg l−1 . Samples (1.0 ml) from each tube were filtered through 0.2-μm black polycarbonate membranes (Poretics Products, Kent, WA), and microscope slides were prepared as previously described . Filters were examined microscopically at magnification ×1,000 using a Zeiss Axioskop 20, and ten fields per filter were counted. Direct Plate Counts from Sediment Biofilm Bacterial Communities Two grams of sediment (wet weight) from each well were suspended in 8 ml of filter-sterilized formation water containing 0.01 M tetrasodium pyrophosphate . Samples were mixed thoroughly for 1 min, held on ice for 30 min, mixed again for 1 min, then sonicated for 20 min as described above. After sonication, samples were vortexmixed for 1 min and allowed to settle for 15 min. Each supernatant was serially diluted using filter-sterilized formation water. Four dilutions (1:1,000–1:1,000,000) were
J. L. Lowell et al.
plated in triplicate on agar made with a filter-sterilized mixture (1:9) of soil leachate and river water and 1% Noble agar (Difco). Briefly, soil leachate was created by continuous agitation of floodplain soil/sediment (excavated from the top meter of the water table) in Milli-Q water for 1 h on a rotary shaker (150 rpm) followed by centrifugation (5,000 rpm for 10 min) and filter sterilization (0.2 μm). The leachate was analyzed for dissolved organic carbon using the above method and diluted to 3 mg/L DOC (about three times over in situ DOC levels). Plates were incubated for 4 weeks at 10°C in the dark (in situ temperatures range from 4 to 12°C). Colonies were counted on plates containing 30–300 total colony-forming units (CFUs). DNA Extraction from Sediment Biofilm Bacterial Communities Bacterial community DNA was recovered from 1.0 g of sediment collected from each well in each season as previously described [16, 53] with modification as follows: Bead beating was performed using the Genogrinder 2000 (Clifton, NJ), followed by freezing in liquid nitrogen, and then boiling for 2 min. The bead-beating and freeze–thaw cycles were repeated four times before proceeding to the extraction and precipitation steps . Before molecular analyses, DNA extractions were further purified using the Promega Wizard® DNA cleanup system (Madison, WI). Samples were quantified using a NanoDrop® ND-1000 UV-Vis spectrophotometer (Wilmington, DE) and stored at −20°C for downstream analyses. PCR Amplification and DGGE Analysis PCR-DGGE was performed on individual DNA samples and also on pooled replicate DNA samples from each location and season for a total of 72 samples. Each 50-µl reaction contained 1× PCR buffer with 1.5 mM MgCl2, 200 μM dNTPs, 0.8 μM of each general bacterial 16S rRNA gene primer 536fc, which included a GC clamp for DGGE  and 907r , 5 μg of bovine serum albumin (Roche, Pleasanton, CA), 1.25 U HotStar Taq polymerase (Qiagen, Valencia, CA), and 100 ng of template DNA. Reactions were performed in a PTC-100 thermal cycler (MJ Research, Waltham, MA) at 95°C for 5 min followed by 30 cycles of 95°C for 1 min, 56°C for 1 min, and 72°C for 1 min and a final elongation step of 72°C for 10 min. PCR product was visualized on 1.5% agarose gels and quantified using a FirstLight UV Illuminator (UVP Inc., Upland, CA) gel image analyzer. PCR amplicons from individual sediment DNA samples were separated by DGGE on a Bio-Rad D-GENE System (Bio-Rad Laboratories, Hercules, CA) to assess microbial community variability across each group of wells (n=5) for
each time-point (season) and for each zone (location). Amplicons from the pooled replicate DNA samples were used for community comparisons on one DGGE gel to avoid introducing gel-to-gel variability issues. A linear gradient of 45–60% denaturant (100% denaturant equals 7 M urea and 40% w/v formamide) in an 8% poylacrylamide gel matrix was used. Gels were run at 70 V for 17 h at 60°C. Following electrophoresis, gels were stained for 1 h at 37°C with SYBR® Gold (Molecular Probes, Eugene, OR) and visualized using a FirstLight UV Illuminator (UVP Inc., Upland, CA) as described above. Similarity matrices were constructed from DGGE banding information using the Jaccard similarity coefficient function in Bionumerics v. 4.61 (Applied Maths, Austin, TX). Cloning and Sequencing of DGGE Bands Major DGGE bands were excised from each lane in the pooled sample gel, re-amplified, cloned, and sequenced. Excised bands were placed in 100 µl of sterile water and incubated at 4°C overnight to elute DNA from the gel matrix. Re-amplification was performed according to the above PCR conditions using 1 µl of the eluted DGGE band DNA as template and an unclamped version of primer 536f and 907r. PCR products were purified using the Qiagen PCR purification kit (Valencia, CA), and TA cloned using the pGEM®-T Easy Vector System and kit supplied JM109 cells (Promega, Madison, WI). Five colonies per each cloned amplicon were chosen in order to detect heterogeneity in DGGE bands caused by the co-migration of multiple PCR amplicons . The selected colonies were grown for 20 h at 37°C in Luria–Bertani broth containing 100 µg/ml of ampicillin after which plasmids were purified with the 5-Prime Inc., Perfectprep Plasmid 96 Vac Kit (Gaithersburg, MD) and submitted to Polymorphic DNA Technologies, Inc. (Alameda, CA) for sequencing. Sequencing of DGGE bands yielded products of approximately 371 bp (following removal of primer sequences) spanning the V4–V5 region of the 16S gene. Recent studies have indicated that relatively short 16S fragments of the V4–V5 region are sufficient for detecting taxon relationships with high jackknife support and that, by using this region, it is possible to obtain the same classification outcome (to the genus level) as from full-length sequences [30, 37]. Statistical Analyses Two-way analyses of variance (ANOVAs) were used to detect variation in SRP, NO3−, DOC, DO, and STN levels and bacterial counts. Similarity matrices were constructed from DGGE fingerprints, and non-parametric ANOVAs (Kruskal–Wallis) were used to assess bacterial community
Hyporheic Flow Path Microbial Communities
structure variability within replicate samples at each gravel bar location. Samples from each time and location were then pooled based on failure to detect significant differences in banding patterns within the replicates. Pooled sample analysis also enabled us to run all samples on one gel, thereby eliminating gel-to-gel variability . Correlations between environmental variables, DAPIs, CFUs, and DGGE similarity matrices were examined by Spearman’s rho tests for correlation. The above statistical tests were carried out using SPSS® v. 15.0 for windows (SPSS Inc., Chicago, IL). Bray–Curtis distance measures were generated from DGGE band presence/absence data, as recommended by . We used distance-based redundancy analysis (DBRDA)  to test the hypotheses that sample location contributed to a significant proportion of the total variance in the bacterial community profiles and that bacterial community structure as determined by DGGE varied with nutrient levels. For these analyses, distance calculations and principal coordinate analysis (PCA) were performed using the PrCoord utility packaged with Canoco (ver. 4.5).
combined communities. The weighted UniFrac option was chosen as a quantitative measurement of β diversity, accounting for the relative number of times each sequence was observed in each habitat . The UniFrac lineagespecific analysis tool was then used to determine which phylogenetic lineages contributed significantly to β diversity detected by weighted UniFrac . Additionally, the P test was used to examine whether each location harbored distinct phylogenetic lineages (i.e., covaried with the tree phylogeny) and which taxa spatially delineated bacterial communities . The β-proteobacterial lineage was identified as a significant contributor to community differences. Therefore, β diversity measures and P tests were repeated with all sequences in the tree, the proteobacterial sequences alone, and with each major β-proteobacterial group removed from the tree as follows: the Comamonadaceae and each of two groups of Oxalobactereaceae. The tree was also re-evaluated with the α-proteobacterial group removed. Elimination of each group individually allowed us to further assess how potential lineage differences between the three flow path locations influenced β diversity.
Sequence Analysis Results Sequence electropherograms were examined and trimmed using the SeqMan option in Lasergene v. 5.01. Only unique sequences were kept from each of the five clones sequenced per DGGE band. Sequence data were submitted to the Ribosomal Database Project (RDP) Chimera Check (http:// 220.127.116.11/cgis/chimera.cgi?su=SSU) to identify artifactual sequences, and those that appeared chimeric were excluded from downstream analyses. The RDP Sequence Match tool was used to search for the nearest neighbors of the remaining sequences. The 191 retained sequences were aligned using the nearest alignment space termination aligner , and a maximum parsimony tree was constructed using ARB parsimony . The tree was rooted with Thermotoga maritima (GenBank accession number AE000512). The resulting phylogeny was exported for use in UniFrac to measure β diversity among the infiltration, advection, and exfiltration zone bacterial communities [31, 32]. This tool allows comparisons of β diversity using genetic diversity, phylogenetic distances and the degree of divergence between taxa to determine microbial community structure according to environment or habitat . For this study, the infiltration, advection, and exfiltration zones were considered disparate habitats based on significant chemical profile differences. Each sequence in the phylogenetic tree was assigned to the habitat from which the sample came, and UniFrac was used to measure the phylogenetic distances between communities by comparing the genetic diversity of each community to the total diversity of the
Environmental Variables Measurements were taken for SRP, NO3−, and DOC for each location and season across the gravel bar to assess habitat heterogeneity. SRP levels did not differ significantly. Spatially, mean NO3− levels were higher at infiltration and decreased significantly across the flow path (p