Terrestrialization of sediment bacterial assemblages when temporary rivers run dry

Abstract

Bacterial communities in river sediments are shaped by a trade-off between dispersal from upstream or nearby land and selection by the local environmental conditions. In temporary rivers (i.e. those characterized by long drying periods and subsequent rewetting) seasonal hydrological dynamics shape bacterial communities by connecting or disconnecting different river habitats. In this study, we tracked and compared the temporal and spatial changes in the composition of bacterial communities in streambed sediments and floodplain habitats across both permanent and intermittent river segments. Our findings revealed that environmental selection played a key role in assembling bacterial communities in both segments. We argue that distinct environmental features act as filters at the local scale, favoring specific bacterial taxa in isolated pools and promoting some typically terrestrial taxa in dry areas. Considering the prospective extension of drying intervals due to climate change, our results suggest an emerging trend wherein bacterial assemblages in temporary streams progressively incorporate microorganisms of terrestrial origin, well-adapted to tolerate desiccation phases. This phenomenon may constitute an integral facet of the broader adaptive dynamics of temporary river ecosystems in response to the impacts of climate change.

Introduction

Bacterial communities spatially differ in composition according to the environment type (Lozupone and Knight 2007, Barberán and Casamayor 2010) and over time in response to environmental fluctuations (Gilbert et al. 2012, Jones et al. 2012, Shade et al. 2012). River networks are characterized by diverse habitats and hydrological conditions, associated both with the position of a river segment within the river network and to the current hydrological phase. The combination of these two factors drives bacterial community composition over space and time (Freimann et al. 2015, Malazarte et al. 2022).

In river ecosystems, water flow is the principal vehicle for species dispersal, allowing the movement of both aquatic and terrestrial bacteria when river flow (Stadler and del Giorgio 2022). According to the metacommunity theory, differences in bacterial community composition across different areas of the river network are caused by two main mechanisms: (i) dispersal of microorganisms driven by the physical forces governing the ecosystem, and (ii) selection of species based on their physiological requirements and the local conditions (Leibold et al. 2004, Martiny et al. 2006, Hanson et al. 2012, Lindström and Langenheder 2012, Székely and Langenheder 2013). It is well established that the relative importance of these mechanisms varies both spatially and temporally in bacterial communities thriving in the river water column. Previous studies demonstrated that a massive transport of bacteria from neighboring soils (referred to as mass effects) is the main mechanism driving species composition in the upper parts of the river, whereas, in the downstream sections, the main mechanism is environmental selection (species sorting) (Crump et al. 2012, Ruiz-González et al. 2015). The relevance of these mechanisms may also vary temporally between seasons. For instance, in boreal river networks, mass effects dominate in spring and species sorting in summer due to seasonal fluctuations in water flow and temperature (Stadler and del Giorgio 2022).

Additionally, temporal variations such as hydrological fluctuations in rivers enhance the interaction between terrestrial and aquatic microbial communities, a topic primarily studied in boreal and arctic streams (Crump et al. 2012, Ruiz-González et al. 2015) but less considered in temporary river networks where the hydrological dynamics are crucial. Temporary rivers cease to flow at certain points in space and time (Acuña et al. 2014) causing those previously flowing habitats to become dry for extended periods (Messager et al. 2021). Furthermore, the hydrological sequences of flow contraction, desiccation, and expansion in temporary streams vary throughout the river network (Sabater and Tockner 2009). Under these dynamics, some river segments remain permanent, while others may experience short-term water flow interruption, or become ephemeral within a single river network (Llanos-Paez et al. 2023). During the process of flow contraction, water may become restrained to pools (Fazi et al. 2013), intensifying interactions with neighboring terrestrial habitats (Niño-García et al. 2016).

Within this highly variable context, local environmental conditions are expected to impose significant selective pressures on microbial communities (Huber et al. 2020). However, it remains unclear whether these selective pressures vary between river segments with contrasting flow permanence.

Bacterial communities within river ecosystems play a crucial role in cycling nutrients and organic matter, serving as the basis of the food web (Battin et al. 2016). Sediment microbial communities drive the majority of the metabolic activity in river ecosystems (Romaní and Sabater 2001, Gibbons et al. 2014). Shifts in sediment bacterial composition directly influence biogeochemical patterns and ecosystem functions, thereby impacting overall ecosystem health (Shade et al. 2012, Freixa et al. 2016). Studying the dynamics of bacterial community composition is essential for predicting ecosystem resistance and resilience to human-induced disturbances.

