Comparing alternative methods of modelling cumulative effects of oil and gas footprint on boreal bird abundance

doi:10.21203/rs.3.rs-1550900/v1

Context

Oil and gas activity is increasing in the western boreal forest of North America. To manage cumulative effects of this industry, a better quantification of footprint effects on wildlife is needed.

Objectives

We used point-count surveys to evaluate how well dose-response (amount) and zone-of-impact (distance) models for seismic lines, pipelines, well sites, roads, and energy facilities predicted the abundance of 48 bird species.

Methods

We developed models for each species, evaluating the best functional forms for different footprint effects, then predicted how different model structures influenced estimates of regional population size.

Results

Most species exhibited at least one nonlinear response to footprint amount (79% of species) or distance (88%). Species associated with older coniferous forests decreased more often with footprint amount and proximity to footprint, while species associated with open lands and young forests increased with footprint amount and proximity. One-third of species strongly changed in abundance at a threshold distance from at least one footprint. Zone-of-impact models had better model fit than dose-response models for 29 of 48 species, but both model types produced similar population estimates.

Conclusions

Both zone-of-impact and dose-response models were useful for assessing cumulative effects on wildlife and mechanisms causing change. Although we did not find evidence for edge effects distinct from habitat loss or gain for songbirds in our study, zone-of-impact models can provide evidence of positive or negative edge effects for developing management buffers, while dose-response models provide important information on functional changes in bird habitat with oil and gas footprint.

cumulative effects

thresholds

boreal birds

energy

oil and gas

footprint

Wildlife conservation in landscapes undergoing industrial development requires an understanding of species’ responses to combined or cumulative effects of multiple human footprints (Johnson and St.-Laurent 2011). Many wildlife studies involving human footprint place wildlife surveys within and away from different footprint types (i.e., control-impact field designs) (e.g., MacDonald 2000; Manel et al. 2000; Bart and Earnst 2002; Bayne et al. 2016). The spatial extent of most control-impact studies is typically at small scales and tends to focus on statistical significance of one footprint type rather than quantifying the effect size of different footprints (Johnson and St.-Laurent 2011). Quantitative assessments along gradients of industrial development are needed in cumulative effects assessment (CEA) to: 1) identify nonlinear effects to assess if there are threshold values of environmental stressors at which wildlife responses change (Toms and Villard 2015); 2) determine if socially accepted thresholds or criteria have been surpassed; and 3) determine relative impacts of different land uses (Bender et al. 1996; Roloff and Kernohan 1999; Burgman et al. 2011; Mahon and Pelech 2021). How to model cumulative effects is not well established, however.

One way to quantitatively assess cumulative effects of human footprint on wildlife is a “dose-response” approach (Karr and Chu 1997), in which a response (metric of habitat suitability or quality for wildlife like abundance, occupancy, nest success, or survival) is compared to the numeric amount or intensity (“dose”) of a particular stressor per sample unit. This approach provides a relatively easy way to quantify how species react to loss or change of forest due to different types of human disturbance happening within a particular area. However, vegetated and unvegetated footprint types can serve as barriers for some species or introduce edge effects (Murcia 1995) but dose-response designs do not allow us to explicitly test for the magnitude of habitat loss versus habitat fragmentation

In zone-of-impact or zone-of-influence studies, wildlife responses are measured at various distances from human disturbance(s) and the distance relationship is analyzed to determine the threshold distance where the wildlife response changes (Johnson and St.-Laurent 2011). Threshold distances may then be used to identify a buffer around the human disturbance, or zone of impact, within which the wildlife metric increases or decreases in a predictable way relative to outside the zone of impact. Zone-of-impact studies can be used to identify evidence for the extent of edge effects within habitat patches by identifying distances from the edge of the patch at which environmental variables (e.g., species abundance) within the patch strongly change (Laurance and Yensen 1991; Ewers and Didham 2007). For example, edge effects (separate from edge-related habitat loss or creation) may cause species abundance within habitat fragments to change at distances large enough for survey areas to contain no actual edge. Understanding the depth of an edge effect is an intuitive concept that helps to separate direct habitat loss from habitat fragmentation. However, whether the distance to the nearest disturbance fully quantifies the magnitude of cumulative effects in areas with multiple overlapping disturbances is not well established.

Due to overlap among footprint types within many areas, boreal wildlife is likely to be affected by multiple types of human footprint simultaneously (Johnson and St.-Laurent 2011). A large body of studies has tested for effects of forestry (e.g., Schieck and Song 2006, Huggard et al. 2014) and other human footprint types (e.g., oil and gas) on boreal forest wildlife (e.g., ungulates [Hebblewhite 2011]; small mammals [Darling et al. 2019]; birds [Bayne and Dale 2011; Bayne et al. 2016; Mahon et al. 2019]). Effects on wildlife of energy development include direct mortality; temporary or permanent loss of older forests; increases in early successional, open, and non-vegetated lands; and edge effects (e.g., changes in microclimate, noise, predation) extending from the footprint into the surrounding forest. Landscape changes caused by oil and gas development produce persistent disturbance and forest fragmentation patterns very different from what occurs in nature (Pickell et al. 2015). The different types of oil and gas footprint in boreal forests vary with respect to each other and to forest harvest activity in terms of individual footprint size, intensity of disturbance, cumulative extent, and the relative amount of edge bordering boreal forest wildlife habitats. Processing facilities and other infrastructure have the highest disturbance intensity and largest individual footprint sizes but are limited in number and cumulative extent. In contrast, linear footprint elements like seismic lines, pipelines, and roads have small footprint sizes locally, but are numerous and extensive, creating considerable edge habitat relative to both processing facilities and timber harvest (Jordaan et al. 2009). Well sites are intermediate in size, being wider than most linear features but considerably smaller than the size of most processing facilities (Boxall et al. 2005). Well sites result in lower amounts of edge than linear footprint (Jordaan et al. 2009) but are numerous and may result in considerable habitat loss (Wells et al. 2020).

In the sedimentary basin of western Canada, there have been > 30 years of studies on boreal bird habitat associations (Cumming et al. 2010). Conclusions about energy sector effects from these different study methods broadly agree with each other (e.g., bird species associated with young forests and open lands are more likely to respond positively to oil and gas disturbances [Bayne and Dale 2011; Bayne et al. 2016]). However, we do not know how well different study methods agree quantitatively, or whether they have similar predictive power when applied to landscape-level planning. Thus, our objective was to use a standardized data set collated from past boreal bird studies to directly compare dose-response and zone-of-impact methods of modeling cumulative effects of oil and gas footprint on individual boreal bird species. Dose-response methods have been used to analyze effects of the amount of different oil and gas sector footprints on abundance of boreal birds using territory densities from spot mapping plots, point counts, or point count clusters as the response (Bayne et al. 2005; Lankau et al. 2013; Mahon et al. 2019; Wilson and Bayne 2019; Sólymos et al. 2020).

The objective of our study was to compare dose-response and zone-of-impact methods of modelling energy development impacts on boreal bird abundance. First, we examined if both methods produced similar direction of effects (e.g., negative dose response effects = higher bird abundance in areas with less disturbance; negative zone-of-impact effects = higher bird abundance farther away from the nearest disturbance). Second, we determined if both methods produced similar regional population estimates of species. Finally, we compared the predictive ability of the best dose-response and zone-of-impact models for 48 boreal bird species.

Study Area

Data were collected in boreal forest landscapes north of Edmonton, Alberta (Fig. 1). Dominant vegetation consisted of: lowland bogs and fens where black spruce (Picea mariana) and tamarack (Larix laricina) are the primary tree species; mesic uplands dominated by deciduous and mixed-wood forests of trembling aspen (Populus tremuloides), balsam poplar (P. balsamifera), and white spruce (Picea glauca) with smaller amounts of balsam fir (Abies balsamea) and paper birch (Betula papyrifera); and xeric uplands dominated by jack pine (Pinus banksiana). Open vegetated lands in the study areas were associated with non-treed bogs and fens, wetlands including beaver ponds, and shrublands and grass growing up after forest fires, harvest, or in disturbances associated with oil and gas footprint (Alberta Biodiversity Monitoring Institute 2017).

Point Counts

We used point counts (N = 20,055 visits, 11,569 stations) from multiple bird studies completed between 2000–2019. Point count methods differed among studies in point count duration (3–10 minutes) and whether distance-sampling (0–50 m, 0-100 m, and 0-unlimited distance) was used (Boreal Avian Modelling Project 2015). Overall, 71% of the point counts came from control-impact studies where stations were centered within different forest stands and age classes and/or energy sector disturbances by design (e.g., Bayne et al. 2016; Mahon et al. 2016). Six per-cent of the point counts came from 3-minute Breeding Bird Survey counts along roadsides from 2000–2005, with 1–2 randomly selected visits from different years per point count station. More recent BBS data are available for Alberta but we lacked temporally accurate vegetation and footprint data for surveys after 2005. The remaining 23% of point counts were in 10 x 10 square arrays of point counts spaced 600 m apart: these arrays were located near steam-assisted gravity drainage facilities in northern Alberta and varied in the density of oil and gas footprint. Point counts in these arrays were not necessarily centered on or away from any energy sector disturbances (Fig. 1).

