Print
A new granite cave-dwelling Bent-toed Gecko from Vietnam of the Cyrtodactylus irregularis group (Squamata; Gekkonidae) and a discussion on cave ecomorphology
expand article infoAnh The Nguyen, Tang Van Duong§, L. Lee Grismer|, Nikolay A. Poyarkov#
‡ Vietnambirds Foundation, Ho Chi Minh, Vietnam
§ Vietnam National Museum of Nature, Vietnam Academy of Science and Technology, Hanoi, Vietnam
| La Sierra University, Riverside, United States of America
¶ M. V. Lomonosov Moscow State University, Moscow, Russia
# Joint Russian-Vietnamese Tropical Research and Technological Center, Hanoi, Vietnam
Open Access

Abstract

An integrative analysis of genetic, morphological, and ecological data recovered a new granite cave-adapted species, Cyrtodactylus raglai sp. nov., from the Song Giang River Valley, Khanh Hoa Province, Vietnam. Cyrtodactylus raglai sp. nov. is nested within one of two major clades within the irregularis species group where it forms a monophyletic group with C. cryptus and its sister species C. kingsadai. It differs from its sister species by an uncorrected pairwise sequence divergence of 16.5% and 16.8% based on the mitochondrial genes NADH dehydrogenase subunit 2 and its flanking tRNAs and the first subunit of cytochrome oxidase c (CO1), respectively. Cyrtodactylus raglai sp. nov. is a narrow-range endemic restricted to a riparian, granite cave microhabitat and its overall morphology bears that of other granite cave-dwelling ecomorphs in the genus. The Song Giang River Valley and its associated waterways are currently threatened by the construction of a hydropower station which will negatively impact the surrounding ecosystem. Urgent field surveys in this region are necessary in order to obtain critical data on its biodiversity and importance to conservation efforts in southern Vietnam.

Key words

Ecomorph, Indochina, narrow-range endemic, reptiles, taxonomy

Introduction

The Cyrtodactylus irregularis group is the largest, most taxonomically diverse species group within the genus. It contains 22 nominal species and at least 11 undescribed species (Grismer et al 2021) that range throughout central and southern Vietnam, eastern and southern Laos, and is represented by isolated populations in northern and eastern Cambodia (Fig. 1). The diversity of the irregularis group is centered in Vietnam wherein 18 of the 20 nominal species occur, most of which are narrow-range endemics (Orlov et al 2007; Nazarov et al 2008, 2012; Ngo and Bauer 2008; Geissler et al 2009; Nguyen et al 2013, 2017; Ziegler et al 2013; Luu et al 2016, 2017; Pauwels et al 2018; Neang et al 2020; Ostrowski et al 2020). The taxonomic diversity of the irregularis group is mirrored by its wide range in habitat preferences, having species adapted to karstic and granitic forests, tree trunks, and generalist that range through many different habitats (see above references). Grismer et al (2020) demonstrated that species within the irregularis group bearing the same habitat preference did not form monophyletic subgroups, but were nested in different places throughout the group’s phylogeny. Nonetheless, the statistically inferred ancestral habitat preference for the group was that of a habitat generalist from which all other habitat preferences ultimately evolved (Grismer et al 2020). Reported here is a newly discovered population of the irregularis group from Khanh Hoa Province, Vietnam. Genetic data from the mitochondrial gene NADH dehydrogenase subunit 2 and its flanking tRNAs (ND2) and the first subunit of cytochrome oxidase c (CO1), indicate it is not nested within any other species of the irregularis group. Morphological analyses demonstrate it bears several statistically significant mensural and meristic differences from its closest relatives, and an ecological assessment of its lifestyle indicates it is the first member of the irregularis group with a granite cave habitat preference. Integrating the results of these different analyses leads us to the robust hypothesis that this population is a new species and as such, is described below.

Figure 1. 

Distribution of the species of the Cyrtodactylus irregularis group. Dot in the center of an icon indicates the type locality of a species. The small inset map shows Vietnam and the location of the area of studies detailed in a larger map (red rectangle). For details on species localities see Pauwels et al. (2018); Neang et al. (2020); and Ostrowski et al. (2020).

Materials and methods

Morphological data and analyses

Morphological data included both meristic and mensural characters. Measurements were taken on the left side of the body when possible to the nearest 0.1 mm by NAP using Mitutoyo dial calipers under a Nikon SMZ 1500 dissecting microscope and follow Grismer and Grismer (2017) and Grismer et al (2018). Measurements taken were: snout-vent length (SVL), taken from the tip of the snout to the vent; tail length (TL), taken from the vent to the tip of the tail; tail width (TW), taken at the base of the tail immediately posterior to the postcloacal swelling; sternal width (SW), distance between the ventral edges of the opposing glenoid fossae; forearm length (FL), taken on the ventral surface from the posterior margin of the elbow while flexed 90° to the inflection of the flexed wrist; forelimb width (FLW), measured from the anterior to the posterior margins of the brachium immediately adjacent to their insertion points on the body; tibia length (TBL), taken on the ventral surface from the posterior surface of the knee while flexed 90° to the base of the heel; hind limb width (HLW), measured from a point equidistant between its anterior and posterior insertion points on the body to the tip of the fourth toe; axilla to groin length (AG), taken from the posterior margin of the forelimb at its insertion point on the body to the anterior margin of the hind limb at its insertion point on the body; pelvic width (PW), distance between the dorsal tips of the ilia; pelvic height (PH), distance from the dorsal tip of the ilium to the ventral surface of the pubis; head length (HL), the distance from the posterior margin of the retroarticular process of the lower jaw to the tip of the snout; head width (HW), measured at the angle of the jaws; head depth (HD), the maximum height of head measured from the occiput to base of the lower jaw; eye diameter (ED), the greatest horizontal diameter of the eye-ball; eye to ear distance (EE), measured from the anterior edge of the ear opening to the posterior edge of the bony orbit; snout length (SNT), measured from anterior-most margin of the bony orbit to the tip of snout; eye to nostril distance (EN), measured from the anterior margin of the bony orbit to the posterior margin of the external nares; interorbital distance (IO), measured between the anterior-most edges of the bony orbits; ear length (EL), measured as the greatest vertical distance of the ear opening; and internarial distance (IN), measured between the nares across the rostrum.

Meristic characters evaluated were the number of supralabial scales (SL) counted from the largest scale immediately below the middle of the eyeball to the rostral scale; infralabial scales (IL) counted from the mental to the termination of enlarged scales just after the upturn of the mouth; the number of paravertebral tubercles (PVT) between limb insertions counted in a straight line immediately left of the vertebral column; the number of longitudinal rows of body tubercles (LRT) counted transversely across the center of the dorsum from one ventrolateral fold to the other; the number of longitudinal rows of ventral scales (VS) counted transversely across the center of the abdomen from one ventrolateral fold to the other; the number of expanded subdigital lamellae on the fourth toe (4TL) counted from the base of the first phalanx where it contacts the body of the foot to the claw including the claw sheath (see Murdoch et al 2019: Fig. 2); and the total number of enlarged femoral scales (FS) from each thigh combined as a single metric. In some species, only the distalmost FS are greatly enlarged and the proximal scales are smaller whereas in others, the enlarged scales are continuous with the enlarged precloacal scales and the separation of the two scales rows was determined to be at a point even with the lateral body margin (see Murdoch et al 2019: Fig. 3). The number of enlarged precloacal scales (PS); the number of precloacal pores in (PP) in males; the number of rows of enlarged post-precloacal scales (PPS) on the midline between the enlarged precloacal scales and the vent; the number of postcloacal tubercles (PCT) on either side of the base of the tail; and the number of dark body bands (BB) between the nuchal loop (the dark band running from eye to eye across the nape) and the hind limb insertions. Other morphological characters evaluated were whether or not the enlarged femoral scales and precloacal scales were separated by a diastema at the base of the femora or contiguous and if so, if the proximal femoral scales were less than one-half the size of the distal femorals (see Murdoch et al 2019: Fig. 3); and the presence or absence of a pocket in the skin webbing between the digits of the hind and forefeet (see Murdoch et al 2019: Fig. 2).

