What and where? Predicting invasion hotspots in the Arctic marine realm

Abstract The risk of aquatic invasions in the Arctic is expected to increase with climate warming, greater shipping activity and resource exploitation in the region. Planktonic and benthic marine aquatic invasive species (AIS) with the greatest potential for invasion and impact in the Canadian Arctic were identified and the 23 riskiest species were modelled to predict their potential spatial distributions at pan‐Arctic and global scales. Modelling was conducted under present environmental conditions and two intermediate future (2050 and 2100) global warming scenarios. Invasion hotspots—regions of the Arctic where habitat is predicted to be suitable for a high number of potential AIS—were located in Hudson Bay, Northern Grand Banks/Labrador, Chukchi/Eastern Bering seas and Barents/White seas, suggesting that these regions could be more vulnerable to invasions. Globally, both benthic and planktonic organisms showed a future poleward shift in suitable habitat. At a pan‐Arctic scale, all organisms showed suitable habitat gains under future conditions. However, at the global scale, habitat loss was predicted in more tropical regions for some taxa, particularly most planktonic species. Results from the present study can help prioritize management efforts in the face of climate change in the Arctic marine ecosystem. Moreover, this particular approach provides information to identify present and future high‐risk areas for AIS in response to global warming.

This study aims to identify potential high-risk AIS and predict hotspots of invasion for those under current and projected future environmental conditions in Arctic environments. Cold-tolerant planktonic and benthic AIS were scored to produce a list of the top species with the highest relative likelihood of invasion and impact, since these ecological groups include the dominant known highest risk marine invasive species (Molnar, Gamboa, Revenga, & Spalding, 2008). The distribution of suitable habitats for this set of species was then modelled to predict regions of overlap under current and future projected conditions, thus identifying hotspots of potential biological invasions. Although all analyses were done globally, there was an emphasis at the Canadian and pan-Arctic scales.

| Study region
Distribution model outputs were produced with global coverage to allow evaluation of global and pan-Arctic patterns although analyses were focussed on the potential highest-risk species for the Canadian Arctic. The rationale for this approach is the vast expanse of the Canadian Arctic (Archambault et al., 2010), its vulnerability to invasion (Chan, Bailey, Wiley, & MacIsaac, 2013;Chan, MacIsaac, & Bailey, 2016;Goldsmit, McKindsey, Archambault, & Howland, 2019), the ambiguous status of new species that have potentially arrived through human transport (Goldsmit et al., 2014) and the unprecedented warming that it has experienced over the last decade (Bush & Lemmen, 2019;Niemi et al., 2019). The region is immense, accounting for eight of the 19 Arctic marine ecoregions (Figure 1; Spalding et al., 2007). In fact, Canada has the longest coastline in the world, the majority of which is situated in the Arctic (accounting for almost 2,000,000 km 2 of the territorial sea; Archambault et al., 2010). Even though few AIS have been identified in Canada relative to other Arctic regions, Chan et al. (2019) suggest that this accounts for ca. 20% of all AIS recorded from marine Arctic waters. Recently, additional species have been identified for the first time in the Canadian Arctic using molecular tools Chain et al., 2016;Lacoursière-Roussel et al., 2018) or identified as cryptogenic as their origin is unclear (Goldsmit et al., 2014). Moreover, Canadian Arctic ports (especially those in the Hudson Bay Complex; Figure 1) have been identified as being of moderate to high ecological risk of invasion since they could provide suitable habitat for various AIS Goldsmit, Nudds, et al., 2019) and there is evidence that a number of non-indigenous species are being transported into the region by shipping (Chan et al., 2015;Dispas, 2019;Laget, 2017;Tremblay, 2017).

