Aquatic species shows asymmetric distribution range shifts in native and non-native areas

Costa, Eudriano F. S.; Encarnação, João; Teodósio, Maria Alexandra; Morais, Pedro

doi:10.3389/fmars.2023.1158206

ORIGINAL RESEARCH article

Front. Mar. Sci., 02 May 2023

Sec. Global Change and the Future Ocean

Volume 10 - 2023 | https://doi.org/10.3389/fmars.2023.1158206

This article is part of the Research Topic Biological Invasions in Marine and Estuarine Ecosystems View all 9 articles

Aquatic species shows asymmetric distribution range shifts in native and non-native areas

Eudriano F. S. Costa^1,2*

João Encarnação³

Maria Alexandra Teodósio³

Pedro Morais^4,5

¹Institute of Marine Sciences - OKEANOS, University of the Azores, Horta, Portugal
²Department of Oceanography and Fisheries, Institute of Marine Research (IMAR), University of Azores, Horta, Portugal
³Centre of Marine Sciences (CCMAR), University of Algarve, Campus de Gambelas, Faro, Portugal
⁴Sustainable and Smart Aquaculture Collaborative Laboratory, Olhão, Portugal
⁵Institute of Environment, Florida International University, Miami, FL, United States

Introduction: Climate change is reshaping the geographical distribution of species across the globe. In marine ecosystems, climate change creates novel challenges to an environment impacted by numerous anthropogenic stressors. Forecasting shifts in species distribution, including the expansion of non-indigenous species under climate change scenarios, is a management challenge for today’s world.

Methods: We applied Bayesian Additive Regression Tree (BART) models to investigate the environmental factors modulating the occurrence and habitat preferences of the Atlantic blue crab Callinectes sapidus Rathbun, 1896 in native and non-native areas. With BART models, we also aimed to predict its current and future distributions under two climate change scenarios (RCPs 4.5 and 8.5). BART models were performed using global occurrences of the Atlantic blue crab – recorded from 1830 to 2022 – and several environmental covariates (i.e., water temperature, salinity, current velocity, and rugosity). Model accuracy was assessed with the Area Under the Curve (AUC) and True Skill Statistics (TSS) criteria. Cross-validation experiments were made to balance accuracy prediction and uncertainty model prediction intervals.

Results and discussion: AUC and TSS values indicated that data validation was successful for the BART model. Water temperature was the most critical variable affecting the presence probability of the Atlantic blue crab. The BART model predicts asymmetric distribution range shifts on both sides of the Atlantic Ocean. The non-native populations will experience broader shifts in their distribution range than in the native range, and the RCP 8.5 scenario model outputs predict a wider distribution range by the end of the century. Overall, we anticipate significant ecological changes in native and non-native areas. The range expansion in native areas is often equivalent to the ecological shifts induced by invasive species, so lessons learned by ecologists and managers in non-native areas will provide actionable insights to managers in native areas.

1 Introduction

The abundance and distribution of marine species changed quickly throughout the 20^th century due to multiple anthropogenic stressors, including habitat modification and destruction, overexploitation, pollution, and the introduction of non-indigenous species (Halpern et al., 2015; Corrales et al., 2018). In the 21^st century, global climate steers species distribution in terrestrial and aquatic ecosystems has become more evident (Poloczanska et al., 2013; Molinos et al., 2015; Pecl et al., 2016; Molinos, 2020). For some species, their distribution range contracted (e.g., alpine trees (Engler et al., 2011; Charbonnel et al., 2016), others were extirpated (e.g., Bay checkerspot butterfly Euphydryas editha bayensis) (Sternitsky, 1937), southern sea lions Otaria flavescens (Shaw, 1800), mountain wildflower Androsace septentrionalis Linnaeus 1753)) (McLaughlin et al., 2002; Baylis et al., 2015; Panetta et al., 2018), while many others experienced a degradation of their conservation status because they now live under sub-optimal temperature conditions that disrupt fitness, reproduction, survival, or established biotic interactions (e.g., white-beaked dolphin Lagenorhynchus albirostris (Gray, 1846), reef-buiding coral Diploastrea heliopore (Lamarck, 1816), golden kelp Laminaria ochroleuca Bachelot Pylaie 1824, brook trout Salvelinus fontinalis (Mitchill, 1814), and Chinese giant salamander Andrias davidianus (Blanchard, 1871) – Elliott and Simmonds, 2007; Cantin et al., 2010; Warren et al., 2012; Blois et al., 2013; Franco et al., 2018; Zhang et al., 2020). The poleward expansion and tropicalization of temperate regions is another phenomenon indicative of the global changes we are witnessing in terrestrial and aquatic ecosystems, from terrestrial plants to invertebrates, but also reptiles, fish, and mammals in multiple regions across the globe (Zeidberg and Robison, 2007; Azzurro et al., 2011; Encarnação et al., 2019; Pinsky et al., 2020; Osland et al., 2021).

