We conducted an electrophoretic survey to examine geographic genetic variation in samples from 17 localities of the Japanese clawed salamander, Onychodactylus japonicus. This species was divided into six genetic groups (N-Tohoku, S-Tohoku, Tsukuba, SW-Honshu, Kinki, and Shikoku) that were largely concordant with clades or subclades recognized in our previous mtDNA study. Although the relationships among these six groups were not clarified, genetic distances between them were not small (mean Nel's D=0. 146–0.471). Among these groups, the geographically isolated Tsukuba group was genetically distinct, possibly as a result of population isolation. In a locality of western Honshu, two groups, SW-Honshu and Shikoku, were found to occur sympatrically. Although several presumable hybrid individuals were found, hybridization between these two groups seemed to occur very rarely. These results indicate that the Shikoku group is specifically distinct from the SW-Honshu group, whose range includes the type locality of O. japonicus.
INTRODUCTION
Onychodactylus japonicus is a small, stream-breeding salamander endemic to Japan, occurring in montane regions of mainland Honshu and Shikoku Islands (Yoshikawa et al., 2008). Although the distributional range of this species is extensive, much wider than for any other Japanese hynobiid salamander, information on geographic variation in this species has been very limited.
Recently, Yoshikawa et al. (2008) conducted a molecular phylogenetic study of this species using the mitochondrial cytochrome b (cyt b) gene, and clarified that O. japonicus consists of four major clades (Clades I–IV), with considerable genetic differentiation among them. Of these, Clades I and III and Subclade II–A are distributed parapatrically, but Clades III and IV, occurring mainly in southwestern Honshu and Shikoku, respectively, were found to be sympatric in several parts of Honshu. This fact strongly suggests the presence of at least two species in O. japonicus. Two alternative explanations would be possible for the sympatric occurrences of the mtDNA groups: (1) sympatric occurrence of two reproductively isolated species, and (2) past secondary contact and hybridization between two genetic races within a single species, and retention of two mitochondrial lineages in the population. However, Yoshikawa et al. (2008) was unable to form any definite conclusions on this issue, because mtDNA is maternally inherited, and held that a genetic analysis using nuclear markers was needed to resolve this problem. The use of nuclear markers is helpful not only for drawing taxonomic conclusions but also for providing deeper insights into the history of a species by comparing nuclear variation with the results of mtDNA analyses (Shaw, 2002). Thus, many recent studies have utilized nuclear markers to clarify geographic genetic variation and uncover taxonomic relationships among the taxa in question (e.g., Weisrock et al., 2006; Fu and Zeng, 2008; Matsui et al., 2008; Okamoto and Hikida, 2009).
Here we conducted an allozymic analysis to reveal the pattern of geographic genetic variation and to assess the possibility of presence of cryptic species within O. japonicus.
MATERIALS AND METHODS
In all, 226 specimens, including 60 metamorphs and 166 larvae, were collected from 17 localities (Fig. 1; Table 1) representing all clades and subclades identified by previous mtDNA analyses (Yoshikawa et al., 2008). As an outgroup, we used Salamandrella keyserlingii, which represents another hynobiid genus. We were unable to obtain fresh tissues of the topotypic O. japonicus from Hakone-machi, Kanagawa Prefecture, but we obtained samples from Izu-shi, Shizuoka Prefecture (Sample 10), which is located approximately 40 km south of the type locality. The genetic closeness of the samples from these localities was well supported by the mtDNA study by Yoshikawa et al. (2008).
Fig. 1.
Map of Japan showing the distributional range of Onychodactylus japonicus (shaded) and localities where sampling was conducted in the present study (open squares, N-Tohoku group; closed triangles, S-Tohoku group; open triangle, Tsukuba group; open circles, SW-Honshu group; open inverse triangle, Kinki group; closed inverse triangles, Shikoku group). Sample 15 was from a locality where the SW-Honshu and Shikoku groups were sympatric. For sample numbers, refer to Table 1.
Table 1.
Species, sample numbers, localities, sample sizes, assigned genetic groups, and mean snoutvent length (SVL) of adult males and females.
From salamanders anaesthetized with a saturated acetone-chloroform solution, liver tissues (and skeletal muscles of some small larvae) were taken and stored at -80°C until electrophoresis. Voucher specimens were fixed in 10% formalin and later metamorphs were preserved in 70% and larvae in 50% ethanol. The vouchers are deposited in the Graduate School of Human and Environmental Studies, Kyoto University (KUHE), or in Mr. Tanabe's private collection.
Homogenized tissue extracts were subjected to standard horizontal starch gel electrophoresis and histochemical staining procedures (Shaw and Prasad, 1970; Ayala et al., 1972), using Starch Art (Otto Hiler, Madison, Wisconsin, USA) and potato starch (Fluka Chemie GmbH, CH-9471 Buchs, Switzerland) mixed in a 1:1 ratio and then suspended in buffer at a concentration of 12%. We analyzed 14 loci encoding 11 enzymes (Table 2) using four standard buffer systems [CAPM6, citrate-aminopropylmorpholine, pH 6.0 (Clayton and Tretiak, 1972); TC7, Tris-citrate, pH 7.0 (Shaw and Prasad, 1970); TC8, Tris-citrate, pH 8.0 (Clayton and Tretiak, 1972); and TBE8.7, Tris-borate-EDTA, pH 8.7 (Boyer et al., 1963)]. Genetic interpretations of zymograms were based on criteria developed by Selander et al. (1971). Enzyme nomenclature, E. C. numbers, and the notation of loci, electromorphs, and genotypes mainly follow IUBMB (1992) and Murphy (1996). We designated electromorphs by letters, with ‘a’ representing the most rapidly migrating variant.
For each sample, genetic variation was assessed by calculating the mean number of alleles per locus (A), the proportion of polymorphic loci (P), and the mean heterozygosity by direct count (H). Variable loci were checked using chi-square goodness-of-fit tests to determine whether they were in Hardy-Weinberg equilibrium. The expected number of each genotype was calculated using Levene's (1949) formula for small sample sizes. All these statistics were calculated by using the BIOSYS-1 computer program (Swofford and Selander, 1981). To calculate the genotypic similarity among individuals, we classified specimens as different genotypes if they differed in any alleles across all loci. We then counted the number of allelic differences between each genotype and the standard (most frequent) genotype (see RESULTS section).
To infer overall genetic differentiation, Nei's unbiased genetic distance (Nei, 1978) and Cavalli-Sforza and Edwards' (1967) chord distance were calculated between samples. Using allelic frequencies, we constructed a neighbor-joining (NJ; Saitou and Nei, 1987) tree with Cavalli-Sforza and Edwards (1967) chord distances and a continuous maximum-likelihood (CONTML; Felsenstein, 1973) tree. Confidences of tree topologies were tested by 1000 non-parametric bootstrap pseudoreplicates (Felsenstein, 1985). These analyses were performed with PHYLIP ver. 3.5C (Felsenstein, 1993).
To reconstruct genetic relationships among samples of O. japonicus, we also conducted multidimensional scaling (MDS) using SAS (SAS, 1985). This analysis can detect clinal and complicated relationships that might be overlooked in clustering procedures (Felsenstein, 1982). Because the results of MDS are stable regardless of the genetic distance used (Lessa, 1990), we used only Nei's (1978) D in this study.
Table 2.
Enzymes, loci, and buffer systems used in this study.
Table 3.
Allele frequencies and genetic variation at 14 polymorphic loci among 19 samples of O. japonicus and S. keyserlingii. For sample numbers, refer to Table 1.
RESULTS
Two genetic types found in Sample 15
Within O. japonicus, all 14 loci scored were polymorphic, and 74 alleles were detected (Table 3). The most variable loci were PGDH and PGM-3, each with eight alleles, followed by ATA-2 and MDH-1, each with seven alleles. In 11 samples (2, 4–10, 12, 14, and 15), several loci showed significant deviation from HW expectations in heterozygote frequency. Most of these were heterozygote deficient.
Salamandrella keyserlingii had unique alleles for eight loci and was monomorphic at all loci examined.
For samples that possessed loci with heterozygote deficiency, we estimated the genotypic similarity among individuals by a method similar to that employed by Tominaga et al. (2003). As a result, we found no bimodal distribution of genotypic similarity in any of the samples examined (data not shown), except for Sample 15 from Hatsukaichi-shi, Hiroshima Prefecture, Chugoku District, western Honshu.
In Sample 15, five of 14 loci examined were polymorphic (Table 3), and heterozygote deficiencies were found in two of them (ATA-2 and PGM-3). For calculation of genotypic similarity, we assumed two loci as diagnostic and used only these to exclude noise derived from polymorphisms in non-diagnostic loci. We divided all 59 specimens examined into five genotypes by differences in allelic composition. Of these five, a homozygous genotype possessing dd at ATA-2 and ee at PGM-3 was most frequent (N=27), and we thus regarded this genotype as the standard one. Then we calculated genotypic similarity by simply counting the number of allelic differences between each genotype and the standard genotype. As a result, we recognized a bimodal distribution of the genotypic similarity score in this sample (Fig. 2). Of the 59 total specimens, 27 and 26 specimens, respectively, showed genotypic similarity of 0 (standard genotype) and 4 (most different genotype), with the remaining six specimens in between (genotypic similarity= 1–3).
From this result, we suspected that specimens from Sample 15 actually consist of two distinct genetic types, with several hybrids. We thus divided this sample into two subsamples, 15A and 15B (genotypic similarity=4 and 0, respectively) in later analyses by omitting six intermediate individuals that were most likely to be hybrids. In each of these subsamples, no locus showed significant deviation from the HW expectation, and Nei's (1978) D between them was 0.160.
Table 4.
Nei's unbiased genetic distance (above diagonal) and Cavalli-Sforza and Edwards' chord distance (below diagonal) for O. japonicus (Samples 1–17) and S. keyserlingii (Sample 18).
Within-population variation
Within Onychodactylus, after splitting Sample 15 into two subsamples (15A and 15B; see above), the mean number of alleles per locus (A) varied from 1.14–2.14, the proportion of polymorphic loci (P) from 14.29–71.43, and the mean heterozygosity by direct count (H) from 0.016–0.205. The highest A, P, and H values were observed in Samples 6, 5, and 3, respectively, whereas the lowest A, P, and H values were in Samples 8 and 15A, 8 and 15A, and 15A, respectively (Table 3). We observed a large Fst (Wright, 1965), with a mean of 0.60 (range=0.04–0.89) and large values for the PEPIg (0.89), SOD (0.78), MDH-1 (0.74), ATA-2 (0.71), and PGM-3 (0.63) loci.
Overall population differentiation
The relative magnitudes of the pairwise Nei's (1978) distances and the chord distances of Cavalli-Sforza and Edwards (1967) were similar (Table 4), and only the Nei' Ds are described below. Salamandrella keyserlingii differed from samples/subsamples of O. japonicus with large genetic distances (range=0.980–1.599, mean=1.148). Nei's D varied between samples/subsamples of O. japonicus, with the mean of all pairs being 0.236. The greatest D (0.556) was observed between Samples 8 (Tsukuba Mountains) and 14 (Chugoku Mountains), while the smallest (0.000) was found between Samples 2 and 3 from adjacent localities.
Fig. 2.
Frequency distribution of genotypes in Sample 15, expressed in terms of genotypic similarity relative to the standard (most common) genotype for ATA-2 and PGM-3. The homozygous/heterozygous complement of alleles for the two loci (ATA-2/PGM-3) is shown on the bar for each genotype. N, number of specimens.
As shown in Fig. 3, the topologies of the dendrograms constructed by NJ and CONTML were nearly identical, and six major groups were consistently recognized: northern Tohoku (N-Tohoku: Samples 1–3); southern Tohoku (STohoku: 4 and 5); Tsukuba mountains (8); southwestern Honshu (SW-Honshu: 6, 7, 9–12, 14, and 15A); Kinki (13); Shikoku (15B, 16, and 17).
The N-Tohoku group was supported moderately (74%) in the NJ analysis but not significantly (<50%) in the CONTML analysis. Allele ACOH-1 (e) was relatively frequent in and nearly unique to this group. In addition, the fixation of SOD allele (c) was also characteristic of this group (Table 3; Fig. 4). The Fst value for this group was 0.08. The S-Tohoku group was highly supported in both analyses (90% and 84%), and samples of this group possessed nearly unique alleles MDH-1 (c) and PEPIg (c) (Table 3; Fig. 4). The Fst value for the S-Tohoku group was 0.10.
The SW-Honshu group, distributed over the widest range covering the southwestern half of the mainland Honshu, was moderately supported in both analyses (63% and 60%). This group contained the sample close to the type locality (Sample10), and one of the two subsamples divided according to genotypic similarity (15A; see above) was also included in this group. There were no unique alleles in this group, but the allelic frequencies for ACOH-1 and PGM-3 were characteristic (Table 3; Fig. 4), and the Fst value for this group was 0.35. The degree of genetic differentiation within this group was relatively high (Nei's D=0.015–0.226; Table 4). Although a sample from eastern-most part of the range (Sample 7) tended to diverge from the others with large D values (0.081–0.226), no significant geographic subdivision was found within this group.
The Shikoku group comprised two samples from Shikoku (Samples 16 and 17) and extralimital Subsample 15B from Honshu. This group was moderately supported (66% and 62%) in the NJ and CONTML analyses, and was nearly fixed for ATA-2 (d), which is rare in other groups, and PGM-3 (e), which was rare in the neighboring SW-Honshu group (Table 3; Fig. 4). The Fst value for the Shikoku group was 0.19.
Fig. 3.
(A) Neighbor-joining tree based on the Cavalli-Sforza and Edwards (1967) chord distance. (B) Continuous maximum likelihood tree. Nodal values indicate the bootstrap support (1000 replicates). Scales in (A) and (B) indicate the Cavalli-Sforza and Edwards chord distance. For sample numbers, refer to Fig. 1 and Table 1.
Fig. 4.
Pie graphs for four typical loci in O. japonicus, showing geographical changes in allelic frequencies. Large circles indicate sample sizes over 15; middle circles, sample sizes between 10 and 15; small circles, samples sizes under 10.
The Tsukuba group consisted of only Sample 8 from Tsukuba Mountains. This group was the most divergent from other groups, with mean D of 0.405 (0.225–0.556; Table 4). In this group, MDH-1 was fixed for a unique allele (d), and another unique allele ATA-2 (f) was found at high frequency (Fig. 4). Allele GPI (c), which was rare in other groups, was also observed at high frequency (Table 3).
Finally, the Kinki group was composed of only Sample 13 from Kyoto in Kinki district. This group had no unique alleles, but was characterized by high frequencies of the PEPIgg (b) and SOD (c) alleles, which were rare and absent, respectively, in samples from adjacent localities (Table 3; Fig. 4). For the ATA-2 and PGM-3 loci, samples of this group possessed allele frequencies intermediate between the SW-Honshu and Shikoku groups (Fig. 4).
Among these six groups, mean Nei's (1978) genetic distances ranged from 0.146–0.471, and the smallest and largest values were observed between the N-Tohoku and Kinki, and the Tsukuba and Kinki groups, respectively. The branching order was nearly identical between the NJ and CONTML trees. In these trees, a basal split separated the S-Tohoku group from the others, and then the N-Tohoku and Tsukuba groups formed a cluster and diverged from the remaining three groups (SW-Honshu, Kinki, and Shikoku), among which relationships were nearly trichotomous.
In applying MDS to our data, we repeated the estimation procedure for two-dimensional resolution until the converging values reached less than 0.01. A badness-of-fit standard value of 0.045 was finally achieved. Two-dimensional plots of the MDS (Fig. 5) separated the six groups recognized in the cluster analyses. Samples of the SW-Honshu group were dispersed in the MDS space, whereas those of the N-Tohoku and Shikoku groups were close together. Two samples of the S-Tohoku group were relatively distant from each other, and the Tsukuba group was distant from the other groups.
Fig. 5.
Two-dimensional scattergram of the multidimensional scaling (MDS) based on Nei's (1978) unbiased genetic distances among samples of O. japonicus. The N-Tohoku, S-Tohoku, Tsukuba, SW-Honshu, Kinki, and Shikoku groups are indicated as open squares, closed triangles, open triangle, open circles, open inverse triangle, and closed inverse triangles, respectively. For sample numbers, refer to Fig. 1 and Table 1.
DISCUSSION
Two genetic types in Sample 15
From the bimodal distribution of genotypic similarity for two loci, we split Sample 15 into two subsamples. Both the cluster and the MDS analyses supported this division by clearly showing the inclusion of Subsample 15A in the SW-Honshu group and 15B in the Shikoku group. This finding is also concordant with the study by Yoshikawa et al. (2008) using mtDNA analysis, which clarified the sympatry of mitochondrial Sublade III–C (a Clade–lll subclade distributed in southwestern Honshu) and IV–B (distributed mainly in Shikoku) in the sample from this locality. From the geographic information and their allocations in the trees, it is clear that Subsamples 15A and 15B in the present study correspond to Subclades III–C and IV–B, respectively. These results indicate the sympatric distribution of two genetic types in Sample 15. Although six individuals (10.2% of all individuals examined) were regarded as hybrids, this percentage would be low enough to lead the conclusion that there are some pre- and/or post-mating isolation mechanisms between the two genetic types. If these two types mated freely with each other, the proportion of hybrids would be much higher. This genetic isolation in a sympatric locality strongly suggests different specific status for the two genetic types. Because the type locality of O. japonicus is included in the range of the SW-Honshu group, the Shikoku group should be separated from this name.
Two sympatric, closely related salamanders, Hynobius naevius and H. yatsui (formerly called Hynobius naevius Groups A and B, respectively), are differentiated by body size and some ecological features, such as breeding habits and larval life history, to avoid competition (Tominaga et al., 2003, 2005). However, from the limited number of our adult specimens, such a body size difference seems unlikely in O. japonicus (Table 1). Moreover, we found larvae of both types simultaneously in a single stream (Yoshikawa et al., unpublished data).
To reveal the factors enabling the coexistence of these two genetic types, detailed morphological surveys and further field observations are strongly required. Even so, reproductive isolation between the two genetic types, i.e., the SW-Honshu and Shikoku groups, is almost certain, and they should be treated as separate biological species in future studies.
Genetic differentiation among samples of O. japonicus
Based on the cluster analyses and MDS of the allozyme data, we recognized six genetic groups within O. japonicus. (N-Tohoku, S-Tohoku, Tsukuba, SW-Honshu, Kinki, and Shikoku). These groups do not contradict those recognized in the previous mtDNA analysis (Yoshikawa et al., 2008). Based on the geographic allocations, it is certain that the N-Tohoku, S-Tohoku, Tsukuba, SW-Honshu, Kinki, and Shikoku groups correspond to mitochondrial Clade I, Subclade II–A, Subclade H–B, Clade III, Subclade IV–A, and Subclade IV–B, respectively.
Although we could not clarify the relationships among these groups, genetic divergences among them were large (mean Nei's (1978) D=O. 146–0.471). The Tsukuba group was most divergent from the other groups. This was unexpected, because samples from Tsukuba (mtDNA Subclade II–B) formed a well-supported sister group of II–A (=S-Tohoku group), with only moderate mtDNA differentiation (Yoshikawa et al., 2008). This result is likely due to a rapid rate of evolution of the allozyme genes in the Tsukuba group. This group has a limited distribution in the Tsukuba Mountains, and was estimated from mtDNA data (Yoshikawa et al., 2008) to have separated from Subclade II–A (=S-Tohoku) about 2.8 Ma. If we apply a molecular clock for salamanders of 1D=14 MY (Maxon and Maxon, 1979), 2.8 MY of mitochondrial divergence time would result in a Nei's (1978) D of 0.2, which is about half the actual D value between the Tsukuba and S-Tohoku groups (mean D=0.402) and is even smaller than the distance between the Tsukuba group and any other group. This may suggest that genetic differentiation, in terms of allozymes, of the Tsukuba group had been accelerated by some stochastic factors, as the result of population isolation, or founder or bottleneck effects. Small, isolated populations can be strongly affected by such effects to form a unique allele composition (Frankham et al., 2002). The low A, P, and H values for this sample also seem to support this idea.
Although the ranges of the N-Tohoku and Kinki groups are interrupted by the ranges of the S-Tohoku and SW-Honshu groups, the former were genetically close to each other, with the smallest mean D (0.146). Because the mtDNA results did not support a sister relationship between these two groups, this small distance may also have been caused by stochastic noise, such as parallel retention of ancestral polymorphisms.
As mentioned above, the D value between the two sympatric genetic types (Subsamples 15A and 15B), whose different specific status is almost certain, was 0.160. This value is close to or slightly greater than minimum distances reported between distinct species of other salamander lineages, e.g., Nei's (1972) D=0. 15 for Plethodon (Highton, 1989, 1990, 1999); Nei's (1978) D=0. 11 for Hynobius (Kim et al., 2007). If we use a D value of 0.160 as the criterion for separating species, it is possible that some or all of the groups recognized in this study represent different species. However, further genetic and morphological comparisons are necessary to reach a conclusion on this issue.
The parapatric N-Tohoku, S-Tohoku, and SW-Honshu groups seemed to be genetically discontinuous at several loci, suggesting restricted gene flow between them. Kinki group, represented by a single sample, retained unique allelic compositions for two loci, despite its geographical distribution within the range of the SW-Honshu group. Yoshikawa et al. (2008) found a sympatric distribution of mtDNA Subclade IV–A (=Kinki group) and Clade III (=SW-Honshu group) in the Kinki District. The sample representing the Kinki group (Sample 13) was not from the locality of sympatry, but its maintenance of unique nuclear gene markers against the SW-Honshu group may suggest reproductive isolation between the two groups. More detailed studies with dense sampling around the contact zones are necessary to clarify their taxonomic relationship.
Within the N-Tohoku group, extensive genetic diver-gences among samples had been observed in mtDNA, but allozymic differentiation between them was very low despite inclusion of the northern- and southernmost samples. This conflict may have resulted from the difference in modes of inheritance between nuclear and mitochondrial markers. There was no consistent geographic genetic structure within the SW-Honshu group, although three distinct and moderately differentiated subclades were recognized in mtDNA (Yoshikawa et al., 2008). This conflict between the two markers may have resulted from insufficient variability in allozymes (Murphy et al., 1996).
Genetic distances between samples of the Shikoku group from Shikoku and Honshu were small (Nei's D=0.001–0.024), although these groups are isolated by the Seto Inland Sea. This result is congruent with that of mtDNA and suggests recent dispersal of this group from Shikoku onto mainland Honshu during the last glacial age across the area of the current Seto Inland Sea, which disappeared during low sea-level stands (Yoshikawa et al., 2008).
Using the correlation between D values and time of separation (Maxon and Maxon, 1979; see above), divergence times among groups of O. japonicus were estimated as 4.19 Ma (range=2.72–6.13) for the S-Tohoku group from all other groups; 3.46 Ma (1.81–4.94) for the N-Tohoku group from the SW-Honshu, Kinki, and Shikoku groups; 3.7 Ma (2.59–5.25) for the SW-Honshu group from the Kinki and Shikoku groups; and 2.87 Ma (2.66–3.22) for the Kinki and Shikoku groups. These estimates largely overlap, though the last has a lower maximum value than the others. At the same time, however, these values suggest rapid or simultaneous separation of the N-Tohoku, S-Tohoku, SW-Honshu, and ancestral Kinki-Shikoku groups, and later divergence of the Kinki and Shikoku groups, and this is concordant with the hypothesis by Yoshikawa et al. (2008) derived from mtDNA data, although the absolute divergence times estimated from the allozyme data are slightly younger than those obtained from mtDNA. This general concordance between the allozyme and mitochondrial data strengthens the hypothesis, at least concerning the order of divergence events in O. japonicus.
ACKNOWLEDGMENTS
We thank G. Aoki, T. Hayashi, S. Ikeda, Y. Inoue-Watanabe, J. Naito, I. Niibe, S. Okada, H. Okawa, T. Okayama, T. Shimada, Z. Shimizu, T. Sugahara, T. Sugihara, A. Tominaga, and T. Utsunomiya for help in collecting specimens, and N. Kuraishi, T. Shimada, and A. Tominaga for assisting with experiments. We are also grateful to T. Johnson for linguistic advice, and to two anonymous reviewers for valuable comments on the manuscript. This work was supported partly by grants from the Ministry of Education, Science, and Culture, Japan (Nos. 63540599, 01304001, 11640697, and 20510215) and the Ministry of Environment to MM.
