Parapatric speciation found in Limnonectes (Dicroglossidae) species in Thailand
Received: October 14, 2017
Accepted: November 08, 2017
Published: January 02, 2018
Genet.Mol.Res. 17(1): gmr16039856
Biodiversity stems from speciation, and species identification allows scientists to measure species biodiversity. However, some closely related species, particularly amphibians that choose to inhabit the same environment, are difficult to identify. Juvenile specimens of Limnonectes gyldenstolpei (Lg) and L. taylori (Lt) were used as models to study the speciation process in populations with a parapatric distribution. Samples of the two species were collected from six provinces in the north and northeast of Thailand. Inter-simple sequence repeat (ISSR) data of five loci showed 96 polymorphic bands and species-specific bands of the primers UBC811 (Lg) and UBC824 (Lt). There was limited genetic variation within the same population, 8% in Lg and 2% in Lt, whereas a high genetic variation was observed between populations from different provinces 90% and 97% in Lg and Lt, respectively. An unweighted pair-group method using arithmetic averages (UPGMA) tree was constructed for each population using Jaccard’s similarity efficient. The tree demonstrated genetic differentiation between the two species. Principal coordinates analysis also supported the UPGMA diagram by showing no clear clusters with close genetic distances. The results of the molecular variance analysis within the population were high, indicating that they were genetically similar, and that microhabitat is an important factor influencing genetic distribution. Overall, the present study proves the hypothesis that the distribution of the two species is in a Parapatric form and that the ISSR technique is inappropriate tool to study the differentiation of a species complex.
Speciation is an important cause of diversity (Schluter and Pennell, 2017). Speciationmay occur either (1) among species that coexist in the same microhabitat, which includes sympatric and parapatric speciation, or (2) among species with a mutually exclusiv geographic range or allopatric speciation (Cain, 1953; Mayr, 1947; Mayr et al., 1963; Rivas, 1964). Morphologically similar species are difficult to identify, particularly in their larval and juvenile stages. In the present study, we found two similar-looking juveniles of the genus Limnonectes inhabiting the same microhabitats in northern and north-eastern Thailand, which then develop into adults that are easily distinguishable. This makes it interesting to study the genetic variability between L. gyldenstolpei (Lg) and L. taylori (Lt) found in the same areas. It is possible that the genes of these two species are similar.
Dicroglossine frogs of the genus Limnonectes are found across most Asian and Southeast Asian countries. Their dominant characteristic is the odontoid processes called fangs on their lower jaw (Emerson et al., 2000; Inger and Stuart, 2010). L. kuhlii Tschudi (1838) is a complex of cryptic species (McLeod, 2010a). In 2008, a big-headed species of Limnonectes, L. megastomias, phenotypically similar to L. gyldenstolpei in the L. kuhlii complex, was collected from eastern Thailand by David S. McLeod. Two species of the L. kuhlii complex from Thailand were described by Matsui et al. (2010) based on molecular and morphological evidence: L. taylori, a northern Thai lineage, identified as Lineage 12 by McLeod (2010b), and L. jarujini, a southern Thai lineage discussed by Taylor (1962). Amphibians are a critical component of certain ecosystems. In addition, they play a role in the proper functioning of communities that live in their vicinity and are dependent on other animals, for instance, for food. The loss of amphibians from an ecosystem can change primary production, affecting food webs from aquatic insects up to riparian predators (Hocking and Babbitt, 2014; Meredith et al., 2016; Whiles et al., 2006). The permeable skin % of frogs is sensitive to environmental changes, and frogs can act as bio-indicators (Campos et al., 2014; Hopkins, 2007; Lebboroni et al., 2006). Amphibian populations are declining, and some are on the edge of extinction. Small changes in the environment can cause speciation that can be traced back by means of genetic studies. To study genetic variability, many molecular tools are available, including a Technique based on inter-simple sequence repeats (ISSRs). ISSRs are regions in the genome located between microsatellite sequences (Zietkiewicz et al., 1994). To use ISSR markers, prior knowledge of the ISSR target sequences is not required; moreover, the markers are highly reproducible due to their primer length and the high stringency achieved by the annealing temperature, and they have been found to provide highly polymorphic fingerprints (Bornet and Branchard, 2001; Kojima et al., 1998; Zietkiewicz et al., 1994). Here, we examined parapatric speciation in two frog species of the L. kuhlii complex. We randomly sampled the areas where records of the species were found. To study their genetic variability, we then subjected juvenile specimens of L. gyldenstolpei and L. taylori to the ISSR technique.
Material and Methods
Sampling, DNA extraction, and polymerase chain reaction (PCR) amplification using an ISSR marker in total, 64 L. taylori and 59 L. gyldenstolpei samples were collected from six provinces in the north and northeast of Thailand (Table 1 and Figure 1).
|L. gyldenstopei||Mae Wang NP., Chaing Mai||MWN||4||678||18°39¢26.3¢¢||98°40¢54.6¢¢|
|Ban Pa Sak Ngam, Chaing Mai||BPS||12||402||18°59¢16.0¢¢||99°06¢52.4¢¢|
|Khun Jae NP., Chiang Rai||KJN||9||592||19°09¢10.6¢¢||99°23¢35.5¢¢|
|Khun Yuam, Mae Hong Son||KY||5||673||18°48¢29.4¢¢||97°52¢16.5¢¢|
|Doi Phu Fa, Nan||DPF||3||665||18°58¢51.5¢¢||101°10¢57.3¢¢|
|Phu Rua NP., Loei||PRN||5||939||17°28¢55.0¢¢||101°20¢59.4¢¢|
|Phu Luang, Loei||PL||13||865||17°13¢35.0¢¢||101°30¢00.8¢¢|
|Sakaerat, Nakhon Ratchasima||SKR||12||865||14°29¢36.1¢¢||101°52¢25.1¢¢|
|L. taylori||Pang Num Tune, Chaing Mai||PNT||13||968||19°04¢09.9¢¢||99°21¢34.1¢¢|
|Pang Bong, Chaing Mai||PB||2||888||18°58¢37.1¢¢||99°19¢27.1¢¢|
|Pang, Kampheaeng Hin,||PKH||4||1067||18°57¢07.3¢¢||99°20¢50.9¢¢|
|Chaing Mai, Doi Lang, Chaing Mai||DL||11||1950||20°05¢59.5¢¢||99°16¢24.2¢¢|
|Mae Wang NP, Chaing Mai||MWN||4||678||18°39¢26.3¢¢||98°40¢54.6¢¢|
Table 1: Population locations and sizes of L. taylori and L. gyldenstolpei
A GPS location for each specimen was recorded. Fragments of the liver, muscle, or toe were preserved in 95% ethanol for DNA extraction, which was performed with a Genomic DNA Extraction Kit (RBC, Taiwan) following the manufacturer’s instructions.
Spectrophotometry (BioDrop DUO, UK) and was used to determine the concentration of the DNA obtained, and its quality was assessed using electrophoresis in agarose/TAE (1.5%) gel. A total of 20 primers of the UBC primer set (University of British Columbia) were tested for the ISSR analysis (Table 2). PCR amplifications were performed in a 15 μl reaction volume containing 10 ng of template DNA, 0.8 μM primer, 0.2 mM dNTP, 1x buffer containing 2 μM Mg2+, and 1U Taq polymerase. The cycling condition was as follows: denaturation step at 94°C for 5 min, 35 cycles of denaturation at 94°C for 1 min, annealing temperature (Table 2) for 1 min, and an extension at 72°C for 1.3 min, followed by a final extension at 72°C for 5 min (Machkour-M'Rabet et al., 2009).
|Primer name||Primer sequence (5’-3’)||Annealing temperature (°C)||Polymorphic band||Size range (bp)|
Table 2: Primers and their annealing temperature, polymorphic band and size range
After amplification, samples comprising 7 μL of the PCR reaction mixture were subjected to electrophoresis on 1.5% agarose gels and 3× GelRed (Biotium, USA). Electrophoresis was performed in 1× TAE buffer at 130 V for 34 min. The bands were visualized under ultraviolet light, and the gel was photographed for analysis.
The fragments of ISSR were scored as binary characters: present (1) or absent (0) for each individual and each primer. To describe the genetic structure and variability among and between populations, non-parametric analysis of molecular variance (AMOVA) (Excoffier et al., 1992) was performed and a scatter plot of principal coordinates was constructed using the program principal coordinates analysis (PCoA) and GENALEX version 6.502. Nei’s genetic distances were used to perform an unweighted pair group method using arithmetic averages (UPGMA) cluster analysis, and phylogenetic trees were constructed using the POPTREE2 program (Takezaki et al., 2009). The relationship between genetic isolation and geographic1distance was analyzed using Isolation by Distance (IBD) Web Service version 3.23.
ISSR PCR is a technique that uses microsatellite sequences and randomly amplifies them in a similar manner as in RAPD. Thus, ISSR PCR is sensitive to reaction parameters. A total of 20 ISSR primers were screened with the 123 samples as mentioned above. After comparing the effects of the template DNA concentration and temperature during the annealing stage of amplification, only five primers were able to generate clear and reproducible DNA fragments (Table 3). These were then selected for further analysis. A total of 96 alleles were revealed by five ISSR loci across 123 individuals with a size ranging between 100–2,000 base pairs (Figure 2). The primers UBC845 and UBC811 showed the highest and lowest number of ISSR bands, respectively. We found UBC811 to be species-specific with the size of the amplified DNA fragment being 300 base pairs in Lg and 1500 base pairs in Lt.
|Source of variation||L. gyldenstopei||L. taylori|
|d.f.||Sum of square||Percentage variation||d.f.||Sum of square||Percentage variation|
|d.f.=Degree of freedom|
Table 3: AMOVA results of L. gyldenstopei and L. taylori based on ISSR
The dendrogram obtained from the UPGMA cluster analysis of the ISSR marker is shown in Figure 3. The UPGMA dendrogram revealed a close association between L. gyldenstolpei and L. taylori. According to the tree, L. gyldenstolpei divided from L. taylori at IV. The similarity index was 0.385–1.260 (Table 4) and 0.000–2.585 (Table 5) for L. gyldenstolpei and L. taylori, respectively. The sub-clusters correlated with the height of the sampling locations.
Table 4: Genetic distance of L. gyldenstopei accession using ISSR markers
AMOVA was used to partition the genetic variation using hierarchical analysis from the distance matrix. AMOVA showed a high percentage of variation within populations; 90% and 97% for Lg and Lt, respectively. There was little genetic difference between the populations collected from the six different regions (provinces); Table 5. Genetic distance of L. taylori accession using ISSR markerss 8% and 2% in Lg and Lt, respectively. Likewise, the genetic difference among populations within the same region was low; only 2% and 1% in Lg and Lt, respectively (Table 3). The PCoA analysis did not clearly separate the two-species based on the first two principal component scores, which accounted for 4.62% and 3.84% of the total variations (Figure 4). To study the relationships within the populations of L. gyldenstolpei and L. taylori, PCoA was performed. According to the result, the first group (red), which accounted for 4.62%, did not clearly separate from the second group (gray), which accounted for 3.84% (Figure 4a).
Figure 4: Relationship between geographic data and genetic distance. (A) Eigen values showing the variance and cumulative in the principal coordinate analysis (PCoA) used to classify 19 populations by ISSR markers. Limnonectes gyldenstolpei: gray, L. taylori: red. (B) Box plot showing the distribution of elevation. The lines inside the boxes mark mean values. The boxes indicate the upper and lower quartiles divided by the median, and the whiskers indicate the extreme values. (C) Relationship between genetic distance (Fst) and geographic distance of L. gyldenstolpei. (D) Relationship between genetic distance (Fst) and geographic distance of L. taylori.
Table 5: Genetic distance of L. taylori accession using ISSR markers
Isolation by distance
The box plots that resulted from height geography daHta analysis provided a clear graphical representation of the distribution of L. gyldenstolpei and L. taylori (Figure 4b). It was observed that the L. taylori population has a wider distribution than that of L. gyldenstolpei. The mean elevation was found to be 709.88 m and 891.45 m for L. gyldenstolpei and L. taylori, respectively, and a t-test showed no significant difference (t=1.025, p-value=0.3194) in the elevation at which the two species were present. The IBD test showed a relationship between the genetic distance (Fst) and geographic distance. For L. gyldenstolpei, there was a non-significant correlation (r=−0.0963 and p-value=0.6330) between the genetic and geographical distances (Figure 4c). However, for L. taylori, the correlation was significant (r=0.4257, p-value=0.0020) (p<0.05 is the level of significance) (Figure 4d).
In the present study, we investigated the genetic variability of L. gyldenstolpei and L. Taylori in the north and northeast of Thailand using ISSR markers. We found a clear genetic difference between the two species from the UPGMA analysis (data from the individual is grouped into the population), but the genetic distance between the two-frog species was very low, showing that they were closely related. The UPGMA analysis presented here provides an interesting result, particularly when considering the population of L. taylori in Ban Pa Sak Ngam, Chaing Mai (latitude=18°59¢16.0¢¢, longitude=99°06¢52.4¢¢, elevation=402 km).
The reason for this population being distinct from the others could be due to the presence of the Mae Kuang Udom Thara dam, which serves as a geographical barrier for frogs. We propose that the population of L. taylori in Ban Pa Sak Ngam has an allopatric distribution, similar to that of the Chinese wood frog (Rana chensinensis) whose genetic variation is caused by a mountain ridge barrier (Zhan et al. (2009). For both frog species, our data provides a high genetic variance within the same populations based on the AMOVA results. The PCoA plot showed an overlapped distribution of the two populations suggesting that there might share the genetic information from deep coalescence as the two populations are treated as the member of the L. khulii complex. In the present study, we found no significant correlation between elevation and genetic distance. The height at which the frogs inhabit does not significantly affect their genetic distance.
Among the populations of L. taylori, where inbreeding and a high degree of population subdivision are observed, the IBD test indicated a significant effect of geographical distance on the animals’ genetic distance. This effect, however, was not significant among the specimens of L. gyldenstolpei. The AMOVA and PCoA analyses support the phylogenetics of the two species as proposed by Funk et al. (2008). The genetic distance between the two species is very low, indicating their close genetic relationship. A very low genetic variation within the genus Limnonectes suggests that there is widespread inbreeding within this genus, which increases its risk of extinction. Inbreeding can be explained by the fact that amphibians are conservative when choosing their spawning ponds and have a low migration activity (Chernyshov and Truveller, 2006; Zhan et al., 2009). The inbreeding caused by geographical barrier was also reported in previous researches showing that mountain ridges act as barriers to gene flow in Columbia spotted frogs (Funk et al., 2005) and the isolation-by-distance factor applies to the plains populations (Palo et al., 2003). Our results could be interpreted as being consistent with the speciation theory that depends on distribution ranges and geographic barriers, which could be classified as allopatric, sympatric, or parapatric speciation (Coyne and Orr, 2004; Mayr, 1942; Schluter, 2001). The study of species co-occurrence and parapatric distributions can provide an evolutionary history of amphibian species, for example, in tree frogs (Hyla molleri and H. meridionalis) and midwife toads (Alytes obstericans and A. cisternasii). Their combined fuzzy intersection and niche comparison metrics showed the effect of the environment or niche on the survival until segregation (Reino et al., 2017). Similarly, Posso-Terranova and Andrés (2016) showed the relative contribution of geographical and environmental factors to the diversification of poison frogs (Oophaga) using gene genealogies and the phylogenetic method. The study found a parapatric distribution of sister taxa and no significant effect of elevation on the genetic distance. High levels of genetic variability and differentiation are typical for amphibians in general and exceed that for other groups of vertebrates (Mezhzherin, 1992).
In conclusion, juvenile L. gyldenstolpei and L. taylori inhabit geographically different habitats, but they have a similar morphology. Our study shows that both internal and external dynamics affect spatial genetic subdivision, which could potentially lead to Parapatric and allopatric speciations. The present study demonstrates that the ISSR technique is an appropriate tool to study the differentiation of a species complex. L. gyldenstolpei and L. taylori are genetically distinct. The genetic variation and distance point to the recent separation of the two species of Limnonectes. Further studies in other organisms using techniques mentioned in this paper should be conducted in the future.
We would like to express our sincere gratitude to Mr. Somchai Nimuan and Ms. Patcharawadee Thongkumkoon for their precious advice. This project was funded by the Graduate School, Chiang Mai University, Thailand. The Institute of Animals for Scientific Purpose Development (IAD) provided permission for the use of animals in this study (No. U1-01205-2558).
About the Authors
- Bornet B, Branchard M (2001). Nonanchored inter simple sequence repeat (ISSR) markers: reproducible and specific tools for genome fingerprinting. Plant molecular biology reporter 19, 209-215. https://doi.org/10.1007/bf02772892
- Cain AJ (1953) Geography, ecology, and co-existence in relation to the biological definition of the species. Evolution 7: 76-83. https://doi.org/10.2307/2405573
- Campos FS, Trindade-Filho J, Brito D, Llorente GA (2014). The efficiency of indicator groups for the conservation of amphibians in the Brazilian Atlantic Forest. Ecology and evolution. 4:2505-2514. https://doi.org/10.1002/ece3.1073
- Chernyshov K, Truveller K (2006). Quantitative evaluation of the genetic divergence of amphibian taxa of different rank by the spectra of total proteins and nonspecific esterases.Zhurnal obshchei biologii. 68: 149-160.
- Coyne J, Orr HA (2004). Speciation. Sinauer Associates Sunderland, MA.
- Emerson SB, Inger RF, Iskandar D (2000). Molecular systematics and biogeography of the fanged frogs of Southeast Asia. Molecular phylogenetics and evolution. 16:131-142. https://doi.org/10.1006/mpev.2000.0778
- Excoffier L, Smouse PE, Quattro JM (1992). Analysis of molecular variance inferred from metric distances among DNA haplotypes:application to human mitochondrial DNA restriction data. Genetics 131: 479-491.
- Funk WC, Blouin M.S, Corn PS, Maxell BA, et al. (2005). Population structure of Columbia spotted frogs (Rana luteiventris) is strongly affected by the landscape. Molecular ecology. 14: 483-496. https://doi.org/10.1111/j.1365-294x.2005.02426.x
- Funk WC, Pearl CA, Draheim HM, Adams MJ (2008). Range-wide phylogeographic analysis of the spotted frog complex (Rana luteiventris and Rana pretiosa) in north-western North America. Molecular Phylogenetics and Evolution 49:198-210. https://doi.org/10.1016/j.ympev.2008.05.037
- Hockin DJ, Babbitt KJ (2014). Amphibian contributions to ecosystem services. Herpetological Conservation and Biology. 9:1-17.
- Hopkins WA (2007). Amphibians as models for studying environmental change. Ilar Journal. 48:270-277.https://doi.org/10.1093/ilar.48.3.270
- Inger RF, Stuart BL (2010) Systematics of Limnonectes (Taylorana) Dubois. Current Herpetology 29: 51-68.https://doi.org/10.3105/018.029.0201
- Kojima T, Nagaoka T, Noda K, Ogihara Y (1998). Genetic linkage map of ISSR and RAPD markers in Einkorn wheat in relation to that of RFLP markers. Theoretical and Applied Genetics. 96: 37-45. https://doi.org/10.1007/s001220050706
- Lebboroni M, Ricchiardino G, Bellavita M, Chelazzi G (2006). Potential use of anurans as indicators of biological quality in upstream of central Italy. Amphibia-Reptilia. 27:73-79. https://doi.org/10.1163/156853806776052164
- Machkour-M'Rabet S, Hénaut Y, Dor A, Pérez-Lachaud G (2009). ISSR (Inter Simple Sequence Repeats) as molecular markers to study genetic diversity in tarantulas (Araneae, Mygalomorphae). Journal of Arachnology. 37:10-14. https://doi.org/10.1636/a08-27.1
- Matsui M, Panha S, Khonsue W, Kuraishi N (2010). Two new species of the “kuhlii”complex of the genus Limnonectes from Thailand (Anura: Dicroglossidae). Zootaxa. 2615: 1-22.
- Mayr E (1942). Systematics and the origin of species, from the viewpoint of a zoologist. Harvard University Press, UK. https://doi.org/10.2307/1438628
- Mayr E (1947). Systematics and the origin of species.
- Mayr E (1963). Animal species and evolution. Belknap Press of Harvard University Press Cambridge, Massachusetts, USA.
- McLeod DS (2010). Elucidating the relationships within the Limnonectes kuhlii species complex (Amphibia: Anura: Dicroglossidae).
- McLeod DS (2010). Of least concern? Systematics of a cryptic species complex: Limnonectes kuhlii (Amphibia: Anura: Dicroglossidae).Molecular Phylogenetics and Evolution 56: 991-1000. https://doi.org/10.1016/j.ympev.2010.04.004
- Meredith H, Van Buren C, Antwis RE (2016). Making amphibian conservation more effective. Conservation Evidence. 13: 1-6.
- Palo J, O'Hara RB, Laugen AT, Laurila, A, et al. (2003). Latitudinal divergence of common frog (Rana temporaria) life history traits by natural selection: evidence from a comparison of molecular and quantitative genetic data. Molecular Ecology. 12: 1963-1978. https://doi.org/10.1046/j.1365-294x.2003.01865.x
- Posso-Terranova A, Andrés, JA (2016). Complex niche divergence underlies lineage diversification in Oophaga poison frogs. Journal of Biogeography. 43: 2002-2015. https://doi.org/10.1111/jbi.12799
- Reino L, Ferreira M, Martínez-Solano Í, Segurado P, et al. (2017) Favourable areas for co-occurrence of Parapatric species: nicheconservatism and niche divergence in Iberian tree frogs and midwife toads. Journal of Biogeography. 44: 88-98.https://doi.org/10.1111/jbi.12850
- Rivas LR (1964). A reinterpretation or the concepts “sympatric” and “allopatric” with proposal or the additional terms “syntopic” and “allotopic”. Systematic Zoology 13: 42-43. https://doi.org/10.2307/2411436
- Schluter D (2001). Ecology and the origin of species. Trends in ecology & evolution. 16: 372-380.
- Schluter D, Pennell MW (2017). Speciation gradients and the distribution of biodiversity. Nature. 546: 48-55.https://doi.org/10.1038/nature22897
- Takezaki N, Nei M, Tamura K (2009). POPTREE2: Software for constructing population trees from allele frequency data and computing other population statistics with Windows interface. Molecular biology and evolution. 27:747-752. https://doi.org/10.1093/molbev/msp312
- Whiles MR, Lips KR, Pringle CM, Kilham SS (2006). The effects of amphibian population decline on the structure and function of Neotropical stream ecosystems. Frontiers in Ecology and the Environment. 4: 27-34. https://doi.org/10.1890/1540 9295(2006)004[0027:teoapd]2.0.co;2
- Zhan A, Li C, Fu J (2009). Big mountains but small barriers: Population genetic structure of the Chinese wood frog (Rana chensinensis) inthe Tsinling and Daba Mountain region of northern China. BMC genetics 10: 17. https://doi.org/10.1186/1471-2156-10-17
- Zietkiewicz E, Rafalski A, Labuda D (1994). Genome fingerprinting by simple sequence repeat (SSR)-anchored polymerase chain reaction amplification. Genomics. 20: 176-183. https://doi.org/10.1006/geno.1994.1151
- Share This