DOI: 10.14466/cefasdatahub.170
title: Biodiversity patterns under a shifting baseline: Sensitive fish species core areas (Aim 1) 2024

description: The dataset contains outputs from Aim 1 in Bluemel et al. 2024 - Biodiversity patterns under a shifting baseline: important areas for sensitive fish species and ecosystem functioning to assist marine spatial planning (see attached), Cefas Project Report for Defra, 40 pp. The dataset includes regional core-areas (persistent areas of high fish population density) and projected spatial changes in the distribution of sensitive fish species (currently and under future environmental change) for nine focal species, including the common skate complex (Dipturus spp.), spurdog (Squalus acanthias), tope (Galeorhinus galeus), spotted ray (Raja montagui), undulate ray (Raja undulata), starry ray (Amblyraja radiata), John Dory (Zeus faber), Atlantic wolffish (Anarhichas lupus) and Atlantic halibut (Hippoglossus hippoglossus). Historic core areas / home ranges: vectors of historic core-areas (persistent areas of high fish population density) and home ranges of the nine focal species derived from standardised trawl data (survey catch-rates expressed in terms of biomass per area swept by the trawl, kg.km2) from long-term fishery-independent monitoring programmes (Lynam & Ribeiro 2022a). Where available, core areas were defined for each decade (1985 – 1994; 1995 – 2004; 2005 – 2014) where 50% of the population’s kernel density was concentrated. Similarly, the species range was defined as the areas containing 95% of the population’s kernel density. NB, there are no historic core areas available for tope. Species distributions: contains gridded rasters of the raw predicted habitat suitability projections from the Bayesian Additive Regression Trees (BART, Chipman et al., 2010) modeling algorithm for the nine fish species for the model training/current day period 2005 – 2014 ('Training') and expected habitat suitability for a species given future environmental conditions in 2050 (climate projection representative concentration pathway (RCP) 4.5 (medium emissions, high mitigation)) ('RCP4.5_2050'). Current core areas: contains vectors of current core-areas of the nine focal species derived from the expected habitat, based on a novel thresholding approach tailored to the projections of each species distribution model (SDM) by outcomes from the spatial analysis of survey-derived biomass. Future core areas under rcp4.5 emissions scenario in 2050: contains vectors of future core-areas of the nine focal species derived from the expected habitat, based on a novel thresholding approach tailored to the projections of each SDM by outcomes from the spatial analysis of survey-derived biomass. NB, no future core area is predicted for starry ray under this emissions scenario. The spatial reference for all data is WGS84/UTM zone 30N (EPSG:32630)

Spatial analysis of core areas derived from observations:

Standardised trawl data (survey catch-rates expressed in terms of biomass per area swept by the trawl, kg.km2) from long-term fishery-independent monitoring programmes (Lynam & Ribeiro 2022a) were used to delineate spatial changes in populations of the focal species. Survey catch rates are generally low for sensitive species in comparison to commercial stocks, so data were grouped by decadal time interval to ensure that analyses of change in the spatial distribution and core areas were representative for rarer species. The three decadal time periods investigated were: 1985 – 1994; 1995 – 2004; 2005 – 2014; with the latter decade matching the period for which data were subset prior to modelling species habitats.

Due to the varying spatial coverage and temporal extent of scientific surveys investigated, we selected for analyses those surveys that offered a combination of the longest time series, largest extent, and highest catches of the focal species (given a minimum of 50 records of occurrence between 2005 – 2014). For species occurring predominantly in the North Sea (i.e. wolfish A. lupus, halibut H. hippoglossus and starry ray A. radiata), we combined the ICES co-ordinated international bottom trawl surveys from quarter one (Q1) with those in Q3 (GNSIntOT1: 1983 - 2020 and GNSIntOT3: 1998 - 2020, respectively). For undulate ray R. undulata, which predominantly occurs in the English Channel, the French otter trawl survey conducted in Q4 (GNSFraOT4: 1998 - 2020) was used. For species widely distributed across the Celtic Seas ecoregion, we combined the Scottish otter trawl surveys conducted in Q1 and Q4 (CSScoOT1: 1985 – 2020 and CSScoOT4: 1997 – 2020, respectively), except for tope G. galeus as there were very few records of tope in the Scottish surveys. The Irish otter trawl survey conducted in Q4 (CSIreOT4: 2003 – 2020) had the highest occurrence of tope across all surveys between 2005 – 2014, but the time series was limited to only a single decade. For consistency, data from selected surveys was cropped to the spatial extent of the SDMs (excluding parts of the Skagerrak and Kattegat from the North Sea surveys).

Changes in the ranges and core areas of species were calculated using the spatial kernel density function, weighted by biomass, using gaussian kernels and default bandwidths in the spatialEco package v2.0-2 (Evans et al. 2023) in R v4.3.2 (R Core Team, 2023). Core areas were defined for each decade where 50% of the population’s kernel density was concentrated. Similarly, the species range was defined as the areas containing 95% of the population’s kernel density. Additionally, biomass-weighted population centroids (i.e. the geometric centres of each population) were also calculated with spatialEco to identify any historic shifts in the species distribution, since the mid-1980s and within the survey extent.