In our study, we aim to study the mechanisms accounting for the composition of sediment bacterial assemblages in temporary river networks. For that, we selected two river segments with contrasting hydrological conditions (permanent versus intermittent flow), characterizing the composition of bacterial communities in sediments, including aquatic and nearby terrestrial habitats. Specifically, we aimed to contrast the following hypotheses: (i) bacterial composition across different habitats (aquatic and terrestrial) becomes more uniform during periods of hydrological connection, and more heterogeneous when these habitats are hydrologically disconnected; and (ii) environmental selection predominantly shapes bacterial community composition in sediments during periods of river flow intermittency, while dispersal is the major force during periods of water reconnection.

We assert that while the dispersal and mass transport of bacteria is constrained by the spatial and temporal fluctuations in the water flow, the composition of local communities is shaped by their specific local environmental conditions. Given that temporary streams contribute 51%–60% of global river networks (Messager et al. 2021), our results may illustrate potential scenarios that could arise in river ecosystems submitted to a pronounced reduction in the water flow. Ongoing climate change is leading to longer and more frequent drought episodes in river systems (Arias et al. 2021), resulting in significant impacts on aquatic biodiversity (Sabater et al. 2022). Therefore, thus our findings may contribute to revealing and forecasting the implications associated with changes in precipitation patterns. Increased frequency and intensity of extreme events, such as floods could also alter river bacterial community dynamics (Mishra et al. 2021). Overall, rivers are subject to strong effects on their hydrology, exacerbated by the effects of climate change and human intervention, which could affect the distribution and composition of riverine bacterial communities.

Materials and methods

Study area and sampling strategy

This study was performed in the Algars River, a tributary to the Ebro (NE Iberian Peninsula, 41° 12′ N, 0° 15′ E). The 405-km² basin area is largely forested and encompasses crop fields in the alluvial area of the lower part. The area has a Mediterranean climate, which includes a long dry period (May–September) and episodic floods during autumn or winter. The mean annual discharge in the lower river section ranges from 0.1 to 0.2 m³/s (www.chebro.es). We selected two river segments that differed in their geomorphological settings and hydrological conditions. The upstream segment was in a constrained valley section, and it has a 3–8-m wide channel with permanent flow (P-segment: permanent). The downstream segment was in the lower part of the river, where the channel is wider (20–50 m wide) and crosses a large floodplain where water flows partially or stops flowing (I-segment; intermittent). Between April 2019 and February 2020, the P-segment had a continuous surface water flow while the I-segment remained dry for 125 days (Fig. S1).

The two segments under study differed not only on their water regime but also on the complexity of their river channels. The development of floodplain areas, presence of isolated pools, and the river channel width and depth were more complex in the I-segment than in the P-segment. For comparison, we collected sediments from four different habitats studied in the two segments. These were two aquatic habitats (having river water surface), defined as Riffles (R) and Pools (P), and two predominantly terrestrial habitats within the floodplain, defined as Bank sediment (BS, littoral zone) and Floodplain sediment (FS). We sampled these four habitats in the four different hydrological phases, which occurred throughout the hydrological cycle (Fig. 1). These phases were named Contraction (C, when surface water was at its shrinking phase), Fragmentation (F, when the water flow decreased and the river habitats became disconnected and resulted in the formation of pools), Non-flow [NF, when the river is dry at the I-segment (or substantially decreased at the P segment)], and Expansion (E, when a high flood event restarted water flow along the river network) (Fig. S1). We randomly collected three surface sediment cores (first 5 cm depth) in each of the four habitat types over the four hydrological phases at the two segments, resulting in 95 sediment samples.

Water and sediment characteristics

We measured water conductivity, pH, temperature, and oxygen content in Riffles and Pools using portable probes (WTW, Weilheim, Germany). We collected three replicates of filtered river water [0.2 µm pore size, nylon filters, (Merk Millipore Ltd)] for nutrient analyses. Nitrate and ammonium (NO₃⁻ and NH4⁺) were analysed by ionic chromatography (Dionex, ICS 5000) and phosphate (PO₄³) was analysed spectrophotometrically by the ascorbate-reduced molybdenum blue method (Murphy and Riley 1962). Dissolved organic carbon (DOC) was quantified using a total organic carbon analyser (Shimadzu TOC-V CSH, coupled to a TNM41 module).

Water content in the sediments was measured as the difference between the wet and dry weight of sediment samples after 72 h at 65°C, expressed as the percentage of the wet weight. Organic matter in the sediment was measured by combusting the dry sediments at 450°C for 4 h, and the results expressed as ash-free dry weight (AFDW). Changes in the sediment water-extractable organic matter (WEOM) were assessed following the method described by Muñoz et al. (2018) representing the most mobile and bioreactive fraction of the organic matter. Briefly, 30 g of sediment was extracted by shaking in 90 ml of Milli-Q water (1:3 sediment:water ratio) at 100 r/m at 4°C for 48 h in the dark. The resulting water was filtered through 0.7 µm pore diameter GF/F glass microfiber filters previously combusted (Whatman TM, Life Science) to determine the w-DOC and w-dissolved organic nitrogen (DON) of WEOM. w-DOC and w-DON were quantified using a total organic carbon and nitrogen analyser (Shimadzu TOC-V CSH 230 and Shimadzu TNM-1).

