Life in the Fast Lane: Modeling the Fate of Glass Sponge Larvae in the Gulf Stream

Introduction

Sponges fill a wide variety of functional roles in many marine ecosystems, having significant influence on local substrates, bentho-pelagic coupling and associations with other organisms (Bell, 2008; Maldonado et al., 2017). While sponges are often widely distributed as isolated specimens, their ecological importance is greater where they form dense aggregations, known as “sponge grounds” (Hogg et al., 2010). Those aggregations have drawn substantial attention in recent years, particularly from ecological and conservation perspectives. However, they remain understudied compared to other structure-forming biogenic habitats such as coral reefs (Bell, 2008).

Globally, destructive anthropogenic activities, such as bottom contact fishing and deep-sea mining, and climatic variability have increased concern for the protection of deep-sea sponges (Hogg et al., 2010) and the ecosystem services they provide. The large vase-shaped glass sponge Vazella pourtalesii, endemic to the western North Atlantic, is a long-lived ecologically significant species that has been the focus of recent research and conservation efforts. The dense sponge grounds formed by this species have local ecological importance, as species richness and abundances are three to four times higher in the presence of V. pourtalesii (Hawkes et al., 2019) and these sponges significantly influence biogeochemical cycling of silicon (Maldonado et al., 2020) and other nutrients (Bart et al., 2020). Consequently V. pourtalesii is considered to be an ecosystem engineer. Beazley et al. (2021) predicted that its distribution under future (2046–2065) climate change could increase up to four times its present-day size and shift into deeper waters and higher latitudes, particularly in its northern range.

Effective conservation management of V. pourtalesii and other deep-sea sponges, including design of appropriate marine protected areas, can be supported by predictive species distribution models (Knudby et al., 2013; Beazley et al., 2018, 2021; Murillo et al., 2018). However, that approach can only provide some of the required information and often does not consider the ability of a species to disperse into isolated areas of predicted suitable habitat. More generally, an understanding of connectivity within and between spatially separated populations is important for habitat management in order to design effective protection networks and enable restoration of damaged habitats (Kool et al., 2013), and thus attracts broad interest from both scientists and managers (Simons et al., 2013; Kenchington et al., 2019; Ross et al., 2020; Wang et al., 2020). Knowledge of connectivity can enhance our comprehension of species distributions and is an important supplement to species distribution modeling, particularly if species responses under future climate projections are to be interpreted reliably (Holloway et al., 2016). Beazley et al. (2021) did not model colonization potential for their habitat predictions for V. pourtalesii and so whether the species will be a “climate change winner” capable of colonizing the large expanse of seabed identified in their models, remains an open question.

A primary tool for assessing the connectivity of benthic marine species are larval dispersal models that couple numerical hydrodynamic data with a tracking algorithm that follows released “particles,” representing larvae, as they move along flow fields (e.g., Bartsch and Knust, 1994; Marinone et al., 2008; Ashford et al., 2010; Paris et al., 2013; Lange and Van Sebille, 2017; Delandmeter and Van Sebille, 2019). If known, larval characteristics and behavior, including such processes as swimming, mortality, buoyancy and response to extrinsic variables, can be incorporated into the tracking algorithm (e.g., Bartsch and Knust, 1994; Wolanski et al., 1997; Metaxas and Saunders, 2009; Gary et al., 2020).

In this study, we developed passive 3-D Lagrangian particle tracking models to investigate the physical connectivity of V. pourtalesii in the western North Atlantic. We aimed to further develop our understanding of the published predicted distribution of this species under present and future conditions (Beazley et al., 2021) and determine how different populations may be connected, isolated and/or retain larvae. We do so by identifying physical connections among currently occupied areas, estimating the maximal extent of larval transport, and determining the influences of location and depth on potential dispersal patterns.

V. pourtalesii was first described (as Holtenia pourtalesii) from samples taken from the Florida Keys, United States, between 282 and 593 m water depth (Schmidt, 1870; Reiswig, 1996). However, its range extends from there to the Emerald Basin on the Scotian Shelf, off Nova Scotia, Canada (Figures 1, 2), where the species reaches its highest densities and forms extensive sponge grounds which have been protected from bottom contact fishing (Beazley et al., 2018). Further south it is generally only found in low numbers for reasons that are not well understood (Beazley et al., 2021). This species’ distribution falls within the Temperate North Atlantic realm and both the Cold Temperate and Warm Temperate Northwest Atlantic biogeographic provinces of Spalding et al. (2007), a distribution shared with at least 1664 other sponge species (Van Soest et al., 2012) such as Tetilla laminaris symmetrica and Clathria (Clathria) prolifera. A single observation from the Azores (Tabachnick, 1994) has been discounted in recent studies on this species, after discussion with local scientists and failure to locate the voucher specimen (Beazley et al., 2021). In Emerald Basin, V. pourtalesii has been linked to areas of low topographic relief and warm, saline waters (>5°C, > 34‰), relatively low oxygen (≈ 2.6–3 ml l^–1) and high nutrient concentrations (Beazley et al., 2018). That area is subject to strong inter-annual variability in average bottom temperature and salinity (Beazley et al., 2018), ranging from 4 to 12°C and 34 to 35.2, respectively, between 1950 and 2015, leading to a conclusion that this species may be resilient to changing climatic conditions (Beazley et al., 2021).

FIGURE 1

Figure 1. Place names and locations of confirmed occurrences of V. pourtalesii along the eastern coast of North America. Occurrences are colour-coded by their assigned Group (1–7). The inset shows the area of the Blake Plateau with the location of the Charleston Bump. Red line indicates the exclusive economic zones of the United States, Canada, and France (with respect to Saint Pierre and Miquelon). Additional place names, general position of major currents and water masses, topography, and undersea features in the study area are provided in Figure 2.

FIGURE 2

Figure 2. Geographic place names, general position of major currents and water masses, topography, and undersea features in the study area. The inset area is outlined in gray on the large map. See also Figure 1 and Supplementary Figure 1.

The physical oceanography of this region is very well studied, which allows for a detailed understanding of the conditions governing dispersal and distributions and enhanced interpretation of connectivity models. Two major current and water mass systems dominate the area occupied by V. pourtalesii: the fast moving Gulf Stream with its associated Warm Slope Water, and the continental-shelf waters from the Gulf of St. Lawrence to Cape Hatteras (Figure 2). The Gulf Stream exits from the Gulf of Mexico, transporting very warm, saline water and flows very swiftly along the shelf break, northwards to Cape Hatteras, North Carolina, where it breaks away from the continent, thereafter heading northeast and east to pass south of Grand Bank (Rossby et al., 1985). Before reaching Cape Hatteras, the Gulf Stream flows over the Blake Plateau, which projects from the continental slope at depths of approximately 700 m. At the northern end of the Blake Plateau, the Charleston Bump rises to about 375 m depth, deflecting the Gulf Stream, intensifying bottom currents and causing downstream eddy formation, mixing and upwelling (Popenoe and Manheim, 2001). The resulting Charleston Gyre becomes more prominent when the Gulf Stream is strong and transports Gulf Stream water inshore (McClain and Atkinson, 1985). Beazley et al. (2021) predicted some increase in suitable habitat for V. pourtalesii in this region in their climate projection models.

After passing Cape Hatteras, the Gulf Stream takes the form of a complex of meanders and eddies which entrain continental shelf waters, forming a new water mass, the Warm Slope Water (WSW), between the shelf break and the Gulf Stream’s northern “Cold Wall.” The WSW circulates in a gyre over the continental slope between Cape Hatteras and Cape Cod (Csanady and Hamilton, 1988; Fratantoni and Pickart, 2007; Kang and Curchitser, 2013) but a portion escapes the gyre, passes the flank of Georges Bank and thence heads east to round Grand Bank between the shelf break and the north wall of the Gulf Stream (Figure 2). South of Nova Scotia, the shelf/slope boundary, which delimits the WSW at the surface, usually lies south of the shelf break but some of that water mass extends beneath the shelf waters and, subsurface, reaches to the continental slope. It intermittently floods over a saddle between two offshore banks, Emerald and LaHave, filling the depths of Emerald Basin (Figure 2). There it provides the warm, saline layer within which the V. pourtalesii aggregations are found (Drinkwater et al., 2003; Hanz et al., 2021).

The waters of the continental shelf from Nova Scotia to Cape Hatteras have their principal origin in the Gulf of St. Lawrence, where a three-layer system forms in summer (two-layered in winter: Koutitonsky and Bugden, 1991). Only the very cold, low-salinity Cold Intermediate Layer (CIL) and, in summer, the warm surface layer, can cross the shallow banks of the eastern Scotian Shelf (Han and Loder, 2003), whence they flow southwestwards, with the Nova Scotia Current following the coast (Figure 2). Clockwise gyres circle the various banks; some moving anticlockwise in Emerald Basin over the V. pourtalesii. These water masses form the surface waters over the sponge grounds in Emerald Basin. All three layers in the Gulf of St. Lawrence outflow merge with cold, low salinity water from the north and are carried around Grand Bank by the Labrador Current. Together, they continue southwestwards along the edge of the Scotian Shelf as a shelf-break current, contributing water to the gyres around the outer banks (Loder et al., 1997, 1998; Han et al., 1999). The CIL and summer surface layers from the Scotian Shelf flow into the Gulf of Maine, which also receives a deep inflow of WSW through the Northeast Channel (Townsend et al., 2015), covering the sponge grounds there. Beazley et al. (2021) predicted range expansion into the Gulf of St. Lawrence in their future scenarios, accounting for much of the predicted increase in suitable habitat for V. pourtalesii.

An alternative water mass, Labrador Sea Water (LSW), traverses the upper continental slope off Nova Scotia (Figure 2). LSW, which forms much of the deep water of the North Atlantic, circulates within the Labrador Sea and the Subpolar Gyre. Although some moves southward, contributing to the Atlantic Meridional Overturning Circulation (AMOC), it does not follow the continental margin. Rather, it breaks away at Flemish Cap, flows east to the Mid-Atlantic Ridge and does not return westwards until far to the south (Biló and Johns, 2019). Only a minor amount of LSW rounds the Tail of Grand Bank where, at upper-slope depths, it lies between the WSW and the continental slope. That tongue of LSW sometimes extends only to the western flank of Grand Bank but can penetrate much further, historically far beyond Cape Cod (Figure 2). Along that extent of the continental margin, the WSW and colder, less saline LSW form what has been described as a “coupled slope-water system,” with the front between the two water masses moving along the continental slope in response to basin-scale atmospheric forcing (Petrie and Drinkwater, 1993; Marsh et al., 1999; MERCINA Working Group, 2001; Greene and Pershing, 2003). In 1997–1998, there was a major intrusion of LSW onto the shelf, which extended to beyond Georges Bank. It flushed the depths of Emerald Basin, lowering temperatures there by more than 4°C, and penetrated into the Gulf of Maine through the Northeast Channel (Drinkwater et al., 2003), flooding the V. pourtalesii grounds in both areas before being replaced by WSW once more after a few months. These events explain the strong inter-annual variability in temperature and salinity (Beazley et al., 2018) experienced by the sponge grounds in these areas, as noted above.

The waters in the Gulf of Maine are modified by local freshwater inputs and by tidal mixing but they generally circulate anticlockwise around the Gulf and then clockwise around Georges Bank (Loder et al., 1998; Townsend et al., 2015). The Gulf of Maine was another key area predicted to see large increases in suitable habitat of V. pourtalesii under 2046–2065 climate conditions (Beazley et al., 2021).

By way of a flow around Cape Cod and another off the southern flank of Georges Bank, the shelf waters pass Cape Cod and continue southwards toward Cape Hatteras, where they contribute to the formation of the WSW. While velocities there are generally low, there is a well-developed current over the upper continental slope south of 39°N (Levin et al., 2018). Beazley et al. (2021) identified a very narrow and fragmented margin of this upper slope area, between Cape Cod and Cape Hatteras, as suitable habitat that did not show large changes in their habitat projections for the future.

