Ultraconserved elements-based phylogenomic systematics of the snake superfamily Elapoidea, with the description of a new Afro-Asian family

The highly diverse snake superfamily Elapoidea is considered to be a classic example of ancient, rapid radiation. Such radiations are challenging to fully resolve phylogenetically, with the highly diverse Elapoidea a case in point. Previous attempts at inferring a phylogeny of elapoids produced highly incongruent estimates of their evolutionary relationships, often with very low statistical support. We sought to resolve this situation by sequencing over 4,500 ultraconserved element loci from multiple representatives of every elapoid family/sub-family level taxon and inferring their phylogenetic relationships with multiple methods. Concatenation and multispecies coalescent based species trees yielded largely congruent and well-supported topologies. Hypotheses of a hard polytomy were not retained for any deep branches. Our phylogenies recovered Cyclocoridae and Elapidae as diverging early within Elapoidea. The Afro-Malagasy radiation of elapoid snakes, classified as multiple subfamilies of an inclusive Lamprophiidae by some earlier authors, was found to be monophyletic in all analyses. The genus Micrelaps was consistently recovered as sister to Lamprophiidae. We establish a new family, Micrelapidae fam. nov. , for Micrelaps and assign Brachyophis to this family based on cranial osteological syn-apomorphy. We estimate that Elapoidea originated in the early Eocene and rapidly diversified into all the major lineages during this epoch. Ecological opportunities presented by the post-Cretaceous-Paleogene mass extinction event may have promoted the explosive radiation of elapoid snakes.


Introduction
Elapoidea is ecologically and morphologically-one of the most diverse superfamilies of advanced or caenophidian snakes.The superfamily consists of the cosmopolitan family Elapidae, which contains numerous medically significant venomous snakes, and additional lineages that have been treated as family or subfamily-level taxonomic groups, including venomous, mildly venomous, and non-venomous species (Kelly et al., 2009;O'Shea, 2018).The non-elapid subclades of Elapoidea reach their highest diversity in Africa and Madagascar, but some members occur in Asia and southern Europe and one subclade is endemic to the Philippines (Kelly et al., 2009;Weinell and Brown, 2018).
Non-elapid elapoids have historically been assigned to the 'catch-all' family Colubridae (e.g., Boulenger, 1893Boulenger, , 1894Boulenger, , 1896;;FitzSimons, 1912), except the unusual burrowing asps which have been variably placed either in Viperidae (Boulenger, 1896) or their own family (Günther, 1858).Starting with the landmark work of Bogert (1940), systematic revisions increasingly discovered that many of these snakes, especially the African and Madagascan genera, are more closely related to each other than they are to other colubrids (Bourgeois, 1968;McDowell, 1987).The African and Malagasy genera, then assigned to Colubridae, rendered Colubridae non-monophyletic with respect to Elapidae and Atractaspis in the serum albumin-based phylogenies estimated by Cadle (1994).He was the first to suggest a plausible rapid radiation involving the basal nodes.Early DNA sequence-based phylogenies further resolved the relationships of elapids to the rest of the Afro-Malagasy elapoid radiation (Kraus and Brown, 1998;Gravlund, 2001;Vidal and Hedges, 2002).
Although most previous molecular phylogenetic studies, with dense sampling of elapoid species, supported the monophyly of currently recognized family/subfamily level lineages, relationships among families/subfamilies have been highly conflicting (Lawson et al., 2005;Nagy et al., 2005;Vidal et al., 2007;Kelly et al., 2009;Zaher et al., 2009Zaher et al., , 2019;;Pyron et al., 2011Pyron et al., , 2013Pyron et al., ,2014;;Figueroa, et al., 2016;Zheng and Wiens, 2016).None of these studies produced a phylogeny with a topology that was consistently corroborated in another work.Additionally, some genera such as Micrelaps, Buhoma and Psammodynastes, could not be assigned to any family/subfamily level taxa and were recovered at different positions in the elapoid phylogeny in different studies (Kelly et al., 2009;Pyron et al., 2013;Zaher et al., 2019).A major source of disagreement seems to be the question of reciprocal monophyly of the Afro-Malagasy radiation and the Elapidae.Whereas some authors (e.g., Vidal et al., 2007;Pyron et al., 2011Pyron et al., , 2013Pyron et al., , 2014) ) have recovered a monophyletic Afro-Malagasy radiation (with respect to Elapidae family, within the Elapoidea superfamily), others (e.g., Kelly et al., 2009;Zaher et al., 2019) found Elapidae to be nested within the Afro-Malagasy radiation, rendering the latter paraphyletic.This has led to conflicts in the classification schemes proposedauthors recovering a monophyletic Afro-Malagasy radiation recognised two families: an inclusive Lamprophiidae (for the Afro-Malagasy subclades, treated as subfamilies) and Elapidae (e.g., Vidal et al., 2007;Pyron et al., 2011;Pyron et al., 2013) and authors inferring a paraphyletic Afro-Malagasy radiation (with respect to Elapidae) classified each elapoid subclade as a family (e.g., Kelly et al., 2009;Zaher et al., 2009Zaher et al., , 2019)).A general weakness of some of those conflicting phylogenies has been a lack of statistical support for the resolution of deeper internodes.Except for the work of Pyron et al. (2014) that generated a dataset of 333 loci, these studies used few (mostly less than ten) Sanger-sequenced loci, many of which were mitochondrial genes.
Very short branches at the base of the elapoid phylogeny, closely spaced Paleogene divergence times for major elapoid subclades, and poor branch support led Kelly et al. (2009) to postulate that elapoid diversity has been a result of an ancient, rapid radiation, in line with the hypothesis of Cadle (1994).Many groups of organisms have experienced a rapid phase of cladogenesis in their evolutionary history, with the 'Cambrian explosion' of Metazoa (Lee et al., 2013;Briggs, 2015), the mostly post K/Pg radiations of birds (Reddy et al., 2017), placental mammals (McCormack et al., 2012), and snakes (Klein et al., 2021), being prominent examples.Ancient, rapid radiations have always been a recalcitrant stumbling block to phylogenetic systematists because of a lack of phylogenetic signal on the short, ancestral branches, saturation along the longer descendant branches and pervasive incomplete lineage sorting on the short branches (Rokas and Carroll, 2006;Whitfield and Lockhart, 2007;Whitfield and Kjer, 2008).In recent years, usage of genome-scale data from various reduced representation genome sequencing techniques, especially target capture, to solve ancient rapid radiation has produced promising results (viz. McCormack et al., 2012;Longo et al., 2017;Léveillé-Bourret et al., 2018).
A robust phylogeny is an indispensable necessity for the study of any macroevolutionary (and often also microevolutionary) processes.Elapoidea contains multiple medically significant venomous snakes in different subclades (Spawls & Branch, 2020) and taxonomic stability, itself stemming from a well-supported phylogeny, is desirable for further work towards a better understanding of venom, venom apparatus acquisition during evolution and snakebite epidemiology.With the two principal goals of reconstructing a robust phylogeny, and revising higher-level taxonomy to reflect our phylogeny, in mind, we sequenced ultra-conserved elements of elapoid genomes to infer the phylogeny.Specifically, we aimed to resolve the problematic, deeper divergences.We examined the inferred phylogenies to better understand potential sources of conflict in the data and probable presence of any multifurcating cladogenesis with various exploratory methods.To look for potential diagnostic characters for the recognised subclades and any novel higher taxa, we utilised cranial osteological data from µ-CT scans.A secondary objective was to estimate a time calibrated phylogeny as previous elapoid phylogenies have been highly unstable and therefore a reappraisal of the tempo of elapoid evolution, hypothesised to be a rapid, Eocene radiation, appears necessary.

