Phylogenetic relationships of xenodermid snakes (Squamata: Serpentes: Xenodermidae), with the description of a new genus

Xenodermidae is a generally poorly known lineage of caenophidian snakes found in South, East and Southeast Asia. We report molecular phylogenetic analyses for a multilocus data set comprising all five currently recognised genera and including new mitochondrial and nuclear gene sequence data for the recently described Stoliczkia vanhnuailianai. Our phylogenetic results provide very strong support for the non-monophyly of Stoliczkia, as presently constituted, with S. borneensis being more closely related to Xenodermus than to the Northeast Indian S. vanhnuailianai. Based on phylogenetic relationships and morphological distinctiveness, we transfer Stoliczkia borneensis to a new monotypic genus endemic to Borneo, Paraxenodermus gen. nov. We also present new morphological data for P. borneensis.


Introduction
with a record from southern Cambodia (Fig. 1). The two species of the recently described Parafimbrios are thus far recorded only from northern Vietnam, Laos and Thailand (Fig. 1). The monotypic Xenodermus may be the most widespread xenodermid species, occurring in southernmost Myanmar, Thailand, peninsular Malaysia, Borneo, Sumatra and Java (Fig.1). The genus Stoliczkia currently includes three poorly known species with a particularly disjunct distribution, two occurring in Northeast India (S. khasiensis Jerdon, 1870 andS. vanhnuailianai Lalronunga, Lalhmangaiha, Zosangliana, Lalhmingliani, Gower, Das & and one in northern and western Borneo (S. borneensis) (Das 2021;Stuebing et al. 2014) (Fig.1). Previously, molecular data were available only for S. borneensis (Vidal and Hedges 2002), and few morphological data for the <10 reported specimens of Stoliczkia (sensu lato) had been published (Lalronunga et al. 2021).
In this paper, we report the first molecular data for Northeast Indian Stoliczkia and new morphological data for S. borneensis. We test the monophyly of Stoliczkia, and describe a new genus for the Bornean species.

Phylogeny
We aligned the new sequences for Stoliczkia vanhnuailianai with eight other xenodermids, and an outgroup, the non-xenodermid caenophidian Acrochordus granulatus. We checked for stop codons in unexpected regions by translating nucleotide alignments to amino acids for pro-tein-coding genes (cytb, co1, cmos, nt3) using MEGA 7 (Kumar et al. 2016). We aligned sequences using Clust-alW (Thompson et al.1994) in MEGA 7 (Kumar et al. 2016) with default settings (alignments online from the Natural History Museum data portal: https://doi.org/ 10.5519/gbzyjuli).
First, we built individual gene trees using Maximum Likelihood (ML). Based on availability of sequence data, we selected one species per xenodermid genus (though included both S. borneensis and S. vanhnuailianai for Stoliczkia) and the outgroup. We then aligned and concatenated the six gene sequences into a single dataset (3122 basepairs in length) with ten tips, including the outgroup (Table 1).
We used PartitionFinder2 (Lanfear et al. 2017) to identify the best-fit partition scheme for the concatenated dataset and the best-fit model of sequence evolution for each partition as determined by the Bayesian Information Criterion (BIC), using the default greedy algorithm linked to branch lengths (Lanfear et al. 2012). The best-fit scheme for the concatenated dataset comprises six partitions, by gene and by codon position (Table 2). We performed Maximum Likelihood (ML; Felsenstein 1981) phylogenetic analyses with RAxML GUI Ver. 2.0 (Edler et al. 2021), using the GTRGAMMA model of sequence evolution, which is recommended over GTR+G+I because the 25 rate categories account for potentially invariant sites (Stamatakis 2006). For Bayesian (BI) phylogenetic analyses we used MrBayes 3.2.6 (Ronquist et al. 2012) via the XSEDE portal CIPRES Science Gateway v3.3 (Miller et al. 2010), with default prior settings and with all six partitions assigned their best-fit model as determined by PartitionFinder (Table 2). We set up two separate runs with four Markov chains each, initiated from random trees and allowed to run for one million generations, sampling every 100 generations and discarding the first 25% of trees as "burn-in". We terminated the analyses when the standard deviation of split frequencies was less than 0.005, and then constructed majority rule consensus trees. We checked for effective sample size (ESS) values using Tracer 1.7 (Rambaut et al. 2014), all parameter values had ESS values >200. We quantified support for internal branches in ML and BI trees using bootstrap (500 replicates) and posterior probability, respectively. We assessed levels of support for relationships incompatible with optimal trees by inspecting bipartition tables of ML bootstrap or BI posterior probability trees using PAUP* 4.0a 169 (X86) (Swofford 2003). We rooted the trees with Acrochordus granulatus because it is a non-xenodermid caenophidian snake (Figueroa et al. 2016;Deepak et al. 2018;Zaher et al. 2019).

Molecular dating
We aligned a larger dataset with 68 tips including two scolecophidians (Gerrhopilus mirus and Liotyphlops albi rostris) and representatives of all subfamilies of Alethinophidia, including nine xenodermids (sampling all five currently recognised genera). We aligned this dataset Table 1. GenBank accession numbers and voucher numbers for the sequences used in this study. Sequences used in the ML and BI concatenated phylogeny are indicated with an asterisk. Accession codes for sequences newly generated in this study are in bold text. -separately using the same methods outlined above (alignments available at: https://doi.org/10.5519/gbzyjuli). We applied seven fossil calibrations (Table 3), largely those recommended by Head (2015) and Head et al. (2016) (Miller et al. 2010) under a Yule tree process. We assigned a relaxed log-normal clock for each partition of the concatenated BEAST analysis. We set up two independent runs, each employing the Markov Chain Monte Carlo (MCMC) for 100,000,000 generations, sampling every 5,000 trees. We obtained effective sample size (ESS) values using Tracer 1.7 (Rambaut et al. 2014). ESS values were below 100 for the priors and posteriors employing the best-fit model identified using PartitionFinder. We also repeated the analysis implementing the less-complex HKY model for the partitions but otherwise using the same settings. However, in this second analysis, we recovered ESS values above 200 for all the priors and posteriors for the two independent runs.

Morphology
We provide here morphological and meristic data for two specimens of Stoliczkia borneensis (BMNH 1946. 1.15.58 and UNIMAS 8002) and additional published information on unspecified specimens from Stuebing et al. (2014). Total length, snout-vent length and tail length were measured with thread and a ruler to the nearest 1 mm. Other dimensions were recorded with dial calipers, to the nearest 0.1 mm. Bilateral scale counts separated by a comma are reported in left, right order. Ventral scales were counted following Dowling (1951). Length and width of head scales were measured at the longest and the widest points of the respective scale(s). Eye diameter was measured horizontally.

Phylogeny
The single-gene ML trees are shown in Fig. 2. Depending on taxon sampling (limited by availability of sequence data), generally S. borneensis and S. vanhnuailia nai show close affinities with Xenodermus and with Achalinus, respectively. Although ML bootstrap support for many relationships are not strong (<90), support for Stoliczkia monophyly in the four gene trees for which both species were sampled is negligible, being only 25 for 16S and 0-0.2 for 12S, cmos and nt3. The ML and BI trees derived from the concatenated dataset agree in the set of relationships depicted (Fig. 1), with generally moderate (70-90 ML bootstrap; 0.80-0.90 BI posterior probability) to high support (>90 ML; >0.95 BI). Importantly, there is zero bootstrap or posterior probability support for Stoliczkia monophyly in these latter trees. Instead, the best-supported relationships that are incompatible with this optimal set of relationships for Stoliczkia spp. are for Xenodermus javanicus being more closely related to Fimbrios and Parafimbrios (ML bootstrap = 20; BI posterior probability = 0) and for S. borneensis being more closely related to Fimbrios and Parafimbrios (ML bootstrap = 5; BI posterior probability = 9). Thus, we conclude that the available DNA sequence data provide good to strong support for S. borneensis being more closely related to Xenodermus than to S. vanhnuailianai, and for S. vanhnuailianai being more closely related to Achalinus than to S. borneensis, and very strong support for non-monophyly of Stoliczkia.

Morphology
Previously, extensive data were available for only a single vouchered specimen (the holotype, BMNH 1946.1.15.58) of Stoliczkia borneensis (Lalronunga et al. 2021). Data for an additional specimen (UNIMAS 8002) are presented in Table 4. This specimen agrees with data presented by Lalronunga et al. (2021) corroborating that S. vanhnuailianai resembles the type species of the genus, S. khasiensis much more closely than either does S. borneensis. Notable differences between the Bornean species and the two Northeast Indian species include presence of 4-6 small scales between the frontals and prefrontals in S. borneensis which are absent in the Northeast Indian species; supralabials not contacting the eye in S. borneensis versus contacting the eye S. vanhnuailianai and S. khasiensis; 10 or 11 supralabials versus 8 or 9 supralabials. Although S. borneensis is seemingly most closely related to Xenodermus (Fig. 2), the two taxa differ markedly in external morphology-for example, X. javanicus lacks large scales on the head other than at the snout tip whereas S. borneensis additionally has large parietal and frontal shields. Xenodermus javanicus and S. borneensis share a derived condition of having more small, irregular head scales than are present in other xenodermids.  Anterior temporals 0,0 0,0
Distribution. This genus is restricted to Northeast India (Fig. 1). Stoliczkia khasiensis is thus far known only from Khasi hills, Meghalaya state, India and the recently described Stoliczkia vanhnuailianai is known only from Mizoram state, India.
Etymology. The genus is named after the Moravian-born Ferdinand Stoliczka (1838-1874). A geologist-natural historian, he was appointed as a palaeontologist with the Geological Survey of India in 1863. Stoliczka collected vertebrates and molluscs from northern India, Andaman and Nicobar Islands, Myanmar and the Malay Peninsula. He served as the official Naturalist with the Second Mission to Yarkand, in central Asia. A biography and a list of published works and reports by Stoliczka can be found in Kolmaš (1982).

Paraxenodermus borneensis
(3) posterior one-third of the head and posterior temporal region covered with small scales like those of the anterior of the body, (4) numerous small scales between parietal and supralabial shields immediately behind eye, (5) a row of 4-6 small scales between the frontal and prefrontal shields, (6) 10-11 supralabials, (7) nostril in a large concave nasal, (8) body slender and somewhat laterally compressed, (9) ventrals large, and (10) (Boulenger 1899) and the contiguous Crocker Range, both in Sabah, in the northeastern part of Borneo (Das 2006a), as well as in the isolated Gunung Murud (Das 2006b), in Sarawak State. https://www.inaturalist.org/observations?taxon_id=28573). Information is not available for the holotype, but all other reported individuals were found late at night, moving slowly on rocky banks of streams at elevations of 950-2,100 m above sea level (Das 2006a).
Variation. The two examined specimens of Paraxenodermus borneensis, the holotype BMNH 1946.1.15.58 and UNIMAS 8002, differ slightly in the number of small scales lying between the frontal and prefrontals, being six and four, respectively. We counted six small scales in this position in images of a live individual on the internet (https://www.inaturalist.org/taxa/28573). UNIMAS 8002 also differs from BMNH 1946.1.15.58 in having a two more ventrals (208 versus 206) and five additional subcaudals (128 versus 123), and in being smaller (713 mm versus 789 mm total length).
Etymology. The generic name Paraxenodermus is composed of the modern Latin generic name Xenodermus and the Latin adjective par (paris), meaning, among other possibilities, "similar to".

Discussion
Taken at face value, our phylogenetic results and the distribution of xenodermid genera (Fig. 1) indicate that there are two main radiations within Xenodermidae; one in Northeast India, northern mainland Southeast Asia and Japan (Stoliczkia + Achalinus sensu stricto) and one in eastern mainland Indochina and southeast Sundaland 1) some supralabials in contact with eye in Stoliczkia, separated by circumorbital scales in Paraxenodermus; 2) fewer supra-and infralabials in Stoliczkia than in Paraxenodermus; 3) single prefrontal in Stoliczkia versus 2-3 in Paraxenodermus, 4) fewer scales between parietal and supralabials immediately behind eye in Stoliczkia than in Paraxenodermus, and 5) small row of scales between frontal and prefrontals absent in Stoliczkia, present in Paraxenodermus. Note small scales behind the temporals are indicative rather than precisely accurate. Pale grey coloured areas are bare skin exposed between scales. Illustrations by V. Deepak and Surya Narayanan. Scale bars = 10 mm.   Millions of years ago (Fimbrios, Parafimbrios, Paraxenodermus, Xenodermus). The most parsimonious interpretation is that the most recent common ancestor of these two main xenodermid radiations occurred in mainland Indochina, suggested by our dating analyses to be approximately 66.7-44.6 Ma (Fig. 6). However, this would be better tested in future by undertaking probabilistic biogeographic analyses of a more broadly taxonomically sampled tree. Establishment of a new genus for S. borneensis and a new understanding of phylogenetic relationships removes the exceptional geographic disjunction presented by the previous concept of Stoliczkia. These results also strengthen evidence for endemic radiations within both Borneo (e.g. Blackburn et al. 2010;Wood et al. 2012;Hertwig et al. 2013;Fritz et al. 2014) and Northeast India (e.g., Pawar et al. 2007;Kamei et al. 2012).