Two new syntopic species of wolf snakes (genus Lycodon H. Boie in Fitzinger, 1826) from an imperiled ecosystem in the Song Giang River Valley of southern Vietnam (Squamata: Colubridae)

An integrative taxonomic analysis of species in the colubrid genus Lycodon Fitzinger, 1826 recovered two new syntopic species of the L. rufozonatus complex from the imperiled Song Giang River valley in Khan Hoa Province, of Southern Vietnam. Although L. truongi sp. nov. and L. anakradaya sp. nov. are syntopic, they are not particularly closely related and can be differentiated from each other and all other species in the L. rufozonatus complex on the basis of meristics, morphometrics, color pattern, and uncorrected pairwise genetic distance based on the mitochondrial gene cytochrome b. The discovery of these two new range-restricted species and a previously described range-restricted gekkonid in the genus Cyrtodactylus Gray, 1828 from the same valley, underscores the necessity of continued field work in the Song Giang River valley so as to catalog the unrealized herpetological diversity in this area and establish research-based conservation programs.


Introduction
The colubrid snake genus Lycodon Fitzinger, 1826 comprises a large radiation of at least 71 non-venomous, terrestrial to arboreal, nocturnal species that occupy a range of forested habitats extending from the Caspian Sea to Sulawesi (Wallach et al. 2014;Uetz et al. 2021). Having such a wide distribution that crosses so many well-established biogeographic borders, it is not surprising that this genus manifests a broad array of species with trenchant Vertebrate Zoology 72, 2022, 371-384 | DOI 10.3897/vz.72.e82201 Copyright Anh The Nguyen et al. This is an open access article distributed under the terms of the Creative Commons Attribution License (CC BY 4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited. morphological differences and species groups bearing deep phylogenetic substructuring (Siler et al. 2013;Wang et al. 2020).
Recent field work in the imperiled riparian habitats of the Song Giang River Valley in Khanh Hoa Province of Southern Vietnam (Fig. 1) has already begun to report the discovery of a new range-restricted endemic species (Nguyen et al. 2021) and here we report the discovery of two new syntopic species of wolf snakes, Lycodon. All three specimens have dorsoventrally flattened heads; vertically elliptical pupils; large nostrils; strongly arched upper maxillary bones whose anterior sections angle inwards; curved, enlarged, anterior maxillary teeth followed by a diastema; 17 or 18 anterior, 17 midbody, and 15 posterior rows of dorsal scales bearing weakly keeled vertebral rows; and rounded, weakly notched, ventral scales (Lanza 1999)-characters that place them in the genus Lycodon. Furthermore, sequences of the new specimens were added to a cytochrome b (cyt b) phylogenetic data set from Wang et al. (2021) and were recovered as two different species nested within a larger clade of Indochinese species. Additionally, these specimens differ from each other and all other species of Lycodon in the Indochinese clade by having a unique suite of morphological and color pattern characters. We therefore consider them new species and describe them below.

Species delimitation and concept
The general lineage concept (GLC: de Queiroz 2007) adopted herein proposes that a species constitutes a population of organisms evolving independently from other such populations owing to a lack of, or limited gene flow. By "independently," it is meant that new mutations arising in one species cannot spread readily into another species (Barraclough et al. 2003;de Queiroz 2007). Under the GLC implemented herein, molecular phylogenies recovered monophyletic mitochondrial lineages of individuals (populations) used to develop initial species-level hypotheses-the grouping stage of Hillis (2019). Discrete color pattern data and morphological data were then used to search for unique suites of characters consistent with the tree-designated species-level hypotheses-the construction of boundaries representing the hypothesis-testing step of Hillis (2019)-thus providing independent diagnoses to complement the molecular analyses. It is important to note, that delimiting species (phylogeny) and diagnosing species (taxonomy) are independent but overlapping operations that should not be conflated (Frost and Hillis 1990;Frost and Kluge 1994;Hillis 2019).

Morphological analysis
Morphological and color pattern data were taken from three specimens from Song Giang River Valley (12.37079° N, 108.83643° E; at elevation 500 m a.s.l.), Khanh Trung Commune, Khanh Vinh District, Khanh Hoa Province, Vietnam and from literature on Lycodon banksi Luu, Bonkowski, Nguyen, Le, Calame and Ziegler, 2018;L. cathaya Wang, Qi, Lyu, Zeng, and Wang, 2020;L. chapaensis (Angel and Bourret, 1933); L. flavozonatus (Pope, 1928); L. futsingensis (Pope, 1923);L. meridionalis (Bourret, 1935);L. rufozonatus Cantor, 1842;and L. septentrionalis (Günther, 1875) of the L. rufozonatus complex (Luu et al. 2018(Luu et al. , 2019Jansen et al. 2019;Wang et al. 2020Wang et al. , 2021. Scale counts and scale nomenclature follow Wang et al. (2021). Given the small sample sizes (n=1 and n=2), standard statistical analyses would not be viable. We are well-aware that the acquisition of new material could render some of these diagnostic characters invalid just as it may validate others and/or recover new characters not listed. Nonetheless, as mentioned above, species diagnoses (i.e. morphological differences) and species delimitation (i.e. phylogenetic relationships) are different operations with different endpoints that stand independent of one another and should not be conflated.
All body measurements were made to the nearest millimeter. Morphometric data include snout-vent length (SVL), tail length (TaL), total length (ToL), head length, width, and height (HL, HW, and HH, respectively), eye diameter (ED), snout length (SnL), eye to narial distance (EN), and internarial distance (IND). Meristic data include maxillary teeth (MT); supralabial and infralabial scales (SL and IL, respectively); of SL contacting the eye (SL-E), loreals (LoR); LoR contacting the eye (LoR-E); preoculars (PrO); postoculars (PtO); anterior temporals (aTMP); posterior temporals (pTMP); dorsal scale rows one head length posterior to the head, at midbody, and one head length anterior to the vent presented in that order; ventral scales (VEN); subcaudal scales (SC); light-colored body bands (BB); and light-colored tail bands (TB). Discrete characters evaluated were body scale texture, cloacal plate divided or single, adult head color pattern, presence or absence of a wide nuchal-occipital collar in adults and/or juveniles, dorsal ground color, body bands thin (1-3 scales rows) or wide (> three scales rows), color of body bands, ventrolateral body pattern, and ventral pattern. The institutional acronym SIEZC refers to the Zoological collection of the Southern Institute of Ecology in Hochiminh City, Vietnam.
Genomic DNA was isolated from muscle tissue stored in 95% ethanol following Nguyen et al., 2021. A fragment of the 5'-end of michondrial Cyt b gene was amplified using a double-stranded Polymerase Chain Reaction (PCR) under the following conditions: 1.0 µl genomic DNA (~ 35 µg of DNA), 0.75 µl forward primer (10μM) cytbL14910 5'-GAC CTG TGA TMT GAA AAC CAY CGT TGT-3' or CytbL14919 5'-AAC CAC CGT TGT TAT TCA ACT-3', 0.75 µl reverse primer (10μM) Cyt-bH16064 5'-CTT TGG TTT ACA AGA ACA ATG CTT TA-3 (Burbrink et al. 2000), 12.5 µl of Master Mix 2x (CWBIO, China), and 10 μl ultra-pure H2O. PCR reactions were completed using a Bio-Rad T100™ gradient thermocycler with the following reaction conditions: initial denaturation at 94°C for 5 mins, second denaturation at 94°C for 1 min, annealing at 50°C for min followed by an extension cycle at 72°C for 1 min per cycle for 35 cycles and the final extension step 72 °C for 10 mins. PCR products were electrophoresed on 2.0% agarose gel, and then target band was cut for purification using the GeneJET Gel Extraction Kit (Thermo Fisher Scientific Co., USA). The purified DNA fragments were sent to National Key Laboratory of Gene Technology (Institute of Biotechnology, VAST, Hanoi, Vietnam) for sequencing using the forward primers of the amplification step. The sequences were viewed by DNA Baser v.5, then uncertain nucleotides at both ends were removed and exported under fasta format. Sequences were aligned by MAFFT version online and checked by eye using Bioedit v.7.0.5.2 to ensure the correct amino acid reading frame.

Phylogenetic analysis
Maximum likelihood (ML) and Bayesian inference (BI) were used to estimate the phylogenetic relationships among the sampled species in our sequence alignment. An ML phylogeny was estimated using the IQ-TREE webserver (Nguyen et al. 2015;Trifinopoulos et al. 2016) preceded by the selection of substitution models using the Bayesian Information Criterion (BIC) in Mod-elFinder (Kalyaanamoorthy et al. 2017), which supported TPM2+F+I+G4 as the best fit model for codon position one, TIM3+F+I+G4 for position 2, and TIM2+F+I+G4 for position 3. One-thousand bootstrap pseudoreplicates via the ultrafast bootstrap (UFB; Hoang et al. 2018) approximation algorithm were employed and nodes having UFB values of 95 and above were considered highly supported (Minh et al. 2013). A Bayesian inference (BI) phylogenetic analysis was carried out in MrBayes 3.2.3. (Ronquist et al. 2012) on XSEDE using the CIPRES Science Gateway (Cyberinfrastructure for Phylogenetic Research; Miller et al. 2010) employing default priors and models of evolution that most closely approximated those selected by the BIC and used in the ML analysis. Two independent Markov chain Monte Carlo (MCMC) analyses were performed, each with four chains, three hot and one cold. The MCMC simulation ran for 100 million generations, was sampled every 10,000 generations, and the first 10% of each run were discarded as burn-in. Convergence and stationarity of all parameters from both runs were checked in Tracer v1.6 (Rambaut et al. 2014) to ensure effective sample sizes (ESS) were well-above 200. Post-burn-in sampled trees from both runs were combined and a 50% majority-rule consensus tree was constructed. Nodes with Bayesian posterior probabilities (BPP) of 0.95 and above were considered highly supported (Huelsenbeck et al. 2001;Wilcox et al. 2002). After removing the outgroups and other species not part of the Lycodon rufozonatus complex (Fig. 2), MEGA7 (Kumar et al. 2016) was used to calculate uncorrected pairwise sequence divergences among and within species using the complete deletion option which removes missing data and gaps.

Coloration in life
Etymology. The specific epithet "truongi" is a patronym honoring Professor Dr. Quang Truong Nguyen for his long-standing extensive contributions to the herpetology of Vietnam and his broad international collaborations.
Natural history (Fig. 7). The Song Giang River and its riparian habitat course through a forested river valley formed by east-west tending mountains in the northwestern portion of Khanh Hoa Province in the vicinity of the Song Giang Hydropower Station. The surrounding mountains form the northeastern slopes of Langbian Plateau and are covered with polydominant montane evergreen tropical forest that are dissected by a rich network of small streams and rivulets that feed into the Song Giang River. The specimen was found while foraging on the ground, in a relatively flat forest floor habitat surrounded by small rivulets and was taking refuge in the hollow of a tree. Table 2). Lycodon truongi sp. nov. could not be separated statistically from other species of the L. rufozonatus complex because of its small sample size (n=1). Therefore, the comparisons below are based on discrete differences some of which are color pattern characters tradionally used to separate species (Jansen et al. Luu et al. (2018Luu et al. ( , 2019, Wang et al. (2020Wang et al. ( , 2021, and Song (2021). Lycodon truongi sp. nov. differs from L. chapaensis and L. septentrionalis in the third, fourth, and fifth supralabials contacting the eye as opposed to the fourth and fifth in the former and second or third in the latter. Lycodon truongi sp. nov. differs from L. cathaya, L. chapaensis, and L. anakradaya sp. nov. by having 14 as opposed to 10-12 maxillary teeth and it differs from L. septentrionalis which has eight. Lycodon truongi sp. nov. differs from all other species of the L. rufozonatus complex by having 19 light-colored dorsal bands as opposed to 12 (L. anakradaya sp. nov., see below) or 19-115 collectively in the remaining species. Lycodon truongi sp. nov. is further separated from L. anakradaya sp. nov. by having two as opposed to three posterior temporals. The    Table 3. Meristic, morphometric, and discrete color pattern and scale morphology characters of the type specimens of Lycodon truongi sp. nov. and L. anakradaya sp. nov. All measurements are in millimeters. pil; rostral triangular, hardly visible from above; nasal vertically divided by a furrow along posterior margin of nostril; two square internasals in wide, medial contact, and in contact with two large, square prefrontals posteriorly; followed by a single, azygous, subpentagonal frontal, nearly as wide as long; two large, elongate parietals, contacted laterally by upper anterior and posterior temporals and a larger paraparietal; 1/1 wide, elongate supraocular; 1/1 small preocular, contacting third supralabial; 2/2 postoculars, upper slightly larger than lower; 1/1 elongate loreal not contacting eye, in contact with second and third supralabials ventrally, the prefrontal dorsally, preocular posteriorly, the nasal anteriorly; 8/8 supralabials, second to sixth higher than wide; first and second supralabials in contact with nasal; third, fourth, and fifth, supralabials contacting eye; sixth and seventh supralabials largest; two elongate anterior temporals, lower temporal largest; three square posterior temporals of similar size; 10/10 infralabials; first pair of infralabials contact medially forming a deep, medial groove; first five infralabials in contact with first pair of chinshields; similarly sized anterior and posterior pair of elongate chinshields, bearing deep, medial grooves contiguous with groove separating first pair of infralabials; 12 maxillary teeth. Body elongate, somewhat laterally compressed; SVL 790 mm; TaL 190 mm; ToL 980 mm. 225 ventrals, 87 paired, subcaudals; cloacal shield divided; dorsal scales in 18-17-15 rows with three keeled vertebral rows; vertebral scale row not enlarged; no apical pits. (Figs 4, 5). Dorsal ground color of head, body and tail dark-brown to nearly uniformly black; body bearing 13 wide orangish body bands extending to ventral scales, seven wide orangish caudal bands incompletely encircling tail; anterior half of venter orangish with dark-brown mottling grading into a uniform darkbrown venter and subcaudal region. Table 3). The paratype SIEZC 20248 resembles the holotype in overall color patter except the dorsal and caudal bands are white and the head plates are edged in white. Presumably the 12 dorsal and five caudal bands become orangish and the head will be-come unicolor dark-brown/black with increasing SVL. The paratype lacks keeled vertebral scales which also may become keeled with increasing SVL. Meristic differences are listed in Table 3.

Variation (Figs 4-6;
Etymology. The specific epithet "anakradaya" is given in a reference to the Ede people ("Anak Radaya" in Ede language), an Austronesian ethnic group living in the upland forested areas of southern Vietnam including those that surround the Song Giang River valley.
Natural history (Fig. 7). The holotype was observed between 18:30-22:30 h swimming across a small, clear stream with a sandy bottom in an effort to escape detection. Once crossing the stream, it climbed up a steep river bank. The stream is at 600 m elevation and bordered by moss-covered granite boulders less than 4 m wide beneath a canopy of broad-leaved evergreen forest. The stream supports many small fishes (Poropuntius cf. deauratus, Channa cf. gachua, Schistura sp.) and various species of amphibians (Hylarana sp., Odorrana sp.). The juvenile paratype of L. anakradaya sp. nov. was found on a branch 30 cm above the surface of the water approximately 1 km upstream from the type locality. The description of the new species described here from southern Vietnam continue to underscore the fact that the mountainous areas of southern Vietnam are very likely to harbor additional new species to science. As noted by Nguyen et al. (2021), these recent discoveries along with that of Cyrtodactylus raglai from the same river valley is a clear indication that the true diversity of this area is unrealized. This creates an urgent need for continued field work in this region, and in particular, forested riparian areas where additional range-restricted endemics are likely to be discovered. To date, C. raglai, L. truongi sp. nov. and L. anakradaya sp. nov. are known only from a narrow area within the Song Giang River Valley. In 2014, a significant part of this valley was partially flooded after the construction of the Song Giang Hydropower Station. A second hydropower station on the Song Giang River is currently under construction and is expected to be completed in 2022. Economic growth in Khanh Hoa Province requires additional electricity, and the plans for further development of hydropower stations on the Song Giang River will likely have a deleterious impact on the hydrological regime of the river and the surrounding riparian ecosystems. Currently, this area has no legal protection. Therefore, the need for additional herpetological surveys in the Song Giang River Valley is not only urgent but crucial for estimating biodiversity and importance of this region for nature conservation in southern Vietnam. Given the relatively small estimated range of these three new species, the increasing threats to their habitat requires additional field work to clarify the extent of their distribution, population trends, and conservation status. Because these are the only known specimens of L. truongi sp. nov. and L. anakradaya sp. nov., we suggest they should be categorized as Data Deficient (DD) according to the IUCN's Red List categories (IUCN Standards and Petitions Committee, 2017).