# Highly structured populations of copepods at risk to deep-sea mining: Integration of genomic data with demogenetic and biophysical modelling

Coral Diaz-recio Lorenzo, Adrien Tran Lu Y, Otis Brunner, Pedro Martínez

Arbizu, Didier Jollivet, Stefan Laurent, Sabine Gollner

► To cite this version:

Coral Diaz-recio Lorenzo, Adrien Tran Lu Y, Otis Brunner, Pedro Martínez Arbizu, Didier Jollivet, et al.. Highly structured populations of copepods at risk to deep-sea mining: Integration of genomic data with demogenetic and biophysical modelling. Molecular Ecology, 2024, 33 (9), pp.e17340. 10.1111/mec.17340. hal-04758530

HAL Id: hal-04758530

https://hal.science/hal-04758530v1

Submitted on 29 Oct 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
1

**Title** Highly structured populations of copepods at risk to deep-sea mining: integration

2

of genomic data with demogenetic and biophysical modelling

3

4

**Authors**

5

Coral Diaz-Recio Lorenzo¹*, Adrien Tran Lu Y², Otis Brunner³, Pedro Martínez Arbizu⁴,

6

Didier Jollivet¹, Stefan Laurent⁵, Sabine Gollner⁶,⁷

7

8

**Affiliations**

9

¹ Adaptation et Diversité en Milieu Marin (AD2M), Station Biologique de Roscoff,

10

Sorbonne Université, CNRS, Roscoff, 29680, France

11

² UMR MARBEC, University of Montpellier, IRD, Ifremer, CNRS, Sète, France.

12

³ Okinawa Institute for Science and Technology, 1919-1 Tancha, Onna-son, Kunigami-

13

gun, Okinawa, Japan.

14

⁴ Senckenberg am Meer, German Centre for Marine Biodiversity Research, Südstrand 44,

15

26382 Wilhelmshaven, Germany.

16

⁵ BioNTech, Mainz, Germany

17

⁶ NIOZ Royal Netherlands Institute for Sea Research and Utrecht University. Landsdiep 4,

18

1797 SZ, 't Horntje (Texel), The Netherlands.

19

⁷ Utrecht University, Budapestlaan 4, 3584 CD, Utrecht, The Netherlands.

20

21

**\*Corresponding author:**

22

Email: [cdrrlc@gmail.com](mailto:cdrrlc@gmail.com)

23

Address: Adaptation et Diversité en Milieu Marin (AD2M), Station Biologique de Roscoff,

24

Sorbonne Université, CNRS, Roscoff, 29680, France
25

**Abstract**

26

27

Copepoda is the most abundant taxon in deep-sea hydrothermal vents, where hard substrate is available. Despite the increasing interest for seafloor massive sulphides exploitation, there have been no population genomic studies conducted on vent meiofauna, which are known to contribute over 50% to metazoan biodiversity at vents.

28

29

30

31

To bridge this knowledge gap, restriction site-associated DNA sequencing, specifically 2b-RADseq, was used to retrieve thousands of genome-wide single nucleotide polymorphisms (SNPs) from abundant populations of the vent-obligate copepod *Stygiopontius lauensis* from the Lau Basin. SNPs were used to investigate population structure, demographic histories, genotype-environment associations at a basin scale.

32

33

34

35

36

Genetic analyses also helped to evaluate the suitability of tailored larval dispersal models and the parameterization of life history traits that better fit the population patterns observed in the genomic dataset for the target organism. Highly structured populations were observed on both spatial and temporal scales, with divergence of populations between the north, mid and south of the basin estimated to have occurred after the creation of the major transform fault dividing the Australian and the Niuafo'ou tectonic plate (350 kya), with relatively recent secondary contact events (< 20 kya). Larval dispersal models were able to predict the high levels of structure and the highly asymmetric northward low-level gene flow observed in the genomic data. These results differ from most studies conducted on megafauna in the region, elucidating the need to incorporate smaller size when considering site-prospecting for deep-sea exploitation of seafloor massive sulphides, and the creation of area-based management tools to protect areas at risk of local extinction, should mining occur.

37

38

39

40

41

42

43

44

45

46

47

48
**Keywords:** Hydrothermal vents, copepods, connectivity, demography, larval dispersal modelling, deep-sea mining

51

52 **Author Contributions:** Coral Diaz-Recio Lorenzo: design of the project, DNA extractions, library preparation, data analyses, writing of draft manuscript; Adrien Tran

53 **Lu Y:** demographic modelling, contribution to draft manuscript; Otis Brunner: larval dispersal modelling, contribution to draft manuscript, Stefan Laurent: supervision of population genomic analysis, revision of draft manuscript, Pedro Martinez Arbizu: providing resources (training and consumables) for library construction, revision of draft manuscript, Didier Jollivet: providing of samples, revision of draft manuscript, Sabine Gollner: design of the project, revision of draft manuscript.

60

61

## 1. Introduction

63 Genetic connectivity is essential in maintaining diversity and promoting the evolution of marine organisms, and depends on factors such as dispersal ability, reproductive traits, seafloor topology, and environmental conditions (De Wit et al., 2023; Blanchard and Gollner, 2022; Virtanen et al., 2020; Brunner et al., 2023; Pereira et al., 2023; Portanier et al., 2023; Yahagi et al., 2017; Ingels et al., 2023; Mouchi et al., 2023; Plouviez et al., 2013; 67 2019; Faria et al., 2021; Antich et al., 2022). High dispersal ability promotes gene flow over long distances, while seafloor topology can create physical barriers that impede movement and lead to allopatric speciation (Yahagi et al., 2017; Plouviez et al., 2013; 69 70 2019; Ingels et al., 2023). Additionally, environmental factors such as temperature, salinity, and oxygen concentrations influence species distribution, potentially resulting in 72
local adaptation, and ultimately sympatric speciation (Virtanen et al., 2020; Faria et al., 2021; Antich et al., 2022).

Deep-sea hydrothermal vents, found along Mid Ocean Ridges (MORs) and Back-Arc Basins (BABs), support diverse meta-communities of mega-, macro-, and meiofauna (Mullineaux et al., 2018) that exhibit deeply entwined interactions across organismal groups. The community composition, succession, functional importance, and connectivity of hydrothermal vent mega- and macrofauna has been explored since hydrothermal vents were discovered in 1977. Recently, hydrothermal vent meiofauna has gained interest due to their potentially important role in the food web (limen et al., 2007; 2008; Nomaki et al., 2023), broad horizontal range of vent-endemic genera (Gollner et al., 2011; 2016) and the ability of non-vent-endemic genera to adapt to an unprecedented range of hydrothermal vent conditions (Nomaki et al., 2019; Nakasugi et al., 2021). Understanding the genetic connectivity of smaller size classes in these habitats is crucial due to their ephemeral nature and potential to create isolated populations at risk of local extinctions despite observed gene flow among populations along fast and slow spreading MORs (Craddock et al., 1995; Jollivet et al., 1995; Won et al., 2003; Hurtado et al., 2004; Mullineaux et al., 2010; Bweedessee et al., 2013; Gollner et al., 2011; 2016; Breusing 2015, 2016; Teixeira et al., 2011, 2012; Breusing et al., 2016; Yahagi et al., 2019; Yearsley et al., 2020; Perez et al., 2021). The life history strategies of organisms in these environments are adapted to volcanic eruptions, exhibiting R-strategist survival mechanisms characterized by rapid growth, early maturation, high fecundity, and dispersal capabilities (Van Dover et al., 2014). Recent studies have explored the demographic history of vent megafauna populations in BABs, revealing contrasting dispersal strategies (Thaler et al., 2014;
97 Plouviez et al., 2019; Breusing et al., 2021; Tran Lu Y et al., 2021; Poitrimol et al., 2022).

98 However, there is currently no research on the genomic connectivity of the meiofauna
99 (animals between 32 µm and 1mm in size), in particular the Dirivultidae family of
100 copepods, which thrives on hard substrates and hydrothermal vents (Figure 1, inset).

101
102 The Dirivultidae (Siphonostomatoida) exist in the tens to hundreds of thousands of
103 individuals per vent wherever basalt rock, hard-shelled megafauna assemblages, or
104 patches of dead organisms creating landscapes of exoskeletons exist (Gollner et al.
105 2010;2016; Diaz-Recio Lorenzo et al., 2021). The Dirivultidae (Siphonostomatoida)
106 contains at least 80 described species that have adapted to vent environments,
107 dominating areas of high temperatures with high concentrations of hydrogen sulphide
108 and low oxygen levels where they feed predominantly on chemoautotrophic bacteria
109 (Gollner et al., 2016, Nomaki et al., 2019; 2023; Watanabe et al., 2021; Diaz-Recio Lorenzo
110 et al., 2021; Nakasugi et al., 2021). They have a lecithotrophic developmental mode and
111 likely disperse via hydrothermal plumes, settling as vent-obligate adults. Species of this
112 family exhibit large population sizes, high levels of mitochondrial gene diversity, and
113 demographic analysis of the latter has shown that the vast majority of Dirivultidae species
114 are undergoing population expansion (Gollner et al., 2016; Watanabe et al., 2021, Diaz-
115 Recio Lorenzo et al., 2023).

116
117 During the last decade, Next-Generation Sequencing technology significantly increased
118 the possibility of obtaining large genomic datasets. However, genome-scale population
119 genetic studies have focused mostly on flagship systems and organisms. Restriction-site
120 associated DNA sequencing (RAD-seq) has been widely applied as an alternative to obtain
meaningful genome-wide information in the form of Single Nucleotide Polymorphisms (SNPs) and in particular, 2b-RAD (Wang et al., 2012) has the ability to target a large number of markers throughout the genome of non-model organisms. This information has now widely been applied to non-model organisms to address question about contemporary population structure as well as how populations separated and reconnected through time, by modeling the effect of migration and population growth on genetic divergence (Gutenkunst et al., 2009; Excoffier et al., 2013; Rougeux et al., 2017; De Jode et al., 2022).

The mining of polymetallic seafloor massive sulphides (SMS) poses the most significant threat to hydrothermal vent communities, resulting in short- and long-term damage, including changes in vent fluid strength and chemistry, toxic plumes, habitat removal, fragmentation, and biodiversity loss (Blasiak et al., 2020; Drazen et al., 2020; Gollner et al., 2017, 2021; Levin et al., 2020; Miller et al., 2018; Van Dover et al., 2018; Blanchard and Gollner, 2022). The impacts of this damage may extend beyond the vent realm (Blanchard and Gollner et al., 2022). Back-arc basins in the western Pacific contain high-grade SMS deposits, with historical exploration tenement blocks issued by some island states (Petersen et al., 2016; Du Preez and Fisher, 2018). However, certain Pacific Island states have imposed a 10-year moratorium on mining activities to understand the environmental and societal impacts better (Kakee, 2020). Advancing our understanding of population connectivity across species size fractions, incorporating life-history traits, evolutionary adaptations, larval dispersal modeling, and demographic modeling, can contribute to defining the vent transition zone and determining the size and location of Area-Based Management Tools (ABMTs) (Blachard and Gollner, 2022).
145

146

147

148

149

150

151

152

153

154

155

156

157

158

159

160

161

162

163

## 2. Materials and Methods

164

### 2.1 *Sampling*

165

166

167

168

The resilience of hydrothermal vent meta-communities depends on the dispersal and colonization abilities of organisms, which play a crucial role in their ability to adapt and survive (Mullineaux et al., 2018). Previous studies have focused on modeling the dispersal ability of a subset of vent fauna to identify isolated populations at risk of local extinction from human disturbances (Mitarai et al., 2016; Breusing et al., 2021). However, these models have primarily considered megafaunal dispersal behavior, which differs inherently from the meiofauna due to variations in energy costs associated with their size differences (Gollner et al., 2015). This study aims to investigate the demographic connectivity of the copepod *Stygiopontius lauensis*, the most abundant deep-sea hydrothermal vent animal in the Lau basin. By utilizing genomic markers across four vent sites, we examine contemporary population structure, the relationship between genomic variation and environmental factors, and the demographic history between populations. Additionally, we assess the effectiveness of tailored larval dispersal models in predicting connectivity between nearby and distant vents within the Lau basin. This multifaceted approach allows for the identification of source and sink populations, and those most at risk to the impacts of deep-sea mining.
During expedition FK160407 samples were collected using the remotely operated vehicle (ROV) ROPOS, using scoops and slurp guns from *Alviniconcha* and *Ifremeria* megafaunal assemblages at the vent fields Tui Malila, ABE, and Tahi Moana. During the Chubacarc expedition aboard L'Atalante, samples were obtained using the ROV Victor 6000 using slurp guns from *Alvinichoncha* megafaunal assemblages at the Tu'i Malila, ABE, and Mangatolo vent sites (Figure 1, Table 1). Megafaunal collections were processed on board of the ship as followed: Copepods were washed off from the megafauna and retained on a 32µm sieve. Copepods collected in 2016 were fixed in 96 % Ethanol, samples collected in 2019 were fixed in 99.9 % EtOH. Copepods were identified first to family level (Dirivultidae), and then identified morphologically to species using taxonomic keys described by Gollner et al. (2010). The most abundant species, namely *Stygiopontius lauensis* was chosen for this study.

181

## 2.2 Library preparation and sequencing

All samples were given three consecutive ethanol baths to remove potential external contaminants before DNA extraction. Extractions from whole specimens were carried out using the E.Z.N.A.® Mollusc DNA Kit (omega BIO-TEK). 2µl of DNA extract were used to amplify the cytochrome C oxidase, subunit I (COI) gene, using the universal Folmer primers (Folmer et al., 1994) LCO1490 (5'GGTCAACAAATCATAAAGATATTGG'3), CO2198 (5'TAAACTTCAGGGTGACCAAAAAATCA'3), M13 tails M13F-pUC (-40) (5'GTTTTCCCAGTCACGAC'3), M13R-pUC (-40) (5'CAGGAAACAGCTATGAC'3). PCR cycle conditions were initial denaturation for 3 minutes at 94°C followed by 40 cycles of denaturation at 30 seconds for 94°C, annealing for 45 seconds at 45°C and extension for 45 seconds for 72°C and a final extension for 2 minutes at 72°C. Samples were sequenced
193 with Macrogen Europe. *Cox1* was chosen as a species delimiter due to its conserved nature
194 between species, and specimens that did not match at least 98% with *S. lauensis* when
195 blasting COI sequences against the NCBI database were not included in the library. Out of
196 168 specimens barcoded, 149 matched at least 98% with *S. lauensis* and were used to
197 construct the 2b-RAD libraries.

198

199 The library was constructed using a custom index array, and the enzyme *bgcl*. DNA
200 extraction and library preparation was conducted at the Senckenberg Institute/German
201 Center for Marine Biodiversity Research (DZMB). 2b-RAD was chosen due to the ability of
202 type IIb enzymes to find a large number of restriction sites across the genome that are not
203 error prone and the type IIb restriction enzyme *bgcl* was used due its simplistic
204 laboratory usage. DNA was digested, shearing DNA upstream and downstream of a
205 recognition site, resulting in fragments of a uniform size of 32 bp (Wang et al., 2012).
206 Population-specific and individual-specific barcodes were ligated to the fragments in
207 order to multiplex a large number of individuals, reducing the number of runs and
208 subsequent costs. Adapter-ligated fragments were then amplified using the Fusion2X
209 master mix, with 2X 20 cycles until 2 clear bands were visible, one containing non-ligated
210 adapter, and a second, longer (by 32 bp) fragment containing the RAD sequence with a
211 total fragment length of 182 bp. The longer fragments were cut from the gel and purified
212 using the Monarch® Genomic DNA Purification Kit (New England Biolabs). Paired-end
213 150-bp sequencing (400 million clusters) was conducted on an Illumina Nextseq500
214 platform at the University Medical Centre of Utrecht University (UMC) on samples from
215 2016. After processing samples from 2019, paired-end 150-bp sequencing (1 billion
216

clusters) was conducted on an Illumina Nextseq2000 platform at UMC. Samples were
demultiplexed by UMC.

