Environmental Gradients Shaping the Freshwater Bryophyte Communities of Croatia (Western Balkans)

A comprehensive field survey of 527 sites on 293 watercourses across Croatia revealed 76 sites (14.42%) in which bryophytes were the dominant part of the macrophyte vegetation. Using classification and ordination analyses, we obtained five community types segregated across the gradients of several climatic, physiographic and water chemistry parameters. The Didymodon tophaceus–Apopellia endiviifolia and the Berula erecta-Cratoneuron filicinum communities were mostly confined to the clean and basic karstic rivers of the Dinaric Ecoregion under the influence of the Mediterranean climate, with the Didymodon tophaceus–Apopellia endiviifolia community being a tufa-forming community associated with the seasonally dry watercourses of small catchment areas and cascades along the larger karstic rivers, while the Berula erecta–Cratoneuton filicinum community was mostly associated with rivers with larger catchment areas and permanent flow. On the other hand, the Oxyrrhynchium hians–Chiloscyphus pallescens community and the Fissidens pusillus–Veronica beccabunga community were associated with eutrophic water restricted to small rivers of the Pannonian Ecoregion under the influence of the temperate climate and flowing over silicate bedrock. The most represented and widespread in Croatia was the Cinclidotus community, displaying the widest ecological range in the study. It was mostly associated with the relatively clean karstic rivers of large catchment areas belonging to the Dinaric Ecoregion, with the majority of the sites under the influence of a temperate climate with higher precipitation during the warm period of the year. The geographical patterns of the freshwater bryophyte communities showed that the relatively clean, fast and cold karstic rivers belonging to the Dinaric Ecoregion provide habitats that harbour a greater diversity of bryophyte communities than the watercourses of the Pannonian Ecoregion, where bryophyte-dominated communities are restricted to a small number of small lowland and semi-montane rivers and predominantly occupy periodically flooded microhabitats such as river margins.


Introduction
Bryophytes are an important part of freshwater biodiversity, inhabiting a wide variety of aquatic and riparian habitats and ecological and hydrological niches associated with running and standing waters [1,2]. They are a dominant part of macrophyte vegetation in headwater and mountain streams, where they thrive in oligotrophic, clear, cold water with fast and usually torrential flow over very stable rocky substrates in a harsh environment unsuitable for the majority of other macrophytes [3,4]. Here, a wide variety of adaptations enables them to withstand mechanical stress from high water velocity and associated drag forces and shear stress [5]. The same features enable their complete dominance in waterfalls, where they develop the most luxuriant communities [6]. Namely, in such turbulent habitats, the thin boundary layer positively affects the amount of CO 2 available for photosynthesis, as well as nutrients. On the other hand, in middle and lower river sections dominated by the geographical segregation of the bryophyte communities follows this division and that those communities would show some degree of affinity towards a particular region.
Given that systematic and comprehensive studies on freshwater bryophyte communities have never been conducted in Croatia and do not exist with respect to the Western Balkans, including the corresponding part of the Mediterranean, we aimed to (1) explore the distribution, diversity and species composition of freshwater bryophyte communities; (2) explore the inherent variability of particular ecoregions in terms of bryophyte communities; and (3) identify the environmental gradients that influence bryophyte communities. This will improve the knowledge on this subject on the European level and provide a basis for further monitoring and protection, including the mitigation of the negative impacts of climate change and associated changes in hydrological regimes, as well as of humaninduced eutrophication and changes in both hydrological regimes and morphology of streams and rivers.

Results
Bryophytes were the dominant component of macrophyte vegetation in 76 sites out of the 527 (14.42%) surveyed sites on streams and rivers situated across the whole Croatian territory and encompassing the heterogeneity of Croatian watercourses in terms of the recent typology of the waterbodies developed for WFD implementation. The majority of the sites (61, to be precise) were situated in the Dinaric Ecoregion, accounting for 31.12% of a total 196 sites surveyed within the particular region. The Dinaric-Continental Subecoregion was the richest, with as many as 40.23% of sites with bryophyte communities out of a total 87 sites surveyed, while 23.85% out of 109 sites in the Dinaric-Mediterranean Subecoregion harboured macrophyte vegetation with bryophyte predominance. On the other hand, this proportion was comparatively low in the Pannonian Ecoregion, amounting to only 4.53%.

