INTRODUCTION

Among the aquatic ecosystems, rivers benefit human communities by providing a supply of water, nutrient retention, removal of toxins, microclimate stability, opportunities for tourism, and are valued by local cultures (^{Brismar, 2002}; ^{Dudgeon, 2019}). The main cause of loss of the ecological integrity and deterioration of these ecosystems are human activities (^{Carpenter et al., 2011}; ^{Dudgeon, 2019}), while the major threats for freshwater biodiversity are overexploitation, water pollution and flow modification, the invasion of exotic species, land use change, and climate change (^{Dudgeon, 2019}). At the basin scale, land use changes influences in the stream conditions modifying water characteristics, sediment supply and deposition, affecting bank stability, and consequently the aquatic biota (^{Strayer et al., 2003}; ^{Townsend et al., 2003}; ^{Allan, 2004}).

In order to design adequate proposals for new management and conservation of fluvial ecosystems, it has been proposed to study selected indicator groups and how the ecological elements and processes in these catchments have changed in the long term to reveal their current ecological condition and future threats (^{Ramírez & Gutiérrez-Fonseca, 2014a)}. Historical analyzes of changes in aquatic communities offers information about the current conservation status of aquatic ecosystems, in order to infer factors that have impacted these systems and obtain insight into the changing conditions of the surrounding watershed (^{Karr, 1981}; ^{Fausch et al., 1990}). Aquatic macroinvertebrates have a range of preferences for environmental conditions, so shift in the assemblages may reflect changes in the aquatic ecosystem and human impacts over time (^{Li et al., 2012}). However, the long-term perspectives and historical comparisons in aquatic macroinvertebrate communities remains often short (^{Jackson & Füreder, 2006}).

Aquatic macroinvertebrates are especially useful to evidence changes in river ecosystems due to anthropic impacts (^{Barbour et al., 1999}; ^{Bonada et al., 2006}; ^{Ligeiro et al., 2013}). Among them, aquatic insects are generally the most abundant and diverse, as they are one of the most ecologically important groups (^{Macadam & Stockan, 2015}) especially in tropical and subtropical zones (^{Dudgeon, 2008}). They are the main primary consumers and are responsible for transferring the energy of primary productivity to other trophic levels of food chains, and there are elements within this group that are important predators (^{Hanson et al., 2010}; ^{Macadam & Stockan, 2015}). Aquatic insects can have highly specific functions in the ecosystem, such as filterers, gatherers, shredders, predators, piercers and scrapers (^{Merritt et al., 2008}; ^{Hanson et al., 2010}; ^{Ramírez & Gutierrez-Fonseca, 2014b)}. Importantly, because aquatic insects deploy a wide array of generalist and specialist feeding strategies, occupy several microhabitats, and have different responses and sensitivities to habitat degradation, they are considered highly useful biological indicators of stream ecological condition (^{Karr & Chu, 1999}).

On the other hand, we must also bear in mind that macroinvertebrate assemblages can not only be affected by pollution or degradation. It has been seen that these assemblages vary due to the flow regime and sediment deposition (^{Díaz-Rojas et al., 2020}). So that high areas of a basin may have greater diversity because the variations of flows given by the slopes allows greater heterogeneity in the landscape than in the low sections where the slope decreases as well as flow (^{Mesa, 2010}). Likewise, when the flow increases in the rainy season, the communities are modified (^{Quesada-Alvarado et al., 2020}). Therefore, it is important to consider the seasonality of the samples, being also throughout the year that differences have been seen in the macroinvertebrate assemblages as the flows and chemical composition of the rivers is modified, for example, by the leachate of the soil in rainy season (^{Leal-Bastidas et al., 2021}).

