Functional changes in littoral macroinvertebrate communities in response to watershed-level anthropogenic stress.

Watershed-scale anthropogenic stressors have profound effects on aquatic communities. Although several functional traits of stream macroinvertebrates change predictably in response to land development and urbanization, little is known about macroinvertebrate functional responses in lakes. We assessed functional community structure, functional diversity (Rao's quadratic entropy) and voltinism in macroinvertebrate communities sampled across the full gradient of anthropogenic stress in Laurentian Great Lakes coastal wetlands. Functional diversity and voltinism significantly decreased with increasing development, whereas agriculture had smaller or non-significant effects. Functional community structure was affected by watershed-scale development, as demonstrated by an ordination analysis followed by regression. Because functional community structure affects energy flow and ecosystem function, and functional diversity is known to have important implications for ecosystem resilience to further environmental change, these results highlight the necessity of finding ways to remediate or at least ameliorate these effects.

Introduction

Widespread anthropogenic modification of the landscape is jeopardizing freshwater ecosystems and their services. Watershed-level anthropogenic stress has negative effects on macroinvertebrate communities in lakes, streams and wetlands [1],[2],[3]. Development in the watershed is characterized by reduction/removal of natural vegetation, increased road density and increased proportion of other impervious surfaces, as well as higher human population density. The effects of development-associated stressors and those associated with agricultural land use on macroinvertebrate communities are mediated through increased nutrient loading, point and non-point source pollution, sediment loads, altered hydrologic and temperature regimes, and habitat destruction and fragmentation in riparian zones and littoral areas [1].

Several functional characteristics (e.g., flow, drag or silt adaptations, respiration and locomotion techniques, feeding habits, voltinism, reviewed in [3]) are considered to be important indicators of the state of an ecosystem and its potential resilience to further anthropogenic modification. In particular, functional diversity (FD) is a critical property of a group of organisms at any scale because increased trait space breadth is likely to be associated with a greater diversity of ecosystem processes and nutrient pathways, which in turn increase resistance and resilience to perturbations [4],[5]. FD can be related to reticulation of the food web, which has important implications for the resilience of food webs as demonstrated in a theoretical study [6]. In terrestrial plant assemblages, greater functional diversity has been shown to maximize resource use in heterogeneous environments, and affect energy flow and ecosystem function (reviewed in [4], [7], [8]), and similar patterns were observed in theoretical studies [9]. However, little is known about the effects of reduced FD in littoral systems. Another important functional trait, voltinism, may have important implications for temporal redistribution of nutrient processing and ecosystem stability, and proportion of univoltine and other longer-lived organisms was shown to decrease with increasing land use stress [10]. Furthermore, changes in the relative abundance of different invertebrate functional groups can alter nutrient processing characteristics [11], [12], potentially triggering cascading effects in higher trophic levels, which may also affect littoral-pelagic and aquatic-terrestrial habitat coupling, by virtue of this group’s central position in aquatic food webs.

Although invertebrate functional trait responses to increasing anthropogenic stress have been described in streams-particularly changes in voltinism and respiration type [3] - those findings may not apply to lake and wetland systems because the effects of watershed-scale stressors are no longer mediated primarily through hydrological alteration and canopy clearing [1]. In addition, some of the more sensitive functional attributes, such as the proportion of semi- and merovoltine taxa and anoxia-intolerant taxa, are already under-represented in lentic systems because these systems routinely experience greater variability in parameters such as dissolved oxygen, temperature, and other factors [13]–[15]. It is important to understand whether functional attributes of littoral macroinvertebrate communities change in response to watershed-scale anthropogenic stress in order to design management and conservation strategies specific to freshwater littoral systems.

To assess the functional responses of lentic macroinvertebrates, we used wetland macroinvertebrate community composition data from the Great Lakes Environmental Indicators project [16] collected across a full gradient of watershed development stress (basin-wide minimum to maximum) to 1) test whether FD and the relative abundance of longer-lived (uni-, semi- and merovoltine) organisms are affected by watershed-level development and agriculture, and 2) investigate the functional changes in invertebrate communities. The a priori prediction was that FD and proportion of longer-lived taxa would decline with increasing development and agriculture due to an overall decrease in taxonomic diversity and reduced habitat availability resulting from the combined effects of direct habitat degradation and pollution commonly associated with these anthropogenic stressors [1]. Due to the large number of taxa and their wide range of physiologic tolerances and habitat requirements, macroinvertebrates can serve as an indicator of the changes impacting other assemblages, as well as a warning signal of changes in littoral food webs.

Materials and Methods

