https://doi.org/10.5061/dryad.brv15dvkb

Description of the data and file structure

Files and variables

File: FunctionalTraits.xlsx

Description: Functional trait database adapted from Tachet et al. 2010

Variables

FunTraits excel sheet:

•  Final_assignment: Final taxonomic assignment of macroinvertebrate taxa.
• FinalFun: Final taxonomic assignment for functional traits. This is a tradeoff between the resolution of the Tachet et al. functional trait database and the final taxonomic assignment (Final_assignment). When final assignment is more precise than the functional trait database, taxa were homogenised to the next best taxonomic level. For example, Baetis spp. were homogenised to the genus level.
• TaxLevel: The taxonomic resolution aquired for FinalFun
• NbGene_1 - Vit_4: Fuzzy-coded functional traits for macroinvertebrates. See the trait info sheet for additional information.

Trait info excel sheet:

• Trait name:  Fuzzy-coded functional traits for macroinvertebrates.Corresponds with names from the FunTraits sheet.
• Trait_type: Indicates the type of functional trait, note F refers to fuzzy-coded traits
• Fuzzy_name: Indicates the fuzzy trait, e.g. Nb.gen, number of generations
• Extra_info: Provides additional information for each fuzzy trait

File: CommunityData.xlsx

Description: Invertebrate and periphyton community sampling

Variables

• Sample_ID: Sample ID can be broken down into: site name => Hydrological ambiance (1,2 or 3), Sample (A, B or C), Sampling campaign (C0 to C6; see Table 1) ==> Sampling year
• Site: Sampling site
• Campaign: Sampling campaign, refers to the order of sampling in a given year
• Year: Sampling year
• Date: Sampling date
• The remaining columns provide abundance data for invertebrate (individuals per 0.10m²) and periphyton (cells per 0.01m²) taxa.

File: Microhabitat.xlsx

Description: Description of microhabitat conditions, samples correspond to invertebrate and periphyton sampling

Variables

• Sample_ID: Sample ID can be broken down into: site name => Hydrological ambiance (1,2 or 3), Sample (A, B or C), Sampling campaign (C0 to C6; see Table 1) ==> Sampling year
• Site: Sampling site
• Campaign: Sampling campaign, refers to the order of sampling in a given year
• Year: Sampling year
• Date: Sampling date
• Ambiance: Hydrological ambiance of the sample. Hydrological ambiances are broken down into: (1) 30-50 cm s-1, (2) 70-90 cm s-1 and (3) > 100 cm s-1
• Largest_substrate: Largest substrate size-class in sampling area. Refer to the README panel in the excel file for additional information.
• Smallest_substrate: Largest substrate size-class in sampling area. Refer to the README panel in the excel file for additional information.
• Dominant_substrate: Dominant (most abundant) substrate size-class in sampling area. Refer to the README panel in the excel file for additional information.
• Substrat_accessory1-6: Presence of accessory substrate size-class in sampling area in order of abundance. Up to 6 accessory substrates were recorded if present. Refer to the README panel in the excel file for additional information.
• Substrate_richness: The total number unique substrate size-classes in the sampling area. Refer to the README panel in the excel file for additional information.
• Largest_underlayer: Largest underlayer substrate size-class in sampling area. Refer to the README panel in the excel file for additional information.
• Smallest_underlayer: Largest underlayer substrate size-class in sampling area. Refer to the README panel in the excel file for additional information.
• Dominant_underlayer: Dominant (most abundant) substrate size-class in sampling area. Refer to the README panel in the excel file for additional information.
• Undercover_accessory1-3: Presence of accessory underlayer substrate size-class in sampling area in order of abundance. Up to 3 underlayer accessory substrates were recorded if present. Refer to the README panel in the excel file for additional information.
• Underlayer_richness: The total number unique substrate size-classes in the sampling area. Refer to the README panel in the excel file for additional information.
• Clogging: The degree (0-5) of interstitial substrate clogging by fine sediment and/or algae. Refer to the README panel in the excel file for additional information.
• Depth: River depth (cm)
• Water velocity (cm s-1) and associated depth measurements: | v3cm.s | ~ | Water velocity at 3cm from riverbed | | :------ | :- | :----------------------------------- | | h_0.2h | ~ | Height (from river bed) at 20% depth | | v0.2h | ~ | Water velocity at 20% depth | | h_0.4h | ~ | Height (from river bed) at 40% depth | | v0.4h | ~ | Water velocity at 40% depth | | h_0.8h | ~ | Height (from river bed) at 80% depth | | v0.8h | ~ | Water velocity at 80% depth |

File: Hydrology.xlsx

Description: 

Variables

• README is included in the xlsx file

Code/software

All statistical analyses and data manipulation were performed in R v4.20 (Core Development Team, 2023) running in RStudio (RStudio, 2011).