In the tropics, the influence of hydrological, physical and chemical alterations upon macroinvertebrate communities remains poorly understood (^{Md Rawi et al., 2014}; ^{Ramírez et al., 2015}). Several studies on the ecology of aquatic insects in Latin America have been reported, with emphasis on the relationship with abiotic factors (^{Ramírez & Gutierrez-Fonseca, 2014a}). For example, recently, ^{Kohlmann et al. (2021)} include in their study a analyzes of the relationships between functional feeding groups of aquatic macroinvertebrates with physicochemical such as NO³⁻ K⁺, biochemical oxigen demand (BDO), oxygen saturation, and pH; ^{Díaz-Rojas et al. (2020)} relates depth, flow velocity, channel width and roughness of the substrate with macroinvertebrate assemblages composition and functional traits; ^{Quesada-Alvarado et al. (2020)} describe the relationship between the aquatic macroinvertebrate assemblages with physicochemical and habitat variables, such as NO³⁻, substrate and flow; and, ^{Mosquera-Restrepo & Peña-Salamanca (2019)} explain the relationships between aquatic macroinvertebrates assemblages with dissolved oxygen, BOD, total dissolved solids, and turbidity.

In Mexico, this approach has been used to assess the biotic integrity of rivers in the Río Chiquito basin in the State of Michoacán (^{Piñón-Flores et al., 2014}), variation of macro-invertebrates in the Laguna de Tecocomulco in the State of Hidalgo for one year period (^{Rico-Sánchez et al., 2014}), and impacts of mining activities in three rivers of the Sierra Gorda Biosphere Reserve (^{Rico-Sánchez et al., 2022}). However, these studies involve brief spatial and temporal scale. For these reasons, the present study has the main goal of compare aquatic insect assemblages data of 1997 with data of 2014. Community structure, diversity, functional feeding groups, and the associations with water characteristics were described in two major subtropical river drainages in east-central Mexico, to interpret the ecological impairment indicated by the patterns founded.

MATERIAL AND METHODS

Study area: The study area includes the Pánuco and Lerma-Chapala river drainages, located in east-central Mexico (Fig. 1). It has a subtropical area in the northeast, located in the Eastern Sierra Madre and the Neovolcanic Belt. Central Mexico has the most degraded river drainages in the country (^{Mercado-Silva et al., 2006}). The Lerma-Chapala and Pánuco river drainages are two of the most important basins of this region, and have been highly impacted by loss of vegetation cover (>30%), expansion of cultivated pastures for livestock, increased agricultural activities, combined with expanded industrialization and urbanization (^{Cuevas et al., 2010}). The Lerma-Chapala river drainage shows an evident problem of physical and chemical anthropogenic transformation, and is considered as the most degraded in Mexico (^{Cotler-Avalos et al., 2004}). At the present, the headwaters of both drainages are being considered for special protection status as water reserves in Mexico by the National Commission of Water (^{Comisión Nacional del Agua, 2011}).

Figure 1 Geographic location of study sites. 1 = Fracción Sánchez, 2 = Planta-La Hacienda, 3 = Puente la Plazuela, 4 = Pinihuan, 5 = Canoas, 6 = Quinta Matilde, 7 = El Realito, 8 = Quiotillos, 9 = El Salto, 10 = Presa del Carmen, 11 = Presa de Rayas, 12 = Comonfort, 13 = La Quemada, 14 = Los Galvanes, 15 = El Xote, 16 = El Oasis, 17 = Chuveje, 18 = Carpintero, 19 = Rascón, 20 = Tamasopo, 21 = Jalpan, 22 = Ayutla, 23 = Santa María (before of Adjuntas), 24 = El Carrizal (Santa Maria after of Adjuntas), 25 = Río Grande, 26 = Calvillo, 27 = Sabinolandia (El Salto de los Salados). 

A total of 27 sampling sites were selected in permanent rivers and were sampled in 1997 and 2014. We chose some of the main waterways in the following five states (Fig. 1) which included: 1) Aguascalientes: San Pedro River and Calvillo River; 2) Jalisco: Grande River; 3) Guanajuato: Laja River and Apaseo River; 4) Querétaro: Extóraz River, Huimilpan River, Querétaro River, San Juan River, Jalpan River and Santa María River; and, 5) San Luis Potosí: Verde River. The field work was conducted in the dry season (from February to May), when the conditions of habitat and biological community of rivers are more stable (^{Pérez-Munguía et al., 2007}; ^{Lyons et al., 1995}) and the effect of the human activities are more evident (^{Moncayo-Estrada et al., 2015}).