Macroinvertebrates were sampled in 101 coastal wetlands of the U.S. coastline of the Laurentian Great Lakes in the summer of 2002 and 2003 (Fig. 1). This dataset was previously collected for a different purpose, as part of a multidisciplinary effort to identify indicators of anthropogenic stress in the Great Lakes coastal zone (the Great Lakes Environmental Indicators project). No permits were required for macroinvertebrate sampling by the U.S. at the time of the original study and no recognized endangered or threatened invertebrate taxa occur in those coastal wetlands. All appropriate protected area sampling permits were secured by the original study.

Figure 1

A map of the study sites in the coastal wetlands of the Laurentian Great Lakes, overlaid with the land-use stressor gradient for a) agriculture and b) development.

The Laurentian Great Lakes include five glacial till lakes: Superior, Michigan, Huron, Erie and Ontario, located in the temperate part of North America. Lakes range from oligotrophic to eutrophic, with the surrounding land use spanning a gradient from completely unimpacted to highly impacted by anthropogenic activities including development and agriculture. Wetlands were equally distributed among lakes and across a full gradient of anthropogenic stress, previously defined as a composite of five major classes of anthropogenic pressure: agriculture, atmospheric deposition, land cover, human development, and point source pollution. The site selection procedure ensured representation of four geomorphic wetland types, including riverine, barrier-beach protected, and lacustrine coastal wetlands and embayments [17]. All wetlands were hydrologically connected to a Great Lake, and most wetlands had well-developed submerged and emergent macrophyte communities.

Proportion of human development and agriculture (by area) in a wetland’s watershed was derived from the USGS National Land Cover Dataset (2001) for each catchment delineated using ArcHydro with 10-m Digital Elevation Models (see [18] for details). Development was defined to include all residential, commercial and industrial areas, but did not include road density.

Invertebrate sampling and functional metrics

To ensure the most complete sampling of the different habitats present within a wetland, macroinvertebrates were collected from all representative near shore land-use and shoreline (e.g., sand beach, cobble beach, lawn, etc.) zones in each wetland. In each land-use-shoreline zone, two transects were extended perpendicular to shore, and macroinvertebrates were collected using 30-second D-frame dipnet sweeps at 0.25 and 0.75-m water depths along each transect. Sweeps were done through the water column from the bottom to the surface, in a forward direction parallel to the shore, regardless of vegetation type or presence. Samples were rinsed in a 250 µm sieve net or bucket to remove fine particles, and preserved in Kahle’s solution for laboratory processing. In the laboratory, macroinvertebrates were identified to the highest possible resolution (genus for most insects) using the most current keys available at that time [19], [20]. Data from the two transects in each zone were first averaged by depth, then averaged across zones to achieve site-level data. After taxonomic identification and sample averaging were completed, traits were assigned to each taxon using the latest reviews [14], [19], [21]–[27] and expert judgment (see Table S1 for details).

FD was measured using Rao’s Quadratic entropy (Q), which accounts for the relative abundances of species and for the functional differences between species by measuring differences between two randomly selected individuals with replacement (see [28], [29] for formulas and performance evaluation). It is closely correlated with the index of functional dispersion based on the centroid-distance approach [30], although it has also been demonstrated that this metric can be conservative under some scenarios due to negative covariation with species richness [29]. This analysis was performed with trait variables including trophic status, feeding mechanism, locomotion and primary and secondary functional feeding groups (Table S2), because those were the traits for which information was most consistently available across all encountered taxa. All traits were combined into a single Q-space.

The voltinism measure was expressed as the proportion of a sample comprised of taxa with long-lived aquatic phases (i.e. the proportion of individuals in a sample that were uni-, semi or merovoltine) to all other taxa. Voltinism is defined only for insect taxa; this life-history information was available for 77 taxa and unavailable for the remaining 85 insect taxa (of the total of 222 insect and non-insect taxa). This was mostly due to incomplete knowledge, as voltinism is one of the most difficult traits to describe, and to a lesser extent due to the level of resolution (variable life history at the species level, but identification to genus level) and differences in life history across a large geographic range of a taxon. Although all studies looking at voltinism unavoidably have this limitation, we chose to analyze this trait because it is a very important indicator of the state of macroinvertebrate assemblages and there is no a priori reason to suspect that the effect of this missing information is directional and not conservative on our ability to understand how longer-lived insects are affected by land use. Fifty-eight taxa in our samples were uni-, semi- or merovoltine, primarily belonging to Odonata and several genera of Trichoptera and Ephemeroptera, along with a few rarely occurring groups.

Statistical Analyses