Community Groups
TWINSPAN classification of bryophyte-dominated sampling sites at maximal distance established five groups after three levels of division (Soerensen dissimilarity, max distance = 0. 77). An ANOSIM test confirmed the overall significant difference among the TWINSPAN groups (coded hereafter as 1-5), i.e., the existence of discrete communities among the sampling sites (overall R = 0.50, p(same) < 0.0001) based on the species composition. Furthermore, ANOSIM pairwise comparisons showed that all community groups were significantly different (Table 1).
With respect to the distribution of 76 bryophyte-dominated sites within the different sub-and ecoregions, 80.26% were located in the Dinaric Ecoregion (61 sites), with 46.05% (35 sites) situated in the Continental Subecoregion and 34.21% (26 sites) in the Mediterranean Subecoregion. The Pannonian Ecoregion was comparatively poor, with only 19.74% (15 sites) out of 76 sites ( Figure 1, Table 2). Table 1. R statistics of the pairwise ANOSIMs of bryophyte-dominated communities obtained from TWINSPAN classification (groups 1-5) (overall R = 0.50, p(same) = 0.0001), * p < 0.001, ** p < 0.005. Some main patterns are recognisable when taking into account the distribution of particular TWINSPAN communities between the sub-and ecoregions; however, a certain overlap is present. Groups 1 and 3 are mainly confined to the Pannonian Ecoregion, while all others are more frequent in the Dinaric Ecoregion, with Group 2 being exclusive for this region and equally distributed in both the Continental and the Mediterranean subecoregions. On the other hand, Group 4 was more frequent in the Mediterranean Subecoregion, while Group 5 was more characteristic of the Continental Subecoregion (Figures 1 and 2; Table 2).  Some main patterns are recognisable when taking into account the distribution of particular TWINSPAN communities between the sub-and ecoregions; however, a certain overlap is present. Groups 1 and 3 are mainly confined to the Pannonian Ecoregion, while all others are more frequent in the Dinaric Ecoregion, with Group 2 being exclusive for this region and equally distributed in both the Continental and the Mediterranean subecoregions. On the other hand, Group 4 was more frequent in the Mediterranean Subecoregion, while Group 5 was more characteristic of the Continental Subecoregion (Figures 1 and 2; Table 2). The characteristic species of Group 1 (Oxyrrhynchium hians-Chiloscyphus pallescens community) (Appendix A) were mostly hygrophytic species confined to periodically submerged river margins, such as the mosses Oxyrrhynchium hians, Plagiomnium undulatum, Pohlia melanodon, and the liverworts Chiloscyphus pallescens, Pellia neesiana, Conocephalum salebrosum, along with the moss Dichodontium pellucidum, collected from periodically submerged rocks within the riverbeds. However, the constant species included rheophytes such as Rhynchostegium riparioides and Leptodictyum riparium, as well as an ampyphyte, Cratoneuron filicinum. The hygrophytes Didymodon tophaceus, Apopellia endiviifolia and Funaria hygrometrica, as well as the amphiphyte Eucladium verticillatum and the rheophyte Fissidens crassipes, were the characteristic species of Group 2 (Didymodon tophaceus-Apopellia endiviifolia community). Group 3 (Fissidens pusillus-Veronica beccabunga community) was characterized by bryophytes such as Brachythecium rutabulum, Fissidens pusillus and Oxyrrhynchium speciosum, as well as by the vascular helophytes Veronica beccabunga and Persicaria dubia (Appendix A). In general, Groups 1, 2 and 3 were characterized by a higher The characteristic species of Group 1 (Oxyrrhynchium hians-Chiloscyphus pallescens community) (Appendix A) were mostly hygrophytic species confined to periodically submerged river margins, such as the mosses Oxyrrhynchium hians, Plagiomnium undulatum, Pohlia melanodon, and the liverworts Chiloscyphus pallescens, Pellia neesiana, Conocephalum salebrosum, along with the moss Dichodontium pellucidum, collected from periodically submerged rocks within the riverbeds. However, the constant species included rheophytes such as Rhynchostegium riparioides and Leptodictyum riparium, as well as an ampyphyte, Cratoneuron filicinum. The hygrophytes Didymodon tophaceus, Apopellia endiviifolia and Funaria hygromet-rica, as well as the amphiphyte Eucladium verticillatum and the rheophyte Fissidens crassipes, were the characteristic species of Group 2 (Didymodon tophaceus-Apopellia endiviifolia community). Group 3 (Fissidens pusillus-Veronica beccabunga community) was characterized by bryophytes such as Brachythecium rutabulum, Fissidens pusillus and Oxyrrhynchium speciosum, as well as by the vascular helophytes Veronica beccabunga and Persicaria dubia (Appendix A). In general, Groups 1, 2 and 3 were characterized by a higher frequency of hygrophyte bryophyte species (1-65.9%, 2-44.1%, 3-56.3%), followed by rheophytes (1-18.7%, 2-37.0%, 3-31.3%) (Figure 3).
The characteristic species of Group 4 (Berula erecta-Cratoneuron filicinum community) were the vascular helophytes Berula erecta, Menta aquatica and Sparganium erectum, while the rheophyte mosses Rhynchostegium riparioides and Fontinalis antipyretica and the amphyphyte Cratoneuron filicinum, along with Mentha aquatica, were constant species with high frequencies within the group (Appendix A). This group was characterized by the highest frequency of rheophytes in all the groups (57.1%), followed by hygrophytes (23.2%) and amphyphytes (19.6%), while mesophyte and xerophyte bryophytes were completely absent ( Figure 3). The characteristic species of Group 5 (Cinclidotus community) were Cinclidotus riparius and Cinclidotus aquaticus, both being constant species as well, along with Cinclidotus fontinaloides, Fontinalis antipyretica, Rhynchostegium riparioides and Cratoneuron filicinum. The rheophyte bryophytes were predominant in this group as well, accounting for 55.3% of all bryophyte occurrences within the group, followed by hygrophytes (32.3%) and amipyphytes (9.9%) ( Figure 3). Analysis of the TWINSPAN groups regarding macrophyte taxa richness revealed that sites belonging to the Group 1 had the highest mean value (12.25 ± 1.33), followed by Group 2 (10.56 ± 1.21) and Group 5 (10.24 ± 1.02). The same pattern was observed when considering the bryophyte species alone. Namely, the mean bryophyte species richness of sites belonging to Group 1 was highest (11.38 ± 1.18), again followed by Groups 2 (8.94 ± 0.91) and 5 (8.21 ± 0.83) ( Table 3). The share of the bryophyte species in total number of taxa was the highest in Group 1 (86%) and over 50% within Groups 5 (56.3%) and 2 (59.3%), while in Group 4, it amounted to 33.33%. On the other hand, this group harboured the overall highest number of vascular plant species (28 species) and had a lower mean taxa richness and mean bryophyte species richness (4.4 ± 0.46) than Groups 1, 4, and 5. Finally, Group 3 was the most taxa-poor group, with the lowest mean taxa richness and mean bryophyte species richness (Table 3).
When present, vascular plants were mostly represented by helophyte species. A comparison of species richness and Shannon-Wiener alpha diversity index of bryophytes and vascular plants between the groups revealed that the vascular plant alpha diversity was the highest in Group 4, followed by Group 3, reaching that of bryophytes in some localities within Group 4 ( Figure 4). By contrast, Groups 1 and 5 had a very low vascular plant alpha diversity; it was somewhat higher in Group 2, but still considerably lower than the bryophyte alpha diversity ( Figure 4). The characteristic species of Group 4 (Berula erecta-Cratoneuron filicinum community) were the vascular helophytes Berula erecta, Menta aquatica and Sparganium erectum, while the rheophyte mosses Rhynchostegium riparioides and Fontinalis antipyretica and the amphyphyte Cratoneuron filicinum, along with Mentha aquatica, were constant species with high frequencies within the group (Appendix A). This group was characterized by the highest frequency of rheophytes in all the groups (57.1%), followed by hygrophytes (23.2%) and amphyphytes (19.6%), while mesophyte and xerophyte bryophytes were completely absent ( Figure 3). The characteristic species of Group 5 (Cinclidotus community) were Cinclidotus riparius and Cinclidotus aquaticus, both being constant species as well, along with Cinclidotus fontinaloides, Fontinalis antipyretica, Rhynchostegium riparioides and Cratoneuron filicinum. The rheophyte bryophytes were predominant in this group as well, accounting for 55.3% of all bryophyte occurrences within the group, followed by hygrophytes (32.3%) and amipyphytes (9.9%) ( Figure 3).
Analysis of the TWINSPAN groups regarding macrophyte taxa richness revealed that sites belonging to the Group 1 had the highest mean value (12.25 ± 1.33), followed by Group 2 (10.56 ± 1.21) and Group 5 (10.24 ± 1.02). The same pattern was observed when considering the bryophyte species alone. Namely, the mean bryophyte species richness of sites belonging to Group 1 was highest (11.38 ± 1.18), again followed by Groups 2 (8.94 ± 0.91) and 5 (8.21 ± 0.83) ( Table 3). The share of the bryophyte species in total number of taxa was the highest in Group 1 (86%) and over 50% within Groups 5 (56.3%) and 2 (59.3%), while in Group 4, it amounted to 33.33%. On the other hand, this group harboured the overall highest number of vascular plant species (28 species) and had a lower mean taxa richness and mean bryophyte species richness (4.4 ± 0.46) than Groups 1, 4, and 5. Finally, Group 3 was the most taxa-poor group, with the lowest mean taxa richness and mean bryophyte species richness (Table 3). When present, vascular plants were mostly represented by helophyte species. A comparison of species richness and Shannon-Wiener alpha diversity index of bryophytes and vascular plants between the groups revealed that the vascular plant alpha diversity was the highest in Group 4, followed by Group 3, reaching that of bryophytes in some localities within Group 4 ( Figure 4). By contrast, Groups 1 and 5 had a very low vascular plant alpha diversity; it was somewhat higher in Group 2, but still considerably lower than the bryophyte alpha diversity (Figure 4).
The chorological comparison of TWINSPAN groups based on major biomes revealed large chorotype overlapping, with a dominance of temperate chorotypes; however, some biogeographical differences were highlighted ( Figure 5). The southern-temperate chorotype had the highest frequency in Group 2 (54.2%), while its lowest frequency was in Group 1 (18.7%). Furthermore, Group 4 was characterized by a higher frequency of boreo-temperate (36.5%) and wide-boreal elements (19.0%) than other groups. The Mediterranean-Atlantic chorotype was completely absent from Groups 3 and 4, with the highest, still quite low in proportion (2.2%) in Group 5. Boreal-montane and boreo-arctic montane chorotypes were most represented in Group 1, with 2.2% and 5.5% respectively.   The chorological comparison of TWINSPAN groups based on major biomes revealed large chorotype overlapping, with a dominance of temperate chorotypes; however, some biogeographical differences were highlighted ( Figure 5). The southern-temperate chorotype had the highest frequency in Group 2 (54.2%), while its lowest frequency was in Group 1 (18.7%). Furthermore, Group 4 was characterized by a higher frequency of boreotemperate (36.5%) and wide-boreal elements (19.0%) than other groups. The Mediterranean-Atlantic chorotype was completely absent from Groups 3 and 4, with the highest, still quite low in proportion (2.2%) in Group 5. Boreal-montane and boreo-arctic montane chorotypes were most represented in Group 1, with 2.2% and 5.5% respectively. The chorological comparison based on the eastern limit showed the dominance of circumpolar and European chorotypes in all communities, while other chorotypes were absent or present with frequencies lower than 5%, except for the Eurosiberian chorotype, which accounted for 7.7% of Group 5.
Bryophyte lifeforms were not evenly distributed among the TWINSPAN groups, with the most conspicuous difference in the share of aquatic trailings, rough mats and turfs ( Figure 6). Namely, Groups 4 and 5 had the highest proportion of aquatic trailings, 46.0% and 33.9%, respectively, while this category was absent from Group 3 and represented with low frequency in Group 1 (6.6%). On the other hand, these latter two groups had a higher frequency of rough mats than the other groups, 37.5% for Group 3 and 24.2% for Group 1. Furthermore, Group 1 was characterized by a high frequency of turfs (39.56%), similar to Group 2, where turfs accounted for 30.5%. The chorological comparison based on the eastern limit showed the dominance of circumpolar and European chorotypes in all communities, while other chorotypes were absent or present with frequencies lower than 5%, except for the Eurosiberian chorotype, which accounted for 7.7% of Group 5.
Bryophyte lifeforms were not evenly distributed among the TWINSPAN groups, with the most conspicuous difference in the share of aquatic trailings, rough mats and turfs ( Figure 6). Namely, Groups 4 and 5 had the highest proportion of aquatic trailings, 46.0% and 33.9%, respectively, while this category was absent from Group 3 and represented with low frequency in Group 1 (6.6%). On the other hand, these latter two groups had a higher frequency of rough mats than the other groups, 37.5% for Group 3 and 24.2% for Group 1. Furthermore, Group 1 was characterized by a high frequency of turfs (39.56%), similar to Group 2, where turfs accounted for 30.5%. Regarding the life strategies, all TWINSPAN groups feature the overall dominance of perennials, followed by colonist bryophyte species (Figure 7). The share of perennials was lowest in Group 2 (19.1%) and highest in Group 3 (37.5%), followed by Group 4 (33.3%). A similar pattern was revealed when observing all perennial categories together; they were most represented in Group 4 (68.3%), followed by Group 3 (56.3%) and Group 5 (50.9%). The share of colonists within Groups 2, 4 and 5 was over 35% and was lowest within Group 3 (12.5%). Taking into account all three colonist categories, the highest proportion was recorded in Group 2 (41.2%), followed by Groups 5 (40.9%) and 1 (40.7%). Life strategy spectrum of freshwater bryophytes for TWINSPAN groups and a total sample of 76 bryophyte dominated sites (p-perennials, pc-competitive perennials, ps-stress-tolerant perennials, c-colonist, cp-pioneer colonist, ce-ephemeral colonist, a-annual shuttle, s-shortlived shuttle, l-long-lived shuttle, f-fugitives).