Data collection: Water physical and chemical parameters were measured with a multimeter probe (Hach Hydromet Quanta, Loveland, Colorado, USA), and we included pH, dissolved oxygen (mg/L), and temperature (°C). Aquatic insects were sampled using a D-net (300 mm of diameter and 300 μm of mesh size) in all different types of reachable habitat, with a sample effort of 60 minutes per study site. During 1997, aquatic insects were preserved in alcohol in 125 ml jars and brought back to the laboratory and separated from detritus. During 2014, insects were separated in situ and were deposited into a plastic vial and preserved in 80% alcohol solution for further transport to the laboratory (Biotic Integrity Lab at Universidad Autónoma de Querétaro). The aquatic insects were identified to the taxonomic level of family using specialized keys (e.g., ^{Arce-Pérez & Roughley, 1999}; ^{Merritt et al., 2008}; ^{Bueno-Soria, 2010}; ^{Springer et al., 2010}). We used the family taxonomic level because it has proven to be a good indicator of the level of ecological disturbance in fluvial ecosystems (cf. ^{Marshall et al., 2006}; ^{Serrano-Balderas et al., 2016}; ^{Wright & Ryan, 2016}), allows for the categorization of functional traits for the different families in most cases, and highest taxonomic level still providing sufficient resolution regarding biological traits of the organisms, saving time and effort to reach lower taxonomic categories. The functional feeding groups (FFG) were obtained from ^{Ramírez & Gutiérrez-Fonseca (2014b)}.

Statistical analysis

Water physicochemistry: Paired T-test were made to assess the difference of the values of each water parameter between two years (1997 and 2014). To analyze and elucidate patterns of all physical and chemical parameters in both years, a principal component analysis (PCA) was conducted, and we normalize all variables using division by their standard deviations because the variables are measured in different units.

Aquatic insect assemblages: We calculated the true diversity as proposed by ^{Jost (2006, 2007)}, through assessing of effective numbers of elements at family level, that refers to the numbers of taxa equally probable and necessary to obtain a diversity value (^{Jost, 2007}). This approach is considered logical and works intuitively, unlike other indices such as the Shannon entropy, which measures the uncertainty degree of a species (^{Jost, 2006}). True alpha and gamma diversity of first order was obtained to sensitize the index to the abundant species, because in aquatic insect communities are common to find this pattern of very abundant and rare taxa, which means an inequitable distribution of abundances among all taxa. Jackknife estimator was used because it is appropriate for this group of organisms (^{Basualdo, 2011}; ^{Martínez-Sanz et al., 2010}). In addition, the true beta diversity was obtained as the effective number of elements in the data set (true gamma diversity) divided by the average number of effective elements of the samples (true alpha diversity); where, one is the minimum number which we can obtained, indicating that all communities are exactly the same, and maximum value are equal to the total of communities (N) (^{Jost, 2007}). We applied paired T test to assess the difference of diversity values between 1997 and 2014, and a similarity percentage analysis (SIMPER) using the Bray-Curtis similarity measure (multiplied with 100), based in abundance per family and FFG was used to identify which taxon discriminates among periods (^{Clarke, 1993}).

Responses of aquatic insect assemblages to water physicochemistry: To assess the effect of water physicochemistry on aquatic insect assemblages we made a non-metric multidimensional scaling (NMDS) based on the Bray-Curtis similarity index, which can be used with cero values in data sets (^{Bray & Curtis, 1957}). NMDS was applied using the abundance per family and per FFG. In the NMDS the environmental variables were associated to the axis and represented with vectors in the plot (^{Hammer et al., 2001}). We correlated by the Spearman method the PCA values and NMDS scores to understand the relationships between abiotic variables and aquatic insect assemblages overall, considering the intrinsic relationships, as was used by ^{Escalera-Vázquez & Zambrano (2010)} in a study of the effect of variation in abiotic factors on fish assemblages.

The paired T-test and Spearman rank order correlation analysis (^{Zar, 2014}) were made using the statistical software SPSS version 20 (^{IBM corp., 2011}). The true diversity values (^{Jost, 2006}, ²⁰⁰⁷) were estimated with SPADE software (^{Chao & Shen, 2010}). The multivariate analysis PCA, NMDS and SIMPER (^{Quinn & Keough, 2002}) were obtained using PAST version 3.07 (^{Hammer et al., 2001}).

