Skip to main content
Advertisement
Browse Subject Areas
?

Click through the PLOS taxonomy to find articles in your field.

For more information about PLOS Subject Areas, click here.

  • Loading metrics

Phylogeographic Analyses Strongly Suggest Cryptic Speciation in the Giant Spiny Frog (Dicroglossidae: Paa spinosa) and Interspecies Hybridization in Paa

Abstract

Species identification is one of the most basic yet crucial issues in biology with potentially far-reaching implications for fields such as conservation, population ecology, and epidemiology. The widely distributed but threatened frog Paa spinosa has been speculated to represent a complex of multiple species. In this study, 254 individuals representing species of the genus Paa were investigated along the entire range of P. spinosa: 196 specimens of P. spinosa, 8 specimens of P. jiulongensis, 5 specimens of P. boulengeri, 20 specimens of P. exilispinosa, and 25 specimens of P. shini. Approximately 1333 bp of mtDNA sequence data (genes 12S rRNA and 16S rRNA) were used. Phylogenetic analyses were conducted using maximum parsimony, maximum likelihood and Bayesian inference. BEAST was used to estimate divergence dates of major clades. Results suggest that P. spinosa can be divided into three distinct major lineages. Each major lineage totally corresponds to geographical regions, revealing the presence of three candidate cryptic species. Isolation and differentiation among lineages are further supported by the great genetic distances between the lineages. The bifurcating phylogenetic pattern also suggests an east–west dispersal trend during historic cryptic speciation. Dating analysis estimates that P. spinosa from Western China split from the remaining P. spinosa populations in the Miocene and that P. spinosa from Eastern China diverged from Central China in the Pliocene. We also found that P. exilispinosa from Mainland China and Hong Kong might have a complex of multiple species. After identifying cryptic lineages, we then determine the discrepancy between the mtDNA and the morphotypes in several individuals. This discrepancy may have been caused by introgressive hybridization between P. spinosa and P. shini.

Introduction

Biologists usually categorize species by morphology as the original tool [1], and morphological traits remain immensely important for taxonomic description and identification of species in the field. However speciation is not always accompanied by a morphological change. In the last decade, cryptic species have been revealed in most types of organism from fungi, insects, and mammals to other types by using molecular data [2], [3]. The actual number of biological species is likely to exceed the current tally of nominal species, most of which are delineated on purely morphological grounds [4]. Identification and characterization of the world’s numerous cryptic species complexes have become an ongoing challenge for evolutionary biologists with potentially far-reaching implications for fields such as biodiversity conservation [4], epidemiology [5], [6], and biological control [7], [8].

Molecular phylogenies and DNA sequence data are important for identification and discovery of cryptic biological diversity that usually cannot be distinguished by traditional non-genetic approaches [4], [9]. Genetic markers provide critical information regarding evolutionary relationships even in instances of morphological homogeneity, and even few markers can clearly demarcate species boundaries when groups are strongly diverged. A number of studies have revealed cryptic diversity among species of anurans by molecular genetic approaches in the last decade [10], [11], [12].

Chinese montane forests are among the world’s major centers that harbor large numbers of species relative to their extremely mountainous terrain [13]. Frogs have species recognition and mate choice systems that rely on non-morphological characteristics, which can be storehouses of vast cryptic diversity [12]. With these factors considered, a taxonomic understanding of anurans that live in Chinese mountain area is essential for future conservation and management decisions. The giant spiny frog Paa spinosa is characterized by keratinized skin spines on the chest in males. The males also develop hypertrophied forearms during the breeding season [14]. These characteristics are viewed as adaptations to males grasping females during amplexus [14], [15]. P. spinosa is a widespread native species of amphibians in the southern and southeast part of China. This species inhabits rocky streams in evergreen forests and open countryside on hills and mountains from 500 m to 1500 m above sea level [16], [17]. Since ancient times the thinking of ‘food as medicine’ has existed in Chinese medical theories and Chinese food therapy [18], and this species are used to be considered have both dietary and medicinal functions [19]. The wild populations of this species have been significantly decreasing because of ecological and environmental degradation as well as consumption as food.