Environmental Gradients
DCA analysis showed the separation of discrete groups with some overlapping. The axis 1 eigenvalue was 0.34, and for axis 2, it was 0.29. The lengths of axes 1 and 2 were 3.8 and 2.9, respectively. The nature of the established gradients in the DCA analysis was further assessed with weighted Ellenberg indicator values, passively projected over the Regarding the life strategies, all TWINSPAN groups feature the overall dominance of perennials, followed by colonist bryophyte species (Figure 7). The share of perennials was lowest in Group 2 (19.1%) and highest in Group 3 (37.5%), followed by Group 4 (33.3%). A similar pattern was revealed when observing all perennial categories together; they were most represented in Group 4 (68.3%), followed by Group 3 (56.3%) and Group 5 (50.9%). The share of colonists within Groups 2, 4 and 5 was over 35% and was lowest within Group 3 (12.5%). Taking into account all three colonist categories, the highest proportion was recorded in Group 2 (41.2%), followed by Groups 5 (40.9%) and 1 (40.7%). Regarding the life strategies, all TWINSPAN groups feature the overall dominance of perennials, followed by colonist bryophyte species (Figure 7). The share of perennials was lowest in Group 2 (19.1%) and highest in Group 3 (37.5%), followed by Group 4 (33.3%). A similar pattern was revealed when observing all perennial categories together; they were most represented in Group 4 (68.3%), followed by Group 3 (56.3%) and Group 5 (50.9%). The share of colonists within Groups 2, 4 and 5 was over 35% and was lowest within Group 3 (12.5%). Taking into account all three colonist categories, the highest proportion was recorded in Group 2 (41.2%), followed by Groups 5 (40.9%) and 1 (40.7%). Life strategy spectrum of freshwater bryophytes for TWINSPAN groups and a total sample of 76 bryophyte dominated sites (p-perennials, pc-competitive perennials, ps-stress-tolerant perennials, c-colonist, cp-pioneer colonist, ce-ephemeral colonist, a-annual shuttle, s-shortlived shuttle, l-long-lived shuttle, f-fugitives).