RESULTS

Water physicochemistry: The paired T-test shows a significantly decrease (p<0.001) of dissolved oxygen between 1997 and 2014 considering all study sites (i.e., both basins), from 8.2 ± 3.3 mg/L to 3.6 ± 2.2 mg/L. There are no significant differences between temperature (20.3 ± 4.31 to 20.5 ± 4.8 °C) and pH (8.01 ± 0.43 to 7.8 ± 0.35) of both years (p>0.05). Nevertheless, the PCA analysis showed a subtle tendency gradient of segregation of data between 1997 and 2014. Three of the main components (PC) had moderately related variables (0.75-0.50). Of these PC1 (eigenvalue=1.76) explained 58.73% of the variance, PC2 (eigenvalue=0.75) explained 25% and PC3 (eigenvalue=0.48) 16.26%. The first component (PC1) showed moderate positively association with the three variables with correlation coefficients of 0.63(pH), 0.58 for dissolved oxygen (DO) and 0.51 for temperature. While the second was strongly positively associated (>0.75) with temperature (0.81) and moderate negatively associated with DO (-0.55) and not so with the pH (-0.15). Whereas the PC3 showed a negative association with pH (-0.76) and positive with DO (0.59) and temperature (0.25). The study sites ordination resulted located diagonally from upper left corner to the lower right corner, following a decrease of dissolved oxygen and pH values, and from the lower left corner to the upper right corner following an increase in temperature. The sites located in the upper left corner zone, comprises mainly the sites sampled in 2014 (Fig. 2).

Figure 2 Principal Components Analysis based on physicochemical parameters in rivers of two subtropical river drainages in east-central Mexico (Lerma-Chapala and Pánuco). Data from 1997 (triangles), data from 2014 (circles). Ayutla = Ayu; Calvillo = Cal; Canoas = Can; Carpintero = Car; Chuveje = Chu; Comonfort = Com; El Carrizal (Santa María after of Adjuntas) = SMD; El Oasis = EO; El Realito = ER; El Salto = ES; El Xote = EX; Fracción Sánchez = FS; Jalpan = Jal; La Hacienda = LH; La Quemada = LQ; Los Galvanes = LG; Pinihuan = Pin; Presa de Rayas = PR; Presa del Carmen = PC; Puente la Plazuela = PP; Quinta Matilde = QM; Quiotillos = Qui; Rascón = Ras; Río Grande = RG; Sabinolandia (El Salto de los Salados) = Sab; Santa María (before of Adjuntas) SM; Tamasopo = Tam. 

Aquatic insect assemblages: A total of 71 aquatic insect families were obtained including both drainages. We collected 47 families during 1997 and 61 families during 2014 (Table 1). We found more representativeness of taxa during 2014 and more gamma diversity for both river drainages (Lerma-Chapala, ¹D_a = 8.2 and Pánuco ¹D_a = 7.82 during 1997; and Lerma-Chapala ¹D_a = 9.17 and Pánuco ¹D_a = 12.13 during 2014). The beta diversity was higher in 2014 (¹D_b = 2.8) than in 1997 (¹D_b = 1.49) on Lerma-Chapala River drainage, and lower in 2014 (¹D_b = 1.4) than 1997 (¹D_b = 1.61) on Pánuco River drainage. Global gamma and beta diversity was also higher in 2014 (¹D_g= 13.34; ¹D_b = 2.08) than 1997 (¹D_g= 9.30; ¹D_b = 1.81). However, the Lerma-Chapala river drainage showed higher alpha and beta diversity of families in the 70% of study sites during 1997 (Table 2). These results are consistent with the results of paired T-test that showed not significantly difference of alpha diversity of all sites between years (p=0.181).

