Shell shape and genetic variability of Southeast Asian Box Turtles ( Cuora amboinensis ) from Borneo and Sumatra

Distinguishing between species is an essential aspect of animal research and conservation. For turtles, morphology and genetic analysis are potentially valuable tools for identification. Shell shape is an important component of phenotypic variation in turtles and can be easily described and quantified by geometric morphometrics (GM). Here, we focus on carapace and plastron shape discrimination of immature Southeast Asian box turtles (Cuora amboinensis) from two of the Greater Sunda Islands with partially distinct faunas. GM analysis identified significant differences in carapace and plastron shape between turtles from Borneo and Sumatra. The discrimination success amounted to 90% and 83.7% for carapace and plastron, respectively. The correlations of carapace and plastron shapes were high for Sumatra (0.846), and less pronounced for Borneo (0.560). We detected no differences in the ontogenetic trajectories of the shell shape between the two islands. We conclude that shell shape can be used for reliable geographic assignment of C. amboinensis of unknown origin. In addition to the comparison of shell shapes, turtles from Borneo, Sumatra, Seram, and turtles of unknown origin from two Czech zoos were studied genetically. Analysis of the complete mitochondrial cytochrome b gene confirmed the distinctness of turtles from Borneo and Sumatra, with p-distance 2.68 – 4.09% sequence difference. Moreover, we discovered considerable genetic difference in Seram turtles of previously unknown haplogroup (p-distance 6.00 – 8.68%) revealing the need for the revision of the whole species complex of Cuora amboinensis.