Environmental Gradients
DCA analysis showed the separation of discrete groups with some overlapping. The axis 1 eigenvalue was 0.34, and for axis 2, it was 0.29. The lengths of axes 1 and 2 were 3.8 and 2.9, respectively. The nature of the established gradients in the DCA analysis was further assessed with weighted Ellenberg indicator values, passively projected over the . Life strategy spectrum of freshwater bryophytes for TWINSPAN groups and a total sample of 76 bryophyte dominated sites (p-perennials, pc-competitive perennials, ps-stress-tolerant perennials, c-colonist, cp-pioneer colonist, ce-ephemeral colonist, a-annual shuttle, s-shortlived shuttle, l-long-lived shuttle, f-fugitives).

Environmental Gradients
DCA analysis showed the separation of discrete groups with some overlapping. The axis 1 eigenvalue was 0.34, and for axis 2, it was 0.29. The lengths of axes 1 and 2 were 3.8 and 2.9, respectively. The nature of the established gradients in the DCA analysis was further assessed with weighted Ellenberg indicator values, passively projected over the ordination as vectors. This revealed a gradient from the sites with higher mean indicator values for temperature, light and moisture and low values for continentality (Group 4) compared to those with higher values for continentality but lower temperature and light values (Group 1), with the sites belonging to Group 5 being intermediate across this gradient ( Figure 8). Group 2 included the sites with higher Ellenberg indicator values for reaction. DCA axis 2 was the most strongly correlated with nutrient content, indicating the higher values in Group 3, as well as in some sites of Group 5.  The first axis of the CCA explained 23.65% and the second axis explained 42.39% of the variation in the relationship between vegetation data and environmental factors. Eigenvalues of the first and second axes equalled 0.35 and 0.27, respectively. Ordination was statistically significant (F = 1.41, p = 0.001) according to the Monte Carlo permutation test (999 permutations). As expected, variables that explain the majority of the variance in the data were also highly correlated with axis 1. An ordination plot of the CCA analysis revealed a strong gradient along axis 1 from sites with high values of orthophosphates and biochemical oxygen demand, and elevated total nitrogen values (groups 1 and 3) to the sites with low values of these water chemistry parameters and high values of total alkalinity, as well as dissolved oxygen (Groups 2 and 4) ( Table 4, Figures 9 and 10). However, the Mann-Whitney pairwise test showed that the latter two groups differed significantly in pH values, with Group 2 being associated with more basic water (Appendixes B and C). Sites belonging to Group 5 were distributed along the longest part of the CCA gradient, suggesting, in general, intermediate values for the abovementioned parameters compared to other groups. Total suspended solids (TSS) showed a similar pattern among the groups ( Figure 10, Appendix B), and the Mann-Whitney pairwise test indicated a significant difference between Groups 1 and 3 with high values and the other three groups with low TSS (Appendix C). The first axis of the CCA explained 23.65% and the second axis explained 42.39% of the variation in the relationship between vegetation data and environmental factors. Eigenvalues of the first and second axes equalled 0.35 and 0.27, respectively. Ordination was statistically significant (F = 1.41, p = 0.001) according to the Monte Carlo permutation test (999 permutations). As expected, variables that explain the majority of the variance in the data were also highly correlated with axis 1. An ordination plot of the CCA analysis revealed a strong gradient along axis 1 from sites with high values of orthophosphates and biochemical oxygen demand, and elevated total nitrogen values (groups 1 and 3) to the sites with low values of these water chemistry parameters and high values of total alkalinity, as well as dissolved oxygen (Groups 2 and 4) ( Table 4, Figures 9 and 10). However, the Mann-Whitney pairwise test showed that the latter two groups differed significantly in pH values, with Group 2 being associated with more basic water (Appendices B and C). Sites belonging to Group 5 were distributed along the longest part of the CCA gradient, suggesting, in general, intermediate values for the abovementioned parameters compared to other groups. Total suspended solids (TSS) showed a similar pattern among the groups ( Figure 10, Appendix B), and the Mann-Whitney pairwise test indicated a significant difference between Groups 1 and 3 with high values and the other three groups with low TSS (Appendix C). Occurring mainly in small, semi-montane watercourses with small catchment areas and under the influence of a temperate climate; in waters with high nutrient levels, BOD and TSS, and lower dissolved oxygen concentrations. Species-poor community, characterized by a higher share of hygrophytes and rough mats in lifeform spectrum, which grow on periodically submerged substrates. Group 4-Berula erecta-Cratoneuron filicinum Community Characteristic species: Mentha aquatica, Berula erecta, Sparganium erectum Constant species: Cratoneuron filicinum, Rhynchostegium riparioides, Fontinalis antipyretica, Mentha aquatica Distribution: Mainly Dinaric Ecoregion, Mediterranean Subecoregion Ecology: Transitional community of karstic rivers with large catchment areas and permanent flow, where vascular species start to outcompete bryophytes. Occurring in clean water with low nutrient content and BOD values, where helophytes occupy the river margins and shallower water, while bryophytes are confined to the riverbed. Characterized by a high share of rheophytes and aquatic trailings in lifeform spectrum, but low overall bryophyte species richness. Group 5-Cinclidotus Community Characteristic species: Cinclidotus riparius, C. aquaticus Constant species: Cinclidotus riparius, C. aquaticus, C. fontinaloides, Rhynchostegium riparioides, Cratoneuron filicinum, Fontinalis antipyretica Distribution: Mainly Dinaric Ecoregion, i.e., its Continental Subecoregion Ecology: The most widespread community, with a wide ecological range; in general, in waters with intermediate values of the water quality parameters (orthophosphates, total nitrogen and biochemical oxygen demand) and climatic variables associated with precipitation and water availability. Mostly in permanent karstic rivers with large catchment areas flowing over carbonate bedrock in neutral to basic water. Species-rich community characterized by a high share of rheophyte species and aquatic trailing in lifeform spectrum.
Additionally, CCA indicated the importance of climatic variables in the segregation of the investigated sites along axis 1 as well. This was especially prominent for precipitation in the coldest quarter (bio19) with the highest values associated with sites belonging to Group 2, followed by Group 4. This variable was highly positively correlated with the precipitation of the wettest month (bio13, r s = 0.97, p < 0.001) and the precipitation of the wettest quarter (bio16, r s = 0.94, p < 0.001), the values of which had the same pattern among the groups. Similarly, Groups 2 and 4 were associated with higher values of the mean temperature of the driest quarter (bio9) than the other groups. This variable was highly and significantly (p < 0.001) correlated with other climatic variables (bio5, r s = 0.89; bio6, r s = 0.96; bio10, r s = 0.91; bio11, r s = 0.98; bio14, r s = −0.73; bio 18 r s = −0,86) describing the hot and dry summer conditions characteristic of the Mediterranean climate associated with the sites aggregated on the right side of the CCA plot, in contrast to the sites on the left end of the CCA axis 1, influenced by more temperate climatic conditions (Table 4). On the other hand, higher values of the mean temperature of the wettest quarter (bio8) were associated with the sites on the left part of the CCA ordination plot.
Additionally, CCA indicated the importance of climatic variables in the segregation of the investigated sites along axis 1 as well. This was especially prominent for precipitation in the coldest quarter (bio19) with the highest values associated with sites belonging to Group 2, followed by Group 4. This variable was highly positively correlated with the precipitation of the wettest month (bio13, rs = 0.97, p < 0.001) and the precipitation of the wettest quarter (bio16, rs = 0.94, p < 0.001), the values of which had the same pattern among the groups. Similarly, Groups 2 and 4 were associated with higher values of the mean temperature of the driest quarter (bio9) than the other groups. This variable was highly and significantly (p < 0.001) correlated with other climatic variables (bio5, rs = 0.89; bio6, rs = 0.96; bio10, rs = 0.91; bio11, rs = 0.98; bio14, rs= −0.73; bio 18 rs= −0,86) describing the hot and dry summer conditions characteristic of the Mediterranean climate associated with the sites aggregated on the right side of the CCA plot, in contrast to the sites on the left end of the CCA axis 1, influenced by more temperate climatic conditions (Table 4). On the other hand, higher values of the mean temperature of the wettest quarter (bio8) were associated with the sites on the left part of the CCA ordination plot.
Regarding the physiographic variables, the catchment area was the most important in explaining the patterns in vegetation and environmental data, with groups 1 and 3 having generally smaller catchment areas (Table 4). This was corroborated by Mann-Whitney pairwise test; i.e., Group 3 was significantly different from Groups 2, 4 and 5, while Group 1 differed significantly from Groups 2 and 5 in the Mann-Whitney pairwise test (Appendix C).  Table 5; outliers: o-"out" values, *-"far out" or "extreme values").
Regarding the physiographic variables, the catchment area was the most important in explaining the patterns in vegetation and environmental data, with groups 1 and 3 having generally smaller catchment areas (Table 4). This was corroborated by Mann-Whitney pairwise test; i.e., Group 3 was significantly different from Groups 2, 4 and 5, while Group 1 differed significantly from Groups 2 and 5 in the Mann-Whitney pairwise test (Appendix C).

