Open access

Integrating disparate survey data in species distribution models demonstrate the need for robust model evaluation

Publication: Canadian Journal of Fisheries and Aquatic Sciences
8 November 2023

Abstract

Marine spatial planning and conservation initiatives benefit from an understanding of species distributions across larger geographic areas than are often sampled by any one survey. Here, we test whether the integration of disparate survey data can improve habitat predictions across a region not well sampled by a single survey using Dungeness crab (Metacarcinus magister) from British Columbia as a case study. We assemble data from dive, trawl, and baited-trap surveys to generate six candidate generalized linear mixed-effect models with spatial random fields. To compare single-survey and integrated models, we evaluate predictive performance with spatially buffered leave-one-out cross-validation and independently with two novel approaches using fisheries catch data. We find improved predictive performance and reduced uncertainty when integrating data from surveys that suffer from small sample size, low detectability, or limited spatial coverage. We demonstrate the importance of robust model evaluation when integrating data and predicting to unsampled locations. In addition, we highlight the need for careful consideration of sampling biases and model assumptions when integrating data to reduce the risk of prediction errors.

1. Introduction

Species distribution models (SDMs), which relate species observations to environmental covariates, can be used to generate predictions in unsampled areas to inform a variety of management and conservation decisions (Guillera-Arroita et al. 2015; Péron et al. 2018). When an SDM is needed for management decisions (e.g., where to locate area-based fisheries closures), it is important to consider the quality of the underlying species data (Fletcher et al. 2019). The utility of an SDM can be degraded when species observations are sourced from sampling protocols that are spatially unbalanced, that imperfectly detect the species being targeted, or that cover a truncated environmental range compared to the species’ full niche space or the management area of interest (Dorazio 2014; Guillera-Arroita et al. 2015; Guillera-Arroita 2017; Kéry et al. 2010; Yates et al. 2018).
In the marine realm, sourcing high-quality data for SDMs across large geographic areas and environmental gradients can be especially challenging. For marine spatial planning, often occurring at broad regional scales, SDM predictions that span the entire planning area are beneficial in reducing knowledge gaps in the planning process. However, SDMs are typically generated from individual surveys with predictions limited to within sampled areas (e.g., Anderson et al. 2019; Carrasquilla-Henao et al. 2018). When faced with an application that requires species predictions to be generated across an area not well covered by a single sampling protocol, two options emerge. One option is to fit an SDM with observations from a single sampling protocol that covers only a subset of the area or environmental range and hope it transfers well to the larger area. This is widely accepted as problematic, as transferring models to areas with novel environmental conditions can lead to inaccurate predictions with high uncertainty, or worse, inaccurate predictions with low uncertainty (Perrin 1904; Randin et al. 2006; Wenger and Olden 2012; Eger et al. 2017; Péron et al. 2018).
A second option is to integrate data from multiple surveys that collectively sample the species of interest across a wide geographic and environmental range. This approach can harness the growing number of available survey datasets to produce more robust SDMs (Miller et al. 2019). The aim of data integration is to combine low-quality data (e.g., spatially unbalanced) with higher quality data (e.g., collected with a standardized sampling design and fixed protocol), to allay biases, improve model predictions, and reduce uncertainty (Dorazio 2014; Pacifici et al. 2017). Several methods for combining data exist, the most appropriate will depend on the types of data available (Fletcher et al. 2019; Zipkin et al. 2021). Model-based integration methods require at least a single source of high-quality data to directly inform the latent state being modelled (Pacifici et al. 2017; Isaac et al. 2020; Rufener et al. 2021). High-quality data are data where sampling effort and detection rates are known or can be assumed to be constant across space and time. In some cases, incorporating spatial random effects, which allow for information sharing among data sources when observations are close in space, may improve integrated model predictions (Pacifici et al. 2017; Pinto et al. 2019; Simmonds et al. 2020). When only lower quality data are available, such as opportunistic observations, data integration methods are not appropriate and the simple pooling together of multiple data sources is a commonly used method (e.g., Pirtle et al. 2019).
The application of data integration methods to marine systems is relatively novel (e.g., Grüss and Thorson 2019; Grüss et al. 2018a; Thompson et al. 2022; Runnebaum et al. 2018; Pinto et al. 2019; Watson et al. 2021). Several methods have been applied, using both fixed and random effects to account for sampling differences (e.g., detectability) across the source data (Grüss and Thorson 2019; Grüss et al. 2018a; Thompson et al. 2022; Runnebaum et al. 2018). Previous methods for combining occurrence data (i.e., data pooling) that assume sampling differences between data sources to be minimal (e.g., Nephin et al. 2020; Chu et al. 2019) should be avoided unless there is strong evidence to support the assumption that detectability does not vary across surveys.
Regardless of the method, careful consideration of the quality and limitations of the individual data sources must be taken when combining data from multiple sources (Zipkin et al. 2021). Integration can mask differences between data sources and lead to prediction errors (Fletcher et al. 2019; Simmonds et al. 2020). Given the variety of approaches and potential for prediction errors, there is an even greater need for thorough and robust model evaluation when combining disparate data. Araújo et al. (2019) and Zurell et al. (2020) suggest model evaluation should address both ecological plausibility and predictive performance—ideally with independent data to assess model generality across space and time.
In this paper, we model the habitat of Dungeness crab (Metacarcinus magister), a data-limited coastal marine species, to evaluate the efficacy of data integration when making predictions to geographic areas larger than the area covered by any one data source. In British Columbia, Dungeness crab are sampled regionally and sporadically with a variety of sampling gears and survey protocols, making them an ideal case study to investigate whether the integration of disparate surveys can improve habitat predictions. The available dive, trawl, and baited-trap research surveys each suffer from either a small sample size of detection or non-detection observations, low detectability, or limited spatial coverage. To compare single-survey models with integrated models, we apply a robust spatially buffered leave-one-out cross validation (LOOCV) routine and introduce novel methods for assessing predictive performance across a larger area, outside of the area covered by Dungeness crab surveys, using fisheries catch records. Using these methods, we test the hypothesis that integrated models are more transferable and have lower uncertainty compared to single-survey models.

2. Methods

To facilitate reproducibility and the application of these methods to other study areas and species, we present the methods under two main sections. First, we outline the general modelling approach including the statistical model, assumptions, validation, prediction, and the novel evaluation procedures. Second, we summarise the species and environmental data sources used in this application and detail how the general methods were applied to the Dungeness crab case study.

2.1. Modelling approach

2.1.1. Statistical model

We fit geostatistical models: generalized linear mixed-effect models (GLMMs) with spatial Gaussian random fields using sdmTMB in R (Anderson et al. 2022). We modelled probability of detection Y at location s using a Bernoulli model with a complementary log–log link (following Fithian et al. 2015):
(1)
(2)
(3)
Where μi,s represents the expected probability of detection at location s under survey-type i, Xs is a vector of predictors at location s, βT represents a vector of corresponding parameter coefficients to be estimated, ωs represents the spatial random effects, and αi is the intercept for each survey-type i ∈ {trawl, …, trap}. The integrated models described in Section 2.2.6 allow each αi to be a free parameter.
We also included year to account for temporal variation in occurrence, using a penalized spline smooth term (Wood 2003, 2017). We limited complexity of the smoother by constraining the basis dimension, setting k to 3 or 4, depending on the model. Constraining the complexity helped to separate space–time effects, ensuring that the year term was not flexible enough to absorb some of the spatial or environmental variation when survey data were irregular in their spatial coverage through time.
The spatial random effects were assumed to be drawn from a multivariate normal distribution with a covariance matrix Σω centered on zero and constrained to follow a Matérn covariance function, which controls the rate that spatial correlation decays with distance. We approximated spatial random effects using the Stochastic Partial Differential Equation (SPDE) approach (Lindgren et al. 2011). To implement the SPDE approach, we constructed a mesh using INLA (Lindgren et al. 2011), with a 5 km minimum allowable distance between mesh vertices (Fig. S1). We used a simplified coastline as a barrier over which spatial correlation decayed five times faster across land than water (Bakka et al. 2019). We selected a decay rate of five times faster as this was the middle value tested by Bakka et al. (2019) who found predictive accuracy to be robust to the choice of decay rate within a range of tested values. The barrier was applied to limit spatial correlation across land masses. We fit these models with sdmTMB (Anderson et al. 2022), which maximizes the models marginal likelihood integrating over the random effects with the Laplace approximation using TMB (Kristensen et al. 2016).
We applied model regularization by setting penalties on model terms using priors. We used Normal(0, 1) priors on all fixed effect terms, except for year and survey-type, and used a penalized complexity prior for the Matérn covariance function (Fuglstad et al. 2019) with P(range < a) = 0.05 and P(σ > b) = 0.05. The range is the distance at which two points are effectively independent (correlation ∼ 0.13), σ is the marginal standard deviation of the random field, a was 1/5th an initially estimated range, and b was 2.5 times an initially estimated standard deviation. We ran initial models with each data source without priors or a temporal term to estimate the probabilistic upper and lower bounds for range and standard deviation. See the supplemental materials for Matérn correlation curves, range estimates, and priors (Fig. S2) and spatial random field maps (Fig. S3).

2.1.2. Assumptions

Given the challenges of modelling with survey data that differs in spatial coverage and detectability, accurate model inference will be sensitive to the following model assumptions: (1) the probability of detecting an individual present within the sample area (e.g., bait plume or transect) is constant for each survey (i.e., does not vary with environmental covariates); (2) sampling locations are independent of the distribution of the species being modelled; (3) no unmeasured influential covariate exists with differing distributions across the survey types as this would confound estimates of survey-type intercepts; and (4) all measured covariates are included in a model in correct functional forms. Under these assumptions, sampling biases can be accounted for with a survey- or gear-type fixed effect, and the probability of detection with each survey and the probability of occurrence are proportional, with the constant of proportionality equal to the survey detectability constant.

2.1.3. Validation