Collectively these currents, associated eddies and water masses control to varying degrees the distribution of benthic invertebrate species in this region which rely on physical dispersal of gametes and larvae. For sponges, which are believed to have very short larval dispersal durations of hours to days (Kenchington et al., 2019; see section “Biological Traits of V. pourtalesii for Particle Tracking Modeling”), particle tracking models initiated from sites of known occurrence (Figure 1) can greatly assist in interpreting species distributions models of predicted occurrence by identifying dispersal kernels. Areas of high predicted probability of occurrence that also lie within dispersal ranges of known occurrences are more likely to be populated as both favorable environmental conditions and recruitment trajectories align. Conversely, areas of predicted habitat that cannot be occupied through passive dispersal from known occurrences are less likely to be colonized. Further, particle tracking models can help to identify “source” and “sink” populations which is important information for the construction of protected area networks (Crowder et al., 2000).

Materials and Methods

Distribution of V. pourtalesii

The known distribution of V. pourtalesii was mapped on a 1 km² grid, from presence only records (Figure 1). Occurrence data for waters off Canada were drawn from catches of Fisheries and Oceans, Canada (DFO) research-vessel trawl surveys (2007–2019), DFO and Natural Resources Canada optical benthic surveys (1989–2018) and fishery observer records of commercial fishery bycatches (1997–2007 and 2010–2015). For areas off the United States, occurrence data were obtained from Reiswig (1996), the NOAA Deep-Sea Coral Data Portal (2009–2011) and NOAA Okeanos Explorer research cruises EX1806 and EX1903 (2018–2019). In total, there were known occurrences found in 136, 1 km² grid cells. As those were highly clustered, the occurrences were further classified into one of seven regional groups, numbered sequentially from their southern to northern locations (Figure 1 and Table 1).

TABLE 1

Table 1. Number of Vazella pourtalesii occurrence records, average seabed depth of occurrences, and number of particles released in each simulation model for each of seven Groups, and the six occurrence sites in Group 4.

Biological Traits of V. pourtalesii for Particle Tracking Modeling

Empirical data on the distribution of early life-history stages of benthic marine organisms, and sponges in particular (Leys and Ereskovsky, 2006), is extremely limited, and only a few sponge species have been studied to date (Leys and Ereskovsky, 2006). Therefore, particle tracking models are often the only tool available for assessing larval connectivity (Ospina-Alvarez et al., 2020). Modeling connectivity through larval drift using particle tracking requires input of at least the spawning season and the duration of the planktonic phase (pelagic larval duration; PLD). The timing of V. pourtalesii spawning is not currently known and all Hexactinellida are currently thought to be viviparous (Leys and Ereskovsky, 2006). Preliminary histological evidence suggests larval release in summer (July and/or August; Dr. A. Riesgo, Museo Nacional de Ciencias Naturales, Madrid, Spain, personal communication). Consistent with this initial observation, a sediment trap deployed from September to June within the sponge grounds on the Scotian Shelf contained no V. pourtalesii larvae (Hanz et al., 2021), as has been documented for other invertebrate larvae (Hanna, 1984) which are thought to passively sink through the water column, and the benthic lander supporting the traps was not colonized. We therefore modeled releases during each calendar month but, for some questions, focused solely on July and/or August, the most likely reproductive season.

The PLD of V. pourtalesii is also unknown, but all sponge larvae observed to date from the limited number of observations made, have been observed to be lecithotrophic and thus their planktonic phases are thought to be very short, generally from a few minutes to about 2 weeks (Leys and Ereskovsky, 2006; Maldonado, 2006; Kenchington et al., 2019). To delimit the maximum connectivity for V. pourtalesii, we therefore assumed a PLD of 2 weeks.

Most sponge larvae are ciliated (Maldonado, 2006) and some may be capable of a mean swimming velocity of 5.78 ± 2.25 mm s^–1 (Montgomery et al., 2019, 2020). However, those results have been challenged (Lanna and Riesgo, 2020) and even if valid, there remains a lack of information on the directionality of sponge larval swimming behavior needed to make use of the information in particle tracking models. Vertical distribution has important consequences for larval transport and dispersal trajectories (Edwards et al., 2007; Gary et al., 2020; Wang et al., 2020). However, the vertical movements of V. pourtalesii larvae, like those of many sponge species, remain unknown (Kenchington et al., 2019), although settlement of V. pourtalesii has been observed on artificial substrates about 5 m off the bottom (Busch et al., 2020). In the absence of better information, larval release was modeled at three depths: the average seabed depth of recorded occurrences of V. pourtalesii in the relevant regional group (Table 1), mid-water from the seabed to the surface, and the surface, to capture three vertical strata that might be reached through active swimming behavior (Vaz et al., 2016).

Lagrangian Particle Tracking

Three-dimensional (3-D) passive particle tracking was performed using the Parcels framework version 2.1 (Lange and Van Sebille, 2017; Delandmeter and Van Sebille, 2019) in a spatial domain bounded by the 20°–50°N parallels, the 40° and 85°W meridians and the coastline. Climatological mean and monthly–mean vertical and horizontal currents, averaged over 1990–2015 and for selected years, were obtained from the Bedford Institute of Oceanography North Atlantic Model (BNAM) (Brickman et al., 2016, 2018; Wang et al., 2016, 2018, 2019). BNAM is an eddy-resolving model of the North Atlantic Ocean, based on Nucleus for European Modeling of the Ocean (NEMO) 2.3, with a nominal resolution of 1/12°. BNAM can resolve eddies with a diameter greater than ∼80 km in the Gulf Stream region. For eddies smaller than 80 km and greater than 50 km (meso-scale eddies), BNAM will permit this type of eddy to exist, but cannot resolve them, while eddies smaller than 50 km (small-scale eddies) cannot be resolved. BNAM models the entire water column in 50 levels, their thicknesses increasing from 1 m at the surface to 200 m at a depth of 1,250 and 460 m at the bottom of the deepest basins (Wang et al., 2018). Partial cells were introduced in BNAM for the bottom layer, improving representation of the seabed; bottom depths being taken from the 2019 General Bathymetric Chart of the Oceans (GEBCO) global bathymetric grid.¹ BNAM results have been independently validated using a variety of observational data (Wang et al., 2018).

Ocean Model Data

Our particle tracking model experiments used different extractions from BNAM. The first was long-term data (1990–2015) averaged over the full time series, referred to hereafter as “long-term averaged currents.” Beazley et al. (2018, 2021) used BNAM long-term averaged current velocities, along with other similarly extracted oceanographic variables, as environmental predictors in their species distribution modeling of V. pourtalesii, and use of long-term averages as environmental predictors are widely used in species distribution models (e.g., Morato et al., 2020; Downie et al., 2021). Similar to the long-term averaged currents we also extracted data averaged over individual months across the full time series, referred to hereafter as “long-term averaged monthly currents.” These were needed to examine seasonal differences in connectivity. In that data set all data from a given month, e.g., July, was extracted and averaged. Our use of long-term averaged currents extracted from BNAM (Wang et al., 2018) provides important information on overall mean dispersal patterns (Edwards et al., 2007) and enhances broad-scale understanding of species distributions (Beazley et al., 2018, 2021). However we recognize that anomalous events and eddies may also play a role in dispersal and delimiting distributions, and that averaging over the long-term masks some of this variation. To investigate the roles of inter-annual variability and narrower averaging windows on connectivity, we additionally selected six individual years within the 1990–2015 time period to examine the effects on dispersal during periods of Gulf Stream variability (see section “Selection of Years for Examination of Inter-Annual Variability”). Climatological monthly averaged currents were extracted from BNAM for each of those years individually, and used to compare the dispersal trajectories with models using the long-term averaged currents and long-term averaged monthly currents. For context, Kang and Curchitser (2013) reported that the average duration of all anticyclonic eddies in this region is 31 days, with maximum durations of 165 days. We refer to these data as “monthly averaged currents within selected years” hereafter to distinguish from the long-term averaged monthly currents. Using the monthly averaged currents within selected years data greatly enhanced the eddy fields in the study area over the long-term averaged currents (Supplementary Figure 1). For the selected year 2010 only (see section “Selection of Years for Examination of Inter-Annual Variability”), 5-day averaged current products were available from BNAM. We additionally extracted those data for July 2010 in order to further examine the effect of using averaged data on the connectivity matrices.

Selection of Years for Examination of Inter-Annual Variability

The Gulf Stream is the essential component of the upper limb of the Atlantic Meridional Overturning Circulation (AMOC), transporting heat from low to high latitudes, and modulating climate at both local and global scales. Change in the strength of the Gulf Stream can be a good indicator for the changes in the AMOC. Wind and heat flux are the two dominant driving factors for ocean circulation patterns, and wind is believed to drive the western boundary currents, such as the Gulf Stream. Winter convection in the Labrador Sea resulting from heat flux has been found to be important in driving the Labrador Current, which interacts with the Gulf Stream at the tip of the Grand Banks, potentially impacting the Gulf Stream. Numerous studies have identified the significance of the North Atlantic Oscillation (NAO) variability on the AMOC. Wang et al. (2019) demonstrated that the NAO’s role on AMOC variability is through the wind-driven Ekman transport component of the AMOC, and its impact is shown on changes of the AMOC at high frequency. Winter convection is an essential mechanism for the AMOC variability at low frequency (Wang et al., 2019). Wang et al. (2019) reported that the AMOC was high before the early 2000s, and became weaker afterward.

In this study we chose the years 1993, 1996, 2001, 2008, 2010, and 2015, to investigate single year connectivity and inter-annual variability in connectivity patterns. Each year was selected for the expression of particular features in the BNAM model: 1993 was a year in the high winter NAO years during the early 1990s with a strong AMOC; 1996 saw a large drop in the winter NAO index, however, there was a strong AMOC; 2001 was the year that Wang et al. (2016) reported a shift of the Labrador Current, and in which year the winter NAO action centers shifted westward; Hakkinen and Rhines (2009) reported a shift of the surface current in the North Atlantic Ocean starting from the early 2000s as well; 2008 was the year when the winter deep convection event in the Labrador Sea returned after years of absence (Yashayaev and Loder, 2009); 2010 was a year with a low winter NAO index, and Häkkinen et al. (2013) reported that the principal component (PC1) of the sea surface height in the North Atlantic Ocean had an abrupt change that year; and 2015 was a recent year with a high winter NAO index. The years 2001, 2008, 2010, and 2015 all had weak AMOC compared to the selected years from the 1990s.

Inter-Model Comparisons

A source of uncertainty in our analyses is the choice of hydrodynamic ocean model used in the particle tracking experiments. Differences among model outputs are to be expected (Delhez et al., 2004; Callies et al., 2011; Hufnagl et al., 2017) but in areas where strong currents influence distribution, it is anticipated that results will be similar (see Kenchington et al., 2019 vs. Wang et al., 2020). We conducted an inter-model comparison to examine the robustness of our conclusion of no exchange of V. pourtalesii larvae from known occurrences between areas north and south of Cape Hatteras, with an estimated pelagic larval duration of 2 weeks. The two ocean models compared were BNAM and GLORYS.² The GLORYS12VI model is a global ocean eddy-resolving product defined on a standard regular grid with 1/12° horizontal resolution, and 50 vertical levels (interpolated from an Arakawa C native grid). Products are available as daily mean, monthly mean, and pluri-annual-climatology-mean values from 1993 to 2019 (Fernandez and Lellouche, 2021). Here we used daily mean GLORYS data and the monthly mean BNAM data from July in each of the 6 years for particle released from the surface when k_h = 0 m² s^–1, to conduct the inter-model comparisons. Surface releases are expected to show the largest potential for connectivity because of the greater spatial coverage shown in their particle releases compared to those of mid-depth and seabed releases (Supplementary Figures 11–22, 24). The use of the daily mean values from GLORYS, unavailable in BNAM, further allows exploration of the impact of averaging over different time periods on the connectivity matrices.

