Introduction

Regardless of how an urban area is defined—whether by human population size or density, by built-up area or by biome—with time urbanization has increased (O’Driscoll et al. 2010). Urbanization strongly affects biodiversity, both the community composition and species richness (Hill et all. 2018). Our understanding of the mechanism in which urbanization affects communities is still incomplete. This is partly caused by the lack of research referring to invertebrates: most research has focused on vertebrates while invertebrates remain relatively poorly studied (Magle et al. 2012; Prescott and Eason 2018).

In general, there are a small number of studies related to the impact of urbanization on macroinvertebrates in aquatic ecosystems (Gal et al. 2019). Our study considers the impact of urbanization on macroinvertebrate communities both in lotic and lentic habitats and characteristics of the realized ecological niche for typical urban/rural species, mitigating the evident lack of similar papers (Moreyra and Padovesi-Fonseca 2015). Numerous studies have shown that environmental factors in the catchment area (e.g. percentage of agricultural areas, forest and seminatural areas, presence of artificial substrate, etc.) have huge effects on structuring macroinvertebrate communities in rivers (Savić et al. 2013; Theodoropoulos et al. 2015) and ponds (Hill et al. 2018). We hypothesized that urbanization affect the response of macroinvertebrate taxes and their communities depending on environmental characteristics of the (sub) catchment system. In this study, we tested how communities of aquatic macroinvertebrates in lotic and lentic environments respond to urbanization assets.

It seems that the impact of urbanization on species richness in aquatic ecosystems varies within different taxonomic groups (Prescott and Eason 2018). Some studies have shown that species richness of macroinvertebrates often decreases with the increasing levels of urbanization (Hansen et al. 2005; Prescott and Eason 2018). On the other hand, in some aquatic groups, species richness does not seem to be related to the level of urbanization (Faeth et al. 2011; Jones and Leather 2012). Moreover, some taxa have higher species richness in highly urbanized areas (Kudavidanage et al. 2011).

In many aquatic ecosystems, especially in lotic systems, urbanization generally reduces species richness (Roy et al. 2003; Collier and Clements 2011; Johnson et al. 2013; Hassall 2014). This pattern in lentic ecosystems is less clear: some previous studies have reported a significantly lower diversity in macroinvertebrate communities (Hassall and Anderson 2015), while other studies supported the hypothesis that there was no significant difference in the diversity of urban and non-urban lentic systems (Hill et al. 2017).

The impact of urbanization on aquatic ecosystems cannot be separated from the impact of other factors such as the introduction of allochthonous species (Havel et al. 2015), habitat fragmentation (Fuller et al. 2015) and coastal disturbances (Renöfält and Nilsson 2008). For example, the destruction of canopy cover in urban areas (Somers et al. 2013) and increased water pollution (Johnson et al. 1992) can lead to an increase in temperature. In urban areas, similar chemicals can be used as on agricultural land (e.g., phosphorus and herbicides; Creason and Runge 1992), while runoff from urban areas may contain many different pollutants (Voelz et al. 2005). Ecosystems with a greater degree of urbanization cause their constituent species to face altered habitat conditions. Species which are urban exploiters can be identified as “urban” species, while those which are urban avoiders can be identified as “non-urban” species (Palacio 2020). “Urban” species typically have broad habitat tolerances (McKinney and Lockwood 1999; Palacio 2020) and strong dispersal capabilities (Bierwagen 2007; Vergnes et al. 2013). In this paper, we have tested the hypothesis that “rural” species have narrower habitat tolerances than “urban” species. One way to show habitat tolerance is to use species response curves (SRC). The advantage of using SRC along ecological or time gradients led to their recognition as a useful tool in determining species attributes such as optimums or niche width (Jansen and Oksanen 2013). In hierarchical logistic regression modelling, the best model is usually chosen among the available models by using statistical information criteria, i.e., a balance between a model fit to the data and the simplicity of the model (Jansen and Oksanen 2013). There are several different techniques to model SRC, and their benefits were previously determined by various authors. For example, bell-shaped Gaussian distributions were used to model the response of species to environment gradients, but this approach lacked flexibility (Benavides and Vitt 2014). Generalized Additive Models (GAMs) offer an alternative, but estimations of tolerances and optimums are usually evaluated by visual inspection of the fitted models, so the disadvantage of these models is a high dose of subjectivity (Oksanen and Minchin 2002; Benavides and Vitt 2014). On the other hand, most authors (e.g., Oksanen and Minchin 2002; Jansen and Oksanen 2013; Jenačković et al. 2016, Savić et al. 2020) emphasized the use of Huisman–Olff–Fresco (HOF) models as the best option from the ecological point of view.

