Exploring fish functional trait assemblages in Eastern Mediterranean rivers: a study across multiple scales using network analysis

Functional traits of river fish assemblages have rarely been studied in the Eastern Mediterranean region. A dataset of 254 Greek electrofished near-natural sites (427 samples) was analyzed to explore native fish functional structure across three spatial scales: (a) among six ecoregions, (b) within a single ecoregion, and (c) in a river basin. Overall, 76 native fish species were assigned eight functional categories, involving 24 different traits. Bipartite networks were used to interpret spatial patterns of species trait data. Two major trait assemblage types were detected in all three spatial scales: the upland type, dominated by insectivorous, benthic, rheophilic, intolerant species, and the lowland type, incorporating omnivorous, limnophilic, phytophilic, and tolerant species. In order to provide further insights between traits and key environmental variables, redundancy analysis (RDA) was employed. The RDA depicted a strong upstream–downstream environmental gradient. Upland and lowland assemblage types form two distinct functional zones along the upstream and downstream river stretches, respectively. Yet, no consistent boundary criteria seem to exist between them. Notably, within a region of high fish assemblage heterogeneity such as the Eastern Mediterranean, functional patterns follow well-established phenomena along the upstream–downstream longitudinal profile of temperate lotic waters.


Introduction
The Mediterranean-climate areas of the world hold distinct within-area environmental idiosyncrasies (Cid et al. 2017;Zogaris et al. 2021) as well as a remarkably diverse biome (Mooney 1988). Lotic ecosystems of the European Mediterranean region are primarily defined by considerable seasonal variability (Blondel et al. 2010). Streams of the European Eastern Mediterranean are typically exposed to dramatic and unpredictable climate changes (water flow, temperature, etc.), far more than most temperate rivers (Guasch and Sabater 1994). In addition to these harsh and dynamic conditions, streams and rivers are increasingly subjected to numerous anthropogenic stressors (Haidvogl 2018), placing them among the most disturbed and threatened ecosystems on the planet (Feld et al. 2016). These human interventions alter the physical characteristics of rivers and modify the entire fluvial biotic structure (Leprieur et al. 2008). At the same time, they undermine ecological processes (Taylor et al. 2006;Bellwood et al. 2012) and degrade ecosystem services (Villeger et al. 2017).
To assess degradation as well as to support management and restoration of lotic ecosystems, knowledge and understanding of species ecological functions are essential (McGill et al. 2006;Schleuter et al. 2012;Benoit et al. 2021). Within this context, it is important to comprehend various functional facets such as species traits, functional groups, and functional syndromes (suites of co-evolved traits; for terminology, see Frimpong and Angermeier 2010;Keddy and Laughlin 2021).
Function-based approaches induced revolutionary changes in community ecological research (Dawson et al. 2021). They originate from plant ecology, but quickly expanded to several other taxa including the freshwater ichthyofauna (Aarts and Nienhuis 2003; Violle et al. 2007;Cano-Barbacil et al. 2020). They were primarily employed for practical reasons, e.g., for reducing data dimensionality in fish community studies (Frimpong and Angermeier 2010;Benoit et al. 2021;Wolter et al. 2021). Yet, they pose several critical advantages compared to taxonomy-based techniques. For example, functional study results are directly comparable between different fish species pools (Ibanez et al. 2009). Additionally, they can easily be transferable to other regions and/or different spatial scales, in various community ecology applications (e.g., Pont et al. 2006Pont et al. , 2007Logez and Pont 2011;Vesk et al. 2021). Finally, they are easily implemented as a tool for a mechanistic interpretation of the cause-effect relationship between anthropogenic pressures and fish species responses (Welcomme et al. 2006;Reyjol et al. 2014).
The spatial distribution of fish in lotic ecosystems is characterized by diverse assemblage patterns. These patterns reflect the role played by geographic and historical-geological factors in fish communities (Hoeinghaus et al. 2007;Zogaris et al. 2012;Economou et al. 2016). Bioregionalization and synecology are commonly used in exploring species spatial patterns. In this context, regardless of the scale of analysis, hierarchical clustering techniques are employed in species data (e.g., Unmack 2001;Kreft and Jetz 2010;Procheş and Ramdhani 2012;Holt et al. 2013;Vavalidis et al. 2019). Yet, network analysis has grown popular as a general approach to investigate community structure and biogeographical patterns (Heleno et al. 2014), since visualization of networks improves the interpretation of complex data (Bellard et al. 2017). Networks have been used to study a wide variety of taxa in different spatial scales, e.g., global amphibian species and vascular plants (Vilhena and Antonelli 2015), global indigenous freshwater fish species (Leroy et al. 2019), and regional non-indigenous freshwater fish assemblage types (Koutsikos et al. 2021). Network analysis literally provoked a major shift in exploring bioregions and community patterns, since it appears to be more sensitive than clustering methods (Bloomfield et al. 2018;Leroy et al. 2019).
Typically, every kind of community exploratory analysis is followed by correlating community patterns with the environment. The environment drives biotic evolution in a predictable direction (Ibanez et al. 2009) and defines the fundamental niche of every species (Aarts and Nienhuis 2003). Thus, it is important to consider the gradient of the physicalenvironmental factors in understanding the strategies and dynamics of communities (Vannote et al. 1980). Several studies have explored fish assemblage structure and functions along the environmental gradients (Ibanez et al. 2007(Ibanez et al. , 2009Oliveira et al. 2012;Logez et al. 2010Logez et al. , 2013. Eastern Mediterranean rivers are naturally characterized by diverse ecological processes and considerable fish assemblage heterogeneity, compared with the rest of Europe (Ferreira et al. 2007). But their functional benchmarks as regards river fish communities have not yet been explored. None of the major pan-European relevant reports (e.g., Logez et al. 2013;Schleuter et al. 2012) has ever incorporated sufficient data of this particular area. In the present study, we attempt to investigate Eastern Mediterranean fish assemblage patterns, through an extensive research effort in Greece. The aim of the analysis was threefold: (a) to depict fish functional assemblage types in near-natural (minimally disturbed) lotic ecosystems, at three different spatial scales; (b) to provide insights on the functional structure and trait shifts that take place across spatial scales; and (c) to assess the influence of environmental factors upon the fish assemblage patterns. Ultimately, the study intends to fill in the long information gaps on freshwater fish functional traits, as emphasized by several authors about this region of the Mediterranean Benejam et al. 2015;Cano-Barbacil et al. 2020).

Study area and primary datasets
Greece stands out as a European hotspot for freshwater biodiversity (Reyjol et al. 2008;Darwall et al. 2014). The freshwater ichthyofaunal comprises one of the highest degrees of fish endemism in the Mediterranean basin (Crivelli 1996;Freyhof and Brooks 2011). Fish assemblages also exhibit distinct patterns (Ferreira et al. 2007) and a high degree of heterogeneity across the country's eight discrete freshwater ecoregions (Zogaris and Economou 2017), reflecting diverse geological and unique biogeographical processes (Economidis and Banarescu 1991).
The primary dataset of the study comprises fish community abundance data (Sample-SpecieS matrix), acquired through electrofishing field sampling of fluvial ecosystems (both perennial and intermittent lotic waters). Sampling was mainly conducted during the low-flow period by the Hellenic Centre for Marine Research (HCMR) between 2002 and 2016 (see IMBRIW-HCMR 2013;Zogaris et al. 2018). In brief, at each site, a typical stretch of the river was electrofished, waded, or by boat (single pass, no stop nets; see Pont et al. 2006). Fish were identified into species level on-site and directly released back into the river. Data recording follows standardized procedures (Breine et al. 2005;Schmutz et al. 2007;Economou et al. 2016;Tachos et al. 2016), while species taxonomy and nomenclature follow Barbieri et al. (2015).
Another separate dataset, involving 11 anthropogenic pressure variables across sites/samples (Sample-preSSure matrix; 1586 samples X 11 pressures ), has been obtained through the work of Tachos et al. (2022). In brief, human modifications at each site were assessed following an ecological degradation pressure analysis (Tachos et al. 2022). Key matrix operations are displayed in the Appendix (Appendix Fig. 5).

Data treatment and trait characterization
To explore fish assemblages in near-natural rivers, 435 samples characterized as "unimpacted" or "slightly impacted" (see Tachos et al. 2022) were incorporated in the Sample-preSSure matrix (for relevant applications, see Pont et al. 2006;Logez et al. 2013). Next, the Sample-SpecieS primary matrix was properly updated, discarding all data belonging to "impacted" sites.
The prevalence of non-indigenous species is still rather low in Greek lotic ecosystems (Koutsikos et al. 2019). Nevertheless, all non-indigenous fish species were discarded from the Sample-SpecieS matrix. Alien and translocated fishes were removed and only the autochthonous ichthyofauna was retained in each sample. During this procedure, 427 samples were retained, and eight samples implicating exclusively non-indigenous species were discarded. Ultimately, 254 sampling sites were incorporated for analysis ( Fig. 1), encompassing 34 river basins within six freshwater ecoregions in Greece. Five river basins in the northern part of the country are transboundary, shared with Albania, North Macedonia, Bulgaria, and Turkey.
After screening, the Sample-SpecieS matrix comprised abundance data for 76 autochthonous fish species (427 samples × 76 species ) and a total of 121.187 fish specimens. Abundance data were further transformed to relative abundance (proportion per sample) within Sample-SpecieS matrix. Separately, all fish species were assigned into eight general functional categories (e.g., feeding, reproduction, migration) and characterized according to twenty-four distinct traits (Table 1). This procedure developed the SpecieS-Function matrix (76 species × 8 functions ). Species trait characterization was based on related catalogs provided in several studies (Balon 1975;Aarts and Nienhuis 2003;Angermeier and Davideanu 2004;Pont et al. 2006;Noble et al. 2007;Logez et al. 2013, Barbieri et al. 2015Zogaris et al. 2018) as well as expert judgment. For those species that no reliable information was available, certain trait characterizations were left empty. The final matrix (Sample-trait matrix: 427 samples × 24 traits ) was designed by consolidating Sample-SpecieS and SpecieS-Function matrices (see Logez et al. 2013 for relevant applications). Ultimately, this new matrix (Sampletrait matrix) incorporated relative frequency data for twenty-four (24) distinct traits and for each of the 427 samples of the survey. All data matrix operations are provided in the Appendix (Fig. 5).
For all sites, altitude, slope, catchment area, source altitude, and mean temperature of January were obtained from geographical information systems (ESRI-ArcGIS v. 10.2), as well as from climatic data (years 1950-2000) from WorldClim (http:// www. world clim. org). These five key environmental descriptors encompass important physical characteristics of the sites and have been previously assessed as non-redundant by relevant extensive analysis of the study area Tachos et al. 2022).
Network analysis, assemblage type delineation, and trait classification A limited number of studies concerning either functional traits or functional structure of fish assemblages are available for the Eastern Mediterranean region. Within the context of a fish bioassessment index, Zogaris et al. (2018) provided a list of trait characterizations for much of the Greek freshwater fish fauna. Additionally, Vagenas et al. (2022) studied the trophic guild patterns of freshwater ichthyofauna, across ecoregions falling within the Greek territory. Yet, to our knowledge, no relevant research has ever been conducted in Greece or the wider Eastern Mediterranean to delineate functional patterns of indigenous freshwater fishes.
A bipartite network is created from a series of nodes connected to each other by links, called edges. In the present study, both samples and traits of species represent the nodes of the network. When a trait occurs in a sample, a linkage between the sample (node) and this particular trait (node) is drawn, creating an edge. Several edges are created this way, since several traits are present in most samples. By definition, samples cannot be connected to another sample; likewise, traits cannot be connected to another trait. Network construction and the assessment of the modularity tests were computed with Gephi v. 0.9.2 (Bastian et al. 2009) for three different spatial scales: (a) throughout the six ecoregions (matrix analyzed: 427 samples × 24 traits ), (b) in one selected ecoregion, Ionian ecoregion (matrix analyzed: 244 samples × 22 traits ), and (c) in one river basin within the previously selected ecoregion, Acheloos river basin (matrix analyzed: 96 samples × 16 traits ). The Acheloos river basin as well as the Ionian ecoregion were selected since they incorporate the largest number of sites/samples on the primary dataset. Spatial illustration of the networks was employed under the force-directed ForceAtlas2 algorithm (Jacomy et al. 2014). The algorithm clusters strongly interconnected nodes, i.e., samples that share the same traits, and disperses all other nodes that are not interconnected. The position of a node in the final network cannot be interpreted on its own; it must be compared to the other nodes, according to Jacomy et al. (2014). Therefore, for the interpretation of the network, the relative position of the nodes and the visual approximation of the graph is of great importance.
For the detection of the modules (trait assemblage types), the Louvain modularity optimization algorithm was employed (Blondel et al. 2008), as integrated in Gephi software. Since the algorithm "kicks off" from a random layout of nodes, the outcome is not deterministic. Consecutive runs may result in different outcomes, in terms of the number of modules or in terms of the values of the modularity statistic. For all networks computed in the study, an arbitrary number of 50 iterations of the algorithm was employed and the most prevalent module layout, with the highest modularity value, was finally adopted. In the end, for each spatial scale: (a) we produced the predominant assemblage types, to assist with functional patterns exploration, and (b) we acquired a classification for each trait, according to assemblage types, to gain insights on the functional structure and trait shifts across different scales.

Functional interaction with the environmental gradient
Redundancy analysis was employed to reveal the relation between the five key environmental descriptors and each trait and trait assemblage type (see CANOCO 4.5 for Windows; Ter Braak and Smilauer 2002). Ordination analysis was applied for all sites (i.e., the six ecoregions scale data), together with the corresponding environmental data (matrices: 427 samples × 5 env. descriptors and 427 samples × 24 traits ). The gradient length of the first axis in detrended correspondence analysis, used to specify the unimodal/linear response of the data, was estimated at 2.1 (< 3). This indicates a linear response; therefore, suggesting RDA as the appropriate multivariate regression analysis. Prior to RDA, the Monte Carlo permutation test (999 permutations) was computed, to specify the statistically significant environmental descriptors that take part in the model (p < 0.05).
To avoid intercorrelated environmental variables which have no unique contribution to the model, only those variables with inflation factor < 20 were used (Ter Braak and Smilauer 2002). In this case, all environmental variables were incorporated in the analysis since the inflation factor was below 2 for each of them. All data were log(x + 1) transformed prior to analysis. To further explore the association of the five environmental descriptors in relation to the trait assemblage types revealed from network analysis, five boxplots were constructed. These box-plots display the values of the environmental descriptors for each trait assemblage type category. Welch's independent samples t-test for significance was also performed for all descriptors, among the box-plot trait categories. All box-plots and t-test analyses were computed using R software (R Core Team 2021).
Supplementary statistics for each species (i.e., absolute abundance, mean relative abundance, frequency of occurrence across samples) are provided in the Appendix (Table 3), together with box-plot graphs for the 24 traits analyzed (Fig. 6).

Functional patterns of fish assemblage
For more than 80% of the iterations in all three scales, bipartite network analysis identified two major trait assemblage types (Fig. 2a, b, c). These two functional modules were named as upland type (pink color) and lowland type (light blue color) across all scales, since each type largely reflected the same site localities across all scales. More specifically, in each scale, the lowland type generally comprised sites located in the lower reaches of the rivers (low altitudes), while the upland type incorporated sites in the upper and middle river stretches (higher altitudes).
Since the dataset included only reference sites, most of them were clustered under the upland type, in all three spatial scales (Fig. 2). Typically, this is the type where anthropogenic disturbance is less present. Particularly, in the six ecoregions scale (Fig. 2a), 78% of the sites (n = 198) were classified into the upland, while 22% (n = 56) were grouped into the lowland type. In the one ecoregion scale (Fig. 2b), these site groupings were 79% (n = 180) and 21% (n = 64), while in the river basin scale (Fig. 2c), they were 67% (n = 33) and 33% (n = 16) respectively. Furthermore, for all networks, no strict boundary seems to exist between sites belonging to the lowland and those to the upland type (Fig. 2). Therefore, nodes belonging to, e.g., the upland type may in fact place themselves closer to nodes belonging to the lowland type, rather than to nodes of the same type.
Functional structure and trait shift across spatial scales Within each spatial scale, most of the traits were classified into the lowland type rather than the upland type (Table 2), suggesting that functional richness increases from the upper to the lower river reaches. Furthermore, in the case of the two ecoregional scales analyzed (i.e., six ecoregions and one ecoregion), all distinct traits known to inhabit the upper reaches of the rivers (e.g., INSV, RHEO, Hab. INTOL) were "correctly" classified into the upland type. On the other hand, all traits considered to be present to the lower reaches (e.g., OMNI, EYRY, PHYT, Hab.TOL) were also "properly" arranged into the lowland type. However, in the case of the river basin scale, this concordance was not as evident, since rheophilic (RHEO), rheolithophilic (RH. LITH), and water column (WC) traits were arranged contrariwise to what would be expected.
Between spatial scales, a conformity of classification is apparent between six ecoregions and one ecoregion scale (Table 2). In short, when shifting from the six ecoregions towards the one ecoregion scale, all traits involved are grouped in the same assemblage type, with no exception. Yet, from the one ecoregion scale towards the river basin scale, three distinct traits (BENTH, RH.LITH, RHEO) shift from the upland to the lowland type, and one trait (WC) shifts in reverse.
Functional patterns along the environmental gradient Multivariate RDA analysis was employed for all sites (i.e., the six ecoregions scale data), setting the environmental descriptors as explanatory variables and trait values as response variables (Fig. 3) axes) explained 83.6% (F-ratio = 88.114; p-value = 0.001), which is most of the variation in the model. The second axes explained 92.2%, and all four axes explained 99.8% (F-ratio = 21.989; p-value = 0.001) that is almost all the variation of the data. The main findings of the multivariate ordination analysis were as follows: (a) the theoretical environmental gradient comprises two opposing components, mainly along the first theoretical axes (x-axes). The "downstream" component (on the left), which involves high values of January temperature and high values of catchment area, and the "upstream" component (on the right) exhibiting high values of altitude, slope, and source-altitude; (b) except for HERB, all other traits that were classified into the upland type according to network analysis, arrange along the "upstream" component. On the other hand, all traits belonging to the lowland type align with the "downstream" component; and (c) several distinct traits exhibit quite a strong correlation (large vector size) with the environmental gradient. For example, the proportion of intolerant (O2.INTOL, Hab.INTOL) and insectivorous (INSV) species undoubtedly increases along the "upstream" gradient. In a smaller degree, this involves lithophilic (LITH), rheolithophilic (RH.LITH), and rheophilic (RHEO) species as well. Similarly, it is evident that the proportion of tolerant (O2.MTOL, Hab.TOL) and omnivorous (OMNI) species increases along the "downstream" gradient. Once more, to a smaller degree, this is also true for phytophilic (PHYT), phytophils (Repro. PHYT), and eurytopic (EYRY) species. In contrast with all cases mentioned above, several traits belonging to both assemblage types (e.g., BENTH, HERB, PELA, LIMNO, OSTRA.VIVI) show weaker correlations (small vector size) with the environmental gradient.
The distribution of the sites according to lowland and upland assemblage type classification is illustrated in Fig. 4. The box-plots for the five key environmental descriptors according to lowland and upland type, as well as the t-test results for significance between the two types, are also shown. Corroborating the findings of the RDA, upland sites exhibit high values of altitude, slope, source-altitude, and low values of catchment area and temperature of January. These sites primarily distribute in the upper courses or in the middle course of the rivers. For lowland sites, environmental descriptors act exactly in the reverse direction; these sites are normally found in the lower courses of the rivers, or they belong to the typical xerothermic Mediterranean stream type. All t-tests exhibited significance (at least p < 0.05), suggesting once more that the RDA environmental components are in fact distinct, and that they broadly incorporate the two different trait assemblage types. Upland type arranges in the "upstream," while lowland type places itself along the "downstream" gradient. Supplementary statistics of the environmental descriptors (i.e., standard deviation, mean, median, range) are also provided in the Appendix (Appendix Table 4).

Functional assemblage structure and shifts
In our dataset, network analysis depicted two major trait assemblage types, in each spatial scale. This functional designation of fish communities along the longitudinal profile of the rivers comprised of two distinct broad functional zones (Oliveira et al. 2012) or functional syndromes (Logez et al. 2013). These two functional zones, i.e., upland and lowland, largely incorporate different quantitative and qualitative functional characteristics.
Within all three scales, trait richness was evidently higher in the lowland compared to the upland type. Traits in species known to inhabit the upper parts of the rivers (e.g., insectivorous, benthic, rheophilic, potamodromous, intolerant species) were largely classified into the upland type, while traits considered to dominate the lower reaches (e.g., omnivorous, phytophilic, limnophilic, and tolerant species) were arranged into the lowland type. These findings are consistent with the typical pattern of the functional richness increase that takes place from the upstream to the downstream river reaches (Penáz and Jurajda 1993;Aarts and Nienhuis 2003). They are also in accordance with the expectations provided by general theoretical concepts, such as the River Continuum Concept (Vannote et al. 1980). The only important exception to the above was in the case of the river basin scale, where benthic and habitat rheophilic species arrange into the lowland type. We deem that this inconsistency is due to the spatial scale "resolution" of this analysis (see Da Silva et al. 2019). To elaborate, we postulate that functional structure tends to be more sophisticated in smaller spatial resolutions; therefore, it does not always follow theoretical conceptions. This latter assumption has also been confirmed by Oliveira et al. (2012). The authors of this review propose multi-scale approaches in order to fully assess the factors that govern the functional organization of biotic assemblages in the Mediterranean-climate areas.
Finally, our trait-based analysis indicates that no strict boundary seems to exist between lowland and upland types (low modularity values; see Fig. 2). Along the longitudinal profile of the rivers, a progressive spatial shift between the upland and the lowland type takes place. This result is similar to other broad scale studies (Aarts and Nienhuis 2003;Erős et al. 2017) and in concordance with the concept that natural variability within assemblage types can sometimes be greater than the variability between them (Erős et al. 2017).  Our results confirm that environmental changes are mirrored in the functional structure of fish assemblages in the study area. The environmental gradient involved two contra-wise major components: the "upstream" and the "downstream" component. The former includes sites with high values of altitude, source-altitude, and slope, while the latter comprises sites that exhibit high values of mean January temperature and high values of catchment area.
Notably, our network analysis also highlights two major trait assemblage types; the upland type which similarly to the "upstream" environmental component involves sites in the upper river reaches and the lowland type that mainly includes sites in the lower river stretches, analogously to the "downstream" environmental component. In functional terms, considering only traits with high  Fig. 4 Distribution of the 254 study sites according to upland and lowland classification type (via network analysis). Boxplots of the five key environmental descriptor values (altitude, slope, catchment area, source altitude, and temperature of January) across upland and lowland type. Box-plot outliers have been masked for better illustration; box represents the interquartile range, bold line the median, and X the mean. T-test for significance, for each descriptor values between lowland (n = 106) and upland (n = 321), is also provided loadings with the environmental gradient, this condition could be described as, "the proportion of insectivorous, lithophilic, rheolithophilic, and intolerant (to habitat and oxygen degradation) species exhibit an affinity with the upstream river stretches, while omnivorous, phytophilic, eurytopic, and tolerant (to habitat and oxygen degradation) species dominate the downstream of rivers." These findings are consistent with similar studies conducted on a continental or on a regional level, across Europe. In a broad European scale study, Logez et al. (2013) describe an analogous functional structure along a similar environmental gradient; however, the East Mediterranean region was not included. Briefly, their analysis also uncovers two distinct gradient components (the catchment area-distance from source-slope gradient opposing the temperature gradient) and two functional syndromes (insectivorous-intolerant-rheophilic species opposing omnivorous-tolerant-eurytopic species). The former functional syndrome dominates the upper parts and the latter the lower parts of the European rivers. Once more, in a more regional scale study at three distinct European ecoregions (Spain, France, and Belgium), Logez et al. (2010) state that omnivorous and tolerant species are positively correlated with temperature and negatively correlated with slope and distance from source. Likewise, in a survey concerning three large European rivers (Doubs, Rhine, and Meuse), Aarts and Nienhuis (2003) show that the proportion of rheophilic and lithophilic species tend to increase in the upper zones of large rivers, while the proportion of limnophilic, phytophilic, and eurytopic species typically increase in the lower reaches. Finally, Oliveira et al. (2012) demonstrate a strong functional upstream-downstream pattern, within the fish assemblages in Western Iberian Peninsula.

Limitations and uncertainties
Because traits are essentially human constructs, assigning species to traits can be problematic (Aarts and Nienhuis 2003), due to arbitrary criteria resulting in reduced reliability (Frimpong and Angermeier 2010;Benoit et al. 2021). This is especially true when considering categorical traits (Cano- Barbacil et al. 2020). Alongside, trait characterization commonly ignores intraspecific variability (Violle et al. 2012) and focuses primarily on adult strategies. As a result, traits do not fully capture the complexity of fish communities (Blondel 2003). The spatial scale of analysis (Oliveira et al. 2012;Logez et al. 2013) and the limited ecological information about the Eastern Mediterranean endemic species (Welcomme et al. 2006;Noble et al. 2007;Benejam et al. 2015;Cano-Barbacil et al. 2020) further complicate trait-based ecological approaches.
Reliability gaps also take place when analyzing functional data. For example, when a considerable amount of either not classified, incorrectly classified  or skewed trait data (Zhang et al. 2021) are introduced in the analysis. For example, eurytopic and phytophilic species are underrepresented in this study. These two traits naturally occur at the lower reaches and the floodplains, but both localities are negatively skewed in our data. All these adversities tend to over-simplify the functional structure along the longitudinal gradient of rivers, subsequently affecting all relevant bioassessment analysis.

Implications for ecological applications
The functional structure of fish assemblages shows a relative stability in space (Logez et al. 2013;Cano-Barbacil et al. 2020) and in time (Aarts and Nienhuis 2003). Apparently, this provides a unique feature for functional-based analyses in applied frameworks, e.g., when computing broad modelbased bioassessment indices using traits. Species functional groupings within assemblages have proven to be more effective than species-level approaches in fish-based bioassessment in Europe (Pont et al. 2006. In addition, the use of networks to disentangle complex ecological structures is undisputed (Bellard et al. 2017). Networks are scientifically informative, visually appealing, and may trigger stakeholder engagement (Pocock et al. 2016). Therefore, the application of networks should be encouraged in the functional exploration of freshwater fish communities (Koutsikos et al. 2021).
Despite recent advances, there is still some debate on functional approaches. Several studies criticize the use of specific traits, e.g., Frimpong and Angermeier (2010) disapprove of the use of tolerance concept, as it is not well defined and consequently it is prone to inconsistency. Likewise, feeding traits due to the dietary flexibility of species may not be appropriate correlates of the environment (Pusey et al. 1995;Welcomme et al. 2006). However, tolerance and feeding traits were very strongly correlated with the environmental gradient in this study (e.g.,INSV,OMNI,Hab.TOL,O2.INTOL,etc.). We deem, at this point, that these distinct traits operate just as major components of the two functional syndromes depicted in this study (i.e., lowland and upland assemblage types), rather than functioning as discrete traits. We, therefore, support the multi-trait approach when an integrative understanding of ecosystem functions is required (Raffard et al. 2017), even if specific functional categories or certain traits seem to be redundant or inappropriate.
As regards the Eastern Mediterranean region, it is expected that several fish species have developed specific functional adaptations to compensate for the harsh climate of summer droughts and the consequent temporal food scarcity (Hershkovitz and Gasith 2013;Benejam et al. 2015;Vardakas et al. 2015). For example, fish may become feeding opportunists (i.e., generalists) and grow greater tolerance to environmental variations (Magalhaes et al. 2002;Ferreira et al. 2007). In this analysis, all intermittent sites (see Skoulikidis et al. 2017) and regardless of their environmental characteristics (altitude, slope, catchment area, etc.) cluster within the lowland assemblage type. This type essentially involves all relevant functional adaptations, e.g., OMNI, EURY, Hab.TOL, and O2.TOL. We subsequently confirm the primary hypothesis that this functional "fine-tuning" of fish species takes place in the Mediterranean freshwater ichthyofauna. We consider that this applies unquestionably in species living in the most typical of the Mediterranean ecotype spectrum; the intermittent rivers and streams, including those that exhibit partially intermitted and ephemeral flow reaches.
This study produced three discrete trait sets: common traits, restricted traits, and data variant traits. Common traits are widespread traits shared by many species that represent major river ecosystem processes for the study area. Restricted traits are shared by a few species that exhibit limited occurrence across samples and reflect ecosystem processes constrained to specific river habitats or due to geographic and/or historic factors. However, data variant traits depend on the dataset used for the analysis and in this case mainly represent ecological functions that would normally be recorded in the lentic/downstream habitats of the rivers. However, the global scarcity of unimpacted sites in downstream river reaches (Benejam et al. 2015) remains an overwhelming difficulty, which can produce misleading benchmarks in ecological or bioassessment applications and should be considered during research.
Lastly, the functional stability across the two ecoregional scales is undisputed; same traits involved in each assemblage type. On the other hand, certain traits shift between assemblage types in the river basin scale. Therefore, it is essential to carefully consider research aims and priorities before deciding on the scale of analysis. For instance, broad trait-based bioassessment indices can effectively adopt a multiecoregional or even ecoregional scale approach, while trait-based biodiversity protection and restoration should normally implicate a narrower scale of analysis, such as the river basin scale or even smaller.

Conclusions
Assigning fish species into functional groups is not a straightforward procedure and requires careful review. In the process of exploring the community ecology of stream fishes, this step provides new and revised matrices of functional traits. This can be of particular importance in effectively applying monitoring and bioassessment applications in understudied areas. By using network analysis, we portray patterns at three commonly used spatial scales and identify two discrete fish functional zones (upland and lowland) among different biogeographical ecoregions. These functional zones are not well defined within strict spatial boundaries, but they are clearly incorporated along the upstream-downstream longitudinal gradient tradition, which is prominent in temperate European river studies. Insights gained from the analysis of a large dataset, using the best available least-impacted sites, have helped develop an understanding for trait-based applied ecology in this part of the Eastern Mediterranean.
Ethics approval Fish community monitoring was carried out with the official endorsement of the Greek Ministry of Rural Development & Food as well as the Greek Ministry of Environment, Energy, and Climate Change (permits: 4AM00-KHI, BIHY0-718, and 9ZE24653P8-ΖΟ6). Fish sampling was carried out by trained personnel from HCMR following the established protocols outlined in the Inland Waters Fish Monitoring Operations Manual (IMBRIW-HCMR 2013). These guidelines were developed in-house to address concerns related to conservation and animal welfare during electricity-based sampling procedures. The research was done according to ethical standards.

Competing interests
The authors declare no competing interests.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/. Figure 5 Table 3  Figure 6  Table 4 Fig. 5 Key operations of data matrices to uploading Sample-trait matrix in Gephi  Barbieri et al. (2015) Fig . 6 Box-plots for the 24 traits analyzed (relative frequency of a trait per functional category). Box represents the interquartile range, bold line represents the median, X represents the mean, and circles represent outliers. Abbreviations in Table 1