Bioinformatics pipeline
Raw sequence data were demultiplexed using the process_radtagsfunction of Stacks (Catchen et al., 2013), with checking for intact RAD sites and reads quality-checked and trimmed. Of the 503.8 million total reads, 99.9% were retained after demultiplexing. Trimmomatic v0.39 (Bolger et al., 2014) was used to remove possible adapter contamination, and reads were aligned to the reference genome assembled in this study with BWA v0.7.17-r1188 (Li & Durbin, 2009). All samples mapped well (mean mapping rate = 95.16%; SD = 0.67, range = 92.83–96.46%). Thegstacks function of Stacks was used to call a catalogue of variants, filtered with the populations function by retaining variants in with a minor allele frequency of 0.01; maximum heterozygosity of 0.8; genotyped in 30% of samples; and retaining only one SNP per tag (–write-random-snp). SNPs were filtered using a custom R script to retain SNPs with a minimum average allelic depth of 2.5× per allele; a coverage difference between alleles of ≤ 80%; a genotyping rate per locus ≥ 80%; and a reproducibility of genotype calls between replicates of 100% (Wright et al., 2019). The mean error rate between replicate pairs was 1.17% (SD = 5.58; range = 0.60–2.19%). Further filtering of the SNPs is outlined in Table S2.