The aims of this study were (1) to determine the main driver in structuring lotic and lentic macroinvertebrate communities and the factors that influence them in urban versus rural landscapes, (2) to test whether the diversity of macroinvertebrate communities in urban lotic/lentic ecosystems is lower than that in rural landscape, (3) to identify characteristic macroinvertebrate species for urban and rural lotic/lentic ecosystems, and (4) to determine species tolerance ranges and species optimums, with special emphasis on characteristic “urban” and “rural” species, and compare the widths of their realized ecological niches.

Materials and methods

Study site

The Krapiel is a small (length about 70 km) lowland coastal river situated in northwestern Poland. The Krapiel has its source above Kamienny Most Lake near Ścienne village (53° 27′ 09.3″ N 15° 28′ 09.1″ E) and empties into the Ina River near Stargard town (53° 19′ 07.3″ N 15° 03′ 06.6″ E). It has a diversified character: next to stretches with a rapid current and a bottom of stones and gravel, there are stretches with a slower water flow and a sandy or muddy bottom, as well as large marginal pools. The river flows through a variety of landscape structures: forests (hornbeam and oak, riparian forests and alder carrs), reeds and sedges, meadows, pastures, and urbanized areas. The valley contains floodplains and many oxbow lakes. Also present are numerous springs, predominantly helocrenes. The water bodies of the river valley occur in various landscape structures: dense forest complexes, small, isolated forested areas and open land.

Macroinvertebrate sampling and environmental characterization

The research covered the entire length of the river: 13 sites (K1–K13) were established along the river (Fig. 1A). At each of the 13 sites samples were collected from several subsites. The number of subsites taken into account resulted from the spatial differentiation of the particular site and was as follows: two subsites at localities K4, K9 and K12, three at localities K2, K8 and K11, four at K3, K5, K7, K10 and K13 and five at K1 and K6. Altogether there were 45 subsites distributed in such a way as to cover all habitat types in which macroinvertebrates occurred. Fieldwork was conducted from April to October 2010. Samples were collected once a month, in the middle of each month. Each sampling consisting of ten sweeps was taken using a kick-net sampler with 300 µm mesh size and covered an area of about 0.5 m2. In some months samples from some sites and/or sub-sites were not collected because of hydrological conditions (very high or very low water level). A total of 411 samples were collected (248 lotic samples and 163 lentic samples).

Fig. 1
figure 1

 A Map and localities; B example of urban localities; C example of rural localities

Full size image

The catchment basin of the Krąpiel River was divided into 13 subcatchments in such a way that every two adjacent sites (K1–K2, K2–K3, etc.) formed the boundaries of one subcatchment. Each of the 13 subcatchments corresponded to a particular study site (subcatchment 1 to site K1, subcatchment 2 to study site K2, etc.). The analysis of the spatial structure of the catchments was based on a set of landscape metrics calculated using TNTmips software by MicroImages. The classification was based on data from Landsat TM7 28-05-2003. Land cover classes were determined according to the Corine Land Cover database (EEA 2006). In order to categorize sites as urban/rural three parameters at the landscape level were measured and described for each subcatchment: the percentage of built-up area (patch area in the subcatchments), the distance of built-up area from the localities and the percentage of insolation.

Temperature, pH, electrolytic conductivity and dissolved oxygen content were measured with a CX-401 multifunction meter (Elmetron, Zabrze, Poland), BOD5 by Winkler’s method, insolation with a CEM DT-1309 light meter (Merazet, Poznań, Poland) and the remaining parameters (NO3, NH4, PO4) with a Slandi LF205 photometer (Slandi, Michałowice, Poland). Three measurements were performed each time, and the median was used for further analysis.

