Conservation partners are concerned that oil and gas development in the Prairie Pothole Region may reduce the abundance of breeding duck pairs using associated wetland habitat. We conducted wetland-based surveys for breeding pairs of 5 species of dabbling ducks in the Bakken oil field during 2015–2017 across a gradient of oil and gas development intensity to test the hypothesis that the abundance of breeding duck pairs on survey wetlands would decrease as the development of oil and gas resources increased. We included covariates traditionally used to predict breeding duck pairs (i.e. wetland size and class) and a spatiotemporal index of disturbance when developing zero-inflated Poisson models relating pair abundance to environmental predictors. Similar to past analyses, pair abundance was strongly associated with wetland size. Our results were mixed and suggested that the abundance of early and late nesting species was positively and negatively related, respectively, to an index of disturbance that was largely driven by oil and gas development. Regardless of the direction of the relationship, effect sizes were small and not considered biologically significant. Our findings indicate that in our study area, strategies to conserve wetland resources for breeding duck pairs should not deviate from previous prioritization metrics within the range of oil and gas development we observed. We believe that our findings may have implications to similar landscapes within the Bakken.
LAY SUMMARY
Since 2008, the oil and gas development in the North Dakota and Montana portion of Bakken Oil Formation has increased dramatically.
There is considerable overlap between the Bakken Oil Formation and important Prairie Pothole Region wetlands critical for waterfowl production.
We surveyed breeding Blue-winged Teal, Gadwall, Mallard, Northern Pintail, and Northern Shoveler pairs from 2015 to 2017 to determine if breeding pair abundance was lower in proximity to a gradient of disturbance from oil and gas development.
Our results were mixed but regardless, changes in pair abundance were small and we considered the potential biological effect to be small.
We recommend that existing conservation tools continue to be used to identify important grassland and wetland resources in the region given that we did not observe a biologically significant reduction in breeding duck pairs.
INTRODUCTION
The distribution, density, and upland accessibility by breeding duck pairs are used to prioritize different conservation actions in the US portion of the Prairie Pothole Region. Annual water conditions of small, shallow, wetlands in this landscape affect the resources available to breeding pairs, and the subsequent waterfowl carrying capacity of the Prairie Pothole Region (Batt et al. 1989). Because of the region's high value for breeding waterfowl, it is considered a high priority landscape for resource conservation (NAWMP 2012, PPJV 2017), and the U.S. Fish and Wildlife Service has acquired a considerable conservation portfolio (see USFWS 2004–2019). Since 2004, the U.S. Fish and Wildlife Service has invested >$54 million (US) (USFWS 2004–2019) on wetland and grassland habitat protection within the extent of the Bakken Formation in North Dakota. While the conservation investment has prioritized waterfowl habitat, other grassland- and wetland-dependent wildlife species certainly benefit from those conservation actions (PPJV 2017). However, a continued waterfowl conservation program delivery within a landscape that fails to consider significant changes or disturbances risks incorrect prioritization of habitat for breeding duck pairs if the disturbances alter or diminish habitat use (e.g., duck pair occupancy or abundance). The primary threat to wetlands important for waterfowl production in the Prairie Pothole Region is the drainage for crop production agriculture (Dahl 2014). However, displacement of breeding pairs from sources of energy development, such as wind, has also been observed (Loesch et al. 2013). A recent surge in oil and natural gas development has resulted in concern about similar effects where disturbance may reduce the value of wetlands within developed landscapes for breeding waterfowl.
Beginning in 2008, advances in horizontal drilling and hydraulic fracturing (hereafter fracking) resulted in a substantial increase in oil and gas extraction-related development in northwest North Dakota and northeast Montana (Gaswirth et al. 2013). In 2012, the number of active drilling rigs peaked in North Dakota at 217. As of 2018, there were 8,725 active wells in the US Prairie Pothole Region of North Dakota, and 78% of those wells had been drilled since 2008.
The effect of oil and gas extraction on wildlife and associated habitat is diverse (Northrup and Wittemyer 2013). As with other forms of energy resource development (e.g., wind and coal-bed methane), extraction sites, infrastructure development, and increased heavy equipment and associated vehicle traffic can affect wildlife through mortality and habitat loss (Kuvlesky et al. 2007, Francis et al. 2011, Naugle et al. 2011, Brittingham et al. 2014, Thompson et al. 2015). The loss of grasslands in the northern Great Plains has resulted primarily from conversion to cropland agriculture. The remaining grasslands are important for breeding grassland birds, several of which are experiencing long-term declines (Sauer et al. 2017, Rosenberg et al. 2019). Since 2000, an additional 6,155 ha of grassland conversion in the Williston Basin (Montana, North Dakota, South Dakota, USA; Manitoba, Saskatchewan, CAN), which encompasses the Bakken Formation (hereafter, the Bakken), can be attributed to roads or well pad construction, and an additional 2.7% of the remaining 20,920 ha of grasslands are predicted to be affected in the future (Preston and Kim 2016). Specific to the portion of the Bakken north of the Missouri River in North Dakota, Dyke et al. (2010) estimated that ∼3,100 ha of grassland would be lost by the year 2020 due to direct impacts by oil and gas development. This grassland loss will reduce available areas for nesting and may also negatively affect duck nest survival (Greenwood et al. 1995, Reynolds et al. 2001, Stephens et al. 2005), the density of breeding duck pairs on associated wetlands (Reynolds et al. 2007), and brood abundance (Carrlson et al. 2018).
Changes in animal abundance during all or portions of their annual cycle can occur due to factors other than habitat loss. Disturbance accompanying energy development (e.g., increased vehicle traffic, operating noise, and infrastructure construction activity) may reduce both the wildlife value and use of nearby habitat (Habib et al. 2007, Bayne et al. 2008, Gilbert and Chalfoun 2011, Blickley et al. 2012, Loesch et al. 2013, Thompson et al. 2015, Shaffer and Buhl 2016). While this has been commonly noted in sage-grouse (Blickley et al. 2012; see Patricelli et al. 2013) and passerines (Habib et al. 2007, Francis et al. 2011, Thompson et al. 2015), few studies isolated the mechanisms driving these behaviors. Visual disturbance requires line of sight from the occupied habitat and may change the spatial distribution of breeding duck pairs or reduce their abundance and vital rates as a result of reduced foraging opportunities relating to increased vigilance (Graeme et al. 2016). Anthropogenic noise is difficult to measure, and response varies among species (Bowles 1995). Small passerines are often thought to avoid noisy infrastructure because it obscures their breeding songs (Thompson et al. 2015). In addition to potential changes in abundance, waterfowl mortality occurs in skim pits (Esmoil and Anderson 1995) and centralized oil field wastewater disposal facilities that use open evaporation ponds (Ramirez 2010). This impact may be negligible in the Bakken because most of the maintenance and produced water in the Bakken is injected into wells (Kurz et al. 2016). Soil disturbance events (e.g., well pad construction, spills, and pipeline installation) can result in the establishment of invasive plants (Trombulak and Frissell 2000, Evangelista et al. 2011) that may lower habitat quality for some species (Nelson et al. 2017). Produced water spills during drilling (e.g., brine) result in sodium toxicity, which kills terrestrial and aquatic vegetation (Halvorson and Lang 1989), negatively affects water quality, and is difficult to remediate (USEPA 2011, Gleason and Tangen 2014, Latta et al. 2015, Lauer et al. 2016). Spills that enter wetlands also negatively affect macroinvertebrates (Preston and Ray 2017), which are important protein sources for breeding female ducks and developing ducklings (Krapu and Reinecke 1992, Cox et al. 1998).
Disturbance from oil and gas development could negatively impact a relatively large geographic area (US Prairie Pothole Region and Bakken overlap = 87,410 km2) that provides habitat for some of the highest numbers of breeding dabbling duck pairs in the Prairie Pothole Region of the United States or anywhere else in the North America (Reynolds et al. 2006, PPJV 2017). Our objective was to measure the potential change in abundance from disturbance over a gradient of oil and gas resource extraction in the Bakken to breeding pairs of the 5 most common species of breeding dabbling ducks in the Prairie Pothole Region. These 5 species included Blue-winged Teal (Spatula discors), Gadwall (Mareca strepera), Mallard (Anas platyrhynchos), Northern Pintail (Anas acuta), and Northern Shoveler (Spatula clypeata). We hypothesized that if breeding duck pairs were negatively affected by oil and gas development from disturbance or decreased habitat quality, pair abundance would be lower on wetlands in areas with higher oil and gas development. Our study area and hypotheses paralleled other concomitant investigations of the impacts of oil and gas development on duck nest survival (Skaggs et al. 2020) and duck brood abundance (Kemink et al. 2019) in the Prairie Pothole Region of northwest North Dakota and northeast Montana.
Study Area
The Prairie Pothole Region was historically an area of rolling topography and a mosaic of grassland and abundant depressional wetlands (Bluemle 2000). While many wetlands remain, 70% of the previously extensive grasslands have been converted to other uses, predominantly agriculture (Doherty et al. 2013, Dahl 2014, Gleason and Tangen 2014). The climate of the study area was continental in nature and characterized by a short, hot growing season (mean temperature: 18.1°C [June to August]) and long, cold winters (mean temperature: –7.3°C [November to March]; NOWData 2020). Precipitation varied temporally and averaged 36.4 cm (range: 16.0–55.4, Williston, ND, 1900–2019) and ranged from 30.5 to 39.4 cm during our study (Williston, ND, 2015–2017; NOWData 2020).
Our study area and much of the Bakken overlapped the Missouri Coteau physiographic region (Figure 1). The Missouri Coteau was the terminal moraine of the Wisconsin Glaciation and contained abundant hummocky collapsed basins that, when ponded, provided some of the most important breeding waterfowl habitats in the US portion of the Prairie Pothole Region (Bluemle 2000, Reynolds et al. 2006, Doherty et al. 2015, PPJV 2017). High densities of pothole wetlands attracted breeding and migrating waterfowl, and highly variable local wetland conditions annually influenced the distribution and abundance of breeding duck pairs (Johnson and Grier 1988, Reynolds et al. 2006, Niemuth and Solberg 2003, Niemuth et al. 2010). On average (years 1987–2016), the habitat in the Bakken supported nearly one-third of the breeding duck population in the North Dakota, South Dakota, and northeastern Montana portions of the Prairie Pothole Region using models from Reynolds et al. (2006) and Niemuth et al. (2010).
METHODS
Wetland Sample Selection
Our potential sample wetlands were within a 25,391 km2 area where the Bakken overlapped the Prairie Pothole Region west of the Souris River in North Dakota and north of the Missouri River in North Dakota and eastern Montana. We defined a convex polygon around geospatial well locations obtained from the North Dakota Industrial Commission (2018) and Montana Board of Oil and Gas Commission (2019). We used well locations that were active beginning January 01, 2004, because that date coincided with a year of considerable increase in fracking activity in the Bakken and included >87% of the active oil and gas wells within our project area in 2015.
FIGURE 1.
Location of 81 plots (3.2 × 3.2 km) within the Bakken Formation where wetlands were surveyed for breeding Blue-winged Teal, Gadwall, Mallard, Northern Pintail, and Northern Shoveler pairs during spring, 2015–2017. Duck pair density information originates from models updated from those presented in Reynolds et al. (2006).
We used a GIS to subdivide the sample area into a 3.2 × 3.2 km grid where each cell was considered a potential survey plot. We used this size because it approximated the home range of a breeding pair of mallards (Cowardin et al. 1988) and was also used as the foundation of a long-term breeding waterfowl population and production survey conducted annually by the U.S. Fish and Wildlife Service since 1987 Cowardin et al. 1995). We then removed plots where the U.S. Fish and Wildlife Service conducted its annual population and production survey (N. Wright, personal communication) so that we did not interfere with their survey.
We developed a geospatial perennial cover layer from the US Department of Agriculture National Agriculture Statistics Service Cropland Data Layer (Boryan et al. 2011). We reclassified the grassland/pasture, winter wheat, and shrubland cover classes identified in the Cropland Data Layer ( https://data.nal.usda.gov/dataset/cropscape-cropland-data-layer/resource/035d2848-dacd-48d2-be92-b5379fc8eb4d accessed March 2015, 2016, and 2017) to perennial cover for the calendar year preceding the survey year, and we used a focalsum function to calculate the amount of perennial cover in the 10.4 km2 around the center point of each 3.2 km × 3.2 km grid cell that represented a potential sample plot. We removed plots with <40% cover and <100 wetlands in the 3.2 × 3.2 km plot from consideration as sample plots in order to minimize potential variation in pair abundance relative to upland composition and wetland density from our sample universe for the respective year. Each year, we randomly selected 64–71 of the remaining plots on which to conduct breeding pair surveys (Figure 1). The number varied depending on the available workforce, wetland conditions, and the collection of aerial photography.
We stratified our wetland sample across all sample plots each year by the changes in the level of wetland disturbance within 3.64 km of each wetland's centroid. For the purposes of wetland sampling, we defined wetland disturbance as the number of active oil well pads within 3.64 km and selected 1,000 wetlands in each of the 4 disturbance strata (control: 0 well pads, low: 1 well pad, medium: 2–3 well pads, and high: >3 well pads). These disturbance strata were only used to distribute sampling effort across the gradient of well density and were not used in subsequent analysis. Our wetland sample focused only on temporary, seasonal, and semipermanent wetlands <5 ha (see Johnson and Higgins 1997), which make up 97% and 43% of the wetland basins and hectares, respectively, in the Bakken. We focused on these wetland classes because of their importance for wetland pairs (Reynolds et al. 2006, Loesch et al. 2012) and limited wetland size because of the logistics of surveying large wetlands. We only surveyed wetlands if the extent was completely within the plot.
Breeding Pair Habitat Assessment
Small prairie pothole wetlands that are critical to breeding waterfowl are often dry. Consequently, we used the number of wet basins and the area of ponded surface water within a plot to assess annual wetland conditions. We used high-resolution (1.5 m) photography collected with a fixed-wing aircraft to assess the spring wetland habitat conditions for all sample plots in each of the survey years (Cowardin et al. 1995). We collected images of the plots from May 3 to May 23. We post processed the images to segment related pixels and then manually classified ponded surface water on each plot relative to wetlands that were mapped by the U.S. Fish and Wildlife Service National Wetlands Inventory circa 1980.
Breeding Pair Surveys
We conducted a ground visit to each sampled wetland to survey duck social groups (i.e. one or more interacting conspecifics) once during an early survey period (i.e. May 1–17) and once again during a late survey period (i.e. May 18 to June 5) to account for differences in the breeding chronology of different species of ducks in the region (Stewart and Kantrud 1973, Higgins et al. 1992, Cowardin et al. 1995, Naugle et al. 2000). A key assumption of the survey was that breeding territories had been established, and pairs were identifiable relative to the respective survey periods (e.g., early: Mallard, Northern Pintail; late: Northern Blue-winged Teal, Gadwall, Northern Shoveler; Dzubin 1969, Cowardin et al. 1995). We recorded social groups of interacting conspecifics, but separated by space from other conspecifics, as ratios of males to females. The sex ratios from social groups identified during the relevant survey times (i.e. early and late) were interpreted to 7 categories that were subsequently used to identify indicated and observed breeding pairs for each species (Dzubin 1969, Hammond 1969). Indicated and observed breeding pairs were combined and used as the dependent variable in the analysis.
We approached each survey wetland on foot and attempted to reduce disturbance to avoid artificially influencing survey counts. We used binoculars to scan the sample wetland for the presence of the 5 target species. Ducks arriving on unscanned portions of the wetland were included in the survey, and conversely, we did not include ducks arriving on previously scanned portions of the wetland. Additionally, ducks that were disturbed by the observers and left the wetland were recorded, but we did not include ducks that left the wetland of their own volition before they were surveyed. We suspended surveys when visual cues suggested winds exceeded 48 km hr–1 and pairs were likely to seek cover (e.g., waves on wetlands, bend of trees from wind) or during steady rainfall when visibility was limited. Surveys resumed when weather conditions improved.
The length of time for each wetland survey varied and was influenced by the size and configuration of the wetland, vegetation composition, and the presence and abundance of waterfowl. When the wetland survey was completed, we recorded a global positioning system location at the edge of the ponded water for each wetland to verify that the correct wetland was surveyed. Wetlands containing ponded water during the early visit were surveyed again during the late survey period. If a sample wetland was dry during either survey period, a previously identified replacement wetland within the same sample plot was surveyed. Wetlands that were dry during the first survey period were not visited during the second survey period of the same year.
In addition to waterfowl count data, we recorded several environmental variables during and after each survey. These variables included survey time, wetland UTM zone 14 northing and easting, ocular estimates of the percent of wetland with surface water, wetland class (i.e. temporary, seasonal, and semipermanent) where the class is a function of the duration of water permanence assigned by the National Wetlands Inventory (i.e. regime modifiers A, C, and F; Cowardin et al. 1979, Johnson and Higgins 1997), and dominant emergent vegetation height. We estimated the emergent vegetation height (i.e. <25 or >25 cm) and classified the water and emergent cover distribution according to guidelines described by Stewart and Kantrud (1971). Ocular estimate of the percent of surface water relative to the wetland area mapped by National Wetlands Inventory was only used to inform the water mapping process for the images collected during the spring. If no water was detected during image processing, the field visit was used to confirm the absence of water in sample wetlands.
Well Pad and Extraction Activity and Disturbance Index
We used the presence and density of physical features associated with oil and gas extraction infrastructure that included well pads, transportation infrastructure, and surface gravel mines as a proxy for the combination of auditory and visual disturbance to breeding waterfowl pairs. Potential disturbance related to these features included vehicle traffic and noise, active drilling, well-site equipment and operation, well site construction and maintenance, road construction, and human activity. A single extraction pump well pad typically consisted of a 0.2-ha gravel pad enclosed by a berm to contain spills, an extraction pump, an oil-water separator unit, and oil and produced-water storage tanks (Kemink et al. 2019).
Each year, we conducted ground checks of each survey plot and associated 3.2 km buffer area to confirm the presence of oil and gas-related infrastructure in the available geospatial data (MBOGC 2019, North Dakota Industrial Commission 2018). We also mapped additional well pad locations, the number of extraction pumps present, gravel pads containing obvious oil and gas-related features (e.g., storage tanks), active fracking sites, and surface gravel mines that were not present in existing geospatial data.
Transportation infrastructure was developed by updating available linear, geospatial road data for the study area ( https://datagateway.nrcs.usda.gov, accessed December 2, 2017). We used geographic information system techniques to update the transportation layer to include gravel and paved surface roads by overlaying the existing road data on current aerial imagery collected by National Imagery Program (USDA 2016). Because the geospatial data were linear and our disturbance index required a spatial footprint, we randomly chose 60 gravel and 60 paved road segments to determine the average width for each surface type (Kemink et al. 2019). We observed high variation and similarity between the road surface types (gravel x̄ = 7.8 ± 1.9 m, paved x̄ = 9.9 ± 2.7 m). We buffered the roads using the widest average size in our sample (10 m).
Analysis
Breeding pair abundance. We used zero-inflated Poisson (ZIP) models to develop pair abundance models for each of the 5 most common dabbling duck species in our study area. A generalized linear model (GLM) with Poisson errors typically provides a good starting point for count data (McCullagh and Nelder 1989). However, our data included high amounts of over-dispersion relative to standard Poisson assumptions (excess zeros and infrequent large counts; Zuur et al. 2007). To deal with this over-dispersion, we selected a ZIP model, which allowed us to account for the large proportion of observed zero values (Lambert 1992, Hall 2000). While formal selection processes are available to determine when ZIP models are preferable over traditionally distributed GLM models, these tests are not considered appropriate or necessary when a clear departure (e.g., >50% of the data consists of zeros) from the Poisson distribution is indicated (Arab et al. 2008).
TABLE 1.
Log-scale parameter estimates of the posterior median (PM) and 95% credible intervals (lower credible interval [LCI] and upper credible interval [UCI]) from the best-fitting abundance models of breeding Blue-winged Teal, Gadwall, Mallard, Northern Pintail, and Northern Shoveler. Models were developed using breeding pair data collected during spring, 2015–2017, wetland surveys in the Prairie Pothole Region of the Bakken Formation of North Dakota and Montana, USA. Intercept coefficient is for seasonally ponded wetlands.
We developed an index to represent disturbance caused by oil and gas development for inclusion in our ZIP models. We developed this index relative to each sample wetland based on the geospatial footprint of annually updated transportation infrastructure, oil well pads, and surface gravel mines to characterize oil and gas disturbance. We calculated the proportion of a buffer circle around each sample wetland that was covered by the disturbance features. We used biologically relevant buffer distances such that each of the 5 modeled pair species' buffer circle radii was represented by distances derived from previous studies that report travel distances from core breeding wetlands to upland nest sites (Northern Pintail = 4.03 km, Mallard = 3.62 km, Gadwall and Blue-winged Teal = 1.61 km, and Northern Shoveler = 1.21 km (Reynolds et al. 2006; Table 1). Using a fuzzy algebraic sum, a common operator in geospatial analyses, we then ascertained a final value for the disturbance index (Bonham-Carter 1994, Malczewski 1999, Theobald 2013, Kennedy et al. 2019). The advantages of the disturbance index were 2-fold. First, it addressed comparisons across equal proportions of cover, allowing for unambiguous interpretations (Riitters et al. 2009). Second, it dealt with potential issues of collinearity while still allowing for the assessment of additive effects of existing and new disturbance on the landscape.
Using past analyses of pair abundance within the Prairie Pothole Region, we selected well-established wetland- and landscape-scale covariates to include in our species-specific models of abundance, in addition to the disturbance index (Cowardin et al. 1988, 1995, Loesch et al. 2013). We included wetland class (i.e. temporarily flooded, seasonally flooded, and semipermanently flooded; Johnson and Higgins 1997), ponded water surface area, and a log-transformed wetland mapped wet area variable (Cowardin et al. 1988, 1995, Reynolds et al. 2006, Loesch et al. 2013). We also included a plot-year covariate, the number of wet wetlands during May, and the amount of perennial herbaceous vegetation as landscape covariates in our species-specific global models.
We applied Akaike's information criterion (AIC) to test our hypotheses regarding energy development and pair abundance (Burnham and Anderson 2002). The model selection process involved 2 stages. In the first stage, we used changes in AIC to acquire reduced models (Burnham and Anderson 2002). We removed a single variable from the fully parameterized ZIP model and recorded the change in the AIC value after running the model (Chambers 1992, Crawley 2007). We repeated this process for each covariate in the model and the final, reduced model consisted of the covariates whose removal resulted in increases of >2 AIC units (Arnold 2010). We used this final model in the last stage of our analysis wherein we ran mixed effect models for each of the 5 species in a Bayesian framework. We incorporated wetland and plot-year random effects in the count models to assess the influence of random noise on covariate values. We conducted all statistical analyses in the R environment (R Development Core Team 2016). We generated Bayesian estimates of model parameters using the rjags package to run Markov chain Monte Carlo iterations (Plummer 2016). For each model, we ran 2 Markov chain Monte Carlo chains for 700,000 iterations and discarded the first 200,000 iterations to minimize the influence of starting values and prior distributions. We used minimally informative priors and random starting values for parameters and random effects. We used Bayesian P-values and convergence diagnostics to assess model fit.
We used the resulting abundance models to estimate the potential pair abundance change relative to the range disturbance index values observed for each of the 5 species and 3 wetland class combinations over the 3 years of the study. Zero was the minimum disturbance index value for all species and was considered the control; however, the maximum value varied relative to species and wetland class. We used the average scaled value for May wetland number and percent perennial cover (i.e. mean = 0) and the species- and wetland class-specific abundance models were applied to the sample wetlands for each year.
We present log-scale posterior median estimates of parameters and credible intervals to describe both the direction and strength of relationship of parameter influences on subsequent pair estimates.
RESULTS
We conducted 15,018 ponded wetland visits from April 28 to June 11 during 2015–2017 (i.e. early count = 7,155; late count = 7,863). Most of the wetlands were surveyed 1 yr (n = 7,072); however, 2,731 wetlands were surveyed 2 yr, and 828 wetlands were surveyed all 3 yr. We identified 39,049 total breeding pairs (i.e. 16,497 Blue-winged Teal, 7,493 Gadwall, 5,694 Mallard, 2,120 Northern Pintail, and 9,858 Northern Shoveler) over the 3 survey years.
For the first stage of our analysis, the models describing Mallard abundance were the only models that were improved by removing variables from the global model. Thus, moving forward into the second stage of our analysis, we used the global model for all species except Mallard, for which we removed the wetland regime variable (Table 1). The natural log of the wetland area was a strong predictor of pair abundance for all species (Figure 2).
FIGURE 2.
Log-scale posterior medians based on 500,000 iterations of zero-inflated abundance models for Blue-winged Teal, Gadwall, Mallard, Northern Pintail, and Northern Shoveler pairs in the Missouri Coteau of northwest North Dakota and northeast Montana, USA. Parameter estimates with 95% credible intervals (CI) that do not overlap 0.0 are black. Parameter estimates with 95% CI that overlap 0.0 are gray. B1 is the intercept and represents seasonally flooded wetland; semipermanently flooded wetland (SEMI) and temporarily flooded wetland (TEMP) are intercept adjustments. Other parameter estimate positions to the left of zero indicate a negative relationship to pair abundance and positions to the right of zero indicate a positive relationship to pair abundance.
The range of disturbance index values was 0–0.3921 and varied by species (Blue-winged teal and Gadwall: 0–0.0733, Mallard: 0.0–0.0474, Northern Pintail: 0.0–0.1724, and Northern Shoveler: 0–0.3921). The effect size of the disturbance index was relatively small, and credible limits overlapped with zero for 3 of the 5 species (Figure 2). While a significant positive relationship existed between the disturbance index and abundance of Mallard and Northern Pintail pairs (Figure 3), the increase in pair abundance for the maximum positive disturbance index value for the median-sized seasonally flooded wetland was small (e.g., Mallard: 0.0026 pairs [0.58%] and Northern Pintail: 0.0031 pairs [1.24%]). The decline in predicted pair abundance for the median-size seasonally flooded wetland surveyed relative to the maximum negative disturbance index observed for Blue-winged Teal, Gadwall, and Northern Shoveler (Figure 3) was also small (e.g., Blue-winged Teal: 0.0002 pairs [–0.07%], Gadwall: 0.0028 pairs [–0.04%], and Northern Shoveler: 0.0074 pairs [–0.21%]). The change in annual predicted pair abundance for the sum of the 15 species:wetland class combinations (i.e. species = 5 and wetland classes = 3) for sample wetlands was small for all 3 yr (mean = 0.19%).
The effects of other landscape covariates on pair abundance estimates were similar among species. For all species except Blue-winged Teal, the credible intervals for the perennial cover estimate overlapped zero. We also observed a negative but nonsignificant relationship (95% credible limits did not overlap 0.0) between the number of wetlands containing ponded water and predicted numbers of breeding pairs (Figure 2).
FIGURE 3.
Model-based log-scale median predictions of wetland-level pair abundance in the Prairie Pothole Region of the Bakken Formation (2015–2017) relative to the disturbance index for (A) Mallard, (B) Northern Pintail, (C) Gadwall, (D) Blue-winged Teal, and (E) Northern Shoveler. All other covariates in the model are held constant at their mean values. Categorical variables are held to their intercept values where the wetland flooding duration is seasonal. Dotted lines represent 95% credible intervals.
DISCUSSION
We found little evidence that breeding duck pairs avoided wetlands or were less abundant in proximity to increasing levels of oil and gas development in the overlap of the Bakken and Prairie Pothole Region that we studied. Species-specific responses were mixed relative to the covariates we measured, but the effect size of the disturbance index in all models was small. Changes in pair abundance estimates for all species and sample wetlands from model-based predictions were small across the range of disturbance index values. For 4,680 wetlands in our study plots, we estimated that these changes would result in a total decrease of 0.20% for the 5 species (n = 43 pairs). This was a change that we did not consider biologically significant.
This overall change in abundance is dampened by the differences in species-specific responses (Table 1). We observed a stronger, positive effect for pairs and the disturbance index for early arriving species (i.e. Mallard and Northern Pintail) and a weaker, negative effect for mid-to late-arriving species (Northern Shoveler, Blue-winged Teal, and Gadwall). A noticeable difference between the species groups is the differing size of their home range (Baldassarre 2014). The larger home ranges of early arrivers may have allowed them to distance themselves from oil and gas development and tolerate the additional disturbance, whereas the other 3 species could not. This difference in effect may be further exacerbated by the population contribution in the region by the species where the late arrivers were 87% of the number of pairs observed in our study. Regardless, this underscores the possibility that the impact of disturbance caused by oil and gas development may have differed between the 2 groups of waterfowl and require additional study.
In contrast to our results, other studies have revealed significant negative impacts of oil and gas related to habitat loss, increased traffic, and noise on wildlife abundance across a range of habitats (Bayne et al. 2008, Doherty et al. 2008, Francis et al. 2011, Gilbert and Chalfoun 2011, Thompson et al. 2015, Barton et al. 2016). Blickley et al. (2012) found that lek attendance of male Greater Sage-grouse (Centrocercus urophsianus) was negatively related to road noise associated with energy development. Similarly, Doherty et al. (2010) demonstrated that Greater Sage-grouse avoided habitat that was otherwise suitable when in proximity to energy development, and Loesch et al. (2013) found decreased duck pair abundance in wind energy facilities.
Both grassland and wetland habitat are important for breeding waterfowl. The 5 duck species in this study establish wetland-based territories and females nest predominantly in upland habitats (Baldassarre 2014). While both habitats are crucial for production, wetlands are the predominant driver of breeding pair settling distribution and abundance (Johnson and Grier 1988, Reynolds et al. 2006, Doherty et al. 2015). Most of the habitat loss related to oil and gas development has occurred on uplands, including native prairie grassland. Preston and Kim (2016) document nearly 13,000 ha of well pads and associated infrastructure that comprise 2.3% of the landscape in the Bakken in the United States. However, only 65 ha (0.5%) of the losses were aquatic habitat. Similarly, Howden et al. (2019) rarely documented habitat loss from wetland drainage from oil and gas development. However, the ability to measure surface water depletion as part of lawful or unlawful water use for fracking is difficult to measure and it may be temporary. Because known wetland losses have been small, the amount of breeding duck pair habitat may be relatively unchanged as a result of the recent increase in oil and gas development.
Availability does not imply quality though, and wetlands close to oil and gas development might be creating ecological traps for breeding duck pairs. Spills of oil and flowback-produced water from fracking and extracting oil and gas occur regularly in the Bakken (Lauer et al. 2016). From 2005 to 2019, there were 5,589 uncontained spills documented in the North Dakota oil fields (NDEQ 2019). The high density of wetlands in portions of the Bakken and the proximity of active wells (Gleason and Tangen 2014) result in a large number of wetlands being at risk of contamination from these spills. Reiten and Teschmak (1993) and Preston and Ray (2017) found relatively high proportions of wetlands they sampled to possess indicators of brine contamination, which can negatively impact the availability of food resources important to breeding waterfowl pairs and broods (Blewett et al. 2017, Preston and Ray 2017).
We observed a negative relationship between pair abundance and the number of wetlands with ponded water regardless of proximity to oil and gas development. While this relationship seems counterintuitive, this phenomenon is documented and is likely a product of social interactions with conspecifics (Johnson and Grier 1988, Cowardin et al. 1995). Breeding waterfowl pairs have been observed to “pack” into available habitat during dry years resulting in higher pair densities. Conversely, breeding duck pairs “spread out” during wet years when more ponded water is present, and although pair densities may be lower, the overall population size may be higher because there is more available habitat.
Although we did not find evidence of impact on breeding duck pairs along the gradient of oil and gas development that we studied, the current development intensity in the Bakken is 20% of the anticipated buildout to fully exploit reserves. If the intensity and distribution of future disturbance increase past the levels examined herein, future research will be warranted. Specifically, smaller-scale investigations to monitor the health of wetland resources in the Bakken will be needed (see Brittingham et al. 2014, Gleason and Tangen 2014).
Our research on pair abundance represented the third component of a collaborative demographic investigation that assessed the effects of oil and gas development in the US portion of the Bakken on waterfowl production. The accompanying studies examined other aspects of waterfowl productivity including brood abundance (Kemink et al. 2019) and nest survival (Skaggs et al. 2020) and found that neither was negatively impacted by the intensity of oil and gas development. Collectively, we measured 2 important demographic parameters and 1 important vital rate for waterfowl population growth (Hoekman et al. 2002). While the collaborative studies did not evaluate other vital rates important to waterfowl population growth (e.g., renesting, hen success, hen survival, and brood survival; Hoekman et al. 2002), we have likely avoided pitfalls of errors in ascribing no effect from oil and gas development when effects related to other vital rates may be present (Van Horne 1983) by collecting associated pair and brood information. Given that minimal effect was observed in pair abundance, nest survival, or brood abundance, which was also considered a proxy for recruitment because the analysis in Kemink et al. (2019) used class 2 and 3 broods and survival for these age classes to fledging is high (Rotella and Ratti 1992, Sargeant and Raveling 1992, Guyn and Clark 1999), relative to oil and gas development, it appeared that population level impacts were minimal within the disturbance range, geography, and habitat conditions of these studies.
Management Implications
The conservation of waterfowl in the US Prairie Pothole Region that overlaps the Bakken focuses primarily on the protection of remaining grassland and wetland resources, along with restoration when opportunities arise. The U.S. Fish and Wildlife Service and its partners have protected >630,000 ha of grassland and wetland habitat in the Bakken, and the importance of the habitat for breeding waterfowl should be consistently reevaluated as potential stressors developed in the region (see Loesch et al. 2012). Given the results of this study and other recent studies on duck brood abundance (Kemink et al. 2019), and duck nest survival and density (Skaggs et al. 2020) in the Bakken, we conclude that the breeding waterfowl value of past conservation investments has not been negatively affected at current levels of disturbance in our study area. We recommend that existing conservation prioritization decision support tools (PPJV 2017) continue to be used but acknowledge that future information is necessary to ascertain the continued value of the conserved habitat for waterfowl and other species as the anticipated oil and gas extraction buildout in the US portion of the Bakken is approached.
ACKNOWLEDGMENTS
We thank the landowners in our study area for permission to access their land to conduct surveys and also numerous seasonal and permanent personnel who assisted with many aspects of training and field data collection. Thanks to N. Wright who developed and applied algorithms to automate water classification of spring images, Tanner Zack who conducted the water classification and created additional spatial data used in the analysis, R. Pritchert for reviewing the manuscript, and N. Niemuth and who provided a critical review of the manuscript and assistance evaluating alternative measures of disturbance. I would also like to thank the anonymous reviewers who provided both their valuable time and comments that benefited greatly improved the manuscript. The findings and conclusions in this article are those of the author(s) and do not necessarily represent the views of the U.S. Fish and Wildlife Service.
Funding statement: We received funding, personnel, and/ or logistical support for data collection, project management, lodging, and analysis from the Prairie Pothole Joint Venture, Ducks Unlimited, Inc., the U.S. Fish and Wildlife Service, the North Dakota Game and Fish Department (funded in part by U.S. Fish and Wildlife Service Pittman-Robertson Wildlife Restoration Act [CDFA# 15.611]), the Central Flyway Council, and the Enbridge Ecofootprint Grant Program via the Minnesota Research Council.
Author contributions: K.M.K. initiated discussions about waterfowl related issues in the Bakken that resulted in this project and along with C.R.L. and C.T.G. formulated the questions. C.T.G, R.C.-S., and M.S. collected data and supervised field collection, entry, and quality control of data. K.M.K. analyzed the data, and C.R.L. wrote the paper.
Data depository: Analyses reported in this article can be reproduced using the data provided by Loesch et al. (2021).
