Estimation of hybrid ancestry and population structure
We obtained genomic DNA from 270 trees using DNeasy 96 Plant Kits (Qiagen). These trees included the 30 maternal trees described above along with 154 parental P. strobiformis and 86 parental P. flexili s trees (Menon et al ., 2018; Fig. S2a). Following the ddRADseq library preparation detailed in Parchman et al . (2012), we size selected 300-400 bp regions and performed single-end sequencing (1 \(\times\ 150\) bp) on the Illumina HiSeq 4000 platform (Novogene Corporation) for multiplexed sets of 96 trees. The resulting data were processed using dDocent v.1.0 (Puritz et al ., 2014) as in Menonet al . (2018, 2021) to obtain genotype calls. Following post-filtering steps implemented in VCFTOOLS v.0.1.153 (Danecek et al ., 2011) and custom functions (Methods S2), we obtained 23,261 SNPs. We used this full dataset of 270 individuals and 23,261 SNPs to estimate hybrid ancestry through the genotype likelihood-based clustering algorithm implemented in NGSAdmix (Skotte et al ., 2013), with the number of clusters (K) set to the two hybridizing species (K = 2). We estimated multilocus F ST from 11,432 of the 23,261 SNPs that were polymorphic across the 30 maternal trees using the method of Weir and Cockerham (1984) as implemented in the hierfstat v.0.04-22 package (Goudet, 2005) in R v.4.0.2 (R Core Team, 2020).