Sampling was carried out in the Red River and its estuary on the Disko Island, West Greenland (69°N) during August 2013 (Figure 1). Around the time of sampling the river flow was 5.7 m³s⁻¹. The bedrock consists of iron-rich basalt. As the glacier and the stream erode the bedrock the iron precipitates and gives the marked red color of the river. The river drains directly into the Disko Bay and the freshwater and sediments supplied from the river are mixed with the saline bay waters under varying wave influence. A river plume of high concentrations of sediment is often visible indicating how the supplied sediment is dispersed.

Five locations were sampled in the river with three replicates at each site (Figure 1). The top sample (R1) being just upstream of an outlet from an adjacent proglacial lake and the second sample (R2) at the outlet from the proglacial lake. The third sample (R3) being at another outlet to the river supplying water directly from the glacier and the fourth sample (R4) 100 m downstream of R3. Sample R4 was collected at the eastern bank on the opposite site of the upstream outlet from the glacier stream (R3), while all other river samples were collected at the center of the river. The last and fifth sample (R5) was collected close to the river mouth. The distance of each sample to the river mouth is supplied in Table 1. In the bay, sampling was done along three transects perpendicular to the coast (Figure 1). Each transect consisted of four sampling locations at distances of 100, 300, 700, and 1100 m from the river mouth. At each distance two samples were collected, one surface sample at 0.5 or 1 m from the surface and one deep sample 1 m from the bottom. At water depths above 20 m, the deepest water sample was collected at 20 m depth. In transect 1, 100 m into the estuary (E100) the deep sample is missing so that there are only two replicates (transect 2 and 3) of E100 samples. At E700 one sample, which should have been sampled at 20 m, was sampled at 1 m depth, so that there are 4 replicates of surface samples and 2 replicates of deep samples for E700. Water was sampled by grab sampling using sterile 50 ml syringes (Sarstedt, Germany) either collecting water directly from the river or collected from a 5 L Niskin water sampler (KC Denmark, Denmark) that had been filled at the sampled depth. The 50 ml water samples collected in the syringes were forced through Sterivex™ filters (Merck Millipore, MA, USA) and the filters were afterwards partly dried by forcing air from the syringes through the filters. The filters were frozen and kept at −20°C until analysis.

Temperature, turbidity, and oxygen saturation were measured at all sites using a YSI 6600-V2 CTD sensor with attached probes (YSI, OH, USA). In the bay, the size (in equivalent spherical diameter, ESD), total area and total number of particles were measured in 6 mL water using a laser sheet camera system, the Pcam (Markussen et al., 2016). Individual water samples were taken at all locations and transferred to new 100 ml polyethylene bottles, frozen as quickly as possible and shipped to Copenhagen for further analysis. The total nitrogen (TN) and dissolved organic carbon (DOC) were determined on a Shimadzu TOC-V total organic carbon analyzer (Shimadzu, Japan) with a TNM-1 total nitrogen measuring unit and pH was measured using a Radiometer Analytical SAC90 autosampler (Hach, CO, USA). DOC measurements were based on triplicate measurements. A standard curve using 1000 ppm sodium hydrogen phthalate with concentrations ranging from 0 to 5 ppm were made and a 100 ppm certified Total Organic Carbon (TOC) standard (SCP Science, QC, Canada) was diluted to 1 ppm for use as reference.

DNA was extracted from the Sterivex™ filters using the PowerWater© Sterivex™ DNA extraction kit (MO BIO Laboratories, CA, USA), using the protocol provided by the manufacturer. The extracted DNA was stored at −80°C until library preparation.