We performed model validation qualitatively by examining randomized quantile residuals (Dunn and Smyth 1996) and assessed ecological plausibility by examining model coefficients, the relative influence of covariates, and the marginal effects of covariates using our prior understanding of the species–environment relationships. We estimated randomized quantile residuals using the maximum likelihood estimate for the fixed effects and Markov chain Monte Carlo sampling for the spatial random effects with Stan (Monnahan and Kristensen 2018; Stan Development Team 2017; Rufener et al. 2021).
We estimated the relative influence of covariates using a permutation method, based on the method implemented in MaxEnt software (Phillips et al. 2006). For each covariate we: (1) randomized the covariate with respect to the observations, (2) fit a model with the randomized covariate and all other non-randomized covariates, and (3) accessed model performance with the area under the receiver operating characteristic curve (AUC) metric. We completed steps (1) through (3) 10 times and returned a mean AUC value from the 10 permutations. We then calculated the ΔAUC, the difference between the non-randomized model AUC and the permuted mean AUC from the randomized models. Finally, for each covariate, we divided the ΔAUC by the sum of ΔAUC values from all covariates to obtain the relative influence. A large ΔAUC indicates that the randomized covariate has a large influence, while a small ΔAUC indicates that the covariate has little influence on the model fit. For spatial random fields, we adjusted the procedure as randomization of sampling location (latitude and longitude) was not appropriate. We calculated the influence of the spatial random field by dropping it from the model, then measuring ΔAUC between the model with random fields and the model without.
We calculated marginal effects of the covariates using the evaluation strip method (Elith et al. 2005) and omitting random field uncertainty. When estimating marginal effects for each covariate, we held all other covariates at their mean values, except latitude and longitude, which were held at median values. We calculated mean and median values for marginal effects from the observed values at the location of the sampling events.

2.1.4. Habitat predictions

To quantify differences in predicted habitat area between single-survey and integrated models, we converted the probability of detection predictions to binary predictions of habitat and unsuitable areas using multiple thresholds. We calculated optimal thresholds for each model, instead of using a default 0.5 threshold, to facilitate model comparison between models derived from data with different prevalences (Freeman and Moisen 2008b). We found threshold values to be sensitive to threshold optimization methods. To minimize the influence of any one threshold optimization method on our results, we opted to use three methods: the maximized sum of sensitivity and specificity, percent correctly classified, and kappa. The kappa statistic measures the agreement between observed and predicted values as the percent correctly classified after accounting for the probability of chance agreements (Freeman and Moisen 2008b). We calculated optimized thresholds using the final model predictions generated from all observations, both train and test data, with the PresenceAbsence R package (Freeman and Moisen 2008a) by comparing detection/non-detection observations to the model fitted values. Then we applied each of the three thresholds to classify model predictions, resulting in three binary layers of habitat and unsuitable areas. Lastly, we calculated the mean habitat area (and standard error) from the three habitat layers.

2.1.5. Prediction uncertainty

We approximated spatial prediction uncertainty by simulating 500 draws from the joint parameter precision matrix (sensu Gelman and Hill 2006). Each draw predicted a value for each grid cell within the study area. To visualize uncertainty, we calculated approximate 68% confidence intervals by adding and subtracting one standard error from the mean in link space and applying the inverse link.

2.1.6. Leave-one-out cross validation

A successful integrated model should be able to describe the observed data from each survey at least as well as the corresponding single-survey model. To judge this, we used LOOCV with a spatial buffer (following Roberts et al. 2017) to assess the out-of-sample predictive performance of the models separately on the data from each survey. Under the assumption that no unmeasured influential covariate with vastly differing distributions across survey-type is present, the model judged to have the highest average predictive performance across all surveys should be preferred.
Initial tests with random and spatially blocked cross validation revealed inflated performance. These metrics were particularly biased toward models without a fixed effect to account for sampling biases since the spatial random fields were able to absorb the unexplained variation associated with the sampling bias. Thus, to ensure predictive performance metrics were not inflated by spatial auto-correlation, we used a buffer of 150 km, roughly two times the estimated spatial range derived from the source data observations, around each leave-one-out observation. We fit a model for each observation in the dataset, holding out a single test observation for evaluation as well as all observations within a 150 km radius of the test observation each time.
We assessed predictive performance with two metrics: AUC (Pearce and Ferrier 2000) and Tjur’s R2(Tjur 2009) discrimination metric. Tjur R2 favours models with good discrimination ability that predict very high probability of detection where a species was detected and very low probability of detection where a species was not detected. This is in contrast to AUC that only considers the rank order of probabilities and not their magnitude. To calculate AUC and Tjur R2 metrics requires more than a single test observation. Therefore, to calculate these metrics we generated predictions for each test observation using the LOOCV models, then we combined the predictions into five spatially blocked folds. We created 150 by 150 km grid cells (blocks) that covered the spatial extent of the observations. Then we randomly assigned each spatial block, and the observations that fall within them, to one of five folds. The random assignment procedure was iterative. We repeated the procedure up to 250 times until a solution was reached that balanced, as much as possible, detection and non-detection observations across the folds. We calculated AUC and Tjur R2 five times, once for each fold. We generated the data partitions for the spatially buffered LOOCV and the five-fold spatial blocking using the blockCV R package (Valavi et al. 2019).

2.1.7. Independent evaluation with fisheries catch

Independent data that covers a larger area, outside of the surveyed area, are essential in evaluating the ability of integrated and non-integrated models to predict to unsampled areas. Unlike the LOOCV routine, evaluation scores from independent data are not confounded by potential unmeasured influential covariates that may differ in distribution across the survey data used to build integrated models, so long as the sampling protocols generating the independent data are consistent across space. Independent data can be hard to come by in the under-sampled marine environment, but commercial catch data can be plentiful and often cover a large spatial area. When taking into account sampling effort (i.e., preferential sampling in areas with higher catches of the target species), fisheries data have produced similar spatial patterns of habitat as those derived from survey data (Pennino et al. 2016) and therefore should be applicable for use in model evaluation.
We used a catch per unit effort (CPUE) index layer derived from commercial catch records to assess predictive ability over a larger spatial extent by deriving hotspot areas from both the CPUE index and model predictions and measuring their spatial overlap. We identified hotspot areas using quantiles to ensure the total area of hotspots would be the same for a given study area, which would facilitate comparisons between models. We calculated percent hotspot overlap twice, following the method used by Hardouin and Hargreaves (2023). We divided the hotspot overlap area by both the total prediction hotspot area and the total CPUE index hotspot area, then used the mean of those two percentages. To account for uncertainty in hotspot identification, we used several quantile cut-offs ranging from 0.65 to 0.95 to derive several hotspot layers. We calculated hotspot overlap at each cut-off level. In using hotspot overlap for model evaluation, we assume that catch increases (possibly non-linearly) with occurrence such that regions of high CPUE and high probability of occurrence should align spatially. This is a logical assumption given that the probability of detecting an individual should generally increase with local population density (Royle and Nichols 2003).
In addition, we evaluated whether areas of predicted habitat were correctly classified based on a commercial fishery footprint area. Using a fishery footprint as a benchmark for which to evaluate predicted habitat would be appropriate when the fishery spans a large extent of the study area and under the assumption that the area within the footprint is on average more likely to be high-quality habitat than the area outside the footprint. We generated a fishery footprint from a smoothed CPUE index layer by delineating all areas with greater than zero CPUE as inside the footprint. We classified habitat (above threshold) predictions as correct when within the fishery footprint and incorrect when outside the footprint for each of the three habitat layers (derived from the three threshold optimization methods). Conversely, we classified unsuitable area (below threshold) predictions as correct when outside the footprint and incorrect when inside the footprint for each of the three habitat layers. We then calculated the mean ratio of correctly to incorrectly classified habitat (and standard error) from the three habitat layers.

2.2. Dungeness crab case study

2.2.1. Prior knowledge and environmental data

In British Columbia, Dungeness crab support a large coastwide commercial and recreational fishery, and are relied on by Indigenous communities for food, social, and ceremonial purposes (Ban et al. 2018; Frid et al. 2016). Ecologically, Dungeness crab are an abundant prey and predator species in the coastal marine ecosystem (Dunham et al. 2011) and may be vulnerable to climate change and competition from invasive species (Berger et al. 2021; Jamieson et al. 1998).
Dungeness crab are wide-ranging, occurring from Alaska to southern California across a range of temperatures and salinities (Dunham et al. 2011). Adult Dungeness crab inhabit sandy and muddy areas from the intertidal zone to at least 230 m depth (Rasmuson 2013; Pauley et al. 1989; Dinnel et al. 1986). From lab experiments and field observations, Dungeness crab are known to prefer higher salinities (∼32 Practical Salinity Unit (PSU)) and become stressed at lower salinities (below 24 PSU) (Curtis and McGaw 2012), although Dungeness crab have been found to occupy low salinity areas (below 10 PSU) (Curtis and McGaw 2008). Conversely, high seasonal temperatures have been found to increase Dungeness crab mortality (Sulkin et al. 1996).
To model Dungeness crab distribution, we relied on environmental variables that were known or hypothesised to have an influence on the distribution of Dungeness crab: bottom depth (meters), slope (degrees), fine- and broad-scale bathymetric position index (BPI), mean bottom temperature (celcius), mean bottom salinity (PSU), bottom salinity range (PSU), and substrate indices representing the amount of sandy and muddy substrate within an area. For details on source data and processing steps see the supplemental materials. The distribution of environmental variables at the sample locations and their pairwise correlations can also be found in the supplemental materials (Fig. S4).

2.2.2. Species data