Projecting future core areas:

To project changes in core areas of sensitive fish species under future environmental conditions, we combined the core areas defined from survey biomass with projections from a Species Distribution Model (SDM). These models determine a “climatic envelope” encompassing suitable environmental conditions for a species by examining environmental data in areas where the species is present or absent.

Species occurrences were sourced from fisheries-independent groundfish survey datasets from the Northeast Atlantic (Lynam & Ribeiro 2022a,b). Additional records from the early years (2006 – 2014) of the quarter one southwest beam trawl survey (for ICES area 7.e that were not included in the data products) were sourced from the Cefas Survey System database (FSS code Q1SWBEAM). Additional Spanish groundfish survey data not included in the aggregated dataset was provided by Francisco Velasco at the Spanish Institute of Oceanography. Data were also sourced from the Cefas English and Welsh Observer Programme, which samples catches and discards from English and Welsh registered fishing vessels. For elasmobranchs, the dataset was supplemented with observations from the Ocean Biogeographic Information System (OBIS; Grassle, 2000; OBIS, 2021) and from the Global Biodiversity Information Facility (GBIF, 2021).

The data for each species were plotted and checked visually for inconsistencies, with duplicate records removed (where occurrence records had been included within multiple data sources). Species occurrence data from 2005 to 2014 was used to train the Species Distribution Models (SDM) (matching the 2005 to 2014 time period of the environmental data) and projected onto a 0.125° x 0.125° study grid. Grid cells with distribution data were classified as “presence” sites and the remaining cells “pseudo-absence” sites. Absence records from scientific trawls were not considered true absences due to the likelihood of under-sampling and potential for false negatives for some rarer species.

Models were trained using Bayesian Additive Regression Trees (BART; Chipman et al., 2010). BART is a technique based on an ensemble of classification tree models, each built by sequentially splitting the data in two groups based on the value of the explanatory variables, with the aim to separate presence and absence sites. For further technical details of the BART modelling process, environmental data layers used (including bathymetry, sea surface temperature and salinity, difference between surface and bottom temperature, surface chlorophyll, median grain size, and gravel and sand percentage), and model validation, refer to section 4.8 in Astarloa et al., 2023). Once trained, the BART models were used to project expected habitat suitability for a species given future environmental conditions. The future climate projections used correspond to the “representative concentration pathway” (RCP) 4.5 (medium emissions, high mitigation) projections from the Intergovernmental Panel on Climate Change (IPCC)'s fifth phase of the Coupled Model Intercomparison Project (CMIP5) RCP dataset. They were produced using the POLCOMS-ERSEM coupled model running at a 10 km resolution, underpinned by the global climate model MPI-ESM-LR (Kay et al., 2018; available from Copernicus Climate Data Store, Kay 2020).

Future core areas were defined from the expected habitat, based on a novel thresholding approach tailored to the projections of each SDM by outcomes from the spatial analysis of survey-derived biomass. First, were plotted the biomass-derived core area for the period 2005 – 2014 (hereafter “biomass-CA”) for the selected survey per species overlaid with the SDM computed for the same decade for the same species. In each case we observed that the biomass-CAs overlapped with areas of high projected habitat suitability, but not perfectly so. Biomass-CAs also typically extended beyond modelled areas of (i.e., into areas of lower habitat suitability adjacent to high suitability areas). We sought a threshold value for modelled habitat suitability so that model grid cells with suitability above this value could be described as core for the species in a similar way to the biomass-CA. To maintain consistency with the biomass-CA, we sought to capture ≥ 90% of the cells that occur within the biomass-CA within the SDM-derived core area. This was achieved by determining a threshold value that captured the modelled areas of high habitat suitability with an additional spatial buffer to ensure any adjacent, sub-optimal habitats that were contained within the biomass-CA were also included (buffered until at least 90% of the biomass-CA was incorporated).

To define the species-specific threshold and spatial buffer that encompassed ≥ 90% of the biomass-CA, we first extracted the modelled habitat suitability values at locations with known species occurrence ordered from low to high suitability and calculated the percentiles of suitability values in increments of 0.05. We then used each of these percentiles as model thresholds with a range of buffers (from 50km to 200km in increments of 50) to determine the spatial extent of the resulting core area. We selected the combination that included the highest percentile (i.e., containing the highest possible probability values from the model) with the smallest buffer to ensure ≥ 90% of the biomass-CA was incorporated with minimal spatial smoothing.

To check the consistency of these modelled-derived core areas beyond the selected survey extent, we examined the percentage overlap of biomass-CAs derived from other surveys (for those surveys that contained > 50 records of the species between 2005 – 2014) and when combining data from all surveys that used Grand Ouverture Verticale (GOV) trawls (i.e., not accounting for potential survey and/or vessel effects or spatio-temporal variation) for all quarters between 2005 – 2014.

For each species, the final threshold and buffer were applied to the SDM projections for the 2005 – 2014 model training period and the RCP4.5 projection in 2050 to determine any changes between current and future core areas.