Discussion
The present study is the first comprehensive study dealing with the freshwater bryophyte communities and the environmental gradients underpinning their diversity, composition and distribution in Croatia and the Western Balkans, filling the gap in the existing knowledge on this subject at the European level.
While the presence and cover of bryophytes in freshwater habitats are primarily determined by riverbed stability, substrate size, stream slope and localized flow type [4,6,11], the diversity and community structure are governed by environmental variables operating on a larger scale. Geological, physiographic, and climatic factors, as well as water chemistry parameters, have been identified as essential drivers shaping freshwater bryophyte communities [4,6,8,[11][12][13]39].
Our findings confirm the importance of climatic, physiographic and water chemistry factors as major drivers influencing the diversity and composition of freshwater bryophyte communities, as well as their geographical segregation. Regarding water chemistry parameters, such as pH and alkalinity, which reflect the underlying geology, the distinction between hard-and soft-water bryoflora has been demonstrated by several authors on the European level and beyond [4,8,12,[39][40][41][42]. Water trophy level, i.e., water nutrient content, has also been recognized as an important factor influencing bryophyte cover and diversity, as well as community structure, with anthropogenically influenced eutrophication having a detrimental impact on freshwater communities [13][14][15][16]18,41,43,44]. Finally, climatic factors, especially those related to precipitation and water availability, as well as their distribution over a year, have been proven to regulate freshwater bryophyte communities, especially in highly seasonal Mediterranean rivers [8,39].
In our study, communities confined to karstic rivers with high alkalinity and pH values, reflecting the dominant carbonate bedrock, and clean and oxygenated water (Groups 2 and 4), were characterized by basophilous species. The Didymodon tophaceus-Apopellia endiviifolia community (Group 2) showed a stronger bryophyte dominance in macrophyte species composition, with tufa-forming species such as Didymodon tophaceus and Eucladium verticillatum being among the characteristic species of the community, along with basophilous Apopellia endiviifolia and Fissides crassipes. The cooccurrence of other basophilus species of oligotrophic water, such as Palustriella commutata and P. falcata and the liverwort Jungermannia atrovirens [45], was recorded within the community as well. Similar species assemblages have already been described for calcareous rivers in Europe, especially from the Mediterranean area, and were regarded as typical of neutral to basic clean water with low nutrient content in undisturbed flush flow fed streams with regular or low current conditions, as well as from cascades [46][47][48]. This community showed a prevalence of hygrophyte taxa in our study, and this was corroborated by a high frequency of turfs, species with vertical stems with little or no branching [49], within the lifeform spectrum. These are known to thrive in seasonally flooded habitats, with the strong impact of water [50] suggesting an interplay between flash flows and low water table periods in the Didymodon tophaceus-Apopellia endiviifolia community. Additionally, this community had the highest proportion of colonists in the life-strategy spectrum, indicative of a higher share of microhabitats flooded seasonally with strong discharge [50]. Regarding bryophyte composition, the Didymodon tophaceus-Apopellia endiviifolia community corresponds well with the community described by Vieira et al. [8] as a freshwater bryophyte community that, in Mediterranean Europe, has an extensive predicted presence, according to environmental niche modelling, in highly seasonal rivers in Spain, southern France, Italy and Greece [8]. The same research concluded that this particular community was characteristic of highly seasonal streams at low to moderate altitudes and high values of precipitation in the driest quarter that sustain permanent flows. However, our study revealed that the most important bioclimatic variable influencing both this community and the Berula erecta-Cratoneuron filicinum community (Group 4) within Croatia was the precipitation of the coldest quarter, a good surrogate for the hydrological patterns of Mediterranean rivers, as well as the mean temperature of the driest quarter. These communities were associated with higher values for these bioclimatic parameters, characteristic of a Mediterranean climate with dry and hot summers where higher precipitation occurs during warm winters [51], and were subsequently characterized by the high proportion of southern-temperate chorotypes. On the other hand, they were characterized by overall large catchment areas, i.e., hydrological watersheds. Namely, Groups 2 and 4 were recorded on karstic rivers of the Dinaric Ecoregion, which are most often a part of large and complex hydrographic networks, and receive water from numerous springs, with overground and subterranean courses supplying water from the Dinaric mountains of both Croatia and neighbouring Bosnia and Herzegovina [52]. However, the Didymodon tophaceus-Apopellia endiviifolia community (Group 2) included several sites on intermittent rivers with small catchment areas, while on the sites with large catchment areas, it was mostly confined to cascades. The Berula erecta-Cratoneuron filicinum community (Group 4) had higher mean catchment area values and significantly lower water pH than the Didymodon tophaceus-Apopellia endiviifolia community. Being the transitional community in which vascular species start to outcompete bryophytes, it harboured the highest vascular species number and alpha diversity in comparison to all other groups. Vascular helophytes such as Mentha aquatica, Berula erecta and Sparganium erectum which thrive in the shallow and slower water along the river margins, indicated the gradual transition of completely bryophyte dominated vegetation towards the herbland vegetation of small freshwater streams and the shallow waterbodies of temperate Europe belonging to the alliance of Glycerio-Sparganion Br.-Bl. et Sissingh in Boer 1942 [53]. However, with the greatest proportion of rheophyte bryophytes confined to the riverbed, this was the most truly aquatic bryophyte community within the study. This was supported by the analysis of the lifeform spectrum, which revealed the highest proportion of aquatic trailings and a moderate proportion of smooth mats, lifeforms best adapted to permanent submergence [50]. Aquatic trailings are mostly associated with slower currents, whereas smooth mats prefer the more torrential water zones [2,50], so their ratio in the Berula erecta-Cratoneuron filicinum community reflects a permanent, slower and more streamlined flow. Regarding the life-strategy spectrum, perennial species were the most represented, indicating constant ecological and hydrological conditions in these watercourses.
The Oxyrrhynchium hians-Chiloscyphus pallescens community (Group 1) was mostly restricted to small lowland rivers located in the Pannonian Ecoregion [54], with a quite low mean value for the catchment area. The sites belonging to this community were characterized by eutrophic and turbid water with low alkalinity. The latter is a result of the underlying geology, since the silicate bedrock is dominant within the Pannonian Ecoregion, while the predominant substrates in these localities are sandy and gravelly alluvial deposits of silicate origin. Furthermore, the higher values of the mean temperature of the wettest quarter were associated with this group, indicating more temperate or even continental climatic conditions present in the Pannonian Ecoregion. The overall high proportion of hygrophytes within this community corresponds with the prevalence of turfs and rough mats (creeping pleorocarpous species with lateral branches erect) [49] in the lifeform spectrum, which inhabit periodically submerged margins of river stretches flowing through forests of the Pannonian lowland. Regarding the rheophytes, the most frequent were Rhynchostegium riparioides and Laptodyctium riparium. While L. riparium is unambiguously recognized as a pollution-tolerant aquatic moss [55], with a preference for eutrophic waters [16,17,41,44], R. riparioides was omnipresent in our study, reaching the threshold set for constant species in all communities and being present along the entire gradient of nutrient concentration covered by this study, which suggests that the species is weakly linked to trophic conditions, highlighting its wide ecology range as previously reported by several authors [41,44].
The species-poorest community, Fissidens pusillus-Veronica beccabunga (Group 3), was the least represented in our study and mostly restricted to eutrophic and turbid small semi-montane watercourses within the Pannonian Ecoregion, with a single locality in the Dinaric-Continental Subecoregion on a small, lowland watercourse with a pebbly-gravelly substrate [54]. The vascular helophytes Veronica beccabunga and Persicaria dubia were among the constant and characteristic species of this community, which harboured the second highest vascular plant species richness following that of the Berula erecta-Cratoneuron filicinum community. The high share of hygrophytes, such as the characteristic species Brachythecium rutabulum, Oxyrrhynchium speciosum and Marchantia polymorpha, as well as a complete absence of aquatic trailings and a high share of rough mats (a lifeform with an adaptive advantage in microhabitats occurring above the normal level of maximum floods) [50] in the lifeform spectrum, confirmed that the periodically submerged rocks and alluvial sediments of river margins make the dominant microhabitat available to freshwater bryophytes in the watercourses of the Pannonian Ecoregion.
The Cinclidotus community (Group 5) was dominant in our study and displayed the widest ecological range regarding water quality and climatic variables associated with precipitation. While the majority of the sites were situated in the Dinaric-Continental Subecoregion, with a temperate climate with high values of precipitation in the warmest quarter, some of the sites were recorded in the source areas, as well as in the lower courses of Mediterranean rivers with permanent flow. Temperate chorotypes were dominant within the community, which was the case for all communities in the study, but here, they were represented by a high proportion of both southern-temperate and boreo-temperate ele-ments. Additionally, the presence of the boreal, wide-boreal and Mediterranean-Atlantic chorotypes distinguished this community from the Didymodon tophaceus-Apopellia endiviifolia community and the Berula erecta-Cratoneuron filicinum community, which were dominantly under Mediterranean influence. Furthermore, the Cinclidotus community was characterized by neutral to basic water, related to the dominant carbonate bedrock in the area of its distribution and the slightly lower mean values of alkalinity in the Didymodon tophaceus-Apopellia endiviifolia and Berula erecta-Cratoneuron filicinum communities. The characteristic species of this group, Cinclidotus aquaticus and C. riparius, have already been reported in situations of high alkalinity and are associated with carbonate bedrock [39,56,57]. Cinclidots aquaticus is characteristic of clean, cold, well-oxygenated waters with low nutrient content [17,39] and is a typical species in the fast water of source areas, permanent torrential watercourses, rapids and waterfalls [39,58]. On the other hand, C. riparius was usually found in more sheltered microhabitats, not directly exposed to water dragging forces, while C. fontinaloides, a constant species of the Cinclidotus community, was most often found on periodically submerged rocks or tree stumps. Among constant species, Fontinalis antypyretica was the most truly aquatic species, with the least desiccation tolerance, growing completely submerged and attached either to rocks or to logs in moving water. The occurrence of F. antipyretica has not been closely related to specific physicochemical or trophic conditions in the majority of the available studies [11,14,20,44], although an increase in its frequency was observed with the decreasing concentrations of nitrates and phosphates [17]. Additionally, the constant species Rhynchostegium riparioides and Cratoneuron filicinum, as well as frequent cooccurrence of Fissidens crassipess, Leptodyctium riparium and Apopellia endiviifolia, make this community quite close to the bryophyte community most commonly found in the Mediterranean and predicted to occur in the freshwater streams of its eastern part, as well as in northern Spain and France [8]. This community was regarded as having a high proportion of exclusively aquatic species characteristic of riverbeds and many boreal elements as compared to other communities identified for the Mediterranean in the particular study. Similarly, the Cinclidotus community in our study was a prominently aquatic community with a high proportion of rheophyte species, which was corroborated by the high proportion of aquatic trailing, as well as smooth mats.
The geographical patterns of the freshwater bryophyte communities in Croatia show that the Dinaric Ecoregion provides substantially more suitable habitats than the Pannonian Ecoregion, harbouring all five communities identified in this study, with different bryophyte communities recorded at as many as 31.12% of all surveyed sites. This was expected, since fast, relatively clean and cold karstic rivers are a prominent feature of this ecoregion, and freshwater bryophytes are known to thrive and are the dominant component of the macrophyte vegetation in conditions of fast and turbulent flow, rocky substrates and low temperatures, which vascular plants cannot withstand [3,4,6,46,53,59]. They are a prominent part of the vegetation in highly seasonal Mediterranean rivers as well, where they successfully cope with the interchange of dry periods and periods with flash flows of strong water currents due to their diverse morphological and physiological adaptations. In contrast, the Pannonian Ecoregion harbours a very small number of sites with bryophyte vegetation, with only 4.53% of all surveyed sites having this vegetation type. Watercourses in the Pannonian Ecoregion are mostly slow, eutrophic lowland streams and rivers with unstable sandy and gravelly alluvial sediments unsuitable for bryophytes, which are here additionally subjected to competition with vascular plants [5,7]. Furthermore, the majority of these watercourses are subjected to flow regulation through canalization, riverbed deepening and embankment, as well as considerable changes in land-use practices, including the removal of the riparian vegetation [60], all of which have a negative impact on aquatic vegetation in general.
Our study was limited by the predefined survey sites that were included in the assessment of ecological status, which, according to the WFD, includes waterbodies with catchment areas greater than 10 km 2 , while omitting the majority of the source areas and smaller headwater streams in which bryophyte communities are expected to occur.
With this in mind, future research should focus on these habitats, especially in parts of the Pannonian Ecoregion with mountain areas of high geological and geomorphological heterogeneity. However, we want to emphasise the importance of the WFD, which encouraged our research into freshwater bryophyte communities by including this group in its assessment of the ecological status of waterbodies. So far, substantial progress has been made [3,8,17,39], which is especially important in regions where bryophytes are still generally poorly researched [24], as in the case of Southeast Europe [13,16], with our findings contributing to knowledge with respect to the ecology and vegetation of bryophyte communities of the Mediterranean and providing the first insights into this subject for the Western Balkans.

