Hybridization drives the evolutionary trajectory of many species or local populations, and assessing the geographic extent and genetic impact of interspecific gene flow may provide invaluable clues to understand population divergence or the adaptive relevance of admixture. In North America, hares (Lepus spp.) are key species for ecosystem dynamics and their evolutionary history may have been affected by hybridization. Here we reconstructed the speciation history of the three most widespread hares in North America - the snowshoe hare (Lepus americanus), the white-tailed jackrabbit (L. townsendii) and the black-tailed jackrabbit (L. californicus) - by analyzing sequence variation at eight nuclear markers and one mitochondrial DNA (mtDNA) locus (6 240 bp; 94 specimens). A multilocus-multispecies coalescent-based phylogeny suggests that L. americanus diverged ~2.7 Mya and that L. californicus and L. townsendii split more recently (~1.2 Mya). Within L. americanus a deep history of cryptic divergence (~2.0 Mya) was inferred, which coincides with major speciation events in other North American species. While the isolation-with-migration model suggested that nuclear gene flow was generally rare or absent among species or major genetic groups, coalescent simulations of mtDNA divergence revealed historical mtDNA introgression from L. californicus into the Pacific Northwest populations of L. americanus. This finding marks a history of past reticulation between these species, which may have affected other parts of the genome and influence the adaptive potential of hares during climate change.
SPTBN1_pre-phase
SPTBN1 alignment before phase detemrination (see GenBank accession numbers in Sup. Table 1).
1_SPTBN1_raw.fas
PRKCI_pre-phase
PRKCI alignment before phasing (see GenBank accession numbers in Sup. Table 1).
2_PRKCI_raw.fas
DARC_pre-phase
DARC alignment before phase detemrination (see GenBank accession numbers in Sup. Table 1).
3_DARC_raw_final_final.fas
KITLG_pre-phase
KITLG alignment before phase detemrination (see GenBank accession numbers in Sup. Table 1).
4_KITLG_raw.fas
TF_pre-phase
TF alignment before phase detemrination (see GenBank accession numbers in Sup. Table 1).
5_TF_raw.fas
POLA1_pre-phase
POLA1 alignment before phase detemrination (see GenBank accession numbers in Sup. Table 1).
6_POLA1_raw_final.fas
GRIA3_pre-phase
GRIA3 alignment before phase detemrination (see GenBank accession numbers in Sup. Table 1).
7_GRIA3_raw.fas
SRY
SRY alignment (see GenBank accession numbers in Sup. Table 1).
8_SRY_raw.fas
CYTB
CYTB alignment (see GenBank accession numbers in Sup. Table 1).
9_CYTB_raw.fas
SPTBN1_phased
SPTBN1 phased alignment (see GenBank accession numbers in Sup. Table 1).
1_SPTBN1_phased.fas
PRKCI_phased
PRKCI phased alignment (see GenBank accession numbers in Sup. Table 1).
2_PRKCI_phased.fas
DARC_phased
DARC phased alignment (see GenBank accession numbers in Sup. Table 1).
3_DARC_phased_final.fas
KITLG_phased
KITLG phased alignment (see GenBank accession numbers in Sup. Table 1).
4_KITLG_phased.fas
TF_phased
TF phased alignment (see GenBank accession numbers in Sup. Table 1).
5_TF_phased.fas
POLA1_phased
POLA1 phased alignment (see GenBank accession numbers in Sup. Table 1).
6_POLA1_phased_final.fas
GRIA3_phased
GRIA3 phased alignment (see GenBank accession numbers in Sup. Table 1).
7_GRIA3_phased_final.fas
SPTBN1_lnrb
SPTBN1 phased alignment of longest non-recombining blocks (see GenBank accession numbers in Sup. Table 1).
1_SPTBN1_lnrb.fas
PRKCI_lnrb
PRKCI phased alignment of longest non-recombining blocks (see GenBank accession numbers in Sup. Table 1).
2_PRKCI_lnrb.fas
DARC_lnrb
DARC phased alignment of longest non-recombining blocks (see GenBank accession numbers in Sup. Table 1).
3_DARC_lnrb_final.fas
KITLG_lnrb
KITLG phased alignment of longest non-recombining blocks (see GenBank accession numbers in Sup. Table 1).
4_KITLG_lnrb.fas
TF_lnrb
TF phased alignment of longest non-recombining blocks (see GenBank accession numbers in Sup. Table 1).
5_TF_lnrb.fas
POLA1_lnrb
POLA1 phased alignment of longest non-recombining blocks (see GenBank accession numbers in Sup. Table 1).
6_POLA1_lnrb.fas
GRIA3_lnrb
GRIA3 phased alignment of longest non-recombining blocks (see GenBank accession numbers in Sup. Table 1).
7_GRIA3_lnrb_final.fas
SPTBN1_with_outgroup
SPTBN1 phased alignment of longest non-recombining blocks with outgroup (see GenBank accession numbers in Sup. Table 1).
1SPTBN1_with_outgroup.fas
PRKCI_with_outgroup
PRKCI phased alignment of longest non-recombining blocks with outgroup (see GenBank accession numbers in Sup. Table 1).
2_PRKCI_with_outgroup.fas
DARC_with_outgroup
DARC phased alignment of longest non-recombining blocks with outgroup (see GenBank accession numbers in Sup. Table 1).
3_DARC_with_outgroup_final.fas
KITLG_with_outgroup
KITLG phased alignment of longest non-recombining blocks with outgroup (see GenBank accession numbers in Sup. Table 1).
4_KITLG_with_outgroup.fas
TF_with_outgroup
TF phased alignment of longest non-recombining blocks with outgroup (see GenBank accession numbers in Sup. Table 1).
5_TF_with_outgroup.fasta
POLA1_with_outgroup
POLA1 phased alignment of longest non-recombining blocks with outgroup (see GenBank accession numbers in Sup. Table 1).
6_POLA1_with_outgroup.fas
GRIA3_with_outgroup
GRIA3 phased alignment of longest non-recombining blocks with outgroup (see GenBank accession numbers in Sup. Table 1).
7_GRIA3_with_outgroup_final.fas
SRY_with_outgroup
SRY alignment with outgroup (see GenBank accession numbers in Sup. Table 1).
8_SRY_with_outgroup.fas
CYTB_with_outgroups
CYTB alignment with outgroups (see GenBank accession numbers in Sup. Table 1).
9_CYTB_with_outgroups.fas
CYTB_Lepus_spp
CYTB alignment of haplotypes obtained in this work and other Lepus species available in GenBank (see Sup. Figure 4).