16S rRNA gene amplicon sequencing and data processing

The bacterial community composition in the sediments was characterized through amplicon sequencing of the 16S rRNA gene. First, we extracted DNA from 1 g of each sediment sample using the FastDNA Spin Kit for soil (MP Biomedicals), according to the manufacturer’s instructions. DNA concentrations and quality were measured using Qubit 1.0 fluorometer (Life Technologies) and Nanodrop2000 spectrophotometer (ThermoFisher Scientific), respectively. The resulting DNA extracts were subjected to a high-throughput sequencing of the V4-specific region of the 16S rRNA gene using specific primers with the Illumina MiSeq system (PEX 250) at the Sequencing and Genotyping Unit of the Genomic Facility/SGIker of the University of Basque Country (Leioa, Spain). Briefly, extracts were used as a template in Polymerase Chain Reaction (PCR) reactions using the primer pair 515F/806R (Parada et al. 2016). complemented with Illumina adapters and sample-specific barcodes using the MiSeq reagent kit v2 500 cycles. The thermal PCR conditions consisted of an initial denaturation cycle at 94°C for 3 min, then 94°C for 45 s, 50°C for 60 s, 72°C for 90 s for 35 cycles followed by a final extension at 72°C for 10 min. Amplification was verified on a 1.5% agarose gel. Amplicon libraries were prepared according to the Miseq manufacturer’s protocol (Caporaso et al. 2012).

Demultiplexing, quality filtering, clustering into amplicon-sequence variants (ASV) and construction of the ASV table were performed in QIIME2 (Bolyen et al. 2019). DADA2 algorithm with default parameters was used to denoise and merge paired reads, filter chimaeras, and dereplicate sequences. The feature classifier script implemented in QIIME2 was employed for taxonomic assignment using the SILVA reference database v.138 (Quast et al. 2013) ASV table was filtered to remove ASV appearing in less than 10 samples and having less than 100 counts across all samples. For community analysis, the filtered ASV table was rarefied by a subset of 10 000 sequences per sample to minimize bias due to different sequencing depths between samples. A total of 2990 ASVs were identified in the dataset.

Data analysis

We used a combination of multivariate and network analyses to assess (i) the alpha diversity of bacterial communities in the I-segment and P-segment; (ii) the spatiotemporal changes of the bacterial community composition between and within segments, and (iii) disentangle the community assembly processes (dispersal versus selection) underlying the structure and dynamics of river bacterial meta-communities using a null modeling approach. A detailed description of the data analyses related to our aims and expectations is compiled in Table S1.

Alpha diversity

Bacterial diversity was assessed using the Hill numbers (q = 0, q = 1, and q = 2) (Chao et al. 2014) to obtain a unified framework of the taxonomic and phylogenetic facets of biodiversity. The Hill numbers are a parametric family of diversity indices differing only by the parameter q, which determines the sensitivity to the relative abundance of the taxa in the community (Jost 2007) Studies typically use the first three numbers [q = (0, 1, 2)] as they allow to assess the importance of rare and dominant species in a continuum. Moreover, these three numbers are also equivalent to some of the most used diversity indices. When q = 0, species abundances are not considered, and therefore, the estimation represents the species richness in the community. When q = 1, the species are weighed in proportion to their relative abundances, analogous to the Shannon–Wiener index of diversity. When q = 2, estimates are equivalent to the Simpson index, where the most abundant species determine the diversity values. Taxonomic Hill numbers were computed in R software (R Development Core Team 2011) using the package hilldiv (Alberdi and Gilbert 2019) and considering ASVs as entities equivalent to species.

Spatial and temporal changes in bacterial community composition

We explored the differences in bacterial community composition between segments, periods, and habitats at the level of ASVs, computing a permutational multivariate analysis of variance (PERMANOVA) and using the Bray–Curtis dissimilarity matrix based on 9999 permutations with the adonis function in vegan (Oksanen et al. 2020). When significant differences between segments, habitats or periods were found, Bonferroni post hoc tests were employed to make pairwise P-value comparisons (using pairwise.adonis2 function in R). A principal coordinate analysis (PCoA) was conducted using Bray–Curtis dissimilarity distance to explore differences in bacterial community composition across periods and habitat using the function ordinate in the phyloseq package (McMurdie and Holmes 2013). Bar plots were generated in phyloseq by grouping the data based on taxonomic classification and transforming it to relative abundance using tax_glom and prune_taxa functions. We also compared the presence of each taxonomical unit by plotting an upset chart using the package UpSetR in R (Conway et al. 2017) based on a binary matrix (presence or absence of each ASV) combining each hydrological period and habitat.