The nucleic acid concentrations of all samples were assessed by spectrophotometer (Nanodrop® ND-1000, Saveen Werner, SE) to be within the range of 3–5 ng μl⁻¹. DNA was then amplified in triplicate using universal prokaryotic primers targeting the variable region V4 of the 16S rRNA gene (Caporaso et al., 2011), forward primer 515F (GTGCCAGCMGCCGCGGTAA) and reverse primer 806R (GGACTACHVGGGTWTCTAAT), using the HiFi polymerase (PCR-Biosystem, UK). The primers were supplied with 12 distinct barcode sequences of 4–6 bases each and combined as differential sets, thus labeling the samples with individual differently tagged sequences. All PCR runs included triplicate positive (E. coli) and negative (dd H₂O) controls. The resulting PCR products (350 bp) were quality controlled by quantification of concentrations using the Qubit® 2.0 dsDNA HS Assay Kit (Life Technologies, CA, USA) and visual inspection of band size following gel electrophoresis. The amplified DNA was then purified using the HighPrep™ PCR size selective carboxyl coted magnetic beads (Magbio, MD, USA). The resulting DNA (average concentration 19.3 ng μl⁻¹) was then pooled and ligation of adaptors was performed according to manufacturer's instructions following the Low Sample (LS) Protocol (TruSeq DNA PCR-Free Sample Preparation Guide, Illumina, CA, USA) with minor modifications. Overhangs on the 3′ ends were removed and 5′ ends filled in by end repair, performed as described in the protocol on 1 μg DNA. Size selection was replaced with a clean-up step with magnetic bead based chemistry (HighPrep™ PCR, CleanNA). A volume of 100 μl from the end repair reaction was purified according to manufacturer's instructions and subsequently eluted in 20 μl molecular biology grade water (MO BIO Laboratories, CA, USA). Following this step, the 3′ ends were adenylated (adding an “A” nucleotide) to prevent them from ligating to one another during the ligation reaction. Then adaptors with a “T” overhang were ligated onto the DNA fragments of the two assemblages as described in the protocol. The assemblage was ligated with the AD012 index adaptor (CTTGTA), and then subjected to a clean-up step with purification beads provided in the kit. Quality control of the ligated amplicons (~400 bp) was performed by PCR amplification using primers targeting the index adaptor followed by gel electrophoresis and visual inspection. Finally, the amplicon assemblage was diluted to a concentration of 3.3 ng μl⁻¹ and sequenced with MiSeq 250PE (Illumina), adding 30% PhiX DNA. Demultiplexed merged reads are deposited in the NCBI Sequence Read Archive (SRA) database under SRA accession SRP076603.

The sequencing data was quality checked using FastQC (Patel and Jain, 2012) and read pairs were merged with the paired-end read merger PEAR (Zhang et al., 2014). Only properly merged reads were used for downstream analysis. Merged reads were processed using Qiime version 1.8.0 (Caporaso et al., 2010a). Demultiplexing with split_libraries_fastq.py was performed with quality filtering at phred threshold ≥ 20. Chimeric sequences were removed from demultiplexed data with USEARCH uchime reference-based chimera removal using the Greengenes database from May 2013 as reference (Edgar et al., 2011). Chimera check removed 12.3% of sequences. Operational taxonomic units (OTUs) were subsequently picked based on 97% identity using de novo OTU picking, which also includes taxonomy assignment using PyNAST alignment against the Greengenes core set of 16S rRNA sequences (Caporaso et al., 2010b). Sequences only represented once in the dataset were removed, which reduced the dataset with 13.7%.

Shannon indices (Shannon, 1948), Chao1 richness (Chao, 1984) and rarefaction plots were computed using alpha_rarefaction.py. Chloroplast sequences were removed and samples were rarefied to the shallowest sample depth of 12,180 sequences per sample with R version 3.1.0 (R Development Core Team, 2008) and R package Vegan (Oksanen et al., 2015). BIOENV analyses were used to assess how well the community structure was explained by environmental variables using non-factorial metadata (Table 1; Clarke and Ainsworth, 1993). For BIOENV analysis the Vegan package was used to create distance matrices of environmental data (Euclidean distances) and community composition (Bray-Curtis distances), which were then compared through Spearman's rank coefficients. DOC, TN, temperature, salinity, oxygen saturation, and turbidity were log transformed prior to analysis. Depth was not included for BIOENV analyses including only river samples as depth was constant and particle data was not included for any analyses including river samples, as the data was not available. LabDSV package in R was used for non-metric multidimensional scaling (NMDS) and indicator species analysis. Indicator species are here denoted indicator OTUs and are defined as OTUs having a higher abundance at one site compared to other sites with indicator values d ≥ 0.3 at a significance level of p ≤ 0.05. Indicator values are a product of relative abundance of an OTU in samples from one site (between 0 and 1) and the relative average abundance of that OTU across all sites (Dufrene and Legendre, 1997). The used indicator species and indicator OTU concept in this study are not equal to the Indicator Species concept representing species that are markers for certain environmental variables in an ecosystem. NMDS analyses were conducted using Bray-Curtis distance matrices. NMDS stress values are included in Figure 3.

Dissolved organic carbon (DOC) and total nitrogen (TN) concentrations in the river were in the same range as those in the estuary (Table 1). River site R2 by the lake outlet had the highest DOC and TN concentrations of all samples in the river and estuary. pH values in the river were slightly lower than in the estuary except for R2, which had a higher pH comparable to the estuary samples. Temperature ranged from 1.7°C at 20 m depth 700 m into the estuary to 12.6°C at river site R2. Temperatures were generally lower in the deep water samples from the estuary compared to the surface samples. Salinity in the river samples was 0.01 PSU for all samples except for R2 where it was 0.023 PSU. The higher salinity in the water from the lake can be explained by the accumulation of ions in the lake due to longer WRT in the lake compared to the river allowing for evaporation of water from the lake. The longer WRT may also explain the higher temperature at site R2. In the estuary, salinity was consistently lower in the shallow water samples compared to deep water samples at the same distance from the river mouth. This was expected from the lower density of the freshwater from the river being mixed into the estuary. Turbidity across all samples, excluding river site R2, ranged from 8.0 to 19.8 NTU, while it was remarkably lower at R2 (1.2 NTU). Camera data from the estuary showed that particle mean diameter was generally higher at shallow depths compared to deep water samples.

