Studies of the effects of urbanization on stream ecosystems have usually focused on single metropolitan areas. Synthesis of the results of such studies have been useful in developing general conceptual models of the effects of urbanization, but the strength of such generalizations is enhanced by applying consistent study designs and methods to multiple metropolitan areas across large geographic scales. We summarized the results from studies of the effects of urbanization on stream ecosystems in 9 metropolitan areas across the US (Boston, Massachusetts; Raleigh, North Carolina; Atlanta, Georgia; Birmingham, Alabama; Milwaukee-Green Bay, Wisconsin; Denver, Colorado; Dallas-Fort Worth, Texas; Salt Lake City, Utah; and Portland, Oregon). These studies were conducted as part of the US Geological Survey's National Water-Quality Assessment Program and were based on a common study design and used standard sample-collection and processing methods to facilitate comparisons among study areas. All studies included evaluations of hydrology, physical habitat, water quality, and biota (algae, macroinvertebrates, fish). Four major conclusions emerged from the studies. First, responses of hydrologic, physical-habitat, water-quality, and biotic variables to urbanization varied among metropolitan areas, except that insecticide inputs consistently increased with urbanization. Second, prior land use, primarily forest and agriculture, appeared to be the most important determinant of the response of biota to urbanization in the areas we studied. Third, little evidence was found for resistance to the effects of urbanization by macroinvertebrate assemblages, even at low levels of urbanization. Fourth, benthic macroinvertebrates have important advantages for assessing the effects of urbanization on stream ecosystems relative to algae and fishes. Overall, our results demonstrate regional differences in the effects of urbanization on stream biota and suggest additional studies to elucidate the causes of these underlying differences.
Cities contain an increasing fraction of the global population and have a disproportionate impact on local, regional, and global human and environmental systems (Grimm et al. 2008, Pickett et al. 2008). At local scales, urbanization has been linked to changes in stream habitat (Chin 2006, Roy et al. 2003), water chemistry (Mahler et al. 2005, Gilliom et al. 2006), and hydrology (Poff et al. 2006, USEPA 1997). Altered physical and chemical characteristics of urban streams, in turn, can serve as stressors that affect benthic algae (Coles et al. 2004, Walker and Pan 2006), fish (Meador et al. 2005, Roy et al. 2007), and macroinvertebrates (Roy et al. 2003, Cuffney et al. 2005).
Despite the large and growing body of scientific work documenting the effects of urbanization on stream ecosystems, most urban ecosystem studies have examined the relationship between urbanization and stream ecosystems in a single location. In their reviews, both Walsh et al. (2005b) and Wenger et al. (2009) identified some consistent patterns across studies from single locations. They also noted variability in the relative importance of individual effects common to all locations and the importance of additional effects not observed at all locations. Such variability probably is partially the result of real differences in stream processes among environmental settings. However, different study designs and methods used by different researchers also can affect the results because different suites of environmental variables are measured and responses of different groups of biota assessed. The need to better understand interactions of stressors across urban landscapes (Walsh et al. 2005b, Wenger et al. 2009) requires measurement and analysis of similar suites of environmental variables in different locations. Responses to environmental stressors differ among groups of biota (e.g., algae, macroinvertebrates, fish). Therefore, designers of large monitoring programs have incorporated sampling of multiple taxa (e.g., US Environmental Protection Agency [EPA] Environmental Monitoring and Assessment Program [EMAP], USEPA 2002a; US Geological Survey [USGS] National Water-Quality Assessment [NAWQA] Program, Gurtz 1993).
The USGS initiated a series of urban stream studies as part of the NAWQA program ( http://water.usgs.gov/nawqa/urban/). These studies were based on a common study design, with standardized data-collection and processing methods (Tate et al. 2005). The goal was to compare the effects of urbanization on stream ecosystems in 9 major US metropolitan areas (MAs)—Boston, Massachusetts (BOS); Raleigh, North Carolina (RAL); Atlanta, Georgia (ATL); Birmingham, Alabama (BIR); Milwaukee-Green Bay, Wisconsin (MGB); Denver, Colorado (DEN); Dallas-Fort Worth, Texas (DFW); Salt Lake City, Utah (SLC); and Portland, Oregon (POR) (Fig. 1, Table 1). The least urbanized sites in each watershed were representative of the landscapes being developed in each MA. In each MA, a sampling network was established along a gradient of urbanization intensity to characterize differences in stream environmental conditions and biota relative to urbanization. This design was a space-for-time substitution intended to give insight into the changes in streams as urbanization occurs in a region.
Our study was designed to address 3 questions: 1) How do responses of stream hydrology, geomorphology, and chemistry to urbanization vary based on ecoregion, historic or present land use, or other environmental features of these MAs? 2) How do biotic responses to urbanization vary based on ecoregion, historic or present land use, or other environmental features of these MAs? 3) As watersheds become urbanized, what is the pattern (e.g., linear, threshold) of response by biological communities?
Conceptual framework
We considered our questions within a conceptual framework that explicitly linked human and environmental systems. In each MA, biophysical and human influences shaped a biophysical and human landscape, resulting in a pattern of urbanization characteristic of that MA. Investigators limited the variation in the biophysical factors that might influence stream biotic response by selecting study watersheds from within 1 or 2 level IV ecoregions (Omernik 1987, 2004). Study watersheds were selected to represent the gradient of urban development intensity that existed in each MA (Tate et al. 2005, Cuffney and Falcone 2009).
The effects of urbanization on stream biota in a particular stream reach are mediated by a variety of factors, including hydrology, water quality, geomorphology, temperature, and stream flow (Paul and Meyer 2001, Walsh et al. 2005b, Wenger et al. 2009). The composite effect of these factors on the stream ecosystem has been called the urban stream syndrome (Walsh et al. 2005b). These factors can operate as stressors at several spatial and temporal scales. Sudden changes beyond the normal range of variability in these factors can occur because of processes operating at scales much broader than that of an individual watershed (e.g., regional floods). Such sudden changes also can originate from human decisions, such as removal of dams or changes in water management. Other changes are more chronic (e.g., increased nutrient concentrations, stream flashiness, riparian vegetation loss, stormwater management). The relative importance of individual factors varies in time and space. The responses of fish, invertebrates, and algae to these factors is a function of the exposure and sensitivity of the organisms and communities (Turner et al. 2003).
Methods
We report several major conclusions emerging from the results from the 9 MAs in 3 sections: patterns of urbanization, environmental changes associated with urbanization, and response of biota to urbanization. Readers interested in details of study design and other results not reported here can refer to the original papers on water quality (Sprague et al. 2007, Sprague and Nowell 2008), algae (Coles et al. 2009), and fish (Brown et al. 2009). Results for hydrology, habitat, and benthic macroinvertebrates, which are not yet published, are available in Giddings et al. (2009), as are original data and results of additional statistical tests for all components. Below, we briefly summarize our approach and methods for the results presented in this paper.
Study design
In each MA, 28 to 30 similarly sized study watersheds, typically drained by 2nd- to 3rd-order streams, were selected to represent a gradient of urbanization (Falcone et al. 2007). Variation in natural features (ecoregion, elevation, geology, soil type) in each MA was limited to increase the opportunity to detect a signal of urban influence. The BOS, BIR, and SLC studies were conducted in 1999 to 2000; ATL, DEN, and RAL in 2002 to 2003; and DFW, MGB, and POR in 2003 to 2004. The original studies used indices of urbanization specific to each MA. For the national analysis, we used an urban intensity index based on 3 variables (housing density, % developed land in the basin, and road density) common to all MA-specific indices. The index was scaled to represent urban intensity independently in each MA (MA urban intensity index, MA-NUII) or to allow comparisons of urban intensity across all 9 studies (national urban intensity index, NUII) (Cuffney and Falcone 2009).
The MA-NUII ranged from 0 (least urbanized) to 100 (most urbanized) in each MA. The NUII also ranged from 0 to 100, but was scaled to the NUII value expected at the maximum population density observed over all sites in all MAs. Thus, individual MAs might not have sites with NUII values across the entire range of possible values (0–100). A maximum NUII value <100 in an MA indicates that the maximum level of development in MA was less than the maximum observed across all MAs.
Housing density data were obtained from the 2000 US Census (Geolytics 2001). Infrastructure data were obtained from the Census 2000 Topologically Integrated Geographic Encoding and Referencing database (GeoLytics 2001). Information on land use/land cover was obtained from National Land Cover Data 2001 (NLCD01) data set classification scheme and protocols (USGS 2005) and were compiled with geographic information system (GIS) software.
Hydrology
Continuous stream-stage data were collected with a submersible pressure transducer with an internal data logger and temperature sensor attached. In 6 of the 9 MAs (excluding BOS, BIR, and SLC), a surveyed channel cross section at the location of the pressure transducer or the USGS streamflow gaging station was combined with the continual stage record to calculate cross-sectional area at each stage value to provide a continuous hourly record of cross-sectional area. Cross-sectional area was calculated because it is a more useful measure of hydrologic change than simple stage, and it can be used to calculate a wide range of hydrologic metrics. Continuous discharge was not computed because limited discharge measurements during the short data collection period (<1 y) precluded construction of a stage–discharge rating curve.
Stream stage or cross-sectional area values were used to calculate a number of hydrologic-condition metrics (McMahon et al. 2003) based on the first measurement taken during each hour. Metrics were calculated to summarize overall hydrologic variability; the rate of change-of-flow levels (or flashiness); and the magnitude, frequency, and duration of high and low flows. We present results for 6 selected metrics (Table 2) representative of the patterns observed in our larger data set (see Giddings et al. 2009 for a complete list of variables and data). Calculations were based on a 4- to 12-mo data record that included the biological sampling event, except in DFW. Missing values were not estimated.
Habitat
Stream habitat characteristics were measured at all sites during late-summer low flows with standard NAWQA protocols for wadeable streams (Fitzpatrick et al. 1998). Habitat measurements were collected from the same stream reach as biological information. Reach lengths sampled for habitat and biological assessments (described below) optimally were ≥20× the average wetted stream width and included ≥1 full meander bend (Leopold et al. 1964). A reach length ≥150 m and up to 300 m was preferred for wadeable streams. For some MAs, individual reaches in areas of high urbanization were sometimes <150 m to avoid dams, confluences of other water bodies, and road crossings.
Quantitative and qualitative data on channel geometry, streambed substrate, water velocity, habitat complexity, cover, bank conditions, and riparian conditions were collected at 11 equally spaced transects along the sampling reach (Fitzpatrick et al. 1998). We report 7 variables (Table 3), representative of the patterns observed in our larger data set (Giddings et al. 2009).
Additional information was collected on local channel modifications and natural controls that might affect habitat and geomorphic responses to urbanization. We compiled remarks from field forms, notes on field-drawn reach diagrams and maps, photographs of the reach taken during the time of sampling, topographic maps, and aerial photographs. We then estimated the percentage of the reach with bank stabilization or channelization and number of grade-control structures (weirs, low-head dams, culverts) within the reach. We also estimated, within a distance of 1 reach length upstream or downstream of the reach, the presence or absence of bedrock outcrops in the channel, and presence or absence of depositional features (lateral, mid, or point bars). The degree of confidence in these techniques for quantifying channel modifications and natural controls varied by MA, depending on the level of detail included on field forms and maps and the number of field photographs (FF, written communication). The calculated values represent a minimum for each type of control.
Water chemistry
Stream chemistry conditions were measured by collecting grab water samples once during high flows and once during low flows (Sprague et al. 2007). Water samples were collected at equal-width increments across the stream channel and processed on site with standard USGS protocols (Wilde et al. 1999, 2002). Water samples were analyzed for selected nutrients, ions, and pesticides at the USGS National Water-Quality Laboratory (NWQL) in Denver, Colorado. Most individual pesticides were detected only sporadically (Giddings et al. 2009), so we analyzed several summary variables for pesticide compounds, including the number of detections and total concentrations of insecticides and herbicides.
We deployed semipermeable membrane devices (SPMD; Bryant et al. 2007) at each site (excluding BOS and SLC). SPMDs are passive samplers that concentrate trace levels of hydrophobic organic compounds in the water column (Bryant et al. 2007). At the end of the 6-wk deployment period, compound residues concentrated in the SPMDs were recovered.
Data from 2 assays (an ultraviolet [UV] fluorescence scan and a P450RGS test; see Bryant et al. 2007 for details) were obtained from each site. The UV fluorescence scan provided a semiquantitative screen for polycyclic aromatic hydrocarbons (PAH). The P450RGS test provided a rapid screen for aryl hydrocarbon receptor (AhR) compounds that included polychlorinated biphenyls (PCB), PAHs, dioxins, and furans. A portion of each SPMD sample was sent to the NWQL for identification and quantification of the target compounds (T. Leiker, US Geological Survey, written communication). Results of the bioassays and chemical analyses were normalized for time of exposure because time of exposure directly affects concentrations in the SPMD. The full range of data was not collected in BIR, so those data are not included in our paper (see Giddings et al. 2009 for the available data).
Biological data
Biological data were collected from the stream reaches established for characterizing instream habitat with standard protocols (Moulton et al. 2002), but slightly different protocols for algae (Porter et al. 1993), benthic macroinvertebrates (Cuffney et al. 1993), and fish communities (Meador et al. 1993) were used in BOS, BIR, and SLC. Quantitative samples for the algal and invertebrate assemblages were composed of 5 subsamples from riffles with cobble and gravel substrates or woody snags. These subsamples were then composited into a single composite sample. This sample was called the richest targeted habitat (RTH) sample because, in these streams, riffles or woody debris (when riffles were absent), were presumed to contain the richest assemblage of algae and macroinvertebrates. In ATL and DFW, RTH samples were collected from woody debris and snags by cutting pieces of branches and scraping the surfaces into a collection container. The total area of snag samples was determined by considering each branch as a simple cylinder, calculating the area of each cylinder, and summing the areas. RTH habitat was consistent within an MA, except for a single site in SLC where samples were collected from woody debris because no riffles were available.
Algae
Algal subsamples were collected in the RTH by scraping the upper surface of cobbles collected from 5 riffle areas (3–5 cobbles from each area). Foil templates were made of the sampled surfaces for determination of sample areas. An additional quantitative sample was collected from a composite of 5 samples from depositional targeted habitats (DTH) in areas of the stream that had little or no stream flow. These areas provide a unique habitat for algal assemblages that can include motile or sensitive species not found elsewhere. These samples were collected by inverting a 47-mm-diameter plastic Petri dish, gently pressing it into the sediment surface, sliding a spatula under the Petri dish to trap the sediment, and removing the Petri dish full of sediment. Samples were preserved in 5% buffered formalin and sent to the Academy of Natural Sciences of Philadelphia for identification and enumeration (Charles et al. 2002). Diatoms generally were identified to species. Other algae generally were identified to species or genus.
Invertebrates
Macroinvertebrate subsamples were collected in the RTH with a net (0.5-m wide and 0.25-m deep) with 500-µm mesh, except 425-µm mesh was used in BOS, BIR, and SLC. The composite sample from 5 separate riffle areas had a combined area of 1.25 m2. The mean area sampled from snags in ATL and DFW was 1.4 m2. In addition, a qualitative multihabitat (QMH) sample, which consisted of macroinvertebrates collected from as many habitats in the stream reach as were accessible, was collected. The QMH sample was collected with a dip net (500-µm mesh, except a 210-µm mesh was used in BOS, BIR, and SLC) supplemented with hand-picking of substrates. Sampling effort (measured as time) was apportioned as equally as possible among accessible habitats in the sampling reach. Macroinvertebrate samples were preserved in 10% buffered formalin and sent to the NWQL for taxon identification and enumeration. Invertebrate samples were processed with standard NAWQA protocols (Moulton et al. 2000) for RTH samples (randomized 300-organism count) and QMH samples (fixed processing time designed to maximize the number of taxa enumerated). In general, insects, mollusks, and crustaceans were identified to species or genus. Oligochaetes generally were identified to family. Other groups, such as flatworms and nematodes, were identified at higher taxonomic levels.
Fish
Two electrofishing passes were completed in each stream reach, except in SLC, where a single pass was used. Fishes were identified in the field to species and counted. Species identifications were made or verified by regional taxonomic specialists. Voucher specimens and specimens that could not be definitely identified in the field were preserved in 10% formalin and brought to the laboratory for processing and identification.
Data analysis
We used analysis of similarity (ANOSIM; Clarke and Gorley 2006) to determine if biotic assemblages differed in taxonomic composition among MAs. The ANOSIM statistic R varies from 0 (no difference among groups) to 1 (groups are completely distinct). ANOSIM is conceptually similar to analysis of variance and compares similarity within and among groups.
Biotic assemblages were characterized using nonmetric multidimensional scaling (NMDS). Analyses were based on Bray–Curtis similarities between samples. Ordination axes were rotated (varimax rotation) and relabeled so that the 1st axis (NMDS axis 1) represented the greatest portion of variance and the 2nd axis (NMDS axis 2) represented the 2nd greatest portion of variance. All stress values for 2-dimensional solutions were <0.20, indicating biologically meaningful ordinations (Clarke and Gorley 2006).
Algae NMDS axes were not very responsive to MA-NUII (see Results), so we also examined correlations of MA-NUII with 3 metrics derived from the algal data, including: 1) diatoms tolerant of eutrophic conditions, 2) algae capable of motility, and 3) diatoms sensitive to pollution. Coles et al. (2009) selected these metrics based on responsiveness to environmental gradients, primarily of water quality, in each MA.
We used Spearman rank correlations (rs) as a means of identifying environmental variables or biotic responses that were strongly associated with urbanization. We present results for correlations of MA-NUII with selected environmental variables that demonstrate the relationships observed in the larger data set (see Giddings 2009). Simple linear regression also was used to examine relationships between benthic macroinvertebrate assemblages and urbanization because the correlation results indicated strong and consistent associations.
Because of the large number of correlations involved, we largely limited our interpretation to correlations we considered moderate (0.50 ≤ |rs| < 0.70) or strong (|rs| ≥ 0.70). All such correlations were statistically significant (p ≤ 0.003 and p < 0.001, respectively). These cutoffs were arbitrary, and we used them primarily to focus our interpretations on variables with relatively strong associations with urbanization. These levels corresponded to ∼50% and 25% of the variance explained in linear relations and seemed unlikely to occur by chance. The sign of correlations of NMDS axes with other variables was not important because the ordination only identified gradients in assemblage composition and did not assess assemblage condition against a specific condition, such as a set of reference sites.
Results and Discussion
Patterns of urbanization
Across the 9 studies, the 4 eastern MAs (BOS, ATL, RAL, and BIR) had lower maximum NUII values than the midwestern and western MAs (Fig. 2). The pattern of urbanization in watersheds at the high end of the urban gradient was characterized by lower maximum values for our combined measure of housing density, % developed land in the watershed, and road density in the eastern MAs than in the western MAs. These results can be interpreted in at least 2 ways. First, urbanization might have a more sprawling pattern in the eastern MAs we studied than in the western MAs (Ewing et al. 2002, Harden 2005, Bluestone et al. 2008). Second, stream burial might have occurred more often at high levels of urbanization in the eastern than in the western MAs (Elmore and Kaushal 2008, Roy et al. 2009). If this is true, then the maximum level of urbanization sampled in the eastern MAs would tend to be less than in the western MAs simply because fewer above-ground streams remain to be sampled. More small, above-ground streams might have been preserved in western MAs because these streams are used to distribute irrigation water.
The dominant land cover in watersheds at the low end of the urban gradient (MA-NUII ≤ 10) was presumed to be the land cover being converted to urban uses. Except for the 3 MAs that bordered the Great Plains region (MGB, DEN, and DFW), the dominant land cover in the less-developed watersheds was forested land (Fig. 3). In MGB, DEN, and DFW, the land cover being urbanized was associated primarily with agricultural uses. The presence of forested land in the less-developed watersheds of the eastern and western US does not mean that these areas represented undisturbed conditions because many of these forested lands have been reclaimed from past agriculture (Harding et al. 1998).
Environmental changes
Hydrology
The frequency of high-flow events, measured as the frequency of flow changes >9× the median flow, was moderately or highly positively correlated with MA-NUII in 6 of 9 MAs (Table 2). The only other correlation important in >½ of the MAs was between flow magnitude (measured as skew of the hourly flows) and MA-NUII. All significant correlations were positive, a result indicating that the data became more skewed toward high flows as urbanization increased. Duration of high flows generally decreased with urbanization. Low-flow duration (duration of longest pulse <10th-percentile value) was not correlated with urbanization. Except for BOS, DEN, and POR, each MA had ≥1 strong correlation between a hydrologic variable and MA-NUII (Table 2). DEN and POR had several moderate correlations. BOS had no significant correlations with the hydrologic measures of high or frequent flows.
Our results agree with those of previous studies highlighting increased frequency and magnitude of high-flow events in urbanized areas, primarily because of increased impervious surface associated with urbanization (Schueler 1994, Booth and Jackson 1997). However, our results suggest that this pattern is not universal among MAs. Some of the variability in response might be caused by differences among MAs in stormwater management, which can affect the strength of the relationship between hydrology and urbanization (Booth and Jackson 1997, Walsh et al. 2005b). Reservoirs also can moderate high flows, and the presence of small reservoirs might have been associated with the lack of hydrologic responses in BOS (Coles et al. 2004). Low-flow responses to urbanization are much more variable than high-flow responses. Low-flow responses can depend on a variety of site-specific factors (e.g., landscape irrigation, point-source discharges, leaking pipes, groundwater and surface water withdrawals) (e.g., Lerner 2002, Roy et al. 2009), so the lack of correlations between low-flow responses and urbanization in our studies is not particularly surprising.
Habitat
Few reach-level habitat variables were correlated with MA-NUII (Table 3). The only strong correlation was between channel size and MA-NUII in RAL, a result indicating larger channel size at high levels of urbanization. Responses of physical-habitat characteristics to urbanization have been variable in other studies. Some studies of geomorphic variables documented increased bed erosion, bank erosion, and channel enlargement in urban streams compared to nonurban streams (e.g., Booth and Jackson 1997, Fitzpatrick et al. 2005). Others did not find consistent patterns (e.g., Doyle et al. 2000, Cianfrani et al. 2006).
Geomorphic responses reflect the balance between physical driving and resisting forces over multiple temporal and spatial scales. Geomorphic processes can change over time as sediment sources and transport change during different stages of urban development (e.g., Wolman and Schick 1967, Chin 2006). Incision associated with urbanization in high-gradient headwater streams can cause sedimentation in downstream low-gradient reaches (Doyle et al. 2000). However, local increases in resisting forces, such as outcrops of erosion-resistant bedrock, artificial bank stabilization, culverts at road crossings, and dams might limit expected geomorphic responses (e.g., Booth 1990, Fitzpatrick et al. 2005). In MGB, a strong correlation was observed between amount of urbanization and bank erosion when streams with bank stabilizations were excluded from the analysis (Fitzpatrick and Peppler 2007). Furthermore, stream modifications from centuries past, such as mill dams associated with historical land use, can continue to affect present day stream habitats (Walter and Merrits 2008). Teasing these relationships apart is necessary to understand when and how urbanization affects stream geomorphology.
Water chemistry
Specific conductance was the most consistent indicator of water quality and increased with urbanization at both high and low flows (Table 4). A strong correlation between specific conductance and urbanization was observed in MGB at high flow, but not at low flow, and urbanization did not affect specific conductance in DEN or DFW.
In several MAs where urbanization was occurring on forested land (BOS, RAL, ATL, and POR) total N or NO3+NO2 concentrations increased with increasing MA-NUII (Table 4). No moderate or high correlations were found between urbanization and dissolved N concentrations in the studies where agriculture (MGB and DFW) or shrubland (DEN) were the predominant background land covers. In these MAs, dissolved N concentrations in watersheds with minimal urban development were considerably higher than in RAL, ATL, and POR (Sprague et al. 2007, Giddings et al. 2009). This result suggests that agriculture generates elevated levels of dissolved N that are not affected substantially by urbanization (Sprague et al. 2007). The relationships between dissolved P concentrations and MA-NUII were variable (Table 4). Dissolved P increased with urbanization at high flows in BOS and POR, but increased at low flows only in POR. In MGB, dissolved P declined with urbanization at high flow. Sprague et al. (2007) suggested that the opposite direction of the response in MGB was caused by high inputs of P from agriculture, particularly dairy farming, in areas with little urbanization. Sprague et al. (2007) hypothesized that urban inputs are less than agricultural inputs of P, resulting in a decreasing trend in P concentrations as urbanization of agricultural land occurs. Walsh et al. (2005b) indicated that increases in dissolved nutrient concentrations (N and P) with urbanization were a general pattern in urban streams. Studies have documented declines in N retention (Grimm et al. 2005) and N assimilation (Kaushal et al. 2006) in urban streams. Our results suggest that these relationships are variable and might partially depend on other land uses where urbanization is occurring.
Total herbicide number and concentration generally increased as MA-NUII increased (Table 4) except in MGB and DEN. In MGB, DEN, and DFW, total herbicide concentrations were highest in watersheds with low urban and relatively high agricultural land cover. This result suggests that, similar to its effects on dissolved N, agriculture can elevate herbicide concentrations to levels that overwhelm an urban effect. Conversely, total insecticide number and concentrations increased significantly with increasing urbanization in both forested and agricultural landscapes, particularly at high flows. This result indicates that urban inputs of insecticides were higher than inputs from agricultural or undeveloped land. These results are consistent with observations that herbicides are the most common pesticide in agricultural streams and insecticides are the most common pesticides in urban streams (e.g., Gilliom et al. 2006). Thus, urbanization of forested lands could be expected to increase both herbicides and insecticides but development on agricultural lands would be more likely to affect insecticides because the urban inputs would be additive.
Twenty-five measures derived from the SPMD samples were moderately or highly correlated with MA-NUII in ≥1 of the 6 MAs (no data from BOS, BIR, and SLC; Table 5). These measures included concentrations of 13 PAHs, compounds that are potentially toxic to aquatic life (Eisler 1987). These results are very similar to those obtained using % urban land use as the indicator of urbanization (Bryant et al. 2007, Bryant and Goodbred 2008), and are consistent with previous work on PAHs in urban settings (Hwang and Foster 2006, McCarthy 2006). Urban sources of PAHs include combustion (e.g., vehicle exhaust) and petroleum products (e.g., oil spills, leaking tanks, sealcoats). Mahler et al. (2005) estimated that, in the 4 watersheds they studied in Fort Worth and Austin, Texas, most of the PAHs present could be attributed to runoff from parking-lot sealcoats; however, PAH concentrations can differ greatly depending on the use and type of sealcoat (Van Metre and Mahler 2005). Overall, the results for pesticides and SPMDs support the generalization that toxicant inputs to urban streams increase with urbanization (Walsh et al. 2005b). However, understanding the direct and indirect effects of such inputs on stream biota requires further research.
Changes in biota
Our study was designed to compare responses of biota to urbanization among MAs. One reason for this design was the expectation that species composition of biotic assemblages would vary among MAs. The ANOSIM results verified this expectation for algae (RTH samples, R = 0.72, p < 0.001; DTH samples, R = 0.78, p < 0.001), fish (R = 0.87, p < 0.001), and benthic macroinvertebrates (R = 0.85, p < 0.001). These results are consistent with those of many other studies showing biogeographic differences in species composition across the US (McAllister et al. 1986, Stevenson and Pan 1999).
Algae
Few clear relationships were found between algal assemblage composition and MA-NUII (Table 6, Fig. 4). Algal NMDS axis 1 was moderately correlated with MA-NUII in BOS, RAL, and MGB (Table 6). The data were highly variable and did not provide evidence for resistance of algal assemblage composition to urbanization or any other thresholds (Fig. 4). The initial response (MA-NUII < 50) appeared approximately linear.
At least 2 of the 3 selected algal metrics also were moderately to strongly correlated with urbanization, as characterized by MA-NUII, in these 3 MAs. Coles et al. (2009) analyzed the algal data in more detail in relation to water-chemistry data and found that algal species composition and metrics were more responsive to composite water-chemistry variables derived from principal components analysis than to MA-NUII. In our study, the moderate correlations found in BOS, RAL, and MGB occurred because MA-NUII was strongly correlated with composite water-chemistry variables in those MAs. A moderate correlation between % motile algae and MA-NUII was found in SLC. High percentages of motile algae are common in areas where sedimentation is a problem because motility allows algae to avoid burial.
The limited number of correlations of algae with MA-NUII in our study and the findings of Coles et al. (2009) that algae respond more closely to multiple-variable measures (e.g., axes from principal components analysis) than to individual measures of water chemistry suggest that algae are more responsive to short-term changes in water chemistry from urban runoff or inputs from other land uses, such as agriculture, than to the average conditions of water chemistry that are captured by the measures of land cover and population used to characterize urban intensity. This result is not surprising because algae are considered more responsive than other taxa to short-term perturbations (Lowe and Pan 1996, Duong et al. 2007), and do not integrate environmental effects over seasons or years (Barbour et al. 1999). Sonneman et al. (2001) found that diatom assemblages better reflected water-quality conditions than watershed disturbance. Walker and Pan (2006) noted that diatom assemblages captured a rural–urban landuse gradient, but the gradient also was associated with changes in water quality.
Fishes
Relationships between fish assemblage composition and MA-NUII were variable (Table 6). Strong Spearman rank correlations (rs ≥ 0.70) were found in only BOS, ATL, BIR, and POR (Table 6, Fig. 5). Only 2 correlations were with NMDS axis 1 (BOS and POR). The 2 correlations with NMDS axis 2 (ATL and BIR) indicate that fish species composition can be responsive to urbanization, as characterized by the MA-NUII, even when other environmental gradients are affecting the fish assemblage. Based on species traits, fishes were responding to habitat (pool, riffle, run) and substrate conditions that were independent of MA-NUII in BIR; however, no such relationships were found in ATL (Brown et al. 2009).
The relationship between fish assemblage composition and MA-NUII appeared to be linear in ATL, BIR, and BOS (Fig. 5). In POR, there was a rapid decline at values of MA-NUII <∼20. At MA-NUII < 20, levels of both agriculture and urban development were low; however, at MA-NUII > 20, the combined area of land disturbed by agriculture or urbanization was ≥50% (Waite et al. 2008). Thus, the nonlinearity in the relation might be explained by combined effects of human disturbances rather than a threshold response to urbanization.
The weak correlations between fish assemblages and MA-NUII in MGB, DEN, and DFW might be related to land use before urbanization. In these MAs, urban land uses replaced agricultural land uses, particularly row-crop agriculture and grasslands that support grazing (Fig. 3). Agricultural land uses are associated with changes in fish assemblages (Brown 2000). The effects of prior land use, including agriculture, on aquatic species diversity appear to last for decades (Harding et al. 1998). Thus, urbanization occurring in former agricultural lands might be perceived as having little effect because the fish community already has been altered.
In the MAs with strong correlations between fish assemblages and MA-NUII, sites at the low end of the urban gradient were predominantly forested. RAL did not fit this general pattern; fish assemblages were not correlated with MA-NUII even though sites with low urbanization were forested rather than agricultural. The presence of large reservoirs in the middle of some of the drainage systems sampled in RAL suggested that connectivity might be an issue in this MA (Brown et al. 2009). Reservoirs can affect fish assemblages of associated streams by disrupting migrations of some species and isolating streams from source populations when localized extirpations occur from disturbances, such as drought (Falke and Gido 2006). Reservoirs also can be sources of introduced species that disrupt native assemblages in both upstream and downstream areas.
Fish assemblages also were not correlated with MA-NUII in SLC, where natural vegetation was a mix of shrubland and riparian forest. However, not all sites were included in the fish analyses because of sampling problems or because of doubt regarding whether populations were self-sustaining (Giddings et al. 2006, Brown et al. 2009). Also, interwatershed transfers of water (diversions and augmentation) and water-supply infrastructure were common in SLC and might have affected environmental conditions and dispersal of native and invasive species. Hydrologic alteration can have profound effects on stream fishes and, in combination with species invasions, can severely alter natural assemblages independent of urbanization (e.g., Light and Marchetti 2007, Propst et al. 2008).
Benthic macroinvertebrates
Macroinvertebrate ordination scores (NMDS axis 1) were strongly correlated with MA-NUII in 5 of 9 MAs (Table 6). Results from RTH samples were similar to results from QMH samples in all MAs. Examination of plots of the raw data in MAs with significant relationships indicated that changes in NMDS axis-1 site scores began as soon as the background land cover was removed (Fig. 6). As with algae and fish, no discernible evidence indicated that invertebrate assemblages had any initial resistance to urbanization.
Linear regressions of NMDS axis-1 site scores on MA-NUII were statistically significant for 6 of the 9 MAs based on RTH and 8 of 9 MAs based on QMH samples (Table 7). The strongest regressions were seen in the eastern MAs (BOS, RAL, ATL, BIR) where MA-NUII accounted for 61 to 84% of the variation in the NMDS site scores (Table 7). Invertebrate assemblages in western MAs (SLC, POR) also were significantly related to MA-NUII, but MA-NUII explained only 31 to 49% of the variation in site scores. These simple regression models explained higher proportions of variance in ecological conditions than have models in many other published studies (e.g., Paul et al. 2009). Invertebrate assemblages in the 3 central MAs (MGB, DEN, DFW) showed no significant relationship to MA-NUII based on RTH samples, and assemblages in only 2 MAs (MGB, DFW) showed weak relationships based on QMH samples. The poor relationships between invertebrate responses and MA-NUII in MGB, DEN, and DFW were associated with the conversion of agricultural land (row-crop agriculture and grasslands that support grazing) to urban land.
Key Points
Place matters
The set of hydrologic, physical-habitat, water-quality, and biotic variables responding to urbanization varied among MAs, with the exception that contaminant inputs to urban streams generally increased with urbanization. All 3 biotic assemblages examined in our studies varied significantly in species composition across MAs as determined from the ANOSIM test. Ecoregion or, more generally, regional-scale environmental setting was the primary factor structuring stream biota assemblages across the 9 MAs. The effect of ecoregion on assemblage structure was sufficient to persist even when assemblages were strongly affected by urbanization. These findings, while obvious to ecologists, have implications for managing the impacts of urbanization. Legislators, board members, and planners, whose decisions shape the urban landscape, should not expect the trajectory of physical and biotic responses to urbanization to be similar in all regions. Because physical, chemical, and biological response will vary by region, management approaches that rely on aquatic-life uses and water-quality criteria have to be shaped with an understanding of the physical, chemical, and response gradient of a particular region. Our consistent study design lends strength to this conclusion because variations in results among MAs cannot be attributed to variations in methods or design.
Prior land use matters
At the scale of an individual MA, the most important factor affecting the urban stream syndrome was prior land use. The effects of urbanization were detected least frequently in MAs where agricultural land uses (crops or grazing) were being replaced by urban land uses. In particular, no significant correlations of algae or fish response variables with measures of urbanization were found in DFW, DEN, and MGB. Invertebrate responses to urbanization in these MAs were lacking or weak, depending on the particular response variable used. Analysis of water-quality data suggested little difference in N concentrations with urbanization in MAs with high levels of agricultural land use. Total insecticide concentrations increased with urbanization in all MAs, but herbicide concentrations did not change with urbanization in agricultural MAs. The confounding effect of agriculture on effects of urbanization have been noted in other studies (e.g., Fitzpatrick et al. 2004, Heatherly 2007). Our results do not necessarily mean that urbanization does not affect agricultural landscapes. Rather, they suggest that ongoing or prior disturbance from pre-urban land uses could make effects of additional disturbances from urbanization difficult to detect.
Response patterns
The one general pattern noted across biological assemblages was the absence of any initial threshold effects. That is, there appeared to be no level of urbanization that had no discernible negative effect on fish, algae, or invertebrate assemblages. This pattern was best illustrated by macroinvertebrate responses (Fig. 6), but also was true of algal and fish responses in some MAs (Figs 4, 5). We examined response patterns for linear and threshold forms because recent literature has reported both (Allan 2004). Thresholds have been reported at ∼10% development (e.g., Schueler 1994, Paul and Meyer 2001) or lower (e.g., Ourso and Frenzel 2003, Walsh et al. 2005a). In MAs where invertebrates responded strongly to urbanization, degradation began as soon as the background land cover was disturbed, and assemblages had no apparent ability to resist change during any phase of urbanization.
The apparent absence of a resistance phase indicates either that the assemblages lack the ability to compensate for changes associated with urbanization or that the watersheds we perceived as relatively undisturbed were inadequate reference sites for undisturbed conditions. In other words, our reference watersheds might already have been disturbed beyond the ability of the assemblages to compensate. This possibility raises the question of what condition should constitute the target for urban stream restoration and protection. Should the target reference condition be existing least-disturbed or historical conditions? One approach would be to look for sites outside the immediate study area that might represent the desired least-affected condition. Another approach would be to develop regional models from large data sets that include least-affected sites (Carlisle et al. 2008). Other approaches are also available (Whittier et al. 2007).
The absence of an exhaustion threshold at high levels of urbanization might indicate that a decrease in disturbance intensity could have a positive effect on assemblage condition no matter where it occurs along the urban gradient. Furthermore, remediation of urban effects in areas where agricultural and grazing lands are being converted to urban land could attain conditions better than the best remaining habitat if the effects of both urbanization and the prior land use can be corrected. However, this hypothesis must be tested to determine if such improvements actually are possible.
Indicators of urbanization
Benthic macroinvertebrates appear to have important advantages for assessing the effects of urbanization on stream ecosystems because the response of macroinvertebrate assemblages to a gradient of urbanization was predictable in most MAs (Table 7). Neither algae nor fishes responded as predictably as macroinvertebrates to the urbanization gradient. Algae generally were responsive to reach-scale water-quality conditions and hydrologic conditions that were influenced by urbanization, rather than to watershed-scale measures of urbanization (Coles et al. 2009). Fishes might be particularly vulnerable to loss of connectivity in stream systems arising from disturbances, such as dams and water conveyance systems.
Macroinvertebrates, followed by fishes and algae, are the most commonly monitored taxa in the US (USEPA 2002b). Most programs monitor >1 assemblage. Monitoring taxa other than macroinvertebrates might be important when taxa have important societal value (e.g., anadromous salmonids in urban streams on the US Pacific coast), when a taxon is well suited to monitoring the effects of management interventions focused on reach- rather than watershed-scale environmental factors, or when an endpoint is important in a regulatory or policy context.
Overall
Our results demonstrate the utility of applying a consistent study design in different geographic locations. Comparisons among geographic areas are facilitated when differences in methods and design are minimized. Our studies, which controlled for nonurban gradients, provided results that were much less variable than those from other studies (e.g., Paul et al. 2009, Purcell et al. 2009) and provided enhanced ability to model and interpret responses to urbanization. However, continued progress in understanding the effects of urbanization on stream ecosystems undoubtedly will rely on a combination of studies involving a variety of study designs and approaches and conducted at scales ranging from single MAs to entire countries. The urban-gradient design has limitations (Carter et al. 2009), especially in situations dealing with multiple interacting stressors. Urban streams are complex systems, involving interactions among multiple factors at several scales, and probably cannot be adequately understood by relying solely on pairwise correlations and simple regression analyses to clarify relationships among stressors and response variables. We are currently exploring use of multilevel hierarchical models to understand the relationships between invertebrate responses and watershed-scale (% urban land cover) and regional-scale variables (e.g., climate and background land use) (Nadkami and Shenoy 2004). Further development of this and other modeling approaches (e.g., Bayes Net; Nadkami and Shenoy 2004) will provide improved understanding of these systems.
Acknowledgments
We thank the many past and present biologists and hydrologists with the NAWQA program who diligently collected the data used in our analyses. Peer review from Daren Carlisle (USGS), Anne Brasher (USGS), and several anonymous reviewers greatly improved this manuscript.
Literature Cited
Figure 2.
Comparison of the nationally scaled urban intensity index (NUII) across 9 metropolitan areas and NUII scores standardized to 0 to 100 scale in each metropolitan area (MA-NUII). Metropolitan areas are ordered from east (BOS) to west (POR). See Table 1 for abbreviations.
Figure 3.
Land use/land cover in nonurban watersheds (metropolitan area urban intensity index [MA-NUII] ≤ 10) in the 9 metropolitan areas. Numbers above columns indicate number of watersheds with MA-NUII ≤ 10. Land cover in SLC included a large percentage of shrubland. Percentages do not add up to 100% because minor land covers are not shown (e.g., wetlands, barren land). See Table 1 for abbreviations.
Figure 4.
Plots of site scores from axis 1 of nonmetric multidimensional scaling (NMDS) analysis of algal species composition vs metropolitan area urban intensity index (MA-NUII) for metropolitan areas in which the correlation was moderate (0.50 ≤ |rs| < 0.70) or strong (|rs| ≥ 0.70). See Table 1 for abbreviations.
Figure 5.
Plots of site scores from axes 1 and 2 of nonmetric multidimensional scaling (NMDS) analysis of fish species composition vs metropolitan area urban intensity index (MA-NUII) for metropolitan areas in which the correlation was moderate (0.50 ≤ |rs| < 0.70) or strong (|rs| ≥ 0.70). See Table 1 for abbreviations.
Figure 6.
Plots of site scores from axis 1 of nonmetric multidimensional scaling (NMDS) analysis of macroinvertebrate species composition (qualitative multihabitat sample) vs metropolitan area urban intensity index (MA-NUII) for 9 metropolitan areas. See Table 1 for abbreviations.
Table 1.
Minimum and maximum values for selected watershed characteristics, 2000 population density, mean annual precipitation (1980–1997), and level III ecoregion (Omernik 1987) for 9 metropolitan study areas in the US.
Table 2.
Spearman rank correlations of metropolitan area urban intensity index (MA-NUII) with selected hydrologic characteristics for 9 metropolitan areas in the US (see Table 1 for abbreviations). Moderate correlations (0.50 ≤ |rs| < 0.70) are italicized and strong correlations (|rs| ≥ 0.70) are bolded. All such correlations are statistically significant (p < 0.05).
Table 3.
Spearman rank correlations of metropolitan area urban intensity index (MA-NUII) with selected habitat characteristics for 9 metropolitan areas in the US (see Table 1 for abbreviations). Moderate correlations (0.50 ≤ |rs| < 0.70) are italicized and strong correlations (|rs| ≥ 0.70) are bolded. All such correlations are statistically significant (p < 0.05).
Table 4.
Spearman rank correlations of metropolitan area urban intensity index (MA-NUII) with selected water-quality characteristics for 9 metropolitan areas in the US (see Table 1 for abbreviations). Moderate correlations (0.50 ≤ |rs| < 0.70) are italicized and strong correlations (|rs| ≥ 0.70) are bolded. All such correlations are statistically significant (p < 0.05). NC = not calculated.
Table 5.
Spearman rank correlations of metropolitan area urban intensity index (MA-NUII) with toxicity indices, total number of compounds, and concentrations of compounds obtained from semipermeable membrane devices for 6 metropolitan areas in the US (see Table 1 for abbreviations). Moderate correlations (0.50 ≤ |rs| < 0.70) are italicized and strong correlations (|rs| ≥ 0.70) are bolded. All such correlations are statistically significant (p < 0.05). Correlations were not calculated (NC) for compounds detected at <25% of sites within a metropolitan area or nationally (Bryant et al. 2007). PAH = polycyclic aromatic hydrocarbons, AhR = aryl hydrocarbon receptor.
Table 6.
Spearman rank correlations of metropolitan area urban intensity index (MA-NUII) with the first 2 axes of nonmetric multidimensional scaling ordinations (NMDS) of algal, benthic macroinvertebrate (richest targeted habitat [RTH] and qualitative multihabitat [QMH] samples), and fish assemblage species composition for 9 metropolitan areas in the US (see Table 1 for abbreviations). Three algal metrics (Coles et al. 2009) also were analyzed. Moderate correlations (0.50 ≤ |rs| < 0.70) are italicized and strong correlations (|rs| ≥ 0.70) are bolded. All such correlations are statistically significant (p < 0.05).
Table 7.
Results from linear regressions of nonmetric multidimensional scaling (NMDS) axis 1 scores from richest targeted habitat (RTH) and qualitative multihabitat (QMH) samples against urbanization as characterized by the metropolitan area urban intensity index (MA-NUII) for 9 metropolitan areas in the US (see Table 1 for abbreviations).
[1] 5 E-mail addresses: lrbrown@usgs.gov
[2] 6 tcuffney@usgs.gov
[3] 7 jcoles@usgs.gov
[4] 8 fafitzpa@usgs.gov
[5] 9 gmcmahon@usgs.gov
[6] 10 jjsteuer@usgs.gov
[7] 11 ahbell@usgs.gov
[8] 12 jasonmay@usgs.gov