Statistical analysis

Statistical analyses were performed using PRIMER 7.0 (Clarke and Gorley 2015), MVSP v3.21 (Kovach 2007) SPSS 19, JUICE 7.0 (Tichý 2002) and R statistical package (R Development Core Team 2013).

During the study period (7 months) on some sites one or more urbanization factors were changed. That was the reason why we did not consider our measures (samples) as repeated in statistical procedures.

We categorized the sites as urban (example on Fig. 1B) or rural (example on Fig. 1C) according to following environmental data: the characteristic patch area in the catchments (percentage of built-up area); the parameters in catchments (distance of built-up area from the lotic/lentic ecosystems) (Oudin et al. 2018); the canopy cover characteristic (percentage of insolation) (Berland 2012) using cluster analysis. For cluster analysis based on environmental data, centered and standardized environmental data were classified by the Euclidean distance with between groups average linkage.

Normality of data distribution was tested with the Kolmogorov–Smirnov test. In order to examine diferences in parameters between lotic/lentic and rural/urban samples we used the Mann–Whitney test.

NMDS (non-metric multidimensional scaling) was used for testing the grouping of samples based on the parameters of urbanization, while ANOSIM (analysis of similarities) was used for testing difference between urban and rural samples. In order to determine which physicochemical parameters contributed the most to separation of rural vs. urban samples, we performed PCA (principal component analysis). This analysis was performed on centered and standardized physico-chemical data and parameters in the catchments (percentage of built-up area, distance of built-up area from the lotic/lentic sites, percentage of insolation).

In order to detect the main driver between environmental factors in structuring macroinvertebrate communities in lotic/lentic urban and rural ecosystems we used CCA analysis. CCA (ter Braak 1986) was applied to test the influence of environmental variables on the communities investigated. We first performed a forward selection of environmental variables (Legendre et al. 2011). Also, we used the unrestricted Monte Carlo permutation test (ter Braak and Wiertz 1994) to test the null hypothesis that the selected variables are unrelated to macroinvertebrate communities in the investigated localities.

The response of the macroinvertebrate assemblage diversity to urbanization was determined based on two macroinvertebrate community metrics: average species richness and Shannon index of diversity.

Using the SIMPER (similarity percentages breakdown) procedure, dissimilarities between, and similarities within the above-mentioned groups (lotic/lentic; urban/rural) can be explained with individual species and the composition of macroinvertebrate communities. We conducted SIMPER analysis to identify species with strong preferences for specific habitat types (characteristic “urban” and “rural” species) and thus determine whether certain species typified urban versus rural lentic and lotic sites.