Illumina sequencing of variable region V4 of the 16S SSU rRNA gene from a total of 38 samples resulted in 462,840 individual sequences after rarefaction to 12,180 sequences per sample, which were binned into 63,624 unique OTUs (97% sequence identity). The number of observed OTUs was not exhausted at this level of rarefaction (Supplementary Figure 1). Shannon indices for the river samples ranged from 5.6 to 10.8 and Chao1 richness in the river ranged from 1408 to 19,117 OTUs per sample (Figure 2). The alpha diversity of the bacterial community represented by both Shannon indices and Chao1 richness decreased at R2 and R3, the sites at which the lake and the glacier stream drains into the river (Figures 1, 2). The alpha diversity rose again at R4, ca. 100 m from the glacial input site. At the river mouth (R5), the diversity increased again and reached a similar level to the first river site (R1) upstream of the glacier and lake input sites.

Shannon indices for the estuary samples ranged from 5.4 to 10.6 while Chao1 richness varied from 2589 to 19,021 OTUs (Figure 2). The diversity and richness were higher in the estuary than in the glacier stream and lake input samples, and slightly lower than in the remaining river samples. There was no apparent pattern in the difference in diversity and richness attributed to different depths of the estuary, sample sites or the distance to the river mouth.

The samples from the first site of the river (R1), upstream of the lake and glacier stream outlets to the river, clustered with samples from the bottom site of the river (R5; Figure 3). These two sites also shared a high number of indicator OTUs (Figures 4A,E) and showed similar diversity and richness (Figure 2). The river site by the proglacial lake outlet (R2) clustered with the river site at the glacier stream outlet (R3). These sites, R2 and R3, also had lower diversity than the other river sites (Figure 2).

In the estuary, the bacterial communities clustered according to sample site for the two sites that were farthest into the estuary (E700 and E1100). The samples from the sites closest to the river mouth (E100 and E300) were dispersed across NMDS 1 and 2 (Figure 3). The samples did not cluster according to sample depth. Samples from the sample sites closest to the river mouth (E100 and E300) clustered more closely with river samples than the samples farthest from the river mouth (E700 and E1100).

BIOENV analysis showed that the total community as well as the non-indicator and indicator OTUs in the river correlated significantly with turbidity at p ≤ 0.05. The strongest correlation was found between the river non-indicator OTUs and turbidity with a Spearman's rank correlation coefficient of 0.586.

BIOENV analysis of the estuary community showed no significant correlations with environmental variables (Table 2).

The number of indicator OTUs in the river ranged from 158 at the second-to-last site of the river (R4) to 678 at the input site from the proglacial lake (R2), which also had the highest percentage of top indicator OTUs (Indicator Value = 1; Table 3). There was a high number of shared indicator OTUs between the top and bottom of the river (Figures 4A,E). Two hundred and eight indicator OTUs from R1 were found at R5 while only 14, 31, and 28 indicator OTUs from R1 were found at R2, R3 and R4 respectively. Taxonomic composition of indicator OTUs at Class level showed similar fractions of Flavobacteria and Gammaproteobacteria across river samples. A greater fraction of Actinobacteria were found in R2-R4 while very few Acidobacteria were found in these samples compared to R1 and R5 where also a higher fraction of Unknown were found (Figure 5). A number of indicator OTU sequences at the uppermost river site (R1) showed similarity to members of Rhizobiales isolated from plant roots and soil (Lee et al., 2005) as well as to strict anaerobes such as Caldilinea, Anaerolineaceae (Yamada et al., 2006), and Desulfobacteraceae (Garrity et al., 2006; Figure 5).

Indicator OTUs identified at the outlet from the lake and glacier stream were found in low numbers at the other river sites (Figures 4B,C). The lake outlet site (R2) had the highest number of indicator OTUs and percentage of top indicator OTUs (Table 3). A number of taxa known to be psychrophilic, such as Moritella (Urakawa et al., 1998), Polaribacter (Gosink et al., 1998), Oleispira (Yakimov et al., 2003), Crocinitomix (Bowman et al., 2003), and Psychromonas (Mountfort et al., 1998) were found among the best matches for the indicator OTUs from the lake outlet, unlike at the other river sites.

