Ancestral gene flow at the base of Heliconius.

(A) Posterior probabilities of species trees for major Heliconius clades inferred from BPP analysis of 200-locus blocks (each spanning 1–3 Mb) across the genome under the multispecies coalescent (MSC) model with no gene flow. The y-axis shows the posterior probability of species trees and ranges from 0 to 1. Colors correspond to the nine most common MAP trees, summarized by lumping species in the melpomene-silvaniform clade (Bes, Num, Mel) into a single tip (BNM); see Figure S1 for full trees. Proportions of coding and noncoding blocks with each tree as a MAP tree are shown in parentheses. (B, C) Two scenarios of early divergence of Heliconius: (B) erato-early versus (C) aoede-early. Each tree is a MAP tree from a block having the MAP tree supporting one of these scenarios, with estimates of branch lengths (posterior means). (D, E) Phylogenetic and introgression histories estimated under an MSC model with introgression (MSC-I) corresponding to the two scenarios based on coding loci in chromosome 1 (3,517 coding loci; see Figures S7 and S8). (F) Proportions of trees of scenarios 1 or 2 in each chromosome in order of increasing number of loci (used as a proxy for chromosome length). The Z chromosome (chr. 21) is placed at the right end. Number of blocks is shown on top of each bar. Tal: Eueides tales, Mel: H. melpomene, Bes: H. besckei, Num: H. numata, Bur: H. burneyi, Dor: H. doris, Aoe: H. aoede, Era: H. erato, Sar: H. sara.

Major introgression events in the melpomene-silvaniform clade.

(A) Blockwise estimates of species trees of the melpomene-silvaniform clade inferred from 200-locus blocks across the genome under the MSC model with no gene flow using BPP (see Table S10 for data, Tables S11 and S12 and Figure S10 for full results). MAP trees are labelled in decreasing order of frequency among blocks. Proportions of coding and noncoding blocks with each tree as a MAP tree are shown in parentheses. (B) The MSC-I model with introgression events that can explain the three major groups of trees in (A). Branch lengths are based on posterior means of divergence/introgression times estimated from 5,341 noncoding loci on chromosome 1 (Table S15). Each internal node is given a label, which is used to refer to a population above the node, e.g. the population between nodes r and bn is referred to as branch bn. Each horizontal arrow represents a unidirectional introgression event, e.g. the arrow from tcmeph1 to n3 represents tcmeph→Num introgression at time τtcmeph1 = τn3 with probability φtcmeph1→n3. (C) Continuous migration (IM) model for the pardalinus-hecale clade, allowing bidirectional gene flow among the three species. (D) MSC-M model for the cydno-melpomene clade. For C and D, branch lengths are based on estimates from noncoding loci in chromosome 1 (Figure S15; Tables S17–19), and arrow sizes are proportional to estimated migration rate (M = Nm). Bes: H. besckei, Num: H. numata, Par: H. pardalinus, Ele: H. elevatus, Hec: H. hecale, Tim: H. timareta, Cyd: H. cydno, Mel: H. melpomene.

Introgression history of the chromosome 15 inversion region (15b)

(A) Blockwise estimates of species trees for the inversion (15b) and the remnant flanking regions (15a and 15c) of chromosome 15. Trees were inferred from 200-locus blocks and 100-locus blocks under the MSC model without gene flow using BPP (Table S20). Tree legends are grouped by whether H. numata clusters with H. ismenius (blue), or H. numata with P1 inversion (Num11) clusters with H. pardalinus (red), or other relationships (green). (B–E) Four possible scenarios of the origin and introgression route of the P1 inversion. Star indicates the origin of the P1 inversion. Green lineages have the inversion, and green arrows indicates introgression of the inversion. Ism: H. ismenius, Num00: H. numata homozygous uninverted 15b (H. n. laura and H. n. silvana), Num11: H. numata homozygous for inversion P1 (H. n. bicoloratus), Eth: H. ethilla. See Figure 2 legend for codes for other species. For the direction of melpomene-cydno clade→H. elevatus, see Figure S18.

Revised phylogeny and introgression history of Heliconius.

Phylogeny of the erato-sara clade was previously estimated using a similar approach (Thawornwattana et al. 2022). Arrows represent introgression events. Introgression of chromosomal inversions are region-specific, indicated by a label (13b or 15b) above the arrow. Status of inversions 13b and 15b in each species is indicated in the two rows at the bottom, with + indicating inversion, – indicating standard orientation, and +/– inversion polymorphism. Grey shading represents a period of continuous gene flow. Triangle represents the origin of pollen feeding near the base of Heliconius. Star indicates a possible origin of the P1 inversion on chromosome 15, and green branches indicate lineages with the inversion (polymorphic in Num, fixed in Par). Him: H. himera, Sia: H. hecalesia, Tel: H. telesiphe, Dem: H. demeter; see Figures 1– 3 legends for other species codes.