In this case, combinations with less than 10 shared species were excluded from the analysis. Finally, the relationship between the bacterial community (ASV table) and the environmental sediment variables was analysed using a canonical correlation analysis (CCA) using the cca followed by anova function in vegan. Finally, to study the differences between sampling periods at each site, based on water quality variables and sediment characteristics, a principal correspondence analysis was performed using prcomp function and results were plotted using the fviz_pca_biplot function in R software.

Processes shaping metacommunity structure and dynamics

We followed a null modeling approach to test whether changes in water flow determined the prevailing mechanism (i.e. dispersion versus environmental selection) structuring bacterial communities. Specifically, we estimated the dissimilarity in community composition (beta diversity) between sampling periods using the Raup–Crick index (Raup and Crick 1979) and then tested how much it deviated from randomness. We followed the approach proposed by Chase et al. (2011) in which the Raup–Crick index is used as a probability metric (ranging from –1 to 1) that indicates whether local communities are (i) more similar (close to 1), (ii) as much dissimilar (close to 0), or (iii) less dissimilar (approaching –1), than expected by random (Stegen et al. 2012). Values close to 1 indicate that bacterial communities are more similar than expected by chance and are considered to be a part of the same regional pool. In this case, differences between samples are explained by environmental selection. Values close to –1 indicate that the communities are less dissimilar than expected by chance because of a higher dispersal between local communities. Finally, values close to 0 indicate that the communities are as dissimilar as it might be expected by chance, which suggest that there is not any governing factor but stochasticity to explain the changes. We implemented a null model that controls richness between samples, since differences in species richness among habitats or sampling dates may introduce some bias in estimating beta diversity (Chase et al. 2011). We used the R script developed by Modin et al. (2020), which randomly assembled samples with predefined numbers of ASVs from the regional pool. We included all the ASVs detected in the two segments as configuring the regional pool. As a randomization scheme, we based this on the relative frequency of ASVs. We performed 9999 randomizations and built a null distribution representing the dissimilarity between two random samples. The null and the observed dissimilarities were then compared to determine if the observed dissimilarity was not driven by randomness. If the observed dissimilarity was higher or lower than the null expectation, we assume that likely there were deterministic factors favoring different or similar taxa in the two habitats (Chase et al. 2011). The script with the pipeline to implement the null model is included as Supplementary material.

We built bipartite networks to disentangle which habitats and periods were more (or less) similar in terms of the composition of bacterial species (ASVs). Specifically, we employed the bipartite package in R to build two segment-level networks that described the similarity of communities at the different hydrological phases for bacterial composition in each habitat (Dormann 2023). Interactions between taxa were derived from pairwise comparisons between networks, focusing on the most abundant ASVs with relative abundance higher than 0.1% (representing 10% of the total number of reads) (Stadler and del Giorgio 2022). Networks were built using two matrices (one per segment) with columns representing the name of each habitat and rows holding the name of the most abundant species. Entries in each matrix represent the sum of interactions (links) among habitats representing each bacterial taxa, providing clues about the bacterial origin. Then, we calculated the beta diversity of the interactions between shared species (βOS) to assess the dissimilarities of the interactions established between species common to both networks and the interaction beta diversity (βWN) to indicate dissimilarity between the two networks (Poisot et al. 2012). These calculations were performed using the function betalinkr_multi in R.

Finally, counts of the most abundant ASV—relative abundance higher than 0.1—were used as the input file for source tracking using SourceTracker2, a Bayesian approach program to estimate the proportion of exogenous sequences in a given (sink) community that come from possible (source) environments (Knights et al. 2011). Source tracker compares the community profiles in the “source” to those of the sink using Bayesian methods to identify the extent of the contribution of each potential source to the sink. A key feature of SourceTracker2 is to create an unknown source environment, which encompasses those sequences that do not have a high probability of being derived from the known and specified source environments. We identified the riffles in each site as a sink of bacteria for the analysis, while the sources of riverine bacteria were identified as samples collected from three different habitats: the floodplain, pools, and bank sediment. SourceTracker2 analysis was performed separately for each river segment (P and I) and conducted using default settings following Knights et al. (2011).

Results

River water quality and sediment characteristics

