Employing Predictive Spatial Models to Inform Conservation Planning for Seabirds in the Labrador Sea

Fifield, David A.; Hedd, April; Avery-Gomm, Stephanie; Robertson, Gregory J.; Gjerdrum, Carina; McFarlane Tranquilla, Laura

doi:10.3389/fmars.2017.00149

ORIGINAL RESEARCH article

Front. Mar. Sci., 22 May 2017

Sec. Marine Conservation and Sustainability

Volume 4 - 2017 | https://doi.org/10.3389/fmars.2017.00149

This article is part of the Research Topic 4th International Marine Conservation Congress View all 29 articles

Employing Predictive Spatial Models to Inform Conservation Planning for Seabirds in the Labrador Sea

$\r\nDavid A. Fifield*$ David A. Fifield¹^*

April Hedd¹

Stephanie Avery-Gomm¹

Gregory J. Robertson¹

Carina Gjerdrum²

Laura McFarlane Tranquilla¹^†

¹Wildlife Research Division, Science and Technology Branch, Environment and Climate Change Canada, Mount Pearl, NL, Canada
²Canadian Wildlife Service, Environment and Climate Change Canada, Dartmouth, NS, Canada

Seabirds are vulnerable to incidental harm from human activities in the ocean, and knowledge of their seasonal distribution is required to assess risk and effectively inform marine conservation planning. Significant hydrocarbon discoveries and exploration licenses in the Labrador Sea underscore the need for quantitative information on seabird seasonal distribution and abundance, as this region is known to provide important habitat for seabirds year-round. We explore the utility of density surface modeling (DSM) to improve seabird information available for regional conservation and management decision making. We, (1) develop seasonal density surface models for seabirds in the Labrador Sea using data from vessel-based surveys (2006–2014; 13,783 linear km of surveys), (2) present measures of uncertainty in model predictions, (3) discuss how density surface models can inform conservation and management decision making, and 4) explore challenges and potential pitfalls associated with using these modeling procedures. Models predicted large areas of high seabird density in fall over continental shelf waters (max. ~80 birds·km⁻²) driven largely by the southward migration of murres (Uria spp.) and dovekies (Alle alle) from Arctic breeding colonies. The continental shelf break was also highlighted as an important habitat feature, with predictions of high seabird densities particularly during summer (max. ~70 birds·km⁻²). Notable concentrations of seabirds overlapped with several significant hydrocarbon discoveries on the continental shelf and large areas in the vicinity of the southern shelf break, which are in the early stages of exploration. Some, but not all, areas of high seabird density were within current Ecologically and Biologically Significant Area (EBSA) boundaries. Building predictive spatial models required knowledge of Distance Sampling and GAMs, and significant investments of time and computational power—resource needs that are becoming more common in ecological modeling. Visualization of predictions and their uncertainty needed to be considered for appropriate interpretation by end users. Model uncertainty tended to be greater where survey effort was limited or where predictor covariates exceeded the range of those observed. Predictive spatial models proved useful in generating defensible estimates of seabird densities in many areas of interest to the oil and gas industry in the Labrador Sea, and will have continued use in marine risk assessments and spatial planning activities in the region and beyond.

Introduction

Offshore resource development has led to an increased need for understanding marine living resource abundance and distribution, and the effects of these activities on biota (Kark et al., 2015). These effects range from point sources and direct mortality (i.e., ship strikes, bottom draggers) to more generalized effects such as eutrophication, pollution, and ocean acidification (Jones, 1992; Laist et al., 2001; Shahidul Islam and Tanaka, 2004; Orr et al., 2005; Rabalais et al., 2009). However, obtaining the needed biological information is logistically challenging, especially in remote regions such as the Arctic (Arctic Council, 2009). Among marine organisms, seabirds are particularly vulnerable to a number of effects from marine industries, including being especially vulnerable to the effects of oil at sea (Wiese and Robertson, 2004; Votier et al., 2005) and other threats such as mortality in fisheries as bycatch (Žydelis et al., 2013; Hedd et al., 2016). Specific to oil pollution, determining the effect of a hydrocarbon release in remote marine environments is difficult when seabird densities in a particular area and season are not known (Wilhelm et al., 2007; Haney et al., 2014).