Specific habitat requirements suggest that P. spinosa is likely to be a species complex and high genetic divergences among populations because of isolated population and poor overland dispersers. Che et al. [20]suggested that P. spinosa is likely a species complex rather than a single species according to nuclear and mitochondrial DNA sequence data, but fewer specimens (only seven specimens from four localities) were used to infer the phylogenetic relationships. Elucidating various aspects of this diversity, such as cryptic species and degree of endemism, is critical in the understanding and formation of geographic models of speciation for these montane regions. However, population genetics and the evolutionary history of P. spinosa have not been systematically investigated. In the present study, a large sample of specimens collected across the whole species range of P. spinosa was investigated using two mtDNA genes (12S rRNA, 16S rRNA). This study aimed to (1) test whether cryptic species existed within P. spinosa, (2) examine population genetic structure and phylogeographic patterns, (3) present a hypothesis of population dynamics explaining the current distribution and phylogeographic pattern, and (4) provide genetic information to further develop effective conservation recommendations for this economically important but threatened animal.

Materials and Methods

Ethics Statement

All frogs were released to the wild at their site of capture after tissue sampling. Tissue samples were preserved in 95% ethanol before they were deposited at Zhejiang Normal University under voucher numbers identified by locality-haplotype numbers. Our experimental procedures complied with the current laws on animal welfare and research in China, and were specifically approved by the Animal Research Ethics Committee of Zhejiang Normal University (Permit No. AREC 20061024). The Provincial Forestry Departments of Anhui, Fujian, Guangdong, Guangxi, Hunan, Jiangxi, Yunnan and Zhejiang provided permits for capturing frogs in China.

Taxon Sampling and DNA Sequence Data

A total of 196 P. spinosa individuals were collected from 12 localities across South China (Table 1, Fig. 1). In addition, 8 individuals of P. jiulongensis, 5 individuals of P. boulengeri, 20 individuals of P. exilispinosa, and 25 individuals of P. shini were sampled and used as the ingroup. Moreover, the 12S rRNA and 16S rRNA genes previously reported as orthologous sequences for the Paa species were added to the dataset for subsequent phylogenetic analysis (Details are included in Table S1). Based on the work by Frost et al. [21], this study used the species Fejervarya limnocharis and F. cancrivora as the outgroups. Sequences of these two species were then obtained from GenBank.

thumbnail
Figure 1. Map of sampling localities for this study.

The hatched area represents the distribution of Paa spinosa. Twelve localities are grouped by city, with the following codes: DH, Dehua; GS, Longsheng; HS, Huangshan; JG, Jinggangshan; JH, Jinhua; PB, Pingbian; PJ, Pingjiang; LS, Lishui; LU, Lushan; WY, Wuyishan; YF, Yongfu; YS, Yangshan.

https://doi.org/10.1371/journal.pone.0070403.g001

thumbnail
Table 1. GPS coordinates and sample sizes for each geographic location.

https://doi.org/10.1371/journal.pone.0070403.t001

Genomic DNA was extracted from muscle and toe clips by standard 3-step phenol/chloroform extraction [22]. Polymerase chain reaction (PCR) was conducted to amplify fragments of the mitochondrial genes 12S rRNA and 16S rRNA in a 50 µl volume reaction under the following conditions: initial denaturation for 4 min at 95°C, 35 cycles of denaturation for 1 min at 95°C, annealing for 1 min at 58°C, extension for 1 min at 72°C; final extension for 10 min at 72°C. 16S rRNA gene was amplified with primers 16Sar 5′-CGCCTGTTTACCAAAAACAT-3′ and 16Sbr 5′-CCGGTYTGAACTCAGATCAYGT-3′ [23]; 12S rRNA gene was amplified with primers 12SF 5′- AGTGCTGAAAACGCTAAGAC -3′ and 12SR 5′- AGGGCGACGGGCGGTGTGTAC -3′ [24], [25]. PCR products were purified using the Wizard PCR Preps DNA Purification Kit (Promega). Each fragment was sequenced in both directions by using the PCR primers. To confirm the identity of the amplified targets, sequence were submitted to BLAST searches before deposition in GenBank.

Data Analysis

Sequences were aligned for each gene by using Clustal X 1.83 [26], and alignment was verified by eye with Bioedit 7.0.0 [27]. We calculated the number of polymorphic sites and haplotypes, as well as haplotype and nucleotide diversity for each taxon and each population by using DnaSP 4.0 [28].