Figure 2. 

Majority-rule consensus tree from 1000 ML bootstrap pseudoreplicates of Cyrtodactylus. Phylogeny is based on 1469 bp of ND2 Black circles represent nodes the UB and BPP support > 90 and 0.90, respectively. Grey circles represent nodes with UB support > 90 only. The white circle represents a node with BPP support > 0.90 only. Nodes lacking circles are unsupported. Scale bar denotes uncorrected pairwise sequence divergence. Species group designations follow Grismer et al. (2021).Photograph by Anh The Nguyen.

Figure 3. 

Majority-rule consensus tree from 1000 ML bootstrap pseudoreplicates of the Cyrtodactylus irregularis group. Phylogeny is based on 691 bp of CO1. Black circles represent nodes with UB support > 90. Nodes lacking circles are unsupported. Scale bar denotes uncorrected pairwise sequence divergence. Photograph by Anh The Nguyen.

Color pattern characters evaluated were the nuchal loop being continuous from orbit to orbit across the nape, or separated medially into paravertebral halves, bearing an anterior azygous notch or not, and the shape of its posterior border; presence or absence of a band on the nape; dorsal body bands bearing paired, having paravertebral elements or not; dark dorsal body bands wider than light interspaces, with or without lightened centers, edged with light tubercles or not, irregularly shaped or more regular (straight-edged or even); dark markings present or absent in the dorsal interspaces; top of head bearing combinations of dark, diffuse mottling or dark, distinct blotches overlain with a light-colored reticulating network or not; anterodorsal margin of thighs and brachia whitish due to a lack of dark pigment; light caudal bands bearing dark markings or immaculate; light caudal bands encircling the tail or not; and dark caudal bands wider than light caudal bands.

Analyses of variance (ANOVA) were conducted on selected meristic and mensural characters (see below) with statistically similar variances (i.e. p values ≤ 0.05 in a Levene’s test) to search for the presence of statistically significant mean differences (p < 0.05) among species across a subset of the irregularis group. Character means bearing statistical differences were subjected to a TukeyHSD test to ascertain which species pairs differed significantly from each other for those particular characters. Boxplots were generated in order to visualize the range, mean, median, and degree of differences between pairs of species bearing statistically different mean values. All statistical analyses were performed in R [v3.4.3].

The morphospatial clustering among the sampled individuals from the selected subset of species were visualized using principal component analysis (PCA) from the ADEGENET package in R (Jombart et al 2010) implemented by the prcomp () command. PCA is a dimension reducing algorithm that decreases the complexity of a data set by finding a subset of input variables that contain the most relevant information (i.e. the greatest variance in the data) while de-emphasizing those characters that do not, thus increasing the overall accuracy of the results by eliminating noise and the potential of overfitting (Agarwal et al 2007). PCA is an unsupervised analysis that recovers morphospatial relationships among the sampled individuals (i.e. data points) and how well they form uncoerced clusters that may or may not align themselves with putative species boundaries delimited by phylogenetic analyses and defined by univariate analyses (see below). It is important to note that clusters of conspecific individuals were not pre-defined in the analysis, but simply color-coded in the scatter plot in order to observe their positions and morphospatial relationships with respect to one another. Given that the meristic data ranged from 8–50, these data were log-transformed prior to analysis in order to normalize their distribution so as to ensure characters with very large or very low values could not over-leverage the results owing to intervariable nonlinearity. To remove potential effects of allometry in the mensural characters, size was normalized using the following equation: Xadj=log(X)-β[log(SVL)-log(SVLmean)], where Xadj=adjusted value; X=measured value; β=unstandardized regression coefficient for each population; and SVLmean=overall average SVL of all populations (Thorpe 1975, 1983; Turan 1999; Lleonart et al 2000). The morphometrics of each species were normalized separately and then concatenated so as not to conflate intra- with interspecific variation (Reists 1986). All data were scaled to their standard deviation to insure they were analyzed on the basis of correlation and not covariance.

A discriminant analysis of principal components (DAPC) was performed on both data sets. The DAPC places individuals from each predefined population into separate clusters (i.e., plots of points) bearing the smallest within-group variance that produce linear combinations of centroids having the greatest between-group variance (i.e. linear distance; Jombart et al 2010). DAPC relies on scaled data calculated from its own internal PCA as a prior step to ensure that variables analyzed are not correlated and number fewer than the sample size. Dimension reduction of the DAPC prior to plotting, is accomplished by retaining the first set of PCs that account for approximately 90–95% of the variation in the data set (Jombart and Collins 2015) as determined from a scree plot generated as part of the analysis. Retaining too many PCs forces false structure to appear in the data while retaining too few, runs the risk of missing true structure (Cangelosi and Goriely 2007).

The type material is deposited in the herpetological collection of the Department of Zoology, Southern Institute of Ecology (SIEZC) in Ho Chi Minh City, Vietnam, and the Zoological Museum of Lomonosov Moscow State University, Moscow, Russia (ZMMU).

Genetic data and analyses

A data set composed of the 310 species of Cyrtodactylus from Grismer et al (2021) plus one specimen of Cyrtodactylus sp. from Song Giang River Valley, Khanh Hoa Province, Vietnam (ZMMU R16688), was assembled and aligned in Geneious. A phylogeny was constructed based on 1469 base pairs of ND2. Mediodactylus russowii and Hemidactylus frenatus were used as outgroups to root the tree based on Gamble et al (2012).