The Algars River experienced marked changes in the water quality in the different hydrological phases, and these changes were more accentuated at the I- than in the P-segment (Table S2 and Fig. S2). Water temperature was comparatively higher in the I-segment, attributed to its shallower water column and increased exposure to sunlight. Additionally, temperature variation between hydrological periods were more pronounced in this segment, fluctuating between 8.4°C and 20.2°C. Dissolved oxygen had similar average values in the two segments, but differed between habitats, being lower in pools (ranging between 4.1 and 12.29 mg/l) than in riffles (7.7–11.5 mg/l; Table S2). Water conductivity was lower in the P-segment (524–629 µS/cm) than in the I-segment (738–2080 µS/cm). Nitrate was 0.98 ± 0.99 mg/l at the P-segment and 33.6 ± 8.87 mg/l at the I-segment (mean ± SD), being this high value associated with the agricultural activities in the lower part of the river. DOC was 1.47 ± 0.26 mg/l at the P-segment and slightly higher at the I-segment (2.0 ± 0.35 mg/l; Table S2).

Higher organic matter in sediments was observed at the P-segment than I-segment (Table S3 and Fig. S2). Organic matter and w-DOC substantially accumulated at the P-segment during the fragmentation and non-flow phases. The highest concentration of w-DOC (18.1 ± 10.9 mg C/kg fresh sediment) and w-DON (4.01± 2.07 mg N/kg fresh sediment) were observed at the I-segment during the fragmentation period. Water content in the sediments of this segment declined during the non-flow period (Table S3). Regarding habitats, the lowest water content was measured in the floodplain sediments habitat, and the highest levels of w-DOC and w-DON were measured in the bank sediments (Table S3).

Sediment bacterial composition across sites

The bacterial community composition was significantly different between the two river segments (PERMANOVA test, Segment factor, F: 11.15, P-value = .001). In the P-segment, the bacterial composition was characterized by a greater proportion of Proteobacteria whereas in the I-segment, Actinobacteria and Firmicutes (bacilli) were comparatively more abundant (Fig. 2). However, the taxonomic alpha diversity did not differ significantly (Fig. S3; pairwise t-test, P-value > .05, Bonferroni corrected). The average number of ASVs (Richness, Hill number, q = 0) was 613 ± 118 (mean ± SD), from which 382 ± 66 (q = 1, diversity) and 237 ± 38 were the most abundant (Hill number, q = 2) for the I-segment, and 580 ± 147 (q = 0), 364 ± 105 (q = 1), and 229 ± 77 (q = 2) were for the P-segment (Fig. S3). The I-segment had slightly higher q = 0 and q = 1 than the P-segment but no significant differences were measured.

Sediment bacterial composition: comparison across hydrological periods and habitats

Changes in the bacterial community composition between hydrological periods were significant in the two segments (PERMANOVA, P-segment: period factor, F:2.17, P-value = .002, I-segment: F: 3.2908, P-value = .001). However, the temporal changes were more pronounced at the I-segment (Fig. 3), exhibiting a different community composition in every phase transition, except in that from the Contraction to the Non-flow phases (Fig. 3B, pairwise test, Table S4). At the P-segment, conversely, the Expansion phase was the only period with a differentiated community composition (pairwise tests comparing period by period; Table S4, Fig. 3A). These differences were reflected in variations in the alpha bacterial diversity between hydrological periods. Richness (q = 1) increased during Fragmentation and decreased during the Non-flow phases, especially at the I-segment. However, no statistically significant differences were measured in the Hill numbers between hydrological periods (ANOVA, P > .05). The relative abundance of bacterial taxa differed among the hydrological phases occurring in the Algars river. The expansion phase showed a higher relative abundance of Gammaproteobacteria (26.4% versus 17.4% in other periods) and Cyanobacteria (8.1% versus 3.2%), and a lower proportion of Actinobacteria (5.3% versus 12.8%) (Fig. 2). Then, representatives of the families Comamonadaceae and Oxalobacteraceae (both within the class Gammaproteobacteria) were the most prevalent in samples collected during the Expansion phase. During the Non-flow phase, sequences affiliated to the class Alphaproteobacteria and the class Thermoleophilia increased and the relative abundance of sequences affiliated to phylum Cyanobacteria decreased, reaching nearly 0 values at the I-segment (Fig. 2).

The bacterial community also differed across habitats in each of the segments (PERMANOVA, P-segment: F: 4.29, P-value = .001; I-segment: F: 2.02, P-value = .001) (Fig. 3). Differences between habitats were more important at the P-segment (Fig. 3A), except for the Bank sediment and Floodplain habitats that exhibited similar taxonomical composition. At the I-segment, significant differences were observed between Pool versus Floodplain and between River versus Floodplain (Fig. 3B, Table S4). In both segments, the relative abundance of Thermoleophilia and Actinobacteria was higher at the Floodplain sediments. Cyanobacteria were more prevalent in Pools at the I-segment and in riffles at the P-segment (Fig. 2). Besides, the taxonomical composition in the pool habitats was distinct at each hydrological phase of the I-segment. These changes in composition between habitats were related to differences in their alpha diversity types (ANOVA, P-value = .008 for Hill q = 0, P-value = .02 for q = 1, and P-value = .048 for Hill q = 2). Overall, the effective number of ASVs singularly increased in the bank sediments in both segments, showing higher bacterial richness (q = 0), higher diversity (q = 1), and a higher number of the most abundant taxa (q = 2).

