Microbial functional genes elucidate environmental drivers of biofilm metabolism in glacier-fed streams.

Benthic biofilms in glacier-fed streams harbor diverse microorganisms driving biogeochemical cycles and, consequently, influencing ecosystem-level processes. Benthic biofilms are vulnerable to glacial retreat induced by climate change. To investigate microbial functions of benthic biofilms in glacier-fed streams, we predicted metagenomes from 16s rRNA gene sequence data using PICRUSt and identified functional genes associated with nitrogen and sulfur metabolisms based on KEGG database and explored the relationships between metabolic pathways and abiotic factors in glacier-fed streams in the Tianshan Mountains in Central Asia. Results showed that the distribution of functional genes was mainly associated with glacier area proportion, glacier source proportion, total nitrogen, dissolved organic carbon, and pH. For nitrogen metabolism, the relative abundance of functional genes associated with dissimilatory pathways was higher than those for assimilatory pathways. The relative abundance of functional genes associated with assimilatory sulfate reduction was higher than those involved with the sulfur oxidation system and dissimilatory sulfate reduction. Hydrological factors had more significant correlations with nitrogen metabolism than physicochemical factors and anammox was the most sensitive nitrogen cycling pathway responding to variation of the abiotic environment in these glacial-fed streams. In contrast, sulfur metabolism pathways were not sensitive to variations of abiotic factors in these systems.

Introduction

Glaciers provide vital water resources[1-3] and harbor complex microbial communities[4-6], particularly in arid regions in Central Asia[7-10]. However, as a dynamic interface between the cryosphere, hydrosphere, atmosphere, pedosphere, and biosphere, glacial ecosystems are extremely susceptible to ongoing climate change[11-14]. Glacier-fed streams integrate upstream catchment processes and constitute a prominent geomorphological and ecological component of the alpine landscape[14,15]. In glacier-fed streams, the proportional contribution of ice-melt and snowmelt to overall stream flow diminishes with increasing distance from the glaciers[11,16,17]. This shifting relationship dramatically influences the stream habitat template, including water temperature, channel stability, conductivity, and solute chemistry[16,18-20]. Accordingly, hydrological alteration associated with glacial retreat has the potential to re-structure the habitat template of glacier-fed streams with consequences for the associated biota, including benthic microbial communities and their ecosystem functions.

The highly variable environments of glacier-fed streams impose strong controls on benthic biota[16,21]. Spatial and temporal variation in glacier runoff affect aquatic community structure and biodiversity[11,22], while nutrient availability influences the growth of microorganisms and alters their community structure[23,24]. More specifically, the microbial diversity of benthic biofilms increases with distance from the source glacier, and declines with glacier area proportion, glacier source proportion, and elevation[14,24]. Importantly, this diversity decreases with glacial retreat, illustrating the sensitivity of benthic microbial communities to the effects of climate change[14].

Microorganisms drive the bulk of biogeochemical cycles on the planet and thus influence ecosystem-level functioning[25]. Functional traits are valuable ecological markers to understand microbial community assembly[26] and help elucidate how natural communities and their functions respond to environmental changes[27,28]. In streams, for example, biofilms are hot spots of microbial activity[29], contributing substantially to metabolism and biogeochemical cycles through nutrient uptake, transfer of nutrients to higher trophic levels, and remineralization[30-33]. Unraveling the functional potential of microbial assemblages is critical for gaining insights into their roles in ecosystem processes and their response to environment change[34,35]. For example, microbial metabolism involves a large set of functional genes and biochemical pathways that power biogeochemical cycling in ecosystems. Microbial nitrogen metabolic processes such as nitrification, denitrification, and nitrogen uptake drive nitrogen cycling in streams, controlling nitrogen exports to downstream[36-38]. In anaerobic environments, sulfur compounds are major electron carriers for anaerobic microbial metabolism and sulfate reduction is an important process in organic matter degradation[39-41]. However, our understanding of microbial functions of benthic biofilm in glacier-fed streams remains extremely limited.