Study Area
Data on the distribution of bryophyte-dominated freshwater communities was collected within the national surface water monitoring scheme conducted to assess the ecological status of waterbodies as required by the Water Framework Directive (WFD) [61]. The sampling sites were originally selected so as to encompass the heterogeneity of different waterbody types recognized by the recent typology developed as a basis for the monitoring of surface waters [54]. According to this typology, the land area of Croatia, 56,594 km 2 , is divided into two hydrological and biogeographical regions-the Pannonian and the Dinaric Ecoregion, the latter being subdivided into the Continental and Mediterranean Subecoregion ( Figure 11). In all, 293 watercourses were surveyed during the vegetation seasons from 2016 to 2021. The survey included as many as 527 sampling sites, ultimately covering the whole of Croatian territory ( Figure 11). The Pannonian Ecoregion refers to the alluvial and diluvial plains in the inland part of the country bounded by three large rivers (Sava, Drava and Danube). This area ranges between 80 and 135 m a.s.l., with a small number of rather low, solitary mountain massifs reaching 1000 m a.s.l. The lithological and geological composition is mostly silicate quaternary deposits, while limestone is present only locally, in higher mountain areas. The climate is temperate, with warm summers throughout most of the area (Cfb), hot summers predominately in the eastern part (Cfa) and no dry season [62]. On the contrary, the Dinaric Ecoregion is predominately built from Mesozoic limestone and dolomite bedrock and is characterized by karstic phenomena. This ecoregion includes the Dinarides, the largest uninterrupted karst landscape in Europe, occupying almost 50% of the territory of The Pannonian Ecoregion refers to the alluvial and diluvial plains in the inland part of the country bounded by three large rivers (Sava, Drava and Danube). This area ranges between 80 and 135 m a.s.l., with a small number of rather low, solitary mountain massifs reaching 1000 m a.s.l. The lithological and geological composition is mostly silicate quaternary deposits, while limestone is present only locally, in higher mountain areas. The climate is temperate, with warm summers throughout most of the area (Cfb), hot summers predominately in the eastern part (Cfa) and no dry season [62]. On the contrary, the Dinaric Ecoregion is predominately built from Mesozoic limestone and dolomite bedrock and is characterized by karstic phenomena. This ecoregion includes the Dinarides, the largest uninterrupted karst landscape in Europe, occupying almost 50% of the territory of Croatia. Because of the predominantly calcareous and dolomite bedrock, many rivers in the area have partly subterranean courses, or flow through impressive canyons or complex systems of barrage lakes, participating in the karst relief formation. The Continental Subecoregion is characterized by a temperate climate (Cfb), while the climate of the Mediterranean Subecoregion is mostly Mediterranean, i.e., temperate with dry and hot summer months (Csa) [62]. The Pannonian watercourses and the majority of the watercourses of the Dinaric-Continental Subecoregion belong to the Black Sea Basin, while the watercourses of the Dinaric-Mediterranean Subecoregion belong to the Adriatic Sea Basin.