Multiple regression (MR) analyses were used to test for relationships between watershed stressor variables and voltinism and FD. Percent development was log-transformed, whereas percent agriculture was not transformed to satisfy the assumptions of regression residuals distribution, as tested using Shapiro-Wilk tests. To address the possibility of the confounding effect of sampling across large geographic scales, we conducted additional multiple regression analyses using latitude as a third predictor for each of the response variables. Latitude was significantly, but weakly correlated with agriculture (linear regression permutation p = 0.010, R2 = 0.06), but not with development (p = 0.26). Functional community structure was summarized using Principal Components Analysis (PCA) conducted on the log-transformed relative abundances of 12 functional traits. The resulting site Principal Component scores for each factor were regressed against the predictor variables as described above. Analyses were done in R 2.12.2 (R Development Core Team, Vienna, Austria), using packages vegan and lattice. Package relaimpo was used to calculate the “lmg relative importance metric” [31] for predictors in the multiple regression models, which produces averages of sequential sums of squares over all orderings of regressors (see [31] and package documentation for details). Although computationally intensive, the lmg procedure has been recommended because it decomposes R2 into non-negative contributions and accounts for direct effects as well as adjustments for other regressors in the model [32]. Natgrid, a two-dimensional random data interpolation package, was used to resample data to a regular grid and produce contour plots. Natgrid, based on nngridr package [33], implements a natural neighbour interpolation method and uses a weighted average method. The package was implemented through matplotlib [34], plotting library for the Python programming language.

Results

Watershed development exerted stronger effects on all functional response variables than did watershed agriculture (Fig. 2). Rao’s functional diversity was significantly reduced by both development and agriculture (Fig. 3, MR F2, 97 = 9.60, adjusted R2 = 0.16, p = 0.003 and 0.004 for development and agriculture, respectively; see Table S3 for details and all regression coefficients). Development was the more important driver of this relationship (Fig. 2, Fig. S1a). To illustrate the extent of this effect, sites with more than 10% development had a 25% reduction in functional diversity compared to the less developed sites. The latitudinal predictor did not contribute significantly to the MR (p = 0.58), and its presence had almost no effect on the relative importance of the two stressors, confirming the lack of confounding effects (Table S3).

Figure 2

Relative importance of the two watershed stressors and latitude in explaining variation in the functional attributes of macroinvertebrate communities, as proportional lmg contribution.

Note the overriding importance of development over agricultural stress and the lack of significant geographic confounding. Asterisks indicate predictors that were significant in the multiple regression (p<0.05).

Figure 3

Macroinvertebrate functional diversity (Rao’s Q) as a function of development and agricultural stress in the watershed.

Note that maximum values of the functional trait are observed mostly below a certain proportion of development (around mid-point of the y-axis, corresponding to 10% untransformed % development), indicating the over-riding contribution of that stressor; and if high values are observed at the higher levels of development, those occur only at sites with minimal agriculture. *[figure footnote] Contour lines divide the figure into a region where values are higher than they are on the contour line itself and a region where they are lower; values change across the line but not along the contour line, and the gradient is larger where contour lines are packed closer together.

Relative abundance of those longer-lived taxa declined with increasing development in the watershed (MR F2, 97 = 11.17, adjusted R2 = 0.17; Fig. 4, p<0.001 for development, also Fig. S1b). The weak (Fig. 2) but significant negative effect of agriculture observed in the two-factor model (p = 0.017) was no longer significant (p = 0.062) when latitude was included as a predictor (see Table S3 for model details and regression coefficients).

Figure 4

Relative abundance of long-lived taxa (uni-, semi- and merovoltine) as a function of development and agricultural stress in the watershed.

The overall functional community structure was weakly but significantly affected by the extent of development in the watershed. The first PC axis, explaining 32% of the variation in functional community structure (eigenvalue = 3.74), was driven primarily by burrowers and filterers on one end vs. scrapers and clingers on the other end (see Table S4 for loadings). This axis was significantly positively correlated with percent development (MR F2, 97 = 5.82, adjusted R2 = 0.11; p = 0.004) but not with agriculture. In other words, the proportions of positively-loading groups, including burrowers and filterers, were positively correlated with greater amounts of watershed development, whereas scrapers and clingers were negatively correlated. This axis score tended to be negatively correlated with agriculture, although not significantly (p = 0.053), with this trend disappearing when we added latitude as a predictor, similar to the pattern observed for other response variables. The second PCA axis (eigenvalue = 3.10, percent variation explained = 19%) was mostly driven by omnivores (Table S4). This axis was weakly negatively correlated with percent agriculture (p = 0.011, R2 = 0.06).

Discussion