We sourced Dungeness crab observations from four different research surveys conducted by Fisheries and Oceans Canada and their partners: SCUBA dive multi-species surveys (dive), small-mesh trawl surveys targeting shrimp (trawl), baited-trap surveys targeting Dungeness crab (trap DC), and baited-trap surveys targeting Carcinus maenas, Green crab (trap GC). See Table 1 for a summary of observations by survey and Fig. 1 for the spatial distribution of detection/non-detection observations by survey. For details on source data and processing steps, see the supplemental materials.
Fig. 1.
Fig. 1. Dungeness crab observations from dive, trawl, trap GC, and trap DC surveys. Coastline basemap provided by the Canadian Hydrographic Service. Map projection is NAD83 BC Albers (EPSG:3005). Trap DC: baited-trap surveys targeting Dungeness crab; trap GC: baited-trap surveys targeting Green crab.
Table 1.
Table 1. Summary of Dungeness crab data sources.
Dive, trawl, trap-GC, and trap-DC surveys varied considerably in spatial coverage across the region (Fig. 1) and their distribution across environmental space (Fig. 2). The distribution of depth, salinity, temperature, and broad BPI were most dissimilar across surveys (Fig. 2). Trawl surveys covered a narrow range of covariate space compared to other surveys, with the exception of depth and substrate variables. Conversely, dive and trap-GC surveys were limited to shallow coastal areas and thus had truncated depth distributions compared to trawl and trap-DC surveys. While trap-GC and trap-DC surveys appear to have good spatial overlap at a broad scale (as seen in Fig. 1), there is limited overlap at a finer spatial scale as they cover different depth and broad BPI ranges (Fig. 2). For the temporal distribution of sampling by survey, see the supplemental materials (Fig. S5).
Fig. 2.
Fig. 2. Distributions of environmental variables at the sampling locations. Points indicate the median value. All variables were extracted from modelled environmental data layers with the exception of depth and year. Bottom depths were recorded at the time of sampling. Trap DC: baited-trap surveys targeting Dungeness crab; trap GC: baited-trap surveys targeting Green crab.
At a broad scale, most of the surveys had some spatial overlap (Table 2). Trap-GC surveys had the greatest overlap with other surveys, with 15% of samples co-occurring with both dive and trap-DC surveys within a 10 by 10 km grid. Dive and trawl surveys were the only surveys to have no spatially co-occurring samples at that scale. See the supplemental materials for a map of the spatial co-occurrence of surveys (Fig. S6).
Table 2.
Table 2. Pairwise percent spatial co-occurrence between surveys. Spatial co-occurrence was calculated within a 10 by 10 km grid.
The total number of detection and non-detection observations varied between surveys (Table 1). Most sampling protocols did not include repeat sampling, thus occupancy across the surveys was confounded with detectability. However, given our knowledge of the sampling protocols, we suggest there are likely to be differences in Dungeness crab detectability between dive, trawl, and trap surveys. Both DC and GC baited-trap surveys were presumed to have high detection rates of Dungeness crab given that baited-traps are the gear type used by the commercial industry. On dive surveys, Dungeness crab were likely to be detected when present within a quadrat, but Dungeness crab could have been disturbed by the placement of the transect line prior to sampling and may have been missed within the quadrat if algae was dense or Dungeness crab were buried in the sediment. No repeat sampling was completed, therefore the rate of false absences and how it varied among divers is unknown. On trawl surveys, detectability may have been low, especially when compared to baited-traps. Dungeness crab are speculated to avoid capture by trawl nets when they are buried in the sediment and may also be able to move out of the net’s path (Pers. Comm. Dan Curtis, Fisheries and Oceans Canada). Therefore, we expect the difference between the modelled fixed effects for trap-DC and trap-GC surveys to be smaller in magnitude than the difference between trap, dive, and trawl fixed effects.
Differences in the detection rates between the trap-DC and trap-GC survey data (Table 1) may be due to differences in trap type or differences in the environmental conditions within which they sampled. There were only small differences in sampling protocols between trap-DC and trap-GC survey (e.g., type of trap used, see supplemental materials for details). However, their sampling locations differed. Trap-DC surveys sampled in known Dungeness crab habitat while trap-GC surveys target Green crab habitat. The sampling differences between trap-DC and trap-GC surveys could bias model inference and prediction if the density of Green crab has an effect on the distribution of Dungeness crab (e.g., through predation or competition mechanisms). Unfortunately, the coastwide distribution of Green crab density is currently unknown.
In addition to detectability, differences in sample area among surveys could also account for prevalence differences. Traps measured occupancy within a neighbourhood surrounding the traps that were deployed in a line (along a ground-line or with individual floats). The bait plume area is estimated to be a 160 000 m2 area (Dunham et al. 2011). Trawl nets measured occupancy within an estimated 35 000 m2 area, calculated from the trawl mouth width and length of tows. Divers measured occupancy within quadrats covering a 25–125 m2 area depending on slope along a transect line. These sample areas represent an estimated 200% difference between the surveys.
Selectivity for juvenile Dungeness crab with less than 100 mm notch width (Butler 1960) may have varied between sampling protocols. Trap-GC surveys sampled a much greater proportion of juvenile crabs, compared to trap-DC surveys, likely a result of both the smaller mesh size traps and sampling within shallow, nearshore areas where juvenile Dungeness crab reside (Gunderson et al. 1990). For both trap-DC and trap-GC surveys, we removed any samples that caught only juvenile-sized crabs, retaining those that caught at least one adult. Unfortunately, lengths were not measured on trawl or dive surveys. On trawl surveys, it is unlikely that juvenile Dungeness crabs are caught given the depth of the sampling. From our knowledge of dive surveys, we know both juvenile and adult crabs are recorded. However, we assume adults are more likely to be detected by divers than juveniles given their larger size.

2.2.3. Model covariates

We extracted covariate values from modelled environmental data layers (20 m in the Salish Sea and 100 m resolution elsewhere in the study area) at the sampling locations for each survey, with the exception of bottom depth and year, which were recorded during sampling. Environmental model covariates represented static long-term spatial patterns as the source environmental data layers did not cover all sampling years across all locations. We included environmental covariates in the models to improve predictive performance outside of sampled areas.
To account for the difference in sample area size between surveys, we created line segments between the start and end positions of sampling and calculated the mean across all grid cell values that intersected the line segments. The estimated maximum widths of the sampling events for all gear types were less than the resolution of the coarser resolution environmental data layers (100 m). Thus, accounting for their widths using polygons to represent sampling events would have little effect on the values of model covariates. We centered and scaled environmental variables by subtracting their mean and dividing by their standard deviation to ensure model parameter coefficients were on the same scale and that the Normal(0, 1) priors would operate as appropriate penalties.

2.2.4. Coastwide predictions

We delineated the coastwide study area based on the depth range of detection and non-detection observations from all surveys sampling Dungeness crab, limiting the extent from the intertidal to 300 m depth. We used all observations for model fitting to generate coastwide predictions of the probability of Dungeness crab detection. For coastwide predictions, we mosaicked environmental data layers from the Shelf and Salish Sea and aggregated them to a 500 m resolution using the R pacakge terra (Hijmans 2022). When making coastwide model predictions, we made gear-effect model predictions using the trap gear-effect level and survey-effect model predictions using the trap-DC survey-effect level. We selected the trap and trap-DC effect levels for prediction to facilitate more direct comparison with the single-survey trap-DC model predictions. We generated predictions for all years between 2007 and 2017, and presented the mean coastwide prediction from that time period to align with the temporal range of commercial fishery records.

2.2.5. Model evaluation

Following the methods in the Modelling Approach section, we applied two methods for model evaluation and comparison. We used spatially buffered LOOCV to compare models based on their out-of-sample predictive performance evaluated with the survey data and a fisheries catch index and fishery footprint to compare spatial patterns of Dungeness crab habitat over a larger geographic area. Figure 3A depicts an example of the spatially buffered leave-out-one data partitioning for a single observation.
Fig. 3.
Fig. 3. Leave-one-out spatially buffered cross validation example (A) and spatially blocked cross validation folds (B). Blocks and leave-one-out buffers were specified with a 150 km range, depicted by a circle in A and colours in B. Coastline basemap provided by the Canadian Hydrographic Service. Map projection is NAD83 BC Albers (EPSG:3005).
We calculated AUC and Tjur’s R2 with the LOOCV test predictions from each spatially blocked fold (Fig. 3B). As a result of unbalanced sampling across the region, some folds did not have observations from all surveys; trawl observations were only present in three folds and dive observations in four folds. To estimate the uncertainty around the mean of the predictive performance metrics, we calculated confidence intervals as ±2 standard errors with sample size as the number of folds with at least one observation.
We derived a mean CPUE index layer from coastwide commercial trap catch records sourced from 2007 to 2017. See the supplemental materials for details on data processing. We assumed the sampling protocols generating the CPUE data were relatively consistent across the coast given: (1) the reasons behind fishers choosing their set locations can be reasonably assumed to be constant across the British Columbia coast (i.e., preferential sampling occurs everywhere as all fishers want to maximize their chances of success); (2) the skill levels and gear qualities of the fishers across regions can be assumed to be homogeneous; and (3) fishers catch Dungeness crab along the entire coast with the exception of some recent small-scale spatial closures (i.e., there are few areas where fishers were not active over the last two decades).
While there are factors other than Dungeness crab abundance that may influence the selection of fishing grounds (e.g., proximity to population centres), we assume that on average the area within the footprint is more likely to be high-quality Dungeness crab habitat than the area outside the footprint. We justify this strong assumption based on multiple lines of evidence. First, the Dungeness crab fishery, one of the largest fisheries in British Columbia, is coastwide with license areas covering the entire coastline (DFO 2000). Second, all areas with commercial abundance are considered to be fully exploited (DFO 2000). Third, high fishing effort is observed away from major population centers (e.g., Hecate Strait) and in exposed areas (e.g., West Coast Vancouver Island), evidence that fishers are not just fishing close to ports and in calmer, inside waters. Last, previous research that examined why commercial Dungeness crab fishers were absent from areas in the central Strait of Georgia concluded fishers were avoiding these area because higher seasonal temperatures were leading to high mortality of juvenile Dungeness crab in the area (Sulkin et al. 1996). Thus, it can be assumed that an absence of commercial fishery catch, over the long term, is evidence of low-quality habitat or unsuitable areas.

2.2.6. Candidate models

We compared six Dungeness crab models summarised in Table 3. We generated single-survey models from dive, trawl, trap-GC, and trap-DC survey observations, and two integrated models from all surveys. Gear-effect and survey-effect integrated models accounted for differences between surveys with a fixed effect for sampling gear (dive, trawl, trap) or survey (dive, trawl, trap GC, and trap DC). We included the survey-effect model to assess whether sampling protocol differences between trap-GC and trap-DC surveys were large enough to merit a unique fixed effect for each trap survey. Other marine SDM studies that integrated survey data from a variety of sources accounted for survey differences using fixed or random effects (Grüss and Thorson 2019; Grüss et al. 2018a; Thompson et al. 2022). We opted to represent surveys and sampling gear types with a fixed effect as it allowed us to estimate detectability coefficients for each of the survey and gear types, and a random effect was not warranted as we were not considering any other survey or gear types.
Table 3.
Table 3. Description of candidate models. All models included depth, salinity, salinity range, temperature, muddiness index, sandiness index, slope, broad-scale BPI, and fine-scale BPI fixed effects unless indicated. All models included spatial random effects.
We included the same environmental covariates in all models except for depth, which was excluded from the dive and trap-GC models as those surveys had little variation in depth. For variables thought to have an optimum or quadratic relationship with Dungeness crab occurrence (e.g., bathymetry and salinity), we first included them as penalized smooth terms (Wood 2003, 2017), but replaced them with linear terms to simplify the models as we did not observe nonlinear relationships in the initial response curves. We did not include interaction terms in the models.