The distribution of estuary indicator OTUs showed a different pattern than the river indicator OTUs (Figures 4F–I). The number of indicator OTUs in the estuary was generally lower than at the river sites. An exception to this was the outermost estuary sample (E1100), which had a number of indicator OTUs comparable to the river sites (Table 3). No top indicator OTUs were found in any of the estuary samples, meaning that no OTUs from the estuary were unique to any of the sample sites. The indicator OTUs for each sample site in the estuary were found only in low numbers at the other sites both in the river and the estuary (Figure 4) and the taxonomic composition at Class level was less similar among the estuary samples than among the river samples (Figure 5).

On average, the bacterial communities in the estuary were made up of 23% river indicator OTU sequences (Figure 6). There was an overall decreasing contribution of river indicator OTUs in the estuary sites with 26–27% river indicator OTUs closest to the river mouth at E100 sites, 17–25% at E300 sites and 8–10% at E1100 sites. E700 sites were exceptions with 22–52% of the community being river indicator OTUs (Figure 6). Closer to the river mouth at sites E100 and E300 there was a larger fraction of the indicator OTUs from the top of the river (R1) and river mouth (R5), except for the deep sample at E100, where the distribution of river indicator OTUs was similar to the E700 estuary sites (Figure 6). In E700 both in the deep and surface samples R2 (lake outlet) indicator OTUs were more abundant than in the other estuary samples, and were more abundant than indicator OTUs from any other river sites (Figure 6). At the estuary sites farthest from the river mouth (E1100) the samples had the highest fraction of non-river indicator OTU sequences (Figure 6).

Shannon indices for the river samples at the lake outlet (R2) and glacier stream outlet (R3) were comparable to a recent study of 87 small streams and rivers in the La Côte-Nord region of Québec, Canada, where OTUs were clustered with the same method as in the present study (Ruiz-González et al., 2015). The remaining river samples had slightly higher diversity than found in previous studies (Galand et al., 2006, 2008; Crump et al., 2009). The difference from less recent studies is likely due to the difference in the technologies applied and the resulting lower number of sequences in the previous studies. Together with sequencing technologies, which have changed dramatically in the last decade, OTU clustering has shown to have a great impact on the detected alpha diversity (Sinclair et al., 2015). Therefore, the comparison of alpha diversity among studies should be interpreted with care. Our results were obtained with the use of the Qiime pipeline (Caporaso et al., 2010a), which has shown to create a larger number of OTUs when compared to other popular clustering methods (Sinclair et al., 2015). Consequently, we might detect a higher diversity because of the clustering method used.

In an extensive study of freshwater networks, small streams were shown to have higher Shannon indices than larger rivers (Ruiz-González et al., 2015). This was attributed to the common terrestrial origin of the microbial community resulting in an initially high diversity in the small streams originating from the surrounding soil (Ruiz-González et al., 2015). This is in contrast with our results showing less diversity in the glacier stream compared to the main river (Figure 2). The lower diversity might be due to the origin of this stream in a glacier rather than subsurface groundwater and surface runoff as described in the above-mentioned study (Ruiz-González et al., 2015). This notion is supported by the fact that the diversity of the glacier stream outlet (R3) is comparable to those recently described for proglacial lakes (Peter and Sommaruga, 2016). Also, the indicator OTUs from the glacier stream were similar to taxa commonly found in freshwater and marine environments (details not shown). These results emphasize the importance of high spatial resolution for assessing the origin of the metacommunity in a complex freshwater network. In this study, the diversity along the network does not follow a uni-directional pattern (Figure 2). Our results illustrate that the origin and structuring of the microbial community might be very different from one network to another. How glaciers and glacier streams affect the metacommunity of freshwater networks is a highly relevant topic yet to be investigated.

The Chao1 richness in the river samples was higher than in the previous papers focused on large Arctic rivers (Galand et al., 2006, 2008) but comparable to that described in a recent paper using the same sequencing platform (Niño-García et al., 2016).

The drop in alpha diversity at the input sites from the lake and the glacier (R2 and R3, Figure 2) shows that the lake and the glacier stream input less diverse bacterial communities into the main river. Lower diversity in lakes compared with the connected rivers has been attributed to longer WRT in lakes (Crump et al., 2012; Ruiz-González et al., 2015; Niño-García et al., 2016). This is especially pronounced in small streams and rivers, where WRT is too short to allow for local sorting of the bacterial community (Crump et al., 2012; Ruiz-González et al., 2015; Niño-García et al., 2016). Downstream of the input sites, the alpha diversity rises again and by the river mouth reaches a level similar to the first river site (R1) upstream of the input sites (R2 and R3). This shows that the volume of water from the lake and glacier outlets does not dilute the downstream river community. Importantly, it also suggests that the less diverse communities from the lake and the glacier stream are concealed downstream of the input sites by the higher diversity of the main river.

