Research Article

Groundwater declines are linked to changes in Great Plains stream fish assemblages

, Keith B. Gido, Jeffrey A. Falke, Kurt D. Fausch, Harry Crockett, Eric R. Johnson, and John Sanderson

  1. aDepartment of Biology, Tennessee Technological University, Cookeville, TN 38505;
  2. bDivision of Biology, Kansas State University, Manhattan, KS 66506;
  3. cUS Geological Survey, Alaska Cooperative Fish and Wildlife Research Unit, University of Alaska, Fairbanks, AK 99775;
  4. dDepartment of Fish, Wildlife, and Conservation Biology, Colorado State University, Fort Collins, CO 80523;
  5. e Colorado Parks and Wildlife, Fort Collins, CO 80526;
  6. fBiology and Conservation Programs, Westar Energy, Topeka, KS 66601;
  7. g The Nature Conservancy, Boulder, CO 80302

See all Hide authors and affiliations

  1. Edited by B. L. Turner, Arizona State University, Tempe, AZ, and approved May 22, 2017 (received for review November 23, 2016)

Loading

Significance

Nature and society depend on groundwater to sustain aquatic ecosystems and human livelihoods, but local and regional groundwater supplies are dwindling where human water extraction exceeds aquifer recharge. Although groundwater depletion is a global problem, ecological consequences for aquatic species such as fishes are rarely examined. We demonstrate that more than half a century of groundwater pumping from the United States High Plains Aquifer has been associated with collapses of large-stream fishes and expansion of small-stream fishes where hydrologic conditions were altered most. Projections indicate that these habitats will continue to shrink over the next half-century if groundwater pumping practices are not modified. Our findings highlight a mechanism for biotic homogenization with global implications given the worldwide extraction of groundwater.

Abstract

Groundwater pumping for agriculture is a major driver causing declines of global freshwater ecosystems, yet the ecological consequences for stream fish assemblages are rarely quantified. We combined retrospective (1950–2010) and prospective (2011–2060) modeling approaches within a multiscale framework to predict change in Great Plains stream fish assemblages associated with groundwater pumping from the United States High Plains Aquifer. We modeled the relationship between the length of stream receiving water from the High Plains Aquifer and the occurrence of fishes characteristic of small and large streams in the western Great Plains at a regional scale and for six subwatersheds nested within the region. Water development at the regional scale was associated with construction of 154 barriers that fragment stream habitats, increased depth to groundwater and loss of 558 km of stream, and transformation of fish assemblage structure from dominance by large-stream to small-stream fishes. Scaling down to subwatersheds revealed consistent transformations in fish assemblage structure among western subwatersheds with increasing depths to groundwater. Although transformations occurred in the absence of barriers, barriers along mainstem rivers isolate depauperate western fish assemblages from relatively intact eastern fish assemblages. Projections to 2060 indicate loss of an additional 286 km of stream across the region, as well as continued replacement of large-stream fishes by small-stream fishes where groundwater pumping has increased depth to groundwater. Our work illustrates the shrinking of streams and homogenization of Great Plains stream fish assemblages related to groundwater pumping, and we predict similar transformations worldwide where local and regional aquifer depletions occur.

  • ecology
  • conservation
  • freshwater
  • Great Plains
  • fishes

Worldwide, irrigation accounts for 90% of human water use and is sustained by the annual pumping of 545 km3 of water from global groundwater sources (1). In North America, major aquifers are important sources of water for 60% of land equipped for irrigation. One of these aquifers, the High Plains Aquifer in the Great Plains, is the single greatest source of groundwater and supports $35 billion [2007 US dollars (USD)] in US market value of agricultural products (2). The total area irrigated with groundwater from the High Plains Aquifer was 8,500 km2 in 1949 when large-scale pumping began, increased to 55,000 km2 by 1980, and reached 63,000 km2 by 2005 (3). Sustaining this level of agricultural productivity will depend on continued extraction of water from the High Plains Aquifer, but the aquifer is experiencing substantial declines in storage (4, 5). Groundwater extraction from the aquifer is occurring faster than recharge (6) and has resulted in the depletion of 410 km3 of stored groundwater, a volume equal to 85% of the water in Lake Erie, North America (5). Long-term depletion of the High Plains Aquifer has caused water tables to drop by more than 50 m in some portions of the Great Plains whereas redistribution of water through surface canals and lateral subsurface flows has contributed to rising water tables in other locations (7, 8). Rates of groundwater depletion caused by pumping from the High Plains Aquifer are similar to those measured on portions of every continent except Antarctica, suggesting that the Great Plains is a microcosm for the effects of global groundwater pumping on the hydrologic cycle (9).