3. Results

3.1. Model comparison

The majority of the candidate models were able to correctly characterize the effects of the covariates based on our prior knowledge. Depth and substrate were, as expected, important environmental drivers of Dungeness crab detection and on average the probability of Dungeness crab detection was found to decrease with increasing depth, salinity range, slope, and fine-scale BPI and increase in muddier and sandier areas, as well as elevated areas (high BPI) at a broad 1–10 km scale (Fig. 4). However, the ecological plausibility of the candidate models did vary. The trawl, trap-DC, and survey-effect models had weaker effects (i.e., coefficients closer to zero) for almost all fixed effects except depth, where they showed a stronger effect. In contrast, dive- and gear-effect models tended to show stronger effects of most covariates except depth. Given that Dungeness crab are known to inhabit sandy and muddy areas, trawl, trap-DC, and survey-effect models appeared to be under-estimating the effect of substrate as both their sandiness and muddiness index parameter estimates were close to or overlapping zero (Fig. 4). The effect of the Normal(0,1) priors, that operated as penalties on parameter coefficients, is shown in the supplemental materials (Fig. S7).
Fig. 4.
Fig. 4. Fixed effect model coefficients from each model. Bars represent the 95% confidence intervals around parameter estimates (Wald confidence intervals: ± 1.96 the standard error). Survey panel showing both survey- and gear-type effects. Depth not included in dive and trap-GC models. BPI: bathymetric position index; trap DC: baited-trap surveys targeting Dungeness crab; trap GC: baited-trap surveys targeting Green crab.
Salinity was also an important environmental driver of Dungeness crab detection; however, in this case, only the trap-GC model characterized the expected effect. Most models predicted high probability of Dungeness crab detection in low salinity areas. This mischaracterization could be driven by data gaps. There were few observations in low salinity areas (below 26 PSU) across all data sources (Fig. 2) and no observations below 20 PSU, thus the models are extrapolating in those areas. Additionally, the tight correlation between bottom salinity and bottom temperature that occurs within the sampled area could be obscuring the true salinity effect, instead showing a combination of the confounded salinity and temperature effects.
The survey-effect and gear-effect integrated models differed in their detectability coefficient estimates. The survey-effect model estimated trap-DC surveys to be 24 times more likely to detect a crab compared to trap-GC surveys and six times more likely to detect a crab compared to trawl surveys. In contrast, the gear-effect model assumes both trap surveys have the same detectability and estimated trap surveys to be 1.4 times more likely to detect a crab compared to trawl surveys (Fig. 4).
Parameter uncertainty, estimated from the mean standard error of fixed effect coefficients, was larger for the single-survey models on average, and confidence intervals around model coefficients were more likely to overlap zero in single-survey models (Fig. 4). The trawl model had the highest mean parameter uncertainty (mean SE = 0.84), followed by the trap-DC model (SE = 0.50). The survey-effect and gear-effect integrated models showed a 62% and 48% decrease in parameter uncertainty (SE = 0.19 and SE = 0.26) compared to the trap-DC model, respectively.
The marginal effects of model covariates were largely consistent between models (Fig. S8). The greatest variability among the marginal effect curves was seen between Trawl and trap-GC single-survey models and the integrated models predicting to those surveys. Less variability among marginal effect curves was observed between trap-DC and survey-effect models. However, there were still differences in magnitude: the survey-effect model predicted lower probability of Dungeness crab detection at deeper depths, greater salinity ranges, lower broad-scale BPI, and higher fine-scale BPI compared to the trap-DC model.
The relative influence of the fixed and random effects varied between models (Fig. S9). Latent spatial effects, year, salinity, and survey- or gear-type covariates were the most influential overall. Spatial random effects accounted for more than 37% of the relative influence on average. The trap-GC model had the highest relative influence of latent spatial effects at 86% and the survey-effect had the lowest at 11%. Year was highly influential in both trawl and trap-DC models, but accounted for very little relative influence in all other models. Survey type was the most influential predictor for the survey-effect model, accounting for 71% of the relative influence, followed by depth with 12%. The only other model where depth was a relatively influential predictor was the trap-DC model. Gear type was much less influential (13%) in the gear-effect model than was survey type in the survey-effect model (71%).
Predictive performance, evaluated using spatially buffered LOOCV, varied considerably between candidate models with AUC ranging from 0.30 to 0.90 (Fig. 5). Dive, trawl, and trap-GC single-survey models showed low predictive performance when evaluated on their respective datasets (mean AUC < 0.61 and mean Tjur’s R2 < 0.10). Conversely, the trap-DC single-survey model had high predictive performance (mean AUC =0.9 and mean R2 = 0.35) when evaluated with the trap-DC data.
Fig. 5.
Fig. 5. Mean predictive performance for all candidate models from spatially buffered LOOCV as measured by AUC and Tjur’s R2 statistics. Error bars represent ±2 standard errors calculated across spatially blocked folds. AUC: area under the receiver operating characteristic curve; LOOCV: leave-one-out cross validation; trap DC: baited-trap surveys targeting Dungeness crab; trap GC: baited-trap surveys targeting Green crab.
When comparing each single-survey model to each integrated model using their respective single-survey datasets, the integrated models performed at least as well as the single-survey models based on both metrics, as indicated by their overlapping confidence intervals (Fig. 5), with the exception of the gear-effect model that did not perform as well as the trap-DC model when evaluated with trap-DC data. Based on the LOOCV routine, the survey-effect integrated model outperformed the gear-effect integrated model. However, the LOOCV routine was unable to differentiate between the highest performing models: the trap-DC model and the survey-effect model (mean AUC > 0.89 and mean Tjur’s R2 > 0.34), which had overlapping confidence intervals when evaluated with trap-DC data.
For validation plots comparing theoretical and sample randomized quantile residuals for each model, see the supplemental materials (Fig. S10).

3.2. Coastwide model predictions and fisheries catch

Here, we compare the coastwide predictions of Dungeness crab detection from the best performing single-survey model with the integrated models, using fisheries catch records as means of independent evaluation. Although it is unlikely that the trap-DC model will perform well when predicting the coastwide occurrence of Dungeness crab, given that the model was trained only on trap-DC data, it is included here both to confirm that hypothesis and serve as a control from which to measure any potential improvements in the integrated models.
Independent data evaluation with the fisheries CPUE index yielded different results to the LOOCV evaluation. In this case, gear-effect predictions were the best performing, with the largest mean hotspot overlap (54.7 ± 4.5 % overlap, ± 2SE) calculated across all quantiles from 0.65 to 0.95 (Fig. 6). Survey-effect predictions had the second highest mean hotspot overlap (49.9 ± 4.9), followed by trap-DC predictions (44.3 ± 4.4).
Fig. 6.
Fig. 6. Dungeness crab commercial CPUE index (top left), percent of overlap between CPUE and probability of Dungeness crab detection hotspots (bottom left), and coastwide mean probability of Dungeness crab detection predictions between 2007 and 2017 from trap-DC, gear-effect, and survey-effect models with 90th percentile CPUE hotspots overlaid in light blue (right). The CPUE index represents the mean catch (kg per trap) from 2007 to 2017, smoothed with inverse distance weighting over a 15 km neighbourhood. Percent of overlap between CPUE and probability of Dungeness crab detection hotspots were estimated across a range of quantile cut-offs from 0.65 to 0.95. The dashed line denotes the 90th percentile cut-off. CPUE hotspot polygons less than 100 km2 were removed for visualization. The .gear-effect model predictions were made with the trap gear-effect level and the survey-effect model predictions were made with the trap-DC survey-effect level. Coastline basemap provided by the Canadian Hydrographic Service. Map projection is NAD83 BC Albers (EPSG:3005). CPUE: catch per unit effort; Trap-DC: model derived from baited-trap surveys targeting Dungeness crab.
Qualitatively, large differences between the trap-DC and the integrated models were observed in their coastwide predictions of probability of Dungeness crab detection (Fig. 6). The trap-DC model predicted high probabilities throughout most of the region, with only deeper areas of the shelf and the inlets showing lower, but still relatively high probability of Dungeness crab detection. In comparison, gear-effect and survey-effect models had much patchier predictions, with more areas of the coast with low probabilities. At a broad scale, the gear-effect and survey-effect coastwide predictions were similar and areas predicted to have high probability of Dungeness crab detection aligned well with the 90th percentile CPUE hotspots (Fig. 6 light blue polygons). All three models predicted high probability of Dungeness crab detection around CPUE hotspots near Dogfish Bank, Chatham Sound, Clayoquot Sound, the southern Gulf Islands, and the banks off Vancouver (Fig. 6). However, the survey-effect predictions differed from the gear-effect predictions along the North Central Coast where the survey-effect model predicted higher probability of Dungeness crab detection in rocky coastal areas.
When comparing coastwide predictions of habitat—areas with greater than threshold probability of Dungeness crab detection—the differences between model predictions were even more pronounced (Fig. 7). The trap-DC model predicted the largest area of Dungeness crab habitat within the study region (82 217 ± 11 008 km2, ± 2SE), double the habitat area predicted by the survey-effect model (43 226 ± 12 445 km2), and more than four times the habitat area predicted by the gear-effect model (18 0292 ± 5391 km2). Standard error values around the mean coastwide habitat area represent the variability across habitat predictions from the three optimized thresholds used to convert probability into binary habitat predictions.
Fig. 7.
Fig. 7. Coastwide predictions of Dungeness crab habitat from each model. Habitat or unsuitable areas were derived from mean probability of Dungeness crab detection predictions between 2007 and 2017 by applying a mean threshold, indicated within the brackets. Thresholds were selected that maximized kappa, the percent correctly classified, and the sum of sensitivity and specificity. The gear-effect model predictions were made with the trap gear-effect level and the survey-effect model predictions were made with the trap-DC survey-effect level. Coastline basemap provided by the Canadian Hydrographic Service. Map projection is NAD83 BC Albers (EPSG:3005). Trap-DC: model derived from baited-trap surveys targeting Dungeness crab; Trap-GC: model derived from baited-trap surveys targeting Green crab.
When comparing the distribution of the probability of Dungeness crab detection inside and outside the fishery footprint (the area exploited by the Dungeness crab fishery between years 2007 and 2017), we found similar results to the hotspot overlap analysis (Fig. 8). Both integrated models had roughly three times the area of correctly than incorrectly classified habitat based on the fishery footprint (gear effect = 2.97 ± 0.27 and survey effect 2.70 ± 0.60, ±2SE). In contrast, the trap-DC model predictions had roughly the same area of correctly and incorrectly classified habitat (0.93 ± 0.27).
To further compare the trap-DC, gear-effect, and survey-effect coastwide predictions, we considered their spatial prediction uncertainty. Prediction confidence intervals differed more between trap-DC and integrated models (Fig. 9). For the trap-DC model, even the lower bound of the coastwide prediction showed high probability of detection predictions across most of the area, exhibiting a similar spatial pattern to the higher bound prediction from the survey-effect model. Conversely, the lower bound of the gear-effect model showed lower probability of detection in most areas. While all three models had wide prediction confidence intervals in some areas of the coast, they all also exhibited lower uncertainty (i.e., narrow confidence interval range) in the high probability of Dungeness crab detection areas around Dogfish Bank—the area with the largest commercial CPUE hotspot—and hotspot areas near Chatham Sound, Clayoquot Sound, Vancouver, and the southern Gulf Islands (see Fig. 6 for the locations of these regions).
Fig. 8.
Fig. 8. Dungeness crab fishery footprint within the study area (left-hand side) and the density distributions of the mean probability of Dungeness crab detection between 2007 and 2017 inside and outside the footprint for trap-DC, gear-effect, and survey-effect models (right-hand side). In the density distributions panels, the blue and orange shaded regions represent the density of probability of Dungeness crab detection values from areas inside and outside the fishery footprint shown on the map in the left panel, respectively. The dashed vertical lines show the three optimized threshold values that define Dungeness crab habitat for each model (only two lines are visible when two of the thresholds have the same values). The gear-effect model predictions were made with the trap gear-effect level and the survey-effect model predictions were made with the trap-DC survey-effect level. Coastline basemap provided by the Canadian Hydrographic Service. Map projection is NAD83 BC Albers (EPS G:3005). Trap-DC: model derived from baited-trap surveys targeting Dungeness crab.
Fig. 9.
Fig. 9. Prediction confidence intervals for trap-DC, gear-effect and survey-effect models. Simulated prediction uncertainty (standard error) was subtracted and added to model predictions in link space to generate the low and high prediction bounds, shown on the response scale. Range represents the difference between the low and high confidence intervals. The gear-effect model predictions were made with the trap gear-effect level and the survey-effect model predictions were made with the trap-DC survey-effect level. Coastline basemap provided by the Canadian Hydrographic Service. Map projection is NAD83 BC Albers (EPSG:3005). Trap-DC: model derived from baited-trap surveys targeting Dungeness crab.