Phylogenetic relationships among haplotypes were estimated using maximum parsimony (MP), maximum likelihood (ML), and Bayesian methods. Heuristic MP search with random-addition sequences and tree bisection and reconnection branch swapping in PAUP*4.0 [29] were conducted. The robustness of these analyses was assessed by non-parametric bootstrap replications using PAUP*4.0 with 100 replicates. ML analyses were performed by heuristic search methods conducted using PAUP*4.0 [29] with the general time-reversible (GTR+I+G) model of DNA evolution having the shape parameter of a gamma distribution (G). The GTR+I+G model was identified as the best-fitting model under the Akaike information criterion (AIC) as implemented in Modeltest 3.7 [30]; In the Bayesian analysis, we used MrBayes 3.1.2 [31] under the GTR+I+G model. BI used four simultaneous Metropolis-coupled Markov Chain Monte Carlo (MCMC) for 6,000,000 generations, we sampled a tree every 100 generations and calculated a consensus topology for 30,000 trees after discarding the first 30,000 trees (burn-in = 3,000,000). We regarded bootstrap values of 70% or greater for MP and ML, and Bayesian posterior probability of 0.95 as indicative of strong support [32], [33].

Species delimitation was guided by phylogenetic relationships based on molecular evidence [34], [35], [36]. The net between and ingroup mean distances between all major mtDNA clades were calculated following the Tamura and Nei (1993) distance module in MEGA 4.0 [37].

We used the BEAST version 1.6.1 [38] to construct an ultrametric tree and erect 95% confidence intervals for node heights. A molecular clock and previously published dates for Paa were used to estimate divergence times. We assumed a substitution rate ranging from 0.2% to 0.5% per Myr for the 12S and 16S genes, which are common in frogs [39], [40], [41], [42]. To estimate divergence times based on previously suggested dates, we calibrated our phylogeny by using one published divergence time to the most recent common ancestor (TMRCA) between the P. boulengeri and P. shini [43]. Fossil records and biogeographic inference were used to obtain the TMRCA of P. boulengeri and P. shini at about 27 Myr. Two MCMC runs were performed, with each running 20 million generations and logging parameters every 2000 generations. A burn-in of 10% was selected after likelihood score parameters were inspected in TRACER for stability.

Results

Sequence Characteristics

The 254 specimens were sequenced for 12S rRNA and 16S rRNA fragments,and new sequences were deposited in NCBI (Accession nos. JX989290-989797). Total alignment consisted of 1333 bp (12S = 805 bp; 16S = 528 bp) for the dataset in which 430 sites were variable and 372 sites were parsimony-informative. Sequences of the 12S rRNA fragment from previous studies [20], [21], [44] were shorter than those obtained in out study; thus, one region with 51 bp of 12S rRNA was excluded for further analysis. A total of 87 haplotypes were identified among the 254 ingroup sequences (Table S1).

Genetic Distances

Uncorrected distances for 12S rRNA to 16S rRNA between all clades ranged from 1.7% to 9.2% (Table 2). All haplotypes from P. spinosa were considered collectively. The mean genetic distance between P. spinosa and other species is 4.7%, and that between clades C (P. spinosa) and D (P. exilispinosa) is 1.7%. Notably, the genetic distances between the clades of P. spinosa ranged from 3% to 8.7%.

thumbnail
Table 2. Numbers below the diagonal are genetic diversity within and between clades from mtDNA data, assessed by K2P distances (median, %).

https://doi.org/10.1371/journal.pone.0070403.t002

Phylogenetic Relationships

MP, BI, and ML analyses generated almost identical tree topology and we present only the topology of ML in Figure 2. Phylogeny exhibits two quite well supported major lineages, and the clade composed of the subgenus Quasipaa (Clades A–G) is the sister taxon to the other species in Paa. Most relationships within these deep clades were resolved (Fig. 2).

thumbnail
Figure 2. Phylogenetic relationships of haplotypes, geographic distribution, and frequency of each Paa haplotype.

Phylogenetic relationships of haplotypes based on the 12S and 16S gene sequences, as determined by MP, ML, and Bayesian inference. Numbers above and below the branches are Bootstrap support values and Bayesian posterior probabilities ≥50 are indicated by the nodes. Legend: asterisk (“*”) indicates 100% ML and MP bootstrap support and 1.0 Bayesian posterior probabilities. “B,” haplotypes of P. boulengeri; “H,” haplotypes of P. spinosa; “J,” haplotypes of P. jiulongensis; “E,” haplotypes of P. exilispinosa; “S,” haplotypes of P. shini. P. spinosa YN and P. spinosa JX, sequences of P. spinosa obtained from GenBank; P. exilispinosa HK and P. exilispinosa FJ, sequence of P. exilispinosa obtained from GenBank; P. jiulongensis FJ, sequence of P. jiulongensis obtained from GenBank; P. shini GX, sequence of P. shini obtained from GenBank; P. taihangnica, P. liebigii, P. conaensis, P. medogensis, P. maculosa, P. arnoldi, P. chayuensis, P. robertingeri, P. yunanensis, P. bourreti, sequence of P. taihangnica, P. liebigii, P. conaensis, P. medogensis, P. maculosa, P. arnoldi, P. chayuensis, P. robertingeri, P. yunanensis, P. bourreti obtained from GenBank respectively.