Both BNAM and GLORYS12V1 are eddy-resolving models and have the ability to nowcast deep ocean mesoscale variability such as eddies and movement of currents and fronts (Hurlburt et al., 2009). However, the two ocean models differ in a number of key aspects: (1) BNAM is a basin-scale North Atlantic Model which requires global ocean model data as open boundary conditions. In BNAM the open boundary data are from the GLORYS reanalysis product (Global Ocean Reanalyses and Simulations, version1; Ferry et al., 2010). GLORYS12V1 is a global model and so no boundary conditions are required; (2) BNAM is a prognostic model, with outputs obtained from governing equations directly, which accords with the physical processes. In contrast, GLORYS12V1 is a data-assimilative model, created by assimilating model output from the NEMO 3.1 ocean model and the LIM2 EVP sea ice model with observational data (Reijnders, 2020); (3) The surface forcing of BNAM is taken from a combination of CORE (Coordinated Ocean-ice Reference Experiments) and NCEP/NCAR reanalysis forcing. Model forcing variables include air temperature, wind velocities and humidity; daily short- and long-wave radiation, and total precipitation (rain plus snow) (Wang et al., 2019). In GLORYS12V1, atmospheric forcings are provided with 3- and 24-hourly frequencies by the ERA-interim dataset, including precipitation and radiative fluxes (shortwave and longwave) corrections (Fernandez and Lellouche, 2021). Many parameter settings are also different between the two models, so discrepancies between trajectories can arise from many aspects.

Lagrangian Transport Modeling

Particle tracking results are sensitive to the resolution of the spatial grid of the underlying model (Putman and He, 2013) and BNAM is limited to scales of 5–10 km. Such a coarse grid is necessary to cover a large spatial domain (McVeigh et al., 2017) such as that covering V. pourtalesii distribution. Oceanographic models with a finer grid could be employed to better understand particle movements but they require extensive ground-truthing, which BNAM has been given (Brickman et al., 2018; Wang et al., 2018). Further, higher resolution models over this spatial extent would require computer resources that are not currently available to most researchers; a situation that will likely improve as access to high-speed computers increases.

Horizontal diffusivity was alternatively modeled as zero or k_h = 250 m² s^–1, the latter value being appropriate for mixing across the Gulf Stream (Bower et al., 1985; Anderson et al., 2000). Use of a non-zero value for k_h is expected to have a greater effect on maximum, than on mean, transport distances and end positions (Edwards et al., 2007), partially simulates small-scale larval swimming behavior through the introduction of a random walk (Codling et al., 2004), and mimics some of the variance in the currents which was mostly removed in the averaging process.

Simulated larvae were advected within the Parcels model using the fourth-order Runge-Kutta method as the integration scheme, with a time step of 10 min. The latter ensures that particles cannot cross a model grid cell in a single step, allowing for more accurate estimates of particle retention. To determine the initial locations for larval release, occurrence records of V. pourtalesii were superimposed on a 0.1° spatial grid. Particles were seeded uniformly over each occupied cell, at 0.002° intervals, resulting in 2,601 particles for each fully populated grid cell (Table 1). Particles overlying the average depth of the occurrence records (Table 1) were retained and used for the simulations. Any particles occurring over deeper or shallower water, such as in regions with high bathymetric variability, were removed prior to modeling in order to relate the model scenarios as closely as possible to known occurrences. The same numbers of particles were used for each model release depth, to facilitate comparisons. The number of particles released in each regional group are provided in Table 1 and ranged from 580 (Group 5) to 75,463 (Group 7), reflecting the number of occurrences in each of the release areas and the slope of the seabed in the 0.1° grid cells. The same procedure was used for particle releases in experiment four (Table 2, section “Model Simulations”) using each of the six occurrence records in that region, each filtered for their specific seabed depth (Table 1).

TABLE 2

Table 2. Summary of sequential particle tracking simulation experiments utilizing long-term averaged and long-term monthly averaged BNAM ocean data, averaged from 1990 to 2015.

The results of simulated larval dispersal experiments were evaluated using connectivity matrices, showing the percentage of particles from each source occurrence cell within a Group that either passed over or terminated in another receiving occurrence cell within a Group. Transit time distributions (TTD; McGuire and McDonnell, 2006) were calculated for particle movements from source occurrence cells within Groups to receiving occurrence cells within Groups. Individual-based dispersal was displayed as histograms of the flow distances traveled from each release point (Kool et al., 2013). Those distances were calculated from individual particle trajectories as the sum across the 2-week PLD of straight-line distances between the start and end positions of each 10-min time step. The median, 95th percentile and maxima of the distances were determined. Dispersal kernels were mapped using individual trajectories and/or averaged end positions.

Model Simulations

Particle tracking models are capable of generating many potential outputs. To clarify their interpretation, this study undertook four discrete experiments (Table 2) all using the long-term averaged currents and/or long-term averaged monthly currents. The first experiment aimed to explore both dispersal (connectivity and maximum flow distances) during each calendar month, and the effect of introducing horizontal mixing by setting k_h at either zero or 250 m² s^–1. The effect of horizontal mixing was pronounced and maintained for all further steps. Although Parcels is a 3-D tracking model, in this region the horizontal velocities created by the Gulf Stream greatly exceed any vertical movement (upwelling or downwelling) and the vertical velocities are generally very small (∼10e^–4 m s^–1). As larval vertical migration and swimming ability is unknown for V. pourtalesii as for most other sponges (section “Biological Traits of V. pourtalesii for Particle Tracking Modeling”), we have created three vertical strata for particle releases in all of our experiments (seabed, mid-water and surface) to allow for swimming-mediated vertical migration, although the seabed scenarios are at this time considered to be the most representative of larval behavior. Results of the first experiment showed little seasonal variation. Therefore, the second experiment assessed the connectivity between occurrence cells among Groups only during the suspected spawning months of July and August using the long-term averaged monthly currents, which were compared to the connectivity simulated using long-term averaged currents (averaged over all months from 1990 to 2015). The connectivity matrices generated in that experiment showed consistent connections in both July and August. Therefore, the remaining steps only considered particle releases in July. In the third experiment, data exported from the second experiment were used to generate frequency histograms of flow distances, for July only. The median, 95th percentile and maximum flow distances were determined. In the fourth and final experiment, larval dispersal from Group 4 was studied in detail as that area is strongly impacted by the Gulf Stream and particles may disperse in various directions (Edwards et al., 2007). Three release sites were located north of the Charleston Bump, and three more to the south, based on the six occurrence records from that area (Table 1). Additional release depths (surface, 100, 200, 400, 600, and 800 m) were used, to better capture variations in swimming-mediated vertical larval transport through the water column. At two release points, bottom depth was insufficient for the deepest releases. Those locations were not used in compiling statistics for the depths concerned.

For the six individual years selected for more detailed study of annual variability in the time series, connectivity between Groups was assessed for releases in July and in January using the respective monthly–averaged currents for each year evaluated (1993, 1996, 2001, 2008, 2010, and 2015) and setting k_h at either zero or 250 m² s^–1. These simulations were also performed for the 5-day averaged currents for July, 2010. Median, 95th percentile and maximum flow distances were calculated for July of each year under both diffusivity settings using the monthly–averaged currents. July simulations reflect the suspected spawning time of V. pourtalesii and allow for comparisons with results obtained from the use of the long-term averaged currents for July. At this season anticyclonic eddies are both large and intense (Kang and Curchitser, 2013). January was chosen as a contrasting month when anticyclonic eddies are fewer in number and on average smaller than in July (Kang and Curchitser, 2013).

Results

Dispersal Characteristics Throughout the Year

The median, 95th percentile and maximum flow distances of released particles (Table 3) and their mean end positions (Figure 3 and Supplementary Figure 2) for each month show that seasonal differences in median flow distances were small for surface releases and negligible for those at the seabed. The 95th percentile and maximum flow distances were an order of magnitude greater than the median in all release depth simulations (Table 3). Adding horizontal mixing (k_h = 250 m² s^–1) greatly increased flow distances (Table 3) but had only a small effect on direction (Figure 3 and Supplementary Figure 2). Depth of release had a small influence on mean end position for Groups 4–7 (Figures 3D–F and Supplementary Figures 2D–F). For Groups 1–3 the influence was larger; particles released at the seabed moved much less than those released at the surface or mid-water (Figures 3A–C and Supplementary Figures 2A–C). This pattern held under both horizontal mixing scenarios of k_h at either zero or 250 m² s^–1.

TABLE 3

Table 3. Median, 95th percentile and maximum flow distances (km) of particles released at each of three depths, by calendar month and alternative values of diffusivity using long-term averaged monthly currents (1990–2015).

FIGURE 3

Figure 3. Mean ending positions of particles released at each of three depths in January, April, July, or October using long-term averaged currents in release sites from (A) Group 1, (B) Group 2, (C) Group 3, (D) Group 4, (E) Group 5, (F) Group 6, and (G) Group 7, modeled with k_h = 0. (For results from modeling with k_h = 250 m² s^–1, see Supplementary Figure 2.) Filled green circles denote mean release position of particles in the Group.

Connectivity Between Groups

Connectivity Based on Long-Term Averaged Currents

The connectivity matrices calculated from tracking particles released in July (Figures 4A,B) were similar to those for August releases, and those produced using annual averages (Supplementary Figures 3, 4). Incorporation of horizontal mixing in the model increased the number of connections (Figures 4A,B) and the area of the dispersal kernels more than their maximum extent (Figures 5, 6 and Supplementary Figures 5, 6).

FIGURE 4

Figure 4. Results of particle tracking simulations using long-term averaged ocean currents from July after a 2-week drift duration for releases from the surface, mid-water depths and the seabed. Connectivity matrices showing the proportion of modeled particles released from each of the occupied grid cells in each of the seven Groups (Source) passing over or terminating in another occupied grid cell in another Group (Receiving Area). Diagonal cells indicate particle retention within the occupied cells within each Group. (A) when k_h = 0 and (B) when k_h = 250 m² s^–1.

FIGURE 5

Figure 5. Trajectories of particles released from occupied grid cells in Groups 1, 2, 3, and 4 [Groups 1 (A,E,I), 2 (B,F,J), 3 (C,G,K), and 4 (D,H,L)] at the surface (A–D), mid-water (E–H) and the seabed (I–L) in July using long-term averaged currents, modeled with k_h = 250 m² s^–1. (Inset shows a close up of trajectories from Group 4 occurrences. Open circles indicate the locations of occurrence records. For results from modeling with k_h = 0, see Supplementary Figure 5).

FIGURE 6

Figure 6. Trajectories of particles released from occupied grid cells in Groups 5, 6, and 7 [Groups 5 (A,D,G), 6 (B,E,H), and 7 (C,F,I)] at the surface (A–C), mid-water (D–F), and the seabed (G–I) in July using long-term averaged currents, modeled with k_h = 250 m² s^–1. (Open circles indicate the locations of occurrence records. For results from modeling with k_h = 0, see Supplementary Figure 6).

Releases in July identified connections among the southerly Groups 1, 2, 3, and 4 but only limited retention within any of them, when modeled with a 2-week PLD (Figure 5 and Supplementary Figure 5). In contrast, retention was identified within each of the northerly Groups 5, 6, and 7 but there was only little connection between the Groups. Group 7, however, was linked unidrectionally to Group 6 (Figure 6 and Supplementary Figure 6). Our simulations showed no exchanges between the four Groups south of Cape Hatteras and the three northern Groups, regardless of release depth or horizontal diffusivity values used. Fewer connections and more retentions were observed with releases at the seabed than with those at the surface, while mid-water releases returned intermediate results (Figures 4A,B). The choice of 2-week PLD was to delimit the maximum possible connectivity in one spawning event. Lesser flow distances are not so much possible as almost certain. If V. pourtalesii larvae settle after a day or so, there would be close to 100% retention within each Group.

Particles released from Group 1 were carried by the Gulf Stream, with the surface releases reaching beyond Cape Hatteras, while the addition of horizontal diffusivity broadened the dispersion of particles, with less effect on flow distance (Figures 5A,E,I and Supplementary Figures 5A,E,I). Particles released from Groups 2 and 3 showed similar end points to those released from Group 1, despite their more northerly start positions (Figure 5 and Supplementary Figure 5). Group 4 releases revealed two different flow modes, some particles being retained in the area, while others moved northward along the continental slope (Figures 5D,H,L and Supplementary Figures 5D,H,L) and is further examined in section “Dispersal from Group 4 on the Blake Plateau using Long-Term Averaged Currents.” Because of the swift transport of particles in the southern portion of the Gulf Stream and the slower movement further north, combined with its own large size, Group 4 received particles from each of Groups 1, 2, and 3, in addition to local retention, when horizontal diffusivity was included in the model. No other Group had equivalent redundancy of larval sources (Figure 4B), potentially conferring increased resilience to Group 4.