217

218

219

## 2.3 De novo 2b-RAD-tag assembly, filtering, and SNP-calling

220

The *STACKS2.6* software module *process_radtags* (Rochette et al., 2019) was used to further demultiplex samples by the adapters used to delimit the populations, remove inline barcodes, and remove PCR duplicates. All reads had a sequence length of 32bp, a characteristic of 2b-RAD sequencing using the enzyme *bcgl*, that removes bias introduced from having to truncate reads to a uniform length. A total of 149 samples were successfully sequenced, with a mean coverage of 25 X. Samples with read depth of <5 and >80 were removed. Initial exploration of the data was conducted to observe outliers, and overall structure. This was done by comparing classic multivariate analyses that visualize data under reduced dimensions such as Principal Component Analysis (PCA), non-metric Multidimensional Scaling (nMDS) and Uniform Manifold Approximation and projection (UMAP), combined with spectral clustering, to visualize population structure.

221

222

223

224

225

226

227

228

229

230

231

232

*S. lauensis* lives in vent environments containing a high abundance of chemoautotrophic bacteria that they also feed on (Limen et al. 2007;2008; Nomaki et al. 2019). Due to the small size of copepods, the entire specimen was used during the extraction to ensure sufficient quantities of pure DNA for digestion and adapter ligation during library preparation. It is therefore likely that some of the RAD tags retrieved may not belong to the copepod genome. RADs created with 2b-RAD are too small to be blasted against existing databases as a quality control method to eradicate non-copepod DNA. Although still under debate, significant correlations have been found between GC content and

233

234

235

236

237

238

239
growth temperature in thermophilic prokaryotes (Hu et al. 2022). GC pairs have an additional hydrogen bond than AT pairs, and GC-rich genomes are expected to provide additional stability to the DNA double-helix under high temperatures. For this reason, RADs belonging to those samples were filtered for GC content higher than 50 % as an extra quality control step.

We used the SNP calling software *DiscoSnpRAD* (Gauthier et al. 2020) for SNP discovery and genotyping for RAD-seq data, which uses a k-mer based approach, implementing the construction of de Bruijn graphs to find bubbles of SNPs without the loss of too much data and without the need for a reference genome. In *DiscoSnpRAD*, the variable parameter is the length of the k-mer. This was tested with lengths between 5 and 31. An optimal k of 15 was chosen, as this resulted in a Variant Calling File (VCF) with the highest number of called SNPs. SNPs were then clustered by loci using the provided post-processing scripts from *DiscoSnpRAD*. Missing data is a common characteristic of RADseq data however studies have highlighted that removing this data could lead to loss of valuable information (Chattopadhyay et al., 2014; Eaton et al. 2016; Tripp et al. 2017; Crotti et al. 2018). Such studies have also reported that setting too conservative parameters leads to detrimental loss of information and that *de novo* SNP calling should be less stringent. Missing data cutoffs were therefore set to 0.5. This method was compared to the classic SNP calling STACKS2 software and results are available in supplementary information (S1, figures 1 and 2).
264 **2.4 Outlier detection for structure analyses**

265 Loci under selection are expected to behave differently from neutral loci in terms of
266 population-related patterns of diversification and are thus typically removed from the
267 global dataset in population genomics studies for investigation of geographic structure
268 (Beaumont and Nichols, 1996; Luikart et al. 2003; Holderegger and Wagner, 2006; Foll
269 and Gaggiotti, 2008; Whitlock and Lotterhos, 2015; Benestan et al. 2016; François et al.
270 2016; Gagnaire and Gaggiotti, 2016; Dalongeville et al. 2018). However, some studies have
271 used outliers to investigate subtle geographic breaks in populations when local
272 adaptation is not overcome by migration (Milano et al., 2013; Araneda et al., 2016; Dorant
273 et al., 2022; Tran Lu Y et al., 2021). Outliers were identified using four different methods.
274 Firstly, the population based *Bayescan v.2.1* (Foll and Gaggiotti, 2008) was used to identify
275 markers deviating significantly from distributions of *F*<sub>st</sub> values across all loci (alpha =
276 0.05) corrected for false positives using *q*-values (qvalue v.2.3.0, Storey et al. 2003).
277 Secondly, *Arlequin v.3.5.2.2* (Excoffier & Lischer, 2010), was used to identify outlier SNPs
278 based on *F*<sub>st</sub> values calculated over all the loci, the *p*-values of which were significantly
279 different from the mean (alpha = 0.05). Thirdly, the individual based *PCAdapt v.4.3.3* (Luu
280 et al. 2017) was used to look for loci that are atypically (i.e., not geographically) related to
281 population structure as measured by latent factors also corrected for with *q*-value cutoffs.
282 A Venn diagram was created to visualize the intersection of outliers detected by all
283 methods software.
284

284

285

286

287
288

## 2.5 Population structure and diversity

289

290

291

292

293

294

295

296

297

298

299

300

301

302

303

304

305

306

307

308

309

310

311
## 2.7 Demogenetic connectivity