In Central Asia, glacial ecosystems are experiencing unprecedented global warming, hydrological variation, and glacier retreat rapidly[42-45]. Glacier-fed stream ecosystems in this semi-arid and cold region are particularly vulnerable. Understanding ecosystem functions and predicting the impacts of climate change on glacier-fed stream ecosystems calls for much better knowledge about microbial processes. High throughput molecular technologies and advanced bioinformatics tools make it possible to predict the functional capacity of microbial communities. Following up on our previous study of microbial communities in glacier-fed streams[24], this research explores the predicted functions of microbial assemblages in benthic biofilms and their relationships to hydrological and physicochemical factors in glacier-fed streams. We predicted metagenomes from 16s rRNA gene sequence data using PICRUSt (Phylogenetic Investigation of Communities by Reconstruction of Unobserved States)[46], identified functional genes associated with nitrogen and sulfur metabolisms based on KEGG (Kyoto Encyclopedia of Genes and Genomes) database[47], and explored the relationships between metabolic pathways and abiotic factors. We hypothesized that the functional genes of biofilm community in glacier-fed streams are closely associated with community taxonomic compositions and abiotic environments, especially hydrological factors. Although metagenomics reflects potential rather than realized functional capacity, our data offer a window into the poorly known metabolism of benthic biofilms in glacier-fed streams, providing insights into the potential impacts of glacial melting on the ecosystem function.

Methods

Study Area

The Tianshan Mountains are located in Central Asia. Based on the dataset from the Second Glacier Inventory of China[48], there are 7934 glaciers distributed in the China’s Tianshan Mountain with a total area of 7179 km2 and an ice volume of 707 km3, which account for 16.3%, 13.9%, and 15.7% of the total number, total area, and total volume of glaciers in China, respectively[49]. With global warming, these glaciers have been shrinking rapidly[42-45]. For example, the Urumqi Glacier No. 1 (GN1, 43°06′N, 86°49′E) in eastern Tianshan Mountain had an area of 1.95 km2 and a terminus altitude of 3730 m in 1962, but by 2010 had shrunk to 1.65 km2 and the terminus retreated to the altitude of 3743 m[50].

Sampling

In June 2016, we collected water and benthic biofilm samples from 11 sample sites distributed along two glacier-fed streams in the Tianshan Mountains (Fig. 1). The samples from the sites S1-S5 (elevations ranged from 3548 to 2840 m) were collected on June 19, 2016. The samples from the sites N1-N6 (elevations ranged from 3828 to 2646 m) were collected on June 20, 2016. At each sample site, we randomly sampled 6–9 submerged rocks from the river cross section at depths of 10 to 30 cm. The benthic biofilm was removed by rigorously brushing an area of 4.5 cm diameter from the upper surface of each stone using a sterilized nylon brush (changed between samples) and rinsing the slurry with sterile water. From the mixed slurry, approximately 10 ml was filtered through a 0.2-μm membrane filter that was immediately frozen in liquid nitrogen in the field and transported frozen to the lab after all the sampling work. In the lab, benthic biofilm samples were stored at −80 °C until DNA extraction.

Figure 1

