The Rough Greensnake, Opheodrys aestivus, is a moderately-sized, semi-arboreal snake broadly distributed throughout eastern North America. While numerous taxa with similar distributions have been shown to be comprised of multiple species, O. aestivus has yet to be examined in a detailed phylogeographic context. Here, we use Sanger-sequence data of one mitochondrial and three nuclear loci for samples from throughout the distribution of O. aestivus to elucidate phylogeographic patterns in this species. We combine this with ddRADseq data for a subset of samples to test patterns on a more genomically comprehensive scale. In both datasets, we find strong support for three deeply divergent clades within O. aestivus: peninsular Florida, central Texas, and a main clade comprising the rest of the distribution, with the Florida clade the earliest diverging lineage of the three. Estimates of divergence time suggest that the central Texas and main clades diverged approximately 1.34 million years ago (Mya), while the peninsular Florida clade diverged from other lineages approximately 2.94 Mya, and these lineages diverged from the sister taxon, O. vernalis, approximately 6.43 Mya. These results also suggest that the historically recognized Florida subspecies, O. a. carinatus, could be elevated to species status. While the divergence of peninsular Florida or central Texas populations is not unique among squamates, nor is low levels of divergence from the Atlantic coast to eastern Texas, this combination of patterns is unusual, and yields important insight into the biogeography of North American biota. Further, our approach helps illustrate how dense geographic sampling with limited genomic sequencing can be used as a guide for the selection of samples to test phylogeographic patterns comprehensively.

We obtained 47 samples of Opheodrys aestivus from throughout the distribution of the species, four samples from the sister taxon O. vernalis, and one sample from each of the more distantly related colubrid species Cemophora coccinea, Lampropeltis getula, Pantherophis guttatus, Pituophis melanoleucus, and the natricid Nerodia sipedon (Appendix 1). Whole genomic DNA was extracted for all samples using a Qiagen DNeasy Blood and Tissue kit and quantified using a Qubit 2.0 fluorometer with the High Sensitivity kit (Invitrogen, Inc.). we targeted one mitochondrial (cytochrome b, cytb) and three nuclear loci (LAT clone, LAT; neurotrophin 3, NT3; and prolactin receptor, PRLR). We used maximum likelihood (ML) and Bayesian inference (BI) methods to infer the phylogenetic relationships among populations of O. aestivus. To further evaluate putative species boundaries and clades recovered in the Sanger sequencing analyses in a more genomically comprehensive context, we also used a double-digest restriction associated DNA sequencing (ddRADseq) for a small subset of the samples, representing all clades identified from Sanger-sequence analyses. A subsample of 13 samples of O. aestivus, spanning the distribution of the species and all lineages recovered in Sanger-data analyses, along with a single sample the sister species, O. vernalis.

Refer to the manuscript for descriptions of analyses leading to the supplementary figures here; figure captions for each supplementary figure are included (see read me file for additional details). In addition, the data files include a zipped folder with the raw ddRAD files, the processed ddRAD files, and a text file of the ddRAD stats. A nexus files of each locus, as described in the manuscript is included (as both a phased and unphased nexus for each nuclear locus).