Pumping from the High Plains Aquifer has also caused surface water to decline in streams of the western Great Plains (10, 11). In the Republican River basin in southwestern Nebraska, depleting groundwater has caused streamflow declines at 70% of US Geological Survey stream gauges (12). In the Arikaree River in eastern Colorado, mean annual discharge declined 60% from the period 1932–1965 to the period 1966–2006 (13). In the western third of Kansas, historically perennial streams are now ephemeral or permanently dry, and groundwater extraction combined with surface diversions permanently dry the Arkansas River just downstream of the Colorado border (14). These losses of surface stream flow are caused by declining water tables below river beds, so much so that they no longer supply groundwater to channels. Although rainstorms can produce floods in this region, base flows are sustained by groundwater (15) and are thus sensitive to water table fluctuations close to the ground surface (16–18). As depth to groundwater increases beneath streambeds, streams become decoupled from the aquifer, and aquatic habitat becomes intermittent or completely dry (6). However, changes in depth to groundwater are not consistent across space and time so broad-scale ecological consequences of decoupling streams from the High Plains Aquifer are sparsely studied despite widely accepted linkages between groundwater and the ecology of surface streams (6, 19, 20).

Although stream drying is a natural and pervasive process in Great Plains streams, the combination of groundwater pumping and habitat fragmentation by diversion dams, reservoirs, and other anthropogenic barriers now prevent many fishes from finding refuge from increased drying (21–23). Across the western Great Plains, stream fish biodiversity declined during 1950–2010 as the area of irrigated land and pumping for irrigation increased, diversion dams and reservoirs fragmented surface habitats, and stream flows diminished due to reduced groundwater input (13, 24–26). The combined effects of increased drying and surface fragmentation "ratchet down" Great Plains fish diversity because barriers prevent recolonization of dry stream segments once precipitation restores flow (27). However, despite the implied connections between groundwater levels and stream fish communities, few studies have explicitly linked groundwater depletion with stream fish community change (6, 13). Conservation of stream fishes in groundwater-dependent ecosystems in the Great Plains and worldwide requires a better understanding of how groundwater pumping reduces streamflow, and how reduced streamflow affects stream fish assemblages (9, 19, 28). Ecological change associated with groundwater depletion represents a rarely studied mechanism of biotic homogenization that ultimately results in loss of freshwater natural resources (29). The ecological and evolutionary consequences of reduced taxonomic, functional, and molecular diversity affect global human livelihoods because biotic homogenization ultimately compromises ecosystem function, services provided, and resiliency to disturbance while simplifying food-web structures and increasing community susceptibility to species invasions (30). Collectively, these consequences of freshwater resource loss represent human diminishment of the very resource required for long-term persistence (31).

We used a retrospective approach to assess the effects of 60 y of historical depletion of the High Plains Aquifer on fish assemblages inhabiting large and small streams in the western Great Plains and combined this approach with prospective modeling to develop baseline expectations for fish assemblage change through 2060. We used historical groundwater-level data from 3,431 observation wells distributed across the upper Kansas River Basin of Colorado, Kansas, and Nebraska (Fig. 1) to create annual interpolated surfaces representing depth to groundwater during 1950–2010. These interpolated depths were subtracted from a digital elevation model for the region (with a resolution of 0.5-km cells) to create spatially continuous estimates of depth to groundwater across the study area. We then used stream locations derived from the National Hydrography Dataset Version 2 Plus (32) to classify stream segments into two categories regarding coupling with the aquifer: "coupled" (i.e., depth to groundwater <1 m under stream location) or "uncoupled" (depth to groundwater >1 m under stream location) (SI Methods). We classified all segments across the region for each year and summed the length of coupled segments to quantify spatiotemporal change in habitat for stream fishes. We used historical fish collections (Table S1–S5 and Fig. S1) made across the region during 1950–2010 to establish relationships between length of coupled stream segments and species occurrences using generalized regression models. Based on these retrospective relationships, we forecasted (2011–2060) change in fish assemblage structure for large streams (order ≥4) (33) and small streams (order <4). Our analysis focused on two spatial scales: the regional scale combining 24 subwatersheds in the upper Kansas River basin and the subwatershed scale that focused on 6 subwatersheds arranged across a gradient of groundwater decline.

Fig. 1.

  • Download figure
  • Open in new tab
  • Download powerpoint

Fig. 1.

Study region in the Great Plains of Colorado, Kansas, and Nebraska, United States showing streams over the High Plains Aquifer, groundwater observation wells, dams, and US Geological Survey streamflow gauges. Hydrologic units (shaded areas) define six subwatersheds in which detailed analyses of hydrology and fish communities were conducted. The Inset illustrates the upper Kansas River Basin boundary (dark shading) within the broader extent of the High Plains Aquifer (blue shading). Data are from refs. 32, 53, and 56.

View this table:

  • View inline

Table S1.

Species, common name, stream size classification, native status, and parameter estimates (SE), Z value, and P value for generalized linear models predicting species capture probability as a function of length of stream coupled with the High Plains Aquifer for the region and the Arikaree River subwatershed

View this table:

  • View inline

Table S2.

Species, common name, stream size classification, native status, and parameter estimates (SE), Z value, and P value for generalized linear models predicting species capture probability as a function of length of stream coupled with the High Plains Aquifer for the North Fork Republican River and South Fork Republican River subwatersheds