Taxon sampling and molecular laboratory work
We sampled 37 genera representing all the eight major elapoid subclades that have been accorded family/subfamily rank, namely Atractaspididae [including representatives of Aparallactinae and Atractaspidinae], Cyclocoridae/Cyclocorinae, Elapidae, Lamprophiidae, Prosymnidae/Prosymninae, Psammophiidae/Psammophiinae, Pseudaspidae/Pseudaspinae and Pseudoxyrhophiidae/ Pseudoxyrhophiinae, and one currently unnamed clade containing the genus Micrelaps.For several genera we included multiple congeneric species.In total, we sampled 45 species assigned to Elapoidea.We used two colubrid species (representing two genera), and one viperid species, as outgroups.Tissue samples originated from natural history collections and tissue collections of university researchers (Supplementary material).
Total genomic DNA was extracted with Macherey Nagel NucleoSpin Tissue kit (Düren, Germany) and extract quality was checked with 1 % agarose gel electrophoresis.Quantifications and further quality checks were carried out with Nanodrop and Qubit.

Bioinformatic processing of raw sequencing data
Phyluce pipeline (Faircloth, 2016) was used to process raw target capture sequencing data in their analysable form.Illumiprocessor (Faircloth, 2013), a tool in the pipeline that utilises Trimmomatic (Bolger et al., 2014), was used for adapter trimming.In same pipeline, we assembled the trimmed reads into contigs with SPAdes assembler (Bankevich et al., 2012).Tetrapods-UCE-5Kv1 probe set sequences were used in the phyluce pipeline first to find the UCE loci in the assembled contigs and eventually to extract those loci.Alignment was performed in the pipeline with MAFFT (Katoh and Standley, 2013).We assembled three datasets of 50, 75 and 95 % levels of completeness.

