The origin and remolding of genomic islands of differentiation in the European sea bass

Abstract

Speciation is a complex process that leads to the progressive establishment of reproductive isolation barriers between diverging populations. Genome-wide comparisons between closely related species have revealed the existence of heterogeneous divergence patterns, dominated by genomic islands of increased divergence supposed to contain reproductive isolation loci. However, this divergence landscape only provides a static picture of the dynamic process of speciation, during which confounding mechanisms unrelated to speciation can interfere. Here we use haplotype-resolved whole-genome sequences to identify the mechanisms responsible for the formation of genomic islands between Atlantic and Mediterranean sea bass lineages. Local ancestry patterns show that genomic islands first emerged in allopatry through linked selection acting on a heterogeneous recombination landscape. Then, upon secondary contact, preexisting islands were strongly remolded by differential introgression, revealing variable fitness effects among regions involved in reproductive isolation. Interestingly, we find that divergent regions containing ancient polymorphisms conferred the strongest resistance to introgression.

Fig. 1

Spatial population structure and local ancestry patterns. a Geographical location of samples, including the newly sequenced genomes (colored circles) and additional reference samples from the Atlantic (dark gray), W-MED (gray), E-MED (light gray), and the Algerian admixture zone (gray crosses). b Principal Component Analysis of newly sequenced genomes combined with 112 individuals genotyped at 13,094 common SNPs (MAF > 0.1). The first PCA axis distinguishes Atlantic and Mediterranean populations while the second axis reveals a subtle population structure between W-MED and E-MED. Some individuals from the Algerian coast represent admixed genotypes between Atlantic and Mediterranean populations. c Schematic representation of a 2 Mb region within chromosome 7, showing the mosaic of ancestry blocks derived from Atlantic (red) and Mediterranean (yellow) populations. For simplification, we only display two individual haplotypes from Atlantic samples, four from W-MED, four from E-MED, and two reconstructed ancestral Mediterranean haplotypes from E-MED samples. d Migrant tract length distribution obtained for the Atlantic (red, showing tracts of Mediterranean origin), W-MED (orange), and E-MED (yellow) populations (showing tracts of Atlantic origin) using four individuals per population and a total of 2,628,725 phased SNPs

Fig. 2

Observed and simulated migrant tract length distributions in the Atlantic and Mediterranean populations. Observed distributions (gray) of migrant tract length are compared with simulated distributions (colored) under the post-glacial secondary contact scenario illustrated in the top-right corner. The abundance of introgressed tracts as a function of their length is represented for observed vs. simulated data in the Atlantic (dark gray vs. red circles, showing tracts of Mediterranean origin) and Mediterranean populations (light gray vs. yellow circles, showing tracts of Atlantic origin)

Fig. 3

Demographic history inferred from the length distribution of IBS tracts. a Flexible demographic model accounting for multiple equal-length episodes of divergence and gene flow between Atlantic and Mediterranean sea bass populations. An ancestral population of size N splits into two populations of size N₁ and N₂, experiencing one to several (m) cycles of interrupted gene flow during T_diff generations followed by migration during T_c generations. Each contact episode contains one to several (n) pulses of admixture (black arrows), replacing Mediterranean and Atlantic populations by a proportion f₁ and f₂ of migrants, respectively. The most recent admixture pulse occurs at time T_c/100 before the end of each contact episode, and preceding pulses are homogeneously distributed every (T_c−T_c/100)/n time interval. b Modeling scenarios fall into three different categories according to the distribution of admixture pulses over time: continuous migration, secondary contact, and periodic pulses. The illustrated example shows these three categories with nine pulses represented by vertical bars along the timeline. c The log-likelihood values obtained for each of the three modeling scenarios from two to ten admixture pulses. Two possible configurations of the periodic pulses scenario exist for models including a total of six, eight, and ten pulses (Supplementary Table 2), only the best of which is represented here. The periodic pulses scenario is not defined for two, three, five, and seven pulses. d Goodness-of-fit of the best model (m = 1, n = 10), showing the length distributions of IBS tracts from observed data (colored lines) compared to model prediction (gray lines). Upper-left: the three distributions observed within ATL and W-MED and between ATL and W-MED populations. Upper-right: model fit within ATL. Lower-left: model fit within W-MED. Lower-right: model fit between ATL and W-MED

Fig. 4

Population genetic statistics calculated in non-overlapping 100 kb windows along chromosomes 7 and 17. a F_ST measured between the Atlantic and contemporary (purple) or ancestral reconstructed (mauve) W-MED or E-MED (b) population. c d_XY calculated between the Atlantic and the W-MED (black) or E-MED (blue) population. d Fraction of introgressed tracts in the W-MED (orange) or E-MED (yellow) population. e RND_min measured between the Atlantic and W-MED (black) or E-MED (red) population. f Average length of introgressed tracts in the W-MED (dark blue) or E-MED (light blue) population. g Population-scaled recombination rate (ρ = 4N_er) averaged between Atlantic and Mediterranean population

Fig. 5

Relationships between divergence (F_ST and d_XY), the population-scaled recombination rate (ρ = 4N_er), and introgression statistics (RND_min and F_intro). a–d The density of points appears in color scale from low (blue) to high (red) densities. e–h The color scale indicates the value of RND_min (e, f) or the frequency of introgression (g, h) in the corresponding window from low (blue) to high (yellow) introgression rate values

Fig. 6

Simulations under different modes of selection to understand the mechanisms underlying genomic islands. a Comparison among simulated distributions of 100 kb window-averaged F_ST (top), d_XY (middle), and RND_min (bottom) obtained under different versions of the secondary contact scenario, including neutral divergence and introgression (Neutral), varying strengths of background selection (BGS, from 0.5 to 0.1 × N_e), and BGS with reproductive isolation (BGS+RI). b Comparison among observed distributions of F_ST (top), d_XY (middle), and RND_min (bottom) for real genomic windows of 100 kb that were either assigned to the Neutral model (7% of windows), the BGS model (76% of windows), or the BGS+RI model (17% of windows) based on Euclidean distances to simulated data