View this table:

  • View inline

Table S3.

Species, common name, stream size classification, native status, and parameter estimates (SE), Z value, and P value for generalized linear models predicting species capture probability as a function of length of stream coupled with the High Plains Aquifer for the Upper Republican River and Frenchman Creek subwatersheds

View this table:

  • View inline

Table S4.

Species, common name, stream size classification, native status, and parameter estimates (SE), Z value, and P value for generalized linear models predicting species capture probability as a function of length of stream coupled with the High Plains Aquifer for the Harlan County Reservoir subwatershed

View this table:

  • View inline

Table S5.

Species, common name, stream size classification, native status, number of occurrences, and the last year reported for fishes excluded from analysis because of rare occurrences (i.e., <30) or because identification was not to the species rank

Fig. S1.

  • Download figure
  • Open in new tab
  • Download powerpoint

Fig. S1.

Mean (SD) stream order for species occurrences in the Kansas Aquatic GAP Database for fishes characteristic of small (<fourth-order mean, blue) and large (≥fourth-order mean, red) streams ("large" and "small" stream classifications based on ref. 37).

SI Methods

In this section, we provide detail for the spatial analysis, fish assemblage analysis, and validation used to characterize changes in depth to groundwater across the upper Kansas River Basin in Colorado, Kansas, and Nebraska.

Spatial Analysis.

We used data from observation wells archived in the US Geological Survey National Water Information System (56) to analyze changes in depth to groundwater through space and time. Methods to measure distance from ground surface to the water table are in ref. 57 and are assumed to be accurate to 0.30 m (8). We selected wells with at least 10 measurements of depth to groundwater during November through February to provide sufficient data for modeling (4) and to avoid intraannual fluctuations caused by pumping during the summer growing season (7). We used geographic information systems to analyze changes in groundwater elevation (i.e., groundwater elevation = ground surface elevation – depth to groundwater) during 1950–2010. This period included an estimate of conditions before pumping began (i.e., 1950) (58) and continued through the most recent decade (2010). Our geographic information system modeling framework tracked changes in groundwater elevation by interpolating among wells across the region to generate a raster layer (see ref. 3 for detailed methods). Not all wells were measured annually so we used the linear change in depth to groundwater through time to predict missing years (4, 36). We then used the "topo to raster" function in the Spatial Analyst Extension of ArcMap 10.0 (Environmental Systems Research Institute) to convert point data for wells into a raster layer for each year. Interpolated groundwater elevations were subtracted from a digital elevation model (1/3 arc-second or ∼10 m resolution; vertical accuracy ±1.55 m) (59) from the National Elevation Dataset (60). For spatial analyses, we used an Albers equal-area conic projection, North American Datum of 1983 (NAD 1983), and raster grid cell size of 0.5-by-0.5 km to match previous analyses of groundwater change in the Great Plains (8).

We used our raster model of depth to groundwater to estimate where groundwater from the High Plains Aquifer reached the surface and was coupled with streams (hereafter, "coupled" streams), versus those places where groundwater did not reach the surface and streams were dry ("decoupled"). We subtracted the modeled groundwater elevations from the ground surface elevation in the digital elevation model ("raster calculator" function) to create a 0.5-km grid cell raster map of depth to groundwater across the study area. We then overlaid streams from the National Hydrography Dataset plus version 2 (32), excluding those over perched or other regional aquifers ("extract by mask" function) (8).

Within the Geospatial Modeling Environment (61), we calculated the length-weighted mean depth to groundwater along each stream segment for each calendar year. We assumed that groundwater was coupled with stream channels if the water table was within 1 m of the land surface estimated from the digital elevation model at the exact location of streams. Although previous studies suggest aquifers exchange water with the surface when depth to groundwater is 0–5 m (16–18), our threshold of 1 m represents an estimate of aquifer–stream coupling within the margins of error for measurements of aquifer depth and the digital elevation model. We then summed the length of coupled streams across the study area for each year to estimate the lengths through time for all streams, as well as for streams of orders first to sixth (33). For this analysis, we held stream order constant according to segment attributes from the National Hydrography dataset rather than reduce stream order as headwaters and tributaries became decoupled from the aquifer. We predicted future groundwater depths, groundwater elevation, and stream coupling for 2011–2060 by assuming a linear relationship between time and depth to groundwater within observation wells (4, 36). We used data from 1980 to 2010 only because the rate of groundwater decline had slowed compared with 1950–1979 (13, 36). The period to 2060 marks the estimated usable life of portions of the aquifer (5) and the life of regional reservoirs given projected drying (36), and is similar to periods modeled for other ecological analyses (13).

Fish Assemblage Analysis.