Shannon indices in the estuary samples were higher than those previously described for Arctic estuaries (Galand et al., 2006, 2008). Previous studies of large Arctic rivers show that bacterial diversity and abundance decrease from rivers to estuaries probably due to upstream input from terrestrial sources (Meon and Amon, 2004; Galand et al., 2006, 2008). Our results from a small river support this conclusion by showing a slight decrease in diversity from the main river sites (R1 and R5) to the estuary (Figure 2). The diversity in the estuary sites closer to the river and in the shallow samples could be expected to be higher than more distant and deep estuary samples due to a higher concentration of the river bacterial community, which is not evident from our results (Figures 2, 6). This indicates that although the community structure from the river to the estuary aligns with previous results by showing a directional decrease in diversity, this directional structure cannot be detected further down the network, in the estuary transects. The lack of pattern in diversity and richness attributed to different depths of the estuary samples or the distance to the river mouth could partly be attributed to an insufficient resolution in sample depth. The low resolution might not allow for the detection of a clear plume and different depth zones. While distinct bacterial communities have been found to be associated with the plume and different oceanic zones in an estuary, these results are from sampled oceanic zones several kilometers farther into the ocean than our samples (Fortunato et al., 2012). Estuary samples were previously discussed as harboring a mix of bacterial communities from the river and the coastal ocean with no distinct autochthonous estuary-community (Fortunato et al., 2012), which resembles our results. The high variability in diversity found among the estuary samples suggests a highly heterogeneous community, also expected in such a region where waters of very different chemistry and origin meet. The large variability in diversity and richness seen at the different estuarine sites (Figure 2) could be explained by the sampling of the different water masses, since sequences from the three different transects were pooled together for each estuarine site and depth. The variability of diversity as well as environmental variables seem to lessen farther into the estuary, which could be expected as homogeneity increases as a greater fraction of the estuary is made up of marine waters (Figure 2 and Supplementary Table 1).

Samples from the first site of the river (R1), upstream of the lake and glacier stream outlets, clustered with samples from the bottom site of the river (R5) while the river site by the proglacial lake outlet (R2) clustered with the river site at the glacier stream outlet (R3) as shown by the NMDS plot (Figure 3). R1 and R5 also shared a high number of indicator OTUs (Figures 4A,E) as well as similar diversity and richness (Figure 2). In contrast, R2 and R3 had lower diversity than the other river sites. The NMDS plots, indicator OTU analysis and alpha diversity results imply that waters sourced from the lake and the glacier stream carry different bacterial communities than that of the main river. The larger volume of the main river community then probably masks the lake and glacier stream communities, thus resulting in the close similarity between the sites R1 and R5. The isolation of the R4 samples from other river samples in the NMDS plots may be explained by the difference in sampling at this site, which was closer to the river bank compared to the other river samples. Another explanation might be the imperfect mixing of water from the upstream lake and glacier outlets with that of the main river at this site. The latter seems to be the best explanation since the comparably low number of indicator OTUs found at R4 suggests that this site contains a mixture of the upstream communities rather than a distinct community from the sampling site (Table 3). This agrees with results from a study on an Arctic tundra catchment, showing that streams leaving lakes have decreasing similarity to the lake microbial community as a function of distance (Crump et al., 2007).

The dispersal of estuary samples on the NMDS plots was in accordance with the diversity measures, which were similar within individual sampling sites independent of sample depth (Figures 2, 3). Samples from the sites closest to the river mouth (E100 and E300) clustered more closely with river samples than the samples farthest from the river mouth (E700 and E1100), consistent with a gradual mixing of the river community with a marine community within the estuary environment. Remarkably, R4 river samples clustered more closely with estuary samples than with the other river samples. Indicator OTUs from R4 are present throughout the estuary transects (Figure 6) and the taxonomic composition of samples from R4 has the largest resemblance to that of estuary site E1100 (Figure 5), which might explain the NMDS results (Figure 3). Furthermore, the bacterial community at this site seems to be a mixture of the different river communities as suggested by the indicator OTU results (Table 3). Therefore, the clustering of R4 samples with estuary samples might also reflect the resemblance to the estuary, in which the river communities are also mixed (Figures 4F–I, 6).

Samples from E100 and E300 were more widely dispersed across the NMDS plots than the samples from farther into the estuary, indicating greater heterogeneity of the bacterial communities. This is not unexpected from a region of mixing of largely different water bodies both in terms of physical and chemical variables as well as origin. The NMDS plot did not indicate that the bacterial communities were stratified according to sample depth. The low resolution of samples through the water column might be part of the explanation. However, the results might also indicate a high degree of mixing through the water column of the Red River estuary at the time of sampling. As previously discussed, this may also be explained by the proximity to the river of the estuary samples in this study compared to other studies, where bacterial communities in the estuary were shown to be stratified according to depth (Fortunato et al., 2012).