Macrophyte Vegetation Sampling
A survey of macrophyte vegetation was performed according to the national methodology for macrophyte sampling [54] from June to September, when macrophyte vegetation is optimally developed, and during the lowest water discharge levels. Watercourses were surveyed for macrophytes along 100 m-long transects from the banks and, if the water depth was low enough, by zigzagging across the channel. In less accessible areas, the river bottom was raked to reach the macrophytes, with the rake either on a long pole or at the end of a rope.
Bryophytes were collected from various substrates (e.g., rocks, boulders, pebbles, xylal) within the riverbed, as well as from the periodically flooded river margins. Other macrophytes were collected for identification in the laboratory only if it was not possible to make accurate identification in the field. The collected material was deposited in herbarium ZA [66]. The nomenclature follows Hodgetts et. al. [67] for bryophytes, Euro + Med [68] for vascular plants and AlgaeBase [69] for algae.

Environmental Data Sampling and Collection
All localities were also sampled for basic water physicochemical and chemical parameters once a month throughout the year. Conductivity, pH, water temperature and dissolved oxygen were measured in situ with a Hach HQ40D Portable Multi Meter under standard conditions. Water samples were collected and analysed in an accredited laboratory (Central Water Management Laboratory, Zagreb, Croatia) for nine additional water chemistry parameters (Table 5).
Furthermore, physiographic and climatic environmental variables were obtained from several data sets using ArcGIS 10.5 software. Altitude was obtained from the digital elevation model EU-DEM v1.0 [70], distance from the source from topographic maps 1:25,000 [71], catchment area from the database of Hrvatske vode-the legal entity for water management and bioclimatic variables from CHELSA climatological datasets [72] (Table 5). Table 5. Environmental variables and abbreviations used.