4. Discussion

With Dungeness crab as a case study, we demonstrate that integrating data from disparate sources can reduce uncertainty and generate more transferable models. Before integrating disparate survey data, we critically evaluated all potential data sources for fitness of purpose and considered the observation process that generated the data from the various surveys as well as the assumptions required for sound model inference. We produced multiple geostatistical spatial GLMMs using different combinations of survey data, allowing us to compare single-survey and integrated model predictions at a coastwide scale. As hypothesised, integrated models outperformed single-survey models. Single-survey models had low predictive performance or over-predicted Dungeness crab habitat, had higher uncertainty, and less overlap with fisheries CPUE hotspots. Coastwide predictions from single-survey models required greater levels of extrapolation, thus these findings largely align with previous studies that have shown greater levels of extrapolation to lead to less accurate and more uncertain predictions (Perrin 1904; Randin et al. 2006; Wenger and Olden 2012; Eger et al. 2017; Elith et al. 2010).

4.1. Comparing single-survey and integrated models

When comparing models using LOOCV, we were expecting data integration to lead to the greatest increases in predictive performance for data-limited surveys, as was observed by Thompson et al. (2022). While we did observe increased performance metrics with integrated models for the more data-limited dive, trawl, and trap-GC models, those increases were not consistent across all metrics (e.g., AUC and Tjur’s R2) and all models. For the higher-performing, less data-limited trap-DC model, we only observed a performance increase from data integration with the survey-effect model in Tjur’s R2 metric. That increase is in line with Grüss et al. (2018b) who reported that integrating data from multiple monitoring surveys led to improved predictions when compared to models built from single surveys. The observed increase in Tjur’s R2 is an indication that the survey-effect integrated model is able to better differentiate, and show more realistic contrasts, between high-detection and low-detection areas, compared to the trap-DC model.
In addition to the increase in predictive performance, we also observed the expected decrease in uncertainty with integrated models. Narrower confidence intervals around parameter estimates were also reported by Rufener et al. (2021) and Thompson et al. (2022) when integrating data from multiple sources and sampling gear types. To understand what was driving these decreases in uncertainty, we compared parameter confidence intervals from the integrated models to a model where we pooled all survey data without accounting for sampling protocol or survey differences with a fixed effect (results not shown). We found parameter confidence intervals to be narrower across all models that included all survey data, including the pooled model. This suggests that the addition of a fixed effect was not solely responsible for reducing uncertainty in the other model parameters. Instead, reductions are likely a result of the larger sample sizes and greater range of environmental covariate space in the integrated models (Wisz et al. 2008; Simmonds et al. 2020).

4.2. Selecting a preferred model

The two evaluation methods, LOOCV and the comparison with fisheries CPUE patterns, produced different results. The LOOCV routine evaluated how well each model could predict the probability of detecting a crab on a dive, trawl or trap survey, while the comparison with fisheries CPUE evaluated how closely model predictions aligned with the spatial patterns in the CPUE index over a larger geographic area. The LOOCV performance metrics suggested that the trap-DC and survey-effect models were preferred, while the comparison with fisheries CPUE indicated that the gear-effect model was preferred.
The LOOCV statistics indicated that the trap-DC model was highly skilled at predicting the probably of detecting a crab on a trap-DC survey. In contrast, the fisheries footprint evaluation showed a high proportion of trap-DC model habitat predictions fell outside the fishery footprint, an indication that the trap-DC model was over-predicting coastwide Dungeness crab habitat. This result was not surprising given the extremely high prevalence and restricted spatial coverage of the trap-DC survey data. Trap-DC surveys also failed to meet the assumption that sampling locations are independent of the distribution of the species being modelled as they were targeting Dungeness crab habitat, thus predictions were biased toward high probability of detection.
Using the LOOCV statistics to compare the performance of the two integrated models, we found the survey-effect model had higher Tjur’s R2 compared to the gear-effect model, but there was no significant difference in AUC. Tjur’s R2 is more sensitive to differences in observer bias (i.e., detectability) compared to AUC, which only takes into account the rank order of predicted probabilities. Thus, higher Tjur’s R2 scores could be an indication that the small protocol differences between trap-GC and trap-DC surveys (i.e., trap mesh size) affected Dungeness crab detectability, justifying the need for a survey-type fixed effect. Alternatively, the LOOCV statistics could be erroneously favouring the survey-effect model if our assumption of no missing influential covariate with differing distributions across survey types is violated. If there is an influential covariate confounding trap-DC and trap-GC survey types, the LOOCV procedure would be unable to accurately rank the survey-effect and gear-effect integrated models. In addition, with our prior knowledge of the sampling protocols (e.g., similar soak times and gear type) we expected trap-GC and trap-DC estimates of detection to be similar in magnitude. Instead, the survey-effect model estimated a large difference in detectability. Thus, a missing confounding variable would need to exhibit a strong effect on Dungeness crab occurrence probability to explain the large estimated difference in detectability between trap surveys.
We hypothesise that Green crab density could be a highly influential confounder. Green crab density is likely to vary between trap-DC and trap-GC surveyed areas since trap-GC surveys are targeting Green crab habitat. In addition, Green crab may compete with and prey on Dungeness crab. Competition-driven niche separation has been theorized to occur between Green crab and Dungeness crab (Gillespie et al. 2015). In laboratory experiments, Green crab displaced Dungeness crab of similar sizes from desired locations and food resources (McDonald et al. 2001). Both species share prey, such as mussels and oysters, but Green crab may be more capable bivalve predators (Behrens Yamada et al. 2010). Additionally, adult Green crab may prey on juvenile Dungeness crab when they settle in shallow estuary and eelgrass areas (Gillespie et al. 2015). These interactions are likely to have a negative affect on Dungeness crab recruitment and survival (McDonald et al. 2001), and could lead to lower probability of Dungeness crab occurrence in areas with high Green crab density and lower detection of Dungeness crab in baited-traps containing Green crab.
If Green crab density is indeed a strong confounder of trap-DC and trap-GC survey-type, then the gear-effect model would be the preferred model. There are several lines of evidence to support this hypothesis: (1) trap-GC surveys targeted areas of high Green crab density; (2) Green crab are thought to be a predator of Dungeness crab juveniles and a competitor of Dungeness crab adults; (3) the LOOCV statistics would be highly sensitive to Green crab density-effects if present; (4) differences in detectability across the two trap surveys were estimated by the survey-effect model to be large despite the protocol differences being known to be small; (5) less realistic species-environmental relationships were estimated by the survey-effect model (e.g., low importance of substrate); and (6) fisheries-dependent data, believed to satisfy the necessary assumptions to be insensitive to Green crab density effects, favoured the gear-effect model.
Conversely, if trap survey types have large differences in detectability and Green crab density is only a weak confounder, then the survey-effect model is preferred. Alternatively, if Green crab density and detectability differences both have an effect, then an average of their model predictions may outperform either model. Future experiments designed to estimate detectability differences between trap-GC and trap-DC surveys are warranted to fully test the missing confounder hypothesis. Fortunately, for end-users interested in these predictions for decision making or subsequent analyses, both the survey-effect and gear-effect models produce relatively similar predictions of the probability of Dungeness crab detection at a coastwide scale. Thus, our results are not overly sensitive to either scenario. We recommend that end-users make use of model averaged coastwide predictions from survey-effect and gear-effect models and treat areas with larger differences between their predictions as regions of greater uncertainty (Fig. 10).
Fig. 10.
Fig. 10. Mean coastwide probability of Dungeness crab detection predictions calculated from gear-effect and survey-effect model predictions (A) and the difference between the coastwide survey-effect and gear-effect model predictions (B). The gear-effect model predictions were made with the trap gear-effect level and the survey-effect model predictions were made with the trap-DC survey-effect level. Coastline basemap provided by the Canadian Hydrographic Service. Map projection is NAD83 BC Albers (EPSG:3005).