Finally, we further explored the temporal and spatial changes in the number of shared ASV combining periods and habitats at each segment (Fig. S4). Interestingly, a high number of taxonomical units were exclusively located in pools. Specifically, at the P-segment, many ASVs were observed in most of the habitats and periods, indicating a high level of homogeneity in bacterial assemblages that promotes a greater number of shared ASVs. In contrast, at the I-segment, a high number of ASVs were exclusively detected in pools, with different ASVs observed in each hydrological period (E:P, F:P, and C:P sets) (Fig. S4).

Environmental drivers of sediment bacterial composition

The environmental variables of the sediments included in the CCA explained a respective 15.2% and 16.3% (variance observed) of the variability in the composition of bacterial communities for each river segment (Fig. S5). Sediment organic matter (AFDW), w-DOC and water content exhibited significant effects in the I-segment (P-value = .003, P-value = .018, and P-value = .003, respectively) whereas in the P-segment, sediment AFDW and also water content were found to have significant P-values, explaining variations in the bacterial community composition (P-value = .024 and P-value = .001, respectively) (Fig. S5).

Dispersal versus selection: results of the null model

The null modeling approach revealed that shifts in bacterial composition between hydrological phases were primarily attributable to the interplay between selection (deterministic process) and stochastic forces (random influence). The Raup–Crick index used in Fig. 4, serves as a probability metric (ranging from −1 to 1) that indicates whether local communities are more similar (close to 1), or less dissimilar (approaching −1) than expected by random chance. Analysis of the pairwise dissimilarities among the three Hill numbers showed varying degrees of consistency, with minimal differences between river segments and hydrological phases (Fig. 4). Particularity, values closer to 1 for q = 1 and q = 2 in the Raup–Crick index (central and right panels in Fig. 4) suggest that environmental conditions predominantly governed the temporal occurrence of the diversity and most abundant species in the two segments. When values approach 1, it indicates that bacterial communities are more similar than expected by chance, implying that differences between samples are primarily due to environmental selection. Conversely, values closer to 0 for q = 0 (Raup–Crick index, Fig. 4) suggested that the temporal occurrence of species richness may be more subjected to stochastic forces. Values close to −1 indicate that the communities are more dissimilar than expected by chance, likely due to higher dispersal between local communities. This implies that fluctuation in species richness (Raup–Crick index q = 0), was more strongly influenced by random processes rather than deterministic environmental factors, such as those affecting the most abundant species (Raup–Crick index q = 2, Fig. 4).

Co-occurrence networks

We built a bipartite network to calculate the dissimilarity (beta diversity) among taxa, which considered both the spatial and temporal dimensions. This analysis calculates the dissimilarity in each of the subwebs but considers only the taxa observed in the paired subwebs (component β_OS) (Fig. 5A and B). The average beta diversity for interactions (β_OS), was 0.04 for segment I and 0.14 for segment P. These values indicate that the interactions among co-occurring ASVs were relatively low, showing a low number of shared taxa between each node of the paired network. This was particularly evident in the I-segment (lower β_OS values), where each habitat and period contained a substantial number of unique ASVs but exhibited a low number of shared ASVs. The β_OS component was higher at the P-segment, particularly between the NF and E periods, suggesting a higher number of common bacterial taxa between habitats in these periods (Fig. 5B).

The pairwise comparison of the networks revealed significant dissimilarities between periods at the I-segment (component β_WN), with values ranging between 0.09 and 0.67 (Fig. 5C and D). This observation unveils the singularity of the NF period against the others. This dissimilarity increased as much as the river contracted at the I-segment (Fig. 5C). Contrastingly, the dissimilarity between networks at the P-segment reached, ranged from 0.12 to 0.23, with higher values during the Expansion period (Fig. 5D).

Sources of bacterial taxa in river sediments

