Origins of the Greenland shark (Somniosus microcephalus): Impacts of ice‐olation and introgression

Abstract Herein, we use genetic data from 277 sleeper sharks to perform coalescent‐based modeling to test the hypothesis of early Quaternary emergence of the Greenland shark (Somniosus microcephalus) from ancestral sleeper sharks in the Canadian Arctic‐Subarctic region. Our results show that morphologically cryptic somniosids S. microcephalus and Somniosus pacificus can be genetically distinguished using combined mitochondrial and nuclear DNA markers. Our data confirm the presence of genetically admixed individuals in the Canadian Arctic and sub‐Arctic, and temperate Eastern Atlantic regions, suggesting introgressive hybridization upon secondary contact following the initial species divergence. Conservative substitution rates fitted to an Isolation with Migration (IM) model indicate a likely species divergence time of 2.34 Ma, using the mitochondrial sequence DNA, which in conjunction with the geographic distribution of admixtures and Pacific signatures likely indicates speciation associated with processes other than the closing of the Isthmus of Panama. This time span coincides with further planetary cooling in the early Quaternary period followed by the onset of oscillating glacial‐interglacial cycles. We propose that the initial S. microcephalus–S. pacificus split, and subsequent hybridization events, were likely associated with the onset of Pleistocene glacial oscillations, whereby fluctuating sea levels constrained connectivity among Arctic oceanic basins, Arctic marginal seas, and the North Atlantic Ocean. Our data demonstrates support for the evolutionary consequences of oscillatory vicariance via transient oceanic isolation with subsequent secondary contact associated with fluctuating sea levels throughout the Quaternary period—which may serve as a model for the origins of Arctic marine fauna on a broad taxonomic scale.