4.3. The benefit of independent evaluation

Although cross validation is a common and effective method for model evaluation (Hijmans 2012; Roberts et al. 2017), there are limitations to its utility when evaluating integrated models with disparate data sources (Zipkin et al. 2021; Isaac et al. 2020). When data used for model evaluation have the same biases as the data used to train the model, predictive performance metrics are likely to be inflated (Roberts et al. 2017). When integrating survey data, each with their own unique biases, selecting a single survey to use for model evaluation can be problematic. Both single-survey and integrated models predicted the probability of detecting a crab with a specific survey protocol. Thus, it is not appropriate to fit a model with observations from one survey and evaluate model predictive performance using observations from another survey, unless the detection rates from each survey are known and can be accounted for.
While the spatially buffered LOOCV routine did allow for model predictive performance to be tested in new spatial areas, the fisheries catch methods facilitated the evaluation of predictive performance across a larger geographic area, outside of the areas sampled with the research surveys. Although evaluating SDMs with independent data is not common practice, many studies have demonstrated its importance (Rooper et al. 2016; Abecasis et al. 2014; Anderson et al. 2016; Elith et al. 2006; Phillips and Dudík 2008; Pinto et al. 2016). Here, fisheries catch data were essential for discriminating between the trap-DC single-survey and integrated model coastwide predictions. In addition, by considering the fisheries footprint as a type of range map, we were able to gather evidence in support of our hypothesis that the trap-DC model would produce poor predictions of Dungeness crab occurrence coastwide. We believe the use of a fisheries footprint area in this manner to be a novel application of the data. Although, expert derived range maps have been similarly applied to validate SDMs in terrestrial ecosystems (Mainali et al. 2020; Pineda and Lobo 2012).
For the Dungeness crab case study, we found the fisheries CPUE hotspots and fisheries footprint to be informative and essential for differentiating between the top performing models. However, fisheries-dependent data may not always be appropriate for model evaluation. For unbiased evaluation, fisheries-dependent data should be generated by similar sampling protocols (e.g., traps and bait) across the study area, and fisheries license areas should span the entire study area and be accessible to all fishers. If fishers are known to avoid productive fishing grounds because of a higher chance of gear loss or difficulty in access, then those biases must be corrected for, or those areas removed, prior to model evaluation.

4.4. Limitations

Working with imperfect, spatially unbalanced data limited our predictions to relative measures of the probability of Dungeness crab detection under a given survey protocol. In addition, we required several assumptions, principally that Dungeness crab detectability within each of the sampling protocols did not vary across environmental gradients. If that assumption does not hold true, accounting for detectability and sample-area differences between surveys using fixed effects would not be sufficient to preserve meaningful inference (Warton et al. 2013). Even so, detectability and scale issues may still affect both inference and prediction (Pacifici et al. 2019; Guillera-Arroita 2017).
Another key limitation in this study was the mismatched spatial and temporal resolution of the environmental covariates and Dungeness crab observations that also varied between data sources. With the exception of bottom depth that was measured at the time of sampling, all other environmental variables were sourced from models. The accuracy, for example, of the substrate index layers at the locations of Dungeness crab sampling is unknown; therefore, it is possible that the substrate layers were not accurately resolving substrate at the sample locations, or that the substrate index layers were not at an appropriate scale to model Dungeness crab distribution.
The environmental data layers were intended to be static, resolving physical environmental properties or long-term oceanographic conditions. Thus, we assume temporal scale mismatch is having a negligible impact on our results. However, future SDMs may still benefit from better temporal alignment between species and environmental data, as well as the propagation of uncertainty from environmental layers into SDM prediction uncertainty (Rocchini et al. 2011). We were also unable to fully account for the sample-area differences between surveys when extracting the model covariates from the environmental layers since their resolutions were much larger than the minimum spatial width of sampling (1 m). However, scale mismatches between species data and environmental covariates is a common and fundamental problem in species distribution modelling and ecology (Guisan and Thuiller 2005; Wiens 1989; Zipkin et al. 2021) and not just an issue when integrating data from multiple sources.
The limitations of our environmental data, stemming both from scale issues and missing variables, were buffered by our use of spatial random fields. Spatial random fields were influential in the majority of models, an indication that important covariates, or covariates at the appropriate scale, were missing from the models. Spatial random fields improved predictions within or nearby to sampled areas, but did not improve predictions in unsampled areas outside the range of spatial auto-correlation. Thus, our predictions in unsampled areas were more uncertain and additional evaluation with fisheries data that spanned a larger geographic area was merited.
We opted to convert the catch counts and weights from trawl and trap surveys to detection/non-detection observations to match the dive observations and buffer against any smaller within-survey protocol differences. Though we made every effort to filter each survey dataset for erroneous records and protocols differences, small year-to-year changes in protocols (e.g., bait composition) existed. Thus, by using binary data, we reduced the number of assumptions required when modelling and the impact of potential protocol changes on model predictions. However, future studies could benefit from preserving the count data, as retaining that additional information may further reduce prediction uncertainty. Our choice of complementary log–log link makes it simple to incorporate count data into the integrated model (Pacifici et al. 2019).
Lastly, model inference was constrained by our limited knowledge of Dungeness crab recruitment timing. The intersection between instances of high Dungeness crab recruitment and sampling events is unknown. However, our inclusion of a year smooth term in the models helped to address this by accounting for year-to-year differences in detection.

5. Conclusion

Integrating multiple data sources can be advantageous when modelling the distributions of marine species, especially when a single source does not comprehensively cover the area of interest. Through a varied and unique approach to model evaluation, we demonstrate the benefits of integrating data when only imperfect data, those suffering from small sample size of detection or non-detection observations, low detectability, or limited spatial coverage, are available. By accounting for possible detection differences between survey types with a fixed effect, we improved predictions of the coastwide distribution of Dungeness crab, while reducing the uncertainty. Accounting for sampling differences between survey types, when the detectability of Dungeness crab in those surveys was unknown, was not straightforward. Careful consideration of the sampling biases inherent in each data source, and assumptions that accompany them, was a critical step.
Although we conclude here that data integration improved predictions, this may not be the case for other species and areas with data sources that have different combinations of sampling biases and detectability challenges. We recommend comparing several candidate models, including single-survey and integrated models, and applying a varied approach to model evaluation, including fisheries-dependent data where appropriate or other data sources that span larger or novel geographic areas. When SDM predictions are used to inform management decisions or marine spatial planning, predictive performance and uncertainty should be assessed across the entire management area, which may require the use of non-traditional data sources such as fisheries-dependent data.

Acknowledgements

We thank the editor and anonymous reviewers for their insightful comments that greatly improved this manuscript. We also thank Katie Gale for their edits and Cole Fields for their contribution to the environmental data processing. Thank you to the Canadian Coast Guard, Fisheries and Oceans Canada shellfish and marine spatial ecology research programs as well as their partners for species data collection. Regional ocean modellers at the University of British Columbia for making the SalishSeaCast data publicly available and Diane Masson and Isaac Fine for providing the shelf model data. Natural Resources Canada, the Canadian Hydrographic Service, and National Oceanic and Atmospheric Administration for their role in collecting, processing, and sharing bathymetry data. Finally, we would like to thank the many software developers that have contributed to the open source software used in this study.

References