Population history reconstruction was conducted to investigate demogenetic connectivity. This was done using a modified version of $\partial\alpha\dot{i}$ that uses a dual annealing optimization function in the diffuse approximation to simulate and optimize the joint allele frequency of pairs of populations under different demogenetic scenarios (scripts can be found at (https://github.com/Atranluy/Scripts-Ifremeria/tree/Main/Dadi scripts). Modeling was performed using the observed folded joint allele frequency spectrum (JAFS) on pairs of populations that had an $F_{st}$ value equal to or greater than 0.1. Twenty eight demographic models, described in Rougueux et al. (2017), were simulated ranging from a low to a high complexity of population parameters (S9, figure 9). The models include strict isolation (SI), isolation with migration (IM), ancient migration (AM), and secondary contact (SC), derived from an ancestral population of $N_{anc}$ size that split into two daughter populations of sizes $N_1$ and $N_2$ with or without demographic changes. The size of the ancestral populations is calculated as follows:

$$N_{anc} = \frac{\theta_{anc}}{4 * \mu * L}$$

$$L = \frac{z * y * 32}{x},$$
where x represents the initial total number of SNPs that were detected on y loci of 32 bp; and z represents the number of SNPs used in the $\partial\alpha\partial i$ analyses (see Tran Lu Y et al., 2022). Time (in years) of divergence ($T_{div}$) in the SC model was then calculated as:

337

338

$$T_{div} = \frac{2N_{anc}(T_s + T_{sc})}{33},$$

339

where $T_s$ represents the time (in generations: Time in coalescent units multiplied by $2N_{anc}$) from the secondary contact to the initial split between the populations, $T_{sc}$ represents the time (in generations: Time in coalescent units multiplied by $2N_{anc}$) from present to the secondary contact, and 33 in the number of generations per year.

A mutation rate of $10^{-9}$ was used, based on updated reports of genomic substitution rates vs mutation rates for the crustacea including the Copepoda (Thomas et al., 2020).

Statistics pertaining to each parameter were computed with each run, and those computed for the best model are reported. The uncertainty of the parameters of the best model computed by $\partial\alpha\partial i$ was realized using the Fisher information matrix (FIM), considering our markers independent. FIM uncertainty was achieved with an epsilon = 0.001 from which all standard deviations and subsequent confidence intervals were calculated. Effective Population Size ($N_e$) was calculated using the formula:

352

353

$$N_e = Nu * b * N_{anc},$$

354

355
where *Nu* refers to the daughter population size and *b* represents the growth factor during population size change. Migration rates (Mij: number of migrants exchanged per generation as calculated by $\partial a \partial i$ in coalescent units), were calculated as follows:

$$M_{ij} = (Nu * b * m1 < -2)/2,$$

Where $m1 < -2$, is the migration rate from population 2 to population 1.

## 2.8 Suitability of larval dispersal models in predicting empirical patterns of gene flow in genomic datasets

To investigate the direction and relative strength of migration between different localities, we generated migration networks using the *divMigrate* function (Sundqvist et al., 2016) from the R package *diveRsity* (Keenan et al., 2013) under the null hypothesis of an n-islands model without demographic changes. We used the relative number of migrants (*Nm*) to infer directional gene flow as an index between 0 and 1, where 0 indicates almost no migrants/generation, and 1 indicates one migrant/generation being transported from population A to population B (Sundqvist et al., 2016). Significant asymmetry was calculated by resampling with 1,000 bootstraps. Migration was calculated as a relative value, by comparing either $G_{st}$ of Joost's *D* estimates between pairs of populations and a hypothetical pool of migrants, allowing for the identification of source and sink populations.
Larval dispersal simulations, using Lagrangian particle tracking methods, were conducted using Python version 3.10.2 (http://www.python.org) and PARCELS version 2.3.1 (Delandmeter and Sebillie, 2019). Lagrangian particle tracking, in this case, refers to the advection and diffusion of simulated particles based on the velocity of currents within an Ocean General Circulation Model (OGCM). The OGCM in question was 'MOI GLORYS12_FREE' from MERCATOR, which has a 1/12th of a degree horizontal resolution and daily temporal resolution at 50 depth bands from 0.5 to 5728m depth. Particles were advected by the Northward and Eastward vertical velocities within the OGCM and effected by horizontal diffusion (Smagorinsky, 1963) to approximate sub grid scale currents. Particles were attributed a larval behaviour so that they would approximate the planktonic nauplii of *S. lauensis*. The principal behaviour attributed to these particles was a temperature dependent Larval Dispersal Duration (LDD). Each particle interacted with OGCM temperature data via an equation to calculate LDD, adapted from O'Connor *et al.* (2007) following the methods of Mitarai *et al.* (2016) and Breusing *et al.* (2021). For *S. lauensis*, we assumed a generational turnover of 33 per year (LDD of 11 days) at 20°C and 10 per year (LDD of 36.5 days) at 10°C (Huntley *et al.*, 1991; Gollner *et al.*, 2016) to calculate the magnitude of the relationship between LDD and temperature using the equations of O'Connor *et al.* (2007). Particles were released at the location of the vent fields 100m above the seabed to replicate the entrainment of larvae into vent plumes (Ivanenko and Martinez Arbizu, 2007; Gollner *et al.*, 2015). One hundred particles were simultaneously released from these locations every four hours for six years (2010 - 2016), resulting in a total of 876,000 particles released from each vent field. We used the horizontal distribution of particles at the time they reached their LDD as a measure of larval dispersal potential. We then extrapolated this distribution using a Lagrangian
402 Probability Density Function (PDF) to calculate the horizontal distribution of dispersal probability. The probability densities can be interpreted as the likelihood that a particle settles in a given location where the density values are proportional (Siegel et al., 2023), smoothed across latitude and longitude using a Gaussian distribution. By extracting the value of the PDF from each vent field at the location of all other vents, we determined the probability of dispersal among all pairs of vent fields in both directions following the methods of Mitarai et al. (2009; 2016). Where possible, we implemented the simulation parameters recommended by the authors of PARCELS (Delandmeter and Sebillie, 2019), and provide the detailed simulation parameters in supplementary information (S10, table S11).

412

## 2.9 Genotype-Environment Association Analysis (GEAA)

414 The relationship between genotypes and environmental factors was investigated using *Sambada*, a genome-environment association (GEA) software designed to search for signatures of local adaptation (Duruz et al. 2019). In addition, a distance-based redundancy analysis (db-RDA), combining regression and ordination was used to investigate the relationship between environmental variables and genotypes as per Benestan, Ferchaud, et al. (2016) and Benestan, Quinn, et al. (2016). A spatial eigenfunction analysis was used to evaluate the relative contribution of geographic distance and asymmetric processes driven by deep-water ocean current circulation to genotype frequency variation in the genomic dataset. Firstly, a distance-based Moran's eigenvector map (db-MEMs) (Legendre et al. 2013), was created to decompose the Euclidean (geographic) distances between sites into vectors (db-MEMs) that can be used for subsequent analysis based on a distance matrix. The matrix consists of Euclidean
distances calculated at a depth of 1800m, the average depth at which hydrothermal
plumes were found to travel in the Lau basin (Speer and Thurnherr, 2012), and therefore
the depth at which *S. lauensis* may disperse (Gollner et al. 2015). These additional
variables account for spatial autocorrelation and allow for the identification of non-
random structure. Further, Asymmetric Eigenvector Maps (AEMs) were generated by
translating the presence of a population connection in a larval dispersal probability
matrix (i.e., dispersal probability from the biophysical model > 0) between all pairs of sites
into a sites-by-edge matrix. We attributed a weight to each edge, which was based on the
probability of dispersal between each pair of sites derived from the biophysical model.
When connectivity between a given pair of sites was greater than zero in both directions,
we selected the direction with the highest probability of dispersal (Xuereb et al., 2018).
The db-MEMs therefore describe the geological and bathymetric setting, while the AEMs
describe dispersal based on oceanographic currents. For both db-MEMs and AEMs the
larger eigenvectors model broad-scale spatial structures such as topographic features
 transform faults, between the north and south or differences in the density of vents),
whereas the smaller eigenvectors model fine-scale spatial structures that could describe
other environmental parameters less related to geography (Dray et al., 2012). We used
the package *adespatial* v.0.0-8 (Dray et al., 2017) to generate db-MEMs both types of
spatial variables.

Samβada was then used to identify markers with a significant correlation to the
environmental variables (figure 3, S8, figure 8), based on several statistics: firstly, a G-
score compares the likelihood of modes with and without environmental parameters.
Secondly, the Wald-score checks whether the effect of the environmental parameters (β₂)
is different from zero, i.e., whether an environmental parameter is associated with a significant change in genotype frequencies. Q-values were once again used to correct for false positives. Markers retained after correction were added to the Venn diagram to check if any markers were shared between outlier detection methods and genotype-environment association methods. Lastly, a db-RDA was conducted on a dataset of significant SNPs (vcf was filtered to only contain significant SNPs found by Samβada). During the model selection in a db-RDA, an adjusted R² statistic is used to correct for the number of explanatory variables. Akaike's Information Criterion (AIC) is then used to select the best variables, and a permutation test of F (p-value) using Analysis of Variance (AMOVA) is conducted. This allows for the testing of significance of the global model, individual parameters as well as the individual RDA axis (S8, table 3).

461

## 3. Results

463

### 3.1 *De novo assembly, filtering, and SNP-calling*

A pre-GC-filtered, *de novo* assembled dataset created using discoSnpRAD++ contained 149 individuals represented by 187,737 variant bubbles, and 104,009 SNPs. After filtering for GC, a *de novo* assembled dataset contained 84,674 variant bubbles, and 18,672 SNPs which were subsequently clustered into 9039 loci, using a k-mer length of 15 and a maximum allowed indel size of zero. The dataset contained 28.4% missing data per SNP, and 11 individuals were removed with missing data above 0.5. After filtering for heterozygosity (cut off = 0.6) and 1 SNP per locus, a final SNP dataset of bi-allelic SNPs contained 138 individuals and 5778 high quality SNPs (S4, table 2), with an average
heterozygosity of 0.2. Because high *F*<sub>st</sub> SNPs may also be caused by drift alone (under strong demographic changes), the fact that the target organism is non-model (no existing genomic information), and the very low number of SNPs found across different outlier methods (pointing to little evidence of selection), all SNPs were retained for the structure and demogenetic analysis.

478

### 3.2 Outlier Analysis

Bayesian analysis did not identify any outliers, instead, only putatively balancing, and neutral SNPs. In contrast, Arlequin and PCAadapt identified 123 and 81 outliers, respectively, with 5 shared outliers. Samβada, which was used to explore the relationship between genotype frequencies and environmental variables (db_MEMs and AEMs), identified 671 SNPs that had a significant (alpha = 0.01) correlation with environmental variables. Fifty-two of them were common between Samβada and PCAadapt, and 12 between Samβada and Arlequin, Due to the low number of SNPs identified as outliers across two or more methods, no outliers were removed from the structure or demogenetic analyses (figure 3).

489

### 3.3 Population structure and diversity

Non-metric Multidimensional Scaling identified four populations, three of which were based on geographical location, and one was based on a temporal sub-structure within the ABE hydrothermal field between sampling years 2016 and 2019 (figure 2A, S3, figure 4A, B, and C). Results were compared to PCA and UMAP with MRCC (S3, figure 4). All methods captured the same structure, statistically supported by UMAP with MRCC that
496 was used to account for sample size and variance. The RCSI, entropy, and consensus index,
497 of the latter, suggest that the most significantly stable amount of genetically differentiated
498 groups within the basin is four (S2, figure 2C, E, and F). However, UMAP does not
499 represent genetic distances, and an nMDS approach was therefore chosen to best
500 represent the genetic distances concordantly throughout the paper.

501

502 ADMIXTURE found the most stable number of groups to be three but with almost
503 negligible difference between cross validation error values of 0.321, 0.319, and 0.322 for
504 values of k=2, k=3, and k=4 respectively (figure 2B and S5, figure 5B, C, and D). For k=3,
505 populations are highly structured across the basin from Mangatolo to Tu'i Malila with
506 some (albeit low) level of shared ancestry. We identified three distinct genetic clusters.
507 Cluster 1 corresponds to Tui' Malila, cluster 2 to ABE and Tahi Moana (from which we
508 conclude shared ancestry for the two populations, (which were therefore grouped for
509 downstream analysis) and cluster 3 to Mangatolo. These results were corroborated by
510 STRUCTURE analysis, where the additional testing for K=4 K=5 did not increase the
511 number of ancestral groups. Half of the individuals from Mangatolo, collected in 2019
512 have a low percentage of ancestry belonging to cluster 1, mostly represented by
513 individuals from Tu'i Malila as is the case for the ABE samples from 2019, indicating that
514 copepods sampled from Tu'i Malila in 2016 may have been a source to other populations
515 sampled in 2019. These ancestral pools were therefore termed south (Tu'i Malila), mid
516 (ABE and Tahi Moana), and north (Mangatolo).

517

518 Percentage variation was 17.7% among the populations and 82.3% within the
519 populations, with a global AMOVA Fst of 0.18 (p = <0.001). Pairwise Fst values were highest
between Mangatolo and ABE(2016) ($F_{st}$ = 0.26, *p* = <0.001), followed by Tu'i Malila and ABE(2016) ($F_{st}$ = 0.25, *p* = <0.001), Tahi Moana and Mangatolo ($F_{st}$ = 0.24, *p* = <0.001) and Tahi Moana and Tu'i Malila ($F_{st}$ = 0.24, *p* = < 0.001). $F_{st}$ values between Tahi Moana and ABE (2016) were almost zero ($F_{st}$ = -0.00084, *p* = 0.12). The $F_{st}$ value between Mangatolo and Tu'i Malila was 0.12 (*p* = <0.001). Interestingly, $F_{st}$ values between Mangatolo and ABE(2019) decreased ($F_{st}$ = 0.05, *p* = <0.001), as well as between Tu'i Malila and ABE(2019) ($F_{st}$ = 0.08, *p* = <0.001) and increased between Tahi Moana and ABE(2019) ($F_{st}$ = 0.11, *p* = <0.001), see table 2 for significance codes. Observed heterozygosity was 0.18, 0.12, 0.08, 0.12, and 0.24 for Tu'i Malila, ABE(2016), ABE(2019), Tahi Moana, and Mangatolo, respectively.

530

531

## 3.4 Demogenetic Connectivity

The folded joint allele frequency spectra (JAFS) display the proportion of alleles shared between the north and mid and the north and south populations for the observed data (figure 4A and 5A, respectively) and the simulated spectra for the best model (figure 4B and 5B, respectively). Secondary Contact (SC) appears to be the most appropriate model in each instance and increasing the complexity with the addition of parameters 2N (linked selection), 2m (heterogeneous migration rate across the genome) and G (population size change) improved the AIC values (figure 4E and 5E, respectively). The time of divergence between the north and mid using the best model (SC2N2mG) was calculated to be ~119 kya (UCI = 189 kya, LCI = 49 kya), with a more recent secondary contact event estimated at ~4 kya (UCI = 4.5 kya, LCI = 3.5 kya) while the time of divergence between the north and south populations using the best model (SC2N2mG) was estimated to be ~326 kya (UCI = 795 kya, LCI = 239 kya), with a recent secondary contact event estimated at ~10 kya.

532

533

534

535

536

537

538

539

540

541

542

543
544

kya (UCI = 17 kya, LCI = 3.2 kya). However, absolute estimates are approximately the same order of magnitude with an increase of 2.5 between them, potentially due to uncertainty in the $\partial a / \partial i$ simulations.

545

546

547

548

The proportion (1-P) of the genome that evolved under heterogeneous migration in models with parameter 2m, was 0.26, meaning that almost one third of the genome is subject to constrained gene flow because of barrier loci. The effective population sizes of the daughter populations Nu1 and Nu2, indicate a severe bottleneck in the populations (b1 and b2 < 1) since the time of population split. The local effect of selection at linked sites seems to affect only 20% of loci (Q = 0.17) with a relatively large value of *hrf* (Hill-Roberston factor = ~0.14) (Table 3). Meanwhile, the proportion of the genome that evolved under heterogeneous migration in models with parameter 2m between north and south, was 0.97, meaning that the majority of the genome is subject to constrained gene flow because of barrier loci. The local effect of selection at linked sites seems to affect almost half of the loci (Q = 0.49) with a relatively large value of *hrf* (Hill-Roberston factor = ~0.11) (Table 3).

549

550

551

552

553

554

555

556

557

558

559

560

561

Both comparisons depict asymmetrical migration rate northward, with 5.97 (m21) > 0.68 (m12) for Mid-North and 4.64 (m12) > 1.03 (m21) for North-South (table 3). The number of migrants (estimated by (N1*b1*m12)/2 and (N2*b2*m21)/2) between the north and mid (1=mid & 2=north) is asymmetrical, and predominantly (albeit weakly) northward oriented with (Nm north <- mid = 0.25) than in the opposite direction (Nm mid <- north = 0.12). Nme was used to explain the orientation of gene flow between north and south as this parameter was highly asymmetrical while Nm was very similar between the two

562

563

564

565

566

567
568

569

570

571

572

573

574

575

576

577

578

579

580

581

582

583

584

585

586

587

populations (Table 3). Gene flow was predominantly northward oriented between the
two populations as captured by Nme (Nme$_{north<-south}$ = 2.40, Nme$_{south<-north}$= 0.66). In this
case, Nme does not represent true migrants, but rather the movement of introgressed
alleles that are not impacted by barrier loci.

**3.6 Effective Population Size (Ne)**

Contemporary *N*<sub>e</sub> estimates using parameters from *dadi* were high. When comparing the north and mid populations, *N*<sub>e</sub> was calculated as 299,738 and 202,206 for the north and mid, respectively, showing a sharp decrease in *N*<sub>e</sub> from an *N*<sub>anc</sub> of 401,202 since the T<sub>div</sub> ~123 kya. When comparing the north and south, *N*<sub>e</sub> was calculated as 701,645 and 865, 000 for the north and south, respectively, indicating an increase in *N*<sub>e</sub> from an *N*<sub>anc</sub> of 677,263 since the T<sub>div</sub>~335 kya. Nevertheless, it remains complicated to distinguish the differences between a decrease in the effective population size caused by a severe bottleneck, as this signature in structured populations can also be caused by a reduction of the gene flow over time, and these results are therefore interpreted with caution.

**3.7 Suitability of larval dispersal models in predicting empirical patterns of gene flow in genomic datasets**

**3.7.1 Larval dispersal models**

Larval dispersal probability models were run with a passive vertical migration parameter, a growth curve for copepods based on temperature response, and millions of particles released daily from 2010-2016. Particles dispersed predominantly northward, particularly in the case of the mid and south vent fields. Particles dispersing from the

588

589

590

591
north dispersed further east and west than the other sites. Larval dispersal from the south and mid was constrained by the topography around the Peggy Ridge and the presence of the Niuafo'ou plate. The lowest pairwise probability of dispersal among vent fields was 7.25⁻³² from ABE to Mangatolo (northward), while the highest was 3.26⁻⁷ from Tui Malila to ABE (northward) (figure 6A). Despite the low probability values, all vent fields were connected to one another by very limited dispersal (values above 0) except the north (figure 6B), and all vent fields displayed high levels of self-recruitment, highest so in the mid and south. The average LDD of Lagrangian particles dispersing from all vent fields was ~104 days (S10, figure 11).

### 3.7.2 Migration inferred from genomic data

Using empirical allele frequency data, *DivMigrate* identified significant, asymmetric, northward migration from the southern and mid populations to the north (Figure 6A). From the south to the north, *Nm* was 1.0 per generation, describing unidirectional gene flow in purely northward orientation. From the mid to the north, *Nm* was 0.31 indicative of a weaker, northward gene flow. From the mid to the south, *Nm* was 0.33, indicative of a weak, southward gene flow (Figure 6A). These values were calculated under the strong assumption of a n-islands model with no demographic changes, however, results from *DivMigrate* are comparable to those from *daći*, with strong northward, asymmetric geneflow, and approximately one third of each generation of migrants emigrating from the mid to the north, and over twice the amount of migrants per generation emigrating from the south to the north.

## 3.8 Genotype-Environment Association Analysis
A db-RDA conducted on a subset of 671 SNPs significantly associated with environmental variables revealed that structural variables of intermediate scales (MEM4 and MEM5) appeared to correlate significantly with differences between Tahi and ABE on one hand, and Tui Malila and Mangatolo, on the other hand (S8, figure 8, table 3). An AEM defining high probabilities of dispersal (AEM1) also appeared to correlate significantly with differences between these two groups of populations. The RDA had a global R² value of 0.84 (p-value <0.001). MEM5, MEM4, and AEM1 were highly significant with p-values of <0.001 for all parameters. RDA1, explaining 83.78 % of the genomic variation was also highly significant with a p-value of <0.001, while RDA2 explained 3.21 % of the variation (p-value = <0.001).

626

627

## 4. Discussion

This study examines population connectivity of the copepod *S. lauensis* in the Lau basin, a hydrothermal vent habitat in the Southwest Pacific BABs. Using 5,778 SNPs, we confirmed strong spatial and temporal structuring, with metapopulations originating from three ancestral pools. Population separation likely occurred around ~119 kya (north-mid) and ~326 kya (north-south), with secondary contacts at ~4 kya and ~10 kya, respectively. Weak migration is observed from mid to north and south, while stronger northward migration occurs from south to north. Bespoke larval dispersal patterns captured similar directionality in the gene flow with the exception of gene flow from the mid to the south and supported high levels of structure among all populations.

637

638

### 4.1 Fine-scale population structure
Previous studies have shown that within-basin structure is unlikely for megafaunal species with planktonic larval dispersal, except for *Munidopsis lauensis* in the Manus basin (Mitarai et al., 2016; Plouviez et al., 2019; Breusing et al., 2021; Tran Lu Y et al., 2021; Poitrimol et al., 2022; Thaler et al., 2014). Although habitat preference and interspecific dynamics are known to influence connectivity patterns (Micheli et al., 2002; Lenihan et al., 2008; Beinart et al., 2012; Johnson et al., 2020), *S. lauensis* experiences high levels of self-recruitment due to its short generation time, reducing the number of migrants being exchanged between populations, thus allowing populations to diverge faster through genetic drift in the absence of migration. Particularly in areas of diffuse hydrothermal venting where *S. lauensis* dominates, nauplii will have less chance of being entrained in a plume vortex and dispersed far into the water column as they are not in close proximity to the vent source. Pronounced geographic structure was evident from ADMIXTURE, IBS, and multiple dimension-reducing techniques within this study.

Temporal genetic structure has been reported for copepods in surface waters resulting from demogenetic changes, due to variations in environmental conditions over time (Tepper and Bradley, 1989; Posavi et al., 2014; Brennan et al., 2022). In the present study, such temporal variations were also depicted at the ABE vent field, as seen in the IBS analysis (S7, figure 7). Allelic introgression from Tu'i Malila (as observed from the ADMIXTURE results in Figure 2B) point to Tu'i Malila having acted as a source for most Lau basin populations after 2016. Furthermore, tectonic rearrangements occurred between 2016 and 2019 with the extinction of the site Kilo Moana which was closer to ABE than Tui' Malila is. However, as the effects of genetic drift act more strongly on small populations due to founder effects, it is also possible that a disruption to the ABE field
caused a sporadic bottleneck. This is supported by the reduction in heterozygosity from 0.12 to 0.08 in samples from ABE from 2016 to 2019. A generation time of 33 generations/year, as reported by Gollner et al. (2016) for copepods at 20°C, was also used, reflecting the r-strategy mechanisms displayed by hydrothermal vent fauna in response to the environmental instability of the vent habitat (O'Connor et al., 2007; Tran Lu Y et al., 2021). Drift acting on a population of reduced size would be especially non-negligible if dirivultid copepods are undergoing 33 generations per year, allowing drift to act on a turnover of almost 100 generations in three years. Bottlenecks also appear to have occurred during the secondary contact event between the north and south, as shown by a decrease in the b1 parameter from the *δαδι* results (Table 3). Disruptions to the sea floor of the Lau Basin are evident after the recent Tonga eruption blanketed a large proportion of the Lau Basin, affecting vast areas of seafloor, including ABE (Seabrook et al., 2023). There is evidence of catastrophic eruptions around the Tonga region throughout the last thousand years (Monzier et al., 1994; Lavigne et al., 2021; Brenna et al., 2022), and this activity has increased over the entire Tonga-Kermadec region within this period of time (Watt et al., 2013), indicating that significant changes to the seafloor may have occurred, giving rise to recurrent bottlenecks in populations of benthic organisms.

## 4.2 Demogenetic history of *S. lauensis* within the Lau basin

The demogenetic history of *S. lauensis* was investigated between the southern, mid, and northern populations of the basin. A mutation rate of $1 \times 10^{-9}$, based on relatively recent reports for crustaceans including copepods and amphipods (Thomas et al., 2020), was used to estimate the time of divergence between the populations. This rate is lower than
the typical rate for mitochondrial DNA and aligns with evolutionary rates between mitochondrial and nuclear genomes.

The commencement of spreading in the Mangatolo Triple Junction (MTJ) is estimated to have occurred between 180 and 600 kya (Mensing et al., 2020, Steward et al., 2022), and the birth of the major transform fault, Peggy Ridge, that separates the Australian plate and the Niuafo'ou plate is estimated to be no older than 350 kya (Stewart et al. 2022). The latter appears to act as a major barrier to dispersal at depths below 600 m, and consistently, the time of split between metapopulations is estimated to have occurred after the formation of this fault (327 kya and 119 kya for the north-mid and north-south, respectively). Time of secondary contact events (~4 kya and 10 kya for the north-mid and north-south, respectively) align with Cox1 data from the same species and vents in the region (Diaz-Recio Lorenzo et al., 2023), that point to an increase in effective population size starting at 10 kya in the Lau basin. During this period, eruptions within the Tonga-Kermadec arc went from making up 10% of the global eruptions to 50%, likely giving rise to new vents across the region, and thus more stepping stone features along which dispersal could have occurred (Watt et al., 2013). This would be consistent with the occurrence of an SC event with sufficient genetic exchanges to cause a mixing of cox1 haplotypes from different geographic regions as reported by Diaz-Recio Lorenzo et al. (2023), while genetic incompatibilities display a pronounced barrier with limited gene flow over a non-negligible portion of the nuclear genome due to a long period of isolation. Despite the rather high standard deviation values of theta for the north-mid (theta = 91.82, SD = 23.55) and the north-south (theta = 155.82, SD = 93.6), demogenetic connectivity patterns appear to be roughly consistent with seabed changes over
711 geological timescales and the more nuanced events that recently contributed to the
712 complex geological setting of the Lau basin (Ruellan 2003; Martinez and Taylor, 2006;
713 Baker et al., 2019; Baxter et al., 2020). However, stochastic hydrodynamic events such as
714 deep-sea eddies (Speer and Thunher et al., 2012), which may modulate dispersal
715 sporadically alongside significant El Nino events, and modifications to regional surface
716 and deep ocean currents due to the onset of the Last Glacial Maximum, may also have
717 played a role (Hu and Piotrowski et al., 2020, Tran Lu Y et al., 2021). The values used in
718 this study are estimates and should therefore be interpreted cautiously. Assumptions
719 were made regarding generation turnover and non-overlapping generations and may not
720 have fully captured the complex reproductive and dispersal mode of *S. lauensis*. Further
721 research on female fecundity, ecology, and stressors is needed to validate these
722 assumptions for this specific species and the dirivultid family, in general.
723
724
725
726
### 4.2.1 Estimation of Effective Population Size (Ne)
727 An overestimation of *N*<sub>e</sub> can be expected when there is immigration from populations with
728 limited differentiation to the focal population, resulting in less drift due to the increase in
729 alleles introduced that are already present in the population (Nadachowska-Brzyska et
730 al., 2021). Conversely, we expect *N*<sub>e</sub> to be underestimated when structured populations
731 exhibit migrants, introducing foreign alleles, leaving signatures in the genome that
732 resemble genetic drift, which could be the case for *S. lauensis*. Based on the structure
733 analysis it appears that individuals from ABE and Mangatolo from 2019 are admixed with
734 individuals from Tu'i Malila from 2016, which increase temporal variations in these two
735

populations and would add to the effect of drift in the structured populations, and yet the
estimates using $\bar{\alpha}\bar{\alpha}_i$ were high with values ranging from 192,576-512,415. However, this
may also be a result of pooling the two temporal states with use in $\bar{\alpha}\bar{\alpha}_i$, and results of $N_e$
must therefore be interpreted with caution.

736

737

738

739

740

## 4.3 Suitability of larval dispersal models in predicting empirical patterns of gene flow in genomic datasets

741

742

743

Larval dispersal models tested two main hypotheses: depth restriction to the spatial distribution of larvae and gene flow orientation across the Lau basin. The models incorporated realistic scenarios of potential vertical movement and temperature-dependent growth, considering the complex current dynamics within the Lau basin (Speer and Thunherr, 2012; Simons et al., 2021) and the relationship between LDD and temperature (O'Conner et al., 2007). The average larval dispersal duration was estimated to be around 104 days, but potentially with a much shorter duration for self-recruiting copepods in the local vent environment due to high temperatures.

744

745

746

747

748

749

750

751

752

The hypothesis that copepod nauplii travel over short distances with deep-water hydrothermal plumes was confirmed by genetic differentiation of the *cox1* gene of populations among basins (Diaz-Recio Lorenzo et al., 2023) and larval dispersal simulations which showed that particles began leaving the Lau basin from each vent at depths of 600 m and shallower. The simulations indicated highly structured populations, even with a long larval dispersal duration of 104 days. Given the low level of gene flow inferred from the genomic data, it appears that few migrants may be able to travel for a

753

754

755

756

757

758
759 prolonged period of time, likely in a state of delayed metamorphosis, acting as weak
760 connections between populations. Vertical passive movement allowed particles to cross
761 the transform fault separating the Australian plate and the Niuafo'ou plate, enabling low-
762 level gene flow. However, southward migration probabilities from the mid to the south
763 were zero, despite weak gene flow observed in genomic analyses. Nevertheless,
764 connections may also occur due to stepping-stones and ghost populations aided by the
765 many other vents that exist in the basin (denoted by the white triangles in figure 1) and
766 seafloor fissures forming potential connectivity corridors aided potentially by bottom
767 currents and sub-surface hydrothermal flow, elements not accounted for in the models.
768

769 This study identified the mid population, as a contemporary source in 2016, contributing
770 migrants to the north and south. Vents flanking the mid population, including Tow Cam
771 and Kilo Moana, form a cluster facilitating bidirectional gene flow. Tu'i Malila in the south
772 serves as an important source population for northward migration in 2019, however
773 historical changes in the distribution of source populations due to tectonic and volcanic
774 events may have influenced population structure through time. These findings emphasize
775 the significance of these populations in maintaining gene flow towards the north which
776 appears to act as a sink, and highlight the need to consider them in conservation
777 strategies, however they also highlight the need to incorporate time series data, as
778 populations may act as source or sinks depending on the sampling year. While
779 acknowledging the limitations of larval dispersal models, this study demonstrates their
780 value in understanding complex dispersal dynamics by incorporating organism-specific
781 traits.
782

782
The results of the GEAA analysis showed that intermediate-scale bathymetric structures and fine-scale asymmetric dispersal played a role in structuring *S. lauensis* populations. MEM4, MEM5 (intermediate scale db-MEMs associated with the spatial arrangement of vents), and AEM1 (a fine scale asymmetric vector describing high probabilities of dispersal predicted by Lagrangian simulations) correlated significantly with the observed differences in genotypic frequencies between populations (S8, figure 8, table 3). The db-MEMs used in the analysis likely captured geological features such as vent density, transform faults, depth gradients, and axial ridges in the Lau basin. MEM5 was associated mostly with values of depth, while MEM4 was associated with values of spreading rate. Despite only having one value per vent for depth and spreading rate, these associations appear meaningful, given the depth gradient between Tu'i Malila, ABE, and Tahi Moana (1888 m, 2155 m, and 2273 m, respectively). Meanwhile, the AEM1 was associated with high probabilities of dispersal, suggesting a higher probability of dispersal between the north-south, compared to the mid-south and the mid-north. Depth may play an important role in structuring vent populations, as nauplii from similarly deep vents, would need to cross a smaller depth gradient in order to disperse (Breusing et al., 2021). Spreading rate is an important factor governing both local extinction rate and the rugosity of the seabed (faults and screes) at each vent location. The differences in spreading rates between the north (38 mm/year -1), mid (69mm/year -1), and south (48 mm/year -1) may result in differences in local vent conditions and bacterial communities on which the copepods feed, likely exerting differences among populations, via local adaptation and recurrent bottlenecks. Furthermore, the frequency of vent occurrence is positively correlated with magmatic supply (Baker and German, 2004), which increases with spreading rate (Olive and Dublanchet et al., 2020) implying that fast spreading zones may harbor dense
aggregations of vent sites, providing stepping-stones for copepods to disperse, a theory evaluated for megafauna of the Southwest Pacific vents (Breusing et al., 2021) and also the East Pacific Rise (Jollivet et al., 1999; Vrijenhoek et al., 2010). To the north of ABE and Tahi Moana, spreading rates increase dramatically from 69 mm/ year⁻¹ to 120 mm/ year⁻¹ (Stewart et al., 2022), as do the number of confirmed vents along that gradient. The aggregations can be seen in e.g., figure 1, where the white triangles denoting known vents appear to form hotspots in the northeast, northwest, and mid/south of the basin. These vents form two main stepping-stone corridors along which vent fauna are distributed, one of which could explain preferential dispersal between north-south, but not between north-mid.

## 4.4 Connectivity of hydrothermal vent meiofauna vs megafauna – implications for resilience to deep-sea mining

Relative migration inferred from empirical genetic data, agrees largely with that of Breusing et al., (2021), who show that populations of *Alviniconcha* snails exhibit a predominantly southward migration from Tahi Moana to Tu'i Malila, and a predominantly northward migration from both of those sites to the north. Plouviez et al., (2019) depicted some population structure for the limpet *Lepetodrilus schrolli* using the *cox1* gene and Poitrimol et al. (2022) showed that cryptic species also exist for a range of megafauna in the Southwest Pacific including within the Lau basin based on *cox1*. However, Breusing et al.'s study as well as other recent population genomics studies (Plouviez et al., 2019; Tran Lu Y et al., 2021) and larval dispersal modelling of megafauna in the region (Mitarai et al., 2016) have found little to no structure in the populations within the Lau basin when using
nuclear markers. Our results highlight that pattern of connectivity cannot be generalized to all fauna in vent systems and that the different faunal groups need to be considered when developing robust area-based management tools (ABMTs). The removal of these abundant copepods from the vent environment through mining may have unprecedented consequences for the ecosystem as a whole. Being the most abundant animal where hard substrate is available, dirivultid copepods represent an unexplored but potentially important functional group in these systems, feeding predominantly on chemoautotrophic bacteria (Nomaki et al., 2023) and forming close associations with a range of megafaunal species including *Rimicaris* shrimp (Gollner et al., 2016), the squat lobster *Shinkai crosnieri* (Senokuchi et al., 2018), Bathymodiolin mussels (Gollner 2016), *Ifremeria* and *Alviniconcha* snails (Gollner et al., 2016) and the vent polychaete *Paralvinella sulfincola*, which was hypothesized to utilize the vent copepod *Stygiopontius quadrispinosus* as a primary food source when densities of the copepod were high (1852 ind. per worm) (Limen et al., 2008). Copepods in the global oceans are hugely important bioindicators, and the role of deep-sea copepods at hydrothermal vents as such is largely unexplored. Further research is needed to investigate the trophic role of this family in the vent environment and the potential risks to the ecosystem should populations of abundant copepods suffer major bottlenecks as a result of anthropogenic disturbance.

849

## 5. Conclusions

851

This study reports the first data on population connectivity in the vent meiofauna, specifically for the vent obligate copepod *Stygiopontius lauensis*. Limited gene flow between north and south populations suggests a potential stepping-stone dispersal model

854
and recent secondary contact events associated with vent dynamics. Divergence between populations likely occurred hundreds of thousands of years ago, with secondary contact events occurring in recent geological history. Our results suggest that migration between vents in the Lau basin is not sufficiently high to offset the effects of drift or local selection, likely leaving populations at risk of local extinction, should mining occur. Larval dispersal models incorporating temperature-dependent dispersal duration and passive vertical movement accurately captured population isolation. However, they did not fully explain the observed low-level southward gene flow from populations in the mid and south, and further research is needed to address this discrepancy. We note that despite the weaknesses in using any single analysis alone to draw conclusions about resilience of vent communities, a multifaceted approach combining multiple complementary analyses provides a much more robust picture of the consequences of anthropogenic activity in these ecosystems. These results add much needed data to connectivity studies in the region, informing dispersal models in the framework of understanding populations at risk of extinction, and shed light on a majorly understudied subgroup of organisms in the Southwest Pacific BAB system.

871

872 Acknowledgements

We warmly thank the crew of L'Atalante, the scientific members of the cruise who sorted the copepod samples (Camille Poitrimol and Eve-Julie Arsenault-Pernet) and the chief scientists of the Chubacarc cruise (Didier Jollivet and Stéphane Hourdez) for kindly donating the 2019 samples. We also warmly thank Roxanne Beinart, chief scientist, for the donation of the samples collected during the 2016 FK160407 cruise and ABIMS Bioinformatics platform of the Station Biologique de Roscoff (France) which provided
computing resources for the demographic analysis. Ship time was supported by the
879
French Oceanographic Fleet programme (CHUBACARC cruise https://
880
doi.org/10.17600/18001111 to Didier Jollivet and Stephane Hourdez), INEE (CNRS) and
881
the Agence Nationale de la Recherche (ANR) that funded the project "CERBERUS"
882
 contr
883
(contract number ANR-17CE02-0003 to Stephane Hourdez). Sequencing was funded by
884
the UU-NIOZ under the project "Protecting deep-seabed hydrothermal vents via area-
885
based management tools" (to Sabine Gollner). We would also like to thank the two
886
anonymous reviewers for their suggestions and comments to improve the manuscript
887

887

888 **Data availability**

889 Raw sequence reads (Individual fasta files), individual genotype data (VCF) and metadata
890 are available on DataDryad (https://doi.org/10.5061/dryad.zkh1893g2). Scripts used in
891 this study (R, eton) are available on a public Github repositories:
892 https://github.com/Atranluy/Scripts-Ifremeria,
893 https://github.com/otistwo/hydrothermal_thesis/blob/main/PARCELS/Copepod PARC
894 ELS.ipynb, and
895 https://github.com/otistwo/hydrothermal_thesis/blob/main/PARCELS/Copepod_post
896 hoc.ipynb.

897

898

899 **4. Bibliography**

900

901 Alexander, D.H., Lange, K., 2011a. Enhancements to the ADMIXTURE algorithm for
902 individual ancestry estimation. BMC Bioinformatics 12, 246.
903 https://doi.org/10.1186/1471-2105-12-246

904

905 Araneda, C., Larraín, M.A., Hecht, B. and Narum, S., 2016. Adaptive genetic variation distinguishes Chilean blue mussels (Mytilus chilensis) from different marine environments. Ecol Evol. 6: 3632-3644. https://doi.org/10.1002/ece3.2110

908

909 Baker, E.T., German, C.R., 2013. On the Global Distribution of Hydrothermal Vent Fields, in: German, C.R., Lin, J., Parson, L.M. (Eds.), Geophysical Monograph Series.
910
911 American Geophysical Union, Washington, D. C., pp. 245-266.

912 https://doi.org/10.1029/148GM10

913

914 Baker, E.T., Walker, S.L., Massoth, G.J., Resing, J.A., 2019. The NE Lau Basin: Widespread and Abundant Hydrothermal Venting in the Back-Arc Region Behind a Superfast Subduction Zone. Front. Mar. Sci. 6, 382.

917 https://doi.org/10.3389/fmars.2019.00382

918

919 Bashevkin, S.M., Dibble, C.D., Dunn, R.P., Hollarsmith, J.A., Ng, G., Satterthwaite, E.V., Morgan, S.G., 2020. Larval dispersal in a changing ocean with an emphasis on upwelling regions. Ecosphere 11. https://doi.org/10.1002/ecs2.3015

922

923 Beedessee, G., Watanabe, H., Ogura, T., Nemoto, S., Yahagi, T., Nakagawa, S., Nakamura, K., Takai, K., Koonjul, M., Marie, D.E.P., 2013. High Connectivity of Animal Populations in Deep-Sea Hydrothermal Vent Fields in the Central Indian Ridge Relevant to Its Geological Setting. PLoS ONE 8, e81570.

927 https://doi.org/10.1371/journal.pone.0081570

928

929 Beinart, R.A., Sanders, J.G., Faure, B., Sylva, S.P., Lee, R.W., Becker, E.L., Gartman, A., Luther, G.W., Seewald, J.S., Fisher, C.R., Girguis, P.R., 2012. Evidence for the role of endosymbionts in regional-scale habitat partitioning by hydrothermal vent symbioses. Proc. Natl. Acad. Sci. U.S.A. 109.

933 https://doi.org/10.1073/pnas.1202690109

934

935 Benestan, L., Quinn, B.K., Maaroufi, H., Laporte, M., Clark, F.K., Greenwood, S.J., Rochette, R., Bernatchez, L., 2016. Seascape genomics provides evidence for thermal adaptation and current-mediated population structure in American lobster (Homarus americanus). Mol Ecol 25, 5073-5092.

939 https://doi.org/10.1111/mec.13811

940

941 Benestan, L.M., Ferchaud, A., Hohenlohe, P.A., Garner, B.A., Naylor, G.J.P., Baums, I.B., Schwartz, M.K., Kelley, J.L., Luikart, G., 2016. Conservation genomics of natural and managed populations: building a conceptual and practical framework.

944 Molecular Ecology 25, 2967-2977. https://doi.org/10.1111/mec.13647

945

946 Benestan, L.M., Rougemont, Q., Senay, C., Normandeau, E., Parent, E., Rideout, R., Bernatchez, L., Lambert, Y., Audet, C., Parent, G.J., 2021. Population genomics and history of speciation reveal fishery management gaps in two related redfish species (Sebastes mentella and Sebastes fasciatus). Evol Appl 14, 588-606.

950 https://doi.org/10.1111/eva.13143

951

952 Bernatchez, S., Xuereb, A., Laporte, M., Benestan, L., Steeves, R., Laflamme, M., Bernatchez, L., Mallet, M.A., 2019. Seascape genomics of eastern oyster (Crassostrea virginica) along the Atlantic coast of Canada. Evol Appl 12, 587-609. https://doi.org/10.1111/eva.12741

955

956

957
958

Blanchard, C., Gollner, S., 2022. Area-based management tools to protect unique

959

hydrothermal vents from harmful effects from deep-sea mining: A review of

960

ongoing developments. Front. Polit. Sci. 4, 1033251.

961

https://doi.org/10.3389/fpos.2022.1033251

962

963

Blasiak, R., Wynberg, R., Grorud-Colvert, K., Thambisetty, S., Bandarra, N.M., Canário,

964

A.V.M., Da Silva, J., Duarte, C.M., Jaspars, M., Rogers, A., Sink, K., Wabnitz, C.C.C.,

965

2020. The ocean genome and future prospects for conservation and equity. Nat

966

Sustain 3, 588-596. https://doi.org/10.1038/s41893-020-0522-9

967

968

Brenna, M., Cronin, S. J., Smith, I. E. M., Pntesilli, A., Tost, M., Barker, S. J., Tonga'onevai,

969

S., Kula, T., Vaiomounga R., 2022. Post-caldera volcanism reveals shallow

970

priming of an intra-ocean arc andesitic caldera: Hunga volcano, Tonga, SW

971

Pacific. Lithos 412-413. Lithos 412-413.

972

https://doi.org/10.1016/j.lithos.2022.106614

973

974

Brennan, R.S., deMayo, J.A., Dam, H.G., Finiguerra, M., Baumann, H., Buffalo, V., Pespeni,

975

M.H., 2022a. Experimental evolution reveals the synergistic genomic

976

mechanisms of adaptation to ocean warming and acidification in a marine

977

copepod. Proc. Natl. Acad. Sci. U.S.A. 119, e2201521119.

978

https://doi.org/10.1073/pnas.2201521119

979

980

Brennan, R.S., deMayo, J.A., Dam, H.G., Finiguerra, M., Baumann, H., Buffalo, V., Pespeni,

981

M.H., 2022b. Experimental evolution reveals the synergistic genomic

982

mechanisms of adaptation to ocean warming and acidification in a marine

983

copepod. Proc. Natl. Acad. Sci. U.S.A. 119, e2201521119.

984

https://doi.org/10.1073/pnas.2201521119

985

986

Breusing, C., Biastoch, A., Drews, A., Metaxas, A., Jollivet, D., Vrijenhoek, R.C., Bayer, T.,

987

Melzner, F., Sayavedra, L., Petersen, J.M., Dubilier, N., Schilhabel, M.B., Rosenstiel,

988

P., Reusch, T.B.H., 2016. Biophysical and Population Genetic Models Predict the

989

Presence of “Phantom” Stepping Stones Connecting Mid-Atlantic Ridge Vent

990

Ecosystems. Current Biology 26, 2257-2267.

991

https://doi.org/10.1016/j.cub.2016.06.062

992

993

Breusing, C., Johnson, S.B., Mitarai, S., Beinart, R.A., Tunnicliffe, V., 2023. Differential

994

patterns of connectivity in Western Pacific hydrothermal vent metapopulations:

995

A comparison of biophysical and genetic models. Evolutionary Applications 16,

996

22-35. https://doi.org/10.1111/eva.13326

997

998

Breusing, C., Johnson, S.B., Tunnicliffe, V., Vrijenhoek, R.C., 2015. Population structure

999

and connectivity in Indo-Pacific deep-sea mussels of the Bathymodiolus

1000

septemdierum complex. Conserv Genet 16, 1415-1430.

1001

https://doi.org/10.1007/s10592-015-0750-0

1002

1003

Brunner, O., Methou, P., Mitarai, S., 2023. The Role of Reproductive Periodicity in

1004

Dispersal Among Hydrothermal Vents and its Implications for Regional
1005

Connectivity and Conservation (preprint). Biophysics.

1006

https://doi.org/10.1101/2023.03.07.531641

1007

1008

Casanova, A., Maroso, F., Blanco, A., Hermida, M., Ríos, N., García, G., Manuzzi, A., Zane,

1009

L., Verissimo, A., García-Marín, J.-L., Bouza, C., Vera, M., Martínez, P., 2021. Low

1010

impact of different SNP panels from two building-loci pipelines on RAD-Seq

1011

population genomic metrics: case study on five diverse aquatic species. BMC

1012

Genomics 22, 150. https://doi.org/10.1186/s12864-021-07465-w

1013

1014

Fisher C. R., and Girguis P., 2016. cruise doi: https://doi.org/10.7284/906664

1015

1016

Chattopadhyay, B., Garg, K.M., Ramakrishnan, U., 2014. Effect of diversity and missing

1017

data on genetic assignment with RAD-Seq markers. BMC Res Notes 7, 841.

1018

https://doi.org/10.1186/1756-0500-7-841

1019

1020

Craddock, C., Hoeh, W.R., Gustafson, R.G., Lutz, R.A., Hashimoto, J., Vrijenhoek, R.J.,

1021

1995. Evolutionary relationships among deep-sea mytilids (Bivalvia: Mytilidae)

1022

from hydrothermal vents and cold-water methane/sulfide seeps. Marine

1023

Biology 121, 477-485. https://doi.org/10.1007/BF00349456

1024

1025

Crotti, M., Barratt, C.D., Loader, S.P., Gower, D.J., Streicher, J.W., 2019. Causes and

1026

analytical impacts of missing data in RADseq phylogenetics: Insights from an

1027

African frog (Afrixalus). Zool Scr 48, 157-167.

1028

https://doi.org/10.1111/zsc.12335

1029

1030

Dalongeville, A., Benestan, L., Mouillot, D., Lobreaux, S., Manel, S., 2018. Combining six

1031

genome scan methods to detect candidate genes to salinity in the

1032

Mediterranean striped red mullet (Mullus surmuletus). BMC Genomics 19, 217.

1033

https://doi.org/10.1186/s12864-018-4579-z

1034

1035

De Jode, A., Le Moan, A., Johannesson, K., Faria, R., Stankowski, S., Westram, A.M., Butlin,

1036

R.K., Rafajlović, M., Fraïsse, C., 2023. Ten years of demographic modelling of

1037

divergence and speciation in the sea. Evolutionary Applications 16, 542-559.

1038

https://doi.org/10.1111/eva.13428

1039

1040

De Wit, P., Faust, E., Green, L., Jahnke, M., Pereyra, R.T., Rafajlović, M., 2023. A decade of

1041

progress in marine evolutionary biology. Evolutionary Applications 16, 193-

1042

201. https://doi.org/10.1111/eva.13523

1043

1044

Delandmeter, P., Van Sebille, E., 2019. The Parcels v2.0 Lagrangian framework: new

1045

field interpolation schemes. Geosci. Model Dev. 12, 3571-3584.

1046

https://doi.org/10.5194/gmd-12-3571-2019

1047

1048

Diaz-Papkovich, A., Anderson-Trocmé, L., Ben-Eghan, C., Gravel, S., 2019. UMAP reveals

1049

cryptic population structure and phenotype heterogeneity in large genomic

1050

cohorts. PLoS Genet 15, e1008432.

1051

https://doi.org/10.1371/journal.pgen.1008432
1052

1053

1054

1055

1056

1057

1058

1059

1060

1061

1062

1063

1064

1065

1066

1067

1068

1069

1070

1071

1072

1073

1074

1075

1076

1077

1078

1079

1080

1081

1082

1083

1084

1085

1086

1087

1088

1089

1090

1091

1092

1093

1094

1095

1096

1097

Diaz-Papkovich, A., Anderson-Trocmé, L., Gravel, S., 2021. A review of UMAP in population genetics. J Hum Genet 66, 85-91. https://doi.org/10.1038/s10038-020-00851-4

Diaz-Recio Lorenzo, C., Bruggen, D.T., Luther, G.W., Gartman, A., Gollner, S., 2021a. Copepod assemblages along a hydrothermal stress gradient at diffuse flow habitats within the ABE vent site (Eastern Lau Spreading Center, Southwest Pacific). Deep Sea Research Part I: Oceanographic Research Papers 173, 103532. https://doi.org/10.1016/j.dsr.2021.103532

Diaz-Recio Lorenzo, C., Patel, T., Arsenault-Pernet, E-J., Poitrimol, C., Jollivet, D., Martinez Arbizu., Gollner S., 2023. Highly structured populations of 42eep-sea copepods associated with hydrothermal vents across the Southwest Pacific, despite contrasting life history traits. PLoS ONE 18(11): e0292525. https://doi.org/10.1371/journal.pone.0292525

Do, C., Waples, R.S., Peel, D., Macbeth, G.M., Tillett, B.J. and Ovenden, J.R., 2014. NeEstimator v2: re-implementation of software for the estimation of contemporary effective population size (Ne) from genetic data. Mol Ecol Res. 14: 209-214. https://doi.org/10.1111/1755-0998.12157

Dorant, Y., Laporte, M., Rougemont, Q., Cayuela, H., Rochette, R., & Bernatchez, L. (2022). Landscape genomics of the American lobster (*Homarus americanus*). Molecular Ecology, 31, 5182–5200. https://doi.org/10.1111/mec.16653

Dray, S., Pélissier, R., Couteron, P., Fortin, M.-J., Legendre, P., Peres-Neto, P.R., Bellier, E., Bivand, R., Blanchet, F.G., De Cáceres, M., Dufour, A.-B., Heegaard, E., Jombart, T., Munoz, F., Oksanen, J., Thioulouse, J., Wagner, H.H., 2012. Community ecology in the age of multivariate multiscale spatial analysis. Ecological Monographs 82, 257-275. https://doi.org/10.1890/11-1183.1

Drazen, J.C., Smith, C.R., Gjerde, K.M., Haddock, S.H.D., Carter, G.S., Choy, C.A., Clark, M.R., Dutrieux, P., Goetze, E., Hauton, C., Hatta, M., Koslow, J.A., Leitner, A.B., Pacini, A., Perelman, J.N., Peacock, T., Sutton, T.T., Watling, L., Yamamoto, H., 2020. Midwater ecosystems must be considered when evaluating environmental risks of deep-sea mining. Proc. Natl. Acad. Sci. U.S.A. 117, 17455-17460. https://doi.org/10.1073/pnas.2011914117

Du Preez, C., Fisher, C.R., 2018. Long-Term Stability of Back-Arc Basin Hydrothermal Vents. Front. Mar. Sci. 5, 54. https://doi.org/10.3389/fmars.2018.00054

Duruz, S., Sevane, N., Selmoni, O., Vajana, E., Leempoel, K., Stucki, S., Orozco-terWengel, P., Rochat, E., Dunner, S., The NEXTGEN Consortium, The CLIMGEN Consortium, Bruford, M.W., Joost, S., 2019. Rapid identification and interpretation of gene-environment associations using the new R.SamBada landscape genomics
1098 pipeline. Mol Ecol Resour 19, 1355-1365. https://doi.org/10.1111/1755-0998.13044

1100

1101 Eaton, D.A.R., Spriggs, E.L., Park, B., Donoghue, M.J., 2016. Misconceptions on Missing Data in RAD-seq Phylogenetics with a Deep-scale Example from Flowering Plants. Syst Biol syw092. https://doi.org/10.1093/sysbio/syw092

1104

1105 Eldon, B., Riquet, F., Yearsley, J., Jollivet, D., Broquet, T., 2016. Current hypotheses to explain genetic chaos under the sea. Curr Zool 62, 551-566. https://doi.org/10.1093/cz/zow094

1108

1109 Pereira, E., Mateus, C.S., Alves, M. J., Almeida, R., Pereira, J., Quintella, B. R., Almeida, P. R., 2023. Connectivity patterns and gene flow among Chelon ramada populations. Estuarine, Coastal and Shelf Science. 281: 108209, ISSN 0272-7714. https://doi.org/10.1016/j.ecss.2022.108209

1113

1114 Beaumont M. A., Nichols R. A., 1996. Evaluating loci for use in the genetic analysis of population structure, 1996. Proc. R. Soc. Lond. B 263, 1619-1626. https://doi.org/10.1098/rspb.1996.0237

1117

1118 Excoffier, L., Dupanloup, I., Huerta-Sánchez, E., Sousa, V.C., Foll, M., 2013. Robust Demographic Inference from Genomic and SNP Data. PLoS Genet 9, e1003905. https://doi.org/10.1371/journal.pgen.1003905

1121

1122 Excoffier, L., Lischer H. E. L., 2010. Arlequin suite ver 3.5: a new series of programs to perform population genetics analyses under Linux and Windows. Mol. Ecol. Res. 10(3): 564-567. https://doi.org/10.1111/j.1755-0998.2010.02847.x

1125

1126 Faria, R., Johannesson, K., Stankowski, S., 2021. Speciation in marine environments: Diving under the surface. J of Evolutionary Biology 34, 4-15. https://doi.org/10.1111/jeb.13756

1129

1130 Ferrini, V.L., Tivey, M.K., Carbotte, S.M., Martinez, F., Roman, C., 2008. Variable morphologic expression of volcanic, tectonic, and hydrothermal processes at six hydrothermal vent fields in the Lau back-arc basin: Hydrothermal Processes in Lau Back-arc Basin. Geochem. Geophys. Geosyst. 9, n/a-n/a. https://doi.org/10.1029/2008GC002047

1135

1136 Foll, M., Gaggiotti, O., 2008. A Genome-Scan Method to Identify Selected Loci Appropriate for Both Dominant and Codominant Markers: A Bayesian Perspective. Genetics 180, 977-993. https://doi.org/10.1534/genetics.108.092221

1140

1141 Folmer, O., Black, M., Hoeh, W., Lutz, R., And Vrijenhoek, R., 1994. DNA primers for amplification of mitochondrial cytochrome c oxidase subunit I from diverse metazoan invertebrates. Molecular Marine Biology and Biotechnology. 3294-299.
1144 François, O., Martins, H., Caye, K., Schoville, S.D., 2016. Controlling false discoveries in
1145 genome scans for selection. Mol Ecol 25, 454-469.

1146 https://doi.org/10.1111/mec.13513

1147

1148 Gagnaire, P.-A., Gaggiotti, O.E., 2016. Detecting polygenic selection in marine
1149 populations by combining population genomics and quantitative genetics
1150 approaches. Curr Zool 62, 603-616. https://doi.org/10.1093/cz/zow088

1151

1152 Gartman, A., Yücel, M., Madison, A.S., Chu, D.W., Ma, S., Janzen, C.P., Becker, E.L., Beinart,
1153 R.A., Girguis, P.R., Luther, G.W., 2011. Sulfide Oxidation across Diffuse Flow
1154 Zones of Hydrothermal Vents. Aquat Geochem 17, 583-601.

1155 https://doi.org/10.1007/s10498-011-9136-1

1156

1157 Gauthier, J., Mouden, C., Suchan, T., Alvarez, N., Arrigo, N., Riou, C., Lemaitre, C.,
1158 Peterlongo, P., 2020. DiscoSnp-RAD: de novo detection of small variants for
1159 RAD-Seq population genomics. PeerJ 8, e9291.

1160 https://doi.org/10.7717/peerj.9291

1161

1162 Gautier, M., Foucaud, J., Gharbi, K., Cézard, T., Galan, M., Loiseau, A., Thomson, M., Pudlo,
1163 P., Kerdelhué, C., Estoup, A., 2013. Estimation of population allele frequencies
1164 from next-generation sequencing data: pool-versus individual-based
1165 genotyping. Mol Ecol 22, 3766-3779. https://doi.org/10.1111/mec.12360

1166

1167 Gollner, S., Fontaneto, D., Martínez Arbizu, P., 2011. Molecular taxonomy confirms
1168 morphological classification of deep-sea hydrothermal vent copepods
1169 (Dirivultidae) and suggests broad physiological tolerance of species and
1170 frequent dispersal along ridges. Mar Biol 158, 221-231.

1171 https://doi.org/10.1007/s00227-010-1553-y

1172

1173 Gollner, S., Govenar, B., Fisher, C., Bright, M., 2015. Size matters at deep-sea
1174 hydrothermal vents: different diversity and habitat fidelity patterns of meio-
1175 and macrofauna. Mar. Ecol. Prog. Ser. 520, 57-66.

1176 https://doi.org/10.3354/meps11078

1177

1178 Gollner, S., Govenar, B., Martinez Arbizu, P., Mullineaux, L.S., Mills, S., Le Bris, N.,
1179 Weinbauer, M., Shank, T.M., Bright, M., 2020. Animal Community Dynamics at
1180 Senescent and Active Vents at the 9°N East Pacific Rise After a Volcanic
1181 Eruption. Front. Mar. Sci. 6, 832. https://doi.org/10.3389/fmars.2019.00832

1182

1183 Gollner, S., Ivanenko, V.N., Arbizu, P.M., Bright, M., 2010. Advances in Taxonomy,
1184 Ecology, and Biogeography of Dirivultidae (Copepoda) Associated with
1185 Chemosynthetic Environments in the Deep Sea. PLoS ONE 5, e9801.

1186 https://doi.org/10.1371/journal.pone.0009801

1187

1188 Gollner, S., Kaiser, S., Menzel, L., Jones, D.O.B., Brown, A., Mestre, N.C., Van Oevelen, D.,
1189 Menot, L., Colaço, A., Canals, M., Cuvelier, D., Durden, J.M., Gebruk, A., Egho, G.A.,
1190 Haeckel, M., Marcon, Y., Mevenkamp, L., Morato, T., Pham, C.K., Purser, A.,
1191 Sanchez-Vidal, A., Vanreusel, A., Vink, A., Martinez Arbizu, P., 2017. Resilience of benthic deep-sea fauna to mining activities. Marine Environmental Research 129, 76–101. https://doi.org/10.1016/j.marenvres.2017.04.010

1194

1195 Gollner, S., Stuckas, H., Kihara, T.C., Laurent, S., Kodami, S., Martinez Arbizu, P., 2016. Mitochondrial DNA Analyses Indicate High Diversity, Expansive Population Growth and High Genetic Connectivity of Vent Copepods (Dirivultidae) across Different Oceans. PLoS ONE 11, e0163776. https://doi.org/10.1371/journal.pone.0163776

1200

1201 Gutenkunst, R.N., Hernandez, R.D., Williamson, S.H., Bustamante, C.D., 2009. Inferring the Joint Demographic History of Multiple Populations from Multidimensional SNP Frequency Data. PLoS Genet 5, e1000695. https://doi.org/10.1371/journal.pgen.1000695

1205

1206 Hey, R.N., Massoth, G.J., Vrijenhoek, R.C., Rona, P.A., Lupton, J., Butterfield, D.A., 2006. Hydrothermal Vent Geology and Biology at Earth's Fastest Spreading Rates. Mar Geophys Res 27, 137–153. https://doi.org/10.1007/s11001-005-1887-x

1209

1210

1211 Holderegger, R., Wagner, H.H., 2006. A brief guide to Landscape Genetics. Landscape Ecol 21, 793–796. https://doi.org/10.1007/s10980-005-6058-6

1213

1214 Honza, E., 1995. Spreading mode of backarc basins in the western Pacific. Tectonophysics 251, 139–152. https://doi.org/10.1016/0040-1951(95)00054-2

1217

1218 Hourdez S. And Jollivet D (2019) CHUBACARC cruise, RV L'Atalante. https://doi.org/10.17600/18001111

1219

1220 Hourdez, S., Jollivet, D., 2020. Metazoan adaptation to deep-sea hydrothermal vents, in: Di Prisco, G., Edwards, H.G.M., Elster, J., Huiskes, A.H.L. (Eds.), Life in Extreme Environments. Cambridge University Press, pp. 42–67. https://doi.org/10.1017/9781108683319.004

1224

1225 Hu, E.-Z., Lan, X.-R., Liu, Z.-L., Gao, J., Niu, D.-K., 2022. A positive correlation between GC content and growth temperature in prokaryotes. BMC Genomics 23, 110. https://doi.org/10.1186/s12864-022-08353-7

1228

1229 Hu, R., Piotrowski, A.M., 2018. Neodymium isotope evidence for glacial-interglacial variability of deepwater transit time in the Pacific Ocean. Nat Commun 9, 4709. https://doi.org/10.1038/s41467-018-07079-z

1232

1233 Hurtado, L.A., Lutz, R.A., Vrijenhoek, R.C., 2004. Distinct patterns of genetic differentiation among annelids of eastern Pacific hydrothermal vents. Mol Ecol 13, 2603–2615. https://doi.org/10.1111/j.1365-294X.2004.02287.x

1236
1237 Ingels, J., Leduc, D., Zeppilli, D., Vanreusel, A., 2023. Deep-Sea Meiofauna—A World on Its Own or Deeply Connected? In: Giere, O., Schratzberger, M. (eds) New Horizons in Meiobenthos Research. Springer, Cham.
https://doi.org/10.1007/978-3-031-21622-0_8

1241
1242 Ivanenko, V.N., Martínez Arbizu, P., Stecher, J., 2007. Lecithotrophic nauplius of the
1243 family Dirivultidae (Copepoda; Siphonostomatoida) hatched on board over the
1244 Mid-Atlantic Ridge (5°S): Nauplius of Dirivultidae (Copepoda). Marine Ecology
1245 28, 49–53. https://doi.org/10.1111/j.1439-0485.2007.00142.x

1246
1247 John, C.R., Watson, D., Russ, D., Goldmann, K., Ehrenstein, M., Pitzalis, C., Lewis, M.,
1248 Barnes, M., 2020. M3C: Monte Carlo reference-based consensus clustering. Sci
1249 Rep 10, 1816. https://doi.org/10.1038/s41598-020-58766-1

1250
1251 Johnson, S.B., Warén, A., Tunnicliffe, V., Dover, C.V., Wheat, C.G., Schultz, T.F.,
1252 Vrijenhoek, R.C., 2015. Molecular taxonomy and naming of five cryptic species
1253 of *Alviniconcha* snails (Gastropoda: Abyssochrysoidea) from hydrothermal
1254 vents. Systematics and Biodiversity 13, 278–295.
1255 https://doi.org/10.1080/14772000.2014.970673

1256
1257 Jollivet, D., Chevaldonné, P., Planque, B., 1999. Polychaete dispersal in the eastern
1258 pacific. 2. A metapopulation model based on habitat shifts. Evolution. 53(4):
1259 1128-1141. https://doi.org/10.1111/j.1558-5646.1999.tb04527.x

1260
1261 Jollivet, D., Desbruyères, D., Bonhomme, F., Moraga, D., 1995. Genetic differentiation of
1262 deep-sea hydrothermal vent alvinellid populations (Annelida: Polychaeta) along
1263 the East Pacific Rise. Heredity 74, 376–391.
1264 https://doi.org/10.1038/hdy.1995.56

1265
1266 Kakee, T., 2020. Deep-sea mining legislation in Pacific Island countries: From the
1267 perspective of public participation in approval procedures. Marine Policy 117,
1268 103881. https://doi.org/10.1016/j.marpol.2020.103881

1269
1270 Keenan, K., McGinnity, P., Cross, T.F., Crozier, W.W., Prodöhl, P.A., 2013. diveRsity : An R
1271 package for the estimation and exploration of population genetics parameters
1272 and their associated errors. Methods Ecol Evol 4, 782–788.
1273 https://doi.org/10.1111/2041-210X.12067

1274
1275 Knaus, B.J., Grünwald, N.J., 2017. VCFR: a package to manipulate and visualize variant
1276 call format data in R. Mol Ecol Res, 17(1), 44–53. ISSN 757.
1277 https://dx.doi.org/10.1111/1755-0998.12549

1278
1279 Lavigne, F., Morin, J., Wassmer, P., Weller, O., Kula, T., Maea, A. V., Kelfoun, K., Mokadem,
1280 F., Paris, R., Malawani, M. N., Farai, A., Benbakkar, M., Saulnier-Copard, S., Vidal,
1281 C. M., Tu'l'afitu, T., Kitekei'aho, F., Trautmann, M., Gomez, C., 2021. Bridging
1282 Legends and Science: Field Evidence of a Large Tsunami that Affected the
1283

Kingdom of Tonga in the 15th Century. Front. Earth. Sci. 9:748755.

1284

https://doi.org/10.3389/feart.2021.748755

1285

1286

Leaché, A., Ogilvie, H.A., 2016. Bayes Factor Delimitation of Species (*with genomic data; BFD*): A Tutorial and Worked Example.

1288

1289

Legendre, P., Gauthier, O., 2014. Statistical methods for temporal and space-time analysis of community composition data <sup>sup</sup>/<sub>sup</sub>. Proc. R. Soc. B. 281, 20132728.

1291

https://doi.org/10.1098/rspb.2013.2728

1292

1293

Lelièvre, Y., Sarrazin, J., Marticorena, J., Schaal, G., Day, T., Legendre, P., Hourdez, S., Matabos, M., 2018. Biodiversity and trophic ecology of hydrothermal vent fauna associated with tubeworm assemblages on the Juan de Fuca Ridge.

1294

1295

Biogeosciences 15, 2629–2647. https://doi.org/10.5194/bg-15-2629-2018

1297

1298

Lenihan, H.S., Mills, S.W., Mullineaux, L.S., Peterson, C.H., Fisher, C.R., Micheli, F., 2008. Biotic interactions at hydrothermal vents: Recruitment inhibition by the mussel Bathymodiolus thermophilus. Deep Sea Research Part I: Oceanographic Research Papers 55, 1707–1717. https://doi.org/10.1016/j.dsr.2008.07.007

1302

1303

1304

Levin, L.A., Amon, D.J., Lily, H., 2020. Challenges to the sustainability of deep-seabed mining. Nat Sustain 3, 784–794. https://doi.org/10.1038/s41893-020-0558-x

1305

1306

1307

Limén, H., Stevens, C., Bourass, Z., Juniper, S., 2008. Trophic ecology of siphonostomatoid copepods at deep-sea hydrothermal vents in the northeast Pacific. Mar. Ecol. Prog. Ser. 359, 161–170.

1308

1309

https://doi.org/10.3354/meps07344

1310

1311

1312

Lotterhos, K.E., Whitlock, M.C., 2015. The relative power of genome scans to detect local adaptation depends on sampling design and statistical method. Mol Ecol 24, 1031–1046. https://doi.org/10.1111/mec.13100

1313

1314

1315

1316

Lowry, D.B., Hoban, S., Kelley, J.L., Lotterhos, K.E., Reed, L.K., Antolin, M.F., Storfer, A., 2017. Breaking RAD: an evaluation of the utility of restriction site-associated DNA sequencing for genome scans of adaptation. Mol Ecol Resour 17, 142–152.

1317

1318

https://doi.org/10.1111/1755-0998.12635

1319

1320

1321

Luikart, G., England, P.R., Tallmon, D., Jordan, S., Taberlet, P., 2003. The power and promise of population genomics: from genotyping to genome typing. Nat Rev Genet 4, 981–994. https://doi.org/10.1038/nrg1226

1322

1323

1324

1325

Luu, K., Bazin, E., Blum, M. G. B., 2017. pcadapt: an R package to perform genome scans for selection based on principal component analysis. Mol. Ecol. Res. 17(1): 67-77.

1326

https://doi.org/10.1111/1755-0998.12592

1327

1328

1329
1329 Mackenzie, C.L., Kent, F.E.A., Baxter, J.M., Gormley, K.S.G., Cassidy, A.J., Sanderson, W.G.,
1330 Porter, J.S., 2022. Genetic Connectivity and Diversity of a Protected, Habitat-
1331 Forming Species: Evidence Demonstrating the Need for Wider Environmental
1332 Protection and Integration of the Marine Protected Area Network. Front. Mar.
1333 Sci. 9, 772259. https://doi.org/10.3389/fmars.2022.772259
1334

1335 Martinez, F., Taylor, B., 2006. Modes of crustal accretion in back-arc basins: Inferences
1336 from the Lau Basin, in: Christie, D.M., Fisher, C.R., Lee, S.-M., Givens, S. (Eds.),
1337 Geophysical Monograph Series. American Geophysical Union, Washington, D. C.,
1338 pp. 5-30. https://doi.org/10.1029/166GM03
1339

1340 Mastretta-Yanes, A., Arrigo, N., Alvarez, N., Jorgensen, T.H., Piñero, D., Emerson, B.C.,
1341 2015. Restriction site-associated DNA sequencing, genotyping error estimation
1342 and de novo assembly optimization for population genetic inference. Mol Ecol
1343 Resour 15, 28-41. https://doi.org/10.1111/1755-0998.12291
1344

1345 Mees, J., Seys, J., 2023. Book of abstracts – VLIZ Marine Science Day, 1 March 2023,
1346 Bruges. Title: VLIZ Special Publication. Vlaams Instituut voor de Zee (VLIZ):
1347 Oostende. ISSN 1377-0950 Volume: 90 Issue: Navigation: vi + 112.
1348 https://doi.org/10.48470/41
1349

1350 Mensing, R., Stewart, M., Hannington, M., Baxter, A., and Mertmann, D., 2020. The
tectonic and volcanic evolution of the Mangatolu Triple Junction, EGU General
1351 Assembly 2020, Online, 4-8 May 2020, EGU2020-19897.
1352 https://doi.org/10.5194/egusphere-egu2020-19897
1353

1354 Micheli, F., Peterson, C.H., Mullineaux, L.S., Fisher, C.R., Mills, S.W., Sancho, G., Johnson,
1355 G.A., Lenihan, H.S., 2002b. PREDATION STRUCTURES COMMUNITIES AT DEEP-
1356 SEA HYDROTHERMAL VENTS. Ecological Monographs 72, 365-382.
1357 https://doi.org/10.1890/0012-9615(2002)072[0365:PSCADS]2.0.CO;2
1358

1359 Milano, I., Babbucci, M., Cariani, A., Atanassova, M., Bekkevold, D., Carvalho, G.R.,
1360 Espiñeira, M., Fiorentino, F., Garofalo, G., Geffen, A.J., Hansen, J.H., Helyar, S.J.,
1361 Nielsen, E.E., Ogden, R., Patarnello, T., Stagioni, M., Tinti, F. and Bargelloni, L.,
1362 2014. Outlier SNP markers reveal fine-scale genetic structuring across European
1363 hake populations (Merluccius merluccius). Mol Ecol. 23: 118-135.
1364 https://doi.org/10.1111/mec.12568
1365

1366 Miller, K.A., Thompson, K.F., Johnston, P., Santillo, D., 2018. An Overview of Seabed
1367 Mining Including the Current State of Development, Environmental Impacts,
1368 and Knowledge Gaps. Front. Mar. Sci. 4, 418.
1369 https://doi.org/10.3389/fmars.2017.00418
1370

1371 Mitarai, S., Siegel, D.A., Watson, J.R., Dong, C., McWilliams, J.C., 2009. Quantifying
1372 connectivity in the coastal ocean with application to the Southern California
1373 Bight. J. Geophys. Res. 114, C10026. https://doi.org/10.1029/2008JC005166
1374

1375
1376 Mitarai, S., Watanabe, H., Nakajima, Y., Shchepetkin, A.F., McWilliams, J.C., 2016. Quantifying dispersal from hydrothermal vent fields in the western Pacific Ocean. Proc. Natl. Acad. Sci. U.S.A. 113, 2976-2981.
https://doi.org/10.1073/pnas.1518395113

1380 Monzier, M., Robin, C., Eissen, J-P, 1994. Kuwae (≈ 1425 A.D.): the forgotten caldera. J. volcanol. Geotherm. Res. 59, 207-218.
https://doi.org/10.1016/0377-0273(94)90091-4

1384 Mouchi, V., Pecheyran, C., Claverie, F., Cathalot, C., Matabos, M., Germain, Y., Rouxel, O., Jollivet, D., Broquet, T., Comtet, T., 2023. A step towards measuring connectivity in the deep-sea: elemental fingerprints of mollusk larval shells discriminate hydrothermal vent sites. Biogeosciences, 21, 145-160. BioRxiv.
https://doi.org/10.1101/2023.01.03.522618

1390 Mullineaux, L.S., Adams, D.K., Mills, S.W., Beaulieu, S.E., 2010. Larvae from afar colonize deep-sea hydrothermal vents after a catastrophic eruption. Proc. Natl. Acad. Sci. U.S.A. 107, 7829-7834.
https://doi.org/10.1073/pnas.0913187107

1394 Mullineaux, L.S., Metaxas, A., Beaulieu, S.E., Bright, M., Gollner, S., Grupe, B.M., Herrera, S., Kellner, J.B., Levin, L.A., Mitarai, S., Neubert, M.G., Thurnherr, A.M., Tunnicliffe, V., Watanabe, H.K., Won, Y.-J., 2018. Exploring the Ecology of Deep-Sea Hydrothermal Vents in a Metacommunity Framework. Front. Mar. Sci. 5, 49.
https://doi.org/10.3389/fmars.2018.00049

1400 Murdock, S.A., Tunnicliffe, V., Boschen-Rose, R.E., Juniper, S.K., 2021. Emergent "core communities" of microbes, meiofauna and macrofauna at hydrothermal vents. ISME COMMUN. 1, 27.
https://doi.org/10.1038/s43705-021-00031-1

1404 Nakasugi, F., Shimanaga, M., Nomaki, H., Watanabe, H.K., Kitahashi, T., Motomura, Y., Iseda, K., 2021. Simple harpacticoid composition observed at deep hydrothermal vent sites on sea knoll calderas in the North-west Pacific. J. Mar. Biol. Ass. 101, 947-956.
https://doi.org/10.1017/S0025315421000874

1409 Nomaki, H., Kawatani, K., Motomura, Y., Tame, A., Uyeno, D., Ogawa, N.O., Ohkouchi, N., Shimanaga, M., 2023. Bacterivory of the hydrothermal-vent-specific copepod *Stygiopontius senokuchiae* (Dirivultidae, Siphonostomatoida) from copepodite through adult stages. J. Mar. Biol. Ass. 103, e21.
https://doi.org/10.1017/S0025315423000139

1415 Nomaki, H., Uejima, Y., Ogawa, N., Yamane, M., Watanabe, H., Senokuchi, R., Bernhard, J., Kitahashi, T., Miyairi, Y., Yokoyama, Y., Ohkouchi, N., Shimanaga, M., 2019. Nutritional sources of meio- and macrofauna at hydrothermal vents and adjacent areas: natural-abundance radiocarbon and stable isotope analyses. Mar. Ecol. Prog. Ser. 622, 49-65.
https://doi.org/10.3354/meps13053

1421
1422 O'Connor, M.I., Bruno, J.F., Gaines, S.D., Halpern, B.S., Lester, S.E., Kinlan, B.P., Weiss, J.M., 2007. Temperature control of larval dispersal and the implications for marine ecology, evolution, and conservation. Proc. Natl. Acad. Sci. U.S.A. 104, 1266-1271. https://doi.org/10.1073/pnas.0603422104

1426

1427 Olive, J.-A., Dublanchet, P., 2020. Controls on the magmatic fraction of extension at mid-ocean ridges. Earth and Planetary Science Letters 549, 116541. https://doi.org/10.1016/j.epsl.2020.116541

1430

1431 Paris, J.R., Stevens, J.R., Catchen, J.M., 2017. Lost in parameter space: a road map for STACKS. Methods Ecol Evol 8, 1360-1373. https://doi.org/10.1111/2041-210X.12775

1432

1433

1434 Perez, M., Sun, J., Xu, Q., Qian, P.-Y., 2021. Structure and Connectivity of Hydrothermal Vent Communities Along the Mid-Ocean Ridges in the West Indian Ocean: A Review. Front. Mar. Sci. 8, 744874. https://doi.org/10.3389/fmars.2021.744874

1435

1436

1437

1438

1439 Petersen, S., Krätschell, A., Augustin, N., Jamieson, J., Hein, J.R., Hannington, M.D., 2016. News from the seabed - Geological characteristics and resource potential of deep-sea mineral resources. Marine Policy 70, 175-187. https://doi.org/10.1016/j.marpol.2016.03.012

1440

1441

1442

1443

1444 Pineda, J., Porri, F., Starczak, V., Blythe, J., 2010. Causes of decoupling between larval supply and settlement and consequences for understanding recruitment and population connectivity. Journal of Experimental Marine Biology and Ecology 392, 9-21. https://doi.org/10.1016/j.jembe.2010.04.008

1445

1446

1447

1448

1449 Poitrimol C, Thiébaut É, Daguin-Thiébaut C, Le Port A-S, Ballenghien M, Tran Lu Y A, Jollivet D, Hourdez S, Matabos M. Contrasted phylogeographic patterns of hydrothermal vent gastropods along South West Pacific: Woodlark Basin, a possible contact zone and/or stepping-stone. PLoS ONE. 2022 Oct 5; 17(10): e0275638. https://doi.org/10.1371/journal.pone.0275638

1450

1451

1452

1453

1454

1455 Plouviez, S., LaBella, A.L., Weisrock, D.W., Von Meijenfeldt, F.A.B., Ball, B., Neigel, J.E., Van Dover, C.L., 2019. Amplicon sequencing of 42 nuclear loci supports directional gene flow between South Pacific populations of a hydrothermal vent limpet. Ecol Evol 9, 6568-6580. https://doi.org/10.1002/ece3.5235

1456

1457

1458

1459

1460 Plouviez, S., Schultz, T.F., McGinnis, G., Minshall, H., Rudder, M., Van Dover, C.L., 2013. Genetic diversity of hydrothermal-vent barnacles in Manus Basin. Deep Sea Research Part I: Oceanographic Research Papers 82, 73-79. https://doi.org/10.1016/j.dsr.2013.08.004

1461

1462

1463

1464

1465 Portanier, E., Nicolle, A., Rath, W., Monnet, L., Le Goff, G., Le Port, A.-S., Daguin-Thiébaut, C., Morrison, C.L., Cunha, M.R., Betters, M., Young, C.M., Van Dover, C.L., Biastoch, A., Thiébaut, E., Jollivet, D., 2023. Coupling large-spatial scale larval dispersal modelling with barcoding to refine the amphi-Atlantic connectivity hypothesis

1466

1467

1468
1469 in deep-sea seep mussels. Front. Mar. Sci. 10, 1122124.

1470 https://doi.org/10.3389/fmars.2023.1122124

1471
1472 Posavi, M., Gelembiuk, G.W., Larget, B., Lee, C.E., 2014. Testing for beneficial reversal of dominance during salinity shifts in the invasive copepod *Eurytemora affinis*, and implications for the maintenance of genetic variation. Evolution 68, 3166–3183. https://doi.org/10.1111/evo.12502

1476
1477 Riginos, C., Crandall, E.D., Liggins, L., Bongaerts, P., Treml, E.A., 2016. Navigating the currents of seascape genomics: how spatial analyses can augment population genomic studies. Curr Zool 62, 581–601. https://doi.org/10.1093/cz/zow067

1480
1481 Rochette, N.C., Catchen, J.M., 2017. Deriving genotypes from RAD-seq short-read data using Stacks. Nat Protoc 12, 2640–2659.
1482 https://doi.org/10.1038/nprot.2017.123

1484
1485 Rochette, N.C., Rivera-Colón, A.G., Catchen, J.M., 2019. Stacks 2: Analytical methods for paired-end sequencing improve RADseq-based population genomics. Mol Ecol 28, 4737–4754. https://doi.org/10.1111/mec.15253

1488
1489 Rougeux, C., Bernatchez, L., Gagnaire, P.-A., 2017. Modeling the Multiple Facets of Speciation-with-Gene flow toward Inferring the Divergence History of Lake Whitefish Species Pairs (Coregonus clupeaformis). Genome Biology and Evolution 9, 2057–2074. https://doi.org/10.1093/gbe/evx150

1493
1494 Ruellan, E., Delteil, J., Wright, I., Matsumoto, T., 2003. From rifting to active spreading in the Lau Basin - Havre Trough backarc system (SW Pacific): Locking/unlocking induced by seamount chain subduction: LAU BASIN-HAVRE TROUGH SYSTEM. Geochem. Geophys. Geosyst. 4, n/a-n/a.
1497
1498 https://doi.org/10.1029/2001GC000261

1499
1500 Siberchicot, A., Julien-Laferrière, A., Dufour, A., Thioulouse, J., Dray, S., 2017. adegraphics: An S4 Lattice-Based Package for the Representation of Multivariate Data. The R Journal. 9 (2):198-212. https://doi.org/10.32614/RJ-2017-042

1504
1505 Siegel, D.A., Kinlan, B.P., Gaylord, B. and Gaines, S.D., 2003. Lagrangian descriptions of marine larval dispersion. Marine Ecology Progress Series, 260, pp.83-96.
1507 https://doi.org/10.3354/meps260083

1508
1509 Simons, E., Speer, K., Thurnherr, A.M., 2019. Deep circulation in the Lau Basin and Havre Trough of the western South Pacific Ocean from floats and hydrography. j
1510
1511 mar res 77, 353–374. https://doi.org/10.1357/002224019833406150

1512
1513 Smagorinsky, J., 1963. GENERAL CIRCULATION EXPERIMENTS WITH THE PRIMITIVE EQUATIONS: I. THE BASIC EXPERIMENT*. Mon. Wea. Rev. 91, 99–164.
1514 https://doi.org/10.1175/1520-0493(1963)091<0099:GCEWTP>2.3.CO;2
1515
1516

1517

1518

1519

1520

1521

1522

1523

1524

1525

1526

1527

1528

1529

1530

1531

1532

1533

1534

1535

1536

1537

1538

1539

1540

1541

1542

1543

1544

1545

1546

1547

1548

1549

1550

1551

1552

1553

1554

1555

1556

1557

1558

1559

1560

1561

1562

Speer, K., Thurnherr, A., 2012b. The Lau Basin Float Experiment (LAUB-FLEX). oceanog 25, 284–285. https://doi.org/10.5670/oceanog.2012.27

Stewart, M.S., Hannington, M.D., Emberley, J., Baxter, A.T., Krätschell, A., Petersen, S., Brandl, P.A., Anderson, M.O., Mercier-Langevin, P., Mensing, R., Breker, K., Fassbender, M.L., 2022. A new geological map of the Lau Basin (southwestern Pacific Ocean) reveals crustal growth processes in arc-backarc systems. Geosphere. https://doi.org/10.1130/GES02340.1

Storey, J.D., 2003. The positive false discovery rate: a Bayesian interpretation and the q-value. Ann. Statist. 31. https://doi.org/10.1214/aos/1074290335

Sunde, J., Yıldırım, Y., Tibblin, P., Forsman, A., 2020. Comparing the Performance of Microsatellites and RADseq in Population Genetic Studies: Analysis of Data for Pike (Esox lucius) and a Synthesis of Previous Studies. Front. Genet. 11, 218. https://doi.org/10.3389/fgene.2020.00218

Sundqvist, L., Keenan, K., Zackrisson, M., Prodöhl, P., Kleinhans, D., 2016. Directional genetic differentiation and relative migration. Ecol Evol 6, 3461–3475. https://doi.org/10.1002/ece3.2096

Teixeira, S., Cambon-Bonavita, M.-A., Serrão, E.A., Desbruyéres, D., Arnaud-Haond, S., 2011. Recent population expansion and connectivity in the hydrothermal shrimp Rimicaris exoculata along the Mid-Atlantic Ridge: Genetic diversity of a hydrothermal vent shrimp. Journal of Biogeography 38, 564–574. https://doi.org/10.1111/j.1365-2699.2010.02408.x

Teixeira, S., Serrão, E.A., Arnaud-Haond, S., 2012. Panmixia in a Fragmented and Unstable Environment: The Hydrothermal Shrimp Rimicaris exoculata Disperses Extensively along the Mid-Atlantic Ridge. PLoS ONE 7, e38521. https://doi.org/10.1371/journal.pone.0038521

Tepper, B., Bradley, B.P., 1989. Temporal Changes in a Natural Population of Copepods. The Biological Bulletin 176, 32–40. https://doi.org/10.2307/1541886

Thaler, A.D., Plouviez, S., Saleu, W., Alei, F., Jacobson, A., Boyle, E.A., Schultz, T.F., Carlsson, J., Van Dover, C.L., 2014. Comparative Population Structure of Two Deep-Sea Hydrothermal-Vent-Associated Decapods (Chorocaris sp. 2 and Munidopsis lauensis) from Southwestern Pacific Back-Arc Basins. PLoS ONE 9, e101345. https://doi.org/10.1371/journal.pone.0101345

Thomas, G.W.C., Dohmen, E., Hughes, D.S.T., Murali, S.C., Poelchau, M., Glastad, K., Anstead, C.A., Ayoub, N.A., Batterham, P., Bellair, M., Binford, G.J., Chao, H., Chen, Y.H., Childers, C., Dinh, H., Doddapaneni, H.V., Duan, J.J., Dugan, S., Esposito, L.A., Friedrich, M., Garb, J., Gasser, R.B., Goodisman, M.A.D., Gundersen-Rindal, D.E., Han, Y., Handler, A.M., Hatakeyama, M., Hering, L., Hunter, W.B., Ioannidis, P.,
1563 Jayaseelan, J.C., Kalra, D., Khila, A., Korhonen, P.K., Lee, C.E., Lee, S.L., Li, Y.,
1564 Lindsey, A.R.I., Mayer, G., McGregor, A.P., McKenna, D.D., Misof, B., Munidasa, M.,
1565 Munoz-Torres, M., Muzny, D.M., Niehuis, O., Osuji-Lacy, N., Palli, S.R., Panfilio,
1566 K.A., Pechmann, M., Perry, T., Peters, R.S., Poynton, H.C., Prpic, N.-M., Qu, J.,
1567 Rotenberg, D., Schal, C., Schoville, S.D., Scully, E.D., Skinner, E., Sloan, D.B.,
1568 Stouthamer, R., Strand, M.R., Szucsich, N.U., Wijeratne, A., Young, N.D., Zattara,
1569 E.E., Benoit, J.B., Zdobnov, E.M., Pfeffer, M.E., Hackett, K.J., Werren, J.H.,
1570 Worley, K.C., Gibbs, R.A., Chipman, A.D., Waterhouse, R.M., Bornberg-Bauer, E.,
1571 Hahn, M.W., Richards, S., 2020. Gene content evolution in the arthropods.
1572 Genome Biol 21, 15. https://doi.org/10.1186/s13059-019-1925-7

1573

1574 Tran Lu Y, A., Ruault, S., Daguin-Thiébaut, C., Castel, J., Bierne, N., Broquet, T., Wincker, P., Perdereau, A., Arnaud-Haond, S., Gagnaire, P., Jollivet, D., Hourdez, S., Bonhomme, F., 2022. Subtle limits to connectivity revealed by outlier loci within two divergent metapopulations of the deep-sea hydrothermal gastropod *Ifremeria nautilei*. Molecular Ecology 31, 2796-2813.
1575
1576
1577
1578
1579 https://doi.org/10.1111/mec.16430

1580

1581 Tripp, E.A., Tsai, Y.E., Zhuang, Y., Dexter, K.G., 2017. RAD seq dataset with 90% missing data fully resolves recent radiation of *Petalidium* (Acanthaceae) in the ultra-arid deserts of Namibia. *Ecol Evol* 7, 7920-7936.
1582
1583
1584 https://doi.org/10.1002/ece3.3274

1585

1586 Tsurumi, M., De Graaf, R.C., Tunnicliffe, V., 2003. Distributional and Biological Aspects of Copepods at Hydrothermal Vents on the Juan de Fuca Ridge, north-east Pacific Ocean. J. Mar. Biol. Ass. 83, 469-477.
1587
1588
1589 https://doi.org/10.1017/S0025315403007367h

1590

1591 Van Dover, C.L., 2014. Impacts of anthropogenic disturbances at deep-sea hydrothermal vent ecosystems: A review. Marine Environmental Research 102, 59-72. https://doi.org/10.1016/j.marenvres.2014.03.008

1594

1595 Van Dover, C.L., Arnaud-Haond, S., Gianni, M., Helmreich, S., Huber, J.A., Jaeckel, A.L., Metaxas, A., Pendleton, L.H., Petersen, S., Ramirez-Llodra, E., Steinberg, P.E., Tunnicliffe, V., Yamamoto, H., 2018. Scientific rationale and international obligations for protection of active hydrothermal vent ecosystems from deep-sea mining. Marine Policy 90, 20-28.
1596
1597
1598
1599
1600 https://doi.org/10.1016/j.marpol.2018.01.020

1601

1602 Virtanen, E.A., Moilanen, A., Viitasalo, M., 2020. Marine connectivity in spatial conservation planning: analogues from the terrestrial realm. Landscape Ecol 35, 1021-1034. https://doi.org/10.1007/s10980-020-00997-8

1605

1606 Vrijenhoek, R. C., 2010. Genetic diversity and connectivity of deep-sea hydrothermal vent metapopulations. Mol. Ecol. 19: 4391-4411. https://doi.org/10.1111/j.1365-294X.2010.04789.x

1609
1610 Wang, S., Meyer, E., McKay, J.K., Matz, M.V., 2012. 2b-RAD: a simple and flexible method for genome-wide genotyping. Nat Methods 9, 808–810.
1611 https://doi.org/10.1038/nmeth.2023

1613

1614

1615 Watanabe, H.K., Senokuchi, R., Nomaki, H., Kitahashi, T., Uyeno, D., Shimanaga, M., 2021. Distribution and Genetic Divergence of Deep-Sea Hydrothermal Vent Copepods (Dirivultidae: Siphonostomatoida: Copepoda) in the Northwestern Pacific. Zoological Science 38. https://doi.org/10.2108/zs200153

1619

1620

1621 Won, Y., Hallam, S.J., O'Mullan, G.D., Vrijenhoek, R.C., 2003. Cytonuclear disequilibrium in a hybrid zone involving deep-sea hydrothermal vent mussels of the genus *Bathymodiolus*. Molecular Ecology 12, 3185–3190.
1622 https://doi.org/10.1046/j.1365-294X.2003.01974.x

1625

1626 Xu, T., Wang, Y., Sun, J., Chen, C., Watanabe, H.K., Chen, J., Qian, P.-Y., Qiu, J.-W., 2021. Hidden Historical Habitat-Linked Population Divergence and Contemporary Gene Flow of a Deep-Sea Patellogastropod Limpet. Molecular Biology and Evolution 38, 5640–5654. https://doi.org/10.1093/molbev/msab278

1630

1631 Xuereb, A., Benestan, L., Normandeau, É., Daigle, R.M., Curtis, J.M.R., Bernatchez, L., Fortin, M.-J., 2018. Asymmetric oceanographic processes mediate connectivity and population genetic structure, as revealed by RADseq, in a highly dispersive marine invertebrate (*Parastichopus californicus*). Mol Ecol 27, 2347–2364.
1632
1633
1634
1635 https://doi.org/10.1111/mec.14589

1636

1637

1638 Yahagi, T., Fukumori, H., Warén, A., & Kano, Y., 2019. Population connectivity of hydrothermal-vent limpets along the northern Mid-Atlantic Ridge (Gastropoda: Neritimorpha: Phenacolepadidae). Journal of the Marine Biological Association of the United Kingdom, 99(1), 179-185.
1639
1640
1641
1642 https://doi.org/10.1017/S0025315417001898

1643

1644 Yahagi, T., Kayama Watanabe, H., Kojima, S., Kano, Y., 2017. Do larvae from deep-sea hydrothermal vents disperse in surface waters? Ecology 98, 1524–1534.
1645 https://doi.org/10.1002/ecy.1800

1647

1648 Yahagi, T., Thaler, A.D., Van Dover, C.L., Kano, Y., 2020. Population connectivity of the hydrothermal-vent limpet Shinkailepas tollmanni in the Southwest Pacific (Gastropoda: Neritimorpha: Phenacolepadidae). PLoS ONE 15, e0239784.
1649
1650
1651 https://doi.org/10.1371/journal.pone.0239784

1652

1653 Yearsley, J.M., Salmanidou, D.M., Carlsson, J., Burns, D., Van Dover, C.L., 2020. Biophysical models of persistent connectivity and barriers on the northern Mid-Atlantic Ridge. Deep Sea Research Part II: Topical Studies in Oceanography 180, 104819. https://doi.org/10.1016/j.dsr2.2020.104819

1656
<table><caption>Table 1: Location of sampling sites (latitude and longitude) with corresponding depth (m) at time of sampling. Spreading rate values were taken from Stewart et al. (2022).</caption><thead><tr><th>Basin</th><th>Site</th><th>Latitude</th><th>Longitude</th><th>Depth (m)</th><th>Spreading Rate</th><th>Substrate</th><th>Sample size</th><th>Year sampled</th><th>References</th></tr></thead><tbody><tr><td>ELSC</td><td>Tu'i Malila</td><td>-21.99</td><td>-176.57</td><td>1888</td><td>48</td><td>Basalt</td><td>24</td><td>2019</td><td>Hourdez and Jollivet (2019)<br>Stewart et al. (2022)</td></tr><tr><td>ELSC</td><td>ABE</td><td>-20.76</td><td>-176.19</td><td>2155</td><td>69</td><td>Andesitic</td><td>24</td><td>2016</td><td>Beinart et al. (2018)<br>Stewart et al. (2022)</td></tr><tr><td>ELSC</td><td>ABE</td><td>-20.77</td><td>-176.20</td><td>2153</td><td>69</td><td>Andesitic</td><td>11</td><td>2019</td><td>Hourdez and Jollivet (2019)<br>Stewart et al. (2022)</td></tr><tr><td>ELSC</td><td>Tahi Moana</td><td>-20.05</td><td>-176.13</td><td>2273</td><td>69</td><td>Andesitic</td><td>15</td><td>2016</td><td>Beinart et al. (2018)<br>Stewart et al. (2022)</td></tr><tr><td>ELSC</td><td>Mangatolo</td><td>-15.41</td><td>185.35</td><td>2039</td><td>37</td><td>Andesitic</td><td>64</td><td>2019</td><td>Hourdez and Jollivet (2019)<br>Stewart et al. (2022)</td></tr></tbody></table>

1660
1661
1662
1663

<table><caption>Table 2: Analysis of Molecular Variance results of pairwise differences between populations, indicating the pairs, the Fst value for each pair and the *p-value associated with each Fst calculation. Significance codes: (***p < 0.001, **p < 0.01, *p < 0.05). Observed heterozygosity for each population was calculated (Hetobs).</caption><thead><tr><th>Population</th><th>Het<sub>obs</sub></th><th>Tahi Moana</th><th>Mangatolo</th><th>Tui Malila</th><th>ABE (2016)</th><th>ABE (2019)</th></tr></thead><tbody><tr><td>Tahi Moana</td><td>0.12</td><td>0</td><td></td><td></td><td></td><td></td></tr><tr><td>Mangatolo</td><td>0.24</td><td>0.24*</td><td>0</td><td></td><td></td><td></td></tr><tr><td>Tui Malila</td><td>0.18</td><td>0.24*</td><td>0.12*</td><td>0</td><td></td><td></td></tr><tr><td>ABE (2016)</td><td>0.12</td><td>0.00</td><td>0.26*</td><td>0.25*</td><td>0</td><td></td></tr><tr><td>ABE (2019)</td><td>0.08</td><td>0.11*</td><td>0.05*</td><td>0.08*</td><td>0.07*</td><td>0</td></tr></tbody></table>

1668
1669
1670
1671
1672
1673
1674
1675
1676
1677
Table 3: Parameters estimated from $\partial a_i$ for the best model (SC2N2mG) with their standard deviations (SD) estimated using the Fisher Information Matrix (FIM), and with parameters describing effective population size (2N), migration rate (2m), and population growth (G) between populations within the Lau basin.

<table><thead><tr><th colspan="3">North-mid</th><th colspan="3">North-south</th></tr><tr><th>Parameter</th><th>Value</th><th>SD</th><th>Parameter</th><th>Value</th><th>SD</th></tr></thead><tbody><tr><td>Nu1 (mid)</td><td>2.41</td><td>0.63</td><td>Nu1 (north)</td><td>1.40</td><td>0.76</td></tr><tr><td>Nu2 (north)</td><td>0.80</td><td>0.61</td><td>Nu2 (south)</td><td>1.24</td><td>1.16</td></tr><tr><td>b1</td><td>0.31</td><td>0.07</td><td>b1</td><td>0.74</td><td>0.21</td></tr><tr><td>b2</td><td>0.63</td><td>0.60</td><td>b2</td><td>1.03</td><td>0.19</td></tr><tr><td>hrf</td><td>0.14</td><td>0.07</td><td>hrf</td><td>0.11</td><td>0.06</td></tr><tr><td>T<sub>s</sub></td><td>4.91</td><td>1.47</td><td>T<sub>s</sub></td><td>7.92</td><td>5.80</td></tr><tr><td>T<sub>sc</sub></td><td>0.17</td><td>0.01</td><td>T<sub>sc</sub></td><td>0.26</td><td>0.07</td></tr><tr><td>m1<-2</td><td>0.68</td><td>0.10</td><td>m1<-2</td><td>5.24</td><td>1.16</td></tr><tr><td>m2<-1</td><td>5.97</td><td>0.40</td><td>m2<-1</td><td>6.24</td><td>2.35</td></tr><tr><td>me1<-2</td><td>0.98</td><td>0.09</td><td>me1<-2</td><td>4.64</td><td>1.39</td></tr><tr><td>me2<-1</td><td>0.51</td><td>0.13</td><td>me2<-1</td><td>1.03</td><td>0.64</td></tr><tr><td>P</td><td>0.74</td><td>0.04</td><td>P</td><td>0.03</td><td>0.23</td></tr><tr><td>Q</td><td>0.17</td><td>0.12</td><td>Q</td><td>0.49</td><td>0.28</td></tr><tr><td>Theta</td><td>91.82</td><td>23.55</td><td>Theta</td><td>155.82</td><td>93.60</td></tr><tr><td>Time to population split (years)</td><td>119,422</td><td>35,764</td><td>Time to the population split (years)</td><td>325,661</td><td>239,304</td></tr><tr><td>Time to secondary contact (years)</td><td>4,050</td><td>274</td><td>Time to the secondary contact (years)</td><td>10,762</td><td>3,231</td></tr></tbody></table>
1684

<figure><img src="image_2.png" /></figure>

Figure 1: Map of sampling sites from 2016, and 2019 within the Lau Basin (Mangatolo, Tahi Moana, ABE, and Tu'i Malila). A dorsal confocal microscopic image of an *S. lauensis* female from the ABE hydrothermal vent field in the Lau Basin. Image taken at the German Center for Marine Biodiversity Research by Dr. Jimmy Bernot.

1685

1686

1687

1688

1689

1690
<figure><img src="image_3.png" /></figure>

<figure><img src="image_4.png" /><figcaption>Figure 2: Structure of populations of *S. lauensis* in the Lau Basin. A) non-Metric Multidimensional Scaling (nMDS), where data is coloured by vent site and shapes denote the sampling year (open circles = 2019, hatched square = 2016). Ellipsoids represent 95 % confidence intervals. B) ADMIXTURE plots for each individual coloured by ancestry and grouped by sampling site/population for the optimal K = 3.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Figure 3: Venn diagram of outlier SNPs (unique and shared) found by the different outlier detection software Bayesian, Arlequin, Sambada, and PCAdapt.</figcaption></figure>
1709

<figure><img src="image_6.png" /><figcaption>Figure 4: A) Joint allele frequency spectrum (JAFS) between North (Mangatolo) and Mid (pooled Tahi Moana and ABE) vent sites within the Lau Basin. B) Simulated JAFS under the SC2N2mG model. The log scale indicates the density of SNPs in each frequency class (shared between the mid and north). C) Residuals of the fit of the simulated model on the data. D) Schematic of the fitted model (N, population size; b, population growth factor; hrf, the Hill-Robertson factor, which simulates linked selection; m, unrestricted migration rate; me, restricted migration rate, which simulates barrier loci). E) Heat-map of the best Akaike information criterion (AIC) value for each parameter combination (population expansion or contraction [G], effect of linked selection [2N] and heterogeneous migration [2m]) and demographic model (strict isolation [SI], isolation with migration [IM], ancient migration [AM], and secondary contact [SC]).</figcaption></figure>