The SIMPER analysis showed that the main families with contribution for abundance dissimilarity between 1997 and 2014 were Chironomidae (24.8%), Baetidae (16.5%), Coenagrionidae (8.5%), Veliidae (4.9%), Corixidae (4.2%), Culicidae (3.9%), Caenidae (3.7%) and Heptageniidae (1.3%). The abundance of Chironomidae was 61 ± 185 in 1997 and 116 ± 257.7 in 2014; of Baetidae was 60.7 ± 104.9 in 1997 and 38.9 ± 57.2 in 2014; Coenagrionidae showed 2.44 ± 6 in 1997 and 46 ± 129.3 in 2014; Veliidae 2.26 ± 5.1 in 1997 and 22.2 ± 55 in 2014; Corixidae showed 6.48 ± 30.7 of mean abundance in 1997 and 20 ± 51 during 2014; Culicidae 0.5 ± 1.5 in 1997 and 40.7 ± 195.2 in 2014; Caenidae 25.3 ± 99.3 in 1997 and 2.29 ± 13.5 during 2014; and Heptageniidae 3.22 ± 10.15 in 1997 and 0.07 ± 0.38 in 2014.

We found the six FFG: gatherers, filterers, predators, shredders, piercers, and scrapers. The most abundant FFG in both years was the gatherers, follow by predators, and piercers were the rarest (Table 3). SIMPER analysis based on abundance per FFG showed that gatherers contributed with 55.6% to the dissimilitude, with change from 79% of gatherers in 1997 to 44% in 2014; predators contributed with 32% and the quantity of individuals changed from 23 ± 23 in 1997 to 120 ± 188 in 2014; filterers contributed with 8.9% and changed from 8 ± 33 to 48 ± 194 individuals between 1997 and 2014; scrapers contributed with 2% and changed from 3.5 ± 10.1 in 1997 to 0.8 ± 2.4 in 2014; shredders contributed with 0.7% and piercers with 0.5%.

Relationships between aquatic insect assemblages and water characteristics: The NMDS based on number of individuals per site (Fig. 3), showed a pattern where taxa were ordinated in a gradient of dissolved oxygen and pH decrease from the upper left corner to the lower right corner. Additionally, the ordination analysis showed a gradient of decrease in abundance per site in the same direction (from the upper left to the lower right corner). We found some sites with contrasting differences between years, where Comonfort (Com) had 1008 individuals in 1997 and 21 individuals in 2014; Sabinolandia (Sab) changed from 298 individuals in 1997 to 36 in 2014; Fracción Sanchez (FS) increased the number of individuals from 85 to 2491 in 1997 and 2014, respectively; however, ~40% (1017 individuals) belong to Culicidae family and the insect diversity decrease in time (Table 2). The NMDS analysis shows a relationship between low concentrations of dissolved oxygen and lower pH values with fewer number of individuals per site; however, it shows no tendency of grouping by sampling years.

Figure 3 Non metric multidimensional scaling (NMDS) based on the number of individual of aquatic insects per family, in rivers in the Lerma-Chapala and Pánuco River drainages, with data of 1997 (triangles) and 2014 (circles). Stress: 0.19. Temp = temperature, DO = dissolved oxygen. Ayutla = Ayu; Calvillo = Cal; Canoas = Can; Carpintero = Car; Chuveje = Chu; Comonfort = Com; El Carrizal (Santa María after of Adjuntas) = SMD; El Oasis = EO; El Realito = ER; El Salto = ES; El Xote = EX; Fracción Sánchez = FS; Jalpan = Jal; La Hacienda = LH; La Quemada = LQ; Los Galvanes = LG; Pinihuan = Pin; Presa de Rayas = PR; Presa del Carmen = PC; Puente la Plazuela = PP; Quinta Matilde = QM; Quiotillos = Qui; Rascón = Ras; Río Grande = RG; Sabinolandia (El Salto de los Salados) = Sab; Santa María (before of Adjuntas) SM; Tamasopo = Tam. 

The NMDS based on relative abundance of insects per FFG (Fig. 4) showed that study sites were ordinated on axis one (left to right) in a gradient from low to high number of predators. The second axis (bottom to up) show a gradient of a greater number of filterers, a smaller number of scrapers and gatherers, and lower values of pH and dissolved oxygen. We found important changes in functional feeding groups at some study sites such as Comonfort (Com) where the abundance of gatherers decreased from 82.9% in 1997 to 9.5% in 2014 and the abundance of filterers (from 0.1% to 52.4%) and predators (from 0.2% to 38.1%) increased drastically between 1997 and 2014. At the Ayutla (Ayu) location, the abundance of gatherer decreased from 84.9% to 62.7%, the abundance of filterers and predator increased from 1.9% to 11.6% and from 3.8% to 22.8%, respectively. Fracción Sánchez (FS) showed considerable increase of filterers from 10.6% to 40.8% from 1997 to 2014. This analysis shows a slightly relationship between lower levels of dissolved oxygen and pH with high abundance of filterers and lower abundance of scrapers.