The BIOENV analysis did not show strong correlations between environmental variables and community composition but it did highlight turbidity as a community-shaping factor in the river (Table 2). The bacterial community in a freshwater network fed by glaciers has recently been shown to be structured along the turbidity gradient (Peter and Sommaruga, 2016). The BIOENV results support the idea that the bacterial community in the Red River freshwater network is partly sourced from the glacier. The lower turbidity at R2 (1.2 NTU) compared to an average of 16.3 NTU (SD = 1.2) at the other river sites is noteworthy since proglacial lakes are known to have high turbidity (Peter and Sommaruga, 2016). While turbidity of the proglacial lake outlet (R2) is higher than that shown for a non-glacier fed lake in the Austrian Central Alps, it is remarkably low compared to other glacier-fed lakes (Peter and Sommaruga, 2016). This might indicate that the proglacial lake is losing hydrological connectivity to the glacier (Peter and Sommaruga, 2016). It should be taken into consideration that the samples are not taken from the actual lake but several 100m downstream (Figure 1).

The Red River is small in size compared to large rivers previously described such as the Mackenzie River (Galand et al., 2008; Garneau et al., 2009) and the Columbia River (Fortunato et al., 2013). For comparison, the average water flow from August to November in the Columbia River was 2988 m³s⁻¹ (Fortunato et al., 2013), while the river flow in the Red River around sampling time was 5.7 m³s⁻¹. An estimated time from top sampling site R1 to the river mouth at R5 is 40 min for the moving water body where the samples are taken. It has been shown that at sites with shorter WRT than 10 days the bacterial community composition was predominantly structured by hydrology (Niño-García et al., 2016). Accordingly, we hypothesized that hydrology would be dominant in shaping the bacterial community in the relatively small Red River with short WRT. Consequently, we did not expect strong correlations between community composition and environmental variables in the river samples. Our study represents a single catchment with short WRT and the results of the BIOENV analysis agrees with previous results by showing weak correlations between the bacterial community and environmental variables (Niño-García et al., 2016).

Salinity has previously been highlighted as a community-shaping factor in estuaries and rivers. For example, the abundance of Alphaproteobacteria, Betaproteobacteria, and Actinobacteria correlated strongly with salinity in the Delaware estuary where a strong negative correlation between Betaproteobacteria and Actinobacteria was shown together with a positive correlation between salinity and Alphaproteobacteria (Kirchman et al., 2005). Salinity, together with temperature, explained 45% of the variation in the community composition in a study of the Mackenzie Shelf (Garneau et al., 2009). Salinity was not identified as a significant factor in the BIOENV analysis. The lack of correlation to salinity in our study is also evident from the NMDS analysis (Figure 3), where R4 river samples cluster with estuary samples despite the large difference in salinity between these environments (Table 1). These results suggest that there may be environmental or hydrological factors other than salinity that explain the observed patterns in taxonomic composition in the study site.

BIOENV analysis of the estuary community showed no significant correlations with environmental variables (Table 2). Also no correlation was found between the bacterial community and spatial variables including distance from the river mouth and depth. The results from the BIOENV analyses indicate that the bacterial community in the estuary is not dispersed according to environmental variables or stratified according to distinct water bodies of riverine or oceanic origin, supporting the results from the diversity assessments (Figure 2) as well as the NMDS plots (Figure 3). Our samples represent a very small fraction of the total estuary; a higher resolution of samples in the estuary might result in more conclusive results.

Indicator OTUs identified at the input sites from the lake and glacier stream were found in low numbers at the other river sites and, therefore, seem to be specific to their respective sources (Figures 4B,C). Notably, of the 678 indicator OTUs from the lake outlet (R2, Table 3), 570 OTUs were not found in the upstream river site (R1) and seem to originate from the proglacial lake. The lake outlet site (R2) had a particularly high number of indicator OTUs and percentage of top indicator OTUs (Table 3). Water bodies with longer WRT have been shown to harbor a less diverse and more differentiated community explained by local sorting of the microbial community (Niño-García et al., 2016). Results from the indicator OTU analysis and diversity of the lake outlet site (R2) show that the lake with a longer WRT has a less diverse and more specialized community compared with the river. A number of taxa known to be psychrophilic, such as Moritella (Urakawa et al., 1998), Polaribacter (Gosink et al., 1998), Oleispira (Yakimov et al., 2003), Crocinitomix (Bowman et al., 2003), and Psychromonas (Mountfort et al., 1998), were found among the best matches for the indicator OTUs from the lake outlet, unlike at the other river sites.