We assessed change in fish assemblage structure using historical fish assemblage collections distributed across the region. Data were obtained from the Colorado Division of Parks and Wildlife (collections during 1977–2010, n = 185), the Kansas Aquatic Gap Program (1950–2010, n= 493), and the Nebraska Game and Parks Commission (1950–2010, n = 262). Collections were defined as all fish species captured from unique sites on unique dates, excluding collections with only one species, and retaining species for analysis if they occurred in at least 30 collections. We used presence/absence of each species in all collections taken during a year as the response variable and the annual modeled length of stream coupled with the aquifer for the same year as the predictor variable. We fit binomial logistic curves to occurrence data to estimate capture probability (i.e., the probability of being present and being detected; range 0–1) as a function of the length of stream coupled with the aquifer for the period 1950–2010. We initially included the number of dams as a second predictor variable, but, because the number of dams was constant through time across most subwatersheds and no additional explanatory power was added, we excluded dams from the final models. Given the exploratory nature of our work, we assessed model significance at α = 0.10 and did not adjust for experiment-wise error. Models were fit for each of 28 species that occurred in the region and subwatersheds (Tables S1–S5). We classified fishes according to stream order using occurrence data from the Kansas Aquatic GAP Program. For this analysis, we used the georeferenced spatial location of fish captures and stream order for the segments from which captures were made to calculate the mean and SD of stream orders used by each species (Fig. S1). This approach provided an objective measure of stream size occurrences for all members of the fish assemblage but was not intended to be a rigorous assignment of habitat suitability or preference on a species-by-species basis. We relied on ref. 37 to define small ("headwater") streams as <fourth order and large streams as ≥fourth order.

Long-term change in fish assemblage composition was assessed using average capture probabilities for small-stream and large-stream fishes. We used predicted values from each logistic regression model to estimate capture probability for each species at both the region and subwatershed scales for the period 1950–2060 using predictions of coupled stream length described earlier (Spatial Analysis). We calculated the mean and 95% confidence interval across all small-stream (<fourth-order) and large-stream (≥fourth-order) species for every year between 1950 and 2060 to illustrate long-term assemblage-scale change.

Depth to Groundwater Estimate Validations.

We validated our spatial predictions of stream coupling in three ways. First, we compared our interpolated changes in depth to groundwater between 1950 and 2011 (i.e., 1 y after the last well observations included in our analysis) with values published by the US Geological Survey (8). This comparison was necessary because different data were used to estimate change in depth to groundwater over time. For example, ref. 8 used a series of "pre-development" (circa 1950) and "post-development" (circa 2011) observation wells distributed across the region to estimate change in depth to groundwater. In contrast, our approach used a smaller subset of these wells with 10 repeated measurements through time, and we projected change for 2011–2060 using linear regression and observations made only during 1980–2010. Spatially explicit comparisons of our modeled change and those reported by the US Geological Survey showed geographically similar changes in depth to groundwater across the study region (Fig. S2). Second, we compared the magnitude of change in depth to groundwater across the study region with values reported by ref. 8. For this comparison, we used the geospatial modeling environment to assign depth to groundwater values to each stream segment over the US Geological Survey layer and our layer for 2011. We then regressed the observations from ref. 8 (expected values) against our estimated values (observed values) to assess agreement between the two datasets. This comparison showed strong agreement (r = 0.90), especially at moderate levels of change (Fig. S3). However, our projected changes were less than those predicted by ref. 8 when increases or decreases in depth to groundwater were larger. Finally, we visited streams in the study area during December and January of 2013 and 2014 and surveyed the presence/absence of water at individual stream segments that were not influenced by upstream water releases (i.e., water impoundment releases, municipal effluents). These data were collected outside of the growing season and therefore represented surface water occurrence when depth to groundwater was not drawn down by pumping activity or evapotranspiration. Data from 2013 and 2014 field surveys were compared with predicted occurrence of water in our spatial analysis for each year. We correctly predicted the presence or absence of water for 87% of field-based observations during 2013 and 2014 (Fig. S4). Misclassifications occurred along the mainstem Republican River and toward the eastern portion of the study area.

Results

Rapid and expansive increases in groundwater pumping and dam construction during 1950–2010 drastically modified stream habitats across the study region. Annual groundwater pumping from the High Plains Aquifer in the Kansas portion of the study area increased almost exponentially from 0 in 1950 to a peak of 1.31 km3 in 1980, and the number of diversion dams and reservoirs throughout the study area increased from 37 to 141 during 1950–1980 (Fig. 2A). During this same period, the total length of stream coupled with the High Plains Aquifer decreased from 2,640 to 2,082 km (558 km lost) (Fig. 2B). During 1980–2010, groundwater extraction in Kansas varied annually from 0.55 to 1.30 with an average of 0.91 km3, and the number of diversion dams and reservoirs increased to 154. Meanwhile, the total length of stream coupled with the aquifer varied from 2,046 to 2,194, with an average of 2,117 km. Projections for the period 2011–2060 indicated a further decrease in length of coupled stream to 1,796 km (another 286 km lost). Most of this loss was driven by declines in length of large (primarily fourth-order) streams although small (second- and third-order) stream length declined rapidly during 1950–1980 and then remained relatively constant or increased slightly during 1980–2060 (Fig. 2B). Although the length of fourth-order coupled streams was 1.4 times greater than second-order coupled streams in 1950, by 2060, lengths of both stream orders are projected to be similar (second, 495 km; fourth, 524 km). The coupled stream length of fifth-order streams was reduced by 58 km during 1950–2010 and is projected to decline by a total of 103 km by 2060.

