Whole-genome sequencing of the blue whale and other rorquals finds signatures for introgressive gene flow

Abstract

Reconstructing the evolution of baleen whales (Mysticeti) has been problematic because morphological and genetic analyses have produced different scenarios. This might be caused by genomic admixture that may have taken place among some rorquals. We present the genomes of six whales, including the blue whale (Balaenoptera musculus), to reconstruct a species tree of baleen whales and to identify phylogenetic conflicts. Evolutionary multilocus analyses of 34,192 genome fragments reveal a fast radiation of rorquals at 10.5 to 7.5 million years ago coinciding with oceanic circulation shifts. The evolutionarily enigmatic gray whale (Eschrichtius robustus) is placed among rorquals, and the blue whale genome shows a high degree of heterozygosity. The nearly equal frequency of conflicting gene trees suggests that speciation of rorqual evolution occurred under gene flow, which is best depicted by evolutionary networks. Especially in marine environments, sympatric speciation might be common; our results raise questions about how genetic divergence can be established.

Fig. 1. MSC tree.

(A) An MSC species tree was constructed from 34,192 individual GFs. Internal branches within Balaenopteridae are numbered 1 to 7. All branches receive maximal support (P = 1.0, ASTRAL analysis). Branch lengths were calculated from an ML analysis. Gray whales, family Eschrichtiidae, are placed inside Balaenopteridae as a sister group to fin and humpback whales. (B) ASTRAL quartet-score analyses for branches 1 to 7 (A). Quartet scores were calculated for the three possible arrangements (q1 to q3) for the respective branch. The principal quartet trees are depicted, with q1 representing the species tree. Branch nos. 2 and 3 receive only limited quartet scores, and no quartet can be significantly rejected.

Fig. 2. Median network of 34,192 GF ML trees with 11% threshold.

Conflicting evolutionary signals characterize the center of the network, which is equivalent to branch no. 3 in the species tree (Fig. 1). In addition, placing the minke whale has some conflicting signal, but the elongated rectangle indicates a higher degree of resolution. The number of supporting GFs is shown for selected splits. Colored circles indicate taxonomic classification. Blue, Balaenoptera; red, Megaptera; yellow, Eschrichtius; green, Balaena and Eubalaena.

Fig. 3. Gene flow signals for baleen whales inferred by the D statistic, D_FOIL, and PhyloNet.

(A) The species tree of baleen whales with gene flow signals detected by the D statistic and D_FOIL indicated by dashed lines. Signals I to IV were inferred by the D statistic, and signals V, VI, and VII were detected by D_FOIL and were partially corroborated by the D statistic. Note that D_FOIL cannot infer gene flow involving the minke whale. (B to D) Rooted networks for the Balaenopteridae sensu lato phylogeny with reticulations inferred from PhyloNet based on 34,192 20-kbp GFs. Reticulations are shown as blue arrows with inheritance probability denoted above or below. Log-likelihood scores are shown below the networks. Notably, inheritance probability around 33% resembles the distribution of quartet scores and the phylogenetic signals from GFs (Fig. 1). (B) The three best networks indicated a reticulation originating at the circled three branches to minke whale. Similar likelihood scores do not allow the identification of a single origin of gene flow; therefore, the networks were merged, and a range of inheritance probabilities is given. (C) The fourth best network has only a marginally poorer likelihood score and indicates a reticulation between the ancestor of the fin and humpback whale and that of the minke whale. (D) The fifth best network has the same likelihood as (C) and finds an alternative placement of gray whale (blue branch) and reticulation from the ancestor of the blue and sei whale to that of the minke whale.

Fig. 4. Demographic history and genome-wide heterozygosity.

(A) Genome-wide heterozygosity estimated from genomic 100-kbp windows. (B) Historical N_e using the PSMC analyses for all baleen whale genomes. The x axis shows the time, and the y axis shows N_e. Plots were scaled using a mutation rate (μ) of 1.39 × 10⁻⁸ substitutions nucleotide⁻¹ generation⁻¹ and species-specific generation times (g). Generation times are noted next to the species names. Light brown shading indicates interglacials (IG) in the Pleistocene and Holocene, and gray shading indicates the MPT and the PPT.

Fig. 5. Divergence time tree of Cetancodonta (56) including the newly sequenced baleen whales, estimated from 234,947 amino acid sites (2778 orthologs).

Rorquals diverged in the late Miocene, 10.5 to 7.5 Ma ago. Four other cetartiodactyl species were also included but not shown due to space constraints; the dog (Canis lupus familiaris) was used as an outgroup. Five calibration points were used for dating (table S8) (, –60).