Abecasis D., Afonso P., Erzini K. 2014. Combining multispecies home range and distribution models aids assessment of MPA effectiveness. Mar. Ecol. Prog. Ser. 513:155–169.
Anderson O.F., Guinotte J.M., Rowden A.A., Clark M.R., Mormede S., Davies A.J., Bowden D.A. 2016. Field validation of habitat suitability models for vulnerable marine ecosystems in the South Pacific Ocean: implications for the use of broad-scale models in fisheries management. Ocean Coast. Manag. 120:110–126.
Anderson S.C., Keppel E.A., Edwards A.M. 2019. A reproducible data synopsis for over 100 species of British Columbia groundfsh. Can. Sci. Advis. Secr. Res. Doc. 2019/041:321. Available from https://publications.gc.ca/collections/collection_2020/mpo-dfo/fs70-5/Fs70-5-2019-041-eng.pdf.
Anderson S.C., Ward E.J., English P.A., Barnett L.A.K. 2022. sdmTMB: an R package for fast, flexible, and user-friendly generalized linear mixed effects models with spatial and spatiotemporal random fields. bioRxiv.
Araújo M.B., Anderson R.P., Márcia Barbosa A., Beale C.M., Dormann C.F., Early R., et al. 2019. Standards for distribution models in biodiversity assessments. Sci. Adv. 5(1):eaat4858.
Bakka H., Vanhatalo J., Illian J.B., Simpson D., Rue H. 2019. Non-stationary Gaussian models with physical barriers. Spat. Stat. 29:268–288.
Ban N.C., Frid A., Reid M., Edgar B., Shaw D., Siwallace P. 2018. Incorporate Indigenous perspectives for impactful research and effective management. Nat. Ecol. Evol. 2:1680–1683.
Behrens Yamada S., Davidson T.M., Fisher S. 2010. Claw morphology and feeding rates of introduced European Green Crabs (Carcinus maenas L, 1758) and native Dungeness crabs (Cancer magister Dana, 1852). J. Shellfish Res. 29(2):471–477.
Berger H.M., Siedlecki S.A., Matassa C.M., Alin S.R., Kaplan I.C., Hodgson E.E., et al. 2021. Seasonality and life history complexity determine vulnerability of Dungeness crab to multiple climate stressors. AGU Adv. 2:e2021AV000456.
Butler T.H. 1960. Maturity and breeding of the Pacific edible crab, Cancer magister Dana. J. Fish. Res. Board Canada 17:641–646.
Carrasquilla-Henao M., Yamanaka K.L., Haggarty D., Juanes F. 2018. Predicting important rockfish (Sebastes spp.) habitat from large-scale longline surveys for southern British Columbia, Canada. Can. J. Fish. Aquat. Sci. 76(5):682–694.
Chu J.W., Nephin J., Georgian S., Knudby A., Rooper C., Gale K.S. 2019. Modelling the environmental niche space and distributions of cold-water corals and sponges in the Canadian northeast Pacific Ocean. Deep. Res. Part I Oceanogr. Res. Pap. 151:103063.
Curtis D.L., McGaw I.J. 2008. A year in the life of a Dungeness crab: methodology for determining microhabitat conditions experienced by large decapod crustaceans in estuaries. J. Zool. 274:375–385.
Curtis D.L., McGaw I.J. 2012. Salinity and thermal preference of Dungeness crabs in the lab and in the field: effects of food availability and starvation. J. Exp. Mar. Bio. Ecol. 413:113–120.
DFO 2000 Dungeness Crab, Coastal Fisheries Areas B, E, G, H, I, & J.Stock Status Rep.C6-14Available from https://waves-vagues.dfo-mpo.gc.ca/library-bibliotheque/331795.pdf.
Dinnel P.A., Armstrong D.A., McMillan R.O. 1986. Dungeness crab, Cancer magister, distribution, recruitment, growth, and habitat use in Lummi Bay, Washington. University of Washington, School of Fisheries. FRI-UW-8612. Available from http://hdl.handle.net/1773/4062.
Dorazio R.M. 2014. Accounting for imperfect detection and survey bias in statistical analysis of presence-only data. Glob. Ecol. Biogeogr. 23(12):1472–1484.
Dunham J., Phillips A., Morrison J., Jorgensen G. 2011. A manual for Dungeness crab surveys in British Columbia. Can. Tech. Rep. Fish. Aquat. Sci. 2964:78. https://publications.gc.ca/collections/collection_2012/mpo-dfo/Fs97-6-2964-eng.pdf.
Dunn P.K., Smyth G.K. 1996. Randomized quantile residuals. Journal of Computational and Graphical Statistics 5(3): 236–244. http://www.jstor.org/stable/1390802.
Eger A.M., Curtis J.M., Fortin M.J., Côté I.M., Guichard F. 2017. Transferability and scalability of species distribution models: a test with sedentary marine invertebrates. Can. J. Fish. Aquat. Sci. 74(5):766–778.
Elith J., Ferrier S., Huettmann F., Leathwick J. 2005. The evaluation strip: a new and robust method for plotting predicted responses from species distribution models. Ecol. Model. 186(3):280–289.
Elith J., Graham C.H., Anderson R.P., Dudík M., Ferrier S., Guisan A., et al. 2006. Novel methods improve prediction of species’ distributions from occurrence data. Ecography 29:129–151.
Elith J., Kearney M., Phillips S. 2010. The art of modelling range-shifting species. Methods Ecol. Evol. 1(4):330–342.
Fithian W., Elith J., Hastie T., Keith D.A. 2015. Bias correction in species distribution models: pooling survey and collection data for multiple species. Methods Ecol. Evol. 6(4):424–438.
Flemming R. 2022. Pacific Multispecies Small Mesh Bottom Trawl Survey (v2.4) [Dataset]. Fisheries and Oceans Canada. OBIS Canada. Available from http://ipt.iobis.org/obiscanada/resource?r=obis_pacshrimp.
Fletcher R.J., Hefley T.J., Robertson E.P., Zuckerberg B., McCleery R.A., Dorazio R.M. 2019. A practical guide for combining data to model species distributions. Ecology 100(6): e02710.
Freeman E.A., Moisen G. 2008a. PresenceAbsence: an R package for presence absence analysis. J. Stat. Softw. 23(11):1–31.
Freeman E.A., Moisen G.G. 2008b. A comparison of the performance of threshold criteria for binary classification in terms of predicted prevalence and kappa. Ecol. Model. 217(1-2):48–58.
Frid A., McGreer M., Stevenson A. 2016. Rapid recovery of Dungeness crab within spatial fishery closures declared under indigenous law in British Columbia. Glob. Ecol. Conserv. 6:48–57.
Fuglstad G.A., Simpson D., Lindgren F., Rue H. 2019. Constructing priors that penalize the complexity of Gaussian random fields. J. Am. Stat. Assoc. 114(525):445–452.
Gelman A., Hill J. 2006. Data Analysis Using Regression and Multilevel/Hierarchical Models. Cambridge University Press, Cambridge, UK.
Gillespie G.E., Norgard T.C., Anderson E.D., Haggarty D.R., Phillips A.C. 2015. Distribution and biological characteristics of European green crab, Carcinus maenas, in British Columbia, 2006-2013. Can. Tech. Rep. Fish. Aquat. Sci. 3120:88. Available from https://publications.gc.ca/collections/collection_2016/mpo-dfo/Fs97-6-3120-eng.pdf.
Gregr E.J. 2012. BC_EEZ_100m: A 100 m raster of the Canadian Pacific Exclusive Economic Zone [Dataset]. SciTech Environmental Consulting, Vancouver, BC. https://bcmca.ca/data/eco_physical_bathymetry/.
Grüss A., Thorson J.T. 2019. Developing spatio-temporal models using multiple data types for evaluating population trends and habitat usage. ICES J. Mar. Sci. 76(6):1748–1761.
Grüss A., Perryman H.A., Babcock E.A., Sagarese S.R., Thorson J.T., Ainsworth C.H., et al. 2018a. Monitoring programs of the U.S. Gulf of Mexico: inventory, development and use of a large monitoring database to map fish and invertebrate spatial distributions. Rev. Fish Biol. Fish. 28:667–691.
Grüss A., Thorson J.T., Babcock E.A., Tarnecki J.H. 2018b. Producing distribution maps for informing ecosystem-based fisheries management using a comprehensive survey database and spatio-temporal models. ICES J. Mar. Sci. 75:158–177.
Guillera-Arroita G. 2017. Modelling of species distributions, range dynamics and communities under imperfect detection: advances, challenges and opportunities. Ecography 40(2):281–295.
Guillera-Arroita G., Lahoz-Monfort J.J., Elith J., Gordon A., Kujala H., Lentini P.E., et al. 2015. Is my species distribution model fit for purpose? Matching data and models to applications. Glob. Ecol. Biogeogr. 24(3):276–292.
Guisan A., Thuiller W. 2005. Predicting species distribution: offering more than simple habitat models. Ecol. Lett. 8(9):993–1009.
Gunderson D.R., Armstrong D.A., Shi Y.B., McConnaughey R.A. 1990. Patterns of estuarine use by juvenile English sole (Parophrys vetulus) and Dungeness crab (Cancer magister). Estuaries 13(1): 59–71.
Hardouin M.E., Hargreaves A.L. 2023. Mapping nationally and globally at-risk species to identify hotspots for (and gaps in) conservation. Proc. R. Soc. B 290:20222307.
Hijmans R.J. 2012. Cross-validation of species distribution models: removing spatial sorting bias and calibration with a null model. Ecology 93(3): 679–688.
Hijmans R.J. 2022. terra: Spatial Data Analysis. R package version 1.6-7. Available from https://cran.r-project.org/package=terra.
Isaac N.J., Jarzyna M.A., Keil P., Dambly L.I., Boersch-Supan P.H., Browning E., et al. 2020. Data integration for large-scale models of species distributions. Trends Ecol. Evol. 35(1):56–67.
Jamieson G., Grosholz E., Armstrong D., Elner R. 1998. Potential ecological implications from the introduction of the European green crab, Carcinus maenas (Linneaus), to British Columbia, Canada, and Washington, USA. J. Nat. Hist. 32(10-11):1587–1598.
Kéry M., Gardner B., Monnerat C. 2010. Predicting species distributions from checklist data using site-occupancy models. J. Biogeogr. 37(10):1851–1862.
Kristensen K., Nielsen A., Berg C.W., Skaug H., Bell B.M. 2016. TMB: automatic differentiation and Laplace approximation. J. Stat. Softw. 70(5).
Kung R. 2021. Canada west coast topo-bathymetric digital elevation model [dataset]. Natural Resources Canada, Ottawa. Available from https://open.canada.ca/data/en/dataset/e6e11b99-f0cc-44f7-f5eb-3b995fb1637e.
Lindgren F., Rue H., Lindström J. 2011. An explicit link between Gaussian fields and Gaussian Markov random fields: the stochastic partial differential equation approach. J. R. Stat. Soc. Ser. B Stat. Methodol. 73:423–498.
Mainali K., Hefley T., Ries L., Fagan W.F. 2020. Matching expert range maps with species distribution model predictions. Conserv. Biol. 34:1292–1304.
McDonald P.S., Jensen G.C., Armstrong D.A. 2001. The competitive and predatory impacts of the nonindigenous crab Carcinus maenas (L.) on early benthic phase Dungeness crab Cancer magister Dana. J. Exp. Mar. Bio. Ecol. 258:39–54.
Miller D.A.W., Pacifici K., Sanderlin J.S., Reich B.J. 2019. The recent past and promising future for data integration methods to estimate species’ distributions. Methods Ecol. Evol. 10(1):22–37.
Monnahan C.C., Kristensen K. 2018. No-U-turn sampling for fast Bayesian inference in ADMB and TMB: introducing the adnuts and tmbstan R packages. PLoS ONE 13(5):e0197954.
Nephin J., Gregr E., Germain C., Field C., Finney J. 2020. Development of a species distribution modelling framework and its application to twelve species on Canada’s Pacific Coast. Can. Sci. Advis. Secr. Res. Doc. 2020/004:107. Available from https://publications.gc.ca/collections/collection_2020/mpo-dfo/fs70-5/Fs70-5-2020-004-eng.pdf.
Pacifici K., Reich B.J., Miller D.A.W., Gardner B., Stauffer G., Singh S., et al. 2017. Integrating multiple data sources in species distribution modeling: a framework for data fusion. Ecology 98(3): 840–850.
Pacifici K., Reich B.J., Miller D.A.W., Pease B.S. 2019. Resolving misaligned spatial data with integrated species distribution models. Ecology 100(6): e02709.
Pauley G., Armstrong D., Van Citter R., Thomas G. 1989. Species profiles: life histories and environmental requirements of coastal fishes and invertebrates (Pacific Southwest)—Dungeness crab. U.S. Fish Wildlife Service. Biological Report 82(11.121). Available from https://apps.dtic.mil/sti/pdfs/ADA224837.pdf.
Pearce J., Ferrier S. 2000. Evaluating the predictive performance of habitat models developed using logistic regression. Ecol. Model. 133:225–245.
Pennino M.G., Conesa D., López-Quílez A., Muñoz F., Fernández A., Bellido J.M. 2016. Fishery-dependent and -independent data lead to consistent estimations of essential habitats. ICES J. Mar. Sci. 73(9):2302–2310.
Péron C., Authier M., Grémillet D. 2018. Testing the transferability of track-based habitat models for sound marine spatial planning. Divers. Distrib. 24:1772–1787.
Perrin E. 1904. On some dangers of extrapolation. Biometrika 3:99–103.
Phillips S.J., Dudík M. 2008. Modeling of species distributions with Maxent: new extensions and a comprehensive evaluation. Ecography 31:161–175.
Phillips S.J., Anderson R.P., Schapire R.E. 2006. Maximum entropy modeling of species geographic distributions. Ecol. Model. 190(3):231–259.
Pineda E., Lobo J.M. 2012. The performance of range maps and species distribution models representing the geographic variation of species richness at different resolutions. Glob. Ecol. Biogeogr. 21:935–944.
Pinto C., Thorburn J.A., Neat F., Wright P.J., Wright S., Scott B.E., et al. 2016. Using individual tracking data to validate the predictions of species distribution models. Divers. Distrib. 22:682–693.
Pinto C., Travers-Trolet M., Macdonald J.I., Rivot E., Vermard Y. 2019. Combining multiple data sets to unravel the spatiotemporal dynamics of a data-limited fish stock. Can. J. Fish. Aquat. Sci. 76(8):1338–1349.
Pirtle J.L., Shotwell S.K., Zimmermann M., Reid J.A., Golden N. 2019. Habitat suitability models for groundfish in the Gulf of Alaska. Deep. Res. Part II Top. Stud. Oceanogr. 165:303–321.
Randin C.F., Dirnböck T., Dullinger S., Zimmermann N.E., Zappa M., Guisan A. 2006. Are niche-based species distribution models transferable in space? J. Biogeogr. 33(10):1689–1703.
Rasmuson L.K. 2013. The biology, ecology and fishery of the Dungeness crab, Cancer magister. Adv. Mar. Biol. 65:95–148.
Roberts D.R., Bahn V., Ciuti S., Boyce M.S., Elith J., Guillera-Arroita G., et al. 2017. Cross-validation strategies for data with temporal, spatial, hierarchical, or phylogenetic structure. Ecography 40(8):913–929.
Rocchini D., Hortal J., Lengyel S., Lobo J.M., Jiménez-Valverde A., Ricotta C., et al. 2011. Accounting for uncertainty when mapping species distributions: the need for maps of ignorance. Prog. Phys. Geogr. 35:211–226.
Rooper C.N., Sigler M.F., Goddard P., Malecha P., Towler R., Williams K., et al. 2016. Validation and improvement of species distribution models for structure-forming invertebrates in the eastern Bering Sea with an independent survey. Mar. Ecol. Prog. Ser. 551:117–130.
Royle J.A., Nichols J.D. 2003. Estimating abundance from repeated presence–absence. Ecology 84(3): 777–790.
Rufener M.C., Kristensen K., Nielsen J.R., Bastardie F. 2021. Bridging the gap between commercial fisheries and survey data to model the spatiotemporal dynamics of marine species. Ecol. Appl. 31(8):e02453.
Runnebaum J., Guan L., Cao J., O’Brien L., Chen Y. 2018. Habitat suitability modeling based on a spatiotemporal model: an example for cusk in the Gulf of Maine. Can. J. Fish. Aquat. Sci. 75:1784–1797.
Simmonds E.G., Jarvis S.G., Henrys P.A., Isaac N.J., O’Hara R.B. 2020. Is more data always better? A simulation study of benefits and limitations of integrated distribution models. Ecography 43(10):1413–1422.
2017. Stan modeling language users guide and reference manual. Version 2.17.0. Available from https://mc-stan.org.
Sulkin S.D., Mojica E., McKeen G.L. 1996. Elevated summer temperature effects on megalopal and early juvenile development in the Dungeness crab, Cancer magister. Can. J. Fish. Aquat. Sci. 53(9):2076–2079.
Thompson P.L., Anderson S.C., Nephin J., Robb C.K., Proudfoot B., Park A., et al. 2022. Integrating trawl and longline surveys across British Columbia improves groundfish distribution predictions. Can. J. Fish. Aquat. Sci.
Tjur T. 2009. Coefficients of determination in logistic regression models—a new proposal: the coefficient of discrimination. Am. Stat. 63(4):366–372.
Valavi R., Elith J., Lahoz-Monfort J.J., Guillera-Arroita G. 2019. blockCV: An R package for generating spatially or environmentally separated folds for k-fold cross-validation of species distribution models. Methods Ecol. Evol. 10(2):225–232.
Warton D.I., Renner I.W., Ramp D. 2013. Model-based control of observer bias for the analysis of presence-only data in ecology. PLoS ONE 8(11):e79168.
Watson J., Joy R., Tollit D., Thornton S.J., Auger-Méthé M. 2021. Estimating animal utilization distributions from multiple data types: a joint spatiotemporal point process framework. Ann. Appl. Stat. 15(4):1872–1896.
Wenger S.J., Olden J.D. 2012. Assessing transferability of ecological models: an underappreciated aspect of statistical validation. Methods Ecol. Evol. 3:260–267.
Wiens J.A. 1989. Spatial scaling in ecology. Funct. Ecol. 3(4):385.
Wisz M.S., Hijmans R.J., Li J., Peterson A.T., Graham C.H., Guisan A. 2008. Effects of sample size on the performance of species distribution models. Divers. Distrib. 14:763–773.
Wood S.N. 2003. Thin plate regression splines. J. R. Stat. Soc. Ser. B Stat. Methodol. 65(1):95–114.
Wood S.N. 2017. Generalized Additive Models: An Introduction with R.2nd ed. Chapman and Hall/CRC. Boca Raton.
Yates K.L., Bouchet P.J., Caley M.J., Mengersen K., Randin C.F., Parnell S., et al. 2018. Outstanding challenges in the transferability of ecological models. Trends Ecol. Evol. 33(10):790–802.
Zipkin E.F., Zylstra E.R., Wright A.D., Saunders S.P., Finley A.O., Dietze M.C., et al. 2021. Addressing data integration challenges to link ecological processes across scales. Front. Ecol. Environ. 19(1):30–38.
Zurell D., Franklin J., König C., Bouchet P.J., Dormann C.F., Elith J., et al. 2020. A standard protocol for reporting species distribution models. Ecography 43(9):1261–1277.