Particle releases from Group 5 showed the greatest effect of including horizontal diffusivity in the model. Without it, in seabed releases, 85% of particles remained in the area, whereas addition of horizontal mixing resulted in 94% of particles leaving the Group.

Transit time distribution analysis (Figure 7) for those groups showing connections (Figures 4A,B), suggests that particles released from Group 3 showed strong peaks in the time taken to establish connections. Peaks in connections (with Group 4) were made in 2.5 days in surface models, 5 days in mid-water models and 12 days in models with particles released from the seabed (Figure 7). Similar peaks were seen in the connection times for particles released from Groups 1 and 2 at mid-water and the surface. Group 7 has more protracted and longer transit times compared to the other groups regardless of release position (Figure 7), with no clear peaks forming in the transit time distributions (Figure 7). This is because of the large distance needed to connect to Group 6 (Figure 1) in addition to the speed and direction of the prevailing currents (Figure 2).

FIGURE 7

Figure 7. Transit time distribution for connections between the occupied grid cells in each of the Groups passing over or terminating in another occupied grid cell in another Group for particles released from the four Groups (1, 2, 3, and 7) using long-term averaged currents showing connectivity (Figures 4A,B) with k_h = 250 m² s^–1. (A) Particles released from the surface; (B) Particles released from mid-water depths; and (C) Particles released from the seabed.

Connectivity Based on Currents Within Selected Years

As observed for the connectivity assessment of inter-month variation in release times utilizing long-term averaged monthly currents from 1990 to 2015 (section “Connectivity Based on Long-Term Averaged Currents”), there was little difference between the July and January connectivity matrices in each of the 6 years that were modeled using monthly mean currents extracted for each year individually. In both July and January models with diffusivity k_h = 250 m² s^–1, 191 connections were made, of which only 13 differed between the month/years, mostly between Groups 7 and 6 in the northern portion of the range (Supplementary Figures 7, 8). When diffusivity was set to k_h = 0 m² s^–1 in the July models, 100 connections were made (Supplementary Figure 9), while in January 102 connections were identified (Supplementary Figure 10)—both less than achieved when using k_h = 250 m² s^–1. There were 24 differences in connections between these months, with the only identical patterns showing at the seabed in 2001 and 2008. The greatest difference between months occurred in 2010 at mid-water release depths where 6 differences in connections were observed with more connections in July than in January.

Comparing the number of connections made between the models utilizing July long-term averaged currents (Figure 4) with those made in July of each of the six individual years, with diffusivity k_h = 250 m² s^–1 (Supplementary Figure 7), there was also very little difference (Table 4). At the surface, only one connection did not appear in one of the years (2008; Table 4). Mid-depth releases showed few differences in connectivity (Table 4) with lost connections from Group 7 to Group 6 occurring in 1993 and 2001 but otherwise with the same connectivity patterns. At the seabed, connections were lost between Group 7 and Group 6 on the Scotian Shelf in 1993 and 2001, while in 2008 and 2010 connections were added between Groups 1 and 2 with 4 (Table 4). Of all of these differences, every loss of a connection in 1 year was maintained in another year (Table 4) resulting in no net loss of connections between the models run using data from individual years and those runs using the long-term averaged data. However, the models with data from the individual years identified two new connections at the seabed not seen in the long-term averaged results: between Group 2 and Group 4 (in 2008 and 2010) and Group 1 and Group 4 (in 2010) (Table 4). When k_h was set to 0 m² s^–1 (Table 4) more differences were identified than with k_h = 250 m² s^–1, however most of those were lost connections in individual years; all but one of which was found in another of the years examined. There was only one novel connection added and that was between Group 1 and Group 4 at the surface in 2015.

TABLE 4

Table 4. Comparison of connectivity patterns from three release depths between results utilizing July long-term averaged currents, with those found in the same month in each of six individual years, that is monthly averaged currents within selected years (differences documented relative to the former; differences in retention not shown).

Together these results show that the loss of variability created by using the long-term monthly averaged currents did not strongly affect the connectivity linkages although with diffusivity k_h = 250 m² s^–1, two new seabed connections were identified, and with diffusivity k_h = 0 m² s^–1, one new surface connection was identified, all in the southern part of the range. Novel particle retentions were seen with k_h = 0 m² s^–1 in Groups 1, 5, and 4 not seen in the long-term monthly averaged current simulations. None of the new connections showed connectivity with Group 5, the relatively isolated population between Cape Hatteras and Cape Cod that was thought to be linked to the southern sub-region through Gulf Stream eddies breaking off north of Cape Hatteras.

Examination of the dispersal kernels for these models run with the July or January averaged currents for each of the 6 years with diffusivity k_h = 250 m² s^–1 (Supplementary Figures 11–16, respectively), compared with those produced from the July long-term averaged currents (Figure 8), show enhanced gyre formation at all three release depths in the Blake Plateau area which is most pronounced in 2015. The surface releases show strong offshore dispersal in July and also to a lesser extent in January in many of the individual years, not seen in models using the July long-term averaged currents (Supplementary Figure 1), indicating that loss of larvae to the system through entrainment into the Gulf Stream is likely to be greater than indicated from the models run with long-term averaged currents. None of the scenarios show connectivity from any of the release sites to Group 5, however in January 1996 the models for mid-depth releases show particles from Group 4 reaching a gyre of Group 5 particles (Supplementary Figure 15B). If spawning occurred at that time of the year, this could provide a means of connection to Group 5, assuming that the larval duration was also greater than 14 days. Similarly, in 1993 in the surface models for January (Supplementary Figure 14A), larvae from Group 3 track close to Group 5. Further, these potential connections are only seen when the diffusivity constant is applied, and without it (k_h = 0 m² s^–1) there are no connections at all with Group 5 which remains very isolated with generally high retention in both July and January in each of the years and for all release depths (Supplementary Figures 17–22). July surface releases in 2015 show offshore movement from Group 5 (Supplementary Figure 17F).

FIGURE 8