We evaluated potential sources of bacterial cells in the riffle habitat, the one more aquatic among those existing habitats. Accordingly, we calculated the relative contributions of bacterial sources using Source Tracker analysis (Fig. 6). Most of the bacterial taxa from riffles were allocated to pool habitats (mean value for all sampling periods and both segments = 45% ± 9%). Remarkably, a substantial fraction of bacterial taxa (19% ± 6%) was from an unknown source. While the sources of bacterial taxa remained relatively consistent across the different sampling periods, clear differences occurred between sampling segments. Our findings revealed that, on average, 33% of the taxa in the P-segment and 42% of taxa in the I-segment had a terrestrial origin (from floodplain and bank sediment habitats), indicating an influence of terrestrial environments on the composition of bacterial communities in the aquatic habitats (our riffles). Of particular interest is the relevant contribution of floodplain habitats as a source of exogenous bacteria at the I-segment during the contraction period (around 18%), whereas the contribution of this habitat type was found negligible at the P-segment (Fig. 6).

Discussion

River bacterial communities are shaped by a trade-off between dispersal and environmental selection, and it is well-established that the relative importance of these mechanisms varies both spatially and temporally (Malazarte et al. 2022, Stadler and del Giorgio 2022). In our study, we observed substantial differences in bacterial community composition across distinct hydrological periods. The results were likely attributed to environmental selection rather than randomness or hydrology-influenced dispersal mechanisms, as confirmed by the null modeling approach applied. Primarily, this response was noted in a dominant group of bacteria within each habitat type and hydrological period, reacting to prevailing environmental factors. Notably, we observed distinct selection pressures between permanent and intermittent river segments, with co-occurrence networks analysis indicating a low number of shared taxa between the hydrological phases, particularly evident in the intermittent segment. This observation suggests that each habitat (including riffles, pools, and bank and floodplain sediments) offers unique conditions for community assembly, predominantly influenced by differences in water content and content of organic matter in sediments, as is explained in the CCA results. We observed significant taxonomic variation between pools within the same river segment and hydrological period, implying that distinct environmental conditions shape the taxonomic composition of bacterial communities at each pool through filtering and reducing the number of common taxa. Flow interruption cannot only restrict species dispersal across river sediment habitats but also modify the environmental conditions of the remaining aquatic habitats, thereby promoting higher taxonomic variability in isolated pools. This variability is influenced by both biotic and abiotic interactions (Fazi et al. 2013). We observed a greater prevalence of Cyanobacteria in certain pools within the intermittent segment, where elevated water temperature and low dissolved oxygen levels (levels around 5 mg/l, temperature >15°C) likely favored these taxa (Paerl and Otten 2013). These findings support the notion that species sorting predominantly influences the composition of sediment bacterial communities, consistent with trends previously reported for stream biofilms (Besemer et al. 2009). The sessile nature of river sediments provides an elegant solution to escape of continuous water mixing, reducing losses, and increasing residence time (Meyer-Reil 1994) compared to the planktonic lifestyle (Lindström and Östman 2011, Niño-García et al. 2016). It is already known that the importance of dispersal versus selection mechanisms differs between planktonic and benthic communities (Gweon et al. 2021, Malazarte et al. 2022). While the composition of bacterioplankton communities is mostly influenced by water residence time (Besemer et al. 2012, Borrego et al. 2020), sediment bacteria are more influenced by water–terrestrial interactions. Previous studies have emphasized that the bacterial colonization of benthic habitats, such as microbial communities inhabiting sediments, maintains a consistent taxonomic set across the river network, as these habitats are less influenced by dispersal mechanisms associated with hydrology (Niederdorfer et al. 2016, Wisnoski and Lennon 2021) and resulting in lower intersite variability across the river network (Gweon et al. 2021, Malazarte et al. 2022).

The Algars River, studied herein, exhibits temporary characteristics at least for some parts of its network (Llanos-Paez et al. 2023). During flow contraction, isolated pools serve as refuges for riverine reach biota (Bonada et al. 2020), and foster the growth of specific bacterial taxa (Fazi et al. 2008). Source tracker analysis results underscore the significance of river pools as reservoirs during drying periods, and their influence on restoring the community composition during subsequent rewetting periods. Additionally, this analysis revealed an increased presence (~45%) of terrestrial taxa (originating from floodplains and bank sediment habitats) in the aquatic habitats, specifically in the riffles of the intermittent segment. This is likely a result of enhanced connectivity between habitats driven by greater hydrological variability at small spatial scales. Thermoleophilia, Actinobacteria, and Bacilli, associated with terrestrial habitats, exhibited higher relative abundances at the intermittent segment where the water content was low. In fact, members of both Actinobacteria and Thermoleophilia have been observed in extremely arid soils and are well-known for their resilience to drought (Bouskill et al. 2013, Gionchetta et al. 2020). Particularly, members of the phylum Actinobacteria have been described as exceptionally tolerant to desiccation and solute stress (Barnard et al. 2013, Stevenson and Hallsworth 2014). On the other hand, the high relative abundance of Bacilli (phylum Firmicutes) across all habitats in the intermittent segment suggests their role as primary colonizers following drought periods. This is supported by the fact that this class includes numerous spore-forming members adapted to arid soil environments (Fazi et al. 2008, Aslam et al. 2016). Overall, the high prevalence of desiccation-tolerant taxa in the intermittent segment indicates their remarkable adaptability to hydrological fluctuations including phases of complete dryness.