The short WRT in the Red River network should accordingly result in a low degree of differentiation, which is confirmed by the low number of top indicator OTUs (i.e., OTUs unique for a particular site), which average < 2% in present study (Table 3). Another study of freshwater networks highlighted an average of 11% unique OTUs between different ecosystems as representing a low number (Ruiz-González et al., 2015). These results differ from our study in that they considered many different lakes while our results are obtained from one lake only (Ruiz-González et al., 2015).

Samples from the lake and glacier outlets (R2 and R3) as well as the sample site just after the glacier outlet (R4) had very few Acidobacterial classes compared to the top and bottom site of the river (Figure 5). Acidobacterial classes were shown to be most common in soil compared to the adjacent freshwater network and Acidobacteria in rivers seem to be sourced from the surrounding terrestrial environment (Ruiz-González et al., 2015). The more differentiated lake community thus seems to harbor a lower fraction of organisms from the surrounding soil community compared to the main river. This might partly be explained by the local sorting of the bacterial community in the lake. It could potentially also to be explained by the presence of taxa with different origin than the main river community, as indicated by the large number of indicator OTUs, which are not present in the upstream site. The lower turbidity at the lake outlet site (R2) might also indicate that there is less input of soil to the lake than to the main river, which causes less mass dispersal effect from the surrounding terrestrial environment. The lower turbidity might, however, also be explained by less suspended particles in the lake because of increased sedimentation due to the longer WRT.

A great number of previous studies of riverine microbial communities have suggested and shown that the river communities are influenced by input of microorganisms from surrounding soil environments (Crump and Baross, 2000; Galand et al., 2006, 2008; Crump et al., 2012, 2007; Ruiz-González et al., 2015; Niño-García et al., 2016). Our results support this by showing that potentially soil-related taxa make up a significantly large fraction of the bacterial community making them part of the indicator OTUs of the main riverine bacterial community (Figure 5). Notably, the results also indicate that, along the river, distinct communities may not have the same degree of influence from the terrestrial surroundings. August is a month of high precipitation and increased erosion around the Red River, which would result in a relatively high influence from the surrounding soil community. The influence from soil may be less pronounced in other months as water flow and erosion levels change.

The dominance of soil microbes in freshwater networks has been established in several recent studies, highlighting soil as the origin of the network metacommunities (Ruiz-González et al., 2015; Niño-García et al., 2016). A gradual differentiation of a stream from an upstream lake as a function of distance has been attributed to the origin of the freshwater communities from a terrestrial metacommunity (Crump et al., 2012, 2007). Our results suggest that glaciers may also supply part of the metacommunity resulting in a different structuring pattern of the network. In our case the structuring pattern was not uni-directional throughout the network but rather showed local changes as different bacterial communities were added to the river. This is illustrated in the diversity results (Figure 2) as well as the NMDS plots of the community composition (Figure 3). These results together with the indicator OTU analysis highlight the importance of additional sources of the metacommunity such as glaciers.

The indicator OTUs from the lake-sourced water (R2) can be found in the second highest abundance in the estuary site 700 m into the estuary (E700) (Figures 4B, 6). The taxonomic composition of E700 differs from the other estuary sites and this site contains a high fraction of Alphaproteobacteria (71%) and a relatively small fraction of Gammaproteobacteria (3%) (Figure 5). Of the 678 R2 indicator OTUs 161 are found at the estuary site 700 m into the estuary, where they make up 20% of the sequences at E700 with a higher fraction in the deep samples compared to the surface samples (Figures 4B, 6). This resembles results from the Columbia River, where the estuary samples were comprised of just over 20% riverine community (Fortunato et al., 2012). The distribution of R2 indicator OTUs suggests that although the organisms from the lake do not form a large enough fraction of the community to be notable along the downstream river, they are transported into the estuary where they form a larger fraction of the community. Our results align with the “landscape reservoir” concept proposed for the Toolik lake, Alaska, where rare organisms from the upslope landscape influence downslope bacterial diversity and become dominant in environments with favorable conditions (Crump et al., 2012).

The taxonomy of the nine indicator OTUs from R2 found in high numbers (>100 sequences) in the E700 samples were mostly related to organisms isolated from oceanic environments such as Marinomonas (Van Landschoot and De Ley, 1983), Oleispira (Yakimov et al., 2003), Pseudoalteromonas (Bowman, 2007), Polaribacter (Gosink et al., 1998), and Sulfitobacter (Sorokin, 1995). Of related non-marine organisms were Glaciecola, which was first described as a Gammaproteobacterium isolated from Antarctic sea ice (Shivaji and Reddy, 2014) and Rhodobacteraceae known from aquatic environments (Pujalte et al., 2014). The fact that indicator OTUs from the proglacial lake outlet to the river are similar to known marine organisms suggests that these organisms are commonly found in marine environments and that they are not originally known from terrestrial environments. Since it is unlikely that organisms are transported from the estuary to the proglacial lake over 2 km upstream, these organisms in the estuary more likely originate from the upstream freshwater network. Possibly, they become such common organisms in the estuarine and marine environments, that these are the environments from which they have become known. It is well established that bacterial communities found in freshwater networks can be traced from upstream positions in the network (Crump et al., 2012; Ruiz-González et al., 2015; Niño-García et al., 2016). We show that in the Red River estuary the river community can be found in the estuary with an overall decreasing fraction from the river mouth toward the ocean (Figure 6). Interestingly, communities that are not notable throughout the river are transported to the estuary where they seem to become an equally large fraction of the estuary as the main river community (Figures 4B, 6). Distinct communities from the river seem to influence the estuary to different extend, so that communities from certain parts of the river make up notably larger fractions of the estuary at some sites (Figures 4, 6).