Fig. 2.

  • Download figure
  • Open in new tab
  • Download powerpoint

Fig. 2.

Regional-scale change through time in (A) annual volume of groundwater pumped from the High Plains Aquifer in the Kansas portion of the study area (blue area) and the cumulative number of barriers constructed (dark-red line) during 1950–2010, (B) length of coupled stream for all sizes (blue area) and for small (second- and third-order, blue lines) and large (fourth- and fifth-order, dark-red lines) streams during 1950–2060, and (C) mean (95% confidence interval) capture probabilities for fish species associated with small (blue line and band) or large (red line and band) streams during 1950–2060.

At the regional scale, the fish assemblage was transformed from dominance by large-stream fishes to small-stream fishes during 1950–1980. The assemblage remained relatively stable during 1980–2010, but projections for 2011–2060 indicate additional suppression of large-stream fishes and expansion of small-stream fishes (Fig. 2C). The average capture probability for large-stream fishes decreased by half, from 0.42 [95% confidence interval (CI) = 0.28–0.56] in 1950 to 0.20 (95% CI = 0.11–0.29) in 1980 and is predicted to decrease further to 0.16 (95% CI = 0.09–0.24) by 2060. The average capture probability for small-stream fishes nearly doubled, from 0.20 (95% CI = 0.18–0.29) in 1950 to 0.39 (95% CI = 0.26–0.51) in 1980 and is predicted to reach 0.49 (95% CI = 0.33–0.65) by 2060.

At the subwatershed scale, minimum stream flows, the number of barriers to fish movement, and total length of coupled stream varied among subwatersheds. From 1950 to 2010, the 90-d minimum flow declined in all but one subwatershed (Table S6 and Fig. 3). The number of barriers remained at zero during 1950–2010 for two subwatersheds (Fig. 3 B and C), increased from 3 to 4 in one (Fig. 3D), increased from 3 to 8 in one (Fig. 3A), and increased from <3 to 27 in two subwatersheds (Fig. 3 E and F). During 1950–2010, total length of coupled stream declined in four subwatersheds (Fig. 3 A, C, D, and E) but remained relatively constant in the remaining two (Fig. 3 B and F). Projections for length of coupled stream during 2011–2060 included continued reduction in length for subwatersheds with historically declining lengths but increased length of coupled streams for subwatersheds with historically stable lengths. Declines in coupled stream length for the four subwatersheds were primarily caused by reduced length of fourth-order streams. Increases in coupled stream length for two subwatersheds were owing to increases in length of second- and third-order streams that offset declines in large streams of fourth- (Fig. 3B) and sixth-order (Fig. 3F).

Fig. 3.

  • Download figure
  • Open in new tab
  • Download powerpoint

Fig. 3.

Subwatershed-scale change through time in (Upper in each panel) 90-d minimum stream flow (gray circles) summarized by generalized additive models (purple line and band, fit and 95% confidence interval) and cumulative number of barriers on surface stream channels (red line) during 1950–2010 and in (Lower in each panel) length of coupled streams (blue area), and length for small (second- and third-order, dark-blue lines) and large (fourth-, fifth-, sixth-order, dark-red lines) coupled streams during 1950–2060. Letters represent subwatersheds as follows: (A) Frenchman Creek, (B) North Fork Republican River, (C) Arikaree River, (D) South Fork Republican River, (E) Upper Republican River, and (F) Harlan County Reservoir (Fig. 1 for locations).

View this table:

  • View inline

Table S6.

Generalized additive model results for the relationship between time (year) and annual 90-d minimum flow statistics calculated using US Geological Survey (USGS) flow data for six subwatersheds

Long-term changes in fish assemblage composition at the subwatershed scale varied across the study area in the Kansas River basin (Fig. 4). During 1950–2010, capture probabilities for small-stream fish exceeded those for large-stream fish in two western subwatersheds (Fig. 4 A and B); they increased from a mean of <0.1 to >0.40 for small-stream fish whereas they declined from a mean of >0.3 to <0.2 for large-stream fishes in two other western subwatersheds (Fig. 4 C and D), and capture probabilities for large-stream fishes exceeded those for small-stream fishes in two eastern subwatersheds (Fig. 4 E and F). Projections for 2011–2060 included relatively stable dominance by small-stream fishes where this group was historically dominant (Fig. 4 A and B), included continued increase of small-stream fishes combined with continued suppression of large-stream fishes where historical shifts in dominance occurred (Fig. 4 C and D), included steady increase of small stream fishes lower in the basin (Fig. 4E), and included continued dominance by large-stream fishes until at least 2050 in the farthest downstream subwatershed (Fig. 4F).

Fig. 4.

  • Download figure
  • Open in new tab
  • Download powerpoint

Fig. 4.