The Algars River experienced a significant flood event preceding the Expansion phase (Fig. S1), which led to sediment mobilization, re-established connectivity among river habitats, and influenced the composition of bacterial communities in both segments. Despite this extreme event, the taxonomic composition of bacterial communities in terrestrial habitats remains unaffected. The studied terrestrial habitats (including bank and floodplain sediment habitats) exhibited lower taxonomical variability compared to aquatic habitats, including a group of ubiquitous taxa that demonstrated resilience to hydrological disturbance. These results suggest that while the “terrestrial” pool of taxa from temporary river sediments is resilient to extreme hydrological events, the “aquatic” pool appears more vulnerable to hydrological variability. This observation has important implications for understanding the temporal dynamics of riverine bacterial communities postdisturbance. It suggests that the most abundant taxa in nearby terrestrial habitats may influence the composition of bacterial communities in aquatic habitats through a priority effect. This effect occurs when the prior occupation of a habitat inhibits the colonization by new individuals of that space (Hubbell 2005, Lepori and Malmqvist 2009). Consequently, the dominance of certain taxa in the river (located in terrestrial habitats) might prevent the establishment of other competitors. Thus, priority effects could influence the taxonomical composition in terrestrial habitats of temporary rivers subjected to long-lasting, high hydrological variability, leaving their imprint in the aquatic bacterial communities at any given period (Fukami 2015, Vass and Langenheder 2017). Considering the predicted increase in the duration of drying periods in many arid and semiarid regions (Arias et al. 2021), we anticipate that the decrease in sediment moisture will modify the composition of bacterial communities in temporary streams toward communities more enriched in terrestrial taxa adapted to desiccation and harsh environmental conditions (Shade et al. 2012, Evans and Wallenstein 2014).

Conclusions

Our findings indicate that environmental filtering predominantly shapes the bacterial composition of dominant taxa in temporary stream sediments. This association is linked to the substantial spatial heterogeneity across the floodplain and aquatic habitats, particularly observed by the variation in water and organic matter content in sediments. Priority effects could be an important factor influencing the bacterial composition of sediment temporary rivers, especially affecting the terrestrial pool of taxa from the floodplain habitat. Our findings underscore the adaptability of bacterial communities to hydrological fluctuations, particularly noting the high prevalence of terrestrial taxa in the flow intermittent segments. Our findings also show that the occurrence of isolated pools may act as reservoirs of true aquatic bacterial taxa, underscoring their role in maintaining biodiversity during drying periods. Since drying periods will increase in the current climate change scenario, we predict shifts in bacterial diversity toward communities enriched in taxa adapted to arid soils and dry conditions, driven by the decrease in the moisture retention capacity and the reduction of isolated pools. Lastly, our study emphasizes the importance of integrating temporal variability and small-scale habitat dynamics into the analysis of bacterial metacommunity assembly in temporary streams.

Acknowledgements

The authors acknowledge the field work assistance of Julio López-Doval and Nina Witteveen as well as the help on DNA extraction of Àlex Sanchez.

Author contributions

Anna Freixa (Conceptualization, Data curation, Formal analysis, Investigation, Visualization, Writing – original draft), Juan David González-Trujillo (Data curation, Formal analysis, Software, Visualization, Writing – review & editing), Oriol Sacristán-Soriano (Data curation, Formal analysis, Writing – review & editing), Carles M. Borrego (Data curation, Formal analysis, Writing – review & editing), and Sergi Sabater (Conceptualization, Funding acquisition, Investigation, Writing – review & editing)

Conflict of interest

The authors declare no conflict of interest.

Funding

This research was funded by the project RIVSTRESS (PID2020-115708RB-C22) of the Spanish Ministry of Science and Innovation (MCIN). We acknowledge funding from the CERCA program and the support through the Consolidated Research Group ICRA-ENV 2021 SGR 01282. A.F. acknowledges the Juan de la Cierva program (IJC2019-039181-I) funded by MCIN/AEI.

Data Availability

Data sets used for this research (DNA sequence dataset) are archived in the public repository NCBI Sequence Read Archive database with accession number PRJNA877098 following this link provided for peer review: https://dataview.ncbi.nlm.nih.gov/object/PRJNA877098?reviewer=cqd48k7lvq1sr1p08sol96upjl. The script with the pipeline to implement the null model is included as supplementary material.

References