Spatio-temporal variation in macroinvertebrate and periphyton communities

We performed Principal Components Analysis (PCA) to characterise spatiotemporal variation in macroinvertebrate and periphyton communities. PCAs were based on total invertebrate and periphyton community density (inv/0.10m² & cells/cm² for invertebrates and periphyton, respectively), taxonomic richness (based on the highest possible taxonomic identification), Shannon and Simpson diversity indices, and taxonomic evenness. Lastly, as benthic communities have a tendency to be dominated by a few hyper abundant taxa (e.g. Battle et al., 2007; Minshall et al., 1992), we calculated the percentage contribution of the most abundant taxon to the total sample density (hereafter dominance). Functional diversity indices were also calculated for the macroinvertebrate community. Functional diversity indices were based on fuzzy-coded traits from Tachet et al., 2010, and included ecological (microhabitat, temperature and water velocity preferences, and feeding strategy) and biological (respiration, diet, locomotion, dispersion strategy, lifestyle, reproduction, number of generations per year) characteristics (Supporting information). Differences in taxonomic assignment between the invertebrate and functional datasets were resolved by summarising groups to the highest shared taxonomic-level (genus, tribe or family; taxa identified to order-level or higher were excluded from functional analyses) and averaging fuzzy trait values (Supporting information). Functionally identical taxa were regrouped into functional entities using the function sp.to.fe *from the *mFD *package (Magneville et al., 2022). Fuzzy-traits were then transformed into two continuous variables (axes 1 & 2) using Principal Coordinates Analysis (PCoA) based on Gower distances (i.e. functional dissimilarity between taxa). We then calculated functional richness, evenness and divergence using *fd_fric, fd_feve and fd_fdiv *functions, respectively from the *fundiversity package (Grenié & Gruson, 2023). The fundiversity package calculates functional metrics according to Villéger et al., 2008, wherein functional richness corresponds to the volume of the convex hull of multidimensional functional space in each sample. Functional evenness corresponds to the minimum spanning tree which links all the species in the multidimensional functional space. Functional divergence is quantified by the extent to which taxa diverge in their distance from the centre of functional space, and is weighted by taxa abundance.

Calculating β-deviation