Change in estimated depth to groundwater between 1950 and 2060 for the study region. Change through time in mean (95% confidence interval) capture probability for fishes characteristic of small (<fourth-order, blue line and band) and large (≥fourth-order, red line and band) streams are shown for six subwatersheds as follows: (A) Frenchman Creek, (B) North Fork Republican River, (C) Arikaree River, (D) South Fork Republican River, (E) Upper Republican River, and (F) Harlan County Reservoir.

Discussion

Our study provides evidence that groundwater pumping during the past half-century has caused declines in the length of stream coupled with the High Plains Aquifer and greatly altered stream fish assemblages of the western Great Plains. Rapid increases in pumping in Kansas during 1950–1980, followed by sustained pumping thereafter, increased depth to groundwater in broad expanses of the study area and resulted in declines in the length of stream coupled with the aquifer in many subwatersheds. Declines in stream length were coincident with reduced flows that further diminished the prevalence of larger streams in the region. Similar declines in groundwater, streamflow, and stream length caused by groundwater pumping have been reported elsewhere in the basin (34), across western Kansas (4), throughout the High Plains Aquifer (2), and globally (9).

The disproportionate decline in the prevalence of larger streams, in turn, contributed to disproportional declines in large-stream fishes. However, decline of large streams and large-stream fishes at the regional scale was confounded by a concurrent increase in the number of barriers, which are demonstrated to have negative effects on some large-stream fishes (24, 28). Our approach of scaling down to subwatersheds disentangled the relative effects of barriers and groundwater declines. Specifically, barrier numbers were greatest, and increased rapidly, in eastern subwatersheds where large-stream fish capture probabilities were highest among all subwatersheds through time. Conversely, barriers were rare or absent in western subwatersheds where small-stream fish capture probabilities were greatest or increased most rapidly, and where large-stream fish capture probabilities were least or decreased most rapidly. Large-stream fish capture probabilities were initially low or declined rapidly among upstream, western subwatersheds where extinction and colonization are naturally high for such fishes (35), and long-term declines attenuated in a downstream direction to the east (Fig. 4 C–F). These patterns collectively illustrate accelerated fish assemblage transformation along a natural stream size gradient initiated by the decoupling of surface streams from groundwater sources and ultimately resulting in the permanent loss of large-stream fishes.

Spatial patterns in historical (pre-1950) and more recent (post-1950) hydrologic alterations contributed to the observed patterns in fish assemblage change. The greatest area of increased depth to groundwater occurred in the North Fork Republican River subwatershed, and larger streams of fourth order generally dried over this area. Streams that remained were small second- and third-order channels near the westward extent of the subwatershed in areas where depth to groundwater is projected to approach the surface and create more small-stream habitat. This pattern explains the dominance by small-stream fishes in this subwatershed. In contrast, the rates of decline in coupled stream lengths and minimum flows were similar among Frenchman Creek, the Arikaree River, and the South Fork Republican River subwatersheds although Frenchman Creek did not experience the same fish assemblage transformation as the Arikaree River and South Fork Republican River subwatersheds. However, Frenchman Creek had a large number of barriers that fragmented habitat before 1950 so it is possible that the fish assemblage was altered earlier than 1950. For example, among the large-stream fishes excluded from analysis because of rare occurrences (Table S5), four were last reported from the region by 1940 (Hiodon alosoides, Macrhybopsis gelida, Macrhybopsis hyostoma, and Macrhybopsis storeriana), and declines among three of these (all Macrhybopsis) have been linked to fragmentation (26). Fragmentation of Great Plains riverscapes began before 1950 (23), and recent work suggests that widespread groundwater extraction may have begun about 1930 instead of 1950 in some portions of the Great Plains (5). More historical information might be necessary for a complete understanding of the ecological consequences of groundwater extraction. Furthermore, declines in minimum flows cannot be attributed to groundwater extraction alone because surface diversions and extractions also reduce water availability in the Great Plains (10, 36), suggesting that the effects of surface barriers on assemblage transformation cannot be ignored.

The mechanisms driving differential responses to groundwater pumping by large- and small-stream fishes are reinforced by permanent water loss and barriers that fragment habitats. Larger streams contain more predictable habitats characterized by greater flows, deeper channels, longer longitudinal connectivity, and higher autochthonous energy production (37). Examples of adaptations of large-stream fishes to this environmental template include migratory behavior, flow-induced synchronized spawning, and spawning within the water column (38, 39). As groundwater pumping shortens the length of stream coupled with the aquifer and surface structures store or divert water, flows decline, channels become shallower, longitudinal connectivity is fragmented, and instream (autochthonous) energy production decreases. Collectively, these processes shift large streams into habitat templates that no longer match the evolutionary history of large-stream fishes. Rather, habitats become characteristic of those that shaped the ecology and life history of small-stream fishes, thus providing increased opportunities for expansion (40). In our study area, eastern subwatersheds with relatively stable depths to groundwater maintained large-stream fish assemblages, but western subwatersheds with increasing depths to groundwater experienced permanent reduction in large-stream fish assemblages. Unfortunately, western assemblages of large-stream fishes that decline or collapse have no opportunity for rebounding (e.g., recolonization, rescue effects) because flows are either permanently desiccated or habitats are isolated from eastern subwatersheds behind dozens of instream barriers. Based on these and related mechanisms of assemblage change (e.g., drought, habitat destruction) (28), the shrinking of Great plains fish assemblages characterized by replacement of large-stream fishes with small-stream fishes is reinforced by the existence of surface barriers (27, 41).