One such remote region is the Labrador Sea. The Labrador Sea contains significant oil and gas reserves and has been a focus of resource exploration for decades (AMAP, 2010). However, the demand for better baseline biological data to support regional scale environmental assessments (C-NLOPB, 2008) is more recent. The Labrador Sea is important to marine birds year-round supporting breeding seabird colonies during summer and providing important staging, migration, and wintering habitat for seabirds from colonies in both the northern and southern hemisphere (Brown, 1986; Huettmann and Diamond, 2000; Bakken and Mehlum, 2005). Recent tracking information shows movements of North Atlantic breeding seabirds to staging and wintering areas in eastern Canadian and Greenlandic waters that include the Labrador Sea, highlighting the international importance of these waters (Frederiksen et al., 2012, 2016; Mosbech et al., 2012; Fort et al., 2013; Jessopp et al., 2013; Linnebjerg et al., 2013; McFarlane Tranquilla et al., 2013). However, explicit local-scale spatial information on seabird distribution in the Labrador Sea has been limited by patchy marine survey coverage (Fifield et al., 2009), mostly due to logistical difficulties. Understanding the spatial and temporal extent of overlap of seabird populations with offshore resource activities will be critical to the environmental assessment process (Camphuysen et al., 2004; Fifield et al., 2009), and to understanding potential risks to seabirds.

Previous assessments of seabird distributions at sea were limited to the time and place where survey data were actually collected. Recent statistical advances in density surface modeling (DSM) provide an emerging framework that is capable of overcoming many of the constraints of earlier approaches while accounting for imperfect detection (Miller et al., 2013). Now, at-sea survey data are used directly to examine relationships between environmental covariates (generally collected by remote sensing) and seabird densities. Once these relationships are understood, predictions of seabird densities in areas of comparable environmental conditions are possible—even when there is no survey coverage in the immediate area (Miller et al., 2013).

The purpose of this paper is to explore the utility of DSM to meet the demand for better baseline biological data to inform conservation and management decision-making in the Labrador Sea. Specifically, we (1) develop seasonal DSMs for all seabird species combined in the Labrador Sea, (2) present measures of uncertainty in model predictions, (3) discuss how DSMs can aid in conservation and management decision making in the region, and (4) explore some of the challenges and potential pitfalls associated with using these complex modeling procedures.

Materials and Methods

Study Area

The study area (464,104 km²) is aligned with the Labrador Shelf Strategic Environmental Assessment (SEA) Area (C-NLOPB, 2008), defined using North Atlantic Fisheries Organization (NAFO) regions (2G, 2H, 2J) within the Canadian Exclusive Economic Zone (Figure 1).

FIGURE 1

Figure 1. Labrador Sea study area, which encompasses the area of the Labrador Sea delineated by NAFO regions 2G, 2H, 2J, out to the Canadian EEZ (Exclusive Economic Zone). Significant discovery licenses and active C-NLOPB land issuance sectors in the Labrador Sea are also indicated.

Ship-Based Surveys

Surveys were conducted within the purview of Canadian Wildlife Service's Eastern Canada Seabirds At Sea (ECSAS) program, which follows a standardized protocol that incorporates Distance Sampling methods (Buckland et al., 2001; Fifield et al., 2009; Gjerdrum et al., 2012). Briefly, surveys consisted of a series of sequential observation periods (nominally 5-min duration) along a continuous line-transect, providing segments of mean length 1.64 km (SD: 0.61 km) as the sample unit for the DSM (see below). GPS coordinates were recorded at the beginning and end of each segment. Environmental variables that might be expected to affect bird detection such as wind speed and wave height were collected and updated at the beginning of each segment (Gjerdrum et al., 2012). We recorded observations of birds on the water continuously along each segment, and used a snapshot approach for flying birds (Tasker et al., 1984; Gjerdrum et al., 2012). Distance categories (0–50, 5–100, 100–200, 200–300 m; Gjerdrum et al., 2012) for each observation were assigned by estimating the perpendicular distance to each individual bird (or the centroid of each group of birds) with the help of a pre-marked custom ruler constructed for each observer-vessel combination (Gjerdrum et al., 2012). Data were entered directly into the ECSAS Microsoft Access database using voice recognition software or recorded on datasheets and entered into the database (Fifield et al., 2009).

Data Analysis

Data Extraction and Filtering

Survey data consistent with the methods described in Section Ship-Based Surveys were extracted from ECSAS database version 3.38. Only data collected from moving vessels whose speed exceeded 4 knots were included (Gjerdrum et al., 2012). Certain seabird species are attracted to fishing vessels, which can artificially inflate densities of these birds around such vessels. Therefore, observations of these taxa (i.e., gulls, Northern fulmar, shearwaters, and black-legged kittiwake) collected during Fisheries and Oceans Canada (DFO) scientific trawl surveys were removed prior to analysis (10 of 35 trips, see Supplementary Tables 1, 2).

