Influences of anthropogenic land use on microbial community structure and functional potentials of stream benthic biofilms

Stream ecosystems are the primary receivers of nutrient and organic carbon exported from terrestrial ecosystems and are profoundly influenced by the land use of the surrounding landscape. The aquatic impacts of anthropogenic land use are often first observed in stream benthic biofilms. We studied the benthic biofilms in streams flowing through forest (upstream) and anthropogenic land use (downstream) areas in southwestern China. The results showed that anthropogenic land use increased nutrient and organic carbon in both stream water and benthic biofilms, which are closely related to the differences in the microbial communities. The taxonomic dissimilarity of the communities was significantly correlated with the functional gene dissimilarity, and the upstream sites had more distinct functional genes. Network analysis showed that upstream sites had more highly connected microbial networks. Furthermore, downstream sites had higher relative abundances of anammox and denitrification suggesting stronger nitrogen removal than upstream sites. Increased nutrients in both the stream water and biofilms caused by anthropogenic land use had severe impacts on the nitrogen cycle in stream ecosystems. Downstream sites also had stronger carbon metabolism than upstream sites. This study provides insights into the influences of anthropogenic land use on microbial community structure and functions of stream benthic biofilms.

and rivers and ultimately impacts the suitability of water resources [11][12][13] . In particular, the nitrogen metabolism in streams is one of the most important ecosystem functions to remove nitrogen 24 and consequently protect downstream lakes from eutrophication. Moreover, bacteria can account for a substantial portion of the transformation and use of organic matter in riverine systems 29,30 , with implications for ecosystem processes, food webs, and greenhouse gases 31 . Thus, the major metabolisms in riverine systems are of importance because of the potential for riverine ecosystems to store, mineralize, and transport organic matter and nutrient before reaching the downstream aquatic ecosystems [32][33][34] . Considerable work focused on export, source, transport, retention, and metabolism of nitrogen and carbon in streams that are impacted by anthropogenic activities 27,31,34,35 . However, the influences of anthropogenic land use on microbial community structures and functional potentials of stream biofilms are remain unclear.
Aquatic ecosystems face increasing pressures and threats from disturbances of watersheds associated with land use 36,37 . Understanding the microbial community and the major metabolic pathways is crucial to obtain insights into ecosystem structures and processes 38,39 . In this work, we studied the benthic biofilms in streams flowing through forests and anthropogenic land use (agriculture and urban) areas in the Erhai Lake watershed in southwestern China. The microbial communities of benthic biofilms in upstream and downstream of anthropogenic land use area were compared. We determined the bacterial communities using high-throughput 16S rRNA gene sequencing and investigated carbon and nitrogen metabolic pathways from the PICRUSt (Phylogenetic Investigation of Communities by Reconstruction of Unobserved STates) predicted functional genes encoding the enzymes required in carbon and nitrogen metabolic pathways. Although metagenomics only reflects a potential rather than a realized functional capacity, our data offered an approach to understand and discuss the anthropogenic influences on the microbial community structures and functions in stream benthic biofilms. We hypothesized that (1) anthropogenic land use leads to nutrient and organic carbon increase and alters microbial community structure, and (2) anthropogenic land use alters the nitrogen and carbon metabolic pathways, which are associated with the community structure and nutrient concentrations.