Our projections of future ecological outcomes caused by continued groundwater pumping from the High Plains Aquifer rest on several critical assumptions. First, we assumed that pumping rates were constant and rates of change in depth to groundwater were linear at individual observation wells during 1980–2010. Although this assumption has considerable precedence (4, 36), rates of change at some wells might be more dynamic or proceed faster than estimated. For example, we validated our projected changes in depth to groundwater using published US Geological Survey data and found evidence that the magnitude of changes included in our study were relatively less compared with changes documented by the US Geological Survey (Figs. S2 and S3). Second, a simplifying assumption in our model was that the aquifer was coupled with surface stream dynamics (i.e., streams received water from the aquifer) if interpolated depth to groundwater was within 1 m of the surface. In reality, coupling between aquifers and surface ecosystems occurs at depths of 0 m when water flows from groundwater sources into streams (16–18), but our use of a 1-m threshold allowed for measurement error for both depth to groundwater and the digital elevation model (SI Methods). A related limitation was the lack of a calibrated groundwater flow model covering the period 1950–2060. We used a method consistent with existing applications for measuring groundwater change (SI Methods), and we expect that future analyses at finer scales will benefit from improved hydrologic models. Finally, our fish assemblage data originated from multiple sources that sampled fish using a variety of gear types with species-specific capture biases. We addressed this issue by using occurrence data originating only from collections targeting entire assemblages. However, change in collection gears over time could still affect capture probabilities. We believe the potential effect of changes in sampling gear to be minor compared with environmental forcing because (i) subwatershed-scale analyses did not illustrate uniform changes in capture probabilities (instead temporal change in observed captures varied longitudinally), (ii) large-stream fish capture probabilities generally declined through time during 1950–1980 despite potential for increased capture efficiency with modern gears (e.g., electrofishing), and (iii) the effort required to capture 90 to 100% of species present in small Great Plains stream ecosystems is relatively consistent among commonly applied gear types (42).

Fig. S3.

  • Download figure
  • Open in new tab
  • Download powerpoint

Fig. S3.

Regression plot comparing estimated change in depth to groundwater (m) between 1950 and 2011 at individual stream segments distributed across the study area as estimated by ref. 8 (x axis) and this study (y axis).

Groundwater depletion caused by pumping is a global environmental problem (9, 42), but analyses of ecological consequences for fish assemblages are rare (6, 13). We documented empirical evidence for the relationship between groundwater depletion in the High Plains Aquifer and the decline in Great Plains stream fish assemblages through extirpation of large-stream fishes and expansion by small-stream fishes. Our analysis provides baseline predictions for stream–aquifer coupling and the ecological status of Great Plains fish assemblages over the next 45 y. The predictions for streams could be used as a benchmark against which long-term conservation goals might be set (e.g., slowing rate of coupled stream loss), and the status of fishes provides a metric of environmental change that is embedded within natural ecosystem functioning (43). Evolving scientific understanding will improve projections and expectations for future change, and adapting short-term groundwater policies and practices to accommodate new information will be at the core of adaptive management (44). This point is particularly critical in the context of climate change, given expected regional precipitation regime shifts, increased evapotranspiration, and accelerated stream loss (45, 46). Groundwater pumping for agriculture is a major cause of alteration to global freshwater ecosystems (47), and increasing the efficiency of human water use is imperative for the future well-being of nature and humans (48, 49). Thus, our study has global implications, given consistencies in groundwater depletion from the High Plains Aquifer and other local and regional aquifers worldwide (9, 41, 50).

Methods

Study Area.

We studied the upper Kansas River Basin in Colorado, Kansas, and Nebraska, United States, where streams are underlain by the High Plains Aquifer (Fig. 1). Here, depths to groundwater have increased up to 30 m since 1950 when extensive groundwater pumping began (8). Natural groundwater recharge rates and aquifer confinement vary across this region (2) so pumping varies from sustainable (Nebraska) to unsustainable (Kansas). Nevertheless, even in Nebraska, pumping and diversions have reduced stream base flows by up to 50% (51). Groundwater declines have decoupled more streams from the aquifer to the west whereas streams to the east are more fragmented by diversion dams and impoundments (23). In Kansas and Colorado, groundwater levels in portions of the aquifer are projected to reach the bottom of the aquifer by 2060 (4), and many stream reaches have already dried (13). Given this variation in groundwater levels and stream drying, we selected six subwatersheds (Fig. 1 and Table S6) arrayed from east to west across the region to analyze hydrologic and fish assemblage responses to groundwater pumping and fragmentation by barriers on the surface.