Environmental covariates used in modeling

We modeled the density (or equivalently, abundance) of seabirds in the Labrador Sea relative to environmental variables (Table 1). These were chosen because each has been demonstrated to correlate with the distribution or abundance of seabirds at sea, or could be expected to correlate (Louzao et al., 2006, 2011; Wakefield et al., 2009; Oppel et al., 2012). To extract (or interpolate) values for all of the dynamic and most of the static variables associated with the starting position of each segment, we used Marine Geospatial Ecology Tools (MGET; Roberts et al., 2010) within ArcGIS (Table 1). Values of sea surface temperature (SST; JPL MUR MEaSUREs Project, 2010), and anomalies in both sea surface height (SSH), and eddy kinetic energy (EKE; AVISO, 2016) were extracted directly using MGET Data Products tools. To estimate spatial gradients in SST (SSTG), monthly rasters of SST climatology were first created, and gradients subsequently estimated as a proportional change (PC) within a surrounding 3 × 3 cell moving window following Louzao et al. (2006): PC = [(SST maximum value − SST minimum value) × 100]/(SST maximum value). SST minimum and maximum rasters were created from monthly climatologies using the Focal Statistics tool, with subsequent algebra executed using the Raster Calculator. Static environmental variables including bathymetry and its derivatives were determined from ETOPO2 grids (National Geophysical Data Center, 2006) which is based on satellite altimetry and shipborne ground-truthing measurements. The spatial gradient in bathymetry (BathyG) was estimated as indicated above for the gradient in SST, and distance to the continental shelf break (Dist1000, defined as the 1,000 m contour) was calculated using the Near tool within ArcGIS Analysis toolset. To avoid collinearity among these covariates, we assessed both seasonal and overall variance inflation factors to ensure that they were all below 3 (Zuur et al., 2010), and to ensure that the absolute value of all Pearson correlation coefficients between pairs of explanatory variables were <0.52 which is below the 0.7 threshold where collinearity has been shown to severely effect model predictions (Dormann et al., 2013).

TABLE 1

Table 1. List of environmental variables used to model the density of seabirds within the Labrador Sea.

Environmental data used to predict seabird abundance or density

Each cell within a 2 × 2 km prediction grid covering the study area was populated with the static variables listed in Table 1, along with seasonally-averaged values for the dynamic covariates recorded over the 10-year period between 2005 and 2014. Seasons were defined as Summer (June–August), Fall (September–October), Winter (November–March), and Spring (April–May). To generate these dynamic covariate layers, we first used MGET Data Products tools to create monthly rasters for the period January 2005 to December 2014, and then the ArcGIS Raster calculator tool (Spatial Analyst toolset) to create seasonal averages. Seasonal spatial gradients in SST were calculated as indicated above for the modeling step, but using the seasonally averaged SST rasters as input. Values for the prediction grid were obtained using an MGET Spatial and Temporal Analysis tool (“Project Raster to Template”), which projected each environmental raster to the coordinate system, cell size and extent of the prediction grid (Lambert Azimuthal Equal Area projection).

Modeling Approach

We conducted an analysis for all seabird species combined (Fifield et al., 2016) producing a seasonal spatial predictive model in R 3.2.3 (R Core Team, 2016) following a two-stage DSM approach (Hedley and Buckland, 2004; Miller et al., 2013) using the Distance v. 0.9.4 (Miller, 2015), dsm v. 2.2.13 (Miller et al., 2015), and mgcv v. 1.8–17 R packages.

Stage 1

In Distance Sampling, the process of fitting a detection function accounts for the fact that some birds are unavoidably missed during surveys (Buckland et al., 2001). A detection function was fitted by modeling detectability as a function of distance from the observer and other covariates including observer identity, species (or guild), flock size, wind speed, wave height, season, and bird behavior (i.e., flying or on the water). The best detection function model was selected using Akaike Information Criteria (AIC), and model fit was assessed using plots of detection function fit to the distance histogram, and through χ² goodness-of-fit tests. The best detection function was then used to correct the estimated density of seabirds present in each survey segment (Hedley and Buckland, 2004).

Stage 2

We constructed seasonal Generalized Additive Models (GAMs, Wood, 2006) of density as a function of environmental covariates in surveyed areas. GAMs were fitted to the density data using quassiPoisson, negative binomial and Tweedie response distributions.