Additionally, to compare the newly discovered population of Cyrtodactylus sp. from Khanh Hoa Province with other members of the irregularis species group, we generated a 672 base pair segment of the 5’-end of (CO1) which included 29 sequences of 26 species and compared it with the CO1-barcoding dataset from Pauwels et al (2018) and Neang et al (2020). Cyrtodactylus intermedius was used as an outgroup to root the tree based on Wood et al (2012. Since at present, the CO1 dataset for members of the irregularis species group has better taxon coverage than the ND2 dataset, we herein applied the CO1 fragment primarily as a DNA-barcoding marker (Murphy et al 2013; Brenan et al 2017), and did not rely on it for reconstructing phylogenetic relationships within the genus. Cyrtodactylus intermedius CO1 sequence was used as an outgroup to root the tree based on Grismer et al (2021). Due to a significant lack of overlap of specimens sequenced for both genes, we elected not to concatenate the gene fragments from the different data sets.

Genomic DNA was isolated from liver tissue stored in 95% ethanol using the following DNA extraction protocol. A small amount of liver tissue (about 1–3 mg) was dried and transferred to 500µl of extraction buffer including 1% SDS, 0.4 mg/ml proteinase K, 10mM Tris HCl (pH 8.0), 1mM EDTA, 100mM NaCl, and subsequently incubated at 56 °C for 3 hours, then centrifuged at 12,000 rpm for one minute. The 450 µl of supernatant were transferred into a new 1.5 ml tube, where 50 µl of 3M sodium acetate solution and 500µl of ice-cold isopropanol were added, mixed and refrigerated at -20 °C for 15 minutes, again centrifuged again at 12,000 rpm for one minute, the supernate was discarded, and the DNA sediment was washed with 500 µl of 70% ethanol. After the ethanol was removed, the total DNA was dried at 50 °C for 8 minutes and subsequently resuspended with 50 µl TE.

The ND2 gene, with parts of adjacent tRNAs, was amplified using a double-stranded Polymerase Chain Reaction (PCR) under the following conditions: 1.0 µl genomic DNA (10–30 µg), 1.0 µl light strand primer (concentration 10 µM), 1.0 µl heavy strand primer (concentration 10 µM), 15 µl Master Mix 2x (CWBIO, China), and 12 µl ultra-pure H2O. PCR reactions were executed on Bio-Rad T100™ gradient thermocycler under the following conditions: initial denaturation at 94 °C for 5 min, followed by a second denaturation at 94 °C for 60 s, annealing at 58 °C for 60 s, followed by a cycle extension at 72 °C for 60 s, for 35 cycles with the final extension step 72 °C for 10 min. For amplification of the full ND2 gene with parts of adjacent tRNAs we used the NADH2-Metf6 (5’-AAGCTTTCGGGCCCATACC-3’) and CO1H (5’-AGRGTGCCAATGTCTTTGTGRTT-3’) primers following Macey et al (1997).

We also amplified the 5’-fragment of the CO1 mtDNA gene; the PCR parameters and the amplification protocol and conditions were the same as for the ND2 gene, but the annealing temperature was at 45 °C. Primers used both for PCR and sequencing were the VF1-d (5’-TTCTCAACCAACCACAARGAYATYGG-3’) and the VR1-d (5’-TAGACTTCTGGGTGGCCRAARAAYCA-3’) following Ivanova et al (2006).

All PCR products were visualized using 1.0 % agarose gel electrophoresis. Successful PCR products were sent to National Key Laboratory of Gene Technology (Institute of Biotechnology, VAST, Hanoi, Vietnam) for PCR purification, cycle sequencing, sequencing purification, and sequencing using the same primers as in the amplification step. Sequences were analyzed from both the 3’ and the 5’ ends separately to confirm congruence between reads. Forward and reverse sequences were uploaded and edited in Geneious 2019.0.4 (https://www.geneious.com). Following sequence editing we aligned the protein-coding region of ND2 and the flanking tRNAs using the MAFTT v7.017 (Katoh and Kuma, 2002) plugin under the default settings in Geneious 2019.0.4 (https://www.geneious.com). Mesquite v3.04 (Maddison and Maddison, 2015) was used to calculate the correct amino-acid reading frame and to confirm the lack of premature stop codons in the ND2 portion of the DNA fragment. The 5’-end partial CO1 sequences were aligned in MAFTT v7.017, and similarly checked for the absence of stop-codons in Mesquite v3.04.

The sequence of ND2 gene and flanking tRNAs fragment and the 5’-end partial sequence of CO1 gene of the newly collected specimen ZMMU R16688 from Khanh Hoa Province, Vietnam, were submitted to GenBank with the accession numbers MW675652 and MW675653, respectively [available upon manuscript acceptance]. GenBank accession numbers for the ND2 sequences of the remaining Cyrtodactylus species in the dataset and the outgroups are listed in Grismer et al (2021).

Phylogenetic analyses

A Maximum likelihood (ML) analysis for the ND2 data set was implemented 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 ModelFinder (Kalyaanamoorthy et al 2017) which selected TVM+F+I+G4 for codon position 1, TIM3+F+I+G4 for codon position 2, GTR+F+ASC+G4 for codon positon 3, and GTR+F+I+G4 for the tRNAs. The best-fitting model selected for CO1 dataset was GTR+G for the first and the third codon positions and HKY+I for the second codon position as suggested by the Akaike Information Criterion (BIC). One-thousand bootstrap pseudoreplicates via the ultrafast bootstrap (UB; Hoang et al 2018) approximation algorithm were employed, and nodes having UB values of 95 and above were considered strongly supported (Minh et al 2013). We considered nodes with values of 90–94 as weakly supported. MEGA7 (Kumar et al. 2016) was used to calculate uncorrected pairwise sequence divergence among the ingroup species using the complete deletion option which removes gaps and missing data from the alignment prior to the anlysis.

A Bayesian phylogenetic tree (BI) was estimated using Bayesian Evolutionary Analysis by Sampling Trees (BEAST) version 2.4.6 (Drummond et al 2012) implemented in CIPRES (Cyberinfrastructure for Phylogenetic Research; Miller et al 2010). Input files were constructed in Bayesian Evolutionary Analysis Utility (BEAUti) version 2.4.6 using a lognormal relaxed clock with unlinked site models, linked trees and clock models, and a Yule prior and run in BEAST version 2.4.6 (Drummond et al 2012) on CIPRES (Cyberinfrastructure for Phylogenetic Research; Miller et al 2010). bModelTest, implemented in BEAST, was used to numerically integrate over the uncertainty of substitution models while simultaneously estimating phylogeny using Markov chain Monte Carlo (MCMC). MCMC chains were run for 300,000,000 generations and logged every 30,000 generations. The BEAST log file was visualized in Tracer v. 1.6.0 (Rambaut et al 2014) to ensure effective sample sizes (ESS) were well-above 200 for all parameters. A Maximum clade credibility tree using mean heights at the nodes was generated using TreeAnnotator v.1.8.0 (Rambaut and Drummond 2013) with a burn-in of 1000 trees (10%). Nodes with Bayesian posterior probabilities (BPP) of 0.95 and above were considered strongly supported (Huelsenbeck et al 2001; Wilcox et al 2002). We considered nodes with values of 0.90–0.94 as weakly supported.

Species delimitation

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 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 were used to recover monophyletic mitochondrial lineages of individual(s) (i.e. populations) in order 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 characters and patterns consistent with the previous designations of the species-level hypotheses—the construction of boundaries representing the hypothesis-testing step of Hillis (2019)—thus providing independent diagnoses to complement the molecular analyses.

Results

Phylogenetic data

For the ND2 dataset, the ML and BI analyses recovered trees with very similar topologies, and the ML topology used here (Fig. 2), is the same as that in Grismer et al (2021). The Cyrtodactylus irregularis group is composed of two clades (clade 1 and 2) which are strongly supported in the ML analysis (UB 99 and 96, respectively) but only clade 1 is strongly supported in the BI analysis (BPP 0.96) . The new population from Song Giang River Valley is nested within clade 1 and forms a well-supported monophyletic group in both analyses with C. cryptus and its sister species C. kingsadai (100/1.00). Together, they form the strongly supported (99/1.00) sister lineage to the remaining species of clade1. The Song Giang River Valley population bears a 13.3% uncorrected pairwise sequence divergence from its sister species C. kingsadai and a 16.5% divergence from C. cryptus. From the remaining species in clade 1, its difference ranges from 8.1–17.0% (Table 1).

The CO1 barcoding data set recovered a polytomy with no backbone support and little resemblance to the ND2 tree (Fig. 3). In both trees however, the new population from the Song Giang River Valley formed a monophyletic group with C. cryptus, an undescribed population from Ba Na, Danang Province, Vietnam, and its sister species C. kingsadai, although none of these relationships were supported in the CO1 tree. The uncorrected pairwise sequence divergence using CO1 between the Song Giang River Valley population and its sister species C. kingsadai was 16.8%, clearly indicating these two populations are not conspecific.

Table 1.

Uncorrected pairwise sequence divergence of ND2 and the flanking tRNAs among the individuals of clade 1 of the Cyrtodactylus irregularis species group. Highlighted cells denote sequence divergence of each species from C. raglai sp. nov.

Species 1 2 3 4 5 6 7 8 9 10 11
1. cryptus 0.00
2. culaochamensis 0.162 0.00
3. kingsadai 0.152 0.150 0.00
4. pseudoquadrivirgatus 0.164 0.094 0.143 0.00
5. raglai sp. nov. 0.165 0.165 0.133 0.170 0.00
6. sp. HLM0316 Kon Ka Kinh 0.160 0.115 0.136 0.122 0.159 0.00
7. sp. HLM0354 Kon Tum 0.157 0.082 0.153 0.087 0.169 0.122 0.00
8. sp. HLM0365 Kon Ka Kinh 0.166 0.134 0.152 0.128 0.177 0.132 0.133 0.00
9. sp. HLM0366 Chu Mom Ray 0.172 0.122 0.150 0.119 0.172 0.123 0.118 0.105 0.00
10. sp. NAP08781 Song Thanh 0.150 0.081 0.143 0.081 0.159 0.122 0.040 0.125 0.120 0.00
11. taynguyenensis 0.159 0.088 0.142 0.100 0.158 0.113 0.096 0.121 0.116 0.097 0.00

Meristic data

The PCA, based on the meristic data set, demonstrated that the new population from Song Giang River Valley clustered separately from its closest relatives and overlapped C. kingsadai, along PC2 (Fig. 4A). PC1 accounts for 43.8% of the variation in the data and loads most heavily for longitudinal rows of tubercles, fourth toe lamellae, and enlarged femoral scales (LPT, 4TL, and FS, respectively; Fig. 4C; Table 2). PC2 accounts for an additional 26.1% of variation in the data set and loads most heavily for labial scales (SL and IL) and ventral scales (VS). Similarly, the Song Giang River Valley population is well-separated from all other species in the meristic DAPC with no overlap of its 67% inertia ellipses (Fig. 4B). The ANOVAs and subsequent TukeyHDS tests of the meristic data demonstrated that the Song Giang River Valley population bears statistically different mean values between it and C. cryptus and C. kingsadai across various combinations of characters (Fig. 4D; Table 3).

Figure 4. 

A. PCA of the Cyrtodactylus raglai sp. nov, C. cryptus, and C. kingsdai based on meristic characters. B. DAPC of same. C. Histograms of the factor loadings of the characters contributing the most to the variation along PC1 and PC2. D. Boxplot comparisons of meristic characters. Light-blue circles are means and black horizontal bars are medians. Species pairs above the plots are those that differ significantly (p < 0.05) from each other based on the ANOVAs and TukeyHSDs.

Table 2.

Summary statistics and principal component analysis scores for meristic characters of Cyrtodactylus raglai sp. nov., C. cryptus, and C. kingsadai. Abbreviations are listed in the Materials and methods.

PC1 PC2 PC3 PC4 PC5 PC6
Standard Deviation 1.622433547 1.25259312 0.919763047 0.763032664 0.520878109 0.314996742
Proportion of Variance 0.43872 0.2615 0.14099 0.09704 0.04522 0.01654
Cumulative Proportion 0.43872 0.70021 0.84121 0.93824 0.98346 1
Eigenvalue 2.632290614 1.568989524 0.845964063 0.582218847 0.271314004 0.099222947
SL 0.358949999 -0.471241487 –0.267598779 0.585203794 –0.411445547 –0.256371476
IL 0.123838432 -0.696295734 0.338218543 –0.001681462 0.609470686 0.118267351
LRT 0.564321114 0.230965469 –0.030580096 0.175925709 0.016848143 0.772028351
VS 0.32917155 0.397242696 0.680146485 0.313822976 0.089154163 –0.405970153
FS 0.45560393 –0.211652886 0.229149873 –0.69862878 –0.445675149 –0.091705728
4TL 0.470565374 0.192769188 –0.545852602 -0.200042007 0.502393653 –0.388634814
Table 3.

Summary statistics of selected meristic and mensural characters of Cyrtodactylus raglai sp. nov., C. kingsadai, and C. cryptus. Mean values in bold differ significantly from the corresponding mean value in C. raglai sp. nov. (p < 0.05).

Meristic Characters Cyrtodactylus raglai sp. nov. (n=3) Cyrtodactylus kingsadai (n=6) Cyrtodactylus cryptus (n=4)
supralabials (SL)
Mean 10.7 12.4 10.3
± 1 SD 0.58 1.42 0.50
Range 10 or 11 10–14 10 or 11
infralabials (IL)
Mean 9.3 10.0 9.0
± 1 SD 0.58 0.89 1.55
Range 9 or 10 9–11 8–10
longitudinal rows of tubercles (LRT)
Mean 14.3 19.2 19.0
± 1 SD 0.58 3.5 1.41
Range 14 or 15 14–23 17–20
ventral scales (VS)
Mean 37.3 42.2 48.8
± 1 SD 1.53 4.00 1.23
Range 36–39 39–46 47–50
subdigital lamellae on 4th toe (4TL)
Mean 21.3 22.7 22.3
± 1 SD 0.58 1.37 0.96
Range 20 or 21 21–25 20–23
Adjusted Mensural Characters
head length (HL)
Mean 3.3 3.2 3.1
± 1 SD 0.02 0.031 0.02
Range 3.29–3.34 3.21–3.29 3.05–3.09
head width (HW)
Mean 2.9 2.8 2.7
± 1 SD 0.01 0.07 0.09
Range 2.90–2.91 2.72–2.91 2.57–2.77
head depth (HD)
Mean 2.3 2.3 2.2
± 1 SD 0.00 0.06 0.04
Range 2.31–2.32 2.21–2.41 2.16–2.27
eye diameter (ED)
Mean 1.9 1.8 1.6
± 1 SD 0.02 0.07 0.08
Range 1.88–1.91 1.73–1.91 1.47–1.67
snout length (SNT)
Mean 2.5 2.4 2.2
± 1 SD 0.036 0.18 0.01
Range 2.49–2.56 2.13–2.61 2.22–2.24

Mensural data

The PCA, based on the mensural data set, demonstrated that the Song Giang River Valley population clustered separately from its closest relatives and overlapped with them along PC2 but not PC1 (Fig. 5A). PC1 accounted for 72.9% of the variation in the data set and loaded most heavily for head length, head width, and eye diameter (HL, HW, ED, respectively; Fig. 5C; Table 4). PC2 accounted for an additional 15.8% of the variation and loaded most heavily for head depth and snout length (HD and SNT, respectively). Similarly, the Song Giang River Valley population is well-separated from all species in the mensural DAPC with no overlap of its 67% inertia ellipses (Fig. 5B). The ANOVAs and subsequent TukeyHDS tests of the mensural data demonstrated that the Song Giang River Valley population bears statistically different mean values between it and C. cryptus and C. kingsadai across various combinations of characters (Fig. 5D; Table 3).

Based on the phylogenetic, PCA, DAPC, and ANOVA analyses and the high degree of uncorrected pairwise sequence divergence between the Song Giang River Valley population and the remaining species in clade 1 (Figs 1, 2, 3), we hypothesize that the Song Giang River Valley population is a discretely diagnosable lineage which shows no evidence of reticulation with any other lineage and as such, should be accorded species status.

Figure 5. 

A. PCA of the Cyrtodactylus raglai sp. nov, C. cryptus, and C. kingsdai based on mensural characters. B. DAPC of same. C. Histograms of the factor loadings of the characters contributing the most to the variation along PC1 and PC2. D. Boxplot comparisons of meristic characters. Light-blue circles are means and black horizontal bars are medians. Species pairs above the plots are those that differ significantly (p < 0.05) from each other based on the ANOVAs and TukeyHSDs. Photographs by Anh The Nguyen.

Table 4.

Summary statistics and principal component analysis scores for mensural characters of Cyrtodactylus raglai sp. nov.,C. cryptus, and C. kingsadai. Abbreviations are listed in the Materials and methods.

PC1 PC2 PC3 PC4 PC5
Standard Deviation 1.90913762 0.890016517 0.538548404 0.467906735 0.23257913
Proportion of Variance 0.72896 0.15843 0.05801 0.04379 0.01082
Cumulative Proportion 0.72896 0.88739 0.94539 0.98918 1
Eigenvalue 3.644806451 0.792129401 0.290034383 0.218936713 0.054093052
HL –0.50087707 0.117559638 –0.026802327 0.458852951 0.723904343
HW –0.459584098 –0.154364348 –0.808742489 –0.307413883 –0.128009361
HD –0.422351581 –0.518980378 0.513977677 –0.518283228 0.139599258
ED –0.49614577 –0.095133377 0.203741279 0.525612741 –0.653459833
SE –0.336694424 0.827018977 0.198830734 –0.387381487 –0.114360833

Taxonomy

Cyrtodactylus raglai sp. nov.

Suggested Common Name: Raglai Bent-toed Gecko; Figures 6, 7, 9

Holotype

Adult male, SIEZC 2.0244, collected 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, by Anh The Nguyen on 22 September 2020.

Paratypes

Two adult females, SIEZC 2.0243 and ZMMU R16688 bearing the same data as the holotype except that the latter was collected on 16 August 2020.

Diagnosis

Cyrtodactylus raglai sp. nov. can be separated from all other species of clade1 of the C. irregularis group by having 10 or 11 supralabials; nine or 10 infralabials; 44–47 paravertebral tubercles; 14 or 15 rows of longitudinally arranged tubercles; 36–39 ventrals; 8–10 expanded subdigital lamellae, 12 or 13 unexpanded subdigital lamellae, and 21–22 total subdigital lamellae on the fourth toe; 18 or 19 enlarged femorals; 12 enlarged precloacals; four rows of enlarged post-precloacals; three postcloacal tubercles in males; five precloacal pores in the male; no pitted precloacal scales in females; enlarged femorals and enlarged precloacals not continuous; proximal femorals less than one-half size of distal femorals; enlarged subcaudals; maximum SVL 111.7 mm; small, irregularly shaped dark blotches on top of head; and four irregularly shaped body bands edged with white tubercles wider than the interspaces (Tables 5, 6).

Table 5.

Diagnostic character states among the nominal species of clade 1 of the Cyrtodactylus irregularis group. / = data unavailable.

Cyrtodactylus raglai sp. nov. Cyrtodactylus cryptus Cyrtodactylus culaochamensis Cyrtodactylus kingsadai Cyrtodactylus pseudoquadrivirgatus Cyrtodactylus taynguyenensis
Maximum SVL (mm) 111.7 90.8 79.8 94.0 83.3 104.1
Ventral scales 36–39 47–50 45–50 39–46 41–57 42–49
Enlarged femoral scales 18 or 19 / / 18–24 / 22–26
Femoral pores 0 / 0–7 0 0
Precloacal pores (males) 5 9–11 7 or 8 7–9 6 6
Precloacal pits (females) 0 0 4–8 5–10 0
4th toe lamellae 21 or 22 20–23 20–23 21–25 16–25 17–21
Enlarged subcaudals yes no yes yes no no
Dorsal bands irregular irregular paired-paravertebral irregular irregular irregular to blotched
Table 6.

Meristic, mensural (mm), discrete morphological, and color pattern data from the type series of Cyrtodactylus raglai sp. nov. Abbreviations are listed in the Materials and methods.

Scale characters SIEZC 2.0244 male holotype SIEZC 2.0243 female paratype ZMMU R16688 female paratype
Supralabials (SL) 10 11 11
Infralabials (IL) 9 10 9
Body tubercles low, weakly keeled yes yes yes
Body tubercles raised, moderately to strongly keeled no no no
Paravertebral tubercles (PVT) 46 47 44
Longitudinal rows of body tubercles (LRT) 14 14 15
Tubercles extend beyond base of tail yes yes yes
Ventral scales (VS) 36 39 37
Expanded subdigital lamellae on 4th toe (TLE) 10 9 8
Unmodified subdigital lamellae on 4th toe (TLU) 12 12 13
Total subdigital lamellae on 4th toe (TTL) 22 21 21
Enlarged femoral scales (R/L) 10/9 9/9 10/9
Total femoral scales (FS) 19 18 19
Femoral pores (R/L) 0 0 0
Total femoral pores in males 0 0 0
Enlarged precolacal scales (PS) 12 12 12
Precloacal pores (PP) 5 0 0
Post-precloacal scale rows (PPS) 4 4 4
Enlarged femoral and precloacal scales continuous no no no
Pore-bearing femoral and precloacal scales continuous no no no
Enlarged proximal femoral scales ~1/2 size of distal femorals yes yes yes
Post-cloacal tubercles (PCT) 3 3 3
Medial subcaudals 2 or 3 times wider than long yes yes yes
Medial subcaudals extend up onto lateral surface of tail no no no
Color pattern characters
Nuchal loop divided medially no no no
2 posterior projections from nuchal loop no no no
Nuchal loop with anterior azygous notch no no no
Posterior border of nuchal loop jagged jagged jagged
Band on nape no no no
Dorsal banding with divided paravertebral elements no no no
Number of body bands (not including nuchal loop) 4 4 4
Dorsal body bands wider than interspaces yes yes yes
Dorsal body bands with lightened centers vertebrally vertebrally vertebrally
Dorsal bands edged with white tubercles yes yes yes
Shape of dorsal bands irregular irregular irregular
Dark markings in dorsal interspaces yes yes yes
Top of head diffusely mottled, blotched, or patternless blotched blotched blotched
Well-defined, light-colored reticulum on top of head no no no
Anterodorsal margin of thighs darkly pigmented no no no
Anterodorsal margin of brachia darkly pigmented no no no
White caudal bands with dark markings yes yes yes
White caudal bands encircle tail no no no
Number of light caudal bands 7 9 8
Number of dark caudal bands 8 10 8
Dark caudal bands wider than light caudal bands yes yes yes
SVL 95.0 111.7 87.5
TL 119.0 135.0 113.4
TW 9.0 9.8 5.8
FL 15.5 18.4 14.2
FLW 5.3 7.6 5.1
TBL 19.0 21.5 16.7
HDW 9.8 12.25 9.1
AG 41.0 46.0 39.7
PW 11.1 14.5 10.4
PH 8.6 11.5 8.1
HL 26.2 31.0 26.0
HW 17.6 21.8 15.4
HD 10.0 11.0 9.5
ED 6.5 7.8 5.6
EE 7.2 8.8 5.9
SNT 11.7 14.5 11.6
EN 9.2 11.5 8.3
IO 3.6 4.3 3.3
EL 3.0 3.4 3.0
IN 3.1 3.4 3.2

Description of holotype

(Figs 6, 7). Adult male SVL 95.0 mm; head moderate in length (HL/SVL 0.28) and width (HW/HL 0.66), flattened (HD/HL 0.38), distinct from neck, triangular in dorsal profile; lores concave anteriorly, weakly inflated posteriorly; prefrontal region slightly concave; canthus rostralis rounded; snout elongate (ES/HL 0.42), flat, rounded in dorsal profile; eye large (ED/HL 0.25); ear opening narrow, elliptical, obliquely oriented, moderate in size; eye to ear distance slightly greater than diameter of eye; rostral rectangular, partially divided dorsally by inverted Y-shaped furrow, bordered posteriorly by large left and right supranasals and one small azygous internasal, bordered laterally by first supralabials; external nares bordered anteriorly by rostral, dorsally by large supranasal, posteriorly by two moderately sized postnasals, bordered ventrally by first supralabial; 10(R,L) rectangular supralabials extending to below midpoint of eye, second supralabial slightly larger than first; 9(R,L) infralabials tapering smoothly to slightly past the termination of enlarged supralabials; scales of rostrum and lores flat, same size as granular scales on top of head and occiput; scales of occiput intermixed with small, distinct, tubercles; superciliaries elongate, largest anteriorly; mental triangular, bordered laterally by first infralabials and posteriorly by large left and right trapezoidal postmentals contacting medially for ~60% of their length posterior to mental; one row of slightly enlarged, elongate sublabials extending posteriorly to fourth infralabial; gular and throat scales small, granular, grading posteriorly into slightly larger, flatter, smooth, imbricate, pectoral and ventral scales.

Figure 6. 

Adult male holotype of Cyrtodactylus raglai sp. nov. (SIEZC 2.0244). A. Dorsal view. B. Ventral view of thighs and precloacal region. C. Right lateral view of head. D. Ventral view of gular region. E. Ventral view of left foot. F. Ventral view of left hand. Photographs by Anh The Nguyen.

Body relatively short (AG/SVL 0.43) with well-defined ventrolateral folds; dorsal scales small, granular, interspersed with small, rounded, semi-regularly arranged, smooth tubercles; tubercles extend from occiput onto base of tail forming transverse rows; approximately 14 longitudinal rows of tubercles at midbody; approximately 46 paravertebral tubercles; 36 flat, imbricate, ventral scales much larger than dorsal scales; five large, pore-bearing, precloacal scales; no deep precloacal groove or depression; and four rows of large post-precloacal scales on midline.

Forelimbs thin, relatively long (FL/SVL 0.16); lacking tubercles, granular scales slightly larger than those on body; palmar scales rounded, slightly raised; digits well-developed, inflected at basal interphalangeal joints; digits slightly narrower distal to inflections; subdigital lamellae transversely expanded, those proximal to joint inflections wider than those distal to inflection; claws well-developed, sheathed by a dorsal and ventral scale; hind limbs thin, more robust than forelimbs, long (TBL/SVL=0.20), covered dorsally by granular scales interspersed with slightly larger, weakly keeled tubercles and anteriorly by flat, slightly larger scales; ventral scales of thigh flat, imbricate, larger than dorsals; subtibial scales large, flat, imbricate; one row of 10(R)9(L) enlarged femoral scales terminating distally before knee, not continuous with enlarged precloacal scales; proximal femoral scales much smaller than distal femorals, the latter forming an abrupt union with smaller, rounded, ventral scales of posteroventral margin of thigh; femoral pores absent; plantar scales flat; digits well-developed, inflected at basal interphalangeal joints, slightly narrower distal to inflections; subdigital lamellae transversely expanded, those proximal to joint inflections wider than those distal to inflection, 10(R,L) transversely expanded subdigital lamellae on fourth toe proximal to joint inflection that extends onto the sole; 12(R,L) narrower lamellae distal to inflection; 22 total subdigital lamellae; and claws well-developed, sheathed by a dorsal and ventral scale at base.

Tail long (TL/SVL 1.25), original, 119 mm in length, 9.0 mm in width at base, tapering to a point; dorsal caudals small, generally square; median row of subcaudals transversely expanded, significantly larger than dorsal caudals, not extending up onto lateral side of tail; transverse rows of 2–4 keeled tubercles on anterior one-half of tail, paravertebral tubercles largest; tubercle rows separated by 7–9 rows of dorsal caudals; base of tail bearing small hemipenal swellings with three large postcloacal tubercles on either side; and postcloacal scales flat, imbricate.

Figure 7. 

Type series of Cyrtodactylus raglai sp. nov. A. Adult male holotype SIEZC 2.0244. B. Adult female paratype SIEZC 2.0243. C. Adult female paratype ZMMU R16688. Photographs by Anh The Nguyen.

Coloration in life (Figs 6, 7). Ground color of top of head, limbs, and dorsum gray; top of head heavily mottled with irregularly shaped, small, dark-brown blotches, light-colored reticulate pattern absent; dark-brown, nuchal loop extends from posterior margin of one orbit to posterior margin of other orbit, anterior and posterior margins on nape jagged, no anterior azygous notch; no dark banding on nape; four dark-brown, irregularly shaped, dorsal body bands reaching ventrolateral folds extending from shoulders to pelvis, wider than light-colored interspaces, bearing lightened vertebral areas, edged in light-colored tubercles; interspaces bearing dark markings; limbs dark-brown with irregularly shaped, light-colored markings; eight, wide, dark-brown caudal bands much wider than seven; light-colored caudal bands do not encircle tail and have darkened centers; iris gold bearing thin, black reticulations; edges of pupils orange; venter beige with faint, dark mottling on lateral edges of belly, undersides of limbs; and subcaudal region dark-brown.

Variation

(Fig. 7). The paratypes resemble the holotype very closely in all aspects of coloration and pattern. Specimens ZMMU R16688 and SIEZC 2.0243 have eight and nine light-colored caudal bands, respectively, as opposed to seven in the holotype and SIEZC 2.0243 has 10 dark-brown caudal bands as opposed to eight in the holotype. Meristic and mensural variation is presented in Table 6.

Distribution

(Fig. 1). Cyrtodactylus raglai sp. nov. is known only from the type locality in the Song Giang River valley, Khanh Trung Commune, Khanh Vinh District, Khanh Hoa Province, Vietnam.

Etymology

The new species name “raglai” is given in a reference to the Raglai people, an ethnic group living in the forested mountain areas of Khanh Hoa Province of Vietnam, including the Song Giang River Valley where the new species was found. In Raglai language, the self-designating word “raglai” also means “forest”, stressing the importance of the tropical forest ecosystem for this people. To reflect this polysemy, the new species name is given as a noun in apposition and hence is invariable.

Comparisons

Within the Cyrtodactylus irregularis group, C. raglai sp. nov. is most closely related to C. kingsadai, C. cryptus, and an undescribed species from Ba Na with which it forms a monophyletic group. It differs from both the described species in having a far greater maximum SVL (111.7 mm vs. 90.8–94 mm collectively) and fewer precloacal pores in the single male specimen (five vs. 7–11, collectively). It differs further from its sister species C. kingsadai in having significantly fewer supralabials (10 or 11 vs. 10–14, collectively), longitudinal rows of tubercles (14 or 15 vs. 14–23, collectively), and ventral scales (36–39 vs. 39–46, collectively), and a significantly longer head (adjusted head length 3.29–3.34 vs. 3.21–3.29, collectively) (Table 2; Figs 4, 5). It differs further from C. cryptus by having significantly more supralabials, significantly fewer longitudinal rows of tubercles, ventral scales, enlarged subcaudals, a significantly longer head, longer snout, larger eye, and a significantly narrower head (Table 2: Figs 4, 5). From the other nominal species in clade 1, C. raglai sp. nov. differs having a greater maximum SVL (111.7 mm vs. 79.8–104.1 mm collectively), fewer ventral scales (36–39 vs. 39–57 collectively), fewer precloacal pores in the single male specimen (five vs. 6–11 collectively), no precloacal pits in females, and enlarged subcaudals (also lacking in C. pseudoquadrivirgatus and C. taynguyenensis) (Table 5).

Natural history

Cyrtodactylus raglai sp. nov. was recorded in the forested valley of Song Giang River in the northwestern part of Khanh Hoa Province (Fig. 1), 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 and are dissected by a rich network of small streams and rivulets that feed into the Song Giang River (Fig. 8A). Scattered outcroppings of large granite boulders (ca. 3–5 m in height, 6–8 m in width) occur throughout the river valley and many of these have small permanent streams running through and beneath them. Here, soil erosion has formed small narrow caves among the boulders varying from 3–5 m in depth (Fig. 8B, C) and C. raglai sp. nov. was only observed inside such caves or within crevices between the boulders. Lizards were elusive and difficult to find and only four specimens were recorded during three excursions to the area, and all were restricted only to specific outcroppings. Lizards were active from 19:00 h to 21:00 h, and were usually observed while foraging on the vertical surfaces of the boulders (Fig. 8C) or hiding within the crevices. The two other species of Cyrtodactylus recorded in sympatry with the new species in the surrounding habitat were C. yangbayensis Ngo & Chan and Cyrtodactylus sp., both members of the irregularis group. Both these species were quite numerous throughout the forested habitats of the Song Giang River Valley but occupied different microhabitats than that of C. raglai sp. nov. This type of habitat partitioning was also observed between the granite cave ecomorph C. hontreensis and habitat generalist C. condorensis on Hon Tre Island, Kien Giang Province, and between the granite cave ecomorph C. eisenmanae and C. condorensis on Hon Son Island, Kien Giang Province. In both cases, these species pairs were sympatric but never observed to be syntopic. A similar situation may exist on Ba Den Mountain, Tay Ninh Province among the cave ecomorphs C. bandeennsis and C. nigriocularis and the granite forest species C. thuongae that is not restricted to the caves (Phung et al 2014). Other gecko species recorded in the Song Giang River Valley were Gekko gecko (Linnaeus) and Gekko grossmanni Günther.

Figure 8. 

Habitat (A) and granite cave microhabitat (B and C) of Cyrtodactylus raglai sp. nov. at the type locality. Photographs by Anh The Nguyen.

Discussion

Evolution of cave ecomorphology

Grismer and Grismer (2017) demonstrated that in the Cyrtodactylus condorensis group, the granite cave-adapted species C. eisenmanae and C. grismeri, had a distinctive, gracile overall morphology adapted for climbing on vertical surfaces in low levels of illumination that was very different from the stout, robust morphology of the closely related sympatric habitat generalists C. condorensis and C. leegrismeri. The cave ecomorphs, C. eisenmanae and C. grismeri, have shorter trunks; longer limbs; a longer head and snout; larger eyes; and a flatter pelvis and head. In a genus-wide analysis involving 243 species, Grismer et al (2020) statistically demonstrated that a granite cave habitat preference evolved independently from a general habitat preference at least twice in the Cyrtodactylus badenensiscondorensis group lineage, once in C. hontreensis of the intermedius group (sec. Grismer et al 2021), and at least one more time from a karst-adapted ancestor in C. nigriocularis of the angularis group (sec. Grismer et al 2021). Although no mensural data have yet been taken to bear out the cave-adapted morphology of C. badenensis, C. hontreensis, or C. nigriocularis, it is clear from their overall body stature, they bear many of the characters listed above for C. eisenmanae and C. grismeri (Fig. 9).

Cyrtodactylus raglai sp. nov. is the sixth species of Cyrtodactylus to have a cave habitat preference and the first species to do so in the irregularis group. Morphometric comparisons of C. raglai sp. nov. to other closely related members of the irregularis group could not be evaluated using the complete data set of Grismer and Grismer (2017) because those data for the other irregularis group species do not exist. Therefore, only morphometrics involving head shape and eye size were evaluated. Despite small samples sizes, C. raglai sp. nov. differed statistically from the habitat generalist C. cryptus in having a longer and wider head, longer snout, and larger eyes (Fig. 5D)—characters observed in the cave ecomorphs, C. eisenmanae and C. grismeri. Its granite-adapted sister species C. kingsadai, has an intermediate morphology between C. raglai sp. nov. and C. cryptus and differs statistically from C. raglai sp. nov. only by having a shorter snout. However, given its intermediate position in the PCA, DAPC, and the boxplots (Fig. 5), a larger sample size would likely increase their statistical differences. Nonetheless, based on the limited morphometric data and ecological observations made at the time of collection, we hypothesize that C. raglai sp. nov. is a cave ecomorph similar in morphology to all the other cave ecomorphs (Fig. 9). A study using a complete morphometric data set comparing all the cave ecomorphs to each other and to their closest relatives is currently under way (Grismer et al in progress).

Figure 9. 

General morphology of the six species of cave ecomorphs in the genus Cyrtodactylus from Vietnam. A. C. raglai sp. nov. from the Song Giang River Valley, Khanh Hoa Province. Photograph by Anh The Nguyen. B. C. hontreensis from Hon Tre Island, Kien Giang Province. Photograph by L. Lee Grismer. C. C. grismeri from Tuc Dup Hill, An Giang Province. Photograph by L. Lee Grismer. D. C. badenensis from Ba Den Mountain, Tay Ninh Province. Photograph by Nikolay A. Poyarkov. E. C. nigriocularis from Ba Den Mountain, Tay Ninh Province. Photograph by Nikolay A. Poyarkov. F. C. eisenmanae from Hon Son Island, Kien Giang Province. Photograph by Ngo Van Tri.

Conservation

The description of Cyrtodactylus raglai sp. nov. continues to underscore the fact that the mountainous areas of southern Vietnam harbor more species of Cyrtodactylus than any other area in the country. This recent discovery of yet another narrow-range endemic, is a clear indication that the true diversity of this area is not yet known. This makes it all the more important the urgent need for continued field work in this region, and in particular, forested areas with granitic outcroppings where more narrow-range endemics are likely to be discovered. To date, Cyrtodactylus raglai sp. nov. is 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 strong impact on the hydrological regime of the river and the surrounding forest ecosystems. Currently, this forested area has no legal protection. Therefore, an urgent need for additional herpetological surveys in the Song Giang River Valley is crucial for estimating its biodiversity and importance for nature conservation in southern Vietnam. Given the relatively small estimated range of Cyrtodactylus raglai sp. nov. and the increasing threats to the ecosystem where it was found, further research is urgently required to clarify the extent of its distribution, population trends, and conservation status. Because C. raglai sp. nov. is known only from three specimens, we tentatively suggest it should be categorized as Data Deficient (DD) according to the IUCN’s Red List categories (IUCN Standards and Petitions Committee, 2017).

Acknowledgements

We thank the Bureau of Forestry, Ministry of Agriculture and Rural Development of Vietnam and of administration of Khanh Hoa Province for permitting the fieldwork. Anh The Nguyen thanks Vu Long (Center for Biodiversity conservation and Endangered Species, CBES) and Pham Van Hieu for numerous support and assistance in during the field surveys. Nikolay A. Poyarkov expresses his gratitude to Andrey N. Kuznetsov (Joint Russian-Vietnamese Tropical Research and Technological Center) for supporting his work in Vietnam. Fieldwork, specimen collection, and molecular phylogenetic analysis for this paper were conducted with the financial support of the Russian Science Foundation (RSF Grant No. 19-14-00050 to Nikolay A. Poyarkov). The authors thank two anonymous reviewers for their useful comments on an earlier draft of the manuscript.

References

  • Barraclough TG, Birky CW, Burt A (2003) Diversification in sexual and asexual organisms. Evolution 57: 2166–2172.
  • Brennan IG, Bauer AM, Ngo VT, Wang Y, Wang W, Zhang Y-P, Murphy RW (2017) Barcoding utility in a mega-diverse, cross-continental genus, keeping pace with Cyrtodactylus geckos. Scientific Reports 7: 5592. https://doi.org/10.1038/s41598-017-05261-9
  • De Queiroz K (2007) Species concepts and species delimitation. Systematic Biology 56: 879–886.
  • Drummond AJ, Suchard MA, Xie D, Rambaut A (2012) Bayesian phylogenetics with BEAUti and the BEAST 1.7. Molecular Biology and Evolution 29: 1969–1973. https://doi.org/10.1093/molbev/mss075
  • Geissler P, Nazarov RA, Orlov NL, Phung TM, Nguyen TQ, Ziegler T (2009) A new species of the Cyrtodactylus irregularis complex (Squamata, Gekkonidae) from southern Vietnam. Zootaxa 2161: 20–32. https://doi.org/10.11646/zootaxa.2161.1.2
  • Grismer LL, Grismer JL (2017) A re-evaluation of the phylogenetic relationships of the Cyrtodactylus condorensis group (Squamata; Gekkonidae) and a suggested protocol for the characterization of rock-dwelling ecomorphology in Cyrtodactylus. Zootaxa 4300: 486–504. https://doi.org/10.11646/zootaxa.4300.4.2
  • Grismer LL, Wood PL, Le MD Quah ESH, Grismer JG (2020) Evolution of habitat preference in 243 species of Bent-toed Geckos (genus Cyrtodactylus Gray, 1827) with a discussion of karst habitat conservation. Ecology and Evolution (in press).
  • Grismer LL, Wood PL, Poyarkov NA, Le MD Kraus F, Agarwal I, Oliver PM, Nguyen SN, Nguyen TQ, Karunarathna S, Welton LJ, Stuart BL, Luu VQ, Bauer AM, O’Connell KA, Quah ESH, Chan KO, Ziegler T, Ngo H, Nazarov RA, Aowphol A, Chomdej S, Suwannapoom C, Siler CD, Anuar S, Tri NV, Grismer JL (2021) Phylogenetic partitioning of the third-largest vertebrate genus in the world, Cyrtodactylus Gray, 1827 (Reptilia; Squamata; Gekkonidae) and its relevance to taxonomy and conservation. Vertebrate Zoology 71: 101–154. https://doi.org/10.3897/vz.71.e59307
  • Grismer LL, Wood PL, Thura MK, Thaw Zin, Quah ESH, Murdoch ML, Grismer MS, Lin A, Kyaw H, Lwin N (2018) Twelve new species of Cyrtodactylus Gray (Squamata: Gekkonidae) from isolated limestone habitats in east-central and southern Myanmar demonstrate high localized diversity and unprecedented microendemism. Zoological Journal of the Linnean Society 182: 862–959. https://doi.org/10.1093/zoolinnean/zlx057
  • Hillis DM (2019) Species delimitation in herpetology. Journal of Herpetology 53: 3–12.
  • Hoang DT, Chernomor O, von Haeseler A, Minh BQ, Vinh LS (2018) UFBoot2: Improving the ultrafast bootstrap approximation. Molecular Biology and Evolution 35: 518–522. https://doi.org/10.1093/molbev/msx281
  • Huelsenbeck JP, Ronquist F, Nielsen R, Bollback JP (2001) Bayesian inference of phylogeny and its impact on evolutionary biology. Science 294: 2310–2314. https://doi.org/10.1126/science.1065889
  • Jombart T, Devillard S, Balloux F (2010) Discriminant analysis of principal components: a new method for the analysis of genetically structured populations. BMC Genetics 11: 94. https://doi.org/10.1186/1471-2156-11-94
  • Kalyaanamoorthy S, Minh BQ, Wong TK, von Haeseler A, Jermiin LS (2017) ModelFinder: fast model selection for accurate phylogenetic estimates. Nature Methods 14: 587. https://doi.org/10.1038/nmeth.4285
  • Kumar S, Stecher G, Tamura K (2016) MEGA7: Molecular evolutionary genetics analysis version 7.0 for bigger datasets. Molecular Biology and Evolution 33: 1870–1874. https://dx.doi.org/10.1093/molbev/msw054
  • Lleonart J, Salat J, Torres GJ (2000) Removing allometric effects of body size in morphological analysis. Journal of Theoretical Biology 205: 85–93. https://doi.org/10.1006/jtbi.2000.2043
  • Luu VQ, Dung TV, Nguyen TQ, Le MD, Ziegler T (2017) A new species of the Cyrtodactylus irregularis complex (Squamata: Gekkonidae) from Gia Lai Province, Central Highlands of Vietnam. Zootaxa 4362: 385–404. https://doi.org/10.11646/zootaxa.4362.3.4
  • Luu VQ, Bonkowski M, Nguyen TQ, Le MD, Schneider N, Ngo HT, Ziegler T (2016) Evolution in karst massifs, cryptic diversity among bent-toed geckos along the Truong Son Range with descriptions of three new species and one new country record from Laos. Zootaxa 4107: 101–140. https://doi.org/10.11646/zootaxa.4107.2.1
  • Macey JR, Larson A, Ananjeva NB, Fang Z, Papenfuss TJ (1997) Two novel gene orders and the role of light-strand replication in rearrangement of the vertebrate mitochondrial genome. Molecular Biology and Evolution 14: 91–104.
  • Miller MA, Pfeiffer W, Schwartz T (2010) Creating the CIPRES Science Gateway for inference of large phylogenetic trees. In: Proceedings of the Gateway Computing Environments Workshop (GCE), 14 Nov. 2010, New Orleans, LA, pp. 1-8. https://doi.org/10.1109/GCE.2010.5676129
  • Murdoch ML, Grismer LL, Wood PL, Neang T, Poyarkov NA, Tri NV, Nazarov RA, Aowphol A, Pauwels OSG, Nguyen HN, Grismer JL (2019) Six new species of the Cyrtodactylus intermedius complex (Squamata: Gekkonidae) from the Cardamom Mountains and associated highlands of Southeast Asia. Zootaxa 4554: 1–62. https://doi.org/10.11646/zootaxa.4554.1.1
  • Nazarov RA, Orlov NL, Nguyen SN, Ho CT (2008) Taxonomy of naked-toe geckos Cyrtodactylus irregularis complex of South Vietnam and description of a new species from Chu Yang Sin Natural Park (Krong Bong District, Dac Lac Province), Vietnam. Russian Journal of Herpetology 15: 141–156.
  • Nazarov RA, Poyarkov NA, Orlov NL, Phung TM, Nguyen TT, Hoang DM, Ziegler T (2012) Two new cryptic species of the Cyrtodactylus irregularis complex (Squamata, Gekkonidae) from southern Vietnam. Zootaxa 3302: 1–24. https://doi.org/10.11646/zootaxa.3302.1.1
  • Neang T, Henson A, Stuart BL (2020) A new species of Cyrtodactylus (Squamata, Gekkonidae) from Cambodia’s Prey Lang Wildlife Sanctuary. ZooKeys 926: 133–158. https://doi.org/10.3897/zookeys.926.48671
  • Nguyen SN, Le TNT, Tran DAT, Orlov NL, Lathrop A, Macculloch RD, Le DTT, Jin J-Q, Nguyen LT, Nguyen TT, Hoang DD, Che J, Murphy RW, Zhang Y-P (2013) Phylogeny of the Cyrtodactylus irregularis species complex (Squamata, Gekkonidae) from Vietnam with the description of two new species. Zootaxa 3737: 399–414. https://doi.org/10.11646/zootaxa.3737.4.4
  • Nguyen L-T, Schmidt HA, von Haeseler A, Minh BQ (2015) IQ-TREE: A fast and effective stochastic algorithm for estimating maximum likelihood phylogenies. Molecular Biology and Evolution 32: 268–274. https://doi.org/10.1093/molbev/msu300
  • Nguyen SN, Zhou W-W, Le TNT, Tran DAT, Jin J-Q, Vo BD, Nguyen LT, Nguyen TT, Nguyen TQ, Hoang DD, Orlov NL, Che J, Murphy RW, Zhang Y-P (2017) Cytonuclear discordance, cryptic diversity, complex histories, and conservation needs in Vietnamese bent-toed geckos of the Cyrtodactylus irregularis species complex. Russian Journal of Herpetology 24: 133–154. https://doi.org/10.30906/1026-2296-2019-24-2-133-154
  • Orlov NL, Nguyen TQ, Nazarov RA, Ananjeva NB, Nguyen SN (2007) A new species of the genus Cyrtodactylus Gray, 1827 and redescription of Cyrtodactylus paradoxus (Darevsky et Szczerbak, 1997) [Squamata, Sauria, Gekkonidae] from South Vietnam. Russian Journal of Herpetology 14: 145–152.
  • Ostrowski S, Do DT, Le MD Ngo HT, Pham CT, Nguyen TQ, Nguyen VTH, Ziegler T (2020) A new species of Cyrtodactylus (Squamata: Gekkonidae) from southern Vietnam. Zootaxa 4789: 171–203. https://doi.org/10.11646/zootaxa.4789.1.5
  • Pauwels OSG, Nazarov RA, Bobrov VV, Poyarkov NA (2018) Taxonomic status of two populations of Bent-toed Geckos of the Cyrtodactylus irregularis complex (Squamata, Gekkonidae) with description of a new species from Nui Chua National Park, southern Vietnam. Zootaxa 4403: 307–335. https://doi.org/10.11646/zootaxa.4403.2.5
  • Phung TM, Schingen MV, Ziegler T, Nguyen TQ (2014) A third new Cyrtodactylus (Squamata: Gekkonidae) from Ba Den Mountain, Tay Ninh Province, southern Vietnam. Zootaxa 3764: 347–363.
  • Reist JD (1986) A empirical evaluation of coefficients used in residual and allometric adjustment of size covariation. Canadian Journal of Zoology 64: 1363–1368. https://doi.org/10.1139/z86-203
  • Thorpe RS (1975) Quantitative handling of characters useful in snake systematics with particular reference to interspecific variation in the Ringed Snake Natrix natrix (L.). Biological Journal of the Linnean Society 7: 27–43. https://doi.org/10.1111/j.1095-8312.1975.tb00732.x
  • Thorpe RS (1983) A review of the numerical methods for recognizing and analyzing racial differentiation. In: Felsenstein, J. (Ed.) Numerical Taxonomy. NATO ASI Series (Series G: Ecological Sciences), Vol. 1. Springer-Verlag, Berlin, 404–423. https://doi.org/10.1007/978-3-642-69024-2_43
  • Trifinopoulos J, Nguyen L-T, von Haeseler A, Minh BQ (2016) W-IQ-TREE: a fast online phylogenetic tool for maximum likelihood analysis. Nucleic Acids Research 44: W232–W235. https://doi.org/10.1093/nar/gkw256
  • Turan C (1999) A note on the examination of morphometric differentiation among fish populations: The Truss System. Turkish Journal of Zoology 23: 259–263.
  • Wilcox TP, Zwickl DJ, Heath TA, Hillis DM (2002) Phylogenetic relationships of the Dwarf Boas and a comparison of Bayesian and bootstrap measures of phylogenetic support. Molecular Phylogenetics and Evolution 25: 361–371. https://doi.org/10.1016/S1055-7903(02)00244-0