A genomic test of sex-biased dispersal in white sharks

Abstract

Mitonuclear discordance has been observed in several shark species. Female philopatry has often been invoked to explain such discordance but has never been explicitly tested. Here, we focus on the white shark, for which female philopatry has been previously proposed, and produced a chromosome-level genome, high-coverage whole-genome autosomal, and uniparental datasets to investigate mitonuclear discordance. We first reconstructed the historical population demography of the species based on autosomal data. We show that this species once comprised a single panmictic population, which experienced a steady decline until recent times when it fragmented into at least three main autosomal genetic groups. Mitochondrial data depict a strikingly different picture, inconsistent with the spatial distribution of autosomal diversity. Using the demographic scenario established from autosomal data, we performed coalescent and forward simulations to test for the occurrence of female philopatry. Coalescent simulations showed that the model can reproduce the autosomal variability, confirming its robustness. A forward simulation framework was further built to explicitly account for a sex-biased reproduction model and track both autosomal and uniparental markers (Y chromosome and mitochondrial DNA). While our model generates data that are consistent with the observed Y chromosome variation, the mitochondrial pattern is never reproduced even under extreme female philopatry (no female migration), strongly suggesting that demography alone cannot explain the mitonuclear discordance. Our framework could, and perhaps should, be extended to other shark species where philopatry has been suggested. It is possible that the proposed widespread occurrence of female philopatry in sharks should be revisited.

Keywords: female philopatry; historical demography; mitonuclear discordance; population genomics; white shark.

Fig. 1.

Worldwide population structure of the white shark. (A) Location of sampling sites. The observed (red shading) and possible (orange shading) range of the white shark extracted from the IUCN Red List assessment (27). Circles represent the sampling locations of this study with the number of individuals for each molecular marker indicated in parentheses, in this order: Whole genome resequencing, target gene capture, mtDNA, and Y chromosome. NAO, North Atlantic Ocean; MED, Mediterranean Sea; WIO, Western Indian Ocean; SPO, Southern Pacific Ocean; NPO, Northern Pacific Ocean; JAP, Japan. (B) Principal component analysis of the 14 individuals of the WGS_aut dataset. The proportion of variance explained is indicated in parentheses. (C) Ancestry proportions estimated with sNMF for K = 3. Each bar corresponds to an individual. (D) Individual level of runs of homozygosity (ROH). Each dot represents an individual. The red line depicts the linear relationship (95% CI in gray) between the total length (sum) and the number of ROH. Dot colors represent the sampling location of each individual.

Fig. 2.

Historical demography through time. (A) SMC++ estimates of effective population size (Ne) variation through time. (B) Variation of Ne as reconstructed by the PSMC algorithm for two high-coverage individuals (GN17768, NAO and GN17525, SPO). Blurred lines correspond to the bootstrap interval. (C) Connectivity graph estimated by SNIF for c = 5 for the same individuals as in (B) averaged over 10 replicates for each of 10 values of distance parameter ω. Line colors correspond to sampling sites. Vertical dotted lines in (B) and (C) correspond to the times t of changes in gene flow set for the SNIF analysis.

Fig. 3.

Demographic modeling of white shark evolutionary history. Maximum-likelihood parameters’ value inferred with FastSimCoal2 under the REC-DIV model. CI computed under a parametric bootstrap approach are provided (SI Appendix, Table S5). Rectangle widths scale with the estimated effective population size. Arrows indicate the number of migrants (Nm) exchanged between populations, with thickness proportional to the inferred values. Dotted lines indicate time events.

Fig. 4.

Mitochondrial and Y chromosome structure. (A) Haplotype network for mitochondrial DNA. (B) Haplotype network for Y chromosome. Mutations are represented by bars perpendicular to the branches of the network. Each haplotype is represented by a circle whose size is proportional to the frequency observed in our sampling. Colors correspond to sampling sites.

Fig. 5.

Bayesian phylogenetic tree and molecular dating of mitochondrial sequences. (A) Bayesian phylogenetic tree built using Carcharodon carcharias mtDNA haplotypes and four outgroups (grouped under Isurus/Lamna, full version of the tree available in SI Appendix, Fig. S11). Node ages are expressed in millions of years, with the 95% credible interval (CI) in brackets. Colors representing C. carcharias lineages correspond to the six sampling sites. The asterisk next to WIO indicates the presence of two SPO individuals within this cluster, as shown in Fig. 4. Branch lengths have been rescaled for better visualization. (B) Time to the most recent common ancestor (T_MRCA) of C. carcharias mtDNA haplotypes computed with BEAST (gray) and expected distribution of the T_MRCA after 10,000 coalescent simulations from the REC-DIV model adjusted for a haploid marker and assuming a sex ratio of 1:1 (black).

Fig. 6.

Pairwise-F_ST and ϕst distribution of autosomal and mitochondrial markers respectively under the three scenarios simulated with SliM. (A) Pairwise-F_ST computed on simulated autosomes using the Reynolds’ estimator, on a dataset reproducing the real sample sizes of the WGS_aut dataset, (B) Pairwise-ϕst computed on simulated mtDNA markers, on a sample reproducing the real mtDNA dataset for the four populations. Colors of the boxplot correspond to the simulated scenarios. Each scenario was replicated 100 times. Red dashed lines represent the observed F_ST and ϕst for the WGS_aut or mtDNA datasets respectively.