Spatial Analysis.

We used a combination of retrospective and prospective approaches to estimate the length of stream coupled with the High Plains Aquifer across the study region during 1950–2060. We first estimated the timing and magnitude of groundwater pumping in the region using data from the Kansas Geological Survey (52). We then used a network of observation wells monitored between 1950 and 2010 and distributed across the study region (Fig. 1) to measure spatiotemporal changes in depth to groundwater (SI Methods). Based on the rates of change in depth to groundwater in individual wells during 1980–2010, we projected future depths to groundwater for the period 2011–2060 using linear regression. We then created interpolated surfaces (cell size 0.5 km) of water table elevations by subtracting depth to groundwater from the elevation of the ground surface at observation wells across the region for every year between 1950 and 2060. We then subtracted interpolated water table elevations from the ground surface elevation at the locations of streams (based on ref. 32) and used groundwater depths to estimate the length of stream coupled with the aquifer each year for the 110-y period. Finally, we validated depth to groundwater estimates using existing published data and field reconnaissance (SI Methods and Fig. S4).

Fig. S4.

  • Download figure
  • Open in new tab
  • Download powerpoint

Fig. S4.

Regional map illustrating the spatial locations of stream segments coupled with (blue lines) or decoupled from (gray lines) the High Plains Aquifer for the year 2013 based on depth to groundwater estimates generated during this study. Segments were visited during December and January of 2013 and 2014 and designated as correctly classified (gray symbols) or misclassified (black symbols) according to the presence of surface water.

Surface Stream Conditions.

We analyzed changes in surface stream channel connectivity and flow in six subwatersheds distributed across the region (Fig. 1). Diversion dams and impoundments that fragment stream channels were located using the 2012 National Anthropogenic Barrier Dataset (53), which included dates of construction for most barriers. We illustrate changes in surface channel fragmentation in each subwatershed by plotting the cumulative number of dams through time. We analyzed changes in stream discharge at the downstream extent of each subwatershed using data from US Geological Survey gauges (Fig. 1 and Table S6). Streamflow data for 1950–2010 were analyzed using the program Indicators of Hydrologic Alteration (54) to calculate the mean discharge for the lowest 90 consecutive days each year (hereafter, 90-d minimum flow) to assess change in low flows through time. We fit generalized additive models with a Poisson error distribution and overdispersion parameter to summarize change in minimum flows through time for each subwatershed using the "mgcv" Package in the program R (55).

Fish Assemblage Analysis.

We assessed change in fish assemblage structure at the regional and subwatershed scales using existing historical data from Colorado, Kansas, and Nebraska (SI Methods). We selected species for analysis based on occurrences across 940 collections and retained only those species with at least 30 occurrences for regression analysis (Tables S1–S5). At the regional and subwatershed scales, we used annual estimates for the length of stream connected to the aquifer as the predictor variable and occurrence of each retained species in collections taken during that year as the response variable for the period 1950–2010. We fit binomial logistic regression models to each species and used these models to predict capture probability (range 0–1) for all years between 1950 and 2060 based on estimated annual lengths of stream coupled with the aquifer. We calculated the mean and 95% confidence interval across all species classified as small-stream or large-stream inhabitants (SI Methods) to illustrate fish assemblage response to spatiotemporal variability in stream lengths coupled with groundwater. Although fish distributions in relation to stream size formed a continuum (Fig. S1), we used the designation of first- to third-order streams as "headwaters" (37) to facilitate comparison of how fishes in small and large streams might be affected by groundwater depletion. Fish collection methods were consistent with protocols approved by the Tennessee Technological University Institutional Animal Care and Use Committee (permit TTU-IACUC-14-15-001 to J.S.P.).

Acknowledgments

R. Waters [Kansas Department of Wildlife, Parks, and Tourism (KDWPT)], M. VanScoyoc (KDWPT), and S. Schainost (Nebraska Game and Parks Commission) provided fish assemblage data. Discussions with D. Chandler, D. Steward, and Y. Yang helped develop methods for linking groundwater elevations to surface stream length. We thank R. Arnold, J. Broska, P. Burch, D. Durnford, and two anonymous reviewers for help with hydrologic modeling and critical reviews of the work. Funding was provided by the Great Plains Landscape Conservation Cooperative and managed by the Wildlife Management Institute. Any use of trade, firm, or product names is for descriptive purposes only and does not imply endorsement by the US Government.

Footnotes

  • Author contributions: J.S.P., K.B.G., J.A.F., K.D.F., H.C., E.R.J., and J.S. designed research; J.S.P., K.B.G., and H.C. performed research; J.S.P. analyzed data; and J.S.P., K.B.G., J.A.F., K.D.F., H.C., E.R.J., and J.S. wrote the paper.

  • The authors declare no conflict of interest.

  • This article is a PNAS Direct Submission.

  • This article contains supporting information online at www.pnas.org/lookup/suppl/doi:10.1073/pnas.1618936114/-/DCSupplemental.