For each distribution, an initial full model was fitted containing seasonal thin-plate regression spline smooths (Wood, 2006). Dynamic variables (EKE, SSH, SST, SSTG; Table 1) were fitted with a single smooth to allow the natural variation within and across seasons to inform the fit. We allowed for a seasonal time component for static variables (Bathy, BathyG, Dist1000), which have no natural seasonal variation, in GAM fitting (through use of the mgcv “by = Season” option). A parametric (i.e., non-smooth) term for each season was also initially included in all models. QQ-plots of initial full models using each of the three response distributions indicated the negative binomial distribution provided a superior fit and the other distributions were not considered further.

Model refinement progressed by iterative backwards selection of the full negative binomial model refining all smooth terms first, followed by parametric terms. At each iteration, the remaining smooth term with the lowest estimated degrees of freedom (EDF; or sum of the four individual component EDFs for seasonal smooths) less than 1.5 was considered for removal, or replacement by a linear parametric term. The term was replaced by a parametric one if the slope of the coefficient was significantly different from 0 (i.e., if its p ≤ 0.01) or removed otherwise. Seasonal smooth terms were replaced by parametric terms with an interaction with season to allow the slope of the linear relationship to vary seasonally. Refinement of parametric terms consisted of removing those whose p ≥ 0.01, and interactions were removed before main effects.

Autocorrelation was assessed via a plot of correlations in the residuals at lags of 1–20 consecutive segments, by inspection of variograms of the residuals at scales of 2, 5, 10, and 100 km, and visual inspection of spatial plots of residuals (bubble plots from R package sp v. 1.2-4).

We used the final fitted model to predict and map seasonal density in a 2 × 2 km grid-cell matrix over the entire study area using the same environmental covariates. The spatial resolution of this grid was chosen to match the segment length. Mapping of model uncertainty in individual grid cells, as measured by the coefficient of variation (CV), was computed at a lesser resolution of 6 × 6 km due to computer memory constraints.

Seasonal Densities and Population Estimates

