Predicting Macroinvertebrate Responses to Water Abstraction in Alpine Streams

Exploitation of hydropower potential in alpine areas undermines the ecological integrity of rivers. Damming and water abstraction substantially alter the physical habitat template of rivers, with strong repercussions on aquatic communities and their resources. Tools are needed to predict and manage the consequences of these alterations on the structure and functioning of macroinvertebrate communities and resource availability in alpine streams. We developed habitat preference models for taxa, functional feeding guilds, and organic resources to quantify the effects of discharge alteration on macroinvertebrate communities in two alpine streams. Our physical habitat model related an indirect measure of bottom hydraulic forces (FST hemispheres) to the distribution of macroinvertebrate taxa and their resources. We observed that flow-dependent habitat availability for macroinvertebrate communities generally decreased with increasing water abstraction. We were able to relate these changes to near-bed hydraulic conditions. Our results suggest, however, the existence of upper discharge thresholds delimiting optimal habitat conditions for taxa. In contrast, we found weak effects of near-bed hydraulic conditions on resource distribution. Overall, our findings contribute towards predicting the impacts of water abstraction on macroinvertebrate communities in small alpine streams and the benefits of baseflow restoration.


Introduction
The hydropower potential of alpine streams is widely exploited for electricity generation. This involves the construction of dams for high-head storage and the development of a network of conduits for water abstraction, trans-catchment transfer and water return to downstream rivers after turbination. The resulting impacts range from baseflow modification to reductions in flow variability, often causing interruption of the river continuum [1,2] and alterations to natural sediment regimes (e.g., retention, flushing) [3,4]. The mechanisms and magnitude of downstream impacts depend on factors related to sediment and water availability, which can cause drastic changes in habitat conditions relevant for local biota [1]. The ongoing intensification of hydropower utilization as an alternative to fossil fuel use, in combination with the effects of climate change [5], will further increase pressures on stream ecosystems globally. In the coming years, the strategic role of hydropower will be consolidated, new hydropower infrastructure will be created, and concessions for hydropower operators will be re-negotiated [6]. These developments not only concern large and conspicuous hydropower schemes, but also, increasingly, those on small streams [7]. At each measuring spot along the transects, we used FST (FliesswasserStammTisch) hemispheres developed by Statzner and Müller [40] to estimate bottom shear stress [40,41], a key feature of near-bed hydraulic conditions. Shear stress controls abiotic characteristics of benthic habitats in streams and exerts selective pressure on flow-related traits of macroinvertebrates. Benthic macroinvertebrate distributions respond to FST range [14], indicating the suitability of FST hemispheres to assess flow-related macroinvertebrate habitat preferences. FST-hemispheres have identical shape and volume but a gradient in weight, due to differences in density of the material used [40]. They are thus labelled from no. 1 (the lightest) to no. 24 (the heaviest). Spot-scale shear stress was estimated from the number of the heaviest FST hemisphere that was not moved by flow. The operational procedure consists in horizontally inserting a base plate (13 × 18 × 0.8 cm) into the top layer of sediment. Hemispheres are then placed on the plate, one at a time, For each stream, one~70-m long reach was selected as study reach. The length of the study reaches was defined to obtain hydro-morphologically uniform reaches of sufficient length that provide enough evenly spaced replicate sampling spots needed for the habitat models (described below). The reach on the Lesgiüna (Figure 1) was a long run with a streambed mostly composed of medium-large cobbles (a-axis between 32 and 128 mm). The study site was located~300 m below the exit from a narrow gorge, and~200 m before entering the Brenno floodplain. This river section had a relatively homogeneous width (~9.5 m), laterally controlled by rock bank protections. The reach on the Orino (Figure 1) Water 2021, 13, 2121 4 of 22 had more complex morphology and featured two riffle-pool sequences. The study site was located~200 m before the confluence with the Brenno river and~200 m downstream from the last weir. The average channel width was 5 m; the streambed was similar to the Lesgiüna, although larger particles were more embedded (authors' observation).

Data Collection and Treatment
Physical habitat features (flow velocity, depth, sediment size, and FST) were measured along 17 transects (~4 m apart) perpendicular to flow at each site, at three discharges. Discharges differed roughly by a factor of 4: low winter discharge (baseflow), medium early spring discharge, and relatively high snowmelt discharge in late spring (respectively, Orino: 0.02, 0.09, 0.43 m 3 ·s −1 ; Lesgiüna: 0.06, 0.26, 1.15 m 3 ·s −1 ). Seventeen transects were the maximum feasible for sampling in one day with a team of three people. Doing so ensured that all measurements and samples were taken under constant discharge conditions. On each sampling day, discharge was calculated from measurements taken with a portable electromagnetic water flow meter (OTT MF pro, OTT, Kempten, Germany) at~30 cm intervals, at three to four depths per point along the same cross-section. Crosssections had a regular profile and were located immediately upstream of the study reaches. To maximize the accuracy of measurements, the streambed was smoothened by removing larger cobbles. The same instrument was used to measure depth and bottom flow velocity (30 s averages) at five evenly spaced measuring spots per transect, with the first and last spot of each transect ca. 25-30 cm from shore. Relative abundances of nine substrate size classes (<0.06; 0.06−2; 2−8; 8−16; 16−32; 32−64; 64−128; 128−256; >256 mm) were visually measured at each spot using a 20 × 25 cm grid and a weighted mean of sediment size was calculated (S1).
At each measuring spot along the transects, we used FST (FliesswasserStammTisch) hemispheres developed by Statzner and Müller [40] to estimate bottom shear stress [40,41], a key feature of near-bed hydraulic conditions. Shear stress controls abiotic characteristics of benthic habitats in streams and exerts selective pressure on flow-related traits of macroinvertebrates. Benthic macroinvertebrate distributions respond to FST range [14], indicating the suitability of FST hemispheres to assess flow-related macroinvertebrate habitat preferences. FST-hemispheres have identical shape and volume but a gradient in weight, due to differences in density of the material used [40]. They are thus labelled from no. 1 (the lightest) to no. 24 (the heaviest). Spot-scale shear stress was estimated from the number of the heaviest FST hemisphere that was not moved by flow. The operational procedure consists in horizontally inserting a base plate (13 × 18 × 0.8 cm) into the top layer of sediment. Hemispheres are then placed on the plate, one at a time, starting from the higher numbers and noting the number of the densest hemisphere that is moved by the force of the current. When even the lightest hemisphere was not moved, the measurement was ascribed to the fictional "no. 0" hemisphere. Accurate hemisphere placement was challenging with the original FST hemisphere setup [40]: hemispheres strongly adhered to the plate due to the smooth contact surfaces, thus impeding accurate measurements [42,43]. This effect was avoided by placing a paint-coated grid (3 mm wire, 1 cm mesh size) on the base plate to create a structured surface [44].
Sampling of benthic macroinvertebrates and particulate organic matter (fine POM (FPOM): 0.25−1 mm; coarse POM (CPOM): > 1 mm) was carried out in early spring 2016 during a period of stable discharge conditions. Sampling in this period minimized the risk of recent disturbance by flood events. In each reach, three replicate samples were taken from spots with FST numbers between 0 and 9 (three samples per number). Because FSTs were measured last in the sampling spot to minimize disturbance of the macroinvertebrate community and resources, we knew the exact FST number only after collecting the benthic sample. We thus completed the three replicates for each FST number 0 to 9 by visually selecting remaining sampling spots in the study reaches. Further samples were collected from spots with FST numbers up to no. 14 to improve preference-curve fitting (described below). This procedure resulted in 60 benthic samples (Lesgiüna (N = 32): three samples for FST hemispheres from no. 0 to 7, two for no. 8, three for no. 9, one each for nos. 10, 12, 14; Orino (N = 28): three for nos. 0 to 8 and one for no. 12). Benthic samples were collected using a 25 × 25 cm Surber sampler (mesh size 0.25 mm). Sediment composition was noted, and flow velocity and depth were measured as for the physical habitat survey. In front of the Surber net, larger cobbles were brushed, and the substrate homogenously disturbed with a screwdriver to a depth of 14 cm for 60 s. Samples were stored in 70% ethanol. FST hemisphere numbers were measured after macroinvertebrate sampling because previous substrate disturbance does not seem to affect the measurement [40].
Water physico-chemistry data ( Table 1; collected four times between February and May 2016) was obtained with a handheld multiparameter probe (HACH HQ40d, HACH, Loveland, CO, USA), and from water samples analyzed in the laboratory (AUA laboratory, Eawag, Switzerland). In the laboratory, macroinvertebrates were separated from organic matter and sediment. Living roots were removed before CPOM and FPOM separation. Some samples from the Orino had a large amount of living root fragments (i.e., not part of particulate organic matter). These, as well as seeds and large wood fragments, were removed from all samples before drying (50 • C, 48 h) and incineration (550 • C, 3 h) to measure ash-free dry mass (AFDM).
Macroinvertebrates were identified to the lowest practical taxonomic level (i.e., mostly species or genus for EPTs and Coleoptera, family, or tribe for Diptera, and variable for the remaining taxa) using a dissecting microscope (M5, Wild Heerbrugg AG, Heerbrugg, Switzerland) and dichotomous identification keys [45][46][47][48][49][50]. Species-specific average individual dry-mass (DM) of macroinvertebrates was estimated from measurements (±0.05 mm) of head capsule width or body length (depending on the type of allometric equation available from the literature [51][52][53][54][55]) for a total of 40 individuals per taxa from each stream. Taxa for which allometric equations were unavailable accounted for less than 3% of the total abundance and were excluded from community biomass calculations. Dry-weight estimates were used to calculate total macroinvertebrate biomass per sample [24] and the relative biomass of dominant functional feeding groups (FFG) as follows: shredders, grazers or scrapers, passive filter feeders, detritus collector-gatherers, and predators. The relative FFG biomass was calculated by multiplying the biomass of a given taxon with an affinity score (1−10) as described in Moog [24]. In addition, the following indices were calculated per sample based on the relative abundance of macroinvertebrate taxa: taxa richness, Shannon (H), and Simpson (D) index: where N is the total number of taxa, and n i is the number of individuals of taxa i.

Data Analysis
Preference models were generated for a total of 14 taxa that met the following criteria: at least 100 individuals in total and occurrence in >10% of the samples at each site. Preference models were also fitted for FFGs, resources (CPOM, FPOM), and community metrics. The models chosen in this study were quadratic polynomial regressions [34]. Predicted values were standardized to obtain a Suitability Index (SI) [11] ranging from 0 (non-suitability) to 1 (maximum suitability or optimal habitat). Before fitting, response variables were log+1 transformed and physical habitat variables centered on the mean. These analyses were carried out using the statistical software R (version 3.6.0) [56].
To estimate discharge-dependent changes in near-bed hydraulic conditions in the two streams, we used the hydraulic model embedded in the FST-based CASiMiR module Benthos (Computer Aided Simulation Model for Instream Flow Requirements [57][58][59]). This model predicts the statistical distribution of FST numbers for the study reach under different discharge conditions expressed as the respective area of each FST number (m 2 ). The prediction is based on a density function of FST numbers that was calibrated with field measurements, by fitting the observations to a log-normal distribution.
The model also coupled the predicted FST distribution with the preference models (see above) to calculate the Weighted Usable Area (WUA) and the Hydraulic Habitat Suitability (HHS) for each response variable and each discharge [11].
where N is the total number of inundated spatial units (measuring spots) within the reach [-]; A i is the surface area of spatial unit i [m 2 ]; and HSI i is the habitat suitability index of spatial unit i [-]. The models allowed predicting taxonomic and functional community composition in response to discharge conditions. Preference models were fitted to data collected in spring, during intermediate discharge conditions (Table 2). The hydraulic models were fitted using data measured in winter (low discharge), early (medium discharge) and late spring (high discharge). The difference between HHS and WUA modelled at the low and mid discharge (∆ 1 ), (Table 2) and the cumulative discharge increase ∆ 1 +∆ 2 (∆ 2 being the difference between mid and high discharge) were calculated to highlight changes in HHS and WUA with changing discharge. Table 2. Measured values of wetted area; FST µ , mean FST number approximated to the closest unit; V µ , mean velocity; D µ , mean depth (mean ± standard deviation); DGS, dominant grain size class under different discharge conditions (season).

Water Physico-Chemistry and Physical Habitat Assessment
Nitrogen, phosphorus, and dissolved organic carbon concentration levels were low and varied little between streams and among sampling dates (Table 1). Dissolved oxygen levels were always near saturation. Water temperature was on average two degrees higher on the Orino, where also a slightly higher fluctuation was observed.
FST mean number, water velocity, and depth varied consistently between the two streams at the three different discharges ( Table 2). The dominant sediment grain size class was constant within sites except that it was higher in late spring at the Lesgiüna. The two smaller sediment size classes (<0.06 mm; 0.06−2 mm) were relatively more abundant in the Orino than the Lesgiüna (12% versus 2%). By contrast, Lesgiüna had substantially more boulders and large cobbles (36%) than Orino (11%).

Hydraulic Model
At both sites, increasing discharge reduced the area of 0 numbers, increased low-to-high numbers (1-10), and slightly increased very high numbers (11)(12)(13)(14)(15)(16). FST numbers between 4 and 10 showed the greatest relative increase in surface, which means that, as expected, streambed shear stress increased overall with increasing discharge. Measured FST distribution during transect surveys, and a comparison between measured and modelled frequency distributions of FST numbers at both sites are reported in Appendix A ( Figures A1 and A2).

Preference Model
The curves fitted using the preference models varied among the response variables. For some variables, the optimal suitability value (SI max = 1) was outside the range of modelled FST numbers (0−16). Three main patterns of responses were observed (Table 3): (a) a bell-shaped curve showing marked increasing/decreasing preference and an optimum peak (e.g., B. rhodani, B. alpinus); (b) a bowl-shaped curve, close to linear, with rapidly increasing suitability starting from~0 and no detected optimum (e.g., Protonemura spp.,

Preference Model
The curves fitted using the preference models varied among the response variables. For some variables, the optimal suitability value (SImax = 1) was outside the range of modelled FST numbers (0−16). Three main patterns of responses were observed (Table 3) Models for taxa including Alainites muticus, Tanypodinae, Leuctra spp., Orthocladiinae, Ecdyonurus spp. and Isoperla spp. had a low fit (R 2 < 0.2). For the FFGs, models for grazers and passive filter feeders had the best fit, whereas model fit was poor for those of shredders, detritus collector-gatherers, and predators. Of the calculated community metrics, only total biomass and total abundance had R 2 ≥ 0.2. Models for resources (CPOM and FPOM) performed poorly (R 2 < 0.2). Models for taxa including Alainites muticus, Tanypodinae, Leuctra spp., Orthocladiinae, Ecdyonurus spp. and Isoperla spp. had a low fit (R 2 < 0.2). For the FFGs, models for grazers and passive filter feeders had the best fit, whereas model fit was poor for those of shredders, detritus collector-gatherers, and predators. Of the calculated community metrics, only total biomass and total abundance had R 2 ≥ 0.2. Models for resources (CPOM and FPOM) performed poorly (R 2 < 0.2).
We combined the hydraulic habitat models with the preference models to estimate HHS and WUA for taxa, FFGs, organic resources, and community metrics at different discharge conditions (Appendix B, Tables A1 and A2). Based on these results, ∆HHS and ∆WUA were calculated to highlight changes with different discharges (Figure 3). Both sites showed similar response patterns. Because of differences in stream size, ∆WUA results varied in magnitude between the two sites. Most of the taxa showed an increase in HHS for both sites as discharge increased from low to intermediate (0.09 m 3 s −1 for Orino and 0.26 m 3 s −1 for Lesgiüna), except for Ecdyonurus spp. and Tanypodinae. With further increases in discharge, B. rhodani and A. muticus had decreasing HHS, while it remained constant (Orthocladiinae) or increased for other taxa. Rheophilic taxa such as R. semicolorata, Rhithrogena spp. and Isoperla spp. showed the largest positive ∆HHS with discharge increase. For FFGs, passive filter feeders showed the strongest response in terms of ∆HHS (Figure 3). For other FFGs, the cumulative increase in HHS from low to high discharge (i.e., ∆ 1 +∆ 2 ) remained below 0.1. Similar results were found for community metrics and CPOM, while FPOM was the only variable showing a decrease in HHS for both discharge increments (i.e., ∆ 1 and ∆ 1 +∆ 2 .). semicolorata, Rhithrogena spp. and Isoperla spp. showed the largest positive ΔHHS with discharge increase. For FFGs, passive filter feeders showed the strongest response in terms of ΔHHS (Figure 3). For other FFGs, the cumulative increase in HHS from low to high discharge (i.e., Δ1+Δ2) remained below 0.1. Similar results were found for community metrics and CPOM, while FPOM was the only variable showing a decrease in HHS for both discharge increments (i.e., Δ1 and Δ1+Δ2.).  The slope of the line connecting Δ1 with Δ1+Δ2 reveals whether a further increase in discharge has positive, neutral, or negative effects on HHS or WUA compared to Δ1 alone. Positive slopes mean that further increase in discharge from mid-to high-flows further increased HHS or WUA. Negative slopes mean that the further discharge increase reduced the positive effects of Δ1 (if Δ1 > 0; e.g., B. rhodani and A. muticus) or intensified the negative effects of Δ1 (if Δ1 < 0; e.g., Ecdyonurus sp. and Tanypodinae). Zero slopes mean that Δ2 does not affect HHS or WUA.Y-axis not in scale among panels.  ) indicates that an increase in discharge from low to mid has positive effects on HHS or WUA, ∆ 1 < 0 indicates a negative effect of the discharge increase. The slope of the line connecting ∆ 1 with ∆ 1 +∆ 2 reveals whether a further increase in discharge has positive, neutral, or negative effects on HHS or WUA compared to ∆ 1 alone. Positive slopes mean that further increase in discharge from mid-to high-flows further increased HHS or WUA. Negative slopes mean that the further discharge increase reduced the positive effects of ∆ 1 (if ∆ 1 > 0; e.g., B. rhodani and A. muticus) or intensified the negative effects of ∆ 1 (if ∆ 1 < 0; e.g., Ecdyonurus sp. and Tanypodinae). Zero slopes mean that ∆ 2 does not affect HHS or WUA.Y-axis not in scale among panels.

Hydraulic Models
The physical habitat models showed a predictable decrease in high-shear stress area under low discharge conditions (Figure 2), with a shift towards very low FST-hemisphere numbers at the lowest discharge (up to 50% of spots had the lowest measurable FST number). However, shear-stress distribution (measured as FST) did not always uniformly change with discharge (Appendix A, Figure A1). At the Lesgiüna, frequencies of FST numbers showed a bimodal distribution, where the first maximum remained constant at 0 and the second moved gradually from 2-4 to 6-8 with increasing discharge. At the Orino, the FST hemispheres frequency curve peaked at 0 at the lowest discharge, was inversely proportional to FST number at the moderate discharge and resembled that of the Lesgiüna at the highest discharge. Differences in FST frequency distribution with discharge could be attributed to differences in channel morphology [60,61].
Due to gentle lateral slopes of the Lesgiüna reach, areas with minimal shear stress (FST numbers 0) were maintained on the newly inundated margins when discharge increased. In contrast, the morphology of the Orino reach was characterized by a pronounced lateral discontinuity, owing to the presence of higher vegetated banks, dry secondary channels, and small, vegetated gravel bars. This morphology was possibly due to channel incision caused by limited sediment supply below the dam. Consequently, the banks were not inundated at intermediate discharge as the flow was constrained in the main channel. This led to a loss of low-flow habitats close to the banks with increasing discharge. A return to a bimodal distribution (including some 0 FST values) at the highest measured discharge follows the re-occurrence of low-flow habitat, as some areas on the banks became inundated (data not shown).
Discharge was also related to the overall hydraulic heterogeneity and the availability of key habitats of macroinvertebrates. At very low discharge (winter baseflow), high FST values (nos. [12][13][14][15][16] were rare in both streams. Nevertheless, the few spots with high shear stress present during baseflow (Figure 2) may be crucial to maintain rheophilic taxa and long-lived taxa whose habitat preference varies with life stage [31] (e.g., Perla grandis [62]). In comparison, at high discharge, flow refugia created by large wood debris and (in our study) boulders are crucial to protect macroinvertebrates from being dislodged [63,64].

Preference Models
Generally, the modelled hydraulic preferences of taxa covered in this study were consistent with published preference models [14]. Overall, spots with intermediate shear stress (FST range 4-11) were optimal for most taxa (Table 3). However, the predictive accuracy of the preference models was low to moderate. The variance not explained by the models can have several sources, including species-specific differences in habitat preference, unmeasured abiotic and biotic factors, tolerance to a range of measured habitat conditions, and stochastic distributions. In our study, coarse taxonomic identification might explain the apparent lack of hydraulic preferences for some taxa (e.g., Orthocladiinae: R 2 = 0.01), which confirms the large variation in hydraulic preferences in this taxon [65]. In addition to physical habitat conditions, other ecological mechanisms define the distribution of taxa within stream reaches, including different habitat requirements of different developmental stages (larger individuals generally occupying spots with higher bottom hydraulic forces [66]), colonization history [67], and biotic interactions (e.g., presence of predators, competitors, and resources) [68].
Furthermore, the measured shear stress distribution was in the range of 0−14 FST numbers, thus not covering the upper portion of the FST range (the complete hemisphere set covers numbers 1 to 24). Although hydraulic preference generally declines at very high shear stress [14,32], the inclusion of spots with higher hemisphere numbers would have been important to increase the accuracy of the biological models, i.e., to capture decreasing habitat suitability with increasing shear stress for rheophilic taxa.
Preference models for FFGs had low predictive power. Nevertheless, our results suggest a preference of passive filter feeders for increasing shear stress (estimated as FST; R 2 = 0.46), probably as a consequence of their feeding mode and ability to resist relatively high shear stress [32,69]. Filter feeders in the study reaches were mostly represented by Simulium, Prosimulium and Hydropsyche sp., which strongly depend on flow for feeding. Conversely, many taxa had a variable grade of affinity to grazing, resulting in a moderate fit, but a relatively flat curve for this FFG, a consequence of the integration across taxa with variable flow-dependent habitat preferences.

Changes in Habitat Suitability in Response to Discharge
Based on our models, water abstraction (e.g., due to hydropower production or climate change) would result in less-suitable hydraulic conditions for most taxa, except Ecdyonurus spp., and Tanypodinae ( Figure 3). Conversely, increased baseflow would improve habitat suitability for most taxa, although some, such as B. rhodani and A. muticus, would suffer a slight decrease in habitat suitability at high discharge. Nevertheless, the effect of an increased wetted area with increasing discharge would in most cases compensate for reductions in HHS (reflected in WUA). Exceptions include Ecdyonurus spp. and Tanypodinae in the Orino, which prefer low discharge conditions that would not be compensated by the increase in wetted area at higher discharge in that stream ( Figure 3). The effect of increased discharge on habitat suitability was not linear. For most taxa, a change from low to medium discharge (winter to spring conditions) was predicted to have a greater beneficial effect on HHS than the increase in discharge from mid to high discharge (summer conditions in our study; Figure 3, Tables A1 and A2). This suggests that some factors (e.g., decreased suitability for taxa preferring intermediate FST) may attenuate the overall positive effects of further increasing discharge, but it also highlights the paramount importance of baseflow and residual flow management for macroinvertebrate communities in alpine streams.

Differences between Streams
In spring, the composition of macroinvertebrate communities was similar between both study sites, i.e., dominated by Baetidae and Chironomidae. However, relative abundances of respective taxa showed marked differences. As water physico-chemistry revealed no signs of anthropogenic impacts at both sites, these differences were probably mainly caused by hydro-morphological dissimilarities. The study reaches had similar average near-bed hydraulic conditions (measured as frequencies of FST numbers) during our study period despite substantial differences in discharge (Table 2). However, the frequency of spots with low shear stress was higher in the Orino. As a result, specialists for slowflow areas (i.e., limnephilic taxa) with fine sediment (e.g., Ephemera danica, Odonata and Oligochaeta) were found almost solely in the Orino. In contrast, rheophilic taxa with preference for higher flow velocities and coarser substrate (e.g., E. alpicola, Rhithrogena spp.) were more abundant in the Lesgiüna.
The more natural flow regime in the Lesgiüna, with seasonal discharge variability and the associated regular streambed disturbance, probably selects for taxa with specific life history traits that allow them to benefit from floods and low flows between floods (see [70]). On the contrary, flow regulation and overall more stable conditions at the Orino would, in the long term, favour taxa which benefit from low-disturbance conditions [71,72]. Community composition and resource availability is also driven by variability in habitat characteristics on longer timescales. Hydrological variability of the stream might be particularly relevant for alpine streams that naturally have high variability but are often constrained by water abstraction and damming. The dam regulating the Orino reduced its hydrological variability and the frequency of floods with consequences on habitat conditions at the study reach, especially when compared to Lesgiüna, as shown in the occurrence of fine sediment cover and vegetated patches particularly in the Orino. The long-term (seasonal) effects of reduced discharge variability in the Orino were evident in the presence of fine sediment patches and the occurrence of macroinvertebrate taxa that prefer these patches, despite similar flow conditions in both streams during the study.

Management Implications and Perspectives
Alpine river management in Europe is subject to two main bodies of law, the Water Framework Directive (WFD) [73] and the Swiss Water Protection Act (814.20). The former determines residual flows as a function of minimum flows (Q 347 ), without reference to ecological conditions. Good ecological status of rivers as defined in the WFD is determined by the assessment of the status of some biotic indicators, including macroinvertebrates. However, despite the well-acknowledged influence of flow regime on ecological status, the hydrological condition is solely intended as a supporting element of these biotic indicators [74]. Our results highlight two main findings that should be considered when assessing the status of macroinvertebrate communities based on flow-dependent habitat conditions and predicting the effect of flow alteration/restoration. We observed that: (1) a decline in baseflow below a certain stream-specific value can substantially alter the composition of the community, and (2) discharge-dependent thresholds could determine optimal conditions for macroinvertebrate taxa. These findings, merging biotic with abiotic conditions, could support the definition of ecologically relevant residual flow targets in streams where hydropower utilization is the dominant impact.
Taxonomic resolution influences the predictive accuracy of preference models. For instance, the large functional diversity within the subfamily Orthocladiinae could explain the poor fit of its preference model. Life history and seasonality probably further contribute to low predictive accuracy. Preference models based on sampling in a specific season are not representative for all larval stages, in particular if only one stage is present during sampling ("survey effect" [14]). In our study, the univoltine Ecdyonurus spp. was present as early instars (head capsule width 1.38 ± 0.88 mm), which were sensitive to strong hydraulic bottom forces (Table 3). This contrasts with its characterization as a rheo-to limnophilic taxa based on other larval stages [75], and with published hydraulic preference models [14,76]. This issue might be less relevant for multivoltine species, where the presence of several cohorts during sampling might result in preference models that are more representative for the species. However, it shows the value of study-specific preference models that can then be compared to published ones for interpretation. Seasonal differences are likely to play a role also for community-level and functional responses. For instance, Merigoux and Dolédec [32] observed that taxa richness was negatively correlated with FST numbers in spring and vice-versa in autumn. Seasonal differences in ecological conditions can modify the intensity of biotic interactions (in addition to the occurrence of different larval stages discussed above), and in turn, alter physical habitat preferences of taxa.
Applying the FFG concept to the same preference models based on flow-dependent physical habitat conditions resulted in high model fits only when feeding strategy was strongly related to flow (i.e., for passive filter feeders; R 2 = 0.46). Models were less accurate for other FFG, probably because of weak direct relationships between the distribution of their respective resources and near-bed hydraulic conditions (the basis of our hypothesis between FFG and flow conditions), as found in this and other studies [21,28,77]. This conclusion is supported by the distribution of shredders, which was strongly related to that of CPOM, whereas neither had a clear relationship with near-bed hydraulic conditions. Evidently, generalization of our results requires further spatial and temporal validation. For instance, higher altitude portions of the catchments should be considered to cover a more representative gradient of conditions of alpine streams. Survey size, i.e., the range of hydraulic conditions and number of sites assessed, also strongly influences the definition of optimal habitat conditions of taxa and its precision [32,[78][79][80]. Developing representative physical habitat models and preference curves for different stream types [60,61,81] would probably facilitate their generalization while also reducing the amount of field data required.
Community composition and relative abundances of taxa are influenced by present environmental conditions, including the effects of anthropogenic stressors such as water abstraction and climate change. Community adaptation to new conditions is a non-linear, multi-faceted ecological response involving several levels of change associated with the physical habitat that influences macroinvertebrate communities but also indirect effects through changes in resource availability and biotic interactions [82,83]. Indeed, it is challenging to predict long-term community shifts and ecosystem adjustments when a stressor is removed or added to an ecosystem (e.g., baseflow restoration or increase in water abstraction). We found near-bed hydraulic preference of macroinvertebrate taxa that are consistent with those reported from other studies. However, further development and validation of habitat models for alpine catchments and increasing seasonal and altitudinal range will be necessary to better describe hydraulic preferences of macroinvertebrate taxa and the distribution of their resources in alpine streams. This is needed to support policy makers, environmental managers, and hydropower operators in predicting the impacts of water withdrawal on taxonomic and functional macroinvertebrate community composition, as well as the ecosystem processes and services that depend on it. Our results point towards the existence of flow-related thresholds for the optimal near-bed hydraulic conditions to safeguard diverse and well-functioning alpine macroinvertebrate communities. However, adaptive water management should not only focus on base-flow levels but also draw on the recent advances in environmental flow design. These account for the natural flow regime (e.g., seasonal flow variability [84,85]) to sustain hydrogeomorphic processes that control ecological dynamics in rivers [86][87][88].     Appendix B
x Odontoceridae x x Enoycila pusilla x Psychodidae x Stratiomidae x Thaumaleidae x Clinocera spp. x