Figure 4 Non metric multidimensional scaling (NMDS) based on relative abundance of aquatic insects per functional feeding group in rivers of two Sub-tropical river drainages including data from 1997 (triangles) and 2014 (circles) Stress: 0.09. Temp = temperature, DO = dissolved oxygen. Ayutla = Ayu; Calvillo = Cal; Canoas = Can; Carpintero = Car; Chuveje = Chu; Comonfort = Com; El Carrizal (Santa María after of Adjuntas) = SMD; El Oasis = EO; El Realito = ER; El Salto = ES; El Xote = EX; Fracción Sánchez = FS; Jalpan = Jal; La Hacienda = LH; La Quemada = LQ; Los Galvanes = LG; Pinihuan = Pin; Presa de Rayas = PR; Presa del Carmen = PC; Puente la Plazuela = PP; Quinta Matilde = QM; Quiotillos = Qui; Rascón = Ras; Río Grande = RG; Sabinolandia (El Salto de los Salados) = Sab; Santa María (before of Adjuntas) SM; Tamasopo = Tam. 

Most of the correlations between the values of PCA with the values from Axis of the NMDS were not significant. The only significant correlation was based on abundance per FFG, using the axis 2 (r_xy = -0.28, p = 0.04). The correlation shows that the increase in water temperature and the decrease in dissolved oxygen is related with more abundance of filterers and less abundance of scrapers; however, in this case the correlation coefficient is very low showing that this pattern is not consistent.

DISCUSSION

The rivers in the Lerma-Chapala and Pánuco river drainages showed symptoms of biological and environmental degradation based on differences in aquatic insect diversity and taxa abundance, FFG and water quality parameters such as dissolved oxygen and pH. The aquatic insect structure and the relationship with water physiochemical variables through the space and time were difficult to interpret at the basin scale. However, our analyzes provided general patterns such as the condition of the aquatic insect fauna and water characteristics in two major subtropical river drainages in east-central Mexico.

It has been shown that the loss of vegetation cover in watershed, mainly riparian vegetation, is a contributing factor in the increase of temperature in freshwater bodies (^{Allan, 2004}; ^{Quinn et al., 1997}). The processes of urbanization are related with degradation symptoms such as the reduction of dissolved oxygen in the water (^{De Jesús-Crespo & Ramírez, 2011}; ^{Ometo et al., 2000}). These patterns are consistent with the decrease of both dissolved oxygen concentrations and pH values, and an increase in water temperature finding in our study and whit the fact that ^{Cuevas et al. (2010)} estimates that the Lerma-Chapala and Pánuco River drainages have lost about 30% and 50% of the vegetation cover respectively, due to the expansion of cultivated pastures, increased agriculture, and urbanization.

Biodiversity measurement has been considered has good indicator of ecosystem stability (^{Maclaurin & Sterelny, 2008}). However, some authors argue that alpha diversity often do not present systematic patterns among habitats, which does not always make them as good indicators of the severity of human impacts (cf. ^{Magurran, 2016}; ^{Pandolfi & Lovelock, 2014}). Another constraint is the fact that all diversity metrics are limited by the ability of researchers to measure them in field, i.e., the community is rarely perfectly measured varying across taxonomic groups, environments, and traits (^{Jarzyna & Jetz, 2016}). The absence of significant differences in the diversity of aquatic insects and the inconsistent patterns in this biological measure in our study, comparing our data from 1997 with data of 2014, are similar to other research, where different gradients of urbanization or river ecosystems degradation were analyzed at the basin scale with no clear responses and patterns in the richness and evenness of aquatic invertebrates (^{Bonada et al., 2006}; ^{Quinn et al., 1997}).

