Atractaspidines are poorly studied, fossorial snakes that are found throughout Africa and western Asia, including the Middle East. We employed concatenated gene-tree analyses and divergence dating approaches to investigate evolutionary relationships and biogeographic patterns of atractaspidines with a multi-locus data set consisting of three mitochondrial (16S, cyt b, and ND4) and two nuclear genes (c-mos and RAG1). We sampled 91 individuals from both atractaspidine genera (Atractaspis and Homoroselaps). Additionally, we used ancestral-state reconstructions to investigate fang and diet evolution within Atractaspidinae and its sister lineage (Aparallactinae). Our results indicated that current classification of atractaspidines underestimates diversity within the group. Diversification occurred predominantly between the Miocene and Pliocene. Ancestral-state reconstructions suggest that snake dentition in these taxa might be highly plastic within relatively short periods of time to facilitate adaptations to dynamic foraging and life-history strategies.

Funding: This work was supported by the Percy Sladen Memorial Fund, an IUCN/SSC Amphibian Specialist Group Seed Grant, K. Reed, M.D., research funds from the Department of Biology at Villanova University, and UTEP (all to EG), National Science Foundation grant DEB-1145459 (to EG and KJ), National Geographic Research and Exploration Grant no. 8556-08 (to EG), National Geographic Okavango Wilderness Project no. EC0715-15 (to WC), Belgian National Focal Point to the Global Taxonomy Initiative (to ZTN). MOR is supported by the Gorongosa Restoration Project and the Mozambican Departamento dos Serviços Cientificos (PNG/DSCi/C12/2013; PNG/DSCi/C12/2014; PNG/DSCi/C28/2015). The University of Texas at El Paso (UTEP) Border Biomedical Research Center (BBRC) Genomic Analysis Core Facility is acknowledged for services and facilities provided. This core facility is supported by grant 5G12MD007592 to the (BBRC) from the National Institutes on Minority Health and Health Disparities (NIMHD), a component of the National Institutes of Health (NIH). MFB received payment as a herpetologist consultant with Flora Fauna & Man, Ecological Services Ltd. (FFMES). The funder provided support in the form of salaries for one author [MFB], but did not have any additional role in the study design, data collection and analysis, decision to publish, or preparation of the manuscript. The specific roles of these authors are articulated in the ‘author contributions’ section. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

Herein, we employ phylogenetic hypotheses in conjunction with temporal biogeographic information to gain a more comprehensive understanding of the evolutionary history of Atractaspidinae. Specifically, we evaluate the following questions: Are currently recognized genera and species monophyletic? Are Atractaspis and Homoroselaps sister taxa? Are Atractaspis genetically partitioned into the ‘bibronii’ and ‘microlepidota’ groups as Underwood and Kochva [ 18 ] suggested? Additionally, we investigate patterns of diversification regarding character traits, including prey selection and fang morphology, within atractaspidines and aparallactines.

The diversification of burrowing asps is particularly interesting because of their unique front fangs, which are starkly different from other lamprophiids [ 21 , 29 – 32 ]. It has been hypothesized that foraging for nestling mammalian prey was a major driver in the evolution of front fangs and “side-stabbing,” which are unique to Atractaspis [ 31 , 33 ]. Both Atractaspis and Homoroselaps have front fangs, which differs from the rear-fang morphology that is common in their aparallactine sister group. Although Atractaspsis and Homoroselaps both contain front fangs, Atractaspis fang morphology is more similar to viperids (Atractaspis was previously and erroneously classified in the Viperidae), whereas Homoroselaps fang morphology is more similar to elapids [ 25 , 31 ]. Underwood and Kochva [ 18 ] suggested a Macrelaps-like ancestor for aparallactines and atractaspidines, which may have foraged above ground and fed on a wide variety of prey items. Specialization on elongated prey items (e.g., squamates and invertebrates) may have taken different evolutionary routes within aparallactines and atractaspidines, which involved morphological changes that facilitated foraging, capture, and envenomation of prey items [ 31 ]. Burrowing asps and their sister group Aparallactinae are ideal groups to study fang evolution, because they possess many fang types (i.e., rear fang, fixed front fang, and moveable front fang) [ 25 , 29 – 32 ]. Additionally, collared snakes (aparallactines) and burrowing asps make interesting models to study fang evolution because of their dietary specializations, especially prevalent within the Aparallactinae, which feed on prey ranging from earthworms to blind snakes [ 25 , 31 ].