Phylogenomic analyses
Individual Maximum Likelihood (ML) gene trees were inferred for all individual UCE loci, in IQ-TREE (Minh et al., 2020b).Models of sequence evolution at each UCE locus were chosen with ModelFinder (Kalyaanamoorthy et al., 2017) as implemented in IQ-TREE.ASTRAL species tree inference has been shown to perform very well with IQ-TREE estimated gene trees, with automated selection of model (Bossert et al., 2021).Ultrafast bootstrapping (UFBoot; Hoang et al., 2018) was performed with 1000 replicates on each gene tree.
Multispecies coalescent phylogeny estimation was performed with the ASTRAL algorithm (Mirarab et al., 2014) which accounts for deep coalescence by using unrooted quartet subtrees induced by the gene trees to find the species tree that agrees with the maximum number of such quartets.We inferred one set of species trees with ASTRAL 5.7.8 using the ML gene trees which were not modified in any way.Another set of species trees were computed with ASTRAL on ML gene trees whose very low support branches (≤ 20 % ultrafast bootstrap support) were contracted into a polytomy with Newick utilities (Junier and Zdobnov, 2010).Contraction of poorly supported nodes into inclusive polytomies has been shown to improve accuracy in phylogenetic inference (Zhang et al., 2017).Branch support was determined with local posterior probability (localPP).
Recently, Zhang and Mirarab (2022) have introduced a weighted ASTRAL (wASTRAL) approach that weights quartets of taxa by branch support, branch length or both (hybrid).This approach has been demonstrated to consistently outperform ASTRAL, with the hybrid wASTRAL being the most accurate.We inferred support weighted (wASTRAL-support), branch length weighted (wASTRAL-length) and hybrid (wASTRAL-hybrid) wASTRAL phylogenies for 50, 75 and 95 % complete datasets.IQ-TREE gene tree sets with UFBoot values and branch lengths were used as input.Branch support was determined with localPP for wASTRAL-hybrid and wASTRAL-support species trees.
A concatenated dataset species tree was inferred with IQ-TREE, with concomitant determination of the appropriate partitioning scheme and models suitable for each partition with ModelFinder (for partioning schemes and models chosen by ModelFinder, see supplementary material).Branch support was determined with 1000 replicates each of UFBoot and Shimodaira Hasegawa-like approximate likelihood ratio test (SH-alrt; Guindon et al., 2010) on IQ-TREE.

Concordance analyses and quartet support
Rapid radiations are very likely to have instances of incomplete lineage sorting (ILS) and hence, considerable disagreement among gene trees for the basal, short branches (Whitfield and Lockhart, 2007).To better understand the potential impact of this challenge for our data/ analyses, we computed gene (gCF) and site concordance (sCF) factors with IQ-TREE on ASTRAL, wASTRAL-hybrid and IQ-TREE topologies.The gCF is the proportion of gene trees containing a particular split to the number of gene trees that could have theoretically contained the split (Minh et al., 2020a).The sCF, on the other hand, is the proportion of sites in a taxon quartet subalignment supporting a bipartition among the sites that are present for all the taxa (in the quartet) and are parsimony-informative (Minh et al., 2020a).
Another metric we employed was the percentage of quartets of taxa supporting a particular relationship and two alternative resolutions of that quartet.This metric is implemented in ASTRAL and it is the one used for estimation of local posterior probability support (Sayyari and Mirarab, 2016).

Test of hard polytomy
To investigate whether the rapid radiation of elapoid snakes had any multifurcations, we employed a test implemented in ASTRAL that is designed to assess hard polytomy in multi-locus datasets while taking incomplete lineage sorting into account.This analysis tests the null hypothesis of gene tree quartets lending equal support to alternative topologies around a branch (Sayyari and Mirarab, 2016).The test was carried out on both unfiltered and filtered gene trees and their corresponding ASTRAL species trees.

Timetree inference
For our time-calibrated the phylogeny, we used Maximum Likelihood implementation of RelTime (Tamura et al., 2012;Tamura et al., 2018) in MEGA 11 (Tamura et al., 2021).RelTime can be run on datasets of the size used here, in tractable timeframes (something that is not feasible with Bayesian methods), with accuracy equivalent to Bayesian methods.Recent versions of RelTime are capable of computing a confidence interval around the estimated divergence time, containing correct times comparable to Bayesian Highest Posterior Density intervals (Tao et al., 2020).The method also allows users to incorporate fossil calibrations as a probability density.We used our weighted ASTRAL hybrid topology (which was identical across datasets) for these analyses.We estimated timetrees with 50 %, 75 % and 95 % complete datasets.Fossil calibrations were set following Head et al. (2016).The oldest calibration was the Colubridae-Elapoidea divergence, based on Coluber cadurci from the Oligocene, which was set as a lognormal distribution with an offset of 30.9, mean 2.0 and standard deviation of 0.7.The other two calibrations, also set as lognormal densities, were Naja (Naja melanoleuca-Walterinnesia aegyptia node; offset 17.0, mean 2.0, standard deviation 0.55) and Laticauda-Suta (equivalent to Laticauda + Oxyuraninae of Head et al. [2016]; offset 10.0, mean 2.0, standard deviation 0.7), based, respectively, on Naja romani and Incongruelaps iteratus from the Miocene.GTR + G + I was set as the substitution model for the divergence time estimation.To check if the lognormal calibration density parameters chosen by us affected the estimates for the time of divergence, we ran the same RelTime analyses on 50 %, 75 % and 95 % complete datasets with calibration densities set as uniform distribution, with a lower (the fossil age) and an upper (the maximum age from lognormal density confidence interval) bound.