Introduction
The Southeast Asian box turtle Cuora amboinensis (Riche in Daudin, 1801), belongs to the most diversified and widespread taxon of the genus Cuora with a distribution range including a major part of Southeast Asia (Iverson 1992).Unfortunately, it is also the most abundant hard-shelled turtle in Chinese markets and frequently used in traditional Chinese medicine (Cheung & DuDgeon 2006;Chen et al. 2009).Thus, it is ex-ploited in huge numbers, especially from Indonesia and Malaysia, despite export quotas and even a total export ban in some regions.As a result, its numbers are rapidly declining and some populations are already extinct (Ives et al. 2008;sChoppe 2008sChoppe , 2009)).Cuora amboinensis is listed in Appendix II of CITES and globally red-listed as 'Vulnerable' (IUCN 2013).In the face of the current Asian Turtle Crisis (Cheung & DuDgeon 2006) and the unsuccessful protection of the species in the wild, ex situ captive breeding programs inside and outside its distribution range are gaining increasing importance.For such captive breeding efforts, the correct identification of subspecies and the geographic provenance of turtles is of paramount importance.Although adults of currently recognized subspecies can be distinguished by standard morphometrics and coloration (rummler & FrItz 1991;mCCorD & phIlIppen 1998), the accurate determination of juveniles still poses serious problems.
Currently, there are four subspecies classified according to morphology and coloration (rhoDIn et al. 2010).Turtles from Sumatra and Java are considered to belong to C. a. couro (Schweigger, 1812).The subspecies C. a. kamaroma (Rummler and Fritz 1991) occurs in Borneo, the Malay Peninsula, Cambodia, Laos, Thailand, and Vietnam.Turtles from Myanmar are identified as C. a. lineata (McCord and Philippen, 1998) 2006;le et al. 2007;prasChag et al. 2006prasChag et al. , 2007;;FrItz et al. 2008;tIeDemann et al. 2014), and several species and subspecies have been described as new to science or resurrected from synonymy (e.g.BlanCk et al. 2006, prasChag et al. 2007, 2009;spInks et al. 2012;Ihlow et al. 2016).A comprehensive genetic study is still lacking for the wide-ranging and polytypic C. amboinensis.The situation is further complicated by the frequent hybridization of geoemydid turtles (wInk et al. 2001;BuskIrk et al. 2005;stuart & parham 2006;Fong et al. 2007;shI et al. 2008 Measuring external morphology using geometric morphometrics (GM) is a practical tool for assessing phenotypic variation of shell shape.This approach is easily applied and yields immediate results, independent from any laboratory work, thus making it highly suitable for taxonomic determination in the field (zelDItCh et al. 2004).
We therefore used GM to analyse the shell shapes of immature C. amboinensis box turtles from Borneo and Sumatra.In addition, we used the mitochondrial cytochrome b gene to genetically investigate the turtles from these islands, and specimens from other locations, in order to gather more information about these species and to compare the morphological results with the genetic findings.

Geometric morphometrics
A total of 195 photographs of C. amboinensis were examined (69 turtles from Borneo and 126 from Sumatra) and 132 (69 Borneo, 63 Sumatra) were chosen for further study.These included only immature individuals of unknown sex, with carapace lengths between 70 and 120 mm.Photographs of carapaces or plastra with abnormalities were discarded as well as photos of closed plastra to avoid perspective bias leaving 130 carapaces (68 Borneo, 62 Sumatra) and 98 plastra (69 Borneo, 29 Sumatra) for analysis.
For each turtle, standard dimensions of the shell (carapace length, carapace width, plastron length, plastron width) were measured using a calliper (0.1 mm precision).The digital images of carapace and plastron of each individual were obtained using a digital camera (Canon EOS 30D with Canon 50/1.8 lens) mounted on a tripod.Twenty-one anatomical landmarks of type 1 on plastron and twenty-five landmarks of type 1 and one of type 3 on carapace following the classification of BooksteIn (1997) were recorded (Fig. 1.) using TPSdig software (rohlF 2008).Each set was then symmetrised and one half was removed using the BigFix6 program (sheets 2003).Statistical examination was performed on half of the landmark sets.We employed the Procrustes superimposition method (zelDItCh et al. 2004) using the CoordGen6 program (sheets 2003) to remove the effects of position, orientation and scale, employing sets of x, y coordinates of landmarks from each specimen.We used the standardization on mean carapace length (for each population separately) to remove the size related shell shape differences in the program Standard6 (sheets 2003).Visualization was performed using CVAGen6 software (sheets 2003).The vectors of the shell shape ontogeny between turtles from Borneo and Sumatra were compared to the variability of the ontogeny vector inside these two samples using the VecCompare6 program (sheets 2003) and 400 permutations.When the vector between the samples is bigger than the 95 th percentile of the ranges of within-sample angles, we can assume that it is not expected that the samples significantly differ in the shell shape vector of the ontogeny randomly.The correlation between carapace and plastron shape was examined using PLSMaker6 software (sheets 2003).The partial warp scores for the further statistical analysis were generated using PCAGen6 software (sheets 2003).The differences in shell shape between turtles from Borneo and Sumatra were tested in the program Statistica 6 (weIss 2007) using Discriminant Analysis.

DNA samples and mitochondrial DNA sequencing
Nine turtles of known geographical provenance (three from Borneo, four from Sumatra, and two from Seram) were studied genetically.Additionally, nine individuals of unknown geographical provenance from zoological gardens (six samples from Zoological Garden Prague, Czech Republic, three samples from Zoological Garden Ústí nad Labem, Czech Republic) were included in this analysis (Table 1).
For each turtle a claw tip was removed and stored in an Eppendorf tube with 96% ethanol prior to DNA extraction.Total genomic DNA was then isolated using the DNAeasy Tissue Kit (Qiagen, Hilden, Germany) following the manufacturer's guidelines.
The DNA amplification was performed with the primers suggested by spInks et al. (2004) for a total length of

Phylogenetic analyses
We used our sequence data to construct a bootstrapped maximum likelihood (ML) tree using RAxML software (version 7.2.8-alpha)(stamatakIs 2006).The relationship between the subspecies was examined using 1,000 bootstrap replicates and the GTRGAMMA model.The average distances between haplotypes from particular groups were calculated in the mega 7.0.18program (kumar et al. 2016) using uncorrected p-distances model.
The correlation for carapace length and centroid size (geometric size) differed between the samples from Sumatra and Borneo (Fig. 4).However, we did not find any significant differences with respect to the growth vectors of carapacial (angle between populations 40, angle within Borneo 43.3, angle within Sumatra 26.6) and plastral shape (angle between populations 32, angle within Borneo 45.3, angle within Sumatra 35.3) between the samples from these two islands.
The analysis revealed a weak correlation (0.5597) between carapacial and plastral shape for turtles from Borneo, and a much stronger correlation (0.8458) for turtles from Sumatra.

DNA analysis
We sequenced the mitochondrial cytochrome b gene (1140 bp) in 18 individuals and found 14 distinct haplotypes.ML revealed two clearly distinct groups.The first group contained C. amboinensis haplotypes from Seram which were deeply divergent from all other haplotypes of C. amboinensis (uncorrected p distances 6.00 -8.68%).Among the remaining haplotypes, uncorrected p distances ranging from 0.00% to 5.36% were observed.The phylogenetic analyses placed the Seram haplotypes as a sister group of a monophyletic group including the remaining C. amboinensis sequences (Fig. 5).The latter group exhibits a clear structure, with sequences from Borneo and Sumatra in distinct parts of the tree.Uncorrected p-distances between haplotypes belonging to the Borneo and Sumatra groups varied within the range of 2.68 -4.09% (Table 2.).The haplotypes from the same island were similar (uncorrected p distances: 0.51-1.53%for Borneo and 0.13 -0.26% for Sumatra) and formed monophyletic groups (Borneo: bootstrap support 87) and, for the sequences from Sumatra (bootstrap support 98), contained additional sequences of unknown geographic origin.The sister relationship between the "Sumatra clade" and the "132 Zoo Ústí nad Labem and 51 Zoo Prague clade" is moderately supported (bootstrap support 70) with uncorrected p-distances varying within the range of 0.64 -1.15%.The relative position of additional sequences from GenBank of C. amboinensis, C. kamaroma, C. cuoro and C. lineata in the tree was not resolved because the phylogenetic relationships within this clade were poorly supported.

Discussion
Shell shape variation shows a clear distinctiveness between Borneo and Sumatra populations of C. amboin ensis, which corresponds with our DNA analyses.Geometric morphometrics therefore provides sufficient   geographical specific population recognition according to shell shape.Thus, it is a useful method for accurate turtle classification even in early ontogenetic stages.Similar studies of the other subspecies, populations and ontogenetic stages are needed to map all phenotypic variability in this species.GM can then serve as a simple yet effective practical tool for the determination of specimens of unknown geographical origin, with particular application in fieldwork.
From our genetic analyses, we discovered that the mitochondrial cytochrome b sequences of Cuora amboin ensis turtles from Borneo and Sumatra are clearly distinct, which corresponds with the currently accepted taxonomy of this species (rummler & FrItz 1991;rhoDIn et al. 2010), with Bornean turtles being identified with another subspecies (C.a. kamaroma) rather than with Sumatran turtles (C.a. couro).Yet, in our phylogenetic analyses, previously published sequences identified with these subspecies (spInks 2004) were distinct from ours, perhaps due to misidentification or geographic variation of these taxa.
A notable finding of our study was that the haplotypes from Seram, geographically corresponding to the subspecies of C. a. amboinensis, are highly divergent.However, it was not possible to prove the congruence between morphology and genetic divergences due to the low number of individuals (n=2) for morphometric analysis.We are therefore not yet able to draw taxonomic conclusions for turtles from the Seram provenience.More detailed DNA analyses, including the use of nuclear markers, are needed to clarify the affinities of this population.Moreover, Cuora amboinensis could in fact be a species complex with similar patterns as in other turtle species with vast geographical distribution which have been studied recently (FrItz et al. 2012;kInDler et al. 2012;kInDler et al. 2016;vargas-ramírez et al. 2010;vargas-ramírez et al. 2013;petzolD et al. 2014;eDwarDs et al. 2016).
The shell-shape and genetic differences between Borneo and Sumatra samples are in compliance with the general differences and separation of the majority of these island faunas.The isolation of the Greater Sunda Islands in the Tertiary caused deep genetic divergences between the reptilian fauna (keogh et al. 2001) and even in mammals (thInh et al. 2010;nater et al. 2011) inhabiting Borneo, Sumatra and the mainland.The similar pattern of distribution in C. amboinensis could be the result of a low dispersion ability due to the semiaquatic lifestyle in marshes, ponds and small streams, as seen in the ecologically close genus Cyclemys (FrItz et al. 2008)   tains in Borneo may act as an effective natural barrier between the suitable lowland habitats of C. amboinensis.A low correlation of the carapace and plastron shape in the Borneo turtles could be a consequence of geographically isolated populations.A similar pattern was observed in the genetic variation found in the Mediterranean turtle (Mauremys leprosa) (FrItz et al. 2006).
Our results revealing distinct morphological and genetic differences between Sumatra and Borneo box turtles are congruent with the recent findings of ernst et al. (2016), who performed classical morphometric analysis on C. amboinensis throughout its distribution range.Their research supported the validity of only two subspecies, amboinensis and kamaroma.According to these authors, amboinensis occupy a range that includes Sumatra, while kamaroma turtles inhabit Borneo.Placement of Borneo populations into continental subspecies C. a. kamaroma sensu lato is confirmed by the clustering of the examined continental haplotype into the Bornean cluster.Nevertheless, the positions of the additional haplotypes from Genbank in our phylogenetic tree suggest that preliminary taxonomic conclusions made by ernst et al. (2016) need to be confirmed by genetic analyses covering the whole range of the species.In light of our results, the use of the name C. a. amboinensis also for Sumatra populations, formerly classified as C. a. couro, is especially problematic.Although we analysed just two samples from the Molucca archipelago, where the typical habitat of C. amboinensis is found, these haplotypes are strongly different not only from those collected in Sumatra, but also from all other examined ones.
In conclusion, the high identification success of immature C. amboinensis specimens based on the phenotypic variation of shell shape using GM clearly demonstrates the usefulness of this method.It could be used on other forms of this species as a practical and effective tool for the determination of C. amboinensis of unknown origin thus contributing to the conservation strategies for this taxon as well as benefitting general scientific research.
From our genetic analyses we have discovered not only that the congruence between morphology and genetic distinctness for Borneo and Sumatra box turtles support deep divergent lineages, but moreover, we uncovered previously unknown haplotypes from Seram suggesting that the species status of this population should be reconsidered.

Fig. 1 .
Fig. 1.Carapace and plastron of Cuora amboinensis showing the landmarks used in this study.

Fig. 2 .
Fig. 2. Canonical plot for carapacial shape.Each carapace is represented by a symbol: × -Sumatra, • -Borneo (A).The first axis accounts for 79.05% of the total between group variation.The thin plate spline diagram shows a change in carapacial shape along the first axis (in direction of arrows) (B).

Fig. 3 .
Fig. 3. Canonical plot for plastral shape.Each plastron is represented by a symbol: × -Sumatra, • -Borneo (A).The first axis accounts for 55.25% of the total between group variation.The thin plate spline diagram shows a change in plastral shape along the first axis (in direction of arrows) (B).
(vorIs 2000)r Cuora species(spInks & shaFFer 2006)with limited distribution ranges.In contrast, a big river species like the Malaysian giant turtle (Orlitia borneensis) with a similar distribution range as C. amboinensis could have benefited from geomorphologic conditions during the glacial periods when the Malay Peninsula, Sumatra and Borneo were connected by a system of rivers i.e.the Siam and West Sunda River drainages(vorIs 2000).It allowed continuous gene flow between the now isolated areas (palupCIkova et al. 2012).The central location of moun-

Table 2 .
Uncorrected p distances (percentages) based on the mitochondrial cytochrome b gene (1140 bp).