The continuous change of world ecosystems due to global climate change, in tandem with other anthropogenic impacts, reshapes communities while opening the door to non-indigenous species to settle or become invasive (Willis et al., 2010; King et al., 2021; Wang et al., 2021). The factors involved in the expansion of non-indigenous species go beyond habitat suitability since biotic resistance and stochastic events may play a significant part in halting or slowing colonization and invasiveness (Brown and Barney, 2021). The Atlantic blue crab Callinectes sapidus Rathbun, 1896 is native to the western Atlantic Ocean, from the coasts of Massachusetts (United States) to central Argentina (Mancinelli et al., 2021). In Europe, it has become one of the fastest-spreading marine invasive species (Mancinelli et al., 2017; Mancinelli et al., 2021; Fuentes et al., 2019; Encarnação et al., 2021a; Encarnação et al., 2022), the reason why it is on the list of the 100 worst invasive species in the Mediterranean Sea (Streftaris and Zenetos, 2006). The oldest records in Europe date from 1900 in the Bay of Biscay (Atlantic Ocean), 1935 in Greece (Mediterranean Sea), 1951 in the Baltic Sea, and 1967 in the Black Sea and Portugal (Atlantic Ocean) (Nehring, 2011; Encarnação et al., 2022). On the western coast of Portugal, the 15 scattered records reported between 1967 and 2019 suggest that the species had not established a viable population (Encarnação et al., 2022), but this changed in 2016 after the first specimens were found in southern Portugal (Morais et al., 2019; Vasconcelos et al., 2019). Soon after, the Atlantic blue crab rapidly colonized the entire southern coast and adjacent estuarine ecosystems (Encarnação et al., 2021a) and then expanded into the west coast to locations where previous introductions have failed (Encarnação et al., 2022).

In the native range, the Atlantic blue crab stocks are overexploited along the northwest Atlantic coast (Aguilar et al., 2008; Semmler et al., 2021), while biotic interactions with the invasive European green crab Carcinus maenas (Linnaeus, 1758) exert additional pressure, mainly to juveniles (MacDonald et al., 2007). Nevertheless, at higher abundances, the Atlantic blue crab may be able to limit the invasion of the European green crab in North America (DeRivera et al., 2005). In Europe, and particularly in the Gulf of Cadiz, the native European green crab does not seem to represent a threat to the Atlantic blue crab because it has become a successful invader in locations where native crabs are abundant (Rodríguez et al., 1997; Oliveira et al., 2000; Sprung, 2001; Encarnação et al., 2013; Morais et al., 2019; González-Ortegón et al., 2022). The Atlantic blue crab has also become a profitable fishery in local fishing communities in Portugal, Spain, Greece, Turkey, and Tunisia (Nehring, 2011; Mancinelli et al., 2017; Vasconcelos et al., 2019; Encarnação et al., 2021a; Clavero et al., 2022).

Predicting changes in the distribution of a species in native and non-native regions has conservation as a common goal. Predicting where and when a species will expand/contract its distribution, or become more or less abundant, adds valuable information about species under fishery management, like the Atlantic blue crab in the native range. The prediction of species distribution allows for the prioritization of conservation resources and increases fishers’ preparedness for changes in fisheries. Failure to do so may result in the collapse of regional fisheries, as happened with the Atlantic cod Gadus morhua (Linnaeus, 1758) in the Gulf of Maine (North America) (Pershing et al., 2015). The poleward range shifts of other commercially important marine species have been documented for the common snook Centropomus undecimalis (Bloch, 1792) in North America and predicted for others, such as the dolphinfish Coryphaena hippurus Linnaeus, 1758 (Salvadeo et al., 2020; Purtlebaugh et al., 2020). Fishing communities must then adapt to changes in species distribution, either by following the fishery or focusing on new resources in their fishing grounds (Rogers et al., 2019; Paukert et al., 2021; Galappaththi et al., 2022). In non-native areas, predicting the invasion pathway of a species is one of the most valuable tools for improving its management and minimizing/mitigating future impacts (Barbet-Massin et al., 2018; Liu et al., 2020; Lyons et al., 2020). The commercial exploitation of invasive species has been used as a mitigation measure when eradication measures are unfeasible (e.g., lionfish Pterois miles (Bennett, 1828) in the Caribean Sea (Barbour et al., 2011; Morris et al., 2011) and Mediterranean Sea (Ulman et al., 2022), weakfish Cynoscion regalis in estuaries of central Portugal (Cerveira et al., 2021), the Atlantic blue crab in the Mediterranean Sea (Mancinelli et al., 2017).

Thus, our objective was to investigate the environmental factors modulating the occurrence and habitat preferences of the Atlantic blue crab in native and non-native areas using Bayesian Additive Regression Trees (BART model), as well as predict the current and future distribution of the species under two climate change scenarios: RCPs 4.5 and 8.5. We hypothesize that 1) water temperature will be closely related to the distribution of the Atlantic blue crab in native and non-native habitats and 2) rising ocean temperature will reshape the current distribution. To test our hypotheses, we used Species Distribution Modeling (SDM) because it allows for addressing issues in conservation planning (Franklin, 2013; Villero et al., 2017), evolution (Pinsky et al., 2020), biogeography (Torn et al., 2020; Pickens et al., 2021), invasiveness potential of non-native species (Jiménez-Valverde et al., 2011; Mainali et al., 2015), and even inferring the impacts of climate change on species distribution (Lezama-Ochoa et al., 2016; Zhang et al., 2019; Guerra et al., 2021). We chose BARTs because it is a flexible SDM tool offering a wide range of options to extract relationships between species and abiotic conditions, elucidating the roles of predictive variables in determining species distributions, environmental condition preferences (two-dimensional range of environmental conditions where the species showed the highest probabilities of presence), posterior distributions on predictions as a measure of uncertainty, and optimal spatial locations for species expansion or invasion (Chipman et al., 2010; Carlson, 2020; Carlson et al., 2022). With this information, we also aim to contribute to more efficiently managing this commercially valuable species in the native and non-native ranges.