Environmental Variable Abbreviation
Water physicochemical parameters

Data Analysis
Seventy-six sites in which bryophytes were the dominant component of macrophyte vegetation were selected for further analysis out of 527 surveyed sites. TWINSPAN analysis, a polythetic divisive classification method [73] modified by Roleček et al. [74], was conducted on vegetation relevés using Soerensen dissimilarity to assess whether discrete bryophyte communities occurred in any of the watercourses surveyed. TWINSPAN analysis was performed in Juice 7.1 [75]. The groups established by TWINSPAN were then tested for a significant difference based on species composition with the nonparametric test ANOSIM (using Bray-Curtis distance and 9999 permutations) [76] run in Past 4.9 software [77].
The groups resulting from these analyses were further analysed with respect to species composition. The synoptic table (Appendix A) was compiled in Juice 7.1. software. For all species within a particular group, e.g., community, ϕ-coefficients were calculated and tested for significance with the Fischer test. This statistical fidelity measure of particular species belonging to the previously defined groups was used to define characteristic species (ϕ ≥ 0.40, p < 0.05) [78]. Constant species were defined as those having the frequency within a particular group equal to or higher than 50%.
Bryophyte communities were analysed based on their affinity to water [10], chorotype composition [79], lifeform [49] and life strategies spectra [10] of bryophyte species. Furthermore, the alpha diversity (species richness and Simpson index) of the communities was analysed and visualized through boxplot graphs in SPSS 22.0 software.
Community structure was assessed using the indirect ordination method, Detrended correspondence analysis (DCA), run in R software (Vegan package) through Juice 7.1.-R connection. In DCA, weighted Ellenberg indicator values (EIV) for continentality, light, moisture, nutrients, temperature and pH reaction [49] were passively projected as vectors over the ordination to assess the possible environmental gradients. DCA revealed that the data were compositional with gradient longer than 3.0 SD units, indicating that a constrained analysis based on a unimodal model was most suitable for describing the data [80]. Subsequently, canonical correspondence analysis (CCA) was used as a direct ordination method to assess the relationship between environmental variables and patterns in species composition [81], while the statistical significance of the relationship between vegetation and environmental variables was tested by Monte Carlo permutation test using 999 permutations. Environmental data were log-transformed using the base-10 logarithm. Preselection of environmental variables was based on a correlation matrix between the variables and explorative DCA. Among highly positively correlated variables, those represented with the longest vectors in DCA were retained. This was performed in Past 4.9 software [77].
Additionally, a basic descriptive statistic (mean ± SE and min-max) of all environmental variables for groups derived from TWINSPAN analysis was calculated (Appendix B). A significant difference between groups was tested for each environmental variable with the nonparametric Kruskal-Wallace test, followed by Mann-Whitney pairwise post hoc tests (Appendix C). Nonparametric tests were selected since the majority of variables (32 out of 35) did not have a normal distribution, which was determined with the Shapiro-Wilk normality test in Past 4.9 [77].

Conclusions
The present study confirmed the importance of the climatic, physiographic and water chemistry gradients as major drivers shaping the diversity, composition and distribution of freshwater bryophyte communities. Comprehensive research revealed five community types and the patterns in their distribution across Croatia. Relatively clean and cold karstic rivers in the Dinaric Ecoregion that flow over carbonate bedrock and stable substrates represent far more suitable habitats, harbouring greater diversity of freshwater bryophyte communities in comparison to the rivers of the Pannonian Ecoregion. The Didymodon tophaceus-Apopellia endiviifolia and the Berula erecta-Cratoneuton filicinum are mostly confined to the karstic rivers of the Dinaric Ecoregion under the influence of the Mediterranean climate. The Didymodon tophaceus-Apopellia endiviifolia community is a tufa-forming community associated with the seasonally dry watercourses of small catchment areas and cascades along the larger karstic rivers, while the Berula erecta-Cratoneuton filicinum community is associated with rivers with a permanent and more streamlined flow. The Cinclidotus community is the most widespread in Croatia, having a wide ecological range and with the centre of its distribution being in the Dinaric-Continental Subecoregin, i.e., in fast and cold karstic rivers with permanent flow due to large catchment areas and high precipitation that are characteristic of this subecoregion. The species-rich Oxyrrhynchium hians-Chiloscyphus pallescens community and the species-poor Fissidens pusillus-Veronica beccabunga community are quite rare and are restricted to small eutrophic and turbid rivers of low alkalinity situated mainly in the Pannonian Ecoregion. These communities are characterized by a low share of rheophytes and a high share of species inhabiting the periodically flooded river margins, which is mainly related to low riverbed stability limiting the development of truly aquatic bryophyte communities in Pannonian rivers.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.

Acknowledgments:
We would like to thank Hrvatske vode-the legal entity for water managementfor providing the data on the water chemistry and catchment areas of the surveyed sites.

Conflicts of Interest:
The authors declare no conflict of interest. Table A1. Characteristic species (ϕ ≥ 40.00) are marked in bold, and constant species (frequency ≥ 50.00) are provided in grey cells.   Oenanthe fistulosa L. .    Table A2. Basic descriptive statistic (mean ± SE and range (min-max)) of 35 environmental variables (for abbreviations, see Table 5). Appendix C Table A3. Results of Mann-Whitney pairwise post hoc tests; significant differences between TWINSPAN group pairs are marked with an asterisk (* p < 0.05, ** p < 0.001) (for abbreviations, see Table 5).