Micro-computed tomographic (µ-CT) scanning
In order look for diagnostic characters of clades, we utilised cranial osteological information available from databases and from µ-CT scans specifically for this study.Heads of two specimens of Micrelaps muelleri scanned for this study were from Steinhardt Museum of Natural History, Tel Aviv University.The specimens were scanned at the University of Helsinki, with a GE Phoenix Nanotom S (GE Measurement and Control Solutions).The scanning was performed with a 1 mm aluminium filter, X-ray beam set at operating voltage 80 kV and current 150 µA, exposure time of 4 X 250 ms per projections and a total of 2000 projections.Voxel size was 15 µm.Volume rendering was performed with Dragonfly (Object Research Systems Inc., Montreal, Canada).CT scans of the cranium of members of most other elapoid family/subfamily level taxa were obtained from MorphoSource (https://www.morphosource.org/) and Digimorph (https://digimorph.org/index.phtml).Many of these scans were generated in the previous projects of the authors of the present study.For those scans obtained from the databases that were not already volume rendered, volume rendering was done with Dragonfly.All the volumes were visualised with Meshlab (Cignoni et al., 2008).The scans were visually examined to look for qualitative osteological characters that can diagnose major elapoid subclades.The examined scans were not limited to the species included in the present phylogenomic study.All the specimens examined in the form of µ-CT scans, with their Morpho-Source, Digimorph and collection identifiers (if any), are listed in the supplementary material.Anatomical terminology follows that of Das et al. (2022).
S. Das et al.

Phylogeny
The 50, 75 and 95 % complete ultraconserved elements (UCE) dataset contain 4561 (5,634,032 bp in total of which 1,030,492 are parsimony informative), 4372 (5,477,173 bp; 1,005,959 bp parsimony informative) and 3600 (4,534,172; 808,398 bp parsimony informative) loci, respectively.One sample, Compsophis infralineatus, yielded too few UCE loci and was discarded during the processing of the dataset in phyluce and therefore, our datasets had 47 taxa of which three were outgroups.
The ASTRAL (both unfiltered and filtered gene tree sets-based), wASTRAL-hybrid, wASTRAL-length and wASTRAL-support multispecies coalescent species trees from 50, 75 and 95 % complete datasets were fully congruent with each other with respect to the phylogenetic relationships of the family/subfamily level elapoid subclades (Fig. 1A, B, Supplementary material Figs.1-15).LocalPP branch support was generally very high for the deeper divergences, being mostly > 0.95 to the maximum of 1 in ASTRAL species trees from unfiltered gene tree sets and always ≥ 0.95 to 1 in the remaining multispecies coalescent trees.In the multispecies coalescent phylogenies, the earliest cladogenesis event is the one between the Cyclocoridae and the rest of the Elapoidea.Elapidae was sister to a major clade containing the whole principally Afro-Malagasy radiation of elapoids (corresponding to the large, inclusive Lamprophiidae sensu Vidal et al. [2007], Pyron et al. [2011], Pyron et al. [2013] etc), and Micrelaps.Micrelaps, a taxonomically problematic genus, was recovered as sister to the Afro-Malagasy radiation.The relationships between the major subclades within the Afro-Malagasy radiation ('Lamprophiidae') were also consistently inferred.The Malagasy Pseudoxyrhophiinae was the first to split from other subclades within this radiation.The other groups are divided into two larger clades, one consisting of Psammophiinae and Atractaspidinae (including genera allocated to Aparallactinae) and another comprising of Prosymninae, Pseudaspidinae (sister to Prosymninae), and Lamprophiinae (sensu stricto).In none of these phylogenies did Atractaspis and Homoroselaps, traditionally considered to comprise the subfamily Atractaspinae (e.g., Portillo et al., 2019), have a sister taxon relationship.The three aparallactine genera, namely Aparallactus, Polemon and Xenocalamus, did not form a clade either, and Xenocalamus was recovered in a node basal to Atractaspis.Every genus for which we could sample two or more species was found to be monophyletic with high localPP in our ASTRAL and wASTRAL species trees.The sole disagreement between these phylogenies was the position of Boaedon olivaceus within the genus Boaedon but the different positions were not well supported by localPP.
Concatenation based Maximum Likelihood (ML) species trees, inferred with IQ-TREE, recovered topologies bearing an overall similarity to the multispecies coalescent ones, albeit with some interesting differences (Fig. 2, Supplementary material fig.16-18).The splits that differed from the coalescent species trees received moderate to low SHalrt and UFBoot support.Most of the other branches, including deeper ones, received very high statistical support.The 75 % and 50 % complete concatenated datasets yielded phylogenies in which Pseudaspidinae was sister to Lamprophiinae and Prosymninae but the (Lamprophiinae, Prosymninae) topology did not receive high support (Fig. 2, Supplementary material fig.17).In the 95 % dataset tree, Pseudoxyrhophiinae had a more nested position within the Afro-Malagasy clade, branching basal to the (Lamprophiinae, (Prosymninae, Pseudaspidinae)) clade, but this was not strongly supported (Supplementary material fig.18).

Gene and site concordance
The percentage of quartets supporting the split representing the common ancestor of elapoids was high, at 71 %.Internal branches below the Afro-Malgasy clade varied from 43 to 46 % quartets, with alternative topologies receiving a quartet support of 29 % or less (Supplementary material fig.19-28).Within the Afro-Malagasy radiation, the quartet support tended to be lower and the gap between the chosen topology Gene and site concordance (gCF and sCF) factors tended to be lower for the deeper branches within the Afro-Malagasy radiation (Supplementary material fig.29-36).Most relationships within this part of the phylogeny were supported by only a little over 33 % sites in the concatenated alignment, indicating the presence of conflicting signal in the dataset (Minh et al., 2020a).For relationships differing among multispecies coalescent and concatenation species trees, the former tended to receive higher gCF than the latter, even if only marginally so.

Hard polytomy
At α = 0.05, the null hypothesis of hard polytomy was rejected for every relationship in both the filtered and unfiltered gene tree set based ASTRAL phylogenies (Supplementary material fig.37-41).At α = 0.01 as well, the null hypothesis could be rejected in all but one branch in 50 and 75 % dataset phylogenies.This branch is the one connecting Homoroselaps with the rest of the atractaspidines and aparallactines and was a very short one, being 0.03 in coalescent units (CU) and > 2 million years in absolute time.There are, however, a few other branches with 0.02 -0.03 CU for which the p-value was zero and therefore, the null hypothesis of an elapoid multifurcation was rejected.

Timetree
RelTime analysis on the 50 % complete dataset, with lognormal node calibration densities, produced a timetree (Fig. 3) that placed the origin of Elapoidea at 54.6 million years (MYA hereafter), in the Ypresian age (56.0 -47.8 MYA) of the Eocene, with a confidence interval (44.9 -64.7 ).Timetrees estimated from 75 and 95 % complete datasets, with lognormal node calibration densities, yielded similar divergence times (Supplementary material fig.43-46).When node calibrations were set to be a uniform distribution with a minimum and a maximum age, estimated divergence times tended to be slightly younger (Supplementary material fig.47-52).For example, the origin of Elapoidea was 50 MYA in the uniform node calibration timetree versus 54.6 MYA in the uniform node calibration Timetree from 50 % complete dataset.However, the estimated divergence times in these timetrees were close to those in lognormal node calibration timetrees, almost always placed at the same geological stage/age and the confidence intervals around node ages in both sets of timetrees were broadly overlapping.

Cranial osteology
Our visual examination of µ-CT scan volumes and survey of literature did not reveal any diagnostic feature, or a combination thereof, for Lamprophiinae and Pseudoxyrhophiinae, two subclades containing an ecologically (and thus morphologically) diverse set of taxa.
Prosymnines are characterised by a combination of cranial characters, namely blade-like posterior maxillary teeth, premaxilla with long, posteriorly directed, long transverse processes (except Prosymna visseri), edentulous pterygoids, well developed nasal horizontal lamina almost or fully reaching frontals, a wide, short parietal and various degrees of cranial fusion (Bourgeois, 1968;Heinicke et al., 2020;present study).
Pseudaspidinae, as currently understood, did not yield any combination of diagnostic features.We examined the µ-CT scanned skull of Psammodynastes, controversially classified into this subfamily (see the review in Zaher et al., 2019), but it did not reveal any character that allows to associate this genus with any other pseudaspidine genus.
Atractaspidines and aparallactines can be diagnosed by a set of characters, that upon closer examination shows to form a transformation series (e.g., gradual loss of the anterior part of the maxilla [Das et al., 2022]), namely fossorial adaptations such well-developed nasal horizontal lamina which closely approaches frontals (except some Aparallactus spp.), premaxilla adapted variously for burrowing, presence of a pseudocoronoid process in various levels of prominence, from a faint ridge to almost a keel (often not developed in the paedomorphic mandible of many Atractaspis spp.), enlarged, grooved rear fang (except Aparallactus modestus) below orbit or prefrontal, with ≤ 8 aglyphous teeth anterior to them or a front fang (a homolog of the rear fang), maxillary ascending process present, in various levels of prominence, and braincase with fossorial adaptations (Portillo et al., 2019;Das et al., 2022;present study).
The cranium of Micrelaps muelleri (Fig. 4, Supplementary material fig.54) is very similar to skulls of snakes belonging to the atractaspidine and aparallactine and this has resulted in Micrelaps being classified with the latter (McDowell, 1987).Micrelaps muelleri has a robust premaxilla with an anteriorly concave ascending process, well developed vertical and horizontal nasal laminae, supraorbital processes of the prefrontal and the parietal bracing the frontal, a long and tubular parietal with a sagittal adductor ridge, a prominent maxillary ascending process, and a grooved rear fang and dentary pseudocoronoid process, similar to many aparallactines (Das et al., 2022.But one difference diagnoses Micrelaps from all the aparallactines and atractaspidines studied by usthe shaft of the ectopterygoid bone in Micrelaps has a prominent, posterolaterally directed protuberance (not homologous to the anterolateral lobe of the ectopterygoid, as the latter is also present) opposite to the bone's articulation to the pterygoid (Fig. 4), perhaps for the attachment of the pterygomandibularis muscle (Das and Pramanick, 2019), which is absent in aparallactines.Although the atractaspidine Homoroselaps has a small, lateral protuberance on its ectopterygoid, it is much rostrad to the ectopterygoid-pterygoid articulation.The monotypic Brachyophis is a very poorly studied genus and has never been included in a molecular phylogenetic analysis.Underwood and Kochva (1993) found Brachyophis to be the sister of Micrelaps based on soft anatomical characters.Although we could not sample Brachyophis revoili for our phylogenomic analyses, examination of the µ-CT scanned crania (Supplementary figure 53) lends support to Underwood and Kochva's hypothesis.The ectopterygoid bone of Brachyophis has the posterolateral protuberance.In this taxon, the protuberance, along with a ridge from the pterygoid, forms a Psammophiines have a 'typical', opisthoglyphous colubroid cranium, often diagnosable within Elapoidea by the combination of a very large optic foramen (especially in Psammophis and Malpolon), high intertrabecular crest on the parabasisphenoid, a grooved rear fang and lack of any fossorial adaptation.
Cyclocorid skulls have been investigated in detail by Weinell et al. (2020).Even though most of the species assigned to the family shows some degree of fossorial adaptation, no synapomorphy or a combination of characters diagnosing this family has been reported.

Phylogenomic resolution of the rapid radiation of elapoid snakes
The principal objective of the present study was to resolve deep-time phylogenetic relationships of the superfamily Elapoidea-the region of elapoid phylogeny where previous molecular phylogenetic studies have inferred highly incongruent topologies, often with poor statistical branch support.In contrast, our ASTRAL (estimated from both filtered and unfiltered gene trees), wASTRAL (branch length and support weighted and hybrid weighted) and concatenated dataset Maximum Likelihood phylogenies, from UCE datasets of different levels of completeness, yielded highly congruent topologies.The multispecies coalescent methods were especially consistent in this regard.Branch support received by deeper splits in the tree were mostly very high (≥0.95localPP and ≥ 95 UFBoot and SH-alrt).
For relationships differing in coalescent and concatenation-based phylogenies, genome and site concordance metrices in IQ-TREE and ASTRAL quartet support usually favoured the coalescent phylogeny resolution, even if only marginally so.Rapid radiations are known to be subject to extensive incomplete lineage sorting, which results in anomaly zones in the phylogeny (Degnan and Rosenberg, 2006;Linkem et al., 2016), and this is known to mislead concatenation-based analyses (Mirarab et al., 2016;Edwards et al., 2016).Therefore, it seems advisable to treat multispecies coalescent phylogenies as the best estimates for ancient, rapid radiations for a given dataset.We treat the wASTRAL topologies as our best estimate of elapoid phylogenetic relationships as wASTRAL has been shown to outperform ASTRAL (Zhang and Mirarab, 2022).
All phylogenies inferred in our study recovered Asian cyclocorids as the first-branching subclade within Elapoidea.Elapidae was inferred as the sister lineage to the Micrelaps + Afro-Malagasy radiation, also diverging early.Interestingly, species-rich phylogenies usually recover Asian coral snakes (Calliophis) as derived from the most basal node within Elapidae (Pyron et al., 2013;Zaher et al., 2019).An Asian origin, therefore, cannot seem to be ruled out for Elapoidea.However, most major elapoid subclades are either endemic to the Afro-Malagasy region or at least, reach their maximum species diversity there.Therefore, the rapid Eocene cladogenesis of Elapoidea, especially that of the Afro-Malagasy radiation, might have occurred in this region.Many other reptile groups, such as chameleons (Tolley et al., 2013), and southern African testudinids (Hofmeyr et al., 2017), diversified during the Eocene in Africa.The Palaeogene presented a window of ecological opportunity for snakes and the early Cenozoic was indeed a period of diversification of many snake clades (Grundler and Rabosky, 2021).Relatively less competition owing to K/Pg mass extinction event may explain the Eocene diversification of elapoids in Africa and this hypothesis may be tested in future studies on the probable factors promoting the rapid cladogenesis in this superfamily.
Phylogenetic relations among genera within a particular elapoid family or subfamily, except Atractaspididae (or Atractaspinae and Aparallactinae), largely corroborated the findings of some previous Sanger-sequencing based phylogenetic studies and/or phylogenomic works with smaller dataset (Kelly et al., 2008(Kelly et al., , 2009;;Pyron et al., 2013;Ruane et al., 2015;Broadley et al., 2018;Burbrink et al., 2019;Heinicke et al., 2020) but there were some disagreements too.For example, all the phylogenies recovered Bothrophthalmus as splitting early from other genera belonging to Lamprophiinae with strong support.This differs from some previous phylogenies (Kelly et al., 2011;Pyron et al., 2013;Zaher et al., 2019 etc).
A surprising disagreement between our work and a large body of previous molecular studies (e.g., Kelly et al., 2009;Pyron et al., 2011;Pyron et al., 2013;Zheng and Wiens, 2016;Portillo et al., 2018Portillo et al., , 2019;;Zaher et al., 2019) was the non-monophyly of Atractaspidinae and Aparallactinae with respect to each other.However, these studies, despite having dense taxon sampling, have almost always used the same markers (usually < 5 loci per species) for taxa assigned to these two subfamilies, which likely explains the similar findings from those studies.It should be mentioned here that the base of the atractaspidinesaparallactine phylogeny is also suspected to constitute a rapid radiation (Portillo et al., 2018(Portillo et al., , 2019;;present study).This part of the phylogeny has very short branches, genomic and site discordance and for the branch uniting Homoroselaps with the rest of the genera, the null hypothesis of a polytomy could not be rejected at α = 0.01.This at least shows that we cannot place high confidence in the reciprocal monophyly of Atractaspidinae and Aparallactinae.It is noteworthy to mention that morphological phylogenies never recovered reciprocal monophyly of these two groups either (Underwood and Kochva, 1993;Das et al., 2022).
Recent phylogenomic studies (e.g., Tilic et al., 2020) have demonstrated that some evolutionary scenarios can only be resolved with a phylogenomic dataset consisting of thousands of genes.Sampling of loci for elapoids has not kept up with the pace of increasing taxon sampling over the last decade (Kelly et al., 2009;Pyron et al., 2011Pyron et al., , 2013;;Figueroa et al., 2016;Zheng and Wiens, 2016;Zaher et al., 2019).These studies, while contributing greatly to alpha-and beta-level relationships, produced highly conflicting, poorly supported topologies when it came to the basal branching order within Elapoidea.This, and the example of Atractaspidinae and Aparallactinae discussed above, seem to indicate that dense taxon sampling alone may not suffice to resolve rapid cladogenesis scenarios.Interestingly, only the > 300 loci phylogenomic study of Pyron et al. (2014), albeit with lesser taxon sampling, has recovered deeper relations within Elapoidea somewhat close to our Maximum Likelihood estimates.While dense sampling of both taxa and loci is desirable, for various practical constrains (viz.the rarity of the samples, logistical and computational limitations), this ideal scenario is often not realised.This necessitates a trade-off between the number of species (or specimens per species) sampled and the number of loci.Systematics has seen multiple debates on the merits of dense sampling of data versus taxa (Rosenberg andKumar, 2001, 2003;Zwickl and Hillis, 2002;Heath et al., 2008).A future study addressing this exact debate in the specific context of rapid radiations seems warranted.

Higher-level systematics of Elapoidea
Disagreement over whether a Pan Afro-Malagasy Lamprophiidae is monophyletic with respect to Elapidae has resulted in two different proposed classifications of these subclades.Kelly et al. (2009) and Zaher et al. (2019) did not recover a monophyletic Afro-Malagasy radiation and classified every subclade at the family level (namely, Atractaspididae [Atractaspidinae and Aparallactinae], Lamprophiidae, Prosymnidae, Pseudaspididae, Psammophiidae and Pseudoxyrhophiidae).Other authors (Pyron et al., 2013;Zheng and Wiens, 2016) recovered a monophyletic Afro-Malagasy clade and treated these subclades as subfamilies under an inclusive Lamprophiidae.Our phylogenies recover the pan Afro-Malagasy Lamprophiidae as monophyletic with high support.We therefore treat the subclades within the Afro-Malagasy radiation as subfamilies.Given the lack of reciprocal monophyly between Aparallactinae and Atractaspidinae, we treat the whole group as Atractaspidinae, without any further division.We did not find any qualitative cranial characters diagnosing Lamprophiidae.However, as noted in the Results section, Prosymninae, Atractaspidinae, and to some extent Psammophiinae, are anatomically diagnosable.
The status of Elapidae as a family has been uncontroversial and we treat it as the same.Dentition serves as a reliable osteological diagnostic character.
Cyclocoridae was initially erected as a subfamily by Weinell and Brown (2018) but has subsequently been treated at a family rank (Weinell et al., 2020).As this clade was recovered as sister to the rest of the Elapoidea, we also treat it as a family.Anatomical synapomorphies for this family are not known.
We could not sample Buhoma and Psammodynastes, two genera with uncertain phylogenetic affinities, for UCE sequencing.We treat Psammodynastes within Pseudaspidinae, along with Pseudaspis and Pythonodipsas, following Pyron et al. (2013), Zheng and Wiens (2016) and Zaher et al. (2019).We did not find any osteological synapomorphy supporting the placement Psammodynastes within Pseudaspidinae.Therefore, we reiterate the caution of Zaher et al. (2019) that a high degree of confidence cannot be placed in the contents of Pseudaspidinae, in its current form, except the placement of Pseudaspis and Pythonodipsas.Pyron et al. (2013) recovered Buhoma within Pseudaspidinae and classified it therein, in disagreement with Kelly et al. (2009), Pyron et al. (2011), Figueroa et al. (2016) and Zaher et al. (2019).Hence, we treat Buhoma as incertae sedis within Elapoidea pending a phylogenomic analysis.Pyron et al. (2013) opined that Micrelaps warrants its own family because it could not be assigned to any elapoid family or subfamily.Our phylogenies showed Micrelaps is not nested within any elapoid subclade.Instead, it was recovered consistently as sister to Lamprophiidae with high support in all phylogenetic analyses.Therefore, we establish a new family to accommodate Micrelaps.We also assign Brachyophis revoili to this family based on shared features in visceral anatomy (Underwood and Kochva, 1993) and cranial osteology.

Taxonomic revision
Micrelapidae new family.Type genus: Micrelaps Boettger, 1880.Type species: Micrelaps muelleri Boettger, 1880.Etymology: Boettger (Böttger) did not give the etymology for the generic nomen but was almost certainly from the Latin adjective micro-, derived from the Greek mikros (small), and elaps, the Latinised form of the Greek noun éllops or élaps (literally sea-fish or serpent, but here in reference to the snake genus Elaps, now a synonym of Homoroselaps).Micrelapidae fam.nov. is derived from Micrelaps by the taking the stem elap-of the root word of the nomen.
Content: Micrelaps muelleri Boettger, 1880, Micrelaps bicoloratus Sternfeld, 1908, Micrelaps vaillanti Mocquard, 1888, Brachyophis revoili Mocquard, 1888.Diagnosis and definition: In the crania of Micrelaps and Brachyophis we examined the ectopterygoid was laterally and medially expanded at the point of contact with the pterygoid, with this expansion not being contiguous with the ectopterygoid anterolateral and anteromedial lobes (Fig. 4, Supplementary material fig.53, 54).The lateral expansion is a posterolaterally and somewhat ventrally directed, very prominent protuberance continuous with a ridge on the ventral surface of the pterygoid.This character state was not present in any other cranium we examined and is very likely a synapomorphy of the family.
Other common cranial characters include a premaxilla adapted for a fossorial lifestyle, premaxillary transverse processes closely approaching the maxilla, a short maxilla with ascending processes abutting the prefrontal, well-developed, grooved fangs below the orbit, preceded by a diastema and 2 -3 teeth, an ectopterygoid deeply forked into anterolateral and anteromedial lobes that articulate with maxillary ectopterygoid processes leaving a foramen in the middle, prefrontal and parietal supraorbital processes laterally bordering the frontal and almost meeting each other, a tendency towards fusion of cranial bones (especially because the supratemporal is absent, very likely fused to the quadrate in Brachyophis and to posterior chondrocranial elements in Micrelaps), and a short quadrate.Brachyophis, however, differs from the type genus in possessing a postorbital (versus postorbital absent in S. Das et al.

Fig. 2 .
Fig. 2. Maximum Likelihood species tree from the concatenated 50 % complete dataset consisting of 4561 loci.Values on the branch indicate Shimodaira Hasegawalike approximate likelihood ratio test and ultrafast bootstrap.Abbreviations as in Fig. 1.

Fig. 3 .
Fig. 3. Time calibrated phylogeny (50 % complete dataset) of elapoid snakes, estimated with the Maximum Likelihood implementation of the RelTime method (with lognormal node calibration densities).Values on the branches indicate the estimated divergence times.The blue bar represents the 95 % confidence intervals around the estimated divergence times.

Fig. 4 .
Fig. 4. A. Dorsal, B, lateral and C. ventral view of the skull of Micrelaps muelleri (SMNH.R 17777).D. Ectopterygoid of the same specimen in ventral view.E. Ventral and F. lateral views of the palatomaxillary arch.