In Arctic marine environments, the evolutionary histories and population demographics of numerous endemic species appear to be influenced strongly by fluctuations in the extent of sea ice associated with cooling glacial and warming interglacial cycles since the onset of the geologic Quaternary period roughly 2.6 Ma ago (e.g., O' Corry-Crowe et al., 2010). Genetic structure in beluga whales Delphinapterus leucas (O'Corry-Crowe, Suydam, Rosenberg, Frost, & Dizon, 1997;O'Corry-Crowe et al., 2010) and walruses Odobenus rosmarus (Andersen & Born, 2000) provides strong indicators that ice effects during the Pleistocene (0.0117-2.588 Ma ago), coupled with additional ecological factors such as large-scale temporal variations in levels of primary productivity, can leave lasting genetic imprints in a species' evolutionary history. Range expansions and contractions have also been shown in other Arctic and sub-Arctic species (Andersen et al., 1998;O'Corry-Crowe, 2008;O'Corry-Crowe et al., 1997;O'Corry-Crowe et al., 2010;Palsbøll, Heide-Jorgensen, & Dietz, 1997;Palo, 2003), associated with sequential cooling and warming periods, respectively. Specifically, for beluga whales, late Pleistocene oscillations (7-36 ka) in sea ice extent may explain recent population subdivisions, with genetic support for recurrent episodes of gene flow that likely coincided with warming periods (O'Corry-Crowe et al., 2010). Species divergence, occurring with such episodic gene flow, can have important impacts on population demographic parameters such as connectivity and size (Hey, 2010b;Moore, Gow, Taylor, & Hendry, 2007;Nosil, 2008), which can then lead to substantial alteration of a species' adaptive potential (Eizaguirre & Baltazar-Soares, 2014;Hellberg, 2009;Pelletier, Garant, & Hendry, 2009). Estimating these parameters is, therefore, an important step toward a better understanding of the evolutionary histories and potential trajectories of different species, especially considering the dynamic environmental conditions that likely govern potential gene flow among them (Seehausen et al., 2008). Disentangling the interplay among demographic parameters, environmental fluctuations, and/or geographic distribution of admixtures and Pacific signatures likely indicates speciation associated with processes other than the closing of the Isthmus of Panama. This time span coincides with further planetary cooling in the early Quaternary period followed by the onset of oscillating glacial-interglacial cycles. We propose that the initial S. microcephalus-S. pacificus split, and subsequent hybridization events, were likely associated with the onset of Pleistocene glacial oscillations, whereby fluctuating sea levels constrained connectivity among Arctic oceanic basins, Arctic marginal seas, and the North Atlantic Ocean. Our data demonstrates support for the evolutionary consequences of oscillatory vicariance via transient oceanic isolation with subsequent secondary contact associated with fluctuating sea levels throughout the Quaternary period-which may serve as a model for the origins of Arctic marine fauna on a broad taxonomic scale.

K E Y W O R D S
elasmobranch, Greenland shark, interspecific gene flow, introgressive hybridization, isolation with migration, Somniosus microcephalus the potential for gene flow, however, can be especially complicated among closely related cryptic sister species whose taxonomic differentiation remains nebulous (Griffiths et al., 2010). Such information is nevertheless urgently needed, especially for species whose main distribution includes large Arctic marine ecosystems, which are currently experiencing major climate-related transformations (Wassmann, Duarte, Agusti, & Sejr, 2011).
The number of currently recognized Somniosus species varies according to sources, with species regarded as geographically isolated from one another despite little differences in physical appearance having occurred among sister taxa over their evolutionary histories (Yano, Stevens, & Compagno, 2004). Greenland sharks (Somniosus microcephalus) are known from the Canadian Arctic eastward, including Greenland and Svalbard through to the Laptev Sea (Chernova et al. 2015), south to Iceland and Nova Scotia on the eastern seaboard of Canada (MacNeil et al., 2012) in the western Atlantic and south to the United Kingdom and southernmost Norway (in deepwater areas). Pacific sleeper sharks (Somniosus pacificus), on the other hand, range from the Bering Sea south throughout the deep waters of the Pacific and into the Southern Ocean (Yano et al., 2004). The smaller Somniosus rostratus is reported from the Mediterranean Sea and the Northeastern temperate Atlantic Ocean, while a fourth species, Somniosus antarcticus has been described from the Southern Ocean (Yano et al., 2004). Current species delimitations, notwithstanding accurate species identification and population delineation for the larger somniosids, (S. microcephlaus, S. pacificus, and S. antarcticus) remain challenging given their overall morphological similarities and the fact that even the few "distinguishing" characteristics for each species often considerably overlap (Benz et al., 2007;MacNeil et al., 2012).
Using partial mitochondrial (mtDNA) cytochrome b (cyt b) gene sequences, Murray et al. (2008) demonstrated that Pacific sleeper sharks and Greenland sharks can be distinguished from one another, with the former exhibiting no geographic structure among samples collected across the Pacific Ocean or more southerly latitudes. More recently, genetic characterization of juvenile Greenland sharks from the Canadian Arctic revealed individuals carrying Pacific sleeper shark mtDNA haplotypes (Hussey et al., 2015). This finding, taken with a growing list of cytonuclear discordance and/or genetic admixture among closely related elasmobranchs in the literature, suggests that hybridization among Somniosus species is possible (Hinojosa-Alvarez, Walter, Diaz-Jaimes, Galván-Magaña, & Paig-Tran, 2016;Kashiwagi, Marshall, Bennett, & Ovenden, 2012;Morgan et al., 2012;Walter et al., 2014). Questions remain whether this potential hybridization is recently developed among once more distinguishable species, whether it results from long-term gene flow between sister species, or whether the shared genetic signatures are a product of incomplete lineage sorting (retained ancestral polymorphism).
In this work, we focus on genetic variation in the Greenland shark across its current known distribution. Our objective was to compare patterns of DNA sequence polymorphism and divergence within and among Greenland sharks sampled throughout the Arctic, and sub-Arctic-temperate regions, with those of the closely related Pacific sleeper shark, and other deepwater somniosids in order to develop a more comprehensive picture of the origins of the Greenland shark.

| Tissue samples
Tissue samples from 247 somniosid sharks were collected across a wide distribution including the Canadian Arctic and sub-Arctic, Western Greenland, Svalbard, Iceland, and the North Atlantic ( Fig. 1) that were identified as Greenland shark (i.e., S. microcephalus) during sample collection. Pacific sleeper shark (S. pacificus; n = 6) tissue samples were collected from Alaska. Additional putative S. microcephalus individuals (e.g., identified as S. microcephalus in the field) were collected from the Mid-Atlantic Ridge (MAR-Azores; n = 1) and the Gulf of Mexico (GM; n = 1), and published Somniosus mtDNA sequences (n = 22) from Murray et al. (2008) were retrieved from GenBank (EF090943-EF090963) and included for comparative purposes (n = 277 Somniosus spp. samples total).

| Molecular methods
Genomic DNA was extracted from all samples using either the standard phenol-chloroform method or the Promega Wizard Extraction kit following the manufacturer's instructions (Promega Corp, Madison WI, USA). We partially sequenced the mtDNA cytochrome b (cyt b) gene, the nuclear recombination activation gene (RAG1), and the nuclear internal transcribed spacer 2 (ITS2) with the flanking 5.8 and 28 ribosomal sequences (5.8S and 28S) herein referred to as the ITS2 locus, for each sampled individual. PCR amplification and sequencing reaction conditions for the markers used are available in the Appendix S1.

| Data analysis
All derived mtDNA sequences were aligned with published GenBank sequences in SEQUENCHER 5.0 (GeneCodes, Ann Arbor, MI) and trimmed to 702 bp. Nuclear sequences were also separately aligned in SEQUENCHER and trimmed to 659 bp and 1,085 bp for RAG1 and ITS2, respectively. To estimate genealogical relationships among haplotypes, a 95% statistical parsimony network was constructed in TCS 1.21 (Clement, Posada, & Crandall, 2000) for each locus. Locus-specific nucleotide models were determined in jModelTest 2.1.10 based on AIC/BIC model comparisons (Darriba, Taboada, Doallo, & Posada, 2012). Haplotype diversity (Hd) and nucleotide diversity (π), Fu's F statistic, Tajima's D, genetic distance (D a ), and tests for recent population expansion from mismatch distributions using cyt b sequences were performed in DnaSP 5.10.01 (Librado & Rozas, 2009). Estimates of genetic differentiation, F ST (frequency-based) and Φ ST (Tamura-Nei) among S. microcephalus sampling sites were calculated in Arlequin 3.5 (Excoffier, Laval, & Schneider, 2005) using the mtDNA data. A demographic history for S. microcephalus based on the cyt b dataset (using HKY+I model derived from jModelTest results) was reconstructed using the Bayesian Skyline Plot (BSP) model (Drummond, Rambaut, Shapiro, & Pybus, 2005) as implemented in the BEAST 1.8.0 package (Drummond, Suchard, Xie, & Rambaut, 2012). Four independent runs were performed, each using 100 mol/L generations and a sampling frequency of 5,000. Each run also used a piecewise-constant skyline model, with 20 groups, and an uncorrelated lognormal relaxed-clock mean rate (ucld.mean) was set to 0.0007 = 7 × 10 −10 Ma −1 as per Martin, Naylor, and Palumbi (1992). Replicates were combined using LogCombiner (BEAST package), applying a burn-in of 50% for each replicate. Finally, the BSP was visualized in the program Tracer 1.5 (Rambaut & Drummond, 2007) to detect changes in population size over time.

| Coalescent analysis and divergence estimation
To estimate long-term demographic parameters, namely N e , migration rates (M; i.e., measures of gene flow) and the time since divergence (t 0 ) between the Greenland shark and the Pacific sleeper shark, we used Isolation with Migration (IM) models as implemented in IMa2 (Hey, 2007(Hey, , 2010a. IM models search parameter space for the most likely genetically based genealogy describing the split between two taxonomic groups, using a Bayesian framework. Estimated genealogies are generated assuming random mating within groups and that neither group exchanges genetic materials with other nonsampled groups (Hey, 2010a;Hey & Nielsen, 2007). Generated genealogies specify branch lengths and group-specific genetic variation (estimated as θ) from which demographic parameters can be computed. IM models were run using the mtDNA cyt b data, with the HKY model of sequence evolution. An initial hypothesized phylogeny between the two species was provided, and uniform priors were set regarding the distribution of each parameter of interest using preliminary runs. Three different substitution rates, each characterizing divergence among elasmobranchs: a conservative estimate (7 × 10 −10 Ma −1 ; Martin et al., 1992), Squalus (9.90 × 10 −10 Ma −1 ; Winchell, Martin, & Mallatt, 2004), and Squatina (1.67 × 10 −9 Ma −1 ; Winchell et al., 2004) at cyt b were used to estimate cyt b substitution rates in this work. These substitution rates respectively range from very conservative (Martin & Palumbi, 1993;Martin et al., 1992) to those possibly more reflective of somniosids (e.g., rates derived for Squalus and Squatina), respectively. Demographic parameter searches were performed with 20 Markov chain Monte Carlo (MCMC) simulations using a geometric heating scheme to ensure adequate mixing among them. Initial searches were set to run between 1 x 10 6 and 12 x 10 6 iterations in order for the genealogies to reach a stationary distribution from which to sample. Once stationarity was reached, genealogies were permutated 10 x 10 6 times and sampled every 100th step for a total of 100,000 recorded genealogies from which demographic parameters were calculated (see IMa2 manual; Hey, 2010a). The different searches were repeated three times, and all derived parameters were compared among all three iterations to evaluate consistency (Hey, 2010a). All genealogies recovered from each search were used to compare the information content and thus likelihood (L-mode) of different IM models. To check for consistency in estimated parameters, model comparison with searches omitting some parameters was also run by specifically setting migration rate priors to zero, hereafter referred to as Isolation Only (ISO). ISO models prevent any migration rate influences on estimated demographic parameters and generate parameter estimates consistent with an incomplete lineage sorting (ILS) scenario.
Due to low levels of polymorphism at nuclear loci, IM analyses were not performed for the RAG1 and ITS2 datasets. Nuclear nucleotide substitution rates for RAG1 (4.85 × 10 −10 Ma −1 ; Eytan & Hellberg, 2010) and ITS2 (1.50 × 10 −9 Ma −1 ; Montoya-Burgos, 2003) were also used to estimate RAG1 and ITS2 divergence times using a simple formula whereby the number of observed mutations at a locus is divided by the estimated gene-wide percent divergence rate per Ma.

| RESULTS
A total of 61 mtDNA haplotypes were recovered across all sleeper shark samples forming two distinct and well-differentiated mtDNA haplogroups (Fig. 2, Tables S1 and S2). The Greenland shark samples consisted of 45 haplotypes, including 30 that are novel to this study (green circles with black text Fig. 2; see Table S2). These samples showed a haplotype diversity (Hd) of 0.7100 and nucleotide diversity (π) of 0.0028, while those recovered for the Pacific sleeper sharks were 0.9700 and 0.0053, for Hd and π, respectively (Table S2). A total of 18 haplotypes were among the Pacific sleeper shark samples, including four novel to this study. Site-by-site statistics are reported in Table 1. DNA sequences are deposited in GenBank under the accession numbers MF537353-MF537393 and MF555591-MF555594.
Only two alleles were recovered from each of the nuclear markers, with these alleles characterized by mostly species-specific single nucleotide polymorphisms (SNPs). For the RAG1 marker, a "C" was F I G U R E 2 Haplotype variation at cyt b, ITS2, and RAG1 loci among sampled somniosids as depicted through a 95% statistically parsimony network. For cyt b, the size of the circles is proportional to the frequency of each haplotype in the 277 individual dataset, ITS2 and RAG1 haplotypes are not scaled for frequency. Haplotype numbers correspond to Table S1, white numbers are previously reported haplotypes, black numbers are novel haplotypes to this study. Individuals identified in the field as Somniosus microcephalus possess haplotypes in green, while those identified as Somniosus pacificus possess blue haplotypes for the Pacific ITS2 allele, but also displayed heterozygosity for the Greenland shark and Pacific sleeper shark alleles at RAG1. Thus, although genetic sequences at the sampled loci tended to separate the somniosid species, evidence of genetic admixture between the two species was clear.  (Fig. 3a), while the BSP revealed a long period of population stability followed by a recent demographic expansion beginning around 0.5 Ma ago (Fig. 3b).

| Demographic expansion
There was no significant genetic differentiation of S. microcephalus among geographic sampling sites with the exception of Cumberland

| Coalescent IM versus ILS analysis and divergence estimates
Across all isolation with migration (IM) models assessing coalescentbased demographic parameter, strongest support was recovered from models with unequal migration rates between the species. The best models were characterized by coalescent migration rates equal to zero from Pacific sleeper sharks into Greenland sharks, while Greenland sharks into Pacific sleeper sharks migration rates were greater than zero. The best model also included varying long-term effective population sizes for extant and ancestral populations, following AIC model selection tests. Following conversion of coalescent rates to "forwardtime," results indicate asymmetrical migration whereby gene flow likely only occurs from Pacific sleeper sharks into Greenland sharks since the species diverged from one another.

| mtDNA (cyt b) analysis
Simulations of IM models (m prior upper bound = 0.2) and ISO (m = 0) models produced similar maximum-likelihood estimates (hereafter MLEs) of species-split time (t 0 ) of 1.190 (Fig. 4a) and 1.130, respectively (see Table 3 for comparative split times in years). For the IM analyses, a second older and less probable peak was identified for t 0 at 7.550 followed by a flat nonzero right-end tail (Fig. 4a). Although overall migration rate parameters (forward-time) were not significant following the Nielsen and Wakeley (2001) log-likelihood ratio test, the population-specific migration rate (2Nm; also forward-time) for Pacific sleeper into Greenland shark was significantly different from zero (p < .05). Marginal MLE of migration were asymmetrical with m PS>GS 0.01330 and m GS>PS 0.00010 (Fig. 4b). Marginal unscaled estimates of q parameters and scaled N e are presented in Table 3. Although generation times have not been reported for either S. microcephalus or S. pacificus, we applied a relatively short estimate based on a 20year age at first reproduction hypothesized for another squaliform species Centroselachus crepidater (Irvine, Stevens, & Laurenson, 2006;Veríssimo, McDowell, & Graves, 2011), and an approximated 150- year generation length following the ~156 age at first reproduction hypothesized by Nielsen et al. (2016). The longer generation time scaled N e estimates were approximately 13% smaller than the shortgeneration time estimates (Table 3).

| Nuclear analysis
Nuclear divergence estimates between species were comparable to the cyt b coalescent analysis. An estimated split time of 3.16 Ma was found based on the single mutation difference between Greenland sharks and Pacific sleeper shark RAG1 alleles; while the ITS2 locus estimate was slightly younger at 2.52 Ma, based on the four mutation steps between both shark species.

| DISCUSSION
Our genetic analyses clearly demonstrate the existence of two ge-  The mixed genetic signatures among individuals identified in the field as Greenland shark could result from two possible processes: (1) ongoing interspecific hybridization (introgression upon secondary contact), or (2) the retention of ancestral DNA polymorphisms at the observed loci (incomplete lineage sorting-ILS). To attempt to distinguish between these two processes we applied IMa2 analyses on the mtDNA data. All IMa2 model-fitting analysis indicated higher support for models including gene flow as the initial species split, and in all cases, this gene flow was asymmetrical with alleles originating in Pacific sleeper shark and migrating into Greenland shark.
These Estimates of species-split (divergence) time for Greenland shark and Pacific sleeper shark using mtDNA did not show a consistent time line for IM and ILS models across a range of substitution rates. The observed bimodal t 0 curve, however, is consistent with a "two attractors" hypothesis described as "a recent signal of no migration, and an older signal with higher migration" (J. Hey as cited in Choleva et al., 2014). Application of a broadly based conservative substitution rate (Martin et al., 1992) and more probable substitution rates (see Table 3) identified that an IM-based estimated species-split time using the cyt

| Ice-olation with migration hypothesis
We introduce the hypothesis of oscillatory vicariance with secondary contact-herein referred to as Ice-olation with Migration. If we consider that the ancestral sleeper shark was a deepwater species that was constrained to this habitat via temperature regimes, this species likely occupied a pan-oceanic distribution prior to the late Miocene, but only in the deep ocean. This is also consistent with the occurrence of somniosids in tropical and temperate regions where they are typically observed at depths greater than 1,000 m and reportedly down to 2,200 m (Benz et al., 2007;Campana, Fisk, & Klimley, 2015;Herdendorf & Berra, 1995), and in contrast to somniosids sampled in the Arctic and northern Quebec, where satellite telemetry data have shown sharks occurring from the surface down to approximately 1,500 m (Campana et al., 2015;Stokesbury, Harvey-Clark, Gallant, Block, & Myers, 2005).
From the late Miocene, and into the Quaternary, a drastic reduction in global temperatures occurred. This decline in sea temperature resulted in sea ice formation in polar regions, including the formation of ice shelves up to 1 km in thickness (Fastook & Hughes, 2013;Grosswald & Hughes, 1999;Polyak, Edwards, Coakley, & Jakobsson, 2001). This ice thickness coupled with the submerged mountainous terrain of the Arctic, may have effectively constrained population connectivity of ancestral somniosids between the Arctic sub-basins/ marginal seas and the remaining pan-oceanic distribution. In support of the Arctic sub-basin isolation, glacial scouring patterns along the Lomonosov Ridge are consistent with ice shelf formations (Polyak et al., 2001;Stein et al., 2016), and it is likely that similar constraints to connectivity occurred among the Arctic and Atlantic Ocean at the East Greenland Rift Basin (Vigdorchik, 1979;Vigdorchik & Vinkovetsky, 1976

| Population expansion and adaptation in the Arctic
Ecological specialization has been proposed as a driver for marine speciation (Corrigan & Beheregaray, 2009;Kashiwagi et al., 2012;Morin et al., 2010;Wang, Chou, & White,1999) and is a plausible explanation for the demographic expansion of sleeper sharks in the Arctic region. Greenland sharks are considered generalists and opportunistic predators feeding on a broad range of benthic and pelagic teleost and cephalopod prey (MacNeil et al., 2012). Marine mammals including ringed-seal (Pusa hispida) have been recovered from the stomach contents of Arctic Greenland sharks (Fisk, Tittlemier, Pranschke, & Norstrom, 2002;Leclerc et al., 2012;Lydersen, Fisk, & Kovacs, 2016;McMeans, Svavarsson, Dennard, & Fisk, 2010;Nielsen, Hedeholm, Simon, & Steffensen, 2014;Ridoux, Hall, Steingrimsson, & Olafsson, 1998), as well as in the stomachs of Pacific sleeper sharks at high latitudes (Bright, 1959;Horning & Mellish, 2014;Hulbert, Sigler, & Lunsford, 2001), identifying both species as high trophic level feeders T A B L E 3 (a) Time-splitting estimates (t 0 ) obtained from IMa2 for Isolation with Migration (IM) and Isolation only (ISO) models unscaled and scaled in millions years among Somniosus microcephalus and Somniosus pacificus, using published mtDNA substitution rates: mtDNA Conservative (Martin et al., 1992), mtDNA Squalus (Winchell et al., 2004), and mtDNA Squatina (Winchell et al., 2004  . When considering the global distribution of shark and pinniped diversity, an inverse relationship occurs, with the highest shark diversity at tropical latitudes where pinniped diversity is low, and few large shark species, including the Greenland and Pacific sleeper shark, occurring at high latitudes where pinniped diversity is highest (Ferguson, Higdon, Tallman, Fisk, & Hussey, 2014). The diversity and abundance of marine mammals in the Arctic might have served as novel food sources for sleeper sharks during their colonization of the Arctic environment.

| CONCLUSION
Our results show that the nearly morphologically identical sister species Greenland shark and Pacific sleeper shark can be genetically distinguished using a combination of mtDNA and nuclear markers, but interspecific hybridization (introgression) between these two species is probable and has likely been ongoing since initial divergence. Our results confirm the fourth reported instance of interspecific hybridization among elasmobranchs (Hussey et al., 2015; see also Morgan et al., 2012;Walter et al., 2014;Cruz et al., 2015), which is likely more common among this group of fishes than previously thought. Our findings, considering a range of substitution rates for the genetic markers employed here, are consistent with a more recent Greenland shark and Pacific sleeper shark speciation event occurring between 1 and 2.34 Ma, more recently than the closure of the Isthmus of Panama.
Our estimates of long-term N e s also indicate slightly higher sizes for Greenland sharks compared to Pacific sleeper sharks. However, recent radiocarbon analysis of S. microcephalus eye lenses has concluded that Greenland sharks may not reach maturity until 150 years of age, and possibly live as long as 400 years (Nielsen et al., 2016). If accurate, this would mean that effective population sizes are but a fraction of our 20-year estimates, and using the N e estimates from the 150-year generation time of Nielsen et al. (2016) is likely more appropriate (Table 3). The lack of geographic structuring for Greenland sharks across the Canadian Arctic and beyond, but with the recurring hybridization with Pacific sleeper sharks at the edge of its distribution in warmer waters, suggests that a continuing warming Arctic may further reduce the genetic integrity of both species. Future research focused on broad-scale population-level connectivity across the range of S. microcephalus is required to clarify the consequences of past climatic and present-day environmental changes on the hybridization rates and evolutionary trajectory of these species.

ACKNOWLEDGMENTS
We thank the crew of the Nunavut Government Research vessel,