Based on scale patterns and counts, Laurent [ 28 ] assigned the known species of Atractaspis into five groups (Sections A–E). Decades later, Underwood and Kochva [ 18 ] partitioned Atractaspis into two groups based on venom gland morphology and geographic distribution: the ‘bibronii’ group and the ‘microlepidota’ group. These authors defined the ‘bibronii’ group as having normal-sized venom glands and a sub-Saharan distribution, and it included the following species: A. aterrima, A. bibronii, A. boulengeri, A. congica, A. corpulenta, A. dahomeyensis, A. duerdeni, A. irregularis, and A. reticulata. The 2nd ‘microlepidota’ group has relatively elongated venom glands and is found in western, central and eastern Africa, including the distinctive horn of Africa, the Sinai Peninsula, and much of Arabia, Israel, and the Levant. This latter group consisted of the following species: A. engaddensis, A. engdahli, A. leucomelas, A. microlepidota, A. micropholis, and A. scorteccii. Moyer and Jackson [ 10 ] reconstructed phylogenetic relationships among 14 species of Atractaspis with morphological data, incorporating Macrelaps and Homoroselaps as outgroups, based on previous studies [ 18 ]. However, the two groups of Underwood and Kochva [ 18 ] were not supported [ 10 ]. More recent molecular phylogenetic studies suggest that Homoroselaps is sister to Atractaspis, whereas Macrelaps is closely related to Amblyodipsas and Xenocalamus [ 8 – 9 , 27 ].

Recently, several studies generated phylogenies of advanced African snakes, including colubrids, lamprophiids, elapids, and viperids [ 1 – 9 ]. In contrast, there has been only one morphology-based, phylogenetic study that focused on atractaspidines [ 10 ]. The Family Atractaspididae was originally erected by Günther [ 11 ] for species of Atractaspis, renowned for their unique and exceptionally long and mobile fangs [ 12 ]. Based on skull morphology, Bourgeois [ 13 ] created the subfamily Aparallactinae (within Colubridae) to accommodate Atractaspis, Aparallactus, and other closely related fossorial snakes. This grouping was supported by jaw musculature studies of Heymans [ 14 – 15 ], who transferred Atractaspis to the Subfamily Atractaspidinae (Atractaspininae, sensu Kelly et al. [ 16 ]). Several recent molecular [ 7 – 9 ] and morphological studies [ 17 – 18 ] recovered a monophyletic group containing both aparallactines and atractaspidines, and with few exceptions [ 19 – 21 ], current classification recognizes Aparallactinae and Atractaspidinae as sister taxa in the Family Lamprophiidae [ 2 , 7 – 9 , 22 – 25 ]. Phylogenetic relationships within atractaspidines are not well known, because many phylogenetic studies that included atractaspidines were limited by low sample sizes [ 2 , 8 – 10 , 21 – 23 , 26 – 27 ].