1710

1711

1712

1713

1714

1715

1716

1717

1718

1719

<figure><img src="image_7.png" /><figcaption>Figure 5: A) Joint allele frequency spectrum (JAFS) between North (Mangatolo) and South (Tu'i Malila) vent sites within the Lau Basin. B) Simulated JAFS under the SC2N model. The log scale indicates the density of SNPs in each frequency class. C) Residuals of the fit of the simulated model on the data. D) Schematic of the fitted model (N, population size; b, population growth factor; hrf, the Hill-Robertson factor, which simulates linked selection; m, unrestricted migration rate; me, restricted migration rate, which simulates barrier loci). E) Heat-map of the best Akaike information criterion (AIC) value for each parameter combination (population expansion or contraction [G], effect of linked selection [2N] and heterogeneous migration [2m]) and demographic model (strict isolation [SI], isolation with migration [IM], ancient migration [AM], and secondary contact [SC]).</figcaption></figure>

1720

1721

1722

1723

1724

1725

1726

1727

1728

1729

1730

1731

1732
<figure><img src="image_8.png" /></figure>

Figure 6: A) Significant relative migration (effective number of migrants, nM) calculated using the Nm statistic and 1000 bootstraps in *DivMigrate* (empirical data = black arrows), and larval dispersal model (LDM) simulations (blue arrows). Blue dashed line denotes a 100 m radium (not to scale around each vent). The thickness of the arrows and size of the arrow heads in each instance denotes the strength of directionality. SR = self-recruitment. Gaussian point distribution probabilities (probability density) for B) Mangatolo, C) Tahi Moana, D) ABE, and E) Tu'i Malila. Blue dashed White triangles denote existing (and possibly stepping stone) vents within the basin from which data was not collected.

1733

1734

1735

1736

1737

1738

1739

1740

1741

1742

1743