The increase in abundance in Chironomidae, Coenagrionidae, Veliidae, Corixidae and Culicidae families, can be related with the environmental degradation. A positive relationship has been reported between the increase of Chironomidae density with land use changes, such as induced grassland and urban sprawl (^{Jones & Clark, 1987}; ^{Quinn et al., 1997}). These land use changes are generally associated with an increase in water temperature and sedimentation, and low dissolved oxygen concentrations (^{Miserendino et al., 2011}; ^{Walsh et al., 2005}). Chironomids are found in a range of conditions more extensive than any other aquatic insect family; it can exploit an almost complete range of gradient in temperature, pH and oxygen (^{Ferrington et al., 2008}). For this reason, is not surprising that this diverse and opportunistic family showed greater relative abundance in the 2014 when compared to 1997, which was correlated with a decrease of dissolved oxygen and lower values of pH. The larvae of some Odonata are also tolerant and often survive relatively low values of dissolved oxygen and subsist better than many other invertebrates in acidic waters (^{Suhling et al., 2015}), which could explain the increment in individuals of the Coenagrionidae family. Some aquatic Heteroptera, especially Gerromorpha, are good indicators of human disturbance having a high tolerance to eutrophication and acidic waters. Corixidae present a great variation among nutrient and pH tolerance (^{Lytle, 2015}). Accordingly, the increment in these two families of hemipterans, especially Corixidae could be a response of lower pH values. The Culicidae increase in 2014 also can be related with anthropogenic stressors. In this sense, ^{Ribeiro et al. (2012)} reported that environmental change, such as the increase in agricultural areas, irrigation ponds, and the reduction in vegetation cover, tends to increase the abundance of opportunistic species of Culicidae, mainly those species that are considered vectors of human diseases (^{Juliano & Lounibos, 2005}).

Although Baetidae is an Ephemeroptera family very common and dominant in tropical and subtropical rivers, in this study showed an abundance decrease (together with Heptageniidae and Caenidae families) from the sampling performed in 1997 compared with the year 2014. Baetidae and Heptageniidae families are reported to be sensitive to land use changes, such as urban and cropland increase (^{Jones & Clark, 1987}; ^{Li et al., 2012}; ^{Quinn et al., 1997}) because many live attached to boulders and feed on the periphyton (^{Flowers & De la Rosa, 2010}). In general, land use changes can result in an increase in fine sediment deposition, reducing available habitat for benthic organisms (^{Wood & Armitage, 1997}) and resulting in a decrease in periphyton (^{Yamada & Nakamura, 2002}) affecting the establishment and development of families such as Heptageniidae.

In terms of the functional feeding groups (FFG), the increase in gatherers, filterers and predators, and the decrease of scrapers are similar to other studies where a reduction in the of river ecosystem integrity is related with agriculture activities and urbanization processes. In this way, ^{Quinn et al. (1997)} and ^{Friberg et al. (2009)} registered an increment of filterers densities and ^{Md Rawi et al. (2014)} document an increase of predators, filterers, and gatherers in association with environmental degradation. This pattern of increase in collectors (filterers and gatherers) can be an indicator of environmental degradation, because filterers have more availability of suspended particles, and gatherers too with the increase in sediment deposition, which implies, in many cases, more fine particulate organic matter as available feeding resources for these groups. The land use change at basin scale, reduction of riparian vegetation cover, wastewater and pollutants discharges, can cause cumulative and additive effects, which impinges on the river community, changing the habitat, water quality and nutrient amount (^{Allan, 2004}). These disturbances provide favorable conditions to some opportunistic groups such as filterers and gatherers. Sedimentation, for example, restricts the suitability for periphyton and biofilm production (^{Wood & Armitage, 1997}; ^{Yamada & Nakamura, 2002}) limiting the success of scrapers that feed on it.

The pattern of increase in temperature shown in the PCA (Fig. 2), were not accurately related with the aquatic insect assemblages. It coincides with ^{Friberg et al. (2009)} and ^{Buendia et al. (2014)} who found no correlation or strong effect among water temperature and aquatic macroinvertebrate diversity. On the other hand, ^{Jacobsen et al. (1997)}, report a positive relationship between water temperature increase and aquatic invertebrate richness; however, their study focuses on Ecuadorian mountain streams and temperate lowland streams, and refers to a lower range of temperatures, contrasting with those of this study. ^{Li et al. (2012)} report that an increase in water temperature was correlated with the scarcity of Ephemeroptera, Plecoptera and Trichoptera. We also found a correlation among lower pH levels and lower dissolved oxygen concentrations with decrease in aquatic insect abundance, more presence of filterers and a decrement on scrapers. These results appear to differ with ^{Townsend et al. (1983)} who found more abundance of filter feeders in sites with higher levels of pH in temperate streams that are commonly acid, and they attribute this pattern to the greater range of resources in fewer acid streams.