Figure 8. Dispersal kernels from particle releases in July using long-term averaged currents from each of the occupied grid cells in each of the seven Groups (Table 1) at three release depths (surface, (A,D); mid-water, (B,E); seabed, (C,F) when k_h = 0 (A–C), and when k_h = 250 m² s^–1 (D–F).

Simulations run using the 5-day average BNAM data for July, 2010 with and without the diffusivity constant (Supplementary Figures 23, 24) showed retention in Group 5 at the seabed as seen with the July monthly averaged data for that year, and no connections to Group 5 from or to any of the other Groups. The dispersal kernels produced from the 5-day averaged data covered a greater area than that of the July monthly averaged data and eddies and meanders were more pronounced in the former, particularly with k_h = 0 m² s^–1. There was also less impact of using the diffusivity constant as expected with the smaller averaging window, and the July monthly averaged results with k_h = 250 m² s^–1 are more similar to those produced with the 5-day averages with and without diffusivity factored in. The greatest difference between the trajectories in the dispersal kernels produced using the 5-day average simulations are seen at the surface (Supplementary Figures 11, 17, 24) with the dispersal kernels becoming more similar at mid-depth (Supplementary Figures 12, 18, 24) and the seabed (Supplementary Figures 13, 19, 24). At the surface, with k_h = 0 m² s^–1, well-developed meanders show in the dispersal kernels as the particles track in the Gulf Stream. The Group 4 particles are much more dispersed over the Blake Plateau and some move directly offshore there.

Flow Distances

For any release depth, the greatest distances traveled were by particles introduced into the Gulf Stream where it is strongest in the south, and were progressively less for those released into it further north. Those released in the three Groups north of Cape Hatteras had flow distances ranked in accord with the strengths of their local currents (Supplementary Table 1 and Figure 9). Except for those originating in Group 4, flow distances were greatest for particles released at the surface and least for those released at the seabed when horizontal diffusivity was incorporated in the model but especially when it was not. For Group 4, flow distances of mid-water releases were generally greatest and, with k_h = 0, the maximum movement amongst particles released at the seabed exceeded that of surface releases. Moreover, the distance interval between 50 and 95% of particle connectivity is 487 km at mid-water and 286 km at the seabed, both of which are the largest spread among all groups (Figure 9). This large distance interval disappears when the diffusivity constant is applied, with the distance interval between 50 and 95% of larval connectivity dropping to 95 km in mid-water particle releases and 54 km in seabed releases (Supplementary Table 1 and Figure 9). The same patterns were seen in individual years (Supplementary Table 2) and flow distances were high similar to those produced using the long-term average data. Group 4 shows a strong bimodality in particle distribution at all three depths (Figure 9), reflective of the physical oceanography of the area of the Blake Plateau.

FIGURE 9

Figure 9. Frequency distributions of maximum flow distances of particles released in July using long-term averaged currents from occupied grid cells in each Group at each release depth, modeled with k_h = 0. (A) Group 1, (B) Group 2, (C) Group 3, (D) Group 4, (E) Group 5, (F) Group 6, and (G) Group 7. Bars show the dispersal probability values over the specified distance interval from the release sites, and red solid lines show the cumulative probability. The red dashed and dotted lines indicate the median and 95th percentile of the distances.

Throughout, the addition of horizontal diffusivity greatly increased the flow distances. Group 5 releases were particularly affected (Supplementary Tables 1, 2). Despite the great increases in flow distance produced when a diffusivity constant is applied, the lack of connectivity between areas north and south of Cape Hatteras, North Carolina at surface, mid-water and seabed release depths is maintained (Figure 8).

Dispersal From Group 4 on the Blake Plateau Using Long-Term Averaged Currents

Within Group 4, the depth and especially the location of particle release has a strong influence on mean end positions (Figure 10). The mean end positions and directions were almost unaffected by introduction of horizontal diffusivity. Particles released south of the Charleston Bump (Group 4 #1, #2) moved only short distances to the south or southeast, regardless of release depth. Releases north of the Bump (Group 4 #4, #5, #6) were transported greater distances to the northeast. Surface releases from Group 4 #3 moved like the particles from Group 4 #1 and #2, but those released at 400 m or deeper moved to the northeast.

FIGURE 10

Figure 10. Release positions (black dots), mean ending points (red dots) and net displacements of particles released within Group 4, north and south of Charleston Bump, at various release depths, using long-term averaged currents modeled with k_h = 0. (A) Surface, (B) 100 m, (C) 200 m, (D) 400 m, (E) 600 m, and (F) 800 m. “#” denotes the individual occurrence record in Group 4.

Inter-Model Comparison

Results using both ocean model products support the main conclusion from our analyses, i.e., that documented populations of V. pourtalesii are not physically connected between areas north and south of Cape Hatteras (Groups 1–4), nor between the Scotian Shelf populations (Groups 6 and 7) and those north of Cape Hatteras (Group 5) (Figure 1). The connectivity matrices produced from the daily averaged GLORYS12VI surface products (Supplementary Figure 25) show no connections between Group 5 and any other Group, with no particle retention. Connections were made between the Groups in the southern sub-region ranging from two Groups connecting in 2008 to four Groups connecting in 1996, with Group 4 showing the only retention and that in 2010 only. However, our results showed differences between the particle trajectories run using the surface July monthly averaged BNAM and daily averaged GLORYS12VI for each of the six selected years, as was expected (compare Supplementary Figures 11, 26, respectively where k_h = 250 m² s^–1 in the former and k_h = 0 m² s^–1 in the later, allowing for the greatest opportunity for connections). Gulf Stream transport shows differences in the shape of the dispersal kernels in particular. Figure 11 compares the dispersal kernels and connectivity matrices from particle releases from each of the occupied grid cells in each of the seven Groups (Table 1) from the surface in July, 2010, using BNAM monthly averaged currents (k_h = 250 m² s^–1), BNAM 5-day averaged currents (k_h = 0 m² s^–1 and k_h = 250 m² s^–1), and GLORYS12VI daily averaged currents (k_h = 0 m² s^–1). The BNAM models with k_h = 250 m² s^–1 but differing in the number of days used to average the products shows a much more diffusive pattern in the dispersal kernels with the smaller averaging period, as would be expected. However there are more connections and retentions in the southern sub-region using the monthly averaged results compared with the 5-day averaged currents. Therefore, while the latter disperses particles over a larger area, they are not tracking over the locations of the known occurrences. Common results across all four models are retention in Groups 7 and 6 in the northern region, and connection between Groups 1 and 2 in the southern sub-region (Figure 11).

FIGURE 11

Figure 11. Dispersal kernels and connectivity matrices from particle releases from each of the occupied grid cells in each of the seven Groups (Table 1) from the surface in July 2010, using BNAM monthly averaged currents (k_h = 250 m² s^–1), BNAM 5-day averaged currents (k_h = 0 m² s^–1 and k_h = 250 m² s^–1) and GLORYS12VI daily averaged currents (k_h = 0 m² s^–1) in the particle tracking models.

Discussion

Regional-Scale Connectivity

Particle tracking has shown that connectivity among V. pourtalesii populations may be very different in the south and the north of its range. South of Cape Hatteras, larvae of this deep-living species are released into the swift-flowing Gulf Stream, which will carry them far enough to provide strong, though unidirectional, connections between the Straits of Florida and the continental slope south of Cape Hatteras, including our Groups 1–4. These connections vary based on the degree to which larvae rise in the water column and their persistence there (Figure 5). Similar rapid dispersal over long distances by the Gulf Stream has been seen in studies of other deep-sea species in this region (Young et al., 2012). In contrast, the tracking of particles released in Group 4 demonstrated the complexity that can arise from the interaction of a major current with seabed topography, specifically the Charleston Bump and its associated gyre (McClain and Atkinson, 1985; Popenoe and Manheim, 2001). That gyre may retain V. pourtalesii larvae south of the Bump, including those from Groups 1, 2 and 3, creating redundancy in recruitment in that region. Redundancy in connectivity networks occurs when there is more than one connectivity path leading to a population. This is a positive feature as it helps to mitigate against the risk of recruitment failure by providing alternate connectivity routes. An outcome of such redundancy may be enhanced genetic diversity conferring greater resilience to climate change to the population (cf. Ramos et al., 2018).

Given the weaker currents of the continental shelf and the absence of the Gulf Stream, particles released from V. pourtalesii areas north of Cape Hatteras had lower dispersal distances (Supplementary Table 2) and greater retention within their source areas than those to the south (Figure 4 and Supplementary Figure 7), especially for those released near the seabed or at mid-water depths. However, unidirectional larval transport from Emerald Basin (Group 7) to the Northeast Channel (Group 6) is possible, if the species PLD is as long as 2 weeks (Figure 6). While that connection may be valuable in maintaining the downstream V. pourtalesii populations over the long term and enhancing its genetic diversity, most larvae dispersed away from their parental areas are likely doomed, whereas those retained confer local resilience through high levels of larval supply in suitable habitats. Retention is a particularly important process for reef-forming organisms (Swearer et al., 1999; Jones et al., 2009), leading to enhanced recruitment and establishment of populations that dominate the community, as with the formation of sponge ground habitats in Groups 6 and 7. One of the outcomes from our assessment of inter-annual variation in connectivity was increased retention over that detected using the long-term averaged current data in our models. Groups 1, 3, and 4 in the southern part of the range showed retention in individual years at all depths, while Group 5 between Cape Hatteras and Cape Cod, showed retention in surface releases not found in the models using the long-term averages. These results reinforce the role of retention in maintaining these populations and the need to protect local populations as a first priority.

No connections were found between the V. pourtalesii in Group 5 occurring on the continental slope between Cape Hatteras and Cape Cod, and those elsewhere in any of the models run (using long-term averaged currents, monthly averaged currents in individual years, 5-day averaged currents in July, 2010 and GLORYS12V1 daily averaged currents in July, 2010). However, warm-core eddies are shed from the Gulf Stream after it passes Cape Hatteras south of the Group 5 occurrence (Kang and Curchitser, 2013). While no connections with Group 5 were found in our models of monthly averaged data in individual years, there was some suggestion that larvae from the south (Groups 3 and 4) could connect with a dispersal gyre in the vicinity of Group 5 in some years in January and with the GLORYS12V1 daily averaged currents in July, 2010. Therefore, there is potential for larvae from southern release sites to reach Group 5, particularly at the surface and mid-water depths with diffusivity included, although the 2-week PLD limits the likelihood that such larvae would be viable both in this area and in eddies cut off the Gulf Stream further north. Connectivity from the Northeast Channel (Group 6) to Group 5 is less likely (Figures 9, 8). In individual years particles that reached the continental slope from Group 6 were entrained in currents that brought them further eastward, rather that westward toward Group 5 (e.g., Supplementary Figures 11, 12F, 13F). Further study is needed to investigate the roles of eddies, particularly smaller and short-lived eddies, for the connectivity in this region as they remain the most likely means of explaining the establishment of the Group 5 population. Beazley et al. (2021) found that suitable habitat for this species is very limited in the bight between Cape Hatteras and Cape Cod. Nevertheless, there may be unknown upstream populations that link to Group 5 providing larval sources, and/or creating links between known occurrences which could also explain the establishment of that population. Its persistence is likely maintained through retention (Figure 4).

This summary of regional-scale connectivity appears little influenced by any seasonal cycle and is consistent in the face of uncertainty over the timing of larval release in V. pourtalesii. It is also largely independent of the effects of horizontal mixing. Introduction of a non-zero value of the horizontal diffusion constant (k_h), a controversial alternative in particle tracking studies (Hunter et al., 1993; Kelly et al., 2020), did have a large effect on the modeled transport trajectories but had little impact on the connectivity matrices presented here, only extending the connections among the southern Groups and adding some retention pathways for those groups. It also had little effect on the mean ending positions, with most of the increased flow distance therefore attributed to the random walks at each of the 2,016 time steps in an individual particle trajectory. Therefore including the diffusivity constant increased the area of the dispersal kernel more than the distance from the source (Figure 8), and comparisons with its effect using different averaging periods demonstrated that it did successfully replace some of the variability lost in using the averaged data.

The effects of inter-annual variation on connectivity in July currents were not as profound as expected. This is likely because the eddy fields were largely offshore from the sponge populations. The most significant difference observed was that the surface currents showed strong offshore dispersal in many of the individual years that was not seen in models using the long-term averaged currents, indicating that loss of larvae to the system through entrainment into the Gulf Stream may be even greater than indicated from the models run with long-term averaged currents. As noted, a few novel connections and retentions were revealed in the models using shorter averaging periods, however, those did not change the general pattern of connectivity in the sub-regional populations north and south of the Gulf Stream and in the area between Cape Hatteras and Cape Cod.

Our results showed differences between the dispersal kernels produced with BNAM and GLORYS12VI, as was expected, however, both models support the main conclusion from our analyses, that documented populations of V. pourtalesii are not connected between areas south of Cape Hatteras (Groups 1–4), and to its north [between the Scotian Shelf populations (Groups 6 and 7) and those on the continental slope north of Cape Hatteras (Group 5)]. This result was produced despite the use of different averaging periods and ocean models, and is likely due to the overwhelming influence of the Gulf Stream, allowing for little opportunity for particles to escape from its offshore trajectory, once entrained. Opportunities for entrainment increased with shorter averaging periods used in the models and with the application of the diffusivity constant. Thus while the sub-regional connectivity and retention patterns varied, the main conclusions regarding the connectivity between northern and southern populations, the downstream connectivity of populations within those regions, and the importance of larval retention appear robust.

Effects of Larval Vertical Distribution and Pelagic Duration

The depth of particle release had substantial effects on modeled movements, particularly when horizontal diffusion was incorporated. Releases at the surface and mid-depths (other than those in Group 4) were transported further and faster, whilst particles released at the seabed were more likely to be retained within their source Groups. Similar results have been obtained in other particle-tracking studies and show that the extent of upward movement from the seafloor is a key determinant in dispersal of larvae of benthic species away from their source sites (Edwards et al., 2007; Young et al., 2012; Kool et al., 2015; McVeigh et al., 2017). Most recently, Gary et al. (2020) investigated the influence of vertical swimming behaviors on the dispersal of deep-sea larvae and found that vertical distribution, especially residence time near the surface, had a large impact on theoretical dispersal.

With only limited available information on larval growth and movement in V. pourtalesii, we did not attempt to further integrate their swimming behavior into the 3-D models, nor to examine the effect of vertical position on the advective pathways directly. Thus, potential effects of larval swimming behavior can only be estimated through comparison of the movements of particles released at various depths. There were clear differences amongst those, notably with releases in the southern portion of Group 4, but the connectivity matrices show a concordance in broad patterns across depths, especially with non-zero horizontal diffusivity included in the model. Thus, the conclusion that there is some retention of larvae, along with unidirectional connections in the south (Groups 1–4) and from Group 7 to Group 6 (Figure 4), may be robust to assumptions about larval vertical distribution.

The more-critical unknown is the PLD. If that of V. pourtalesii is short by the standards of other sponges, perhaps no more than a few minutes (Maldonado, 2006; Kenchington et al., 2019), then almost all larval settlement will be in the immediate area of the adults. Longer PLDs would mean greater connectivity, up to the 2-week planktonic drifts modeled here.

Predicted Distributions and Potential for Larval Supply

Beazley et al. (2021) modeled the distribution of V. pourtalesii, based on the predicted occurrence probability from random forest and generalized additive models, under present and future climate conditions, using Representative Concentration Pathway intermediate (RCP 4.5) and worst-case (RCP 8.5) CO₂ concentration scenarios for the period 2046–2065. Their present-day species distribution models showed a discontinuous distribution of V. pourtalesii with two centers of high probability of occurrence, respectively, the Scotian Shelf and Northeast Channel, and the Blake Plateau in the south, with a low probability of occurrence in between. In random forest models both surface and bottom currents were influential in determining the predicted distributions in the south (Beazley et al., 2021). The present study has shown that there is no direct connection between those centers through larval transport, supporting the species distribution model conclusions.

Our particle tracking suggests that the center of occurrence on the Blake Plateau noted by Beazley et al. (2021) may be supported by recruitment of larvae released into the Gulf Stream further south, from sites as far away as the Straits of Florida (either directly or through settlement at intermediate locations), and aided by retention in the gyre south of the Charleston Bump. V. pourtalesii is not found to the west of the known faunal boundary between the Gulf of Mexico and the Atlantic, which is marked by genetic discontinuities even in species with trans-boundary distributions, indicating contemporary obstructions to gene flow (Avise, 1992). With the powerful flow of the Gulf Stream and limited suitable habitat to the westward, it is unlikely that V. pourtalesii could extend its distribution beyond the reef-edge off the Florida Keys. Our particle tracking modeling found that retentive settlement there, in Group 1, only occurred with particles released at the seabed and only if the model included horizontal diffusivity representing, in this case, the presence of local eddies. Cyclonic eddies have been shown to be very important to recruitment of larval damselfish on mesophotic coral reef ecosystems in the Florida Keys area where Group 1 is located (Vaz et al., 2016) and could also influence sponge dispersal at this southern range boundary.

In the north, the situation is very different. During the Pleistocene, the Scotian Shelf was covered in ice until 20 ka BP. By 12 ka BP, an extensive coastal plain with large islands and shallow seas had formed (Pielou, 1991; Shaw et al., 2002, 2006; Paul and Schäfer-Neth, 2003; Schäfer-Neth and Paul, 2003). Present conditions have developed gradually through the millennia since, presumably allowing V. pourtalesii to colonize the Northeast Channel and Emerald Basin through post-Pleistocene range expansion as hypothesized for other species (Kenchington et al., 2009).

Beazley et al. (2021) predicted that the area of suitable habitat for the species will enlarge in the future in both the south and north under the RCP scenarios tested, with the greatest increases in the deeper parts of the Gulf of Maine and Gulf of St. Lawrence. However, the latter is unlikely to be colonized quickly as it lies upstream of contemporary occurrences of V. pourtalesii and projected climate change is not expected to alter the current directionality in this sub-region (Brickman et al., 2018), leaving distributional expansion dependent on fortuitous flow reversals during extreme events. Conversely, our study suggests that the sponge populations in Emerald Basin and the Northeast Channel may already benefit from a moderate degree of larval retention (enhanced in the models with shorter averaging periods), possibly supplemented by a unidirectional connection from the Basin to the Channel. From the latter, future suitable habitat in the Gulf of Maine should be colonized readily, provided that the existing sponge grounds can be protected from anthropogenic harm.

Conclusion

This study used 3-D biophysical particle dispersal models to simulate hydrodynamic connectivity among known concentrations of the deep-sea hexactinellid sponge Vazella pourtalesii throughout its distribution on the east coast of North America. We have shown that intense ocean currents, primarily the Gulf Stream, shape the contemporary distribution of this species, creating two sub-regions roughly separated at Cape Hatteras, United States, and concurs with present day and future predicted distribution models for this species. The addition of connectivity modeling to complement published distribution models informs both modeling approaches for this species. We found no exchange with a 2-week time period (estimated maximum pelagic larval duration for V. pourtalesii) between areas north and south of Cape Hatteras, at surface, mid-water and seabed release depths, irrespective of month of release, application of a horizontal diffusion constant specific to cross-Gulf Stream diffusivity, or use of different averaging periods for ocean currents. This coincides with a gap in suitable habitat in the area north of Cape Hatteras and south of Cape Cod in both present and future predicted distributions (Beazley et al., 2021). Monthly averaged models for selected years suggest that it may be possible that eddy transport of larvae from the Gulf Stream could supply V. pourtalesii larvae from the south to the area northwest of Cape Hatteras (Group 5), connecting the two sub-regions. Such eddies are not fully resolved in our models and smaller eddies could be the means of larval transport from the south to Group 5. However, the viability of such larvae would require a longer PLD than 2 weeks, which is already thought to be a maximum duration from what is known of sponge larval biology. Furthermore, our models do not include larval mortality which would act to also limit connectivity. At sub-regional scales, the influence of the Charleston Bump, a topographic high rising from the Blake Plateau, southeastern United States, was profound. Particles seeded north of the Bump were transported greater distances than those seeded to the south, which were caught in an associated gyre, promoting retention in seabed releases. Connectivity among populations in the southern sub-region was shown to be highly probable, with the Blake Plateau south of the Bump an area with considerable redundancy from downstream populations, potentially conferring increased resilience through genetic mixing from multiple source populations. This area is also likely a sink, with larvae either being retained or transported by the Gulf Stream into deep waters. Despite weaker currents in the northern sub-region, with higher particle retention, our dynamic models show downstream connectivity between the Emerald Basin and Northeast Channel areas. Population genetic studies, focusing on the two sub-regions, could help to establish whether these patterns also represent genetic subdivisions. Our connectivity modeling could be used in combination with the species distribution models (Beazley et al., 2021) to construct explicit spatial hypotheses for population genetic studies of V. pourtalesii following a seascape genetic framework (Kenchington et al., 2006, 2015; Zeng et al., 2020). Long-term conservation of this ecologically important species should focus on the protection of local populations, given that retention appears to be a prominent factor in our dispersal models.

Data Availability Statement

The 3-D Lagrangian particle simulator (Parcels) is freely available at https://www.oceanparcels.org. Monthly mean currents from the Bedford Institute of Oceanography North Atlantic Model (BNAM) results averaged over the 1990–2015 period for the Northwest Atlantic Ocean are available online on the Federal Government of Canada Open Data Portal (https://open.canada.ca/data/en/dataset/dd5ef0d6-d588-4f4f-99e5-27c4bdfddf6a). The GLORYS12V1 ocean products with daily means for the surface are publically available from the Copernicus Marine Service data portal (see text footnote 2).

Author Contributions

EK and AD acquired the funding. SW performed the data analyses and produced the tables and figures. All authors contributed to the writing of the manuscript, take responsibility for the content, approved the submitted version, and designed the experiments.

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.

Funding

This research was funded through the Fisheries and Oceans Canada (DFO) Strategic Program for Ecosystem-Based Research and Advice (SPERA) project “Evaluation of the Effectiveness of Two Sponge Conservation Areas in the Maritimes Region: Identifying Patterns of Dispersal, Connectivity, and Recovery Potential of the Russian Hat Sponge Vazella pourtalesii,” DFO’s International Governance Strategy Science Program through project “Marine Biological Diversity Beyond Areas of National Jurisdiction (BBNJ): 3-Tiers of Diversity (Genes-Species-Communities),” and the H2020 EU Framework Programme for Research and Innovation Project SponGES (Deep-sea Sponge Grounds Ecosystems of the North Atlantic: an integrated approach toward their preservation and sustainable exploitation) (Grant Agreement no. 679849). This document reflects only the authors’ views, and the Executive Agency for Small and Medium-sized Enterprises (EASME) is not responsible for any use that may be made of the information it contains.

Acknowledgments

We thank Lindsay Beazley and Cam Lirette (DFO) for their support in the early stages of this work. Drs. T. J. Kenchington and F. J. Murillo provided thoughtful reviews as part of the DFO internal review process, which helped to improve the manuscript.

Supplementary Material

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

Footnotes

1. ^ https://www.gebco.net/data_and_products/gridded_bathymetry_data/gebco_ 2019/gebco_2019_info.html; accessed 4/12/2021.
2. ^ https://resources.marine.copernicus.eu/?option=com_csw&view=details&product_id=GLOBAL_REANALYSIS_PHY_001_030

References

Anderson, L. A., Robinson, A. R., and Lozano, C. J. (2000). Physical and biological modeling in the Gulf Stream region: I. Data assimilation methodology. Deep Sea Res. Part I Oceanogr. Res. Pap. 47, 1787–1827. doi: 10.1016/S0967-0637(00)00019-4

CrossRef Full Text | Google Scholar

Ashford, J., La Mesa, M., Fach, B. A., Jones, C., and Everson, I. (2010). Testing early life connectivity using otolith chemistry and particle-tracking simulations. Can. J. Fish. Aquat. Sci. 67, 1303–1315. doi: 10.1139/F10-065

PubMed Abstract | CrossRef Full Text | Google Scholar

Avise, J. C. (1992). Molecular population structure and the biogeographical history of a regional fauna: a case history with lessons for conservation biology. Oikos 63, 62–76. doi: 10.2307/3545516

CrossRef Full Text | Google Scholar

Bart, M. C., De Kluijver, A., Hoetjes, S., Absalah, S., Mueller, B., Kenchington, E., et al. (2020). Differential processing of dissolved and particulate organic matter by deep-sea sponges and their microbial symbionts. Sci. Rep. 10:17515. doi: 10.1038/s41598-020-74670-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Bartsch, J., and Knust, R. (1994). Simulating the dispersion of vertically migrating sprat larvae (Sprattus sprattus (L.)) in the German Bight with a circulation and transport model system. Fish. Oceanogr. 3, 92–105. doi: 10.1111/j.1365-2419.1994.tb00052.x

CrossRef Full Text | Google Scholar

Beazley, L., Kenchington, E., Murillo, F. J., Brickman, D., Wang, Z., Davies, A. J., et al. (2021). Climate change winner in the deep sea? Predicting the impacts of climate change on the distribution of the glass sponge Vazella pourtalesii. Mar. Ecol. Prog. Ser. 657, 1–23. doi: 10.3354/meps13566

CrossRef Full Text | Google Scholar

Beazley, L., Wang, Z., Kenchington, E., Yashayaev, I., Rapp, H. T., Xavier, J. R., et al. (2018). Predicted distribution of the glass sponge Vazella pourtalesi on the Scotian Shelf and its persistence in the face of climatic variability. PLoS One 13:e0205505. doi: 10.1371/journal.pone.0205505

PubMed Abstract | CrossRef Full Text | Google Scholar

Bell, J. J. (2008). The functional roles of marine sponges. Estuar. Coast. Shelf Sci. 79, 341–353. doi: 10.1016/j.ecss.2008.05.002

CrossRef Full Text | Google Scholar

Biló, T. C., and Johns, W. E. (2019). Interior pathways of Labrador Sea Water in the North Atlantic from the Argo perspective. Geophys. Res. Lett. 46, 3340–3348. doi: 10.1029/2018GL081439

CrossRef Full Text | Google Scholar

Bower, A. S., Rossby, H. T., and Lillibridge, J. L. (1985). The Gulf stream—barrier or blender? J. Phys. Oceanogr. 15, 24–32. doi: 10.1175/1520-0485(1985)015<0024:TGSOB>2.0.CO;2

CrossRef Full Text | Google Scholar

Brickman, D., Hebert, D., and Wang, Z. (2018). Mechanism for the recent ocean warming events on the Scotian Shelf of eastern Canada. Cont. Shelf Res. 156, 11–22. doi: 10.1016/j.csr.2018.01.001

CrossRef Full Text | Google Scholar

Brickman, D., Wang, Z., and Detracey, B. (2016). Variability of current streams in Atlantic Canadian waters: a model study. Atmos. Ocean 54, 218–229. doi: 10.1080/07055900.2015.1094026

CrossRef Full Text | Google Scholar

Busch, K., Beazley, L., Kenchington, E., Whoriskey, F., Slaby, B. M., and Hentschel, U. (2020). Microbial diversity of the glass sponge Vazella pourtalesii in response to anthropogenic activities. Conserv. Genet. 21, 1001–1010. doi: 10.1007/s10592-020-01305-2

CrossRef Full Text | Google Scholar

Callies, U., Plüß, A., Kappenberg, J., and Kapitza, H. (2011). Particle tracking in the vicinity of Helgoland, North Sea: a model comparison. Ocean Dyn. 61, 2121–2139. doi: 10.1007/s10236-011-0474-8

CrossRef Full Text | Google Scholar

Codling, E. A., Hill, N. A., Pitchford, J. W., and Simpson, S. D. (2004). Random walk models for the movement and recruitment of reef fish larvae. Mar. Ecol. Prog. Ser. 279, 215–224. doi: 10.3354/meps279215

CrossRef Full Text | Google Scholar

Crowder, L. B., Lyman, S. J., Figueira, W. F., and Priddy, J. (2000). Source-sink population dynamics and the problem of siting marine reserves. Bull. Mar. Sci. 66, 799–820.

Google Scholar

Csanady, G. T., and Hamilton, P. (1988). Circulation of slopewater. Cont. Shelf Res. 8, 565–624. doi: 10.1016/0278-4343(88)90068-4

CrossRef Full Text | Google Scholar

Delandmeter, P., and Van Sebille, E. (2019). The Parcels v2.0 Lagrangian framework: new field interpolation schemes. Geosci. Model Dev. 12, 3571–3584. doi: 10.5194/gmd-12-3571-2019

CrossRef Full Text | Google Scholar

Delhez, E., Damm, J. M. P., de Goede, E., de Kok, J. M., Dumas, F., Jones, J. E., et al. (2004). What can we expect from shelf seas models: the NOMADS2 Project. J. Mar. Syst. 45, 39–53. doi: 10.1016/j.jmarsys.2003.09.003

CrossRef Full Text | Google Scholar

Downie, A.-L., Piechaud, N., Howell, K., Barrio Froján, C., Sacau, M., and Kenny, A. (2021). Reconstructing baselines: use of habitat suitability modelling to predict pre-fishing condition of a vulnerable marine ecosystem. ICES J. Mar. Sci. 0, 1–13. doi: 10.1093/icesjms/fsab154

CrossRef Full Text | Google Scholar

Drinkwater, K., Petrie, B., and Smith, P. (2003). Climate variability on the Scotian Shelf during the 1990s. ICES Mar. Sci. Symp. 219, 40–49.

Google Scholar

Edwards, K., Hare, J., Werner, F., and Seim, H. (2007). Using 2-dimensional dispersal kernels to identify the dominant influences on larval dispersal on continental shelves. Mar. Ecol. Prog. Ser. 352, 77–87. doi: 10.3354/meps07169

CrossRef Full Text | Google Scholar

Fernandez, E., and Lellouche, J. M. (2021). Product User Manual For the Global Ocean Physical Reanalysis product GLOBAL_REANALYSIS_ PHY_001_030. Issue 1.2. EU Copernicus Marine Service – Public Ref: CMEMS-GLO-PUM-001-030.

Google Scholar

Ferry, N., Parent, L., Garric, G., Barnier, B., and Jourdain, N. C. (2010). Mercator global eddy permitting ocean reanalysis GLORYS1V1: description and results. Merc. Ocean Q. Newslett0 36, 15–27.

Google Scholar

Fratantoni, P. S., and Pickart, R. S. (2007). The western North Atlantic shelfbreak current system in summer. J. Phys. Oceanogr. 37, 2509–2533. doi: 10.1175/JPO3123.1

CrossRef Full Text | Google Scholar

Gary, S. F., Fox, A. D., Biastoch, A., Roberts, J. M., and Cunningham, S. A. (2020). Larval behaviour, dispersal and population connectivity in the deep sea. Sci. Rep. 10:10675. doi: 10.1038/s41598-020-67503-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Greene, C. H., and Pershing, A. J. (2003). The flip-side of the North Atlantic Oscillation and modal shifts in slope-water circulation patterns. Limnol. Oceanogr. 48, 39–322. doi: 10.4319/lo.2003.48.1.0319

CrossRef Full Text | Google Scholar

Hakkinen, S., and Rhines, P. B. (2009). Shifting surface currents of the northern North Atlantic Ocean. J. Geophys. Res. Oceans 114:C04005. doi: 10.1029/2008JC004883

CrossRef Full Text | Google Scholar

Häkkinen, S., Rhines, P. B., and Worthen, D. L. (2013). Northern North Atlantic sea surface height and ocean heat content variability. J. Geophys. Res. Oceans 118, 3670–3678. doi: 10.1002/jgrc.20268

CrossRef Full Text | Google Scholar

Han, G., and Loder, J. W. (2003). Three-dimensional seasonal-mean circulation and hydrography on the eastern Scotian Shelf. J. Geophys. Res 108:3136. doi: 10.1029/2002JC001463

CrossRef Full Text | Google Scholar

Han, G., Loder, J. W., and Smith, P. C. (1999). Seasonal-mean hydrography and circulation in the Gulf of St. Lawrence and on the eastern Scotian and southern Newfoundland shelves. J. Phys. Ocean. 29, 1279–1301. doi: 10.1175/1520-0485(1999)029<1279:SMHACI>2.0.CO;2

CrossRef Full Text | Google Scholar

Hanna, C. A. (1984). Initial Settlement of Marine Invertebrate Larvae : The Role of Passive Sinking in a Near-Bottom Turbulent Flow Environment. Doctoral thesis. Dartmouth MA: Massachusetts Institute of Technology. doi: 10.1575/1912/2210

CrossRef Full Text | Google Scholar

Hanz, U., Beazley, L., Kenchington, E., Duineveld, G., Rapp, H. T., and Mienis, F. (2021). Seasonal variability in near-bed environmental conditions in the Vazella pourtalesii glass sponge grounds of the Scotian Shelf. Front. Mar. Sci. 08:2021. doi: 10.3389/fmars.2020.597682

CrossRef Full Text | Google Scholar

Hawkes, N., Korabik, M., Beazley, L., Rapp, H., Xavier, J., and Kenchington, E. (2019). Glass sponge grounds on the Scotian Shelf and their associated biodiversity. Mar. Ecol. Prog. Ser. 614, 91–109. doi: 10.3354/meps12903

CrossRef Full Text | Google Scholar

Hogg, M., Tendal, O., Conway, K., Pomponi, S., van Soest, R., Gutt, J., et al. (2010). Deep-sea Sponge Grounds: Reservoirs of Biodiversity. UNEP-WCMC Biodiversity Series 32. Cambridge: UNEP-WCMC.

Google Scholar

Holloway, P., Miller, J. A., and Gillings, S. (2016). Incorporating movement in species distribution models: how do simulations of dispersal affect the accuracy and uncertainty of projections? Int. J. Geogr. Inf. Sci. 30, 2050–2074. doi: 10.1080/13658816.2016.1158823

CrossRef Full Text | Google Scholar

Hufnagl, M., Payne, M., Lacroix, G., Bolle, L. J., Daewel, U., Dickey-Collas, M., et al. (2017). Variation that can be expected when using particle tracking models in connectivity studies. J. Sea Res. 127, 133–149. doi: 10.1016/j.seares.2017.04.009

CrossRef Full Text | Google Scholar

Hunter, J. R., Craig, P. D., and Phillips, H. E. (1993). On the use of random walk models with spatially variable diffusivity. J. Comput. Phys. 106, 366–376. doi: 10.1016/S0021-9991(83)71114-9

CrossRef Full Text | Google Scholar

Hurlburt, H. E., Brassington, G. B., Drillet, Y., Kamachi, M., Benkiran, M., Bourdallé-Badie, R., et al. (2009). High-resolution global and basin-scale ocean analyses and forecasts. Oceanography 22, 110–127. doi: 10.5670/oceanog.2009.70

CrossRef Full Text | Google Scholar

Jones, G. P., Almany, G. R., Russ, G. R., Sale, P. F., Steneck, R. S., Van Oppen, M. J. H., et al. (2009). Larval retention and connectivity among populations of corals and reef fishes: history, advances and challenges. Coral Reefs 28, 307–325. doi: 10.1007/s00338-009-0469-9

CrossRef Full Text | Google Scholar

Kang, D., and Curchitser, E. N. (2013). Gulf Stream eddy characteristics in a high-resolution ocean model. J. Geophys. Res. Oceans 118, 4474–4487. doi: 10.1002/jgrc.20318

CrossRef Full Text | Google Scholar

Kelly, S. J., Popova, E., Aksenov, Y., Marsh, R., and Yool, A. (2020). They came from the Pacific: how changing Arctic currents could contribute to an ecological regime shift in the Atlantic Ocean. Earths Future 8:e2019EF001394. doi: 10.1029/2019EF001394

CrossRef Full Text | Google Scholar

Kenchington, E. L., Nakashima, B. S., Taggart, C. T., and Hamilton, L. C. (2015). Genetic structure of capelin (Mallotus villosus) in the northwest Atlantic Ocean. PLoS One 10:e0122315. doi: 10.1371/journal.pone.0122315

PubMed Abstract | CrossRef Full Text | Google Scholar

Kenchington, E., Harding, G., Jones, M. W., and Prodöhl, P. A. (2009). Pleistocene glaciation events shape genetic structure across the range of the American lobster, Homarus americanus. Mol. Ecol. 18, 1654–1667. doi: 10.1111/j.1365-294X.2009.04118.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Kenchington, E., Patwary, M. U., Zouros, E., and Bird, C. J. (2006). Genetic differentiation in relation to marine landscape in a broadcast-spawning bivalve mollusc (Placopecten magellanicus). Mol. Ecol. 15, 1781–1796. doi: 10.1111/j.1365-294X.2006.02915.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Kenchington, E., Wang, Z., Lirette, C., Murillo, F. J., Guijarro, J., Yashayaev, I., et al. (2019). Connectivity modelling of areas closed to protect vulnerable marine ecosystems in the northwest Atlantic. Deep Sea Res. Part I Oceanogr. Res. Pap. 143, 85–103. doi: 10.1016/j.dsr.2018.11.007

CrossRef Full Text | Google Scholar

Knudby, A., Kenchington, E., and Murillo, F. J. (2013). Modeling the distribution of Geodia sponges and sponge grounds in the northwest Atlantic. PLoS One 8:e82306. doi: 10.1371/journal.pone.0082306

PubMed Abstract | CrossRef Full Text | Google Scholar

Kool, J. T., Huang, Z., and Nichol, S. L. (2015). Simulated larval connectivity among Australia’s southwest submarine canyons. Mar. Ecol. Progr. Ser. 539, 77–91. doi: 10.3354/meps11477

CrossRef Full Text | Google Scholar

Kool, J. T., Moilanen, A., and Treml, E. A. (2013). Population connectivity: recent advances and new perspectives. Landsc. Ecol. 28, 165–185. doi: 10.1007/s10980-012-9819-z

CrossRef Full Text | Google Scholar

Koutitonsky, V. G., and Bugden, G. L. (1991). The physical oceanography of the Gulf of St. Lawrence: a review with emphasis on the synoptic variability of the motion. Can. Spec. Publ. Fish. Aquat. Sci. 113, 57–90.

Google Scholar

Lange, M., and Van Sebille, E. (2017). Parcels v0.9: prototyping a Lagrangian ocean analysis framework for the petascale age. Geosci. Model Dev. 10, 4175–4186. doi: 10.5194/gmd-10-4175-2017

CrossRef Full Text | Google Scholar

Lanna, E., and Riesgo, A. (2020). Sponge larvae do not swim that fast: a reply to Montgomery et al. (2019). J. Mar. Biol. Assoc. U.K. 100, 181–183. doi: 10.1017/S0025315419000870

CrossRef Full Text | Google Scholar

Levin, J., Wilkin, J., Fleming, N., and Zavala-Garay, J. (2018). Mean circulation of the Mid-Atlantic Bight from a climatological data assimilative model. Ocean Model. 128, 1–14. doi: 10.1016/j.ocemod.2018.05.003

CrossRef Full Text | Google Scholar

Leys, S. P., and Ereskovsky, A. V. (2006). Embryogenesis and larval differentiation in sponges. Can. J. Zool. 84, 262–287. doi: 10.1139/z05-170

PubMed Abstract | CrossRef Full Text | Google Scholar

Loder, J. W., Han, G., Hannah, C. G., Greenberg, D. A., and Smith, P. C. (1997). Hydrography and baroclinic circulation in the Scotian shelf region: winter versus summer. Can. J. Fish. Aquat. Sci. 54(Suppl. 1), 40–56. doi: 10.1139/f96-153

PubMed Abstract | CrossRef Full Text | Google Scholar

Loder, J. W., Petrie, B., and Gawarkiewicz, G. G. (1998). “The coastal ocean off northeastern North America: a large-scale view,” in The Sea, Volume 11: The Global Coastal Ocean - Regional Studies and Syntheses, eds A. R. Robinson and K. Brink (New York, NY: John Wiley & Sons), 105–133.

Google Scholar

Maldonado, M. (2006). The ecology of the sponge larva. Can. J. Zool. 84, 175–194. doi: 10.1139/z05-177

PubMed Abstract | CrossRef Full Text | Google Scholar

Maldonado, M., Aguilar, R., Bannister, R. J., Bell, J. J., Conway, K. W., Dayton, P. K., et al. (2017). “Sponge grounds as key marine habitats: a synthetic review of types, structure, functional roles, and conservation concerns,” in Marine Animal Forests, eds S. Rossi, L. Bramanti, A. Gori, and C. Orejas (Cham: Springer International Publishing), 1–39. doi: 10.1007/978-3-319-21012-4_24

CrossRef Full Text | Google Scholar

Maldonado, M., Beazley, L., Lopez-Acosta, M., Kenchington, E., Casault, B., Hanz, U., et al. (2020). Massive silicate utilization facilitated by a benthic-pelagic coupled feedback sustains deep-sea sponge aggregations. Limnol. Oceanogr. 66, 366–391. doi: 10.1002/lno.11610

CrossRef Full Text | Google Scholar

Marinone, S. G., Ullo, M. J., Parés-Sierra, A., Lavín, M. F., and Cudney-Bueno, R. (2008). Connectivity in the northern Gulf of California from particle tracking in a three-dimensional numerical model. J. Mar. Syst. 71, 149–158. doi: 10.1016/j.jmarsys.2007.06.005

CrossRef Full Text | Google Scholar

Marsh, R., Petrie, B., Weidman, C. R., Dickson, R. R., Loder, J. W., Hannah, C. G., et al. (1999). The 1882 tilefish kill – a cold event in shelf waters off the north-eastern United States? Fish. Oceanog. 8, 39–49. doi: 10.1046/j.1365-2419.1999.00092.x

CrossRef Full Text | Google Scholar

McClain, C. R., and Atkinson, L. P. (1985). A note on the Charleston Gyre. J. Geophys. Res. 90, 11857–11861. doi: 10.1029/JC090iC06p11857

CrossRef Full Text | Google Scholar

McGuire, K. J., and McDonnell, J. J. (2006). A review and evaluation of catchment transit time modeling. J. Hydrol. 330, 543–563. doi: 10.1016/j.jhydrol.2006.04.020

CrossRef Full Text | Google Scholar

McVeigh, D. M., Eggleston, D. B., Todd, A. C., Young, C. M., and He, R. (2017). The influence of larval migration and dispersal depth on potential larval trajectories of a deep-sea bivalve. Deep Sea Res. Part I Oceanogr. Res. Pap. 127, 57–64. doi: 10.1016/j.dsr.2017.08.002

CrossRef Full Text | Google Scholar

MERCINA Working Group. (2001). Oceanographic responses to climate in the northwest Atlantic. Oceanography 14, 76–82. doi: 10.5670/oceanog.2001.25

CrossRef Full Text | Google Scholar

Metaxas, A., and Saunders, M. (2009). Quantifying the “bio-” components in biophysical models of larval transport in marine benthic invertebrates: advances and pitfalls. Biol. Bull. 216, 257–272. doi: 10.1086/BBLv216n3p257

PubMed Abstract | CrossRef Full Text | Google Scholar

Montgomery, E. M., Hamel, J. F., and Mercier, A. (2019). Larval nutritional mode and swimming behaviour in ciliated marine larvae. J. Mar. Biol. Assoc. U.K. 99, 1027–1032. doi: 10.1017/S0025315418001091

CrossRef Full Text | Google Scholar

Montgomery, E. M., Hamel, J. F., and Mercier, A. (2020). Larval nutritional mode and swimming behaviour in ciliated marine larvae – CORRIGENDUM. J. Mar. Biol. Assoc. U.K. 100, 189–191. doi: 10.1017/S0025315419001176

CrossRef Full Text | Google Scholar

Morato, T., González-Irusta, J.-M., Dominguez-Carrió, C., Wei, C.-L., Davies, A., Sweetman, A. K., et al. (2020). Climate-induced changes in the suitable habitat of cold-water corals and commercially important deep-sea fishes in the North Atlantic. Glob. Change Biol. 26, 2181–2202. doi: 10.1111/gcb.14996

PubMed Abstract | CrossRef Full Text | Google Scholar

Murillo, F. J., Kenchington, E., Tompkins, G., Beazley, L., Baker, E., Knudby, A., et al. (2018). Sponge assemblages and predicted archetypes in the eastern Canadian Arctic. Mar. Ecol. Prog. Ser. 597, 115–135. doi: 10.3354/meps12589

CrossRef Full Text | Google Scholar

Ospina-Alvarez, A., de Juan Mohan, S., Alós, J., Basterretxea, G., Alonso-Fernández, A., Follana-Berná, G., et al. (2020). MPA network design based on graph theory and emergent properties of larval dispersal. Mar. Ecol. Prog. Ser. 650, 309–326. doi: 10.3354/meps13399

CrossRef Full Text | Google Scholar

Paris, C. B., Helgers, J., van Sebille, E., and Srinivasan, A. (2013). Connectivity modeling system: A probabilistic modeling tool for the multi-scale tracking of biotic and abiotic variability in the ocean. Environ. Modell. Softw. 42, 47–54. doi: 10.1016/j.envsoft.2012.12.006

CrossRef Full Text | Google Scholar

Paul, A., and Schäfer-Neth, C. (2003). Modeling the water masses of the Atlantic Ocean at the last glacial maximum. Paleoceanogr 18, 1058–1084. doi: 10.1029/2002PA000783

CrossRef Full Text | Google Scholar

Petrie, B., and Drinkwater, K. (1993). Temperature and salinity variability on the Scotian Shelf and in the Gulf of Maine 1945–1990. J. Geophys. Res 98, 20079–20089. doi: 10.1029/93JC02191

CrossRef Full Text | Google Scholar

Pielou, E. C. (1991). After the Ice Age: The Return of Life to Glaciated North America. Chicago, IL: The University of Chicago Press. doi: 10.7208/chicago/9780226668093.001.0001

CrossRef Full Text | Google Scholar

Popenoe, P., and Manheim, F. T. (2001). Origin and history of the Charleston Bump – Geological formations, currents, bottom conditions, and their relationship to wreckfish habitats on the Blake Plateau. Am. Fish. Soc. Symp. 25, 43–94.

Google Scholar

Putman, N. F., and He, R. (2013). Tracking the long-distance dispersal of marine organisms: sensitivity to ocean model resolution. J. R. Soc. Interface 10:20120979. doi: 10.1098/rsif.2012.0979

PubMed Abstract | CrossRef Full Text | Google Scholar

Ramos, J. E., Pecl, G. T., Moltschaniwskyj, N. A., Semmens, J. M., Souza, C. A., and Strugnell, J. M. (2018). Population genetic signatures of a climate change driven marine range extension. Sci. Rep. 8:9558. doi: 10.1038/s41598-018-27351-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Reijnders, B. J. H. R. (2020). Assessing Ocean Surface Connectivity in the Arctic: Capabilities and Caveats of Community Detection in Lagrangian Flow Networks. Master thesis. Utrecht: Utrecht University. doi: 10.5194/egusphere-egu21-201

CrossRef Full Text | Google Scholar

Reiswig, H. M. (1996). Redescription and placement of the rossellid genus Vazella Gray (Hexactinellida : Lyssacinosida). Bull. Instr. Sci. Nat. Belg. 66, 135–141.

Google Scholar

Ross, R. E., Nimmo-Smith, W. A. M., Torres, R., and Howell, K. L. (2020). Comparing deep-sea larval dispersal models: A cautionary tale for ecology and conservation. Front. Mar. Sci. 12:2020. doi: 10.3389/fmars.2020.00431

CrossRef Full Text | Google Scholar

Rossby, T., Bower, A. S., and Shaw, P. T. (1985). Particle pathways in the gulf stream. Bull. Am. Meteorol. Soc. 66, 1106–1110. doi: 10.1175/1520-0477(1985)066<1106:PPITGS>2.0.CO;2

CrossRef Full Text | Google Scholar

Schäfer-Neth, C., and Paul, A. (2003). “The Atlantic ocean at the last glacial maximum: 1. objective mapping of the GLAMAP sea-surface conditions,” in The South Atlantic in the Late Quatenary: Reconstruction of Material Budgets and Current Systems, eds G. Wefer, S. Mulitza, and V. Ratmeyer (Berlin: Springer-Verlag), 531–548. doi: 10.1007/978-3-642-18917-3_23

CrossRef Full Text | Google Scholar

Schmidt, O. (1870). Grundzüge einer Spongien-Fauna des Atlantischen Gebietes. Leipzig: Wilhelm Engelmann.

Google Scholar

Shaw, J., Gareau, P., and Courtney, R. C. (2002). Palaeogeography of Atlantic Canada 13–0 kyr. Q. Sci. Rev. 21, 1861–1878. doi: 10.1016/S0277-3791(02)00004-5

CrossRef Full Text | Google Scholar

Shaw, J., Piper, D. J. W., Fader, G. B. J., King, E. L., Todd, B. J., Bell, T., et al. (2006). A conceptual model ofthe deglaciation of Atlantic Canada. Q. Sci. Rev. 25, 2059–2081. doi: 10.1016/j.quascirev.2006.03.002

CrossRef Full Text | Google Scholar

Simons, R. D., Siegel, D. A., and Brown, K. S. (2013). Model sensitivity and robustness in the estimation of larval transport: a study of particle tracking parameters. J. Marine Syst. 119-120, 19–29. doi: 10.1016/j.jmarsys.2013.03.004

CrossRef Full Text | Google Scholar

Spalding, M. D., Fox, H. E., Allen, G. R., Davidson, N., Ferdaña, Z. A., Finlayson, M., et al. (2007). Marine ecoregions of the world: a bioregionalization of coastal and shelf areas. BioScience 57, 573–583. doi: 10.1641/B570707

CrossRef Full Text | Google Scholar

Swearer, S. E., Caselle, J. E., Lea, D. W., and Warner, R. R. (1999). Larval retention and recruitment in an island population of a coral-reef fish. Nature 402, 799–802. doi: 10.1038/45533

CrossRef Full Text | Google Scholar

Tabachnick, K. (1994). “Distribution of recent hexactinellida,” in Sponges in Time and Space: Biology, Chemistry, Paleontology, eds R. V. K. Van Soest, T. M. G. van Kempen, and J.-C. Braekamn (Rotterdam: Balkema), 225–232.

Google Scholar

Townsend, D. W., Pettigrew, N. R., Thomas, M. A., Neary, M. G., McGillicuddy, D. J. Jr., and O’Connell, J. (2015). Water masses and nutrient sources to the Gulf of Maine. J. Mar. Res. 73, 93–122. doi: 10.1357/002224015815848811

PubMed Abstract | CrossRef Full Text | Google Scholar

Van Soest, R. W. M., Boury-Esnault, N., Vacelet, J., Dohrmann, M., Erpenbeck, D., De Voogd, N. J., et al. (2012). Global diversity of sponges (Porifera). PLoS One 7:e35105. doi: 10.1371/journal.pone.0035105

PubMed Abstract | CrossRef Full Text | Google Scholar

Vaz, A. C., Paris, C. B., Olascoaga, M. J., Kourafalou, V. H., Kang, H., and Reed, J. K. (2016). The perfect storm: match-mismatch of bio-physical events drives larval reef fish connectivity between Pulley Ridge mesophotic reef and the Florida Keys. Cont. Shelf Res. 125, 136–146. doi: 10.1016/j.csr.2016.06.012

CrossRef Full Text | Google Scholar

Wang, S., Kenchington, E. L., Wang, Z., Yashayaev, I., and Davies, A. J. (2020). 3-D Ocean particle tracking modeling reveals extensive vertical movement and downstream interdependence of closed areas in the northwest Atlantic. Sci. Rep. 10:21421. doi: 10.1038/s41598-020-76617-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, Z., Brickman, D., and Greenan, B. J. W. (2019). Characteristic evolution of the Atlantic Meridional Overturning circulation from 1990 to 2015: An eddy-resolving ocean model study. Deep Sea Res. Part I Oceanogr. Res. Pap. 149:103056. doi: 10.1016/j.dsr.2019.06.002

CrossRef Full Text | Google Scholar

Wang, Z., Brickman, D., Greenan, B. J. W., and Yashayaev, I. (2016). An abrupt shift in the Labrador current system in relation to winter NAO events. J. Geophys. Res. Oceans 121, 5338–5349. doi: 10.1002/2016JC011721

CrossRef Full Text | Google Scholar

Wang, Z., Lu, Y., Greenan, B., Brickman, D., and DeTracey, B. (2018). BNAM: an eddy-resolving North Atlantic ocean model to support ocean monitoring. Can. Tech. Rep. Hydrogr. Ocean. Sci. 327, vii+18.

Google Scholar

Wolanski, E., Doherty, P., and Carleton, J. (1997). Directional swimming of fish larvae determines connectivity of fish populations on the Great Barrier Reef. Nat. Sci. 84, 262–268. doi: 10.1007/s001140050394

CrossRef Full Text | Google Scholar

Yashayaev, I., and Loder, J. W. (2009). Enhanced production of Labrador Sea Water in 2008. Geophys. Res. Lett. 36:L01606. doi: 10.1029/2008GL036162

CrossRef Full Text | Google Scholar

Young, C. M., He, R., Emlet, R. B., Li, Y., Qian, H., Arellano, S. M., et al. (2012). Dispersal of deep-sea larvae from the Intra-American Seas: simulations of trajectories using ocean models. Integr. Comp. Biol. 52, 483–496. doi: 10.1093/icb/ics090

PubMed Abstract | CrossRef Full Text | Google Scholar

Zeng, C., Rowden, A. A., Clark, M. R., and Gardner, P. A. (2020). Species-specific genetic variation in response to deep-sea environmental variation amongst vulnerable marine ecosystem indicator taxa. Sci. Rep. 10:2844. doi: 10.1038/s41598-020-59210-0

PubMed Abstract | CrossRef Full Text | Google Scholar