2 Materials and methods

2.1 Species data

The use of Species Distribution Models (SDMs) to predict current and future distribution ranges of marine organisms has increased substantially since the mid-2010s (Robinson et al., 2017; Guerra et al., 2021; Hernández-Ucera et al., 2021; Pickens et al., 2021). SDMs rely on the relationship between a species’ georeferenced record (e.g., presence/absence data) and environmental variables to predict the probability of being found elsewhere (Robinson et al., 2017). So, we downloaded global spatial occurrences data for the Atlantic blue crab from 1830 to 2020 (data made available by Mancinelli et al. (2021) – https://doi.org/10.6084/m9.figshare.12896309). These authors compiled data from multiple databases/sources of published and unpublished records: Global Biodiversity Information System (GBIF, www.gbif.org), Ocean Biogeographic Information System (OBIS, www.obis.org), Biodiversity Information Serving Our Nation (BISON, https://bison.usgs.gov), iNaturalist (www.inaturalist.org). Additional records made between 2019 and March 2022 were retrieved from GBIF, scientific literature (Encarnação et al., 2021a; Chaouti et al., 2022; Encarnação et al., 2022; Encarnação et al., 2023), and unpublished records submitted to our citizen science campaign New Marine Species of the Algarve (NEMA – Novas Espécies Marinhas do Algarve, www.NEMAlgarve.com). The datasets were merged, and then the duplicates, doubtful identifications, and mismatched geographic coordinates (i.e., coordinates outside the marine environment) were removed. In the end, we were left with 27,664 valid records for subsequent analyses (Figure 1).

FIGURE 1

Figure 1 Location of the 2130 Atlantic blue crab Callinectes sapidus Rathbun, 1896 presence points used for Bayesian Additive Regression Tree modeling. Dotted lines represent the maximum depth of 200 m. The red dots represent the records collected in the non-native area, while the green dots represent the records collected in native area.

We used spatial filtering to reduce overfitting and improve model quality (Boria et al., 2014; Tourne et al., 2019) and kept one record per grid cell (spatial resolution of 0.083° × 0.083°) to reduce record clustering, which often results from sampling bias (Boria et al., 2014; Sillero and Barbosa, 2021). Thus, we were left with 2130 presence points to use in model training.

2.2 Environmental predictor variables

The environmental predictor variables, i.e., the covariates, were selected based on data availability and their physiological and ecological importance to marine species. So, we included water temperature (°C), salinity, current velocity (m.s^-1), and depth (m) as covariates. The first three variables were downloaded as mean for the surface marine realm and the latter as mean depth. We used depth to estimate rugosity, an index often used to predict species distribution and suitable habitats, as well as an indicator of hard-bottom habitat (seafloor heterogeneity) (Dunn and Halpin, 2009; McArthur et al., 2010; Fonseca et al., 2017; Guerra et al., 2021). The rugosity index was derived from the depth layer using the ‘Terrain Ruggedness Index’ of the raster package in the R software (R Core Team, 2021). The rugosity index was calculated as the mean of absolute differences between a cell’s value and its eight surrounding cells. Low rugosity values indicate unconsolidated substrate (i.e., mud and sand), whereas high rugosity values are associated with rocky substrates (Riley et al., 1999; Fonseca et al., 2017). All covariate layers were downloaded from the public database Bio-Oracle v2.0 which provides environmental variables under current conditions (2000-2014) and projections for the 2040s (2040-2049) and 2090s (2090-2099), considering the concentration of greenhouse gas emissions in the Representative Concentration Pathway (RCP) trajectories of climate change adopted by the IPCC (Intergovernmental Panel on Climate Change) (Tyberghein et al., 2012; Assis et al., 2018). We used the moderate (RCP 4.5) and the most pessimistic trajectory projections (RCP 8.5) to predict areas of occurrence of the Atlantic blue crab for the 2040s (mid-term) and 2090s (long-term) (Supplementary Table 1; Supplementary Figure 1).

Pearson’s correlation coefficient (r) and Variation Inflation Factor (VIF) were applied to detect correlation and collinearity between covariates by using the R packages corrplot and usdm, respectively (Zuur et al., 2010; Naimi et al., 2014; Wei and Simko, 2021). The lowest and highest correlation recorded were 0.18 and 0.39, respectively (Supplementary Figure 2). VIF values were lower than 2, ranging from 1.13 to 1.28 (Supplementary Table 2). Therefore, no significant correlations and collinearity were detected, and no covariate was excluded before data modeling. The Bio-ORACLE dataset was retrieved on March 19, 2022, from the Bio-ORACLE’s website (Bio-ORACLE, 2022).

2.3 BART modeling

Bayesian Additive Regression Trees (BARTs) were performed to predict the probability of the Atlantic blue crab’s presence as a function of a set of environmental covariates. BART is a Bayesian approach to non-parametric function estimation with the ability to capture complex interactions and nonlinear relationships among predictors. This is accomplished by using a sum-of-trees framework (Chipman et al., 2010). BART model can be expressed as:

Y = f (x) = \sum_{j = 1}^{m} g (x; T_{j}, M_{j}) + ϵ

where $Y$ is the $n \times 1$ vector of response variable, $x$ is the $n \times p$ design matrix (predictors); $m$ represents the distinct regression trees, each composed of a tree structure $T$ , and parameters of the terminal nodes (leaves), as an ensemble $(T_{1}, M_{1})$ , …, $(T_{m}, M_{m})$ , denoted by $M$ ; $g (x)$ is each prediction of a set of n trees; and $ϵ$ = $n \times 1$ is vector noise, representing independent mean-zero normal error with variance σ. A set of posterior draws $f (x)$ of sum-of-trees functions, generated by the Markov Chain Monte Carlo process, create the posterior probability for $P (f / y) = P (t r e e s / d a t a)$ . For binary response $Y$ , like absence (0) presence (1) data, BART uses a logit link: $P (Y = 1 / x) = Φ (f (x))$ , where Φ denotes the standard normal cumulative distribution function. A complete description of BART models can be found in Chipman et al. (2010).

The absence data were randomly generated using the function ‘randomPoints()’ (R package dismo), as pseudo-absences. The total number of pseudo-absences was the same as the number of presences, considering known species distribution and the statistical modeling approach, as recommended by Babet-Massin et al. (2012) for tree-based models (Supplementary Figure 3). Then, both presence and absence data were merged into the BART model.

BART was performed with the R package embarcadero, an SDM wrapper for the dbarts package (Chipman et al., 2010; Carlson, 2020; Dorie, 2022). The BART model was performed with all covariates (current scenario from both native and non-native areas), followed by an automated variable selection procedure using the ‘step.model()’ function. This model was trained on current data and fitted 100 times for six different settings of ensemble size (n= 10, 20, 50, 100, 150, and 200 trees) using the ‘bart()’ function (Carlson et al., 2022). This procedure was conducted to select the model with the most informative covariates and the lowest average model root mean square error (RMSE). Then, the final model was performed with 200 trees, 1000 posterior draws with a burn-in of 100 draws, and values for hyperparameters of power equal to 2.0 and base equal to 0.95 for tree regularization prior, as recommended by Chipman et al. (2010). The final BART model was used to pedict the probability of occurrence in future scenarios using the function ‘predict()’.

2.4 Model validation and evaluation

The three hyperparameters, which govern priors on Tj, Mj, and σ, were tuned in for a full cross-validation experiment to balance accuracy prediction and uncertainty model prediction intervals. The ‘retune()’ function was used to run the experiment across the k hyperprior (values of 1, 2, and 3), power parameter (1.5, 1.6, 1.7, 1.8, 1.9, and 2.0), and the base parameter (0.75, 0.80, 0.85, 0.90, and 0.95). Model accuracy was assessed by using the Area Under the Curve (AUC) of the Receiver Operational Characteristic (ROC) approach and True Skill Statistics (TSS), which considered commission (false positives: type I error) and omission errors (false negatives: type II error) (Allouche et al., 2006). Both AUC and TSS methods measure the ability of the model to correctly discriminate presence and absence data; however, TSS also accounts for both model’s sensitivity (i.e., the proportion of presences accurately predicted) and specificity (i.e., the proportion of absences accurately predicted). AUC values between 0.5 and 0.7 indicate poor prediction, from 0.7 to 0.8 is an acceptable prediction, from 0.8 to 0.9 is excellent, and greater than 0.9 is outstanding (Hosmer and Lemeshow, 2000). TSS is defined as ‘sensitivity + specificity - 1’ and ranges from -1 to +1, where +1 indicates perfect agreement and values of zero or less indicate a performance no better than random (Allouche et al., 2006).

The raster maps of the predicted posterior mean, first (0.025), and third (0.975) probability quantiles were generated from the BART model. Spatial uncertainty was calculated by mapping the differences between the two quantile rasters that generate a credible interval width (Carlson, 2020). Future layers of the predicted values were subtracted from the current predicted layer to obtain habitat change trends. The distribution shifting map under each climate change scenario was generated with positive, zero, and negative grid values that represent habitat range expansion, no change, and habitat range retraction, respectively.

All raster layers had a spatial resolution of 0.083° × 0.083° cell size (approximately 9.2 km at the equator) and were cropped to depths up to 200 m in the coastal zone of the Atlantic Ocean, encompassing the depth range of the Atlantic blue crab. This species occupies shallow estuarine and coastal waters, down to 100 m depth in offshore habitats (Hill et al., 1989; WoRMS, 2022). However, the maximum depth of the calibration area also considered that the species may shift its distribution towards deeper waters under climate change scenarios (Dulvy et al., 2008; Hastings et al., 2020). All statistical analyses were performed using R (version 4.0.0) (R Core Team, 2021). All maps were created using the QGIS Geographic Information System, version 3.24 (Open-Source Geospatial Foundation Project, http://qgis.osgeo.org).

3 Results

3.1 Presence probability and covariates

In general, the TSS (0.69) and training AUC (0.92) values indicated that the BART models performed very well on data validation with cut-off thresholds of 0.54 (Supplementary Figure 4). The final BART model included temperature as the most informative variable with the highest average importance value, followed by salinity and current velocity (Supplementary Figure 5A). (Supplementary Figure 6) On the other hand, rugosity was dropped during the covariates selection process due to the lowest RMSE value (~0.3315) (Supplementary Figure 5B).

The relationship between the probability of presence and temperature showed that the Atlantic blue crab occurred in areas with water temperatures higher than 10°C, with a high probability of being found in marine habitats between ~20 and 25°C (Figure 2A). The species also occurred in a wide range of salinity, but its presence probability decreased in waters with values higher than 31 (Figure 2B). The probability of presence decreased with increasing current velocity since over 70% of occurrences happened in calm waters (Figure 2C). The relationships between temperature – the most important covariate – and other covariates revealed the two-dimensional range of environmental conditions where the species showed the highest presence probabilities (Figure 3). Generally, high presence probabilities were associated with temperatures between ~20 to 25 °C, salinity between 25 and 32 (Figure 3A), and current velocity slower than 0.20 m.s^-1 (Figure 3B).

FIGURE 2

Figure 2 Smoothed fits of covariates modeling the spatial distribution of the Atlantic blue crab Callinectes sapidus Rathbun, 1896 for temperature (A), salinity (B), and current velocity (C). Gray shading indicates a 95% credible interval from the posterior distribution of the current scenario BART model.

FIGURE 3

Figure 3 Two-dimensional partial dependence plots for the Atlantic blue crab Callinectes sapidus Rathbun, 1896 distribution model. The contour plots show the relationship between temperature and the other covariates in predicting the species distribution, namely salinity (A), and current velocity (B). The contour lines represent the posterior median of the $f (x)$ draws.

3.2 Prediction maps, climate change scenarios, and uncertainty

The BART model for current conditions identified areas of importance for the Atlantic blue crab. The highest presence probabilities were found in the Gulf of Mexico, the southeast coast of the United States, and along the southeast Brazilian coast to the Rio de la Plata estuary (Uruguay and Argentina) in the native area, and in the Mediterranean Sea and off the Mauritanian coast in non-native areas. Sites with a medium to high presence probability were also identified in native and non-native habitats (Figure 4A).

FIGURE 4

Figure 4 Predicted distribution of the Atlantic blue crab Callinectes sapidus Rathbun, 1896 (1832-2022) for posterior predictive mean for the BART model in the current scenario (A) and habitat range changes for the RCP 4.5 2050 (B), RCP 8.5 2050 (C), RCP 4.5 2100 (D) and RCP 8.5 2100 (E) scenarios.

For both climate change scenarios, the distribution shifts of the Atlantic blue crab will be less pronounced in native than in non-native areas. Overall, changes in habitat suitability (retraction or expansion) affect the probability presence under the extreme RCP 8.5 in a higher degree scenario than in the moderate RCP 4.5 scenario. Habitat changes in the 2090s were greater than in the 2040s as well. For both scenarios, the model predicted pronounced habitat retraction in the Gulf of Mexico and expansion from the Caribbean Sea to the northern Brazilian coast and eastern Canada (Figures 4B, C). In the non-native area, and particularly by the end of the century, extensive habitat expansion is predicted to occur from the coast of Mauritania to Morocco and from southern Portugal towards northern Atlantic regions, including the North and Baltic Seas, as well as a retraction in the western areas of the Mediterranean basin and expansion in the eastern Mediterranean basin, including the Baltic Sea (Figures 4D, E) (Supplementary Figure 7).

The spatial analysis, which considered the temperature effect as an indicator of thermal habitat suitability on the predicted probability of the Atlantic blue crab over time (current, the 2040s, and the 2090s), also disclosed that temperature plays a decisive role in the spread of the species along the south and north Atlantic Ocean, and Mediterranean Sea (Figures 5A–E). Additionally, thermal habitat suitability showed a similar pattern for the 2040s and 2090s in the native range, while in the non-native range, it increased drastically in the 2090s (Figures 4D, E).

FIGURE 5

Figure 5 Partial spatial dependence plots for the Atlantic blue crab Callinectes sapidus Rathbun, 1896 showing the average contribution of temperature to the BART model predictions: current scenario (A), RCP 4.5 2050 (B), RCP 8.5 2050 (C), RCP 4.5 2100 (D), and RCP 8.5 2100 (E).

Model uncertainty also varied over time and between native and non-native habitats. Uncertainty was higher in the 2090s for both RCP 4.5 and 8.5 scenarios, and it was generally associated with deeper waters. In the native range, uncertainty predictions were higher in areas located in the Caribbean Sea (e.g., coasts of Honduras, Nicaragua, and Costa Rica), northern Brazil, and eastern Canada. In the non-native range, the highest uncertainties areas were mainly found in western Africa (Guinea-Bisau, Cameron), Morocco, and the Baltic Sea (Figure 6). The 2.5% and 97.5% predicted intervals of distribution probability used to estimate uncertainty across models and climate change scenarios are shown in Figure 6.

FIGURE 6

Figure 6 Prediction of the Atlantic blue crab Callinectes sapidus Rathbun, 1896 distribution between 1832 and 2022 (A–C) and under the influence of RCP 4.5 2050 (D–F), RCP 4.5 2100 (G–I), RCP 8.5 2050 (J–L), and RCP 8.5 2100 (M–O). The uncertainty areas are shown in the left column, the 2.5% quantile in the center column, and the 97.5% quantile in the right column.

4 Discussion

4.1 Advantages and limitations of BART models

This study shows that BARTs amplify our understanding of how environmental factors affect the occurrence and habitat preferences of the Atlantic blue crab in native and non-native areas. Equally important, BARTs projected the species distribution under two climate change scenarios, and showcases the applicability of a high-performance machine learning algorithm to study the distribution of an aquatic species.

Our results about the current and future distributions of the Atlantic blue crab were based on the assumption of changes in thermal boundaries. However, other factors affect species distribution, including biotic interactions, adaptative genetic drift, evolutionary changes, migration and dispersal, and anthropogenic stressors (e.g., fishing pressure, species introductions) (Hines et al., 2010; Gamliel et al., 2020; Pinsky et al., 2020; Seaborn et al., 2020; Encarnação et al., 2021a; Mancinelli et al., 2021). Therefore, as modeling and computing capabilities improve the inclusion of such data in future studies will improve prediction accuracy and reduce uncertainty. Yet, current results already provide timely information to managers to elaborate data-driven management policies.

We must also stress that we used environmental data representing the realized niche in the native range and that recently introduced non-indigenous species are often in non-equilibrium with the dynamics of the ecosystem. Invasive ecology theory reports numerous situations when non-indigenous species became invasive and then saw the population numbers collapse (Hui and Richardson, 2017). Therefore, projections about non-indigenous species distribution shifts must always be cautious because of a mosaic of introduction timelines along the non-native area (Hui, 2022). This is the case of the Atlantic blue crab, first reported in Europe in 1900 and with dissimilar invasion timelines between the North Sea, the eastern Atlantic coasts, and the Mediterranean Sea (Welk, 2004; Gallien et al., 2012; Pili et al., 2020; Freeman et al., 2022). Conscious of this limitation, we cannot guarantee that the niche-environment equilibrium assumption of SDM models have always been met.

We opted to analyze the distribution range shifts of the Atlantic blue crab at an ocean scale. At such scale, we are unable to capture the fine spatial resolution required to evaluate the expansion or retreat into and from estuarine ecosystems, and even freshwater ecosystems further inland, despite their importance for the species (Mangum and Amende, 1972; Taylor and Fehon, 2021; Scalici et al., 2022). Another factor precluding us from a finer resolution is the poor global database of environmental information covering inland ecosystems, despite the inclusion of data from several large estuarine ecosystems in the Bio-ORACLE database (e.g., Black Bay, Corpus Christi Bay, Chesapeake, Bay, Mississippi Sound estuary, and Lake Borgner in the US, Rio de La Plata in Uruguay and Argentina, Términos Lagoon in Mexico, and Maracaibo Lake in Venezuela) (Tyberghein et al., 2012; Assis et al., 2018). Therefore, the absence of environmental information on inland ecosystems in the raster layers does not invalidate the main results and conclusions of our study.

4.2 Model performance and effect of covariate

Overall, the BART model performed very well and the robust results illustrated the significance of climate change on the distribution shifts of an aquatic species. Both the AUC (0.92) and TSS (0.62) derived from the BART model achieved excellent performance, showing a high degree of discrimination between the sites where the Atlantic blue crab was present and those where it was absent. In addition, TSS values indicated the models’ high ability to correctly predict true negative (specificity) and true positive rates (sensitivity). The high performance of machine learning techniques has been reported for other tree-based models used to predict the distribution of other marine species (e.g., boosted regression trees, random forest) (Pittman and Brown, 2011; Mainali et al., 2015; Hernández-Ucera et al., 2021; Rahmanian et al., 2022).

The distribution of marine species can be influenced by several environmental predictors. However, only four to five variables can accurately predict the current and future distributions of marine species (Lezama-Ochoa et al., 2016; Zhang et al., 2019; Guerra et al., 2021; Hernández-Ucera et al., 2021). Temperature, salinity, current velocity, primary productivity, and seafloor parameters (e.g., slope, roughness, rugosity) are amongst the most informative predictors used in species distribution models in marine ecosystems (Bosch et al., 2018; Guerra et al., 2021; Pickens et al., 2021). Water temperature was the most informative environmental variable in controlling the potential distribution of the Atlantic blue crab, which aligns with the species’ temperature-dependent trait (Brylawski and Miller, 2006; Hines, 2007; Hines et al., 2010; Azra et al., 2020; Rogers et al., 2022). The life history traits of the Atlantic blue crab vary latitudinally because temperature strongly affects growth, molting, reproduction (e.g., maturation, fecundity, larval development time, frequency of spawning), and settlement. The drastic decrease in the species occurrence probability below 10°C and above 25°C in the BART model coincides with the species’ optimal water temperature range. Indeed, the Atlantic blue crab enters a period of dormancy in waters below 10°C (Brylawski and Miller, 2006; Hines et al., 2010), and their mean thermal preference for the Gulf of Mexico population ranges from 26 to 33°C with an optimum temperature at 24°C (García-Rueda et al., 2021; Marchessaux et al., 2022). Yet, we should always be aware that local adaptations reshape populations’ optimal temperature range.

Like many decapod species, the Atlantic blue crab requires various environmental conditions during its life cycle, which coincides with different habitat use strategies (Harris, 1979; Hines et al., 2010). For instance, the growth phase occurs in low-salinity estuarine waters, while ovigerous females migrate toward higher-salinity waters to release the eggs (Hines, 2003; Hines, 2007). The occurrence probability was closely related to temperature but also influenced by the other environmental predictors considered in this study. Thus, the variability in occurrence probabilities (between 25% and 80%) concerning salinity may be explained by the species’ migratory behavior during which they inhabit ecosystems with striking different water characteristics. Our analysis also shows that the highest occurrence probabilities were registered in calm water areas. Although the BART model did not indicate rugosity as an informative critical variable, these calm waters may be showing unconsolidated substrate. The Atlantic blue crab prefers nearshore aquatic submerged vegetation beds, also making use of unvegetated sandy and muddy bottoms for shelter, refuge from predators, and feeding grounds (Wilson et al., 1990; Hines et al., 2010).

4.3 Current distribution

The current prediction for the Atlantic blue crab distribution in the native range indicated three zones with the highest occurrence probabilities that match with areas with a high concentration of observations (Figure 1). These zones consist of diverse habitats, including bays, estuaries, coastal lagoons, tidal marshes, shallow inshore waters with soft bottoms, rocky bottoms, and mangrove wetlands (Miloslavich et al., 2011; Love et al., 2013; NMFS, 2015). The abundance in the Gulf of Mexico is considerably higher in salt marshes than in open waters (Minello et al., 2008), similar to observations on the Virginia coast, United States, where they are associated with seagrass meadows, salt marshes, channels, and oceanic inlets (Cheng et al., 2022). Therefore, vegetation zones likely enhance habitat suitability for the Atlantic blue crab because it provides refuge from predators and a bounty of food.

Our BART model also predicted high presence probabilities for true absence areas off the coasts of west France, Mauritania, Senegal, Gambia, and South Africa. The environmental conditions in these areas, with bays, salt marshes, mangroves, estuaries, and lagoons, seem suitable for the Atlantic blue crab despite the absence of records. These results should be evaluated with caution due to the lack of comprehensive studies in these regions while underscoring the high invasion probability of these areas under the current climate scenario.

Invasion success is often dictated by propagule pressure (Holle and Simberloff, 2005) to circumvent or minimize biological resistance (Brown and Barney, 2021). On the northwestern coast of Africa, propagule pressure is increasing from northern populations along the Atlantic coast of Morocco (Chaouti et al., 2022; Chairi and González-Ortegón, 2022). Our present-day habitat suitability predictions already identified the coast of Mauritania with a high presence probability, meaning that the species may colonize this area in the upcoming years. This prediction is of particular concern for the Banc d’Arguin National Park in Mauritania, a highly productive ecosystem, the largest marine protected area in western Africa, and a UNESCO World Heritage Site (Trégarot et al., 2020). At these latitudes, non-i Atlantic blue crab populations may compete with three native congeners Callinectes amnicola (de Rochebrune, 1883), Callinectes marginatus (A. Milne-Edwards, 1861), and Callinectes pallidus (Rochebrune, 1883) (Manning and Holthuis, 1981).

4.4 Future distribution and climate change

The shifts in aquatic species distribution due to climate change are multidirectional, either towards the poles or the equator or between shallower and deeper waters (Dulvy et al., 2008; Hastings et al., 2020). The Atlantic blue crab distribution is limited by cold water at higher latitudes (Molina et al., 2021), yet climate change will expand its distribution polewards to escape warmer waters at lower latitudes (Azra et al., 2020; Hastings et al., 2020). Our data shows that the temperature-driven expansion/contraction of the species distribution will vary between native and non-native areas as a function of climate change scenarios.

The asymmetric nature of the future distribution on both sides of the Atlantic Ocean owes to the populations’ intrinsic status, native or non-native. In native habitats – especially in the Gulf of Mexico, the Caribbean Sea, and off South America – and under the RCP 4.5 and 8.5 scenarios, our models predicted an expansion towards the equator (e.g., Caribbean Sea, and from the coasts of Trinidad and Tobago to the North of Brazil. The population retraction predicted in the Gulf of Mexico may be associated with water temperatures above 30°C, disrupting the species’ metabolism when it reaches the upper thermal tolerance limits (Marchessaux et al., 2022). The poleward expansion will be more significant later in the century and particularly in eastern Canada, as supported by recent observations in the Gulf of Maine (Johnson, 2015). In non-native ecosystems, the most significant distribution shifts are expected to occur in the 2090s under the RCP 8.5 scenario. The species will be widespread along most of the Mediterranean Sea, between Mauritania and the Bay of Biscay, the North Sea, and even reaching the Baltic Sea. Indeed, there have been sporadic records between the Bay of Biscay and the North Sea (Cabal et al., 2006; Nehring, 2011; Mancinelli et al., 2021) and some already identified established populations off Belgium and the Netherlands since the 1990s (Nehring, 2011). The coastal areas in the North Sea and Baltic Sea are geographically attractive with their numerous estuaries and salt marsh areas, which may favor a steep increase in abundance when coupled with more amenable sea temperatures (Dijkema, 1990).

The environmental selection for traits that enhance reproduction and survival under climate change scenarios may benefit species’ range expansion (Finch et al., 2021). For the Atlantic blue crab, we expect that two mechanistic phenomena will act in tandem to favor the poleward-expansion along the eastern Canadian waters and northern European seas by the end of the 21^st century. First, the poleward-increased water temperature is conducive to a fundamental niche expansion. Second, decreasing water temperature-induced mortality during winter allows consolidation of the expansion obtained during warmer months (Molina et al., 2021). Additionally, its plasticity to adjust the reproductive strategy and behavior according to temperature changes allows the Atlantic blue crab to tolerate, adapt, or move to deal with environmental shifts in native and non-native habitats.

The expansion of another invasive species to northern European estuaries will disturb ecosystems that already endure significant anthropogenic impacts (e.g., excessive nutrient load, habitat degradation, river flow regularization, urbanization) (Damme et al., 2005; Nehring, 2006; Courrat et al., 2009; Dauvin et al., 2009). The range expansion of the Atlantic blue crab in the native area will set new ecological dynamics often equivalent to a biological invasion process (Tulis et al., 2023).

Overall, environmental and fishery managers should be proactively vigilant because the expansion of the Atlantic blue crab in native and non-native areas will likely shift ecological dynamics and increase ecological stress. Managers from the native area will undoubtedly obtain precious insights from lessons learned in the non-native area. Finally, managers should promote or sponsor citizen science campaigns, not only because of their significant role in reporting the distribution of non-indigenous species and tracking their expansion (Encarnação et al., 2021b) but also because they generate data with great value to improve model validation and applicability (Matutini et al., 2021).

5 Conclusions

This study anticipates that climate change will lead the Atlantic blue crab to expand their distribution in native and non-native areas. An asymmetric range expansion will occur on both sides of the Atlantic Ocean, with temperature playing a decisive role in setting such asymmetry. In native habitats and considering the RCP 4.5 and 8.5 scenarios, we predict an expansion toward the equator and deeper areas during the 2040s. A poleward range expansion will be particularly noticeable by the end of the 21^st century off eastern Canada. In non-native ecosystems and considering the RCP 8.5 scenario, the most significant distribution shifts will occur in the 2090s when the species become common even in the North Sea. The range expansion of native species is often equivalent to the dynamics and impacts of invasive species, so lessons learned by ecologists and environmental managers about the invasive Atlantic blue crab will provide precious insights to colleagues in the native areas.

Data availability statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: All occurrence data of the Atlantic blue crab Callinectes sapidus used in the study are available for download from the Global Biodiversity Information System (GBIF- www.gbif.org/), iNaturalist (https://www.inaturalist.org/). Occurrences from 1832 to 2020 can be downloaded from the link: https://www.nature.com/articles/s41597-021-00888-w#additional-information. Additional records made between 2019 and March 2022 are available at https://zenodo.org/record/7421549#.Y6mYdnbx4mI. GIS rasters providing environmental data for current and future scenarios are available at Bio-ORACLE website (https://www.bio-oracle.org/).

Author contributions

EC: conceptualization, methodology, data curation, performed the analyses, writing - original draft, writing - review & editing. JE: investigation, data curation, methodology, writing - review & editing. MT: writing - review & editing. PM: conceptualization, writing - review & editing. All authors contributed to the article and approved the submitted version.

Funding

This research received no specific funding agency in the public, commercial, or not-for-profit sectors. However, the NEMA is funded and supported by the Portuguese national funds from Foundation for Science and Technology (FCT) (Projects: UIDB/04326/2020, UIDP/04326/2020 and LA/P/0101/2020), ATLAZUL (Poctep/Interreg 0755/ATLAZUL/6/E) and ALFCORAZUL projects (Fundo Ambiental, candidature n° 48). JE has a Ph.D. scholarship (SFRH/BD/140556/2018) funded by Foundation for Science and Technology (FCT, Portugal).

Acknowledgments

The authors express their gratitude to the referees for their thoughtful comments and effort toward improving our article.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmars.2023.1158206/full#supplementary-material

References

Aguilar R., Johnson E., Hines A. H., Kramer M. A., Goodison M. R. (2008). Importance of blue crab life history for stock enhancement and spatial management of the fishery in Chesapeake bay. Rev. Fish. Sci. 13, 117–124. doi: 10.1080/10641260701681599