The negative relationship among dissolved oxygen and aquatic insect abundance found in this study, coincides with other studies where it is suggested that the availability of dissolved oxygen restricts macroinvertebrate diversity (^{Jacobsen et al., 1997}; ^{García-Alzate et al., 2010}; ^{Md Rawi et al., 2014}). Moreover, it has reported that with a decrease in dissolved oxygen there is an increase in insect predators (^{Md Rawi et al., 2014}), which is in accordance with our data increment of the mean abundance of Coenagrionidae, Veliidae and Corixidae, a predator’s groups. Only one location (Fracción Sánchez) showed a substantial increase in the number of individuals, explained by the addition of 1017 individuals (>40% of abundance) of the family Culicidae, but a decrease in insect diversity from 1997 to 2014. The positive relationship between reduced water dissolved oxygen and an increase in the filterers, is mainly explained by the great abundance of individuals of Culicidae (categorized as filterers) in 2014, whose members are independent of water dissolved oxygen as they can obtain this resource directly from the atmosphere (^{Clements, 1992}; ^{Wallace & Walker, 2008}).

Across the ecosystems of the world, freshwaters are the most endangered (^{Nel et al., 2009}), and subtropical streams and rivers are especially threatened because are greatly diverse ecosystems usually more than temperate waters (^{Dudgeon, 2008}). Additionally, some pressures are increasing in developing countries (^{Strayer & Dudgeon, 2010}), most located within tropical and subtropical zones (cf. ^{Sachs, 2001}). The rivers of the Lerma-Chapala and Pánuco River drainages have been affected by the combined and cumulative negative effects of human activities (^{Cuevas et al., 2010}).

The results of this study show symptoms of both chemical and biological degradation of these subtropical rivers. The evidence is supported by an increase in water temperature, and a decrease in dissolved oxygen concentrations and lower pH water levels along the space and time. For the aquatic insects, there was an increase in opportunistic and tolerant taxa with a corresponding decrease in sensitive groups. Also, the patterns in the FFG included an increase in the collectors (filterers and gatherers) and a decrease in scrapers. These symptoms reflected loss of river functional processes including energy transformation, nutrient turnover, storage and processing of organic matter, retention and cycling of nutrients, and transportation and deposition of sediments in both river drainages, with the most severe changes occurring in the Lerma-Chapala. In this basin the degradation condition was evident, because the anthropic impacts, including loss of habitat, contamination of waters, increase in sediment deposition and loss of riparian vegetation cover caused by human population growth and agricultural and livestock activities has been noticeably greater than in the Pánuco River basin (^{Cotler-Avalos et al., 2004}; ^{Cuevas et al., 2010}).

Our data reflected the generalized degradation of rivers of two subtropical river drainages in east-central Mexico, which continues unabated and it is evidence of a deficiency in ecosystems conservation strategies in the country. This degradation it’s a risk for the support ecological systems and the public health, because generates the conditions for the proliferation of mosquitoes, capable of transmitting viral infectious diseases such as Dengue (^{Secretaría de Salud, 2001}; ^{Instituto de Diagnóstico y Referencia Epidemiológicos, 2016}), Chikungunya (^{Staples & Fischer, 2014}) and Zika (^{Secretaría de Salud, 2016}), which have occurred in Mexico and other tropical and subtropical zones.

This is the first analysis in Mexico that explores the relationship between aquatic insect assemblages and water quality variables, using this information to indicate degradation levels in two major river drainages through a long-term time scale comparison. This research contributes to the understanding of the trends in the responses of the aquatic biota related to water parameters, providing a framework for the application of historical comparison studies for evaluating the ecological conditions in rivers and to interpret the surrounding landscape impairment in other similar subtropical zones.