Study area and sample sites. 11 samples were collected from two glacier-fed streams in June 2016 in the Tianshan Mountains, Central Asian. The map is from Ren et al. (2017)[24] and was created in ArcGIS 14.0 (http://desktop.arcgis.com/en/arcmap/) using DEM image download from USGS (https://earthexplorer.usgs.gov/).

DNA Extraction, PCR, and Sequencing

Bacterial 16S rRNA genes were analyzed to determine the predicted functional genes and metabolic pathways of the microbial community in biofilms. DNA extraction, PCR, and sequencing were described in our previous study[24].

Physicochemical factors

At each sample site, water temperature (Temp), dissolved oxygen (DO), pH, and conductivity (Cond) were measured in situ using a YSI handheld meter (Model 80, Yellow Springs Instruments, Yellow Springs, Ohio). Elevation was measured using a GPS unit (Triton 500, Magellan, Santa Clara, California). Water samples were collected for nutrients and dissolved organic carbon (DOC) analysis. Total nitrogen (TN) was analyzed by ion chromatography after persulfate oxidation (EPA 300.0)[51]. Nitrate (NO3 −) was determined by ion chromatography (EPA 300.0)[51]. Ammonium (NH4 −) was analyzed using the indophenol colorimetric method (EPA 350.1)[52]. Total phosphorus (TP) was quantified using the ammonium molybdate method after oxidation (EPA 365.3)[53]. Soluble reactive phosphorus (SRP) was analyzed using the ammonium molybdate method (EPA 365.3)[53]. DOC was analyzed using a Shimadzu TOC Analyzer (TOC-VCPH, Shimadzu Scientific Instruments, Columbia, Maryland). The results of these basic physicochemical factors are shown in our previous study[24].

Hydrological Factors

The glacier distribution data were obtained from the Second Glacier Inventory of China[48]. The river channel network and the catchment boundaries were derived from topographic data, the Digital Elevation Model (DEM) at 90 m resolution (https://earthexplorer.usgs.gov/). Hydrological factors were calculated using the landscape-based hydrological model of Gao et al. (2016, 2017)[10,17], including the glaciated area proportion (GA) for each sample site, runoff proportion generated from glacier source (GS), and distance to glaciers (GD). Methods for calculating these basic hydrological factors are shown in our previous study[24].

Analysis

Bacterial 16s rRNA sequence data (available at National Center for Biotechnology Information, SRP115356) were cleaned using the software package QIIME[54] and then clustered to operational taxonomic units (OTUs) with a complete linkage algorithm on a 97% sequence identity level. The metagenomes were predicted from 16s data using Phylogenetic Investigation of Communities by Reconstruction of Unobserved States (PICRUSt)[46]. The functional genes associated with nitrogen metabolism and sulfur metabolism were identified based on Kyoto Encyclopedia of Genes and Genomes (KEGG) database[47]. Redundancy analysis (RDA)[55] was conducted to reveal the association of microbial communities in relation to abiotic factors based on functional gene abundances using vegan package 2.4[56] in R 3.3.2[57]. A Monte Carlo test was used to assess the significance of abiotic factors in RDA[55]. Mantel tests were run to assess how functional dissimilarity (Bray-Curtis distance calculated from functional gene abundances) correlated with taxonomic dissimilarity (Bray-Curtis distance calculated from taxon abundances), environmental distance (Euclidean distance), and geographic distance. Correlation analysis was conducted to assess the relationships between metabolic pathways and abiotic factors with the P-value adjusted by Bonferroni correction using psych package 1.7.5[58] in R 3.3.2[57].

Results and Discussion

Functional Variation

The 16S rRNA data set consisted of 127,146 sequences clustered in 7545 OTUs and 5771 functional genes. The Bray-Curtis distance calculated from functional gene abundances was significantly correlated with operational taxonomic units (OTUs) dissimilarity (Bray-Curtis distance calculated from taxon abundances, Mantel r = 0.614, P < 0.001) and environmental distance (Euclidean distance, Mantel r = 0.567, P < 0.001), but not with geographic distance (Mantel r = 0.287, P = 0.413, Fig. 2), suggesting that the overall changes in functional genes composition were closely associated with community taxonomic composition and abiotic conditions. Redundancy analysis (RDA) was applied to depict the spatial distribution of microbial communities in relation to effects of physicochemical and hydrological factors in terms of the functional gene composition (Fig. 3). The first two axes accounted for 75.6% of the variance (RDA 1: 50.34%; RDA 2: 25.26%). A Monte Carlo test showed that glacier area proportion (GA), glacier source proportion (GS), total nitrogen (TN), dissolved organic carbon (DOC), and pH were the abiotic factors most closely associated (p < 0.05) with the distribution of microbial communities (Fig. 3). Canonical correspondence analysis based on taxon abundances showed similar results[24].

Figure 2

Plots of Mantel tests showing the relationships between (a) functional genes dissimilarity and OTUs dissimilarity, (b) functional genes dissimilarity and environmental distance, and (c) functional genes dissimilarity and geographic distance. One point represents one sample pair. Blue dotted lines denote the 95% confidence interval.

Figure 3

RDA plot revealing the association of microbial communities (based on functional gene abundances) in relation to abiotic factors. The first and second axes accounted for 75.6% of the variance (RDA 1: 50.34%; RDA 2: 25.26%). Only the abiotic factors that were significantly correlated with the community (forward selection with Monte Carlo test, P < 0.05) are shown as the black solid lines. Blue dots represent the sampling sites and the abbreviations are the names of the sampling sites.

Different biomes clearly harbor distinct microbial assemblages[27,59,60]. It has been demonstrated that metagenomic functional composition is strongly correlated with taxonomic composition in spite of functional redundancy at the genomic level[61,62]. Moreover, recent work also suggests that the metabolic functional potential of microbial communities in ocean and soil is closely related with environmental conditions[63]. Our study also indicates that the functional metabolic potential of glacier-fed streams is strongly influenced by environmental conditions and taxonomic structure.

Potential Nitrogen Metabolism

Nitrogen has various chemical species that are cycled by a suite of coupled biogeochemical processes[64] that are catalyzed by various microbe-derived enzymes[65,66]. To deepen our understanding of glacier-fed stream ecosystems, it is crucial to know which N cycling pathways are operating and how they respond to the surrounding abiotic environment. In this study component, we identified and compared the potential nitrogen cycle pathways that might occur in the glacier-fed streams based on the functional genes encoding the enzymes required for specific nitrogen cycling pathways.

Of the 5,571 detected genes, 48 functional genes were associated with nitrogen metabolism (Figure S1) and the relative abundance of these functional genes ranged from 0.79% to 0.89% of the whole gene reads. Figure S1 summarizes the relative abundance of the genes encoding enzymes that catalyze nitrogen metabolism. The N cycle involves four reduction pathways (assimilatory nitrate reduction, dissimilatory nitrate reduction to ammonia, denitrification, and nitrogen fixation) and two oxidation pathways (nitrification and anammox) (Fig. 4). The relative abundance of genes associated with assimilatory nitrate reduction and dissimilatory nitrate reduction ranged from 0.044% to 0.072% and from 0.066% to 0.151%, respectively (Fig. 5). The assimilatory pathways require energy and start with the reduction of NO3 − to NO2 − and then to NH4 +, which is highly bioavailable and can be readily used by cells for synthesis of amino acids and nucleotides[67,68]. Dissimilatory nitrate reduction to ammonia uses N-compounds to provide energy (ATP) and is another catabolic pathway that can retain the nitrogen in the system in a bioavailable form (NH4 +) for further biological processes[69,70]. Generally, assimilatory pathways were much more prevalent than dissimilatory pathways[71]. In our studied glacier-fed streams, however, genes associated with dissimilatory pathways were much more abundant that those involved in assimilatory pathways, suggesting that the microbial assemblages in glacier-fed streams use inorganic nitrogen more as an energy source more as an N source for biosynthesis. Nitrogen fixation is an energetically demanding process and only a few taxonomic groups genera can carry out this process[72,73]. As a competing process to nitrate reduction pathways, denitrification is the main biological process for removal of N from freshwater systems[70,74]. The relative abundance of genes associated with denitrification and nitrogen fixation ranged from 0.038% to 0.123% and from 0.012% to 0.074%, respectively (Fig. 5). Nitrification is an essential process in nitrogen cycle performed by nitrifiers converting ammonium to nitrate. In our study, functional genes encoding enzymes that catalyze nitrification had a relative abundance ranges from 0.01% to 0.05% (Fig. 5). Anammox plays a significant role in the nitrogen cycle, converting nitrite and ammonia to dinitrogen[75], and is another microbial process that releases fixed nitrogen from the environment as dinitrogen[76-79]. In our study, functional genes encoding enzymes that catalyze anammox had only a very small relative abundance (Fig. 5).

Figure 4

A simplified schematic of the nitrogen cycle constructed with the use of the KEGG nitrogen metabolism pathway. The N cycle consists of interconnected biogeochemical pathways (dissimilatory nitrate reduction, assimilatory nitrate reduction, denitrification, nitrogen fixation, nitrification, and anammox). There are a variety of genes encoding enzymes that catalyze the important transformation reactions of various nitrogen forms ranging from the oxidation states +5 in nitrate to −3 in ammonium. The percentage of associated gene reads were predicted using PICRUSt based on KEGG database.

Figure 5

Relative abundance of major categories of functional genes encoding the enzymes that catalyze nitrogen cycling pathways (dissimilatory nitrate reduction, assimilatory nitrate reduction, denitrification, nitrification, nitrogen fixation, and anammox) based on the KEGG database.

In ecosystems, nitrogen metabolic pathways are particularly susceptible to environmental redox fluctuations because of the large difference between the oxidation state of nitrate (+5) in nitrate and that of ammonia (−3)[80,81]. Moreover, some other environmental factors also strongly impact nitrogen cycling. For example, microorganisms can be stimulated by high concentration of NO3 − to increase populations containing the functional genes encoding the enzymes required for nitrate reduction. Low N:P ratios and low NH4 + can stimulate N2 fixing bacteria[82-84]. It has also been demonstrated that anammox is influenced by physicochemical properties, such as nitrate concentration, C/N ratio, ammonia concentration, and pH[85]. Nitrification can be affected by temperature, salinity, light, organic matter concentration, substrate concentrations, pH, and oxygen concentration[86]. In our study, we explored the relationships between nitrogen cycle pathways and abiotic environment including physicochemical and hydrological factors (Table S1). The results showed that anammox was negatively correlated with GA and GS. Denitrification, dissimilatory nitrate reduction, and nitrification were negatively correlated with GD. In addition, anammox was positively correlated with pH, temperature, and DOC but negatively correlated with TN and NO3. Assimilatory nitrogen reduction was positively correlated with TN and NO3. Nitrification was negatively correlated with pH. Nitrogen fixation was positively correlated with NH4. These results suggest that hydrological factors had more significant influences on nitrogen metabolism than physicochemical factors. Among the nitrogen cycling pathways, anammox was the most sensitive to variation of abiotic environment in the glacial-fed streams we studied (Table S1). Since global warming leads to significant changes of abiotic environment, especially in the hydrology in glacier-fed streams in Central Asia, these results may allow us to predict changes in nitrogen metabolism under climate change.

Potential Sulfur Metabolism

Overall, there were 43 functional genes associated with sulfur metabolism in our data set (Figure S2). The relative abundance of these functional genes ranged from 0.848% to 1.009% of the total gene reads. Figure S2 summaries the relative abundance of the genes encoding sulfur metabolism based on the KEGG database. Assimilatory sulfate reduction occurs in many anoxygenic phototrophic bacteria and leads to the biosynthesis of sulfur-containing compounds[87,88]. Dissimilatory sulfate reduction is a pathway in which sulfate or sulfur serves as the terminal electron acceptor of the respiratory chain, producing inorganic sulfide. In our study, the relative abundance of genes associated with assimilatory sulfate reduction ranged from 0.247% to 0.313%, followed by genes involved in the sulfur oxidation (SOX) system and in dissimilatory sulfate reduction (Fig. 6), suggesting that microbial assemblages in glacier-fed streams regard inorganic sulfur as an S source in biosynthesis more than an electron acceptor. The SOX system, a sulfur oxidation pathway, is a multienzyme pathway found in both photosynthetic and lithoautotrophic sulfur-oxidizing bacteria[39,89,90]. Sulfur oxidation is a critical step in the biogeochemical sulfur cycle[91,92]. Our detection of genes encoding the enzymes catalyzing sulfur oxidation suggests that microbial activity is generating acidity that may also affect mineral formation and stability[88]. Previous studies have suggested that the availability of electron donors[93], light[94], and aerobic conditions[95] are important environmental factors affecting sulfur metabolism. In our study, only assimilatory sulfate reduction was positively correlated with GD (Table S1), suggesting that sulfur metabolism pathways were not sensitive to variations of hydrological and physicochemical factors in glacier-fed streams.

Figure 6

Relative abundance of major categories of functional genes encoding the enzymes that catalyze sulfur metabolism, including assimilatory sulfate reduction, dissimilatory sulfate reduction, and the SOX system.

Conclusions

In glacier-fed streams, microbial assemblages in benthic biofilms change significantly with variation in the abiotic environment and thus appear vulnerable to glacier retreat caused by climate change. As most glaciers on Earth are receding rapidly, the habitat template (especially hydrological factors) of glacier-fed streams will be increasingly disrupted. Community taxonomic composition and the abiotic environment regulate microbe-mediated functions. In particularly, hydrological factors have significant impacts on the distribution of functional genes and nitrogen cycling pathways. Thus, large changes in biogeochemical processes and associated ecosystem functions likely will occur as glacier-fed ecosystems are transformed due to loss of their association with the cryosphere.