The objective of the β-deviation metric was to quantify the degree of non-random invertebrate community structure at each sampling site and date. We calculated β-deviation as the extent to which observed spatial dissimilarity (β-diversity) within macroinvertebrate communities deviated from null expectations (i.e. based on random community reshuffling). β-deviation was estimated for both taxonomic and functional (see section 3.5.1 for trait information) β-diversity indices, separately, with the aim of providing complementary information on the processes that drive community assembly. β-diversities were calculated separately for each sampling date and site (i.e. among the nine Surber samples per campaign) as pairwise Hill numbers (q = 1) based on a Jaccard-type complement. Taxonomic β-diversities were calculated using the pair_dis *function from the *hilldiv *package (Alberdi & Gilbert, 2019) and functional β-diversities were calculated using the *beta_fd_q function (tau = mean, package mFD; *Chao et al., 2019; Magneville et al., 2022) based on Gower distances between functional entities (function *funct.dist, *package *mFD; (Magneville et al., 2022).

The choice of null model is particularly important as it imposes constraints on the null distribution generated by the model, thus modifying the ecological significance of β-deviation (Gotelli & Ulrich, 2012; Harvey et al., 1983; Vellend et al., 2014). We therefore selected the most widely used null model based on the reasoning of Kraft et al., 2011 which randomly assigns individuals to sites while preserving the overall species abundance distributions and the total number of individuals in each site (e.g. Shi et al., 2021; Siqueira et al., 2020; Tucker et al., 2016). Abundance-based null model approaches like the Kraft et al. null model have been shown to be more robust than β-deviation estimates based on presence/absence data (Tucker et al., 2016). It should be mentioned that the Kraft et al. null model has been criticised, notably, for constraining species abundance distributions. Qian et al., 2013 has pointed out that species abundance distributions actually arise largely due to deterministic processes (e.g. Alonso et al., 2008). However, as the abundance distribution of invertebrate communities varies greatly over seasonal timeframes, we justify our use of the Kraft et al. model by our need to explicitly control for such temporal variation in species abundances. β-deviation was calculated as the standardised effect size of a two-sided test comparing the simulated and observed β-diversity values using the oecosimu *function from the *vegan package (n = 1000, method = r2dtable; Oksanen et al., 2008; Patefield, 1981). There are three possible outcomes for this test: (i) the mean β-diversity value falls within 95%CI of simulated values, indicating that dissimilarity could have arisen solely from stochastic processes, (ii) the mean β-diversity value exceeds the 95%CI of simulated values, indicating overdispersal of community variation and (iii) the mean β-diversity is significantly lower than the 95%CI of simulated values, indicating clustered community variation. Positive (i.e. outcome ii) and negative (outcome iii) β-deviation values are expected to arise from a dominant role of heterogenous and homogenous selection processes, respectively. It should be noted that non-significant β-deviation values (outcome i) can also arise due to opposing deterministic forces of community assembly. For example, an equal contribution of heterogenous and homogenous selection would be expected to result in non-significant β-deviation. In the case of a non-significant β-deviation value (i.e. p > 0.05), the value was set to zero. This two-sided test assumes that the null distribution is normally distributed. An assumption that is often invalidated in the case of null model generated metrics (Bernard‐Verdier et al., 2012), including the present study. We therefore rank transformed our β-deviation metrics into 10 distinct classes in order to reduce the magnitude of differences among standardised effect sizes driven by non-normal null distributions.

The trait-environment relationship

We quantified the trait-environment relationship as an additional indicator of deterministic community structure (Cottenie, 2005), and tested its relationship with β-deviation. The trait-environment relationship was quantified as the goodness-of-fit between invertebrate functional β-diversities and environmental β-diversities. Invertebrate and environmental β-diversity were both determined within the context of each sampling date at each site. Invertebrate assemblages were compared using functional β-diversities (see section 3.5.2). Environmental β-diversity was separated into two abiotic components including (i) substrate and underlayer characteristics (dominant substrate, size-class richness, smallest and largest substrate size-class, and clogging class), and (ii) hydraulic characteristics (water velocity mean, max, min, coefficient of variation and at 3cm from the riverbed), hydraulic ambiance (30-50 cm s-1, 70-90 cm s-1 or > 100 cm s-1) and water depth. To calculate abiotic environmental β-diversity we first scaled all microhabitat variables to zero mean. Microhabitat samples were then plotted using Principal Components Analysis (function* dudi.pca*, *package ade4; *Dray & Dufour, 2007) and β-diversity was determined as the Euclidean distance between points. We also quantified (iii) periphyton β-diversity as a biotic dimension of environmental β-diversity. To do so, periphyton community data (log + 1 transformed) were compared using pairwise Bray-Curtis distances. The trait-environment relationship was quantified as the Pearson’s r correlation between functional invertebrate β-diversities and environmental β-diversities (i.e. β-diversityInv. ~ β-diversityenv) for each sampling date at each site and environmental β-diversity component (i.e. substrate, hydraulic and periphyton), separately. This resulted in three trait-environment relationships related to: substrate, hydraulic and periphyton β-diversities.

The relationship between ecological conditions, β-deviation and the trait-environment relationship  

To characterise the relationship between ecological conditions, β-deviation and the trait-environment relationship we fitted causal models. We considered four distinct categories of ecological variables: disturbance, phenology, abiotic conditions and biotic conditions (periphyton and invertebrate communities; see Supporting information for a priori *model structure). However, due to incomplete periphyton data (Table 1) we fitted two causal models (model *i) one using all data types (n = 56) and (model ii) a second excluding periphyton data (n = 80). Ecological variables were all zero-mean standardised. Prior to model selection we verified the normality of taxonomic and functional β-deviation metrics, which were both skewed towards zero. We therefore log + 1 transformed both β-deviation metrics prior to rank transformation (see section 3.5.2) to meet the assumption of normality.

The first step was to select the best performing ecological predictors of β-deviation. Variable performance was evaluated for each ecological compartment according to RandomForest analysis (mtry = 10, ntree = 5000, 100 runs, package randomForest; *Liaw & Wiener, 2002), based on mean square error divided by its standard deviation. The two best performing disturbance, microhabitat, periphyton and macroinvertebrate variables were selected for causal analysis. Causal analysis was performed using the R package *piecewiseSEM (Lefcheck, 2016). The global model was composed of multiple mixed effect models, all of which included a year and site random effect. The final model was obtained by removing all non-significant causal paths. We also tested the significance of dependencies between covariates (function dSep), and included significant dependencies in the final model. We tested the relationship between β-deviation and the trait-environment relationship (i.e. invertebrate β-diversity ~ environmental β-diversity corelation statistic) using the same mixed model structure as indicated above. The β-deviation ~ periphyton trait-environment relationship was tested in model *i *(n = 56), while the β-deviation ~ hydraulic and substrate trait-environment relationships were tested in model *ii *(n = 80). Lastly, to detect site-specific effects, we visualised and tested the relationship between the trait-environment relationships and β-deviation (Pearson’s r correlation and associated *p-*value) by site and type.

Access information

Data was derived from the following sources:

• Functional traits: Tachet, H., Richoux, P., Bournaud, M., & Usseglio-Polatera, P. (2010). Invertébrés d’eau douce: Systématique, biologie, écologie (Vol. 15). CNRS éditions Paris