Two seasonal density (birds·km⁻²) and population estimates were produced: one for the entire study area and a second excluding areas of environmental extrapolation beyond the range of sampled covariates used to build the models. Extrapolations outside these ranges are largely speculative and such areas have been indicated in Figure 2 with crosshatches to appropriately caution the reader (Mannocci et al., 2016). These areas varied seasonally. Further, large portions of the study area may be covered by ice in winter, spring, and summer and therefore be unavailable for seabirds. We removed these areas from density and population estimates by obtaining ArcGIS shapefiles containing weekly climatologies of median ice concentration for the East Coast and Northern Canadian Waters for the period 1981–2010, as compiled by the Canadian Ice Service (http://ec.gc.ca/glaces-ice/). To be conservative, we were interested in the maximum seasonal median ice extent where ice concentrations exceeded 9−9+ (i.e., waters that were >90% ice covered), and so these concentrations were selected and the ArcGIS Merge tool (Data Management Tools) was used to create a single polygon for each season.

FIGURE 2

Figure 2. Seasonal predicted densities of all seabirds based on Generalized Additive Models (GAMs). Areas of extrapolation beyond the range of model covariates are indicated with crosshatching. Seasons are defined as Summer (June–August), Fall (September–October), Winter (November–March), and Spring (April–May).

Sensitivity Analysis

We conducted a sensitivity analysis by individually removing each covariate from the final fitted DSM. For each reduced model, a new set of seasonal maps and population estimates were produced and compared with the results from the full model.

Results

From May 2006 to November 2014, observers surveyed 13,783 linear km, or 4,135 km² (8,392 segments) for seabirds within the Labrador Sea study area, with the most intensive effort occurring between 2012 and 2014. Surveys were conducted during all seasons; however, as survey platforms were largely ships of opportunity, effort was distributed unevenly across the study area by season (Figure 3). The most complete spatial coverage was achieved in summer, when surveys were conducted from Saglek Bank in the north to Hamilton Bank in the south, both on and off the Labrador Shelf, while the poorest spatial coverage occurred in winter when most surveys were concentrated in the south in the vicinity of Hamilton Bank (Figure 3). Effort was most intense in fall, with shelf waters south of and including the Nain Bank being well-surveyed. Seasonal differences in the spatial distribution and intensity of survey effort have implications for the precision of seabird density estimates, as demonstrated below.

FIGURE 3

Figure 3. Seasonal bird observations and coefficient of variation (CV) for predicted densities of all seabirds based on Generalized Additive Models (GAMs). Seasons were defined as Summer (June-August), Fall (September-October), Winter (November-March), and Spring (April-May).

In total, ship surveys detected 33,469 seabirds (12,379 flocks) in 4,638 of 8,392 (55%) segments, with dovekie (Alle alle), northern fulmar (Fulmarus glacialis), black-legged kittiwake (Rissa tridactyla), and murres (Uria spp.) being the taxa most frequently observed (see Supplementary Tables 3, 4 for a complete list of observed species).

Modeling

The best detection function model contained a hazard-rate key function and covariates for taxon, wind, wave height, bird behavior (flying or on water), season, and flock size. Average detection probability was 38% (CV = 0.26), indicating that failure to adopt a distance sampling methodology would have underestimated seabird densities in each surveyed segment by an average of 2.6 times. For the density surface model, a negative binomial distribution with dispersion parameter k = 0.248 provided the best fit to the data explaining 17.6% of the deviance. Retained smooth terms included Bathymetry, BathyG, Dist1000, EKE, SST, and SSH, along with the parametric Season term (see Supplementary Figure 1 model output and plots of regression splines).

Residual autocorrelation at lags of 1–20 segments was relatively low (all r < 0.2), and variograms and spatial plots of residuals showed no obvious patterns.

Predicted Densities and Population Estimates

The Labrador Shelf and adjacent portions of the Labrador Sea are clearly important regions for seabirds, particularly during fall and winter, when densities were 33.9 (CV = 0.26) birds·km⁻² and 21.6 (CV = 0.37) birds·km⁻², respectively (Table 2, Figure 2). During fall and winter, the estimated number of seabirds in the region (to the nearest 100,000) were 15.4 million (95% CI: 9.3–25.5 million) and 5.3 million (95% CI: 2.6–10.9 million), respectively. During fall, relatively high densities were predicted throughout the Labrador Shelf, from the Saglek and Nain Banks south to the Labrador Trough. Predicted densities in a substantial portion of the Saglek Bank and Labrador Trough exceeded 50–75 birds·km⁻² in fall. Areas with high predicted densities in fall overlapped with several significant hydrocarbon discoveries on the Labrador Shelf (Figures 1, 2). Somewhat lower bird densities were predicted in offshore regions in fall; however, the lack of survey coverage in these areas lowers confidence in the predictions. During winter, the model predicted high densities (>25 birds·km⁻²) all along the continental shelf break. In summer, mean density was 15.3 (0.25) birds·km⁻² and particularly high densities were also predicted along the shelf break, with a noticeable peak in the very north. The regional seabird population estimate in summer was 6.4 million (95% CI: 3.9–10.5 million). Predicted densities and numbers of birds were lower overall during spring, and averaged 9.0 (CV = 0.26) birds·km⁻² and 2.7 million (95% CI: 1.6–4.4 million) birds (Table 2). Mean (range) CVs (6 × 6 km cells) were 0.19 (0.13–1.14) in spring, 0.20 (0.14–1.50) in summer, 0.19 (0.11–0.52) in fall, and 0.36 (0.16–0.85) in winter (Figure 3).

TABLE 2

Table 2. Seasonal densities and population estimates (to the nearest 100,000) excluding ice-covered areas for the entire study area and for the area without extrapolation beyond range of model covariates (see Figure 2).

Sensitivity Analysis

Density and abundance estimates were insensitive to removal of individual model terms; confidence intervals for these quantities all overlapped extensively for all models (Supplementary Figure 2). The spatial density pattern was also similar across all models, with the greatest differences evident in fall for models without Dist1000 and Season (Supplementary Figure 3). Likewise, spatial uncertainty patterns were similar across all models except during fall for the model without Dist1000 where uncertainty off the continental shelf (where there was no sampling) was greater than in other models (Supplementary Figure 4).

Modeling and Interpretation Challenges

Areas of extrapolation where the values of environmental predictors were beyond the range of those covariates in the model occurred to varying degrees across seasons (crosshatched in Figure 2). In the majority of cases, this was caused by a lack of vessel-survey coverage in those water depths during those seasons. Predicted densities in these areas must be interpreted with extreme caution. This is of particular interest to the current application since much of the continental slope and deep offshore area, particularly in the southeast, is in the early stages of hydrocarbon exploration. The lack of survey coverage in these regions highlights a priority data gap.

DSM fitting, spatial prediction and variance computation was both CPU and memory intensive. For example, this process typically required 1–2 h for a single run to complete on a Lenovo laptop with an Intel Core i5-3320M processor running at 2.6 GHz with 8 GB ram utilizing all four processors in parallel. During the per-cell prediction and variance calculation steps the R process memory requirements typically grew to more than 4 GB at a 6 × 6 km scale. Attempts to compute the variance at a 2 × 2 km scale failed when R attempted to allocate more than 96 GB of memory.

Discussion

The purpose of this study was to explore the utility of DSM to inform conservation and management needs of seabirds in the Labrador Sea. This is of particular importance given domestic and international commitments to protect 10% of coastal and ocean areas by 2020 (CBD, 2010) and the regional increase in industrial activity. Using at-sea seabird survey data collected onboard ships of opportunity, seasonal models of seabird distribution and abundance were built throughout the region. These will contribute to the baseline information needed for quantitative environmental assessments which precede development of the offshore hydrocarbon industry. While data collection within oil and gas significant discovery license areas themselves was sometimes sparse, by utilizing this relatively recent statistical approach (Miller et al., 2013), models built from data collected in nearby regions enabled plausible predictions of seabird densities throughout much of the area of interest. We expect these outputs will inform future conservation and management decision making in the Labrador Sea and contribute to conservation targets. In addition to marine protected area planning and offshore oil and gas exploration and production, understanding risks from shipping (and associated risk of hydrocarbon releases) and mortality from fishing activities can all be improved by these models.

Our results indicate millions of seabirds use the Labrador Sea at all times of year, with especially high concentrations over the shelf regions in the fall. This is consistent with known migration patterns of seabirds in the North Atlantic. Many of the 3 million breeding thick-billed murres from Canadian colonies (Gaston et al., 2012) migrate along the Labrador Shelf to winter in low arctic regions (McFarlane Tranquilla et al., 2013) and large numbers (~850,000) of thick-billed murres from colonies spanning the eastern and western North Atlantic winter in the Labrador Sea (Frederiksen et al., 2016). Similarly, black-legged kittiwakes from colonies spanning the North Atlantic use the Labrador Sea in winter, especially waters off the shelf edge (Frederiksen et al., 2012). Dovekie are the most numerous seabird in Northwest Atlantic waters, and millions of breeding birds from West Greenland colonies migrate south through the Labrador Sea to wintering grounds off Newfoundland (Fort et al., 2013).

In terms of conservation planning, our models have the ability to inform and refine the boundaries of Ecologically and Biologically Significant Areas (EBSAs) recently defined for the Newfoundland and Labrador Shelves Bioregion (DFO, 2013; Figure 4). The most recent EBSA delineation process for the Labrador Shelf included seabird data collected from pelagic surveys (DFO, 2013), but this data was very coarse and lacked spatial and temporal coverage across much of the Labrador Sea (Fifield et al., 2009). Our analysis shows the importance of the fall migration period for seabirds, especially along the northern parts of the shelf. Portions of this region are covered by the Northern Labrador and the Hopedale Saddle EBSAs (DFO, 2013) but important areas of ocean for seabirds between these two EBSAs are not included. A portion of the important summer concentration of seabirds along the shelf break at the northern part of the study area is included in the Outer Shelf Saglek Bank EBSA, and concentrations of wintering seabirds are within the Labrador Slope EBSA, which covers a large stretch of ocean along the shelf edge in the middle of the study area (DFO, 2013). As Canada progresses with Marine Protected Area planning exercises to meet Convention on Biological Diversity (CBD) targets (Secretariat of the Convention on Biological Diversity the Scientific Technical Advisory Panel—GEF, 2012), our results can be used to further refine EBSA and subsequent Marine Protected Area boundaries. Since these models have a temporal component, boundaries can be created that are dynamic, for example they may change with seasons, or can be changed if future data indicates a significant shift in marine bird concentrations.

FIGURE 4

Figure 4. Ecologically and Biologically Significant Areas (EBSA) within the study area (DFO, 2013).

Our analysis considered all seabirds in the Labrador Sea combined. However, different seabird guilds have different relative risks to specific types of threats in the marine environment (e.g., Wiese and Robertson, 2004; Good et al., 2009; O'Hara et al., 2009; Hedd et al., 2016). Predictive modeling at the guild or species level is important to inform risk and damage assessments for seabirds (Le Rest et al., 2016; Fox et al., 2017). Future studies should address this gap by constructing DSMs on a per guild (or species) basis to better inform management and conservation for marine birds in the Labrador Sea.

A number of static (Bathymetry, BathyG, Dist1000) and dynamic (EKE, SST, SSTG, SSH) oceanographic features we selected were retained in the model as predictors of seabird densities and explained 18% of deviation in the data. Previously, several of these features had been identified as important predictors of at-sea densities for other top predators (Mannocci et al., 2016). This clearly, and not unexpectedly, indicates that the density of seabirds are not uniform in time or space on the Labrador Sea, and smaller scale estimation will provide better estimates of seabird density at scales relevant to marine resource planners and to development of conservation policy. During fall, high predicted densities occurred in relatively shallow waters throughout the continental shelf. In both the summer and winter, however, highest predicted densities occurred along the shelf break and slope, habitat features which have been shown to be important for seabirds in other regions (Catard et al., 2000; Croxall and Wood, 2002). Specific to oil and gas exploration, significant discovery licenses which are on the continental shelf midway along the Labrador coast, overlap with high predicted densities of seabirds, notably in the fall. Predicted seabird densities in areas of interest to the oil and gas sector in the southeastern part of the Labrador Sea identified large areas where birds concentrated during fall and winter, and to a lesser extent during spring and summer. However, in the fall, surveys were restricted to the continental shelf, precluding us from obtaining high confidence in the density estimates in these areas of interest to the oil and gas industry. The best approach for addressing these remaining data gaps would be a designed survey approach that targets high priority areas, as opposed to opportunistic surveys.

Caution and an understanding of the uncertainty associated with model predictions is vital when interpreting DSM-predicted density output. To facilitate appropriate interpretation we produced companion maps of associated uncertainty (in our case, using the CV or coefficient of variation; see also Winiarski et al., 2014; Mannocci et al., 2016), to be viewed alongside fine-scale maps of predicted seabird densities. These maps are to be interpreted in tandem, to highlight seasons and areas where model predictions are associated with higher degrees of uncertainty. In modeling seabirds in the Labrador Sea, the overall range of per-cell CVs was modest with similar precision achieved during spring, summer and fall, and lesser precision in winter with relative differences being driven by seasonal sample sizes. Within each season, the CV tended to increase toward the edges of the study area where some covariates approached the limit of their ranges and so corresponding predictions were less precise.

Large areas of ocean featured environmental extrapolation beyond the range of model covariates, and predictions in these areas must be interpreted with caution (e.g., in deep off-shelf waters during fall and winter). The tendency for GAMs to produce extreme predictions at or beyond the limits of covariate values is well-known (Mannocci et al., 2016). For example, in a previous modeling exercise (Fifield et al., 2016) we included the X and Y coordinates (projected longitude and latitude, respectively) as covariates in our GAMs. However, much of the study area was outside the X and Y range of sampling in one or more seasons leading to extreme predicted densities and low precision. It is thus apparent that the choice of model covariates strongly affects predicted densities and their precision. To draw the reader's attention to this fact in the current analysis, areas of extrapolation have been crosshatched in the prediction surfaces (cf. Mannocci et al., 2016). It can be important to dig deeper and understand which environmental predictor(s) cause the extrapolation, and to what extent they affect predicted densities in such areas. For example, during fall there is a considerable area of moderately high predicted densities (up to 50 birds·km⁻²) in deep off-shelf waters in the northeast where extrapolation is due to a lack of seabird survey effort in deep waters. However, examining the fall regression spline plot for Bathymetry from the GAM output (Supplementary Figure 1) indicates that water depth had no effect on predicted density during this season due to a nearly flat relationship centered at 0. However, Bathymetry was retained in the model due to a significant interaction effect in other seasons. This highlights the tradeoffs inherent between pooling all data in a single model with seasonal interactions vs. modeling each season separately where data may be sparse in some seasons. Instead, these high predicted densities are a result of the Dist1000 model term. Thus an intimate knowledge of modeled relationships may be required to interpret crosshatched plots of environmental extrapolation (cf. Mannocci et al., 2016) and this may guide the degree of caution required in interpreting predictions in such areas. Nonetheless, a qualitative assessment of prediction plausibility should be undertaken in the context of other sources of distribution and abundance information (see below and Mannocci et al., 2016).

Sensitivity analysis indicated that the model was relatively insensitive to single-term deletions with density and population estimate confidence intervals overlapping for all models. Likewise spatial patterns of density and uncertainty (CV) were similar across all models except for in deep off-shelf waters during fall for 2 models; an area where bathymetric extrapolation occurred due to lack of fall deep-water survey coverage. Improving survey effort in this area would likely help to limit model sensitivity.

Residual spatial autocorrelation can bias analyses by inflating Type I error rates, and a variety of statistical methods have been proposed to account for this fact (Dormann et al., 2007; Bivand et al., 2015). Similar to Winiarski et al. (2013), un-modeled residual autocorrelation was relatively low in this study and we feel confident that it did not unduly bias our results. Nonetheless, future work should investigate methods to accommodate autocorrelation.

When deciding whether to use DSM, considerations include the non-trivial commitments of time and computer resources required to execute these complicated techniques. The time required to fit and visualize prediction and uncertainty outputs should not be underestimated, especially when comparing many candidate models across multiple species groups (Fifield et al., 2016). Computing variance for DSM predictions is memory intensive, particularly for large study areas with small grid-cell sizes and may stretch computing resources beyond available limits. Approaches to streamline the process could involve using multiple computers running in parallel, using a centralized high performance computing server environment, or processing the study area in a series of sub-areas and combining the results (although this further complicates the analysis). The statistical approaches involved in DSM are an active area of research (Hedley and Buckland, 2004; Miller et al., 2013), but all are based on advanced techniques and extensive experience and training are required. These are not readily accessible to all and likely help to explain the gap between having data and using data to best inform critical environmental assessment processes, both in the Labrador sea, and in other marine areas of concern where survey coverage may be limited. These caveats notwithstanding, we believe DSMs are a valuable approach to inform seabird conservation over large areas of ocean where survey coverage is incomplete, and would recommend that management agencies develop capacity and expertise in this or similar predictive modeling approaches.

Author Contributions

DF: developed project design, secured funding, contributed to data collection, lead the data analysis, and co-wrote the manuscript. AH: contributed both to data analysis and to the writing of the manuscript. SA: developed project design, contributed to data collection, analysis and editing of the manuscript. GR: developed project design, secured funding, managed project through completion and contributed to the writing of the manuscript. CG: coordinated data collection and observer training, and managed the ECSAS database. LM: contributed to data collection, project management (field logistics), editing the manuscript.

Funding

Funding for this project was received from the Environmental Studies Research Funds (ESRF; project 2010-08S) and Environment and Climate Change Canada.

Conflict of Interest Statement

The funding base for the Environmental Studies Research Funds (ESRF) is provided through levies on frontier lands paid by interested holders such as oil and gas companies. However, the ESRF is directed by a joint government/industry/public Management Board and administered by a secretariat in the Office of Energy Research and Development, Natural Resources Canada, and remains at arm's length from this project and its findings.

Acknowledgments

We thank the Environmental Studies Research Funds (ESRF), Dave Taylor (ESRF Secretariat) for funding, guidance and support. Access to survey vessels was kindly provided in part by DFO chief scientists and CCGS captains, Joey Angnotok, and Putjotik Fisheries Ltd. We thank the many dedicated observers who spent long hours collecting data in often less than ideal weather conditions, and the captains and crews for their hospitality. We thank the following for many insightful conversations regarding study design, data collection, and analysis and mapping approaches: François Bolduc, Steven Duffy, Jack Lawson, Ken Morgan, Pat O'Hara, Rob Ronconi, Christian Roy, Pierre Ryan, Sabina Wilhelm and Sarah Wong. We thank two reviewers for comments that improved the quality of the paper.

Supplementary Material

The Supplementary Material for this article can be found online at: http://journal.frontiersin.org/article/10.3389/fmars.2017.00149/full#supplementary-material

References

AMAP (2010). Assessment 2007: Oil and Gas Activities in the Arctic – Effects and Potential Effects. Vol. 1. Oslo: Arctic Monitoring and Assessment Programme (AMAP), 423.

AVISO (2016). Ssalto/Duacs Altimeter Products. Copernicus Marine and Environment Monitoring Service (CMEMS). Available online at: http://www.marine.copernicus.eu (Accessed November 05, 2015).

Arctic Council (2009). Arctic Marine Shipping Assessment 2009 Report. Protection of the Marine Environment (PAME) and the Arctic Council.

Bakken, V., and Mehlum, F. (2005). Wintering areas and recovery rates of Brunnich's guillemots (Uria lomvia) ringed in the Svalbard archipelago. Arctic 58, 268–275. doi: 10.14430/arctic428