“Seed bank” is a term proposed for the fraction of dormant organisms that may be resuscitated when met with different environmental conditions through e.g., dispersal to other ecosystems (Lennon and Jones, 2011). The concept of seed banks was recently extended to freshwater networks where organisms originating in a soil community were proposed as the seed bank for boreal freshwater networks (Ruiz-González et al., 2015). For freshwater networks it was discussed that shallower sequencing depth might lead to the erroneous conclusion that freshwater communities do not derive from a shared pool of terrestrial microbes (Ruiz-González et al., 2015). This could lead to an incomplete understanding of the mechanisms of assembly and the actual linkages and dispersal of microbes between connected ecosystems (Ruiz-González et al., 2015). We show that sampling resolution not only in terms of sequencing depth but also resolution along the network may result in overlooking distinct microbial communities and how these are distributed and linked to the downstream estuary. Our results indicate that not only the terrestrial surroundings but also upstream glaciers may act as seed banks for freshwater networks. While the uni-directional structure in freshwater networks might be a consequence of the numerical dominance of terrestrial OTUs as shown previously (Crump et al., 2012; Ruiz-González et al., 2015) our results suggest that this does not necessarily imply that the bacterial community in a freshwater network has a common origin from microbes from soil. A higher resolution along the river might reveal distinct bacterial communities of different origin and with different composition, which are introduced downstream in the network. These distinct communities, which might be concealed by the numerically dominant terrestrial community along the river, are able to act as seed banks for downstream environments. The different composition of inputs along the river affects the structure of the community, which is not necessarily uni-directional for all freshwater networks as shown in the present study.

Therefore, sampling with the right resolution, both in terms of sequencing depth and the distance between samples along the network, is crucial for understanding the source of microbial communities found in the estuary. This is especially true at times with high precipitation and erosion. Our study shows that with the right resolution, microbial communities can be valuable in understanding transport pathways of meltwater and matter from source to oceans in that they can serve as both tracers as well as indicators of origin in their adaptation to the environment.

The indicator OTUs for each sample site in the estuary were found only in low numbers at the other sites both in the river and the estuary. This pattern is in contrast with the results from the Columbia River, USA, where indicator OTUs from the upper water-column of the estuary (< 56 m depth) showed generalist taxa qualities by having high relative abundance and occurrence in a high number of samples outside their indicator environment when compared to the river environment (Fortunato et al., 2013). Ocean-specific taxa are not expected to be found upstream in the river, which is also apparent from the distribution of estuary indicator OTUs (Figures 4F–I). The low number of indicator OTUs and lack of top-indicators in the estuary show that the different sites in the estuary do not hold distinct communities. These results are in accordance with the results of the BIOENV analysis and NMDS plots, which suggest that the bacterial communities are not dispersed according to environmental or spatial variables, as well as the highly variable diversity measures in the estuary. As previously discussed, samples in the present study are sampled relatively close to the river mouth and a more distinct stratification of the bacterial communities might become visible farther into the estuary.

We expected the estuary sites to contain a mixture of the communities found in the river and the ocean, with more environmental variability closer to the river mouth due to the mixing of river- and sea-water. This is supported by our data, which show greater variance in environmental data closest to the river mouth (Supplementary Table 1) as well as river indicator OTUs from more of the river sites closer to the river mouth at E100 and E300 (Figure 6). Mixing of river and ocean water may result in an allochthonously dominated community shaped by hydrology rather than by environmental selection, also indicated by the NMDS plot (Figure 3) and BIOENV analysis (Table 2). This is supported by the fact that the outermost sample site of the estuary (E1100) had a higher number of indicator OTUs (Table 3) as variability is expected to decrease with increasing distance from the river mouth and the most distant estuary site is expected to contain a higher number of ocean indicator OTUs. This was supported by the taxonomy of the indicator OTUs that were all similar to marine-related taxa at E1100 (details not shown). This site also had a lower fraction of river indicator OTUs compared to the estuary sites closer to the river mouth (Figures 4, 6).

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.