https://doi.org/10.1371/journal.pone.0070403.g002

Seven major haplotype clades were revealed in our study with strongly support (Fig. 2). Clade A includes all haplotypes of P. boulengeri from Guangdong, Hunan, Jiangxi, and Hubei provinces. Clade D includes P. exilispinosa from Fujian and Guangxi provinces, except three haplotype of P. exilispinosa E6, E7, E8. The haplotype of P. exilispinosa from Hong Kong was diverged far from Clade D. Clade F includes all haplotypes of P. jiulongensis. Clade G includes all haplotypes of P. shini, except P. spinosa H44, which possessed an mtDNA haplotype of P. shini. Within the P. spinosa complex, all haplotypes were clustered into three major groups, each of which corresponded to a geographical region. Clade B was found to be P. spinosa limited to Yunnan Province (Western China) known as the sister taxon of P. boulengeri, which formed a basal branch with a strong support. Clade C includes P. spinosa from Jiangxi, Hunan, Guangdong, and Guangxi provinces as well as Dehua in Fujian Province (Center and Southern China). Clade E includes P. spinosa from Anhui, Zhejiang, and Jiangxi provinces as well as Mt. Wuyi in Fujian Province (Southeast China). Interestingly, the haplotype P. shini S6 was also found in Clade E (Fig. 2).

Diagnostic Sites

Thirty diagnostic sites which could clearly distinguish the three major clades of P. spinosa were found (Fig. 3). For example, the special base composition at sites 25, 92,115, 237, and so on could differentiate the alleles H31–33 from the clade B (western China clade). Diagnostic sites were also found in the haplotypes from clade C at site 265, 443, 548, 593, and the “T” at site 298 only occurred in the haplotypes from clade E (central China clade) (Fig. 3).

thumbnail
Figure 3. Variable sites of 12S and 16S mtDNA fragment in the three clades of P. spinosa, show the diagnostic sites.

Numbers above sequences correspond to base positions of bariable sites. Haplotype names and B, C, E within brackets correspond to those clade in Figure 2. A dot(.) indicates identity with the top sequence. Hyphens(−)indicate gaps. Asterisks(*) indicate diagnostic sites.

https://doi.org/10.1371/journal.pone.0070403.g003

Estimation of Divergence Times

Estimated divergence times among the three clades of P. spinosa are shown in Fig. 4. Our estimates indicated that Clade B (P. spinosa from western China) separated from the remaining P. spinosa populations at approximately 14 Ma. Clades C (P. spinosa from Central China), D (P. exilispinosa from Mainland China), and E (P. spinosa from Eastern China) diverged between 5 and 3 Ma. P. exilispinosa from Hong Kong diverged at 6 Ma.

thumbnail
Figure 4. Time tree for all species in our study.

Divergence time is shown on the node bars. Tree topology is derived from BEAST. Clades A, B, C, D, E, F, G correspond to those in Figure 2.

https://doi.org/10.1371/journal.pone.0070403.g004

Discussion

Cryptic Species Complex

Taxonomic decisions are usually based on recognizable morphological characteristics. However, different species cannot always be recognized by phenotypic changes [4]. In this study, we apply the mtDNA tree-based method as a phylogenetic criterion to propose species boundaries considering geographically separated mtDNA clades as different phylogenetic species. Our data revealed three divergent lineages of P. spinosa, which allows us to propose the existence of three candidate cryptic species within different localities of P. spinosa. Furthermore, each lineage corresponds to specific sampling regions (Figs. 1 and 2). Clade B mainly corresponds to southern China; Clade C mainly corresponds to Central China (Guangdong, Guangxi, Hunan, and Jiangxi provinces and one in Fujian Province); Clade E exclusively corresponds to Southeast China. A previous study showed widely distributed P. spinosa, and Fujian is designated as the place of origin [45]. This study suggests that a higher species diversity is observed in P. spinosa than in one widely distributed species in previous studies.