Supplementary material

Supplementary Material 1 (PDF / 1.38 MB).

Information & Authors

Information

Published In

cover image Canadian Journal of Fisheries and Aquatic Sciences
Canadian Journal of Fisheries and Aquatic Sciences
Volume 80Number 12December 2023
Pages: 1869 - 1889

History

Received: 29 November 2022
Accepted: 13 July 2023
Version of record online: 8 November 2023

Data Availability Statement

Data from the small-mesh trawl surveys are available from Flemming (2022). The trap and dive survey data used in this study are currently only available upon request from data managers in Fisheries and Oceans Canada, Pacific Region. The source data for many of the environmental data layers used in this study are publicly available.
Two bathymetry data sources: Kung (2021); Gregr (2012).
Code used to complete this analysis can be accessed: https://gitlab.com/dfo-msea/dungeness-sdm.

Key Words

  1. species distribution modelling
  2. data integration
  3. Gaussian random fields
  4. Dungeness crab

Authors

Affiliations

Institute of Ocean Sciences, Fisheries and Oceans Canada, Sidney, BC, Canada
Author Contributions: Conceptualization, Data curation, Formal analysis, Methodology, Software, and Writing – original draft.
Institute of Ocean Sciences, Fisheries and Oceans Canada, Sidney, BC, Canada
Department of Zoology, University of British Columbia, Vancouver, BC, Canada
Author Contributions: Conceptualization, Methodology, Software, and Writing – review & editing.
Pacific Biological Station, Fisheries and Oceans Canada, Nanaimo, BC, Canada
Department of Mathematics, Simon Fraser University, Burnaby, BC, Canada
Author Contributions: Methodology, Software, and Writing – review & editing.
Sean Anderson and Christopher Rooper served as Associate Editors at the time of manuscript review and acceptance and did not handle peer review and editorial decisions regarding this manuscript.
Institute of Ocean Sciences, Fisheries and Oceans Canada, Sidney, BC, Canada
Author Contributions: Data curation and Writing – review & editing.
Pacific Biological Station, Fisheries and Oceans Canada, Nanaimo, BC, Canada
Author Contributions: Conceptualization, Methodology, and Writing – review & editing.
Sean Anderson and Christopher Rooper served as Associate Editors at the time of manuscript review and acceptance and did not handle peer review and editorial decisions regarding this manuscript.
Brendan Aulthouse
Pacific Biological Station, Fisheries and Oceans Canada, Nanaimo, BC, Canada
Author Contributions: Data curation and Writing – review & editing.
Pacific Biological Station, Fisheries and Oceans Canada, Nanaimo, BC, Canada
Author Contributions: Methodology and Writing – review & editing.

Author Contributions

Conceptualization: JN, PT, CR
Data curation: JN, AP, BA
Formal analysis: JN
Methodology: JN, PT, SA, CR, JW
Software: JN, PT, SA
Writing – original draft: JN
Writing – review & editing: PT, SA, AP, CR, BA, JW

Competing Interests

The authors declare there are no competing interests.

Funding Information

Fisheries and Oceans Canada: No specific grants or awards
This research was supported by the Government of Canada’s Marine Conservation Targets program.

Metrics & Citations

Metrics

Other Metrics

Citations

Cite As

Export Citations

If you have the appropriate software installed, you can download article citation data to the citation manager of your choice. Simply select your manager software from the list below and click Download.

Cited by

1. Estimating the spatial distribution of the white shark in the Mediterranean Sea via an integrated species distribution model accounting for physical barriers
2. Untangling multi-species fisheries data with species distribution models

View Options

View options

PDF

View PDF

Login options

Check if you access through your login credentials or your institution to get full access on this article.

Subscribe

Click on the button below to subscribe to Canadian Journal of Fisheries and Aquatic Sciences

Purchase options

Purchase this article to get full access to it.

Restore your content access

Enter your email address to restore your content access:

Note: This functionality works only for purchases done as a guest. If you already have an account, log in to access the content to which you are entitled.

Media

Tables

Media

Share Options

Share

Share the article link

Share on social media

Cookies Notification

We use cookies to improve your website experience. To learn about our use of cookies and how you can manage your cookie settings, please see our Cookie Policy.
×