| Species selection
A three-step procedure was used to select potential marine AIS (zoobenthos, phytobenthos and zooplankton) for modelling ( Figure S1).
Most were species that are considered AIS in other regions of the world and may be transported via shipping to the Canadian Arctic; several others had been previously detected in other Arctic environments. These steps included: (a) prescreening analysis and selection; (b) ranking a subset of these species using an invasive species screening tool; and (c) selecting a final list of higher risk species based on tool results.
The first step of the prescreening analysis was to identify potential Arctic AIS. One hundred species were identified using a combination of data sources including published articles, grey literature and webbased global invasive species databases (see details in Figure S1.1).

Species' biological and ecological characteristics for survival in Arctic
conditions (e.g. temperature and salinity tolerances) and their ability to arrive in the region via shipping were evaluated, identifying a total of 31 potential Arctic AIS ( Figure S1.2a). These species were then ranked to evaluate invasion risk using the Canadian Marine Invasive Species Tool (CMIST; Drolet et al., 2016)-a rapid screening-level risk assessment tool to quantify the risk of existing or potential marine invaders in a given area. The semiquantitative tool uses existing information and expert opinion to evaluate the potential for arrival, establishment, spread and impact by a given species and has been applied in a number of eco-regions within Canada (DFO, 2017;Drolet et al., 2016;Moore, Lowen, & DiBacco, 2018;Therriault et al., 2018).
CMIST scores are computed based on responses to 17 questions related to the likelihood and impact of invasion (see Figure S1.2a) according to the information found in the literature and other sources.
Expert assessor knowledge on the risk assessment areas was used to score the potential impact of a given species based on its known impacts observed elsewhere and the availability of suitable habitats and environmental conditions. Uncertainty is scored for each risk question by assigning a qualitative score reflecting the quality of information available to answer each CMIST question. A Monte Carlo randomization procedure is then used to obtain adjusted risk scores that include uncertainty (Drolet et al., 2016).
CMIST-ranked species with either medium or high mean risk scores (N = 18; Table 1) were retained for more detailed assessment using SDM. To this list, five potentially harmful phytoplankton species were added (Table 1; Figure S1.2b). These dinoflagellate species have been found in ballast water tanks and/or in ballast water exchange zones of Canadian domestic ships that discharge their ballast water into Canadian Arctic ports (Laget, 2017). All are known to have the capacity to produce toxins and have been implicated in harmful algal events throughout the world (Harmful Algal Information System, from the Intergovernmental Oceanographic Commission of UNESCO, http://haedat.iode.org/).
The final species list for modelling thus included a total of 23 known marine AIS or harmful algal species from four ecological groups (zoobenthos, phytobenthos, zooplankton and phytoplankton) that pose potential risks for invasion in the Canadian Arctic. To simplify terminology, hereafter when referring to this suite of species they will be termed as 'AIS'. (See Table S1 for information on species' characteristics and impacts). Likewise, the term 'invasion' is used to make reference to the complete process of a species transitioning all invasion stages (transport, arrival, establishment and spread; Lockwood, Hoopes, & Marchetti, 2007).

| Species data
Global scale occurrence data of marine invaders selected for modelling was compiled for both native and invaded ranges using global biodiversity databases such as the Global Biodiversity Information Facility (GBIF-www.gbif.org), Ocean Biogeographic Information System (OBIS-www.obis.org), and invasive species lists with coordinate location information and specific literature (Table S1). An effort was made to maximize the number and quality of occurrence records to best predict potential distributions by doing a vast and complete search of occurrence records and by checking information to find the original source of those records (García-Roselló et al., 2015;Guisan, Graham, Elith, & Huettmann, 2007;Guisan, Zimmermann, et al., 2007;Lobo, 2008). Both native and invaded ranges were used for training and evaluating all models , since invaded areas provide valuable information on species' tolerances to climatic conditions that may not be present in their native range (Marcelino & Verbruggen, 2015) and may improve predictions for the extent of suitable habitat (Broennimann & Guisan, 2008). A single presence record was counted per grid cell to decrease the possibility of overprediction, hence, occurrence points were considered at the same resolution as the corresponding environmental layers (García-Roselló et al., 2015). All points were verified to ensure that they were in sea grids (and not over land). When necessary, points were moved to the closest sea grid using the Near Proximity tool in ArcMap v10.2.2.

| Environmental data
Marine data layers, prepared specifically for ecological modelling, were used as environmental predictors and downloaded at a global scale from Bio-ORACLE v2 (http://www.bio-oracle.org/). These layers were produced with climate data describing monthly averages for the period from 2000 to 2014 representing recent (hereafter referred to as present) conditions (Assis et al., 2018). They were obtained from preprocessed global ocean re-analyses combining satellite and in situ observations at regular two-and three-dimensional spatial grids (Assis et al., 2018). A set of 37 environmental layers was used from this source comprising bottom and sea surface temperature and salinity, ice thickness, chlorophyll, dissolved oxygen, pH, photosynthetically active radiation and minerals and nutrients (iron, calcite, nitrate, phosphate and silicate; Table 1;   Table S2). Long-term maximum, minimum and mean values were used for most predictors when available (Table S2). Resolution of environmental layers was 5 arcmin (approximate 9.2 km at the equator). Bathymetry and land distance layers were obtained from Aquamaps (http://www.aquam aps.org/; Kaschner et al., 2016), but with a resolution of 30 arcmin (approximately 55 km at the equator; Table S2). Bottom substratum type was not included in the analysis due to the lack of availability of a high-resolution global scale database; suitable benthic habitat was assumed to be present within each study region.
Multiple modelling techniques were used so that the results from each approach could be combined into an ensemble model, which has been shown to reduce the biases that single models may have (Araújo & New, 2007).
MaxEnt is a machine learning method based on maximum entropy that predicts the potential geographic distribution of suitable habitat for species using species occurrence data and various combinations TA B L E 1 List of modelled species and the methodology for species selection. Ecological groups were classified as zoobenthos, phytobenthos, zooplankton and phytoplankton. Selection methods were the Canadian Marine Invasive Screening Tool (CMIST; Drolet et al., 2016), and harmful dinoflagellate species found in ballast of vessels discharging in the Canadian Arctic (Laget, 2017)

Kryptoperidinium triquetrum
No common name found Dinoflagellata Phytoplankton Laget (2017) of environmental data layers as input. The method is one of the most widely used SDM algorithms due to its high predictive accuracy and efficiency in modelling range shifts under future climate change scenarios (Bucklin et al., 2015;Elith et al., 2006Elith et al., , 2011Hijmans & Graham, 2006;Pearson, 2007). It has also recently been shown to outperform other modelling techniques in accurately predicting invasive species distributions (Battini, Farías, Giachetti, Schwindt, & Bortolus, 2019). This was complemented with SDM using a suite of other techniques within biomod2, another well-known and widely used SDM tool (Hao, Elith, Guillera-Arroita, & Lahoz-Monfort, 2019) that was selected because it allows multiple models to be run on the same training and testing datasets to better ensure the consistency and accuracy of resultant ensemble models. The four models were run with biomod2: generalized linear model (GLM), random forest (RF), artificial neural network and generalized additive model (GAM) because each represents a fundamentally different modelling technique. Initially, BIOCLIM (a presence-only model) was also included in runs but was later removed as its performance was very low, in agreement with other studies (Elith et al., 2006), and is thus the least frequently used SDM approach within biomod2 (Hao et al., 2019). More detailed information about each of these models included from biomod2 may be found in Thuiller, Lafourcade, Engler, and Araújo (2009).
MaxEnt generates background data to compare with known presence points. This study used the default option, which generates 10,000 random background points. In contrast, selected models from biomod2 require presence and absence data, hence it was necessary to create pseudo-absence values before running the bio-mod2 models (Thuiller et al., 2020). For consistency with MaxEnt models, 10,000 pseudo-absence points were selected randomly using the biomod2 default settings.
Model predictive power for MaxEnt and biomod2 was evaluated using cross-validation with 70% of the occurrence points chosen randomly and used to train the model and the other 30% to test the complete set of models (Araújo, Pearson, Thuiller, & Erhard, 2005). For MaxEnt, data were partitioned by a random process of k = 500 training and validation iterations (Hijmans, 2012), and for biomod2 the data were partitioned randomly. Training points for MaxEnt were selected by random seeding with the convergence threshold set at 0.00001. The hinge feature in MaxEnt was used as it produces complex yet smooth and ecologically meaningful response curves and improves model performance (Merow, Smith, & Silander, 2013;Phillips & Dudík, 2008).
Environmental predictors for each species were selected based on initial models run with MaxEnt, which has been shown to have good discrimination accuracy regarding the identification of important variables compared to other models (Smith & Santos, 2019). Selected predictors included those that are most typically considered as important for related taxa and have been used in other modelling studies (Table 1; Barnes, 1999;Belanger et al., 2012;Chust et al., 2016;Cusson, Archambault, & Aitken, 2007;Gallardo et al., 2015;Jensen, Mousing, & Richardson, 2017;Leidenberger, Obst, et al., 2015;Wagner, 1977), or were identified by taxonomic experts as being important (G. Winkler, personal communication, 2017; A. Rochon, personal communication, 2018). After the first run, predictors with a relative contribution score of <4% were excluded (Jueterbock et al., 2016). Correlations may lead models to produce erroneous response curves to predictors that do not reflect species physiological tolerances (Marcelino & Verbruggen, 2015). Thus, highly correlated environmental predictors (correlation coefficient ≥ 0.7; Dormann et al., 2013), as determined using the SDMtoolbox (Brown, 2014; Table S2), were evaluated and the one with the highest contribution to the modelling exercise retained for model construction. As per , environmental predictors were identified by evaluating the combination of: (a) the response curves for each species-to evaluate if the predictor behaves in a biologically meaningful way in the model (Marcelino & Verbruggen, 2015); (b) a species-specific Jackknife test-to evaluate the contributions of the various parameters and analyse the importance of predictors; and (c) the estimates of the contribution of each environmental layer to model prediction. A minimum of three environmental predictors were included in each model. Each species was then modelled again with the selected environmental predictors (Table S3).
Model performance was evaluated as the area under the curve (AUC) of the receiver operating curve and true skill statistic (TSS). AUC is the probability that the model correctly ranks a random presence site versus a random site from the study area (Phillips et al., 2009). Values close to 1 indicate good prediction in site discrimination, a value of 0.5 indicates a prediction equal to random and values lower than 0.5 indicate a performance that is worse than random predictions. TSS assesses the accuracy of predictions using sensitivity (proportion of correctly predicted presences) and specificity (proportion of correctly predicted absences) in its calculation (TSS = sensitivity + specificity − 1) and is an appropriate evaluation alternative for model predictions converted to binary maps using a threshold (Allouche, Tsoar, & Kadmon, 2006). It ranges from −1 to +1, where values between 0 and −1 indicate performance no better than random, while a statistically reliable model performance is indicated by values >0.4, excellent models by a minimum of 0.7, and 1 indicates perfect agreement with the model (Allouche et al., 2006). Models with a TSS score of <0.7 and sensitivity = 0 were excluded from the final ensemble predictions.
MaxEnt and biomod2 models were each run five times. The continuous values produced by the models were transformed to binary values to identify predicted suitable and unsuitable habitat since continuous model projections may be difficult to interpret.
In addition, binary presence/absence maps are more useful for risk assessment exercises and are required for some model evaluations (Liu, White, & Newell, 2013 Jiménez-Valverde & Lobo, 2007;Liu et al., 2013), and has been used in other studies with similar objectives (e.g. Bellard et al., 2013;Duffy et al., 2017;Wisz et al., 2015). Following transformation, all binary models were delimited using a threshold for the maximum depth each species could inhabit according to their ecological requirements

| Future projections
Once all models were run, tested and validated, predicted habitat suitability was evaluated under projected global change scenarios at global and pan-Arctic scales. The same set of environmental layers that were used for contemporary models were used for projected future models although only temperature, salinity and ice layers were available for projected future scenarios (Table S3). This

| Analysis of environmental predictors and model evaluation
Varying combinations of environmental predictors contributed to explaining the species distribution models of the ecological groups analysed (Figure 2). Sea surface temperature and land distance contributed to models for all ecological groups, although in varying proportions. Bottom temperature and depth were important in explaining models for zoobenthos and phytobenthos, in addition to ice thickness for the latter. Other main predictors for zooplankton included sea surface salinity, ice and depth, whereas important phytoplankton predictors included nutrients and minerals (especially iron; Figure 2; Table S3). All environmental predictors used to construct the final model for each species were within training ranges. in GAM, GLM and Maxent models), but still well above the thresholds for random site prediction (Table S5). A few individual replicated GAM, GLM and RF models (N = 19) for five species were excluded from model-specific binary outputs (maximum 3/species) due to their having TSS < 0.7 and/or sensitivity = 0; which resulted in model averages based on fewer than five replicates for these species (Table S5).

| Hotspots vulnerable to invasions
From the 23 modelled species, up to 20 were predicted to have overlapping suitable habitat in Arctic regions, mostly located near the coasts surrounding Hudson Bay, Northern Grand Banks/Labrador, Chukchi/Eastern Bering seas and Barents/White seas (Figure 3a).
These hotspots are predicted to have potentially greater AIS richness compared to other Arctic regions in the present and both modelled future scenarios. The areal extent of hotspots vulnerable to invasion is predicted to increase over time; this increase is predicted not only in total area, but also in the number of species projected to encounter appropriate habitat there (Figure 3b).
These same regions were predicted to be hotspots for individual ecological groups, particularly zoobenthos and phytoplankton ( Figure 4a,b,d). When evaluated independently, these groups showed the same pattern of future increased areal extent in the predicted suitable habitat, but at varying taxa-dependent scales (absolute area extent: zoobenthos ~6 × 10 6 km 2 , phytobenthos, zooplankton and phytoplankton ~1-2 × 10 6 km 2 ; Figure 4a-d). Despite these differences in magnitude for absolute predicted future suitable habitat, at a relative scale the percentage change in suitable habitat was predicted to be greater through time for some groups (e.g. zoobenthos). Furthermore, the relative change in predicted suitable habitat  Table S1).
At a global scale, total species richness was higher and largely concentrated on the coasts of Northern Europe and the White Sea, Black Sea, Northwest Atlantic, some regions of the Northwest and Northeast Pacific, Iceland and southern regions of Temperate South America and Temperate Australasia ( Figure S2).   (Figures 6 and 7). However, the same analysis predicted an overall suitable habitat loss at a global scale in both future scenarios (−4.3% by 2050 and −3% by 2100; Figures 6 and 7).

| Distribution change
These results show that most modelled species are predicted to find suitable habitats in colder regions, with a trend towards poleward shifts in future distributions, particularly in northern latitudes (the poleward shift to the south is predicted to be lower; Figure 6).
At a global scale, however, planktonic species are predicted to experience habitat loss (except for a small habitat gain predicted for phytoplankton towards the year 2100), while benthic species are predicted to have a positive gain, although much less than that at the pan-Arctic scale (Figure 7; Figure S3).  (Melia et al., 2016;Miller & Ruiz, 2014;Smith & Stephenson, 2013). Likewise, the Antarctic could be affected by native and invasive species from northern regions, given that environmental conditions are similar.

| D ISCUSS I ON
Arctic ports could act as possible sources of AIS to the Antarctic since some ships transit between poles and could act as a potential pathway of transportation, although survival and establishment in such scenarios are unknown (McCarthy et al., 2019).
The use of CMIST as a screening tool to identify high-risk species for specific ecoregions is a recent approach to stay ahead of the invasion process (Drolet et al., 2016;Moore et al., 2018;Therriault et al., 2018). Although the use of CMIST combined with habitat modelling, as done in this study, places more emphasis on the likelihood of invasion than the likelihood of impact component of risk, such a tool may be used as a first step to develop a ranked watch list to guide early detection and monitoring efforts and to prioritize species for detailed risk assessments and potential regulations (Locke, Mandrak, & Therriault, 2011). CMIST provided the opportunity to identify potential AIS of greater risk to the Arctic by considering known species' characteristics and impacts in other regions of the world where they have already invaded. Eighty-seven per cent of the species included in the present modelling assessment belong to groups with the greatest numbers of known invasions in the marine Arctic ecosystem (Chan et al., 2019). The uneven emphasis of organisms retained for modelling from the different ecological groups (maximum of 11 zoobenthic species versus a minimum of 3 zooplankton species) may, in part, be explained by the fact that most known introduced marine organisms are benthic species (Streftaris, Zenetos, & Papathanassiou, 2005).
These proportions are comparable with previous assessment studies (Leidenberger, Obst, et al., 2015), although the uneven coverage does make general trends more difficult to interpret.
Sea surface temperature and land distance were retained in all models, consistent with other studies that found these predictors to contribute significantly to explaining species distributions (Bradie & Leung, 2017;Leidenberger, Obst, et al., 2015;Stelzer, Heyer, Bourlat, & Obst, 2013). Water temperature has been identified as the most relevant predictor of global marine species distribution and land distance of moderate importance for both benthic and planktonic species (Bosch, Tyberghein, Deneudt, Hernandez, & De Clerck, 2018;Bradie & Leung, 2017). Depth was an important predictor, as shown by Bosch et al. (2018) and Snickars et al. (2014) who suggested that bathymetry is of high relevance for modelling various taxa. Ice thickness was moderately important for phytobenthos and zooplankton, perhaps reflecting limits on distribution of the former due to ice abrasion, changes in light exposure and a preference of ice-free waters for the latter (Clark et al., 2013;Kube, Postel, Honnef, & Augustin, 2007;Pascual et al., 2015;Richardson, 1979). Iron was an important predictor for all modelled phytoplankton species, likely as it is known to be important to phytoplankton growth, abundance, dominance and species distributions (Hecky & Kilham, 1988;Moore et al., 2013) and plays a role in the development of harmful algal bloom species (Doucette & Harrison, 1990;Wells, Mayer, & Guillard, 1991). In Arctic regions, meltwater can be a significant bioavailable source of iron to surrounding coastal oceans (Bhatia et al., 2013;Tovar-Sánchez et al., 2010) and evidence suggests a link between ice and iron from glacial meltwater leading to blooms in some phytoplankton taxa (Aguilar-Islas, Rember, Mordy, & Wu, 2008;Joli et al., 2018). This is important given that increased iron due to a shift from an ice-covered to an open water Arctic Ocean (Screen, 2018;Seneviratne et al., 2018) could create favourable conditions for harmful species that may arrive in the region.

Four Arctic regions (i.e. Hudson Bay, Northern Grand Banks/
Labrador, Chukchi/Eastern Bering seas and Barents/White seas) were identified as potential hotspots for invasions. Invasion hotspots could pose even greater risks if they coincide with major shipping routes, biodiversity hotspots or areas of special interest/importance, which is the case for most of the invasion hotspots predicted in the present study. The predicted invasion hotspots overlap with regions that have been highlighted for being special (with regard to their uniqueness, importance for species life histories, threatened species and/ or habitats, biological productivity, diversity, etc.). For example, some Ecologically and Biologically Significant Areas (also known as EBSAs) and Marine Refuges coincide with the predicted invasion hotspots in

Northern Grand Banks/Labrador, Chukchi/Eastern Bering seas and
Barents/White seas (Kenchington et al., 2011;Speer & Laughlin, 2011;Templeman, 2007;www.dfo-mpo.gc.ca/ocean s/oeabc m-amcep z/ refug es/index -eng.html), as well as in the southern Hudson Bay, which has been identified as being an area of high biological importance (Stephenson & Hartwig, 2010). Ecoregions with a high concentration of polynyas could also be at greater risk, given that they may act as important biological hotspots due to their increased productivity and biodiversity (Marchese, 2015). This may be the case for the Chukchi and Eastern Bering seas, which are known to have concentrations of polynyas (CAFF, 2013). This is consistent with previous studies that have identified biodiversity hotspots in marine and terrestrial ecosystems as being at particular risk to future invasions (Bellard et al., 2014;Li, Liu, Kraus, Tingley, & Li, 2016) and climate change (Ramírez, Afán, Davis, & Chiaradia, 2017). The presence of AIS in biological hotspots could endanger native and endemic species as AIS are considered to be one of the leading causes of biodiversity loss and contemporary extinctions (Bellard, Cassey, et al., 2016;Blackburn et al., 2019) with impacts on ecosystem structure and functioning, including changes in food webs, biomass, flux rates, etc. (Ehrenfeld, 2010). Relating where invasion hotspots overlap with biological hotspots could be crucial for prioritizing conservation efforts.
The Barents Sea has been highlighted as a region that already has a substantial number of invasions and environmental conditions that increase the probability of successful establishment (Chan et al., 2019). Two of the species modelled in the present study, C. opilio and P. camtschaticus, are already established there (Chan et al., 2019;Dvoretsky & Dvoretsky, 2009;Hansen, 2015).
The former is thought to have been introduced by ballast water, while the latter was an intentional introduction (Alvsvåg, Agnalt, & Jørstad, 2009;Jørgensen & Nilssen, 2011). The Barents Sea appears to be in transition from a cold Arctic to warm Atlantic climate regime (Lind, Ingvaldsen, & Furevik, 2018), making it particularly vulnerable to invasion. Indeed, it is predicted to suffer one of the largest future habitat losses by endemic species in the Arctic (Renaud et al., 2019), leaving potential niches available for novel species to occupy. Additionally, Arctic ecoregions will be more exposed to potential future arrivals with further ice reduction and increased navigability of the Northern Sea Route and the Northwest Passage, although much greater investment in infrastructure, navigation and communications would be needed to this end (Buixadé Farré et al., 2014). Nonetheless, both human activities and AIS are likely to increase over time in the Arctic, as has been observed over the last few years (Chan et al., 2019;Dawson, Pizzolato, Howell, Copland, & Johnston, 2018).
At a global scale, planktonic organisms generally showed loss of suitable habitat that, in all but one case, exceeded predicted gains.
This may be due to their being dispersive pelagic organisms which have the capability of expanding rapidly and show extensive distribution changes in response to temperature increases related to global warming (Poloczanska et al., 2013). Predicted future changes in sea surface temperatures, in particular, which are expected to outpace shifts in bottom temperature (Levitus et al., 2012), may be driving this pattern. Indeed, models in the present study showed that sea surface temperature was more important for planktonic organisms, whereas bottom temperature was important for benthic organisms. However, the pattern of predicted change by taxa is different at the pan-Arctic scale where all ecological groups showed high positive suitable habitat gains coinciding with other modelling studies (de Rivera et al., 2011;Jueterbock et al., 2013;Townhill et al., 2018;Villarino et al., 2015;Ware et al., 2016). This could be explained by the fact that temperature increases have been shown to be greater A marked poleward shift was predicted, consistent with other forecasting studies for various types of marine organisms, including invertebrates, algae and fish (Byrne et al., 2016;Chust et al., 2013;Jones & Cheung, 2015;Mackey et al., 2012;Townhill et al., 2018;Valle et al., 2014;Wisz et al., 2015). Poleward shifts may serve as an early warning signal of ecosystem change due to climate warming. Multiple studies across different taxa are showing consistency between observed responses of marine organisms to climate change, particularly in high-latitude regions (Poloczanska et al., 2013). Cold-adapted species typically have narrow thermal windows and low-energy demand lifestyles, making them more sensitive to temperature changes (Poloczanska et al., 2016). Predicted poleward shifts in suitable habitat in the present study are disproportionately located in Arctic areas, as has been shown in Cheung et al. (2009), García Molinos et al. (2015 and Jones and Cheung (2015). Invasive species in the Antarctic region have been limited by physiology at cold temperatures rather than geographic limits, but global warming may remove those physiological barriers and alter distributions (Aronson et al., 2007). However, it should be noted that a substantial number of the species modelled here are already distributed in northern/cold regions, which may bias observed patterns. Although, three modelled species (C. intestinalis, B. schlosseri and U. pinnatifida) are among the most likely to become invasive in the Antarctic Peninsula according to a recent horizon scanning study (Hughes et al., 2020).
The locations with greater reductions of the sea ice season over the last 40 years (Stammerjohn, Massom, Rind, & Martinson, 2012) coincide with the regions of predicted invasion hotspots in the present study, suggesting that these areas are already experiencing changing environmental conditions and that such variations are projected to continue in time, including the probability of a complete icefree Arctic during summer (Screen, 2018;Sigmond et al., 2018). The combination of species invasions and predicted reductions in sea-ice cover could alter reproduction/phenology timing, energy pathways and food-web dynamics with subsequent impacts on production at higher trophic levels (Haug et al., 2017;Poloczanska et al., 2013 and references therein). In this context, AIS could take advantage of new habitats and resources and outcompete native species, which are generally expected to be more sensitive to temperature changes given that they live within a narrow low-temperature range (Peck, 2005 provide valuable information to help manage resources in marine ecosystems that will face increasing anthropogenic pressures (Reiss et al., 2014) and has been shown to be a useful predictive tool to assess various taxa concurrently (Gallardo et al., 2015;Leidenberger, Obst, et al., 2015). SDM may also be particularly useful for regions such as the Arctic, where predicting biodiversity changes under global warming effects is challenging due to the paucity of baseline data for most organisms (Renaud et al., 2019;Wassmann, Duarte, Agusti, & Sejr, 2011).
Distributional change of AIS may increase risk in previously unimpacted areas, potentially creating new problems for wildlife and even human health (such as harmful algal blooms). These types of episodes are being facilitated by climate change (Hallegraeff, 2010;Kordas, Harley, & O'Connor, 2011;Poloczanska et al., 2013Poloczanska et al., , 2016. The present study predicts changes over such an extensive area (~6 × 10 6 km 2 , equivalent to almost three times the size of Greenland) that altered community structure may be widespread.
Results presented here provide information on AIS that pose some of the greatest threats to the Arctic, with areas at greatest risk identified as hotspots. This information is valuable given that aquatic invaders are understudied in the context of climate change (Bellard et al., 2018). The Arctic is predicted to be affected by increased habitat suitability for a number of potential AIS and, given its vast area, could be severely impacted by AIS accumulating in specific locations under both current and future environmental conditions.
Identification of hotspots through SDM, predictions of habitat vulnerability and particular areas of concern could guide ballast water practices and other management actions, including prevention, early detection monitoring, rapid response and conservation planning. Information such as that provided by the present study should help guide how to prioritize management efforts in this unique and vast region.

ACK N OWLED G EM ENTS
Special thanks to Gesche Winkler and André Rochon for providing expert knowledge on zooplankton and dinoflagellates respectively. norm.

DATA AVA I L A B I L I T Y S TAT E M E N T
Data sharing is not applicable to this article as no new data were created in this study. The data that support the findings of this study