The isolation and differentiation among the clades are further supported by genetic distances. In this study, the genetic distances among the three lineages were found to range from 3.9% to 8.7%. The distance between Clade C (consisting of P. spinosa) and Clade D (consisting of P. exilispinosa) was markedly larger (Table 2). Vences and Wake [46] indicated that interspecific genetic distances ranged from 1% to 16.5% for 16S rRNA sequences among the species of Mantellidae from Madagascar. The interspecific distances between the recognized genus Engystomops also ranged from 2.9% to 4.1% for the 12S rRNA and 16S rRNA mtDNAs [47]. Notably, genetic distances in P. spinosa were greater than those in other phylogenetic and taxonomic studies of amphibians. Although species delimitation using mtDNA only is controversial, data collected from our research strongly suggest that cryptic species are hidden within P. spinosa. This result also agrees with the suggestion that seven specimens of P. spinosa from four localities do not belong to one lineage, based on nuclear and mtDNA sequence in the study by Che et al. [20]. Though relationship between some P. exilispinosa (E6, E7, E8) and Clade C is still have confusion in this study (Fig. 2) and the relationship between P. spinosa and P. exilispinosa are not resolved in the study of Che et al. [20]. The evidence from genetic data and phylogenetic analysis suggest that Clade E should be designated as the type material for P. spinosa in the results of our study.

Amphibians are considered philopatric, showing a high degree of population subdivision over relatively short geographical distances [48], [49], [50], [51], [52]. The most common mode of amphibian species formation is supposed to be allopatric speciation [46]. Che et al. [43] demonstrated that the genus Paa originated from Indochina, and rapid divergences accompanied by geographical and probably paleoclimatic change occurred. In the present study, a west–east dispersal trend of P. spinosa populations was revealed. The Western China lineage contains samples only distributed in Yunan Province. The Central China lineage contains samples distributed from Guangxi and Hunan provinces to Jiangxi and Fujian provinces. The Eastern China lineage ranges from Jiangxi Province to Zhejiang Province. This high correspondence is not difficult to interpret because these frogs live in isolated montane streams.

Evidence of Hybridization and Cryptic Lineages among Paa

Hybridization is prevalent in rapidly radiating groups, and amphibians and bird species remain adequately compatible to produce viable laboratory hybrids [53]. The clearest sign of introgression is the sympatric sharing of geographically localized mtDNA sequence haplotypes between otherwise genetically and morphologically divergent species [54]. In our study, 1 of 23 P. shini individuals (LUS1, Table S1) is identified as a hybrid that has the mtDNA haplotype of P. spinosa. Moreover, 1 P. spinosa (WY18) has an mtDNA haplotype of P. shini (Fig. 2). Despite the absence of hybridizable pairs of species known in Paini, intraspecific hybridization has occurred among different lineages of P. boulengeri because of secondary contact [55]. Considering that Paa is not an old genus composed through the Oligocene into the Miocene [43], the phylogenetic position of the haplotype S6 within Clade E and the haplotype H44 within Clade G most likely results from reciprocal introgression of mtDNA that occurred between these two amphibian species. In the present study, we also uncovered a population where haplotype Central China lineage and Eastern China lineage of P. spinosa coexisted in the Mt. Lushan area. Hence, Mt. Lushan is likely a rather narrow contact zone between populations of P. spinosa.

The haplotypes of P. exilispinosa HK from Hong Kong, which are sequences obtained from GenBank, are divided far from Mainland China. The separation of P. exilispinosa HK from mainland China started in the Pliocene. P. exilispinosa may also be a species complex rather than a single species; the formation of Hong Kong may explain this species complex. However no samples from Hong Kong was presented. Further systematic sampling from more geographic areas and more molecular markers should be used in future investigations to confirm the present hypothesis.

Implications for Conservation

A recent worldwide effort to evaluate the conservation status of amphibians indicated that more than 400 species (including P. spinosa) are rapidly declining and may be considered threatened in the near future [56], [57]. Failure to diagnose biological diversity can hamper conservation efforts and basic scientific inquiry [58]. Thus, precise and correct species delimitation is essential given that the species are the basic units of analysis in biogeography, ecology, macroevolution, biodiversity assessment, as well as conservation and management. Widespread high-density species are less likely to succumb to extinction compared with small isolated species [59]. Our study suggests that the threatened P. spinosa should at least be divided into three species, with each characterized by the following: a smaller population, more restricted geographic ranges, greater vulnerability to further population decline, and extinction.

Hybrid zones represent a potentially important source of information regarding taxon divergence and population interaction [60]. Protections of the remaining populations have to be prioritized considering the limited resources. Protection of contact zones such as Lushan and Dehua, which have greater gene diversity, should be prioritized considering that various anthropogenic pressures such as commercial overexploitation and habitat destruction have threatened this species complex in the past decades.