Microcomputed tomography (CT) scans of specimens were produced using GE Phoenix V|Tome|X systems at the General Electric Sensing & Inspection Technologies in Scan Carlos, CA and University of Florida’s Nanoscale Research Facility. X-ray tube voltage and current, detector capture time, voxel resolution, and projection number were optimized for each specimen ( S1 File ). The radiographs were converted into tomograms with Phoenix Datos| R, and then rendered in three dimensions with volumetric rendering suite VGStudioMax 3.2 ( http://www.volumegraphics.com ). Tomogram stacks and 3D mesh files for all scans are available on Morphosource.org ( S1 File ).

To understand the evolution of fang morphology and diet selection in atractaspidines, we reconstructed the pattern of character changes on the ML phylogeny herein. For ancestral-state reconstructions, we included all samples of aparallactines and atractaspidines available to us in order to better characterize fang and diet characters. All ancestral-state reconstructions were conducted by tracing characters over trees in Mesquite v3.10 [ 38 ]. We scored taxa using descriptions from the literature [ 25 , 30 – 31 , 45 – 55 ], and from our own data. We evaluated the following characters for fang morphology and diet selection: A. Fang morphology: (0) no fang, (1) rear fang, (2) fixed front fang, (3) moveable front fang, and (4) rear-front fang intermediate (anterior half of the maxilla, but not the anteriormost tooth); B. prey selection (0) rodents, (1) rodents, snakes, fossorial lizards, and amphibians, (2) snakes, (3) amphisbaenians, (4) snakes and fossorial lizards, (5) invertebrates, and (6) fish and amphibians. A ML approach was used for both analyses, because it accounts for and estimates probabilities of all possible character states at each node, thus providing an estimate of uncertainty [ 56 ]. A Markov K-state one-parameter model (Mk-1; [ 57 ]) that considers all changes as equally probable was implemented in our ancestral-state reconstructions. States were assigned to nodes if their probabilities exceeded a decision threshold; otherwise nodes were recovered as equivocal.

The program BEAST v1.8.3 via Cipres Science Gateway [ 40 ] was used to estimate divergence times across atractaspidine phylogenetic estimates. The five-gene data set was used to estimate divergence dates in BEAST. Substitution and clock models were unlinked for all partitions; trees were unlinked across the nuclear loci, but were linked for the two mitochondrial partitions because these evolve as a single unit. We implemented an uncorrelated log-normal relaxed clock model with a Yule tree prior. Two independent analyses were run for 100 million generations, sampling every 10,000 generations. Primary calibration points were obtained from Head et al. [ 44 ] and a secondary calibration point was obtained from Kelly et al. [ 7 ] including: the split between Scolecophidia and all other snakes (120–92 mya); split between Caenophidia and its nearest sister taxon, Booidea (72.1–66 mya); split between Colubroidea and its nearest sister taxon (Acrochordus + Xenodermatidae) (72.1–50.5 mya); the divergence of Colubridae + Elapoidea (30.9 ± 0.1 mya); and the split between Crotalinae and Viperinae (23.8–20.0 mya). All calibrations were constrained with a log-normal mean of 0.01, a normal standard deviation of 2.0 (first calibration point), and 1.0 (the last four calibration points). Parameter values of the samples from the posterior probabilities on the maximum clade credibility tree were summarized using the program TreeAnnotator v1.8.3 via Cipres Science Gateway [ 40 ].

Phylogenetic analyses were conducted for our individual and five-gene concatenated data sets. Data were interpreted using the program SeqMan [ 36 ]. An initial alignment for each gene was produced in MUSCLE [ 37 ] in the program Mesquite v3.10 [ 38 ], and manual adjustments were made in MacClade v4.08 [ 39 ]. The Maximum Likelihood (ML) analyses of single gene and concatenated data sets were conducted using the GTRGAMMA model in RAxML v8.2.9 via the Cipres Science Gateway v3.3 [ 40 ]. All parameters were estimated, and a random starting tree was used. Support values for clades inferred by ML analyses were assessed with the rapid bootstrap algorithm with 1,000 replicates [ 40 ]. We also conducted Bayesian inference (BI) analyses with MrBayes v3.2.6 via the Cipres Science Gateway [ 40 ]. The model included 13 data partitions: independent partitions for each codon position of the protein-coding genes ND4, cyt b, c-mos, and RAG1, and a single partition for the mitochondrial gene 16S. Phylogenies were constructed based on concatenated data, which included 16S and the four protein-coding genes listed above. Concatenated data sets were partitioned identically for ML and BI analyses. The program PartitionFinder v1.1.1 [ 41 – 42 ] was used to find the model of evolution that was most consistent with our data for BI analyses. Bayesian analyses were conducted with random starting trees, run for 20,000,000 generations, and sampled every 1000 generations. Phylogenies were visualized using FigTree v1.3.1 [ 43 ].

Genomic DNA was isolated from alcohol-preserved muscle or liver tissue samples with the Qiagen DNeasy tissue kit (Qiagen Inc., Valencia, CA, USA). Primers used herein are shown in Table 2 . We used 25 μL PCR reactions with gene-specific primers with an initial denaturation step of 95°C for 2 min, followed by denaturation at 95°C for 35 seconds (s), annealing at 50°C for 35 s, and extension at 72°C for 95 s with 4 s added to the extension per cycle for 32 (mitochondrial genes) or 34 (nuclear gene) cycles. Amplification products were visualized on a 1.5% agarose gel stained with SYBR Safe DNA gel stain (Invitrogen Corporation, Carlsbad, CA, USA). Sequencing reactions were purified with CleanSeq magnetic bead solution (Agencourt Bioscience, La Jolla, CA) and sequenced with an ABI 3130xl automated sequencer at the University of Texas at El Paso (UTEP) Genomic Analysis Core Facility.

Specimens from the Subfamily Atractaspidinae were collected from multiple localities in sub-Saharan Africa ( Fig 1 ). We generated sequences of three mitochondrial genes (16S, ND4, and cyt b) and two nuclear genes (c-mos and RAG1) for 91 atractaspidine individuals (Tables 1 and 2 ). This study included sequences from both atractaspidine genera (14/22 species of Atractaspis; 2/2 species of Homoroselaps) [ 24 , 34 ]. Sequences from some of these individuals have been published previously [ 2 , 7 ], and new sequences were deposited in GenBank ( Table 1 ). Concatenated trees were rooted with Acrochordus granulatus (not shown on Fig 2 ). Three genera of Viperidae (Agkistrodon, Atheris, and Crotalus; not shown on Fig 2 ), two genera of Elapidae (Naja and Dendroaspis), six genera of Lamprophiinae (Boaedon, Bothrophthalmus, Bothrolycus, Gonionotophis, Lycodonomorphus, and Lycophidion), Psammophylax, and Micrelaps were used as outgroups for the concatenated analyses ( Table 1 , Fig 2 ). Additionally, we included sequences from six of the eight known aparallactine genera (6/9 species of Amblyodipsas; 7/11 species of Aparallactus; 1/2 species of Chilorhinophis; 1/1 species of Macrelaps; 7/14 species of Polemon; 4/5 species of Xenocalamus) [ 24 , 35 ] for concatenated analyses and ancestral-state reconstructions. For divergence-dating analyses, additional samples from the squamate taxa Scincidae, Leptotyphlopidae, Viperidae, Colubrinae, and Dipsadinae were included ( Table 1 ).

Permission for DFH, MB and EG to collect snakes in Uganda was obtained from the Uganda Wildlife Authority (UWA—permit no. 2888 issued on August 1, 2014, permit no. 29279 issued on August 11, 2015) and the Ministry of Tourism, Wildlife and Antiquities (permit no. GoU/008/2016). Permission for CK, WMM, MMA, and EG to collect snakes in Burundi was granted by the Institut National pour l’Environnement et la Conservation de la Nature (INECN—unnumbered permit from Directeur General de l’INECN dated December 27, 2011). Permission for CK, WMM, MMA, DFH, and EG to collect snakes in Democratic Republic of Congo (DRC) was granted by the Centre de Recherche en Sciences Naturelles (CRSN—LW1/28/BB/MM/BIR/050/07, unnumbered permit from 2008, LWI/27/BBa/MUH.M/BBY/141/09, LWI/27/BBa/MUH.M/BBY/023/10, LWI/27/BBa/MUH.M/BBY/001/011, LWI/27/BBa/CIEL/BBY/003/012, LW1/27/BB/KB/BBY/60/2014, LWI/27/BBa/BBY/146/014), Institut Congolais pour la Conservation de la Nature (ICCN—unnumbered permit by Provincial Director of ICCN, Equateur Province in Mbandaka in August 2013, 004/ICCN/PNKB/2013, 06/ICCN/PNKB/2014, 02/ICCN/PNKB/2015), and Institut Superieur d’Ecologie Pour la Conservation de la Nature (ISEC, Katana—ISEC/DG/SGAC/04/2015, ISEC/DG/SGAC/04/29/2016). The University of Texas at El Paso (UTEP) Institutional Animal Care and Use Committee (IACUC—A-200902-1) approved field and laboratory methods. Permits for WC to collect snakes in South Africa were granted by the Department of Economic Development, Environmental Affairs and Tourism (permit nos. CRO 84/11CR and CRO 85/11CR). Permits for MOR and JP to collect snakes in Mozambique were granted by the Gorongosa Restoration Project and the Mozambican Departamento dos Serviços Cientificos (PNG/DSCi/C12/2013; PNG/DSCi/C12/2014; PNG/DSCi/C28/2015). Additional specimens and samples were obtained from natural history museums and university collections ( Table 1 ) that followed appropriate legal guidelines and regulations for collection and loans of specimens.

For the analyses with diet data, likelihood reconstructions inferred a generalist diet of rodents, reptiles, and amphibians for the ancestral condition of all lamprophiids (99.7%) ( Fig 6[B] ). Several lamprophiines (Lycodonomorphus) subsequently adopted a more specialized diet of amphibians, reptiles, and fish. The common ancestor for aparallactines and atractaspidines was inferred to have a generalist diet of rodents, reptiles, and amphibians (92.4%). Results recovered a more specialized ancestral diet of snakes and lizards (64.5%) for aparallactines, which was favored over a generalist diet (27.7%). The condition of a snake and lizard diet (79.9%) was favored over a generalist diet (16.2%) for the ancestor of Polemon/Chilorhinophis and Amblyodipsas/Macrelaps/Xenocalamus. The latter dietary condition was retained for the ancestor of Polemon/Chilorhinophis (79.4%) and the ancestor of Amblyodipsas/Macrelaps/Xenocalamus (87.6%). Specialized dietary conditions were recovered for the genera Aparallactus (centipedes and other invertebrates, 99.7%), Polemon (snakes, 97.8%), and Xenocalamus (amphisbaenians, 98.8%). Results suggested a generalist diet for Atractaspidinae (92.3%). The ancestor of Homoroselaps was inferred to have a diet consisting of mostly lizards and snakes (99.9%), whereas the ancestor of Atractaspis was inferred to have a broader diet of rodents, reptiles, and amphibians (99.2%).

X-ray computer tomography of collared snakes and burrowing asps can be seen in Figs 3 and 5 . Likelihood reconstructions of atractaspidine ancestral fang morphology inferred a rear fang condition for the ancestral condition of all lamprophiids (96.7%) ( Fig 6[A] ). Subsequently, the Subfamily Lamprophiinae lost a venom delivery fang condition. The common ancestor of aparallactines and atractaspidines was inferred to have a rear fang condition (97.8%). The analyses suggested a rear fang ancestor (72.5%) for the clade containing Homoroselaps and Atractaspis. The ancestor to Atractaspis was inferred to have a moveable front fang condition (97.4%). Results recovered a fixed front fang condition for the ancestor of all Homoroselaps (99.8%). The ancestor to all aparallactines was inferred to have a rear fang condition (99.6%), and this remained consistent throughout most aparallactine nodes with the exception of Polemon (rear/front fang intermediate, 97.8%) and Aparallactus modestus (no specialized fang, 99.7%).

Topologies from the BEAST ( Fig 4 ) analyses were mostly consistent with the results from our concatenated tree analyses (Figs 2 and 3 ). BEAST results recovered A. corpulenta corpulenta/A. corpulenta kivuensis as sister to A. congica/A. dahomeyensis with strong support (Figs 2 – 4 ). Additionally, the relationship between Atractaspis irregularis and A. corpulenta/A. congica/A. dahomeyensis was strongly supported in BEAST analyses ( Fig 4 ). Results from dating analyses suggested atractaspidines split from aparallactines during the early Oligocene around 29 mya (24.8–31.4 mya, 95% highest posterior densities [HPD]) ( Table 3 , Fig 4 ), which is similar to the results (34 mya) of Portillo et al. [ 62 ]. Subsequently, Atractaspis split from Homoroselaps in the mid-Oligocene, and most radiation events within each of the major clades associated with these genera occurred during the mid- to late Miocene and Pliocene ( Fig 4 ). Specific dates with ranges are specified in Table 3 .

The following models of nucleotide substitution were selected by PartitionFinder for BI analyses: 16S (GTR+G), ND4 1 st codon position (GTR+G), ND4 2 nd codon position (TVM+G), and ND4 3 rd codon position (HKY+I+G); cyt b 1 st codon position (TVM+G), cyt b 2 nd codon position (HKY+I+G) and cyt b 3 rd codon position (GTR+G); c-mos and RAG1 1 st , 2 nd and 3 rd codon positions (HKY+I). Preferred topologies for the ML and BI analyses were identical, with similar, strong support values for most clades ( Fig 2 ), and single-gene mtDNA analyses recovered similar topologies (not shown). The ML analysis likelihood score was –46340.867388. The relationships of Elapidae, Lamprophiinae, Micrelaps, and Psammophylax with respect to the ingroup Atractaspidinae, were not strongly supported in ML and BI analyses. However, Atractaspidinae was recovered in a strongly supported clade. Atractaspis and Homoroselaps were strongly supported as sister taxa ( Fig 2 ). The genus Homoroselaps was recovered as a monophyletic group, and H. lacteus was partitioned into several well-supported clades. There were several strongly supported clades within Atractaspis: (1) Atractaspis andersonii, (2) Atractaspis aterrima, (3) A. bibronii, (4) A. bibronii rostrata, (5) A. cf. bibronii rostrata, (6) A. boulengeri, (7) A. congica, (8) A. corpulenta corpulenta, (9) A. corpulenta kivuensis, (10) A. dahomeyensis, (11) A. duerdeni, (12) A. engaddensis, (13) A. irregularis, (14) A. cf. irregularis, (15) A. reticulata heterochilus, and (16) A. microlepidota. There was strong support for a western Asia/Middle East and Africa clade containing A. andersonii, A. engaddensis, A. microlepidota, A. micropholis, A. watsoni, and A. sp. Atractaspis andersonii did not form a monophyletic group, because one of the samples from Oman (AF471127) was recovered as sister to a clade of A. engaddensis with strong support ( Fig 2 ). The western African species A. aterrima was recovered with strong support as sister to a clade containing A. reticulata heterochilus and A. boulengeri. Atractaspis corpulenta kivuensis samples from eastern DRC were strongly supported as sister to A. corpulenta from northwestern Republic of Congo (near Gabon, the type locality). A well-supported clade of Atractaspis irregularis samples was partitioned by strongly supported central (A. cf. irregularis) and western African (A. irregularis) subclades. Atractaspis duerdeni was recovered within a well-supported A. bibronii complex. Atractaspis bibronii rostrata samples were partitioned into two highly divergent clades from southeastern DRC and Tanzania/Mozambique.

Our data set consisted of 3933 base pairs (16S [546 bp], ND4 [679 bp], cyt b [1094 bp], c-mos [605 bp], and RAG1 [1009 bp]). Individuals with missing data were included in the concatenated sequence analyses, because placement of individuals that are missing a significant amount of sequence data can be inferred in a phylogeny, given an appropriate amount of informative characters [ 8 , 58 – 60 ]. Furthermore, Jiang et al. [ 61 ] showed that excluding genes with missing data often decreases accuracy relative to including those same genes, and they found no evidence that missing data consistently bias branch length estimates.

4. Discussion

4.1 Biogeography Atractaspidines are distributed throughout sub-Saharan Africa except for three species of Atractaspis that are found in western Asia/Middle East (Atractaspis andersonii, A. engaddensis, and A. microlepidota) [25, 29–31]. Based on our results, the most likely scenario for Atractaspis is an African origin with a vicariance or dispersal event into the western Asia/Middle East region in the late Miocene (Fig 4). Atractaspis from western Asia/Middle East and Africa last shared a common ancestor during the late Miocene around 12.1 mya (7.8–17.6). Other studies of African-western Asian/Middle Eastern complexes (e.g., Echis and Uromastyx) recovered similar dates during the late Miocene, with the Red Sea proving to be a strong biogeographic barrier [63–69]. However, lineages of Varanus from Africa and the Middle East split from each other 6.9 mya [70], and African and Middle Eastern Bitis arietans last shared a common ancestor around 4 mya [64]. These dating estimates suggest that there were multiple dispersal events, which were taxon specific. Many Middle Eastern amphibians and reptiles have common ancestors in the Horn of Africa [63–71]. Our study lacked multiple Atractaspis species from the Horn of Africa, and future studies should include samples of A. fallax, A. magrettii, A. leucomelas, and A. scorteccii to improve understanding of likely Africa–Asia biogeographic patterns in atractaspidines. Atractaspis began to diversify around the mid-Oligocene simultaneously with many aparallactine genera [62]. Many of the modern species split from recent common ancestors during the mid- to late Miocene (Table 3, Fig 4). The late Miocene was characterized by considerable xeric conditions, which led to the expansion of savannas globally [72–73]. Other studies on Central and East African herpetofauna, including squamates (Adolfus, Atheris, Boaedon, Naja, Kinyongia, and Panaspis) and frogs (Amietia, Leptopelis, and Ptychadena), have shown similar trends of species diversification during the late Miocene [3–5, 62, 74–78]. The diversification of several western and central African Atractaspis was most likely a consequence of increasingly xeric conditions during the Miocene, when forest and other moist habitats were fragmented [72]. These Atractaspis were likely isolated in fragmented patches of forest during the mid- to late Miocene. Atractaspis irregularis is partitioned clearly by western African and central African lineages that diverged in the mid-Miocene, similar to Aparallactus modestus [62]. At this time, southern African and Middle Eastern Atractaspis also diversified. Atractaspis from the Near and Middle East (A. andersonii, A. engaddensis, and A. microlepidota) and southern Africa (A. bibronii and A. duerdeni) are not tropical forest species, and they inhabit deserts or semi-desert savannas and dry woodland [30, 79–80]. This adaptation to more xeric and open habitats would have allowed Near and Middle Eastern, and southern African Atractaspis, to disperse into these habitats during the dry conditions of the mid- to late Miocene. Studies on mammals and birds show most diversification events during the Pliocene [81–84], which is consistent with the timing of diversification for Atractaspis aterrima, A. congica, A. dahomeyensis, and populations of South African A. bibronii (Fig 4). In contrast to Aparallactus jacksonii, Atractaspis bibronii rostrata showed no clear genetic partitioning between populations in the Nguru, Usambara, and Udzungwa Mountains [62]. Aparallactus jacksonii clearly exhibited deep divergence between an extreme northern Tanzanian population, and a population from the Nguru Mountains. These two populations diverged from each other during the late Miocene, suggesting that the habitats of this taxon were fragmented with increased aridity [62]. Other vertebrate taxa that have shown substantial divergences between populations found in extreme northern Tanzania (Usambara, Taita, and Pare Mountains) and those slightly south (Uluguru, Ukaguru, Nguru, and Malundwe Mountains), include the reed frog Hyperolius puncticulatus, the green barbet (Stactolaema olivacea), and the streaky canary (Serinus striolatus) [82, 85]. But like Atractaspis bibronii rostrata, the hyperoliid reed frog Hyperolius spinigularis and the aparallactine Aparallactus guentheri showed no clear biogeographic patterns between populations in different areas of the Eastern Arc Mountains. These results support the hypothesis that the evolutionary history of species from the Eastern Arc Mountains is lineage specific [85]. Atractaspis bibronii rostrata inhabit low-elevation woodlands and grasslands, and transitional habitats, rather than montane forest (i.e., Aparallactus jacksonii) [25]. This would allow taxa such as Atractaspis bibronii rostrata to continuously disperse between the different mountains of the Eastern Arcs, despite increased aridity. Additionally, ecological niche requirements may also explain the different biogeographic patterns seen in Aparallactus jacksonii and Atractaspis bibronii rostrata. Atractaspis bibronii has a generalist diet (mammals, squamates, and amphibians) and could have exploited more habitats than Aparallactus jacksonii, which is a centipede specialist [25].

4.2 Evolutionary relationships and taxonomy of Atractaspidinae Our results indicate that both Atractaspis and Homoroselaps are strongly supported as monophyletic sister taxa. Results from Figueroa et al. [27] recovered a monophyletic group containing aparallactines and atractaspidines, but their results did not recover a monophyletic Atractaspis (A. irregularis was recovered as sister to aparallactines + atractaspidines). This sample was excluded from our analyses, because the only sequence available for this taxon was from BDNF, a gene not used herein. The results from Figueroa et al. [27] may be an artifact of sample size of atractaspidines, or incomplete lineage sorting of the BDNF nuclear gene. Results from our study indicate that A. irregularis is a monophyletic lineage within a strongly supported, monophyletic Atractaspis. Underwood and Kochva [18] recognized two groups within Atractaspis: (1) the ‘bibronii’ group (represented in our study by A. aterrima, A. bibronii, A. boulengeri, A. congica, A. corpulenta, A. dahomeyensis, A. irregularis, and A. reticulata), characterized by a single posterior supralabial, three anterior infralabials, normal-sized venom glands, and a sub-Saharan distribution; and (2) the ‘microlepidota’ group (represented in our study by A. andersonii, A. engaddensis, A. microlepidota, and A. micropholis), characterized by two anterior temporals, highly elongated venom glands, and a North African/Near and Middle Eastern distribution. Whereas our study did not include genetic samples of all known species of Atractaspis, results herein (Fig 2) support partitioning of the genus into two groups sensu Underwood and Kochva [18]. Our results indicated a clear partition between a ‘Middle Eastern + African’ clade (including A. watsoni, a species that was not included by Underwood and Kochva [18]) and a ‘sub-Saharan African’ clade (Figs 2 and 4). These results strengthen the notion that venom gland size and length in Atractaspis are homologous. Our support for the ‘microlepidota’ group is consistent with the “Section A” (A. andersonii, A. fallax, A. leucomelas, A. microlepidota, and A. micropholis) of Laurent [28] and the A. micropholis/A. microlepidota/A. watsoni clade recovered by Moyer and Jackson [10]. However, our phylogeny (Fig 2) contrasts with the remaining “sections” of Laurent [28], most relationships depicted in the morphological phylogeny of Moyer and Jackson [10], and the molecular phylogenies of Pyron et al. [8–9] and Vidal et al. [22]. Based on relatively long branch lengths, several lineages of Atractaspis seem to be cryptic complexes of species. Because of the extensive geographic distribution of A. bibronii in central, eastern and southern Africa, it is unsurprising to find several highly divergent lineages that likely represent cryptic species. Given the proximity (ca. 167–333 km) of our Tanzanian localities of A. bibronii rostrata (Nguru, Usambara, and Udzungwa Mountains) to the insular type locality for this taxon (Zanzibar, Tanzania), the morphological similarity between our voucher specimens and the types [86], and the relatively long branch length and reciprocal monophyly of this clade compared to topotypic South African A. bibronii (Fig 2), it is likely that the former taxon is a valid species. However, additional comparisons to type specimens are needed to clarify the taxonomic status of populations in this clade, including samples from Haut-Katanga Province in southeastern DRC. Our phylogenetic results indicated that several other species, including A. andersonii, A. boulengeri, A. congica, A. corpulenta, A. dahomeyensis, and A. irregularis likely represent more than a single species. For example, topotypic Angolan samples of A. congica are deeply divergent from our eastern DRC sample (Fig 2), which is likely attributable to A. congica orientalis [46]. Like Polemon fulvicollis fulvicollis (Gabon) and P. fulvicollis laurenti (DRC) [62], Gabonese Atractaspis corpulenta and eastern DRC populations of A. corpulenta kivuensis also showed marked genetic divergences between each other (Fig 2). The well-supported clade of A. irregularis from western Africa likely includes topotypic populations, because they straddle the type locality (Accra, Ghana) [87], whereas our Albertine Rift samples are likely attributable to one of the taxon’s many synonyms. One of these, Atractaspis bipostocularis from Mt. Kenya, was named for its two postocular scales, which distinguishes it from the single postocular of topotypic A. irregularis [88]. Because Mt. Kenya is located east of the Kenyan Rift, a major biogeographic barrier to several species of squamates [78], and moreover, all voucher specimens of A. cf. irregularis from the Albertine Rift have a single postocular (EG pers. obs.), A. bipostocularis is likely a distinct species that is endemic to the central Kenyan highlands. Other synonyms of A. irregularis that have one postocular and type localities in or near the Albertine Rift are likely attributable to our well-supported clade of A. cf. irregularis (Fig 2 in [87]), and include the following taxa: A. conradsi Sternfeld, 1908 (type locality: Ukerewe Island, Lake Victoria, Tanzania [89]), A. schoutedeni de Witte, 1930 (type locality: Goma, North Kivu, DRC [90]), A. babaulti Angel, 1934 (type locality: Kadjuju [1500 m elevation] on the western border of Lake Kivu, 15 km north of Katana, DRC [91]), and A. irregularis loveridgei Laurent, 1945 (type locality: Bunia, DRC [46]). Additional sampling and morphological analyses are in progress that will help clarify the correct taxonomy for these lineages. Because of the relative lack of fieldwork in Central Africa in recent decades [92–93] and the relatively rare encounters of these snakes above ground (EG, pers. obs.), it is likely that genetic samples from the above topotypic populations will remain elusive for many years.