To account for varying detection probabilities for each species due to differences in survey methods, time of year and day of survey, and different rates of sound attenuation with distance in different vegetation types, we used the QPAD approach (Solymos et al. 2013). In brief, QPAD creates a statistical offset that adjusts the mean count based on known relationships between bird abundance and the area sampled (i.e., point count radius) and duration of point count. The QPAD approach also enabled us to predict bird abundance in a common unit of density (#/ha) (Sólymos et al. 2013). Code for QPAD is available at https://github.com/psolymos/QPAD/.

Quantifying Vegetation and Human Footprint

Forest age classes (including harvested forests) and non-forest habitats (shrubland, grassland, wetland, seismic lines, pipelines, well sites, roads, and facilities [e.g., power plants, compressor stations, processing facilities]) were quantified within 150 m of each point count. We chose this spatial scale for analysis as it is the approximate radius over which most passerine songbirds are detected using ARUs. We buffered locations in ArcGIS and intersected vegetation and footprint data from the 30-m wall-to-wall Alberta Biodiversity Monitoring Institute Vegetation Inventory (v. 6.1 vegetation/back-filled layer), correcting the intersected data for survey year (Alberta Biodiversity Monitoring Institute 2017, 2018). We summarized forest structure and composition using four variables: proportion of all non-forested lands excluding human footprint apart from harvest areas and including any stands < 10 years old (proportion open); area-weighted mean forest age within 150 m (forest age), proportion of all stands within 150 m that were conifer-dominated (proportion conifer), and proportion of all lands within 150 m that were forested or non-forested wetlands (proportion wetland). We did not model effects of forest harvest per se in our analysis because our primary focus was on comparing alternative methods of measuring oil and gas footprint effects, rather than building the best possible cumulative effects models for birds. We used the “Near” tool in ArcGIS 8.3 to calculate the nearest distance to each of these footprint types from each point count (Environmental Systems Research Institute 2002).

Modelling Each Cumulative Effects Assessment Approach

We ran models for 48 species of songbirds that were detected on at least 1% of point count visits (Mahon et al. 2016) (Table 1). These species varied in their preferences for open land versus complex older forests, coniferous versus deciduous forest, and upland versus wetland habitats (Mahon et al. 2016) and their responses to different types of energy sector development in previous studies (Bayne et al. 2016; Mahon et al. 2019).

Table 1

Number of detections, number of point counts with detections, and proportion of point counts with detections of each species in this study.
Species	Number of Birds Counted	Point Counts With Detections	Proportion of Point Counts
Alder Flycatcher Empidonax alnorum	2694	1911	0.0952
American Redstart Setophaga ruticilla	1208	874	0.0435
American Robin Turdus migratorius	1889	1631	0.0812
Bay-breasted Warbler Setophaga castanea	1080	939	0.0468
Black-and-White Warbler Mniotilta varia	771	705	0.0351
Black-capped Chickadee Poecile atricapillus	963	765	0.0381
Black-throated Green Warbler Setophaga virens	878	740	0.0369
Blue-headed Vireo Vireo solitarius	1131	1019	0.0508
Boreal Chickadee Poecile hudsonicus	761	641	0.0319
Brown Creeper Certhia americana	641	592	0.0295
Canada Jay Perisoreus canadensis	2456	2011	0.1002
Canada Warbler Cardellina canadensis	591	480	0.0239
Cape May Warbler Setophaga tigrina	763	698	0.0348
Chipping Sparrow Spizella passerina	7167	5699	0.2838
Clay-coloured Sparrow Spizella pallida	1091	880	0.0438
Common Yellowthroat Geothlypis trichas	1859	1377	0.0686
Connecticut Warbler Oporornis agilis	934	798	0.0397
Dark-eyed Junco Junco hyemalis	3342	2601	0.1295
Golden-crowned Kinglet Regulus satrapa	814	733	0.0365
Hermit Thrush Catharus guttatus	4554	3405	0.1696
Le Conte's Sparrow Ammodramus lecontei	2410	1679	0.0836
Least Flycatcher Empidonax minimus	1465	1051	0.0523
Lincoln's Sparrow Melospiza lincolnii	4029	3077	0.1533
Magnolia Warbler Setophaga magnolia	1409	1200	0.0598
Mourning Warbler Geothlypis philadelphia	1335	1074	0.0535
Northern Waterthrush Parkesia novaeboracensis	389	353	0.0176
Olive-sided Flycatcher Contopus cooperi	343	325	0.0162
Orange-crowned Warbler Vermivora celata	387	332	0.0165
Ovenbird Seiurus aurocapilla	10229	6062	0.3019
Palm Warbler Setophaga palmarum	1878	1302	0.0648
Philadelphia Vireo Vireo philadelphia	477	434	0.0216
Red-breasted Nuthatch Sitta canadensis	1719	1533	0.0764
Red-eyed Vireo Vireo olivaceus	4217	3171	0.1579
Red-winged Blackbird Agelaius phoeniceus	1050	696	0.0347
Rose-breasted Grosbeak Pheucticus ludovicianus	1660	1471	0.0733
Ruby-crowned Kinglet Regulus calendula	5034	3791	0.1888
Song Sparrow Melospiza melodia	295	239	0.0119
Swainson's Thrush Catharus olivaceus	8864	6262	0.3119
Swamp Sparrow Melospiza georgiana	1484	1096	0.0546
Tennessee Warbler Leiothlypis peregrina	13165	7201	0.3587
Warbling Vireo Vireo gilvus	784	646	0.0322
Western Tanager Piranga ludoviciana	1499	1354	0.0674
Western Wood-pewee Contopus sordidulus	447	392	0.0195
White-throated Sparrow Zonotrichia albicollis	12872	7807	0.3888
Winter Wren Troglodytes hiemalis	2346	2097	0.1044
Yellow Warbler Setophaga petechia	886	695	0.0346
Yellow-bellied Flycatcher Empidonax flaviventris	319	276	0.0137
Yellow-rumped Warbler Setophaga coronata	9929	7252	0.3612

We used generalized linear models (“glm” function in R [R Core Team 2020]) with a Poisson error distribution to create both dose-response and zone-of-impact models. We first built a set of habitat models for each species that also included a linear term for year. We considered but did not use mixed-effects models to account for correlation in point counts clustered in space and time. While some sites had multiple visits within or across years and/or were clustered spatially, the large number of random effects required, the difficulty of assigning random effects in an objective manner, and the large number of single survey point counts within our data were considered prohibitive. Instead, as latitude and longitude were not strongly correlated with other predictors at point counts, habitat models included quadratic functions of these variables to partially account for spatial correlation in the data. Different habitat models for each species included either linear or quadratic functions of the variables forest age, proportion conifer, proportion wetland, and proportion open within 150 m of point counts. Quadratic functions were assessed in case species were more abundant at intermediate values of those variables. We compared models where all four habitat variables were present along with models where one or more habitat variables were excluded and selected the model with the lowest Akaike’s Information Criterion (AIC) (Burnham and Anderson 2002) as the best habitat model for each species. Models were ranked using the “MuMIn” package in R (Barton 2009).

In dose-response models, in addition to all the predictors from the best habitat model, we used the amount of each footprint within 150 m of each point count to measure disturbance effects. We first used AIC to identify the best functional form of each footprint amount. We evaluated quadratic functions to see if species increased with intermediate amounts of footprint, inverse functions to see if species strongly increased or decreased with even small amounts of footprint, and square-root functions to see if rates of species increase or decrease tapered off beyond some value of footprint amount. To model an inverse functional relationship, all footprint proportions < 0.01 were set to 0.01. Once we determined the best functional form for seismic line, pipeline, well site, and road amount, we developed a global dose-response model containing these functional forms.

For zone-of-impact models, we used similar procedures as for testing dose-response models, except that we measured the shortest distance to each footprint type instead of the amount of each footprint. In contrast to dose-response models, negative effects of a given footprint are associated with positive coefficients for distance to that footprint. We evaluated inverse and square-root functions of footprint distance – setting the minimum distance to each footprint to 1 m when modelling inverse functions – but did not consider quadratic functions of distance amount since we did not consider it likely for bird species to reach maximum or minimum abundance at an intermediate distance from footprint. We then developed a global zone-of-impact model containing the best functional forms for each footprint. We identified possible evidence of edge effects of footprint when species were best explained by an inverse function of distance to that footprint, but only if a sharp change in edge response occurred at a distance > 150 m. At that distance, there was none of that footprint within 150 m of a point count and therefore no habitat lost to that footprint. Although habitat could be lost to footprint at larger scales, 150 m was a large enough distance to contain multiple territories of the songbird species we analyzed.

Predictive Ability and Population Consequences of Each Modelling Approach

We compared the ability of different approaches to accurately predict boreal bird responses to cumulative energy sector development. We used AIC and pseudo-R² as measures of explanatory power and relative goodness of model fit of the final habitat, dose-response, and zone-of-impact models chosen for each species.

Using the final dose-response and zone-of-impact model for each species, we predicted total abundance in a landscape of adjacent 300-m cells in a 6-million-hectare area of NE Alberta that aligns with the Alberta Pacific Forest Industries Ltd. Forest Management Area (Al-Pac FMA). We generated 300-m resolution raster layers (i.e., ~equal to the width of the 150-m radius point counts that we modelled) of each model predictor term that potentially could be included in the final models generated for each species. We created raster layers of forest vegetation variables from national 250-m layers clipped to Alberta, resampling via nearest-neighbours methods to 300-m resolution (Beaudoin et al. 2014). We created raster layers from human footprint feature layers, including but not limited to the footprint types within our models (Alberta Biodiversity Monitoring Institute 2018). Although forest harvest was among these human footprint types and is the largest source of human disturbance within the Al-Pac FMA, we assumed that its main effects on forest structure were to reduce mean forest age and increase the proportion of open vegetated lands within the raster layers of forest vegetation. Forest harvest activity also changes the temporal and spatial patterns of forest stand ages and types across the landscape, which will influence the total availability of suitable habitat for boreal birds associated with different forest types, ages, and structural stages. However, these spatial patterns of forest stand ages and types would be captured at least at a coarse spatial scale by the raster layers we used. We also estimated the shortest distance from the centroid of each cell to each of the five footprint types analyzed (seismic lines, pipelines, well sites, roads, facilities), using the “Near” proximity tool in ArcGIS, and created five 300-m distance raster layers, where the value in each cell was the nearest distance to one of the five footprints.

Predicted density layers for each species (#males/ha) was multiplied by 9 (# hectares per 300-m cell), then abundance in each cell was summed for the whole study area, excluding the cells containing human footprint that were not included in our models (e.g., railways, gravel pits, towns, cropland, mines, transmission lines: 116554 out of 789259 cells or 15% of the cells within the Al-Pac FMA, mostly along major highways, and within lakes, settlement areas, and mined areas north of Fort McMurray, Alberta). We calculated correlations and plotted the relationships between predicted population estimates from the best habitat, dose-response, and zone-of-impact models, using each species as an observation.

Direction and Magnitude of Oil and Gas Footprint Effects

A larger proportion of species associated with young, productive ecosystems or with older deciduous forests (Mahon et al. 2016) responded positively to footprint amount within 150 m (34% of species X footprint effects) and closer footprints (42%) than responded negatively to footprint amount (13%) or proximity (22%). In contrast, a larger proportion of species associated with older coniferous upland forests responded negatively to footprint amount (34%) and proximity (36%) than responded positively to amount (8%) or proximity (12%). For species associated with either older deciduous upland forests or with lowland coniferous forests, there were slightly more positive than negative species responses to footprint amount and slightly more negative than positive species responses to proximity to footprint (Tables 2–3, model coefficients and plots in Appendices 1–2).

Table 2

The “best” model function for describing the dose-response effect of the proportion of each footprint type within 150 m, in the global dose-response models predicting 48 species of boreal songbirds in northern Alberta. Type of function: I = inverse; L = linear; Q = quadratic; S = square root; 0 = no effect [footprint not included in the best model] or the overall effect was non-significant (i.e., 95% confidence intervals for effect size overlapped with zero) even if it was included in the best dose-response model; (-) = negative effect of footprint; (+) = positive effect of footprint on species abundance. Habitat associations (from Mahon et al. 2016): YP = young productive ecosystems (young deciduous forests, shrublands, open lands, marsh); OD = older, structurally complex deciduous forests; OC = older, structurally complex upland coniferous forests; LC = lowland coniferous forests. Species in italics were not mentioned in Mahon et al. (2016) but were assigned to those habitat associations based on their ecology.
Species	Habitat	Seismic	Pipeline	Well	Road	Facility
Alder Flycatcher	YP	I(-)	S(+)	S(+)	S(+)	S(+)
American Robin	YP	0	S(+)	S(+)	Q(+)	Q
Clay-coloured Sparrow	YP	0	0	I(+)	S(+)	S(+)
Common Yellowthroat	YP	0	0	I(+)	Q(+)	0
Hermit Thrush	YP	0	0	0	0	L(-)
Le Conte's Sparrow	YP	Q(-)	S(-)	I(+)	I(-)	I(-)
Lincoln's Sparrow	YP	0	I(+)	S(+)	Q(+)	I(+)
Magnolia Warbler	YP	0	0	S(-)	0	0
Northern Waterthrush	YP	0	I(-)	S(+)	Q	L(-)
Orange-crowned Warbler	YP	0	0	S	S(+)	0
Olive-sided Flycatcher	YP	0	0	0	0	S
Red-winged Blackbird	YP	I(-)	Q(+)	I(+)	Q(+)	L
Song Sparrow	YP	I	S	Q	S(+)	S(+)
Swamp Sparrow	YP	0	0	S(+)	Q(+)	0
Western Wood-pewee	YP	Q	0	0	Q(+)	Q
White-throated Sparrow	YP	Q(-)	0	0	Q(+)	0
Yellow-bellied Flycatcher	YP	L	S	I	L(+)	L
American Redstart	OD	0	0	0	L(+)	I(-)
Black-and-White Warbler	OD	0	0	0	0	0
Black-capped Chickadee	OD	I(-)	0	Q(-)	L(+)	Q
Canada Warbler	OD	S(-)	L(-)	L	L	L
Connecticut Warbler	OD	I(-)	0	I	0	I(-)
Least Flycatcher	OD	Q(+)	0	Q(+)	L(+)	Q(+)
Mourning Warbler	OD	0	0	0	0	0
Ovenbird	OD	0	0	0	Q(-)	0
Philadelphia Vireo	OD	0	I(+)	0	0	L
Red-eyed Vireo	OD	0	0	0	L(+)	0
Rose-breasted Grosbeak	OD	L(-)	0	Q(-)	L	I(-)
Warbling Vireo	OD	0	Q(+)	I	I(+)	0
Yellow Warbler	OD	I(-)	Q	I(+)	Q(+)	I
Bay-breasted Warbler	OC	0	0	S(-)	L(-)	I
Black-throated Green Warbler	OC	L(-)	I(-)	I(-)	Q(-)	I
Blue-headed Vireo	OC	0	0	0	L(-)	0
Boreal Chickadee	OC	Q(+)	Q(+)	S	0	0
Brown Creeper	OC	S	I	L(-)	S(-)	L
Cape May Warbler	OC	0	0	0	Q(+)	S(-)
Golden-crowned Kinglet	OC	0	I(-)	L(-)	L(-)	0
Red-breasted Nuthatch	OC	L(-)	0	0	0	0
Swainson's Thrush	OC	0	0	0	0	0
Tennessee Warbler	OC	I(-)	0	0	0	0
Western Tanager	OC	S(-)	L(-)	S(-)	L	I(+)
Winter Wren	OC	0	0	S(-)	Q(-)	L(-)
Yellow-rumped Warbler	OC	I(+)	0	0	Q(-)	0
Canada Jay	LC	L(-)	Q(+)	Q(+)	0	S(-)
Chipping Sparrow	LC	0	I(+)	S(+)	L(+)	0
Dark-eyed Junco	LC	0	Q(+)	Q(+)	0	0
Palm Warbler	LC	S(+)	0	Q(-)	S(-)	Q(-)
Ruby-crowned Kinglet	LC	Q(+)	L(+)	L(+)	L(+)	I(+)

Table 3

The “best” model function for describing the zone-of-impact effect of the nearest distance to different footprint types from each point count, from the global zone-of-impact models predicting 48 species of boreal songbirds in northern Alberta. Type of function: I = inverse; L = linear; S = square root; 0 = no effect [footprint not included in the best model] or the overall effect was non-significant (i.e., 95% confidence intervals for effect size overlapped with zero) even if it was included in the best dose-response model; (-) = negative effect of footprint; (+) = positive effect of footprint on species abundance. In contrast to dose-response models, positive effects result when bird abundance decreases with footprint distance and negative effects result when bird abundance increases with distance from footprint. * = implausible model coefficient estimates from global model, based on the “best” functional forms of each footprint. Habitat associations (from Mahon et al. 2016): YP = young productive ecosystems (young deciduous forests, shrublands, open lands, marsh); OD = older, structurally complex deciduous forests; OC = older, structurally complex upland coniferous forests; LC = lowland coniferous forests. Species in italics were not mentioned in Mahon et al. (2016) but were assigned to those habitat associations based on their ecology.
Species	Habitat Association	Seismic	Pipeline	Well	Road	Facility
Alder Flycatcher	YP	0	0	S(+)	0	S(+)
American Robin	YP	0	0	S(+)	S(+)	S(+)
Clay-coloured Sparrow	YP	0	S(-)	S(+)	S(+)	S(+)
Common Yellowthroat	YP	L(+)	L(-)	S(+)	S(+)	0
Hermit Thrush	YP	S(+)	S(+)	L(+)	0	L(-)
Le Conte's Sparrow	YP	L(+)	0	S(+)	I(+)	L(-)
Lincoln's Sparrow	YP	S(+)	S(-)	S(+)	0	L(+)
Magnolia Warbler	YP	0	L(-)	S(+)	0	L(+)
Northern Waterthrush	YP	0	S(-)	0	I(+)	L(+)
Orange-crowned Warbler	YP	I(-)	0	S(+)	S(+)	0
Olive-sided Flycatcher	YP	L(+)	L	S(+)	S(-)	I(-)
Red-winged Blackbird	YP	0	I(-)	L(-)	S(+)	S(+)
Song Sparrow	YP	I(-)	0	0	S(+)	S(+)
Swamp Sparrow	YP	L(+)	S(-)	0	S(-)	0
Western Wood-pewee	YP	S(-)	0	S(+)	0	L(-)
White-throated Sparrow	YP	S(-)	L(-)	S(+)	0	S(+)
Yellow-bellied Flycatcher	YP	0	0	0	0	0
American Redstart	OD	L(-)	I(-)	L(+)	L(-)	0
Black-and-White Warbler	OD	0	S(-)	0	0	0
Black-capped Chickadee	OD	S(-)	I(+)	0	S(+)	L(-)
Canada Warbler	OD	S(-)	0	0	S(-)	0
Connecticut Warbler	OD	L(+)	0	L(-)	L(+)	0
Least Flycatcher	OD	0	I(-)	S(+)	S(+)	L(-)
Mourning Warbler	OD	S(-)	0	L(-)	L(-)	0
Ovenbird	OD	S(-)	L(-)	L(-)	0	0
Philadelphia Vireo	OD	0	0	0	0	0
Red-eyed Vireo	OD	I(+)	0	L(-)	S(+)	S(+)
Rose-breasted Grosbeak	OD	L(-)	L(-)	S(+)	0	0
Warbling Vireo	OD	L(-)	0	S(+)	S(+)	L(-)
Yellow Warbler	OD	0	S(-)	S(+)	S(+)	I(+)
Bay-breasted Warbler	OC	S(-)	0	L(-)	0	0
Black-throated Green Warbler	OC	0	0	0	I(-)	S(-)
Blue-headed Vireo	OC	0	I(+)	0	0	L(-)
Boreal Chickadee	OC	0	I(+)	0	L(-)	L(-)
Brown Creeper	OC	0	0	S(-)	L(-)	0
Cape May Warbler	OC	S(-)	0	L(-)	0	0
Golden-crowned Kinglet	OC	0	0	0	0	L(+)
Red-breasted Nuthatch	OC	S(-)	0	S(+)	0	L(-)
Swainson's Thrush	OC	0	L(-)	L(-)	0	L(+)
Tennessee Warbler	OC	S(-)	S(-)	S(+)	L(-)	0
Western Tanager	OC	S(-)	0	L(-)	0	S(-)
Winter Wren	OC	0	0	L(-)	0	0
Yellow-rumped Warbler	OC	S(+)	0	L(-)	I(-)	L(+)
Canada Jay	LC	L(-)	S(+)	S(+)	L(-)	L(-)
Chipping Sparrow	LC	S(+)	0	L(-)	L(-)	S(+)
Dark-eyed Junco	LC	S(+)	S(+)	L(-)	0	S(-)
Palm Warbler	LC	S(+)	0	L(-)	S(-)	S(-)
Ruby-crowned Kinglet	LC	S(+)	L(-)	0	0	S(+)

Most species exhibited at least one non-linear response to footprint in both dose-response models (79% or 38 of 48 species) and zone-of-impact models (87.5% or 42 of 48 species). One-third of species varied with an inverse function of distance to at least one footprint type, indicating that abundance of these species strongly changed at a threshold distance from footprint. However, we found that few species significantly changed in abundance (i.e., 95% confidence intervals of effect size excluding zero) at threshold distances greater than 150 m from footprint, suggesting a lack of strong edge effects separate from habitat loss (Tables 3–4, Appendix 2).

Table 4

Approximate thresholds to the nearest 100 m for different footprint types from the global zone-of-impact models for 48 species of boreal songbirds in northern Alberta, for those species that varied as an inverse function of distance from a given footprint. Thresholds indicate where abundance of a species strongly changes with increasing distance. Cells are coloured according to the overall effect of each footprint type on abundance of each species (green = positive [abundance decreases with increasing distance from footprint], yellow = negative [abundance increases with increasing distance from footprint], no colour plus “-“ = inverse function for a particular footprint was not included in the best zone of impact model for a given type of footprint for a particular bird species, so there was no threshold.
Species	Habitat Association	Seismic	Pipeline	Well	Road	Facility
Alder Flycatcher	YP	-	-	-	-	-
American Robin	YP	-	-	-	-	-
Clay-coloured Sparrow	YP	-	-	-	-	-
Common Yellowthroat	YP	-	-	-	-	-
Hermit Thrush	YP	-	-	-	-	-
Le Conte's Sparrow	YP	-	-	-	500	-
Lincoln's Sparrow	YP	-	-	-	-	-
Magnolia Warbler	YP	-	-	-	-	-
Northern Waterthrush	YP	-	-	-	0	-
Orange-crowned Warbler	YP	100	-	-	-	-
Olive-sided Flycatcher	YP	-	-	-	-	1000
Red-winged Blackbird	YP	-	300	-	-	-
Song Sparrow	YP	0	-	-	-	-
Swamp Sparrow	YP	-	-	-	-	-
Western Wood-pewee	YP	-	-	-	-	-
White-throated Sparrow	YP	-	-	-	-	-
Yellow-bellied Flycatcher	YP	-	-	-	-	-
American Redstart	OD	-	200	-	-	-
Black-and-White Warbler	OD	-	-	-	-	-
Black-capped Chickadee	OD	-	0	-	-	-
Canada Warbler	OD	-	-	-	-	-
Connecticut Warbler	OD	-	-	-	-	-
Least Flycatcher	OD	-	100	-	-	-
Mourning Warbler	OD	-	-	-	-	-
Ovenbird	OD	-	-	-	-	-
Philadelphia Vireo	OD	-	-	-	-	-
Red-eyed Vireo	OD	0	-	-	-	-
Rose-breasted Grosbeak	OD	-	-	-	-	-
Warbling Vireo	OD	-	-	-	-	-
Yellow Warbler	OD	-	-	-	-	0
Bay-breasted Warbler	OC	-	-	-	-	-
Black-throated Green Warbler	OC	-	-	-	100	-
Blue-headed Vireo	OC	-	0	-	-	-
Boreal Chickadee	OC	-	0	-	-	-
Brown Creeper	OC	-	-	-	-	-
Cape May Warbler	OC	-	-	-	-	-
Golden-crowned Kinglet	OC	-	-	-	-	-
Red-breasted Nuthatch	OC	-	-	-	-	-
Swainson's Thrush	OC	-	-	-	-	-
Tennessee Warbler	OC	-	-	-	-	-
Western Tanager	OC	-	-	-	-	-
Winter Wren	OC	-	-	-	-	-
Yellow-rumped Warbler	OC	-	-	-	100	-
Canada Jay	LC	-	-	-	-	-
Chipping Sparrow	LC	-	-	-	-	-
Dark-eyed Junco	LC	-	-	-	-	-
Palm Warbler	LC	-	-	-	-	-
Ruby-crowned Kinglet	LC	-	-	-	-	-

Predictive Ability and Population Consequences

The global zone of impact model for each species had a better relative fit (lower AIC) than the habitat or global dose-response model and higher absolute fit as measured by the percentage of deviance explained in all point counts (pseudo-R²), in 29 of 48 species (Table 5). For most species, the relative difference in population size predicted by each type of global model was small (relative ratio ~ 1 [Fig. 3], Table 6). However, there were also some dramatic differences in population estimates for some noteworthy species. Black-throated Green Warbler – a species of special concern in Alberta – was predicted to be 20% more abundant under its best dose-response model and 38% less abundant under the zone-of-impact model relative to the habitat model with no footprint predictors (Table 6). Differences in density distribution maps based on the best habitat, dose-response, or zone-of-impact models were subtle for some species but not for others (e.g., Boreal Chickadee) (Fig. 4, Appendix 3). Nevertheless, there were strong positive correlations between population estimates from the habitat, dose-response and zone-of-impact models for each species overall (Spearman r > 0.9 for all pairwise comparisons).

Table 5

Comparison of AIC statistics and % deviance explained in all point counts from the habitat model, global dose-response model, and global zone-of-impact model for each of 48 boreal songbird species. Cell colour indicates which model had the lowest AIC and highest % deviance explained.
	Habitat Model		Global Dose Response		Global Zone of Impact
Species	AIC	Explained	AIC	Explained	AIC	Explained
Alder Flycatcher	16294	9	16055	11	16144	10
American Redstart	8370	20	8333	20	8340	20
American Robin	12573	8	12110	13	12303	11
Bay-breasted Warbler	7693	15	7627	16	7568	17
Black-and-White Warbler	6293	7	6308	7	6274	8
Black-capped Chickadee	7547	11	7508	11	7489	11
Black-throated Green Warbler	6049	25	5909	28	5985	26
Blue-headed Vireo	8383	8	8376	8	8344	8
Boreal Chickadee	6263	9	6248	9	6213	10
Brown Creeper	5220	13	5181	14	5170	15
Canada Jay	14910	7	14751	8	14116	14
Canada Warbler	4780	18	4711	20	4717	20
Cape May Warbler	5931	13	5895	14	5827	15
Chipping Sparrow	30830	3	30680	4	30761	3
Clay-coloured Sparrow	7748	16	7493	20	7563	19
Common Yellowthroat	12375	9	12253	11	12268	10
Connecticut Warbler	7513	7	7499	7	7436	8
Dark-eyed Junco	18878	7	18864	8	18769	8
Golden-crowned Kinglet	5758	23	5727	24	5720	24
Hermit Thrush	22659	10	22654	10	22525	10
Le Conte's Sparrow	12404	32	12184	34	11967	36
Least Flycatcher	10127	15	10029	16	10091	15
Lincoln's Sparrow	21633	5	21309	7	21429	6
Magnolia Warbler	10322	3	10320	4	10262	4
Mourning Warbler	9179	15	9219	15	9146	15
Northern Waterthrush	3787	4	3738	6	3743	6
Orange-crowned Warbler	3643	10	3619	11	3614	11
Olive-sided Flycatcher	3312	7	3309	7	3259	9
Ovenbird	36139	20	36136	21	35831	21
Palm Warbler	11210	22	11130	23	11155	23
Philadelphia Warbler	4264	10	4255	10	4258	10
Red-breasted Nuthatch	10475	23	10458	23	10420	23
Red-eyed Vireo	21934	10	21729	11	21622	12
Red-winged Blackbird	8254	12	7900	17	8091	14
Rose-breasted Grosbeak	10778	13	10713	14	10720	14
Ruby-crowned Kinglet	22077	18	21859	20	21831	20
Song Sparrow	2759	17	2440	30	2515	27
Swainson's Thrush	33628	13	33665	13	33461	14
Swamp Sparrow	10274	11	10263	11	10252	11
Tennessee Warbler	44000	14	43993	14	43427	16
Warbling Vireo	5063	35	5022	36	5014	36
Western Tanager	9641	17	9674	17	9573	18
Western Wood-pewee	3980	9	3948	11	3940	11
White-throated Sparrow	42609	14	42456	15	42430	14
Winter Wren	13558	14	13456	16	13319	17
Yellow-bellied Flycatcher	2887	16	2878	16	2900	15
Yellow Warbler	6454	20	6289	23	6323	23
Yellow-rumped Warbler	35243	12	35096	13	35109	13

Table 6

Predicted population size and population ratios (estimate + 95% prediction interval) for 48 boreal songbird species within the Alberta-Pacific Forest Industries Ltd. Forest Management Area in northeastern Alberta, Canada. Predictions were generated from the best habitat, dose-response, and zone-of-impact models for each species. DR = dose-response model. ZOI = zone-of-impact model.
Species	Habitat	Dose-response	Zone-of-impact	DR:Habitat	ZOI:Habitat	DR:ZOI
Alder Flycatcher	348763 (306894–397473)	322090 (279263–372993)	380632 (316998–458831)	0.92 (0.91–0.94)	1.09 (1.03–1.15)	0.85 (0.88 − 0.81)
American Redstart	348556 (264070–466952)	356294 (266492–483680)	379464 (280208–522423)	1.02 (1.01–1.04)	1.09 (1.06–1.12)	0.94 (0.95 − 0.93)
American Robin	189630 (164273–219463)	152984 (130195–180413)	158867 (128735–196781)	0.81 (0.79–0.82)	0.84 (0.78–0.9)	0.96 (1.01–0.92)
Bay-breasted Warbler	843837 (700781–1018510)	922871 (759883–1123720)	992065 (790122–1251229)	1.09 (1.08–1.1)	1.18 (1.13–1.23)	0.93 (0.96 − 0.9)
Black-and-white Warbler	252580 (190383–339950)	260525 (192465–358134)	278720 (207747–379506)	1.03 (1.01–1.05)	1.1 (1.09–1.12)	0.93 (0.93–0.94)
Black-capped Chickadee	283349 (219277–377563)	280030 (213339–378321)	175179 (123557–257934)	0.99 (0.97-1)	0.62 (0.56–0.68)	1.6 (1.73–1.47)
Black-throated Green Warbler	129454 (102650–178863)	155021 (118647–213766)	79770 (58070–124506)	1.2 (1.16–1.2)	0.62 (0.57–0.7)	1.94 (2.04–1.72)
Blue-headed Vireo	321854 (270231–387662)	339809 (283087–414020)	212734 (162845–281805)	1.06 (1.05–1.07)	0.66 (0.6–0.73)	1.6 (1.74–1.47)
Boreal Chickadee	609352 (520426–714999)	588880 (494231–703781)	371526 (272198–511881)	0.97 (0.95–0.98)	0.61 (0.52–0.72)	1.59 (1.82–1.37)
Brown Creeper	275219 (202858–371323)	298508 (219532–402358)	305625 (214211–435773)	1.08 (1.08–1.08)	1.11 (1.06–1.17)	0.98 (1.02–0.92)
Canada Jay	916166 (830246–1014202)	885390 (794613–990070)	162022 (121741–218650)	0.97 (0.96–0.98)	0.18 (0.15–0.22)	5.46 (6.53–4.53)
Canada Warbler	126662 (97698–165207)	112808 (85289–150427)	143537 (105952–196607)	0.89 (0.87–0.91)	1.13 (1.08–1.19)	0.79 (0.8 − 0.77)
Cape May Warbler	539213 (434594–672723)	605578 (476938–780982)	559775 (408688–772261)	1.12 (1.1–1.16)	1.04 (0.94–1.15)	1.08 (1.17–1.01)
Chipping Sparrow	1558799 (1466186–1658307)	1497725 (1395418–1608750)	1693056 (1532234–1874281)	0.96 (0.95–0.97)	1.09 (1.05–1.13)	0.88 (0.91 − 0.86)
Clay-coloured Sparrow	167314 (132507–212913)	148240 (114094–194403)	195159 (142743–269060)	0.89 (0.86–0.91)	1.17 (1.08–1.26)	0.76 (0.8 − 0.72)
Common Yellowthroat	174634 (147488–207660)	165929 (137881–200696)	167535 (132349–214305)	0.95 (0.93–0.97)	0.96 (0.9–1.03)	0.99 (1.04–0.94)
Connecticut Warbler	158047 (121549–206123)	172907 (131117–226321)	130360 (92562–185410)	1.09 (1.08–1.1)	0.82 (0.76–0.9)	1.33 (1.42–1.22)
Dark-eyed Junco	904452 (815172–1004766)	851558 (758281–957738)	714335 (620580–823881)	0.94 (0.93–0.95)	0.79 (0.76–0.82)	1.19 (1.22–1.16)
Golden-crowned Kinglet	493581 (395319–617680)	527532 (419481–665167)	808540 (578995–1132844)	1.07 (1.06–1.08)	1.64 (1.46–1.83)	0.65 (0.72 − 0.59)
Hermit Thrush	377904 (348564–410016)	382824 (350579–418406)	393800 (350117–443581)	1.01 (1.01–1.02)	1.04 (1-1.08)	0.97 (1-0.94)
Le Conte's Sparrow	280577 (231126–341326)	313209 (252838–389019)	157696 (124155–206117)	1.12 (1.09–1.14)	0.56 (0.54–0.6)	1.99 (2.04–1.89)
Least Flycatcher	610765 (448063–826712)	604401 (435165–827969)	503097 (348143–726418)	0.99 (0.97-1)	0.82 (0.78–0.88)	1.2 (1.25–1.14)
Lincoln's Sparrow	709624 (647326–778906)	664527 (598757–738833)	878034 (758889–1018157)	0.94 (0.92–0.95)	1.24 (1.17–1.31)	0.76 (0.79 − 0.73)
Magnolia Warbler	350834 (303712–405880)	351708 (300480–412455)	458658 (363757–582311)	1 (0.99–1.02)	1.31 (1.2–1.43)	0.77 (0.83 − 0.71)
Mourning Warbler	250809 (199711–316976)	221277 (174025–283129)	272511 (209808–357496)	0.88 (0.87–0.89)	1.09 (1.05–1.13)	0.81 (0.83 − 0.79)
Northern Waterthrush	46719 (31661–71194)	44896 (29013–72149)	80292 (48946–136048)	0.96 (0.92–1.01)	1.72 (1.55–1.91)	0.56 (0.59 − 0.53)
Olive-sided Flycatcher	29340 (21700–40112)	26015 (18522-1327611)	31624 (20649–50147)	0.89 (0.85–33.1)	1.08 (0.95–1.25)	0.82 (0.9-26.47)
Orange-crowned Warbler	141203 (105658–189755)	124291 (91471–170029)	122513 (84584–180866)	0.88 (0.87–0.9)	0.87 (0.8–0.95)	1.01 (1.08–0.94)
Ovenbird	1256593 (1161004–1359283)	1295086 (1192411–1405841)	1335196 (1217904–1463388)	1.03 (1.03–1.03)	1.06 (1.05–1.08)	0.97 (0.98 − 0.96)
Palm Warbler	528072 (456975–611657)	470871 (399873–556341)	521626 (418034–653393)	0.89 (0.88–0.91)	0.99 (0.91–1.07)	0.9 (0.96 − 0.85)
Philadelphia Vireo	144838 (93454–238525)	149299 (95267–248413)	148239 (94714–248106)	1.03 (1.02–1.04)	1.02 (1.01–1.04)	1.01 (1.01-1)
Red-breasted Nuthatch	332694 (285061–388865)	344542 (294298–403980)	228782 (180454–290905)	1.04 (1.03–1.04)	0.69 (0.63–0.75)	1.51 (1.63–1.39)
Red-eyed Vireo	581469 (510424–664593)	561522 (487617–648881)	736142 (627189–866570)	0.97 (0.96–0.98)	1.27 (1.23–1.3)	0.76 (0.78 − 0.75)
Red-winged Blackbird	70610 (57573–87063)	55261 (43817–70112)	58035 (44247–76615)	0.78 (0.76–0.81)	0.82 (0.77–0.88)	0.95 (0.99 − 0.92)
Rose-breasted Grosbeak	267691 (219063–329376)	259918 (211279–322026)	276673 (210316–367050)	0.97 (0.96–0.98)	1.03 (0.96–1.11)	0.94 (1-0.88)
Ruby-crowned Kinglet	1399546 (1292859–1516002)	1175583 (1074104–1287584)	1763560 (1558252–2000667)	0.84 (0.83–0.85)	1.26 (1.21–1.32)	0.67 (0.69 − 0.64)
Song Sparrow	30378 (20700–44886)	15392 (10013–23885)	39307 (24071–65013)	0.51 (0.48–0.53)	1.29 (1.16–1.45)	0.39 (0.42 − 0.37)
Swainson's Thrush	1418913 (1330182–1516958)	1444864 (1354180–1545054)	1687675 (1544213–1848335)	1.02 (1.02–1.02)	1.19 (1.16–1.22)	0.86 (0.88 − 0.84)
Swamp Sparrow	143853 (122282–169644)	147346 (123812–175885)	167240 (134444–209826)	1.02 (1.01–1.04)	1.16 (1.1–1.24)	0.88 (0.92 − 0.84)
Tennesee Warbler	3464849 (3300262–3633908)	3699161 (3505012–3900618)	3971175 (3676151–4280985)	1.07 (1.06–1.07)	1.15 (1.11–1.18)	0.93 (0.95 − 0.91)
Warbling Vireo	150615 (122381–186283)	140083 (112738–175328)	98572 (68280–143647)	0.93 (0.92–0.94)	0.65 (0.56–0.77)	1.42 (1.65–1.22)
Western Tanager	282166 (235047–330866)	291334 (239577–346007)	279280 (215030–353022)	1.03 (1.02–1.05)	0.99 (0.91–1.07)	1.04 (1.11–0.98)
Western Wood-pewee	65291 (44131–98605)	72540 (48995–109449)	48387 (29661–80610)	1.11 (1.11–1.11)	0.74 (0.67–0.82)	1.5 (1.65–1.36)
White-throated Sparrow	1317210 (1227218–1414679)	1352129 (1256540–1456079)	1540484 (1403722–1691483)	1.03 (1.02–1.03)	1.17 (1.14–1.2)	0.88 (0.9 − 0.86)
Winter Wren	258177 (229086–290453)	285241 (250918–323183)	257014 (220483–299493)	1.1 (1.1–1.11)	1 (0.96–1.03)	1.11 (1.14–1.08)
Yellow Warbler	244717 (170212–357106)	244448 (164086–371528)	220625 (147507–334375)	1 (0.96–1.04)	0.9 (0.87–0.94)	1.11 (1.11–1.11)
Yellow-bellied Flycatcher	91862 (59488–125545)	89734 (57129–125128)	60397 (33366–96351)	0.98 (0.96-1)	0.66 (0.56–0.77)	1.49 (1.71–1.3)
Yellow-rumped Warbler	3551299 (3350832–3766121)	3647058 (3425230–3885540)	3787852 (3493494–4109784)	1.03 (1.02–1.03)	1.07 (1.04–1.09)	0.96 (0.98 − 0.95)

We compared regional-scale predictions of boreal bird abundance as a function of oil and gas footprint based on dose-response versus zone-of-impact models and evaluated the predictive accuracy of both sets of models. Zone-of-impact models usually had higher explanatory power than dose-response models for a given species. As well, footprint impacts were detected for more bird species using zone-of-impact than dose-response models. One reason for the higher explanatory power of zone-of-impact models may be that dose-response models considered amounts of footprint only within 150 m of point counts rather than at larger, potentially influential spatial scales (e.g., Mahon et al. 2019). In contrast, measuring effects of distance to footprint ⸺while not linked to a particular spatial scale⸺probably captured landscape-scale influences of footprint that were not present in dose-response models.

That dose-response and zone-of-impact models generally predicted similar regional population outcomes was not surprising given the reasonably strong correlations between amount and nearest distance to footprint for most stressor types. There were strong negative correlations (Spearman r < -0.5) between the proportion of each footprint within 150 m and nearest distance to that footprint type for most stressors (seismic: -0.70; pipeline: -0.51; well site: -0.51; roadside: -0.50), but there was a weaker negative correlation between the proportion of and distance to facility footprint (-0.35). This weaker correlation for facilities suggests that relationships between footprint amount and distance are conditional on how a study is designed and analyzed. Our collated data set included data from different studies designed to test specific gradients of energy sector footprint, distances to edge, and other objectives. Thus, we did not have a balanced experimental design with similar numbers of surveys across different footprint types, systematically located across our study area. Further, modelling dose-response and zone-of-impact effects at local spatial scales probably weakened the correlation between facility amount and distance to facility. Facility footprint is the least extensive at large spatial scales, being less numerous and more concentrated in smaller areas than other footprint types (Boxall et al. 2005). As a result, there is a much larger range of distances to the nearest facility than any other footprint type. Facilities can be very large and are more likely to cover the entire area sampled for birds (footprint proportions range from 0 to 1 for facilities) than the other footprint types (i.e., footprint proportions range from 0 to 0.25 for seismic). Recognizing that inherently different ranges are possible for footprint types suggests that different approaches to modelling footprint types might be required when creating large-scale regional models of cumulative effects (e.g., Sólymos et al. 2020). Dose-response models may be more useful for examining additive and interactive effects of multiple widely distributed footprint/disturbance types in complex multi-stressor and multi-sector landscapes where footprint types are not spatially or temporally independent (Mahon et al. 2019). In areas with lower amounts of footprint types, greater dispersion of footprint types, or fewer footprint types, zone-of-impact models might be useful because distance effects can be examined independently (spatially independent) and accurately.

The strong correlation between dose-response and zone-of-impact measurements for other footprint types suggests that either model could be used in regulatory decisions. While the improved model fit and greater number of significant zone-of-impact relationships suggest that zone-of-impact approaches are “better”, we caution against relying entirely on one framework over the other. We suggest instead that creating and presenting both kinds of models can provide important complementary insights into the processes generating impacts of human footprint on birds. If a primary goal is to understand cumulative effects of different footprint types, dose-response models may be more appropriate as the total area disturbed by different combinations of footprint at varying spatial scales can be assessed, which can not be done using zone-of-impact models. However, the mechanisms (i.e., habitat loss versus fragmentation) are more challenging to separate with dose-response models.

Zone-of-impact models may be more useful for detecting potential edge effects of footprint on species of wildlife (Laurance and Yensen 1991; Ewers and Didham 2007). There is stronger evidence for an edge effect of footprint if a species’ abundance strongly changes at distances large enough that there is no footprint within the species’ territories within survey areas. For songbirds, evidence of edge effects could be associated with threshold distances greater than the radius of a point count (~ 150 m): a 150-m point count is large enough to contain multiple territories of most forest songbird species. Further, many forest songbird species have song characteristics that make it unlikely to 1) detect those species within forests beyond 150 m using point counts, and 2) detecting conspecifics using habitat including footprint outside of the point count radius (Matsuoka et al. 2012). Cumulative edge effects could be significant given the amount of habitat edge created by well-sites, linear footprint, and forest harvest in boreal forests (Wells et al. 2020) and could be positive or negative for different species (Murcia 1995). However, although we did find strong changes in abundance of species at threshold distances, these distances were generally not large enough to exclude the possibility of habitat loss or creation by footprint. Even if we had identified significantly large zone-of-impact thresholds, the underlying mechanisms of edge effects (e.g., noise, microclimate changes, resources, predators) for different species in this paper have not been identified. A possible alternative explanation is that as there are multiple human footprint types (including harvest) that are close to or overlapping with each other, i.e., impacts which may not be independent. In our study, point count distance from pipelines was positively correlated with distance from roads (r = 0.72); otherwise, correlations between nearest distances to different footprint types were weak (r < 0.5). However, correlations among nearest distances to different footprint types vary spatially, so that potential edge effects of one footprint could vary with the proximity of other footprint types. In addition, the collection of survey locations we used did not have equal numbers at varying distances in all vegetation types. Undoubtedly, the way a species responds to distance-to-edge will be dependent on the suitability of the vegetation type for that species. Future analysis should evaluate whether different thresholds are observed in different vegetation types. Understanding the magnitude of edge effects is essential as it has major implications for conservation or management of boreal bird habitats as estimates of suitable habitat and population size change considerably when negative edge effects are found (Laurance and Yensen 1991; Ewers and Didham 2007), in addition to whatever vegetation has been lost due to oil and gas development (Johnson and St.-Laurent 2011).

Piecewise regression models could be used in place of an inverse function to assess if the relationship between bird abundance and amount of or distance from footprint increased or decreased or changed sharply at a threshold determined from the data as opposed to a predetermined threshold (Toms and Villard 2015). We initially used the “segmented” package in R to run piecewise regression models and estimate distance thresholds at which species abundance changed sharply (Muggeo 2008). In these initial analyses, zone-of-impact models including piecewise regression functions of distance usually performed better AIC-wise than dose-response models and identified distance thresholds greater than 150 m for many songbirds, suggesting that edge effects were important. However, different thresholds were calculated often when piecewise regression models were rerun, making comparisons to other models less reproducible and increasing model prediction uncertainty, with population projections for some species varying by an order of magnitude. For that reason, we have not included piecewise regression results in this paper.

Using both dose-response and zone-of-impact methods, we found more negative effects of oil and gas footprint on species associated with older coniferous upland forests. The dose-response models demonstrate that habitat loss is occurring for these species, while the zone of impact models show thresholds at distances farther than the sampling radius of a point count, suggesting edge effects. These effects seem to be stronger for species that prefer older coniferous and mixed-wood stands than species that rely on older deciduous forest. This difference could be caused by faster vegetation regrowth on abandoned disturbances in deciduous forest (Nijland et al. 2015). Importantly, old coniferous and mixed-wood forest are some of the least common habitats in our study region and at greatest risk of loss/ disturbance from other extensive and intensive industrial activities such as industrial forestry that targets both deciduous and mixed-wood forests for pulp and paper, but also conifer forests for saw timber. Forestry remains the largest human activity in the boreal forest in terms of area disturbed (Wells et al. 2020) but is usually a temporary disturbance. While timber and pulpwood harvest size and distribution can differ strongly from natural disturbance like forest fires, some harvest plans attempt to emulate natural disturbance patterns in long-term effects on wildlife communities (Schieck and Song 2006, Huggard et al. 2014). Nevertheless, it is likely that effects of forestry could compound energy sector effects on forest birds, especially if both industries are operating in the preferred stands and forest ages used by a species.

Future studies that combine dose-response and zone-of-impact models to evaluate cumulative effects of footprint on boreal birds or other boreal forest wildlife can be improved in a few ways. First, as more point count data becomes available for different footprint types, such studies can consider a wider variety of footprints or whether different footprint types can be combined in larger categories (Mahon et al. 2019). Industry regularly asks to understand how each type of disturbance influences species so that it can develop appropriate mitigation but whether this level of detail is needed for cumulative effects management remains uncertain. Second, while we used a single spatial scale (within 150 m of point counts) and the same data to directly compare explanatory power and regional population estimates of dose-response and zone-of-impact models, future studies could consider effects of footprint amount at larger spatial or multiple spatial scales with co-occurring footprints including harvest (e.g., Bayne et al. [2005]; Mahon et al. [2019]). One difference between dose-response and zone of impact methods is that dose-response can be modelled at any spatial scale with any number of survey stations per landscape unit. In contrast, zone of impact models must use point-level data that are independent of scale. Third, since boreal forest landscapes usually include multiple footprint types in proximity (Johnson and St.-Laurent 2011), future studies can consider both additive and interactive effects of different footprint types as in Mahon et al. (2019). For zone-of-impact models, interactive effects might be assessed by determining if the zone-of-impact for a specific footprint varies among landscapes with different amounts of other footprint types, with the suitability of adjacent vegetation for a particular species, with the age or recovery of vegetation in that footprint, given that probability of species using footprint can change over time with vegetation recovery (Lankau et al. 2013, Leston et al. 2018). Yet another alternative would be to combine dose-response and zone-of-impact predictors within the same models.

The differences in approach we discuss are not fundamentally right or wrong when it comes to scientific investigation. However, they demonstrate a key challenge when using wildlife-habitat models for regulatory decision-making and cumulative effects assessment, as well as the complexity of addressing the issue of cumulative effects in a way that identifies the relative impact of different disturbance types. While such studies are useful and needed, they should not be viewed as the last step before management action. Ultimately, what is needed for effective cumulative effects management is a framework of analysis that is standardized and agreed to by all stakeholders. Such an approach allows integration of new data and studies into this framework to improve model predictions by filling key gaps and uncertainties, rather than evaluating how new analytical approaches influence interpretation. To do so requires fulsome scientific evaluation of how model structure, non-linear relationships, and interactions between different footprints influence results so that future data collection can be more effectively targeted.

Acknowledgments

This work was funded under the Oil Sands Monitoring – Terrestrial Biological Monitoring Program (grant ID: AE 22GRRSD18 Bayne). It is independent of any position of the OSM Program. This work is also a contribution of the Boreal Avian Modelling (BAM) Project, an international research collaboration on the ecology, management, and conservation of boreal birds. We acknowledge the BAM Project’s members, data partners, and funding agencies (including the OSM Program, Environment and Climate Change Canada and the U.S. Fish & Wildlife Service), listed in full at <www.borealbirds.ca/index.php/acknowledgements>.

Author Contributions

LL processed and analyzed data, organized results presented in this manuscript, and contributed to writing. EB procured funds, provided data, proposed the analytical design, and contributed to writing. JT and CLM provided data, proposed the analytical design, and contributed to writing. AC, PS, JB, SJS, FKAS, DS, and TD contributed to writing.

Competing Interests

There are no competing interests that would have influenced the results of analysis in this manuscript. Co-authors have no financial stakes in the outcomes from this study and results from this study are independent from the position of the funding agency.

Data

Data used in modelling and scripts used in this analysis are stored at: https://github.com/LionelLeston/Cumulative-effects-models-Alberta. Model outputs can be recreated from this data. Additional data for predicting regional populations is available from the Boreal Avian Modelling Project at http://borealbirds.ca.

Alberta Biodiversity Monitoring Institute (2017). Alberta wall-to-wall vegetation layer including “backfilled” vegetation in human footprints (Version 6). Alberta Biodiversity Monitoring Institute, Alberta, Canada. https://www.abmi.ca/home/data-analytics/da-top/da-product-overview/Data-Archive/Detailed-Vegetation-Maps.html.
Alberta Biodiversity Monitoring Institute (2018). Human footprint inventory 2016 (Version 1). Alberta Biodiversity Monitoring Institute, Geospatial Centre, Alberta, Canada. https://www.abmi.ca/home/data-analytics/da-top/da-product-overview/Human-Footprint-Products/HF-inventory.html.
Ball JR, Bayne EM, Machtans CM (2009) Energy sector edge effects on songbird nest fate and nest productivity in the boreal forest of western Canada: A preliminary analysis. Proceedings of Tundra to Tropics: Connecting Birds, Habitats and People. Proceedings of the 4th international Partners in Flight Conference.
Bart J, Earnst S (2002) Double sampling to estimate density and population trends in birds. Auk 119:36–45. https://doi.org/10.1093/auk/119.1.36.
Barton K. (2009). MuMIn: multi-model inference. R package version 1. 0. 0. http://r-forge.r-project.org/projects/mumin/. Accessed January 1, 2019.
Bayne EM, Van Wilgenburg SL, Boutin S, Hobson, KA (2005) Modeling and field-testing of Ovenbird (Seiurus aurocapillus) responses to boreal forest dissection by energy sector development at multiple spatial scales. Landsc Ecol 20:203–216. https://doi.org/10.1007/s10980-004-2265-9.
Bayne EM, Dale BC (2011). Effects of energy development on songbirds. In: Naugle DE (ed) Energy development and wildlife conservation in western North America, Island Press, Washington, DC, pp95–114. https://doi.org/10.5822/978-1-61091-022-4_6.
Bayne E, Leston L, Mahon CL, Sólymos P, Machtans C, Lankau H, Ball JR, Van Wilgenburg SL, Cumming SG, Fontaine T, Schmiegelow FKA, Song SJ (2016). Boreal bird abundance estimates within different energy sector disturbances vary with point count radius. Condor 118:376–390. https://doi.org/10.1650/CONDOR-15-126.1
Bender LC, Roloff GJ, Haufer JB (1996) Evaluating confidence intervals for habitat suitability models. Wildl Soc Bull 24:347–352. https://www.jstor.org/stable/3783134.
Boreal Avian Modelling Project. 2015. http://www.borealbirds.ca. Accessed March 31, 2016.
Boxall PC, Chan WH, McMillan ML (2005). The impact of oil and natural gas facilities on rural residential property values: a spatial hedonic analysis. Resour Energ Econ 27:248–269. https://doi.org/10.1016/j.reseneeco.2004.11.003.
Burnham KP, Anderson DR (2002) Model selection and multi-model inference: a practical information-theoretic approach. 2nd edition. Springer-Verlag, New York, NY.
Burgman M, Carr A, Godden L, Gregory R, McBride M, Flander L, et al (2011). Redefining expertise and improving ecological judgment. Cons Lett 4:81–87. https://doi.org/10.1111/j.1755-263X.2011.00165.x.
Cumming SG, Lefevre K, Bayne E, Fontaine T, Schmiegelow FKA, Song SJ (2010) Toward conservation of Canada’s boreal forest avifauna: design and application of ecological models at continental extents. Avian Cons Ecol 5:8. http://dx.doi.org/10.5751/ACE-00406-050208.
Darling A, Leston L, Bayne E (2019) Small mammal abundance differs between pipelines, edges, and interior boreal forest habitat. Can J Zool 97:880–894. https://doi.org/10.1139/cjz-2018-0314.
Environmental Systems Research Institute (ESRI). 2002. ArcGIS. Release 8.3. ESRI, Redlands, California, USA.
Ewers RM, Didham RK (2007) The effect of fragment shape and species’ sensitivity to habitat edges on animal population size. Cons Biol 21:926–936. https://doi.org/10.1111/j.1523-1739.2007.00720.x.
Farwell LS, Wood PB, Sheehan J, George GA (2016). Shale gas development effects on the songbird community in a central Appalachian forest. Biol Cons 201:78–91. https://doi.org/10.1016/j.biocon.2016.06.019.
Fischer J, Lindenmayer DB (2007). Landscape modification and habitat fragmentation: a synthesis. Global Ecol Biogeog 16:265–80. https://doi.org/10.1111/j.1466-8238.2007.00287.x.
Hebblewhite M (2011). Effects of energy development on ungulates. In: Naugle DE (ed) Energy development and wildlife conservation in western North America, Island Press, Washington, DC, pp71–94. https://doi.org/10.5822/978-1-61091-022-4_6.
Hobson KA, Bayne E (2000). Breeding bird communities in boreal forest of western Canada: consequences of “unmixing” the mixedwoods. Condor 102:759–769. https://doi.org/10.1093/condor/102.4.759.
Hsu CS, Robinson PR (2019) Production for recovery. In: Hsu CS, Robinson PR (eds) Petroleum science and technology, Springer, Cham, pp97–125. https://doi.org/10.1007/978-3-030-16275-7_7
Huggard DJ, Grover BE, Dzus E, Smith M, Schieck J (2014). Effectiveness monitoring for biodiversity: comparing 15 year old structural retention harvest areas to fires in boreal aspen. Can J For Res 45:153–161. https://doi.org/10.1139/cjfr-2014-0091
Johnson, CJ, St-Laurent MH (2011). Unifying framework for understanding impacts of human developments on wildlife. In: Naugle DE (ed) Energy development and wildlife conservation in Western North America, Island Press, Washington, DC, pp.27–54. https://doi.org/10.5822/978-1-61091-022-4_6.
Jones NF, Pejchar L, Kiesecker JM (2015). The energy footprint: how oil, natural gas, and wind energy affect land for biodiversity and the flow of ecosystem services. BioSci 65:290–301. https://doi.org/10.1093/biosci/biu224.
Jordaan SM, Keith DW, Stelfox B (2009). Quantifying land use of oil sands production: a life cycle perspective. Environ Res Lett 4:024004. https://iopscience.iop.org/article/10.1088/1748-9326/4/2/024004.
Karr JR, Chu EW (1997) Biological monitoring: essential foundation for ecological risk assessment. Hum Ecol Risk Assess 3:993–1004. https://doi.org/10.1080/10807039709383742.
Lankau H, Bayne E, Machtans C (2013). Ovenbird (Seiurus aurocapilla) territory placement near seismic lines is influenced by forest regeneration and conspecific density. Avian Cons Ecol 8. http://dx.doi.org/10.5751/ACE-00596-080105.
Laurance WF, Yensen E (1991) Predicting the impacts of edge effects in fragmented habitats. Biol Cons 55:77–92. https://doi.org/10.1016/0006-3207(91)90006-U.
Leston L, Bayne E, Schmiegelow FKA (2018) Long-term changes in boreal forest bird occupancy within regenerating harvest units. For Ecol Manag 421:40–53. https://doi.org/10.1016/j.foreco.2018.02.029.
MacDonald L (2000) Evaluating and managing cumulative effects: processes and constraints. Environ Manag 26:299–315.
Mahon CL, Holloway G, Sólymos P, Cumming SG, Bayne EM, Schmiegelow FK, Song SJ (2016) Community structure and niche characteristics of upland and lowland western boreal birds at multiple spatial scales. For Ecol Manag 361:99–116. https://doi.org/10.1016/j.foreco.2015.11.007.
Mahon CL, Holloway GL, Bayne EM, Toms JD (2019). Additive and interactive cumulative effects on boreal landbirds: winners and losers in a multi-stressor landscape. Ecol Appl 29:1–18. https://doi.org/10.1002/eap.1895.
Mahon CL, Pelech S (2021). Guidance for analytical methods to cumulative effects assessment for terrestrial species. Environ Rev 29:201–224. https://doi.org/10.1139/er-2020-0037.
Manel S, Buckton ST, Ormerod SJ (2000) Testing large-scale hypotheses using surveys: the effects of land use on the habitats, invertebrates and birds of Himalayan rivers. J Appl Ecol 37:756–770. https://doi.org/10.1046/j.1365-2664.2000.00537.x.
Martensen, A. C., Ribeiro, M. C., Banks-Leite, C., Prado, P. I., & Metzger, J. P. (2012). Associations of forest cover, fragment area, and connectivity with neotropical understory bird species richness and abundance. Conservation Biology, 26(6), 1100–1111. https://doi.org/10.1111/j.1523-1739.2012.01940.x.
Matsuoka SM, Bayne EM, Sólymos P, Fontaine PC, Cumming SG, Schmiegelow FK, Song SJ (2012) Using binomial distance-sampling models to estimate the effective detection radius of point-count surveys across boreal Canada. Auk 129:268–282. https://doi.org/10.1525/auk.2012.11190.
Muggeo VM (2008) Segmented: an R package to fit regression models with broken-line relationships. R news 8:20–25. https://www.r-project.org/doc/Rnews/Rnews_2008-1.pdf. Accessed January 1, 2019.
Murcia, C. (1995). Edge effects in fragmented forests: implications for conservation. Trends in Ecology & Evolution, 10(2), 58–62. https://doi.org/10.1016/S0169-5347(00)88977-6.
Nelitz MA, Beardmore B, Machtans CS, Hall AW, Wedeles C (2015) Addressing complexity and uncertainty: conceptual models and expert judgments applied to migratory birds in the oil sands of Canada. Ecol Soc 20:4. http://dx.doi.org/10.5751/ES-07906-200404.
Nijland W, Coops NC, Macdonald SE, Nielsen SE, Bater CW, White B et al (2015) Remote sensing proxies of productivity and moisture predict forest stand type and recovery rate following experimental harvest. For Ecol Manag 357:239–247. https://doi.org/10.1016/j.foreco.2015.08.027.
Pickell, P. D., Andison, D. W., Coops, N. C., Gergel, S. E., & Marshall, P. L. (2015). The spatial patterns of anthropogenic disturbance in the western Canadian boreal forest following oil and gas development. Canadian Journal of Forest Research, 45(6), 732–743.
R Core Team (2020). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. https://www.R-project.org/.
Roloff GJ, Kernohan BJ (1999). Evaluating reliability of habitat suitability index models. Wildl Soc Bull 27:973–985. https://www.jstor.org/stable/3783655.
Schieck J, Song SJ (2006) Changes in bird communities throughout succession following fire and harvest in boreal forests of western North America: literature review and meta-analyses. Can J For Resear 36:12991318. https://doi.org/10.1139/x06-017.
Sólymos P, Matsuoka SM, Bayne EM, Lele SE, Fontaine P, Cumming SG, Stralberg D, Schmiegelow FKA, Song SJ (2013) Calibrating indices of avian density from non-standardized survey data: making the most of a messy situation. Meth Ecol Evol 4:1047–1058. https://doi.org/10.1111/2041-210X.12106.
Sólymos P, Mahon CL, Bayne EM (2014) Development of predictive models for migratory landbirds and estimation of cumulative effects of human development in the oil sands areas of Alberta. Technical report, Joint Oil Sands Monitoring: Cause-Effects Assessment of Oil Sands Activity on Migratory Landbirds, Edmonton, AB. pp.829.
Sólymos P, Toms JD, Matsuoka SM, Cumming SG, Barker NKS, Thogmartin WE, Stralberg D, Crosby AD, Dénes FV, Haché S, Mahon CL, Schmiegelow FKA, Bayne EM (2020). Lessons learned from comparing spatially explicit models and the Partners in Flight approach to estimate population sizes of boreal birds in Alberta, Canada. Condor 122:1–22. https://doi.org/10.1093/condor/duaa007.
Toms J, Villard MA (2015). Threshold detection: matching statistical methodology to ecological questions and conservation planning objectives. Avian Cons Ecol 10. http://dx.doi.org/10.5751/ACE-00715-100102.
Wells JV, Dawson N, Culver N, Reid FA, Morgan Siegers S (2020). The state of conservation in North America’s boreal forest: issues and opportunities. Front For Global Change 3:90. https://doi.org/10.3389/ffgc.2020.00090.
Wilson SJ, Bayne EM (2019). Songbird community response to regeneration of reclaimed wellsites in the boreal forest of Alberta. J Ecoacoust 3:1–11. https://doi.org/10.22261/JEA.I4B2LF.

No competing interests reported.

Appendix1.zip
href="https://drive.google.com/drive/folders/1voNc4qL9HljYgEGhTVIq-jX2W2-254qz?usp=sharing" rel="noopener noreferrer" target="_blank">Appendix 1. Model coefficients and plots of predicted abundance of 48 boreal songbird species in relation to the proportion of land within 150 m occupied by seismic lines, pipelines, well-pads, and roads. These predictions were derived from the best dose-response model functions for each of these footprint types. Note: to model the inverse functional response to footprint seen in these plots (e.g., proportion of pipeline), all footprint proportion values < 0.01 were set to a minimum value of 0.01, explaining the flat relationship relative abundance of some species and footprint proportions from 0-0.01.
Appendix2.zip
href="https://drive.google.com/drive/folders/1Ev7MIiQScSzSfDL4A-BsN8-13ym3-j-D?usp=sharing" rel="noopener noreferrer" target="_blank">Appendix 2. Model coefficients and plots of predicted abundance of 48 boreal songbird species in relation to distance to the nearest seismic line, pipeline, well-pad, and road from each point count. These predictions were derived from the best zone-of-impact model functions for each of these footprint types.
Appendix3.zip
href="https://drive.google.com/drive/folders/1Sks_czrf8xfB3g7gSKjowO5WyOUEvHmL?usp=sharing" rel="noopener noreferrer" target="_blank">Appendix 3. Predicted density distribution maps of 48 boreal songbird species in the Al-Pac Forest Management Area in northeastern Alberta (Fox Sparrow not included), based on the best habitat-only, dose-response, and zone-of-impact models for each species.

Comparing alternative methods of modelling cumulative effects of oil and gas footprint on boreal bird abundance

Status:

Version 1

Abstract

Figures

Introduction

Methods

Study Area

Point Counts

Quantifying Vegetation and Human Footprint

Modelling Each Cumulative Effects Assessment Approach

Predictive Ability and Population Consequences of Each Modelling Approach

Results

Direction and Magnitude of Oil and Gas Footprint Effects

Predictive Ability and Population Consequences

Discussion

Declarations

References

Competing Interests

Supplementary Files

Status:

Version 1