Supporting Information

Table S1.

Sampling information in Paa including specimens ID, localities, haplotypes and GenBank accession numbers.

https://doi.org/10.1371/journal.pone.0070403.s001

(DOC)

Acknowledgments

We are indebted to the following people who assisted us in collecting samples for this study: Chuntao Liu, Canyang Li, and Xinsheng Tang (Huangshan University). Thanks are also due to Xuming Zhou, Fang Yan, Yuanting Jin, Bing Shen, and Lian Chen for their assistance in the laboratory work.

Author Contributions

Conceived and designed the experiments: RZ JZ. Performed the experiments: HH SY. Analyzed the data: SX GY. Wrote the paper: HH SY.

References

  1. 1. Mayr E (1942) Systematics and the origin of species, from the viewpoint of a zoologist: Harvard University Press.
  2. 2. Hebert PDN, Penton EH, Burns JM, Janzen DH, Hallwachs W (2004) Ten species in one: DNA barcoding reveals cryptic species in the neotropical skipper butterfly Astraptes fulgerator. Proceedings of the National Academy of Sciences of the United States of America 101: 14812.
  3. 3. Murray TE, Fitzpatrick U, Brown MJF, Paxton RJ (2008) Cryptic species diversity in a widespread bumble bee complex revealed using mitochondrial DNA RFLPs. Conservation Genetics 9: 653–666.
  4. 4. Bickford D, Lohman DJ, Sodhi NS, Ng PKL, Meier R, et al. (2007) Cryptic species as a window on diversity and conservation. Trends in Ecology & Evolution 22: 148–155.
  5. 5. Collins F, Paskewitz S (1996) A review of the use of ribosomal DNA (rDNA) to differentiate among cryptic Anopheles species. Insect Molecular Biology 5: 1–9.
  6. 6. Mayr E (1970) Populations, species, and evolution: an abridgment of animal species and evolution: Belknap Press.
  7. 7. Heraty JM, Woolley JB, Hopper KR, Hawks DL, Kim JW, et al. (2007) Molecular phylogenetics and reproductive incompatibility in a complex of cryptic species of aphid parasitoids. Molecular Phylogenetics and Evolution 45: 480–493.
  8. 8. Hafez M, Doutt RL (1954) Biological evidence of sibling species in Aphytis maculicornis (Masi).(Hymenoptera, Aphelinidae). The Canadian Entomologist 86: 90–96.
  9. 9. Velo-Antón G, Burrowes P, Joglar R, Martínez-Solano I, Beard K, et al. (2007) Phylogenetic study of Eleutherodactylus coqui (Anura: Leptodactylidae) reveals deep genetic fragmentation in Puerto Rico and pinpoints origins of Hawaiian populations. Molecular Phylogenetics and Evolution 45: 716–728.
  10. 10. Fouquet A, Recoder R, Teixeira M, Cassimiro J, Amaro RC, et al. (2011) Molecular phylogeny and morphometric analyses reveal deep divergence between Amazonia and Atlantic Forest species of Dendrophryniscus. Molecular Phylogenetics and Evolution. 62: 826–838.
  11. 11. Evans BJ, Brown RM, McGuire JA, Supriatna J, Andayani N, et al. (2003) Phylogenetics of fanged frogs: testing biogeographical hypotheses at the interface of the Asian and Australian faunal zones. Systematic Biology 52: 794–819.
  12. 12. Stuart BL, Inger RF, Voris HK (2006) High level of cryptic species diversity revealed by sympatric lineages of Southeast Asian forest frogs. Biology Letters 2: 470–474.
  13. 13. Harkness J (1998) Recent trends in forestry and conservation of biodiversity in China. The China Quarterly 156: 911–934.
  14. 14. Liu C, Zheng R, Lu Q, Cheng B (2007) Sexual dimorphism and female reproductive characteristics in the giant spiny frog Paa spinosa. Herpetologica Sinica 11: 160–165.
  15. 15. Ohler A, Dubois A (2006) Phylogenetic relationships and generic taxonomy of the tribe Paini (Amphibia, Anura, Ranidae, Dicroglossinae), with diagnoses of two new genera. Zoosystema 28: 769–784.
  16. 16. Hu J, Xie F, Li C, Jiang J (2011) Elevational patterns of species richness, range and body size for spiny frogs. PLoS One 6: e19817.
  17. 17. Zhao E, Wang S, baohuju CGh, Commission ESS (1998) China red data book of endangered animals: Amphibia and reptilia. Amphibia and reptilia/chief compiler Zhao Ermi: Science Press.
  18. 18. Chen J (2000) Regularly control of functional food in China. Scandinavian Journal of Nutrition 22: 130–131.
  19. 19. Yu BG, Ye RH, Zheng RQ, Zhou Y, Liu CT, et al. (2008) The ethogram and activity rhythm of capative Paa spinosa during breeding period. Acta Ecologica Sinica 28: 6371–6378.
  20. 20. Che J, Hu J, Zhou W, Murphy RW, Papenfuss TJ, et al. (2009) Phylogeny of the Asian spiny frog tribe Paini (Family Dicroglossidae) sensu Dubois. Molecular Phylogenetics and Evolution 50: 59–73.
  21. 21. Frost DR, Grant T, Faivovich J, Bain RH, Haas A, et al.. (2006) The amphibian tree of life. Bulletin of the American Museum of natural History: 1–291.
  22. 22. Sambrook J, Russell DW (2001) Molecular cloning: a laboratory manual: CSHL press.
  23. 23. Bossuyt F, Milinkovitch MC (2000) Convergent adaptive radiations in Madagascan and Asian ranid frogs reveal covariation between larval and adult traits. Proceedings of the National Academy of Sciences 97: 6585.
  24. 24. Zhou Y, Zhang JY, Zheng RQ, Yu BG, Yang G (2009) Complete nucleotide sequence and gene organization of the mitochondrial genome of Paa spinosa (Anura: Ranoidae). Gene 447: 86–96.
  25. 25. Kocher TD, Thomas WK, Meyer A, Edwards SV, Pääbo S, et al. (1989) Dynamics of mitochondrial DNA evolution in animals: amplification and sequencing with conserved primers. Proceedings of the National Academy of Sciences 86: 6196–6200.
  26. 26. Thompson JD, Gibson TJ, Plewniak F, Jeanmougin F, Higgins DG (1997) The Clustal X windows interface: flexible strategies for multiple sequence alignment aided by quality analysis tools. Nucleic acids research 25: 4876–4882.
  27. 27. Hall TA. BioEdit: a user-friendly biological sequence alignment editor and analysis program for Windows 95/98/NT 1999: 95–98.
  28. 28. Rozas J, Sánchez-DelBarrio JC, Messeguer X, Rozas R (2003) DnaSP, DNA polymorphism analyses by the coalescent and other methods. Bioinformatics 19: 2496–2497.
  29. 29. Swofford DL (2003) PAUP*: phylogenetic analysis using parsimony, version 4.0 b10.
  30. 30. Posada D, Crandall KA (1998) Modeltest: testing the model of DNA substitution. Bioinformatics 14: 817–818.
  31. 31. Ronquist F, Huelsenbeck JP (2003) MrBayes 3: Bayesian phylogenetic inference under mixed models. Bioinformatics 19: 1572–1574.
  32. 32. Huelsenbeck JP, Hillis DM (1993) Success of phylogenetic methods in the four-taxon case. Systematic Biology 42: 247–264.
  33. 33. Hillis DM, Bull JJ (1993) An empirical test of bootstrapping as a method for assessing confidence in phylogenetic analysis. Systematic Biology 42: 182–192.
  34. 34. Simpson GG (1961) Principles of animal taxonomy: Columbia University Press.
  35. 35. Wiley EO (1978) The evolutionary species concept reconsidered. Systematic Biology 27: 17–26.
  36. 36. De Queiroz K (2007) Species concepts and species delimitation. Systematic Biology 56: 879–886.
  37. 37. Tamura K, Dudley J, Nei M, Kumar S (2007) MEGA4: molecular evolutionary genetics analysis (MEGA) software version 4.0. Molecular biology and evolution 24: 1596–1599.
  38. 38. Drummond AJ, Rambaut A (2007) BEAST: Bayesian evolutionary analysis by sampling trees. BMC evolutionary biology 7: 214.
  39. 39. Evans BJ, Kelley DB, Tinsley RC, Melnick DJ, Cannatella DC (2004) A mitochondrial DNA phylogeny of African clawed frogs: phylogeography and implications for polyploid evolution. Molecular Phylogenetics and Evolution 33: 197–213.
  40. 40. Lemmon EM, Lemmon AR, Cannatella DC (2007) Geological and climatic forces driving speciation in the continentally distributed trilling chorus frogs (Pseudacris). Evolution 61: 2086–2103.
  41. 41. Wiens JJ, Sukumaran J, Pyron RA, Brown RM (2009) Evolutionary and biogeographic origins of high tropical diversity in Old World frogs (Ranidae). Evolution 63: 1217–1231.
  42. 42. Gvoždík V, Moravec J, Klütsch C, Kotlík P (2010) Phylogeography of the Middle Eastern tree frogs (Hyla, Hylidae, Amphibia) as inferred from nuclear and mitochondrial DNA variation, with a description of a new species. Molecular Phylogenetics and Evolution 55: 1146–1166.
  43. 43. Che J, Zhou WW, Hu JS, Yan F, Papenfuss TJ, et al. (2010) Spiny frogs (Paini) illuminate the history of the Himalayan region and Southeast Asia. Proceedings of the National Academy of Sciences 107: 13765–13770.
  44. 44. Jiang J, Dubois A, Ohler A, Tillier A, Chen X, et al. (2005) Phylogenetic relationships of the tribe Paini (Amphibia, Anura, Ranidae) based on partial sequences of mitochondrial 12S and 16S rRNA genes. Zoological science 22: 353–362.
  45. 45. Dubois A (1975) Un nouveau sous-genre (Paa) et trois nouvelles espèces du genre Rana. Remarques sur la phylogénie des Ranidés (Amphibiens, Anoures). Bull Mus Natn Hist nat, Paris (3) 324: 1093–1115.
  46. 46. Vences M, Wake D (2007) Speciation, species boundaries and phylogeography of amphibians. Amphibian biology 7: 2613–2671.
  47. 47. Ron SR, Santos JC, Cannatella DC (2006) Phylogeny of the túngara frog genus Engystomops (Physalaemus pustulosus species group; Anura: Leptodactylidae). Molecular Phylogenetics and Evolution 39: 392–403.
  48. 48. Funk WC, Blouin MS, Corn PS, Maxell BA, Pilliod DS, et al. (2005) Population structure of Columbia spotted frogs (Rana luteiventris) is strongly affected by the landscape. Molecular Ecology 14: 483–496.
  49. 49. Hitchings SP, Beebee TJC (1997) Genetic substructuring as a result of barriers to gene flow in urban Rana temporaria (common frog) populations: implications for biodiversity conservation. Heredity 79: 117–127.
  50. 50. Kraaijeveld-Smit F, Beebee T, Griffiths R, Moore R, Schley L (2005) Low gene flow but high genetic diversity in the threatened Mallorcan midwife toad Alytes muletensis. Molecular Ecology 14: 3307–3315.
  51. 51. Palo JU, Schmeller DS, Laurila A, Primmer CR, Kuzmin SL, et al. (2004) High degree of population subdivision in a widespread amphibian. Molecular Ecology 13: 2631–2644.
  52. 52. Veith M, Vences M, Vieites DR, Nieto-Roman S, Palanca A (2002) Genetic differentiation and population structure within Spanish common frogs (Rana temporaria complex; Ranidae, Amphibia). Folia Zool 51: 307–318.
  53. 53. Mallet J (2005) Hybridization as an invasion of the genome. Trends in Ecology & Evolution 20: 229–237.
  54. 54. Funk DJ, Omland KE (2003) Species-level paraphyly and polyphyly: frequency, causes, and consequences, with insights from animal mitochondrial DNA. Annual Review of Ecology, Evolution, and Systematics: 397–423.
  55. 55. Yan F, Zhou W, Zhao H, Yuan Z, Wang Y, et al. (2012) Geological events play a larger role than Pleistocene climatic fluctuations in driving the genetic structure of Quasipaa boulengeri (Anura: Dicroglossidae). Molecular Ecology. 22: 1120–1133.
  56. 56. Stuart SN, Chanson JS, Cox NA, Young BE, Rodrigues ASL, et al. (2004) Status and trends of amphibian declines and extinctions worldwide. Science 306: 1783–1786.
  57. 57. Xie F, Lau MWN, Stuart SN, Chanson JS, Cox NA, et al. (2007) Conservation needs of amphibians in China: A review. Science in China Series C: Life Sciences 50: 265–276.
  58. 58. Mayden RL, Wood RM. Systematics, species concepts, and the evolutionarily significant unit in biodiversity and conservation biology; 1995.
  59. 59. MacArthur RH, Wilson EO (1963) An equilibrium theory of insular zoogeography. Evolution: 373–387.
  60. 60. Haldane JBS (1948) The theory of a cline. Journal of Genetics 48: 277–284.