This study demonstrates that macroinvertebrate functional diversity and abundance of longer-lived taxa in Laurentian Great Lakes coastal wetlands decreased significantly with greater watershed-level stress. Although macroinvertebrate communities are influenced by factors acting at many spatial scales, as has been demonstrated in a number of stream studies (e.g. [3], [35], [36] but see [37] for lake margins), our results indicate that the negative effects of watershed-level development were sufficiently robust to be detected against the background of a strong geographic gradient, the presence of additional stressors and a wide range of local habitat features. Compin and Céréghino [38] used self-organizing maps to show that, in human-modified landscapes, effects of stream watershed-scale stressors over-rode the influence of natural physical factors. It also means that large-scale stressors can produce detectable changes in ecosystem function, at least with respect to macroinvertebrate functional traits. The effect of watershed development was significant, but not highly predictive, reflecting the hierarchical complexity of the underlying factors affecting macroinvertebrate assemblages and leading to the high variability seen in many such datasets [15], [35].

A reduction in macroinvertebrate FD examined in this study translates into reduced community-level variation in foraging mechanisms and locomotion/substrate relations, which is likely to lead to significant alterations in food web structure and energy flow in those coastal wetland systems. Reduced FD is also likely to translate into decreased asynchrony of taxa responses to environmental perturbations (because more functionally similar taxa would have more synchronized responses to changes in resource abundance), which in turn has been proposed to be one of the key mechanisms reducing ecosystem stability in a theoretical modeling study [39].

The greater relative importance of development vs. agricultural stress is not surprising considering previously reported macroinvertebrate community thresholds at very low levels of development [40]–[44]. Similarly, a previous study found that development in the watershed was the best predictor of several macroinvertebrate metrics, including Aeshna abundance [37]. However, further study is needed to identify the more proximal factors driving this relationship, specifically in littoral systems. For instance, decline in longer-lived organisms could be mediated through the greater cumulative effects of development-associated stressors such as point-source pollution, in particular organic pollution, on these organisms, or through the destruction of biogenic (macrophyte) complexity, which has been shown to support greater abundance of several types of predators [45] and likely sustains more stable, and thereby more diverse, predator-prey populations [46]. It has proven difficult to elucidate repeatable patterns in these coastal wetlands, possibly due to their highly variable hydrology [35], [47], e.g. half a meter or more changes in water level on diel and annual basis due to fetch and climate change. Despite this, there is considerable interest in finding community metrics responsive to anthropogenic stress, and several such metrics have been found including Ephemeroptera-Trichoptera-Odonata richness-based metrics [48].

The observation of shift in functional community composition towards greater relative abundance of burrowers and filterers is consistent with other studies and has been previously related to increased sedimentation in streams [49]. Yet flashier hydrology and increased siltation, which are cited as the most common causes of anthropogenic changes in stream assemblages [1], are less likely to be the same factors responsible for observed changes in littoral assemblages, and lake-specific environmental variables responsible for these effects need to be investigated. Several functional groups (clingers, burrowers and insect filter-gatherers) were previously demonstrated to be affected by the extent of anthropogenic development in this system; however, responses were complex and dependent on several predictor variables as well as land use predictors’ buffer sizes [37].

Trait responses are generally less frequently reported in aquatic studies than trends in diversity and abundance [50]. Effects of watershed land use on FD have garnered even less attention, particularly in non-stream ecosystems. For lotic macroinvertebrates, it has been shown that trait type changed [51] and trait diversity decreased [52] with increasing agricultural land use. The latter study demonstrated that this effect was detectable both at the watershed scale and local patch scale, and was related to increasing sedimentation [52]. Significant effects on these important functional variables (e.g., trait diversity) emphasize the need to further investigate functional responses to anthropogenic stress in lentic macroinvertebrates, the mechanisms potentially underlying these responses and the local factors that may mitigate those effects. Considering additional local factors as well as integrating traits with explicit consideration of trait linkages [53] may increase the predictive capability of future littoral trait models. However, what is of greater interest for future studies is that this uncertainty potentially reflects trait responses to smaller-scale habitat factors, which, if amendable to manipulation or restoration, could be used to ameliorate effects of watershed-scale land use.

In summary, we observed statistically and biologically significant reduction in macroinvertebrate FD and abundance of longer-lived taxa, and noticeable differences in functional community structure associated with increasing proportion of development in the watershed. These findings, along with the previously-observed threshold changes in macroinvertebrate community composition [44], show that lentic fauna exhibit significant functional changes associated with greater levels of watershed-scale anthropogenic stress.

Univariate relationships with the development stressor.

(DOCX)

Click here for additional data file.

List of observed taxa with taxonomic information and trait values.

(XLSX)

Click here for additional data file.

List of possible states for each trait category used in calculating Rao’s Q.

(DOCX)

Click here for additional data file.

Details of multiple regression analyses.

(DOCX)

Click here for additional data file.

Trait loadings for the Principal Component Analysis.

(DOCX)

Click here for additional data file.