Species Response Curves (SRC) were performed with logistic regression models which are also known as Huisman–Olff–Fresco (HOF) models (Huisman et al. 1993) with respect to the physical, chemical and catchment area data. HOF models are a hierarchical set of five models: I—flat with no response, II—monotone decreasing or increasing, III—monotone increasing to a plateau, IV—symmetric unimodal and V—asymmetric unimodal. This statistical routine was developed by David Zelenyand Lubomír Tichy (https://davidzeleny.net/juice-r/doku.php/scripts:species-response-curves) and it was performed using JUICE 7.0 software package (Tichý 2002). All of the response models were fitted with untransformed data on the species abundances and environmental variables.

Results

Relationship between urbanization and faunistic communities

According to cluster analysis based on the percentage of built-up area (patch area in the catchments/subcatchments), the distance of built-up area from localities and the percentage of insolation (Online Appendix A), two clusters were separated: the first cluster (“RURAL”) comprises rural samples, while the second cluster (“URBAN”) includes urban samples. In total 216 lotic (LOTR) and 135 lentic (LENR) samples were assigned as rural, while 32 lotic (LOTU) and 28 lentic (LENU) samples were assigned as urban .

According to the Mann–Whitney U test, differences between urban and rural samples were found to be significant for all three parameters: the percentage of built-up area (p ˂ 0.001), the distance of built-up area from the localities (p ˂ 0.001) and the percentage of insolation (p ˂ 0.001). In comparison to the cluster of urban samples, the group of rural samples had a lower mean percentage of built-up area (1.08 ± 1.1 vs. 15.25 ± 4.8) as well as a lower mean percentage of insolation (46.42 ± 43.5 vs. 87.8 ± 18.1). On the other hand, the group of rural samples had a higher mean distance from built-up areas (493.93 ± 235) compared to the group of urban ones (225.19 ± 53.4).

To visualize the difference between rural/urban samples we used non-metric multidimensional scaling (NMDS) (Fig. 2). NMDS showed a clear spatial trend distinctly separating the urban and rural groups of samples. Moreover, the ANOSIM procedure indicated significant differences between urban and rural samples (R = 0.483 p = 0.001).

Fig. 2
figure 2

Non-metric multidimensional scaling samples according to percentage of built-up area; distance of built-up area from localities and percentage of insolation

Full size image

Significant differences in physicochemical parameters were detected between rural and urban sample clusters in temperature (p ˂ 0.001) and conductivity (p ˂ 0.001). In other physicochemical parameters, no statistically significant differences were found. However, it is worth mentioning that the mean concentration of NO3 was lower in rural (2.52 mg/l) than in urban samples (3.05 mg/l). The situation with the concentration of NH4 was similar. On the other hand, in rural samples, the mean values of BOD5 were higher than in the urban ones (Online Appendix B).

The first and second PCA (Fig. 3) axes explain 29.6% and 18.2% variation of the variables analyzed, respectively. The first PCA axis correlated with the distance of built-up area from the localities (0.49) and conductivity (0.48). The second PCA axis correlated positively with insolation (0.588) and correlated negatively with the concentration of O2 (− 0.425), NH4 (− 0.309) and the values of BOD5 (− 0.353).

Fig. 3
figure 3

Results of PCA showing differentiation urban and rural samples according to environmental variables: percentage of built-up area (PBU), distance of built-up area (DBU), insolation (INS), temperature (T), conductivity (CON), pH (pH), concentration of NH4 (NH4), concentration of NO3 (NO3), concentration of PO3, concentration of O2 (O2), biochemical oxygen demand (BOD5), urban lotic samples (LOTU), urban lentic samples (LENU), rural lentic samples (LENR), rural lotic samples (LOTR). Environmental factors determining the first axis are shown in red; environmental factors determining the second axis are shown in green (color figure online)

Full size image

The results revealed that urban and rural subsamples were clearly separated especially along the first axis which is mostly determined by the distance of built-up area from the localities. There is no clear separation along the second axis, but there is a pronounced tendency for urban subsamples to group toward higher values of insolation and a higher percentage of built-up area. On the other hand, there is no clear separation between lotic/lentic samples (Fig. 3).

CCA was performed with the aim to examine the effects of environmental variables (physico-chemical) and parameters in the catchments (the percentage of built-up area, the distance of built-up area from the samples and the percentage of insolation) on macroinvertebrate community composition. For lotic samples (Fig. 4), the first and second axes explain 31.97% and 19.2% of species variation of the variables analyzed, respectively. Axis 1 is mostly determined by temperature (− 0.475), while axis 2 is determined by the concentration of oxygen (− 0.58), the percentage of built-up areas (− 0.52) and insolation (0.47). The urban lotic sites are mostly grouped along the negative part of the X- and Y-axes (third quadrant) and characterized by higher temperature and a higher percentage of built-up areas (Fig. 4).

Fig. 4
figure 4

Canonical correspondence analysis for lotic communities according to environmental variables: percentage of built-up area (PBU), distance of built-up area (DBU), insolation (INS), temperature (T), conductivity (CON), pH (pH), concentration of NH4 (NH4), concentration of NO3 (NO3), concentration of O2 (O2), urban lotic samples (LOTU), rural lotic samples (LOTR). Key factors: T, O2, PBU and INS. Environmental factors determining the first axis are shown in red; environmental factors determining the second axis are shown in green (color figure online)

Full size image

For lentic samples (Fig. 5), the first and second axes explain 29.02% and 18.8% of species variation of the variables analyzed, respectively. Axis 1 is mostly determined by insolation (− 0.709); BOD5 (0.54), the concentration of oxygen (0.53) and pH (0.51), while axis 2 is determined by temperature (0.9) and the distance of built-up area (0.3). The urban lentic sites are generally grouped along the positive end of the X- axis and characterized by a higher percentage of built-up areas and higher values of NO3 and NH4 concentrations (Fig. 5).

Fig. 5
figure 5

Canonical correspondence analysis for lentic communities according to environmental variables: percentage of built-up area (PBU), distance of built-up area (DBU), insolation (INS), temperature (T), conductivity (CON), pH (pH), concentration of NH4 (NH4), concentration of NO3 (NO3), concentration of PO3, concentration of O2 (O2), biochemical oxygen demand (BOD5), urban lentic samples (LENU), rural lentic samples (LENR). Key factors: INS, BOD5, O2, pH nd T. Environmental factors determining the first axis are shown in red; environmental factors determining the second axis are shown in green (color figure online)

Full size image

Faunistic characterization of rural and urban communities

Overall, 134 macroinvertebrate taxa with 9375 individuals from 411 samples taken from 13 sites were collected. From lotic samples we detected 88 taxa. The most abundant was species Torrenticola amplexa, with 1537 collected individuals. From lentic samples, we detected 111 taxa. The most abundant was Hygrobates longipalpis, (533 ind.). The highest species richness detected in lotic samples was 18 species, while in lentic it was 14. The average Shannon’s diversity in lotic samples (H′ = 2.1) was slightly higher than in lentic ones (H′ = 1.9). The Mann–Whitney U test did not reveal a significant difference in the average number of species between lotic and lentic samples (p = 0.4).

From samples classified as rural we detected 126 taxa with 7814 specimens, while from samples classified as urban, 49 taxa with 1561 specimens were identified. No significant differences in Shannon’s diversity index between urban and rural samples were found. However, it is worth mentioning that the mean values of Shannon’s diversity index for rural samples (H′ = 2.04) were higher than for urban ones (H′ = 1.8). Moreover, no significant differences in regard to Shannon’s index were found between lotic and lentic ecosystems. In regard, to rural samples, the mean values of Shannon’s diversity index were the same for lotic and lentic samples (H′ = 2.0). On the other hand, if we consider only urban samples, a clear difference in diversity values is visible. Shannon’s diversity index is clearly lower in lentic urban samples (H′ = 1.05) than in lotic ones (H′ = 2.7) (Fig. 6).

Fig. 6
figure 6

Mean values of Shanon’s diversity index (a) and number of species (b) in urban lotic ecosystems (LOTU), urban lentic ecosystems (LENU), rural lotic ecosystems (LOTR) and rural lentic ecosystems (LENR)

Full size image

When Shannon’s diversity index was considered separately for lotic samples, significant differences between rural and urban samples were found (p = 0.017). Concerning lentic ecosystem significant differences between rural and urban samples (p ˂ 0.001) were also present.

No significant differences were found in species richness between urban and rural samples (p = 0.9) (when we compared all urban vs. all rural samples), nor between lotic and lentic samples (p = 0.4) (when we compared all lentic vs. all lotic samples). Unexpectedly, the mean number of species in urban lotic samples (6.8) was higher than in rural lotic ones (3.8) (similar to Gadawski et al. 2016). As expected, in rural lentic samples, the mean number of species (4.4) was higher than in urban lentic samples (2.2) (Fig. 6B). Concerning lotic ecosystems significant differences in species richness between urban and rural ecosystems (p = 0.001) were found. Significant differences (p ˂ 0.001) between urban and rural samples were also found for the lentic environment.

Fig. 7
figure 7

The species response curves for ten parameters: p percentage of built-up area, d distance of built up area, i percentage of insolation, t temperature, c conductivity, H pH, b BOD5, O2 concentration of oxygen, PO3 concentration of PO3, NH4 concentration of NH4. See Online Appendix D for detailed data (minimum values, maximum values, tolerance and raw optima values) on species response

Full size image

SIMPER analysis was performed with four previously determined groups of samples: lotic-urban (LOTU), lentic-urban (LENU), lotic-rural (LOTR) and lentic-rural (LENR) (Online Appendix C). According to SIMPER analysis, within LENU group average similarity was 7.92%; within LOTU 13.45%; within LOTR 6.61% and within LENR group 7.71%. The main contributors to structuring the community of macroinvertebrates in LENU group were Piona sp. (36.04%) and Hydrodroma pilosa (29.4%). Torrenticola amplexa (24.81%) and Mideopsis crassipes (23.85%) were characteristic of LOTU assemblages. For rural samples, in LOTR group Torrenticola amplexa (27.9%) and Hygrobates longipalpis (20.48%) were the main contributors; in LENR group, Hygrobates longipalpis (37.55%) and Mideopsis orbicularis (24.19%) were the main contributors. The analysis was conducted for species level, with the exception of the genus Piona, where there were numerous deutonymphs identified only for the genus level. However, since the whole genus of Piona is characteristic of lentic waters, we decided to include it as characteristic of this type of habitats.

Higher average dissimilarities in community’s structure according to SIMPER analysis were detected between urban lentic (LENU) and rural lentic (LENR) then between urban lotic (LOTU) and rural lotic (LOTR) (Online Appendix C).

It was noticed that Torrenticola amplexa prefers lotic ecosystems but without clear preference for either urban or rural type. The presence of Mideopsis crassipes is related to urban ecosystems while the occurrence of M. orbicularis is more linked to the rural environment. Therefore, HOF models for environmental parameters were performed for these two species, in order to test the hypothesis that “urban” species typically have broad habitat tolerances. Additional HOF models were performed for Hygrobates longipalpis, a species characteristic of rural ecosystems and for Piona sp. and Hydrodroma pilosa, two species characteristic of urban lentic ecosystems.

According to HOF models, Hygrobates longipalpis, which was characterized as a typical representative of rural ecosystems, showed narrow tolerance for five out of ten discussed parameters: the distance of built-up area, temperature, conductivity, BOD5 and the concentration of PO3 (Fig. 7). The frequency of encountering this species decreases with the increase of the percentage of built-up areas. It showed narrow optimums for temperature (optimum at lower levels around 16 °C), conductivity (optimum at lower levels around 230 µS/cm), BOD5 (optimum at lower levels around 1.5) and PO3 (optimum around 0.3 mg/l) (Fig. 7). Mideopsis orbicularis showed narrow tolerance levels for four parameters. On the other hand, two species, M. crassipes and Hydrodroma pilosa whose presence was more related to the urban environment, showed a narrow tolerance for just two parameters, while Piona sp., another typical representative of urban ecosystems, did not show a narrow range of tolerance for any of the ten studied parameters (Fig. 7).

A comparison of two species from the same genus showed that M. orbicularis, which is characteristic of the rural environment, had a narrow optimum for temperature (around 18 °C), while its congener, M. crassipes, was indifferent to temperature increase. The frequency of occurrence of the latter species increased with an increase of BOD5, while it tolerated a high level of insolation as well: it remained within the optimum until 50%. The optimum of insolation for M. orbicularis was at somewhat lower levels of 25%.

Discussion

Environmental difference between rural and urban samples

This study aimed to identify the main factors that affect structuring macroinvertebrate communities in urban and rural ecosystems, separately in the lentic and lotic environment. These habitats are distinguished by different combinations of physicochemical factors and also characterized by different faunistic assemblages. Unlike the lotic/lentic differentiation which is relatively easy, the separation of the samples and sites to “urban” and “rural” is not always simple. For the purposes of this paper, samples were classified based on three parameters which are considered as indicators of the degree of urbanization: the percentage of built-up area, the distance of built-up area from the localities and the percentage of insolation.

Our study reveals that urban and rural samples clearly differ in physico-chemical parameters and characteristics in catchment. The PCA analysis showed that distance from the study sites to built-up areas and conductivity were the main factors contributing to their separation. The significance of distance from the study sites to built-up areas and increasing of conductivity was already highlighted by other authors (Chusov et al. 2014; Moreyra and Padovesi-Fonseca 2015; Prescott and Eaton 2018). According to some authors (Novikmec et al. 2016) proportion of urban land was positively related with water conductivity, i.e., water ecosystems with high proportion of urban areas in their catchments were characterized by high values of conductivity. There was a pronounced tendency of urban samples to group toward higher values of insolation and a higher percentage of built-up area. On the other hand, there was no clear separation between lotic/lentic samples. It can be observed that the separation of samples from the urban and rural environment was more pronounced than the separation of samples from the lotic and lentic landscapes.

For the lotic sites, temperature, the percentage of built-up area, insolation, and oxygen concentration were the main factors that were significantly associated with changes in community composition. On the other hand, for the lentic sites, insolation, temperature, and BOD5 were recognized as the main factors were significantly associated with changes in community composition.

The degree of insolation is one of the most important factors in structuring the macroinvertebrate communities both in lotic and lentic ecosystems. Insolation (as well as canopy cover) was used as a surrogate measure of disturbance by agriculture and urbanization (O’Driscoll et al. 2010), emphasizing the importance of the preservation of riparian habitats for the macroinvertebrate community (Rios and Bailey 2006). A separate assessment of this factor during our study yielded data supporting the relationship between riparian habitat disturbance and the macroinvertebrate community in both types of aquatic ecosystems. Such changes in macroinvertebrate communities were registered in the Krąpiel River after hydrotechnical works, which resulted in changes in the vegetation cover (Szlauer-Łukaszewska and Zawal 2014; Stępień et al. 2015, 2019; Zawal et al. 2015a, b, 2016a, b; Buczyński et al. 2016; Dąbkowski et al. 2016; Płaska et al. 2016; Stryjecki et al. 2016).

Temperature is another factor which is important for the macroinvertebrate community of both lotic and lentic ecosystems. Aquatic communities in the urban environment are affected by an additional increase in water temperature due to decreased canopy cover (Somers et al. 2013). Such influence of increasing insolation and temperature has been demonstrated for some species of water mites, mollusks, and ostracods in the Krąpiel River especially in lentic habitats (Stryjecki et al. 2016; Zawal et al. 2016c, 2017, Szlauer-Łukaszewska and Pešić 2020).

It seems that the structure of macroinvertebrate assemblages in the lentic environment was more related to the distance from built-up area and to the percentage of built-up area, while in the lotic environment the percentage of built-up area was a more important factor. The importance of the distance from built-up area has already been recorded in previous studies (e.g., Kinzing et al. 2005; Prescot and Eason 2018). This factor broadly reflects the overall change of anthropogenic disturbance, without identifying specific elements of urbanization that might affect community composition (Kinzing et al. 2005).

Effect of urbanization on macroinvertebrate communities

This study did not reveal statistically significant differences in Shannon diversity index between samples from the rural and urban environment. Moreover, there were no statistically significant differences in Shannon’s diversity between rural samples from both the lotic and lentic environment. On the other hand, Shannon’s diversity of urban samples was statistically significantly lower for the lentic compared to the lotic environment. In lentic habitats, due to a slower water flow, there is a greater accumulation of organic matter than in the lotic environment, which leads to a decrease in oxygen content and, as a result, a reduction in biodiversity. This might be used as a basis for the hypothesis that urbanization has a different impact on lotic and lentic systems, which could lead to the conclusion that threats caused by urbanization are greater for lentic ecosystems or at least for their diversity. This idea is supported by the fact that SIMPER analysis showed higher average dissimilarities in macroinvertebrate community structure between urban and rural samples from the lentic environment than between urban and rural samples from the lotic environment (Online Appendix C).

The mean species richness in samples from the lotic environment in our study was unexpectedly higher in lotic “urban” samples in comparison to “rural” samples. This might be explained by the “intermediate disturbance hypothesis” (Connell 1978), where intermediate disturbance disrupts the succession of the ecosystem and maintain a high equilibrium of biodiversity, restraining the dominant species and allowing less competitive species to coexist. Several studies have stressed that in lotic systems urbanization usually reduces species richness (Roy et al. 2003; Cuffney et al. 2010; Collier and Clements 2011; Hassall 2014). The loss of diversity as well as changes in structuring assemblages vary between individual faunal groups and some authors (e.g. Prescott and Eason 2018) also found that species richness of some particular insect group did not differ between urban and rural sites in lotic systems.

For the lentic environment, the results of our study revealed a different pattern, and the mean number of species in rural samples was statistically significantly higher than in samples from the urban environment. The results of our study support the claim of Hassall and Anderson (2015) that species richness in urban environments is significantly lower than in rural ones.

The number of characteristic species is much higher in “rural” samples. Of the 134 taxa registered by our study, 8 were recorded only in urban samples while 85 were found only in rural samples. SIMPER analysis showed that the characteristic species for the lotic environment was Torrenticola amplexa. This is in agreement with the claim that this species is strongly associated with greater velocity and is represented with the highest number of individuals in mesohabitats with a well-developed lotic zone (Zawal et al. 2017).

On the other hand, the occurrence of M. crassipes was related to urban ecosystems while the presence of its relative M. orbicularis was more related to the rural environment. This is in line with the findings that M. orbicularis is negatively correlated to BOD5 (Zawal et al. 2015a, b), meaning that it prefers “cleaner” ecosystems, in this case the rural ones. The obtained results support the claim that urbanization has very different effects even on closely related taxa (Prescott et Eason 2018).

The characteristic species of rural systems (both lotic and lentic environments) was Hygrobates longipalpis. It has been known for a long time that this species inhabits both standing and running waters throughout Europe, but recently Pešić et al. (2019a) attributed standing water-dwelling populations to H. prosiliens, and populations from slowly flowing rivers to H. longipalpis. Hydrodroma pilosa and Piona sp. are characteristic species for urban lentic ecosystems. As shown by the results of our study, all characteristic species in urban and rural samples belong to water mites. Most water mite species inhabiting rivers are strictly bound to particular microhabitat or stream sectors, which makes them sensitive to climate and human induced environmental changes (Zaval and Pešić 2018). At the local scale, substrate composition and temperature have a decisive impact on the abundance and species composition of water mites, often resulting in each subsite being represented by an assemblage of species which displays distinct habitat preferences (Gerecke 2002). Recent studies hypothesize that biotic metrics that include water mites are more sensitive in assessing the response of macroinvertebrate assemblages to anthropogenic modification (Pešić et al. 2019b).

For species characterizing urban/rural lentic and lotic ecosystems, HOF models were performed according to environmental parameters in order to compare the tolerance ranges and check the hypothesis that “urban” species typically have broad habitat tolerances. Species Hygrobates longipalpis and Piona sp. were at the opposite sides of the spectrum regarding the preference for the degree of urbanization. The first species showed a narrow tolerance range for five out of ten tested parameters, and it was characterized as typical of rural ecosystems. On the other hand, Piona sp., which was characterized as typical of urban ecosystems did not show narrow tolerance ranges for any of the tested parameters. This is supports the hypothesis that highly “urban” species have broad habitat tolerances (at least for some environmental parameters).

The ecological preferences of species belonging to the same genus may be very different, as shown by our study with the species pair M. orbicularis and M. crassipes. The occurrence of the first species was linked with the presence of the organic bottom in lentic habitats, and the second one with the mineral substratum in lotic environments (Zawal et al. 2015a, b). In our study SIMPER analysis revealed that M. orbicularis is characteristic of rural and M. crassipes of urban ecosystems.

The differing responses to urbanization observed in this study were probably a consequence of the ecology of the species themselves as well as the differences between lentic and lotic ecosystems (Prescott and Eaton 2018). The results of our study support the hypothesis that lotic and lentic ecosystems do not react to urbanization in the same way. Different patterns of species richness of lotic and lentic assemblages from urban versus rural landscapes indicate that the impact of urbanization and consequently the implementation of conservation measures should be viewed separately within the lentic and lotic gradient. An increase in urbanization in the future through changes to which riparian cover will be exposed, and, at the same time, through changes in temperature, will lead to additional pressure on aquatic habitats and macroinvertebrate assemblages. We have assumed that the conditions in the macroinvertebrate community may be used to determine many characteristics of the riparian system, including the degree in which it is disrupted, which will eventually require their involvement in the management of the surrounding terrestrial environment. Future studies building up on this work will be directed towards determining the links between macroinvertebrate and riparian areas, with the goal of suggesting detailed measures for both conservation and mitigation.