Methods
Study area. The study streams are located on Cangshan Mountain, Yunnan Province, southwestern China ( Fig. 1). Cangshan Mountain has 19 peaks. The summit of Cangshan Mountain has an altitude of 4122 m, and the other 18 peaks are all above 3500 m. In the valleys between peaks, there are 18 streams originating from the eastern slope of Cangshan Mountain that flow into Erhai Lake, which has a water level of 1964 m and is the 65 th largest lake and the 7 th largest freshwater lake in China. The forest coverage of Cangshan Mountain is above 90%, and the stream reaches flowing through the mountain are difficult for humans to access and are still pristine. However, when the streams flow out of the mountain, they flow through a narrow plain (around 5 km wide) between Cangshan Mountain and Erhai Lake (Fig. 1). The plain has been intensely exploited for agriculture and urban development. In this study, we sampled five streams. In each stream, the upstream site is located at the foot of Cangshan Mountain, where the stream flows outward from the mountain and is not influenced by anthropogenic activities. The downstream site is located at the mouth of the stream, where the stream flows into the lake and is severely influenced by agriculture and urbanization on the plain.
Sampling and chemical analyses. Samples were collected in August 2016. At each sample site, we randomly sampled 6-9 submerged rocks from the river cross-section. The benthic biofilm was removed by rigorously brushing an area of 4.5 cm in diameter from the upper surface of each stone using a sterilized nylon brush (changed between samples) and rinsing the slurry with sterile water. For microbial analysis, 10 mL of mixed slurry was filtered through 0.2 μm membrane filters and immediately frozen in liquid nitrogen in the field and transported to the lab for DNA extraction and subsequent analyses. For the biofilm chemical analysis, 20 mL of mixed slurry was filtered through 0.2 μm GF/F filters that were dried at 60 °C. Biofilm carbon (Bio-C) and biofilm nitrogen (Bio-N) was analyzed using an element analyzer (LECO 628). Biofilm phosphorus (Bio-P) was analyzed spectrophotometrically by ascorbic acid colorimetry following digestion with potassium persulfate (EPA 365.3). At each sample site, water samples were collected for nitrogen and phosphorus analysis. Total nitrogen (TN) was analyzed by ion chromatography after persulfate oxidation (EPA 300.0). Nitrate (NO 3 − ) was determined by ion chromatography (EPA 300.0). Total phosphorus (TP) was analyzed using the ammonium molybdate method after oxidation (EPA 365.3). Soluble reactive phosphorus (SRP) was analyzed using the ammonium molybdate method (EPA 365.3). Ammonium (NH 4 + ) was analyzed using the indophenol colorimetric method (EPA 350.1). Dissolved organic carbon (DOC) was analyzed on a TOC analyzer (TOC-VCPH; Shimadzu Scientific Instruments, Columbia, Maryland).
DNA extraction, PCR, and sequencing. Genomic DNA was extracted using the PowerSoil DNA Isolation Kit (MoBio, Carlsbad, CA, USA) following the manufacturer's protocols. The V3-V4 regions of the 16 S rRNA gene were amplified using 338F-ACTCCTACGGGAGGCAGCA and 806R-GGACTACHVGGGTWTCTAAT (Invitrogen, Vienna, Austria). PCR was performed with a model 2720 thermal cycler (ABI, USA) using the following program: 1-min hot start at 80 °C, 94 °C for 5 min followed by 30 cycles of denaturation at 94 °C for 30 s, followed by annealing at 52 °C for 30 s and 72 °C for 1 min 30 s, with a final extension step at 72 °C for 10 min. Amplified DNA was verified by electrophoresis of the PCR mixtures in 1.0% agarose in 1X TAE buffer and purified using the Gel Extraction Kit (Qiagen, Hilden, Germany). Samples were sent for sequencing on a MiSeq platform.
Analysis. Bacterial 16 S rRNA sequence data were cleaned using the software package QIIME 40 and then clustered into operational taxonomic units (OTUs) with a complete linkage algorithm at a 97% sequence identity level. The metagenomes were predicted from 16 S data using PICRUSt 41 . The functional genes associated with nitrogen and carbon metabolism were identified from KEGG (Kyoto Encyclopedia of Genes and Genomes) database 42 . Redundancy analysis (RDA) was conducted to reveal the association of the microbial communities in relation to nutrient factors based on OTU and functional gene abundances (Vegan package 2.4 in R 3.3.2). Correlation analysis was conducted to assess the relationships between the nutrient factors and nitrogen cycle pathways, and the P-value was adjusted with a Bonferroni correction (psych package 1.7.5 in R 3.3.2). PERMANOVA was used to test whether upstream and downstream sites harbor significantly different microbial communities or metagenomes (using PAST 3.0). Mantel tests were run to assess correlations between functional and taxonomic community dissimilarity matrices based on Bray-Curtis distance. Non-metric multidimensional scaling (NMDS) was applied to reveal differences in community composition between microbial assemblages in upstream and downstream sites (using R 3.3.2 and Vegan package 2.4-1). For network analysis, the relative abundances of the OTUs in each sample were used to construct matrices for visualizing interactions between OTUs in the networks (microbial networks of upstream and downstream sites). A Spearman correlation coefficient R score and a P-value were calculated pairwise between the OTUs (the1000 most abundant OTUs) using the Hmisc package (version 4.0-1) in R (version 3.3.2). Only strong (Spearman's correlation coefficient R > 0.95 or R < −0.95) and strongly significant (P < 0.01) correlations were considered. These correlations were visualized using Cytoscape (version 3.5.0). Each node represents an OTU, and each edge represents a strong and significant correlation. To describe the network topology, a set of node/edge metrics was analyzed using the Network Analyzer plugin within Cytoscape 43 . The modular structures of the networks were analyzed using the ClusterMaker in Cytoscape. Modularity values > 0.4 suggests that the network is modular 44 .

Results
Nutrient changes in stream water and benthic biofilms. In the stream water, downstream sites had significantly higher concentrations of TN and NO 3 − than upstream sites (Table 1, t-test, P < 0.05), and the nitrogen was mainly NO 3 − in both the upstream and downstream sites. Moreover, downstream sites had higher concentrations of TP, SRP, and DOC than upstream sites. Accordingly, downstream sites had significantly higher C, N and P contents in the benthic biofilms (Bio-C, Bio-N and Bio-P) than upstream sites (Table 1, t-test, P < 0.05).  Table 1. Nutrient variations in stream water and benthic biofilms of upstream sites and downstream sites (Means ± SD, abbreviations as in text).
Bio-N was positively correlated with TN, NO 3 − , and TP ( Table 2, P < 0.05). Bio-P was positively correlated with TP and SRP but negatively correlated with DOC (Table 2).

Microbial community variation-taxonomic and functional.
After quality filtering, the valid sequences obtained were assigned to 7869 OTUs. The proportion of detected OTUs that overlapped between upstream and downstream was 57.15%. The proportions of detected OTUs that were unique to upstream and downstream sites were 21.17% and 21.69%, respectively. All the OTUs were associated with 28 phyla. The dominant phyla in upstream sites were Proteobacteria (relative abundance of 45.9%) followed by Cyanobacteria (39.3%) and Bacteroidetes (7.5%) (Fig. 2). The dominant phyla in downstream sites were Proteobacteria (62.9%) and Cyanobacteria (24.8%) (Fig. 2). Downstream sites had significantly higher abundances of Proteobacteria than upstream sites (Fig. 2, t-test, P < 0.05), while upstream sites had a higher abundance of Cyanobacteria than downstream sites. NMDS analysis showed that the upstream and downstream sites had different taxonomic communities in the benthic biofilms (Fig. 3a). The Bray-Curtis distance based dissimilarity tests confirmed that the differences were significant (PERMANOVA, F = 1.630, P = 0.028). RDA indicated that TN, NO 3 − , Bio-C, and Bio-N were the most significant nutrient factors associated with taxonomic composition (Monte Carlo test, P < 0.05, Fig. 4a). The first two axes explained 41.77% of the taxonomic information (RDA 1: 22.56%; RDA 2: 19.21%).
From the PICRUSt-predicted metagenome results, 5994 functional genes were detected. NMDS analysis showed that the upstream and downstream sites had different functional gene compositions in their benthic biofilms (Fig. 3b). However, the Bray-Curtis distance-based dissimilarity tests showed that the differences were nonsignificant (PERMANOVA, F = 2.203, P = 0.104). RDA analysis indicated that TN, NO 3 − , SRP, Bio-C, Bio-N, and Bio-P showed significant relationships with the functional gene composition of the biofilm communities (Monte Carlo test, P < 0.05, Fig. 4b). The first two axes explained 80.24% of the metagenomic variance (RDA 1: 63.92%; RDA 2: 16.32%).
To evaluate whether taxonomic differences between microbial communities affect their functional potential, we conducted Mantel tests to analyses the relationships between Bray-Curtis dissimilarities based on the taxon abundances and the functional gene abundances of microbial communities in stream biofilms. Mantel correlation tests revealed significantly positive correlations between the functional and taxonomic dissimilarities in upstream biofilms (Fig. 5, R 2 = 0.659, P = 0.004) and downstream biofilms (Fig. 5, R 2 = 0.399, P = 0.049). Moreover, the slope of the linear relationships was smaller in the downstream sites than in the upstream sites.
Bacteria co-occurrence. We compared the co-occurrence networks of upstream and downstream microbial communities (Fig. 6). The interactions between OTUs (the top 1000 OTUs) in the upstream and downstream sites were described by the topological parameters of the networks (Table 3). In comparing the topological parameters of these two networks (Table 3), the upstream microbial network exhibited a greater values of network diameter, network centralization, network heterogeneity, and characteristic path length. Both networks had a strongly clustered topology, but the upstream microbial network had a lower modularity than the downstream microbial network (Table 3). Potential metabolisms. The N cycle involves four reduction pathways (assimilatory nitrate reduction, dissimilatory nitrate reduction to ammonia, denitrification, and nitrogen fixation) and two oxidation pathways  Table 2. Correlations between nutrient concentrations in stream water and benthic biofilms (abbreviations as in text). Note: "*" indicates level of statistical significance at P < 0.05, "**" indicates P < 0.01.

Figure 2.
Relative abundance of microorganisms at the phylum level in stream biofilms. Only the phyla with a relative abundance >1% are shown. "Others" represent the unsigned OTUs and the phyla with a relative abundance <1%. Only the relative abundance of proteobacteria was significant between the upstream and downstream sites ("*" indicates level of statistical significance at P < 0.05).
(nitrification and anammox) (Fig. 7). Compared to upstream sites, downstream sites had higher relative abundances of some genes which encoding enzymes associated with denitrification, nitrification, and anammox (Figs 7 and 8), suggesting that downstream sites had higher nitrogen metabolic pathways associated with nitrogen removal. Downstream sites also had a higher relative abundance of nitrogen fixation (Figs 7 and 8). Upstream      Table 3. Topological parameters of microbial networks in the upstream and downstream sites.
sites had a higher relative abundance of assimilatory nitrate reduction than downstream sites (Figs 7 and 8).
Both upstream and downstream sites had a high relative abundance of dissimilatory nitrate reduction (Fig. 8).
Correlation analyses (Table 4) showed that anammox was positively correlated with TN, NO 3 − , Bio-C, Bio-N and Bio-P. Assimilatory nitrate reduction was negatively correlated with TN, NO 3 − , and Bio-C. Denitrification was positively correlated with Bio-C, Bio-N, and Bio-P. Dissimilatory nitrate reduction was negatively correlated with NH 4 + . Nitrification was positively correlated with Bio-C and Bio-P. For the carbon metabolism, we detected 234 functional genes encoding the core carbon metabolism pathways (Fig. 9). There were 44 functional genes that had higher relative abundances in microbial assemblages in downstream sites than in upstream sites, while only 7 functional genes had higher relative abundances in upstream then in downstream microbial assemblages. This result suggested that microbial assemblages in benthic biofilms had a stronger carbon metabolism in downstream sites than in upstream sites (Fig. 9).  42 . The red box indicates the enzyme that had a significantly high relative abundance in upstream sites. The blue box indicates the enzyme that had a significant high relative abundance in downstream sites.

Figure 8.
Relative abundances of functional genes encoding the enzymes that catalyze nitrogen cycling pathways. Only the relative abundances of anammox and assimilatory nitrate reduction were significantly different between upstream and downstream ("*" indicates level of statistical significance at P < 0.05).  Table 4. Relationships between the relative abundance of functional genes in benthic biofilms associated with nitrogen cycle pathways and nutrients in stream water and benthic biofilms (abbreviations as in text). Note: "*" indicates level of statistical significance at P < 0.05, "**" indicates P < 0.01. Figure 9. Carbon metabolism pathways (map01200) performed by KEGG mapper (http://www.genome.jp/ kegg/mapper) 42 . The red line indicates the pathway that had significantly high relative abundance in upstream sites. The blue line indicates the pathway that had significantly high relative abundance in downstream sites.

Discussion
Anthropogenic impacts on nutrients and biofilm microbial communities. Land use significantly controls the export of nutrient and organic matter from catchments to aquatic ecosystems [11][12][13] . Agriculture and urban use increase NO 3 − , NH 4 + , and SRP concentrations as a result of excess fertilizer application 45 , as well as sewage and septic inputs 46,47 . Increased nutrient loading can lead to reduced water quality 48 and altered stream communities 49,50 . Biofilm nitrogen and phosphorus contents increase with increasing nutrient availability in streams 51,52 . In our study, downstream sites had higher nutrient concentrations in both stream water and benthic biofilms than upstream sites. The nutrient contents in benthic biofilms were closely related to stream-water nutrient concentrations. N and P availability may limit the growth of both biofilm autotrophs and heterotrophs 53,54 . It has been demonstrated that the community composition, diversity, biomass, and stoichiometry of biofilms change significantly with available nutrients in streams 23,28,55,56 . The composition of microbial communities can change in streams influenced by catchment land-use features [57][58][59][60] . In our study, downstream sites had a high abundance of proteobacteria while upstream sites had a high abundance of cyanobacteria. The differences in microbial communities were closely correlated with nutrients in both the stream water and biofilms.
Different biomes harbor distinct microbial assemblages [61][62][63] . It has been demonstrated that functional beta diversity is strongly correlated with taxonomic beta diversity across soil microbial communities 61 . In our study, both upstream and downstream sites had significant linear relationships between OTU dissimilarities and functional gene dissimilarities. The high linear regression slope suggested that the microorganisms in upstream sites had more distinct functional genes than in downstream sites. The microorganisms in downstream sites share more functions than upstream sites, indicating high functional gene redundancy. This also indicated that intense anthropogenic impacts can cause the reduction of functional diversity [64][65][66] .
Co-occurrence patterns of organisms have been evaluated to reveal community assembly rules and interaction networks in highly complex systems 67 . Network topology (the node distributions and interactions) can affect the stability of the system 38,68 . In previous macroecological studies, communities with tightly connected species were shown to be more susceptible to disturbance 69,70 . The highly connected microbial network of microorganisms in upstream sites suggested that microorganisms have stronger co-occurrence relationships in upstream sites than in downstream sites. A minor disturbance might cause large impacts on upstream network because of the complex relationships among OTUs. Furthermore, the microbial assemblages in both upstream and downstream sites exhibited a modular structure (Table 3, modularity values > 0.4 suggest that the network is modular, Newman 2006). In a biotic network, highly interconnected species are grouped into a module. In a module, species interactions are more frequent and intensive than in the rest of the community 38,44 . Niche differentiation can lead to higher modularity in a community 38 . In our study, the high modularity of both upstream and downstream microbial networks suggested that the upstream and downstream sites of the stream had heterogeneous habitats offering rich niches for microorganisms, especially in downstream sites which receive more nutrient and organic matter from the catchment.
Anthropogenic influences on biofilm metabolic potentials. Microorganisms in biofilms are major components of the biogeochemical cycles in stream ecosystems 71 . Human activities modify land use and affect quantity and quality of nutrient and organic matter in streams 3,72 . As a result of anthropogenic impacts on catchments, changes in the community composition of benthic biofilms result in the alteration of stream ecosystem functions 6,73 .
Nitrogen has various chemical forms and is cycled by a suite of biogeochemical processes 74 , which are catalyzed by various microbe-excreted enzymes 75,76 . The assimilatory pathways require energy and start with the reduction of NO 3 − to NO 2 − and then to NH 4 + , which is highly bioavailable and can be readily used by the cell for the synthesis of amino acids and nucleotides or can be transformed by nitrification or anammox 77 . Dissimilatory nitrate reduction to ammonia uses N-compounds to provide energy (ATP) to microbes and is another catabolic pathway that can retain the nitrogen in the system in a bioavailable form (NH 4 + ) for further biological processes 78,79 . Nitrogen fixation is an energetically demanding process requiring 16 molecules of ATP to break the triple-N bond. Only a few organisms, including cyanobacteria and certain other genera, can carry out this process 80,81 . Low NH 4 + can stimulate some N 2 -fixing bacteria [82][83][84] . Nitrification is an essential process in the nitrogen cycle performed by nitrifiers converting ammonium to nitrate. As a competing process to nitrate reduction pathways, denitrification is the main biological process for the removal of N from freshwater systems 79,85 . Anammox plays an important role in the nitrogen cycle, turning nitrite and ammonia to dinitrogen 86 , and is another microbial process that releases fixed nitrogen from the environment as dinitrogen [87][88][89][90] . In our study, a higher relative abundance of denitrification, nitrification, and anammox suggested that nitrogen removal was stronger in downstream sites than in upstream sites.
Furthermore, nitrogen metabolic pathways are regulated by environmental redox states 91,92 , nutrient availability, and other abiotic factors. High concentrations of NO 3 − can stimulate microorganisms to increase populations containing the functional genes associated with nitrate reduction. Anammox is influenced by NO 3 − concentration, NH 4 + concentration, C/N ratio, and pH 93 . Low NH 4 + can stimulate functional genes associated with nitrogen fixation [82][83][84] . Nitrification can be affected by temperature, salinity, light, organic matter concentration, substrate concentrations, pH, and oxygen concentration 94 . Our study showed that anammox was positively correlated with TN, NO 3 − , Bio-C, Bio-N and Bio-P. Assimilatory nitrate reduction was negatively correlated with TN, NO 3 − , and Bio-C. Dissimilatory nitrate reduction was negatively correlated with NH 4 + . Denitrification was positively correlated with Bio-C, Bio-N and Bio-P. Nitrification was positively correlated with Bio-C and Bio-P (Table 4). These results suggest that human land use significantly increases nutrient and organic carbon in stream water and biofilms, and the increased nutrient had severe impacts on the nitrogen cycle in streams.
In stream ecosystems, the processing of organic matter is another major function that bacterial communities mediate 29 . Carbon transported by streams provides a major source of energy for aquatic food webs and is a significant component of the global carbon cycle 31,32,95 . Anthropogenic land use contributed significantly to accelerated transport of carbon from terrestrial environments to aquatic ecosystems 31, [96][97][98] . On the other hand, anthropogenic land use can accelerate the transition of streams from transporters to transformers of carbon 31 . In stream ecosystems, a substantial portion of the organic matter are derived from terrestrial environments and serve as a substantial resource fueling riverine bacterial communities 29 . Small streams can mineralize large quantities of the organic carbon load because of the bacterial carbon metabolism, which results in increased losses of CO 2 to the atmosphere and with implications for ecosystem services 8,34,99 . Moreover, the increased nutrient availability can further stimulate whole-stream metabolism of terrestrial carbon 8,100 . In our study, the overall carbon metabolism of the benthic biofilms was enhanced by the intense anthropogenic activities.

Conclusions
Anthropogenic activities have caused the increase of nitrogen and carbon loads in streams and rivers. A portion of increased nitrogen and carbon loads are exported to downstream lake and coastal waters, while a substantial proportion are also transformed along streams and rivers 35,101 . Furthermore, anthropogenic land use accelerates the transition of streams from transporters to transformers of nutrient and organic carbon 31 . In our study, anthropogenic influenced downstream sites had higher relative abundances of functional gene encoding nitrogen and carbon metabolisms compared to upstream sites with relatively pristine catchments. The results allow us to understand how microbial community structures and metabolism response to influences of anthropogenic land use.