The perennial terrestrial orchid Liparis krameri is distributed in mainland China (SW Hubei), the Korean Peninsula, Russian Far East, and Japan (Chen et al., 2009). On the Korean Peninsula, the species grows well on half-shaded slopes under trees (e.g., pine-oak forests) rich in humus with a good drainage. Pollinators of L. krameri are unknown, but the species is visited by small diptera such as fungus gnat in Japan (K. Suetsugu, pers. comm.). In LMR-4, ovaries from induced autogamy were detached, but those pollinated from different individuals set fruit, suggesting self-incompatibility. However, in Japan L. krameri is either self-compatible (but non-autogamous) or self-incompatible (thus, ‘partially self-incompatible’); fruit and seed set of selfed individuals was lower than that of crossed ones (K. Suetsugu, pers. comm.). A similar scenario is found in L. makinoana: it was self-incompatible based on individuals in a location near Mt. Sobaek, South Korea (Oh et al., 2001). However, as L. makinoana in Japan set fruits with viable seeds sometimes produced by induced autogamy, Suetsugu considered it as ‘partially self-incompatible’ (K. Suetsugu, pers. comm.). Whitehead et al. (2018) showed that plant mating systems often vary among conspecific populations based on meta-analyses of 741 populations from 105 species and cautions that estimates of outcrossing rates from single populations are often unreliable indicators of the entire species.

To investigate allozyme variation within- and between-populations of L. krameri, we collected leaf samples (from a total of 177 individuals) from five populations from southern Korean Peninsula (Fig. 1, Table 1). We cut 1 cm from the leaf tip to minimize the damage to plants; we transported leaf samples on an ice box and kept them in a refrigerator once at the corresponding author’s laboratory. We cut and crushed them with a precooled mortar and pestle in a phosphate polyvinylpyrrolidone extraction buffer (Mitton et al., 1979). We absorbed enzyme extracts onto 4 × 6-mm wicks cut from Whatman 3 MM chromatography paper (Whatman International, Maidstone, UK), which were then stored at −70°C until needed.

To compare L. krameri with the previously studied (in the same laboratory) L. kumokiri and L. makinoana (Chung et al., 2005, 2007), we followed the same method as described by Chung et al. (2007). We stained starch gels (12%) for 15 putative loci resolved from nine enzyme systems using three buffer systems: Dia-1, Dia-2, Fdh, Idh-1, Idh-2, Lap-1, Lap-2, Mdh-1, Mdh-2, 6Pgd-1, 6Pgd-2, Pgi-1, Pgi-2, Pgm, and Skdh. We designated putative loci sequentially, with the most anodally migrating isozyme designated as ‘1,’ the next ‘2,’ and so on. We also designated different alleles within each locus sequentially by alphabetical order (‘a’, ‘b’, ‘c’, ‘d’, and ‘e’).

We designated a locus as ‘polymorphic’ when two or more alleles were detected, regardless of their frequencies. We estimated the genetic diversity parameters within populations using the programs POPGENE (Yeh et al., 1999) and FSTAT (Goudet, 1995): percentage of polymorphic loci (%P_P), mean number of alleles per locus (A_P), allelic richness (AR) using a rarefaction method to account for uneven population sample sizes (Hurlbert, 1971; El Mousadik and Petit, 1996), observed heterozygosity (H_oP) and Hardy–Weinberg (H–W) expected heterozygosity or Nei’s (1978) gene diversity (H_eP). Except for AR and H_o, we also estimated these parameters for the total samples as a whole (i.e., at the species level). Hereafter, the subscript “P” indicates population means, while the subscript “S” indicates species’ (or pooled samples) means.

To test for differences between populations of L. makinoana and those of L. krameri for observed statistics, OSx (AR_P, H_oP, and H_eP), we used a permutation scheme (999 replicates) by randomly allocating whole samples to the different groups, keeping the number of samples in each group constant and calculating differences between populations of the two species for randomized statistics, RSx. We then obtained the p-value of the test as the proportion of randomized data sets giving higher values for RSx than for OSx. These calculations were performed using FSTAT (Goudet, 1995).

To test for recent decreases of N_eP (i.e., genetic bottlenecks) in L. krameri, we evaluated differences across loci between the H–W H_e and the equilibrium heterozygosity (H_eq) expected assuming mutation–drift equilibrium. Using the program BOTTLENECK (Cornuet and Luikart, 1996), we evaluated these differences using a sign test and a Wilcoxon sign-rank test under an infinite allele model. Since allelic diversity is generally lost more rapidly than H_eP (Nei et al., 1975), recently bottlenecked populations will exhibit an excess of H–W H_eP relative to H_eq (Cornuet and Luikart, 1996; Luikart et al., 1998).

We used the program SPAGeDi (Hardy and Vekemans, 2002) to calculate population-level F_IS (inbreeding) and its significance level by 999 permutations under the null hypothesis of F_IS = 0. To measure deviations from H–W equilibrium at each polymorphic locus, we calculated averages of Wright’s (1965) F_IS and F_ST (deviations from H–W equilibrium of individuals relative to their local populations, and local populations relative to the total population, respectively) following Weir and Cockerham (1984). Using FSTAT, we constructed 95% bootstrap confidence intervals (CI) (999 replicates) around means of F_IS and F_ST, and considered the observed F_IS and F_ST to be significant when 95% CI did not overlap zero. Statistical significance of differences in F_IS and F_ST between populations of L. makinoana and L. krameri was determined as outlined above for AR, H_oP, and H_eP.

To estimate the levels of genetic divergence among populations, we performed a principal coordinate analysis (PCoA) with GenAlEx (Peakall and Smouse, 2012) based on codominant genotypic distances. In addition, we assessed the genetic structure by means of the Bayesian algorithm implemented in STRUCTURE 2.3.4 (Pritchard et al., 2000). The program estimates the likelihood of the individuals being structured in a given number of genetic clusters (or genetic populations, K). We selected ‘admixture’ and ‘correlated’ as appropriate models for ancestry and allele frequencies, respectively, as events of migration and populations with shared ancestry are presumably to occur in L. krameri. We set the burn-in period and Markov Chain Monte Carlo to 50,000 and 500,000 iterations, respectively, and 20 replicates per K were run. We determined the most likely value of K by the ΔK statistics of Evanno et al. (2005). As the ΔK method tends to identify K = 2 as the top level of hierarchical structure (Janes et al., 2017), we combined it with the method of choosing the smallest K after the log probability of data [ln Pr(X|K)] values reached a plateau (Pritchard et al., 2010). The results of both methods were visualized with the aid of Structure Harvester (Earl and vonHoldt, 2012).

Seven (Dia-1, Dia-2, Pgm, Idh-1, Pgi-2, 6Pgd-1, and 6Pgd-2) out of 15 loci were polymorphic across five populations of L. krameri in southern Korean Peninsula. Variation at allozyme loci within populations was moderate or high (Table 1). At the population level, the average percentage of polymorphic loci (%P_P) was 41.1%, the mean number of alleles per locus (A_P) was 1.59, the AR was 1.55, and the mean observed (H_oP) and expected heterozygosity (H_eP) were 0.136 and 0.145, respectively (Table 1). Higher levels of %P_S (46.7%), A_S (1.93), and H_eS (0.166) were found when samples were treated as a whole (Table 1). Liparis makinoana harbored higher within-population genetic variation than L. krameri (one-sided p-values for AR_P, H_oP, and H_eP were 0.004, 0.010, and 0.010, respectively).

Among the five studied populations of L. krameri, we did not find any population with significant excess of H–W H_e under the infinite allele model (Table 2).

We found no significant deficiency of heterozygotes (at the 0.05 level) relative to H–W expectations in all but one population (LMR-4; F_IS = 0.186) (Table 2). These results, as well as the non-significant multi-population-level F_IS (F_IS = 0.059; 95% CI, −0.047 to 0.146) (Table 2), indicated a general fit to H–W expectations within populations. The value of pooled multi-population F_IS for L. makinoana (F_IS = 0.199) was significantly larger than that of L. krameri (p = 0.025). Deviation from H–W expectations due to allele frequency differences between populations was significantly different from zero for the two species (F_ST = 0.149; 95% CI, 0.095 to 0.216 for L. krameri vs. F_ST = 0.107; 95% CI, 0.067 to 0.148 for L. makinoana). We found no statistically significant differences between the two species (p = 0.587).

In the PCoA (Fig. 2), the first two components accounted for 78.1% (axis 1 = 45.0%; axis 2 = 33.1%) of the total genetic variance. Population LMR-3, which occurs in Gajisan Mt. and is located very close to LMR-1 and LMR-2, clustered, instead, with LMR-4 (Fig. 2). The best clustering scheme of STRUCTURE (K = 2 or K = 3, according to both the ln Pr(X |K) and the ΔK statistic), agreed with the PCoA (Fig. 3).

We anticipated a lower within-population genetic variation and a higher between-population genetic differentiation for L. krameri compared to L. makinoana. Our results support the expectation for within-population genetic variation, but not that for between-population genetic differentiation; although, as predicted, F_ST for L. makinoana was lower than that for L. krameri (0.107 vs. 0.149), these two values were not statistically significantly different.

The levels of genetic variation within populations of L. krameri and within the species (i.e., the samples pooled as a whole) are slightly lower than those compiled for plant species occurring mainly in the BDDG (H_eP = 0.145 vs. 0.159 and H_eS = 0.166 vs. 0.193) (Table 1). However, L. krameri maintains higher levels of genetic variation at population and species levels than those typically found in orchids (H_eP = 0.134 and H_eS = 0.135) (Table 1), in plant species with outcrossing-animal breeding system (H_eP = 0.124 and H_eS = 0.167) (Table 1), and in herbaceous perennials (H_eP = 0.096 and H_eS = 0.125) (Table 1).

Except LMR-4, H_oP is close to H_eP (multi-population-level F_IS = 0.059), which in part is consistent with the life-history traits observed for the species in Korea (self-incompatibility and pollination by unidentified insects). However, a significant deficiency of heterozygotes relative to H–W expectations (multi-population-level F_IS = 0.199) was found in L. makinoana.

Liparis krameri exhibits a relatively low degree of genetic differentiation among populations (F_ST = 0.149). This degree is lower than those compiled for plant species occurring mainly in the BDDG (G_ST = 0.175) (Table 1), orchid species (G_ST = 0.190) (Table 1), plant species with outcrossing-animal breeding system (G_ST = 0.197) (Table 1), and perennial herbs (G_ST = 0.256) (Table 1). Why do populations of L. krameri harbor moderate to high levels of within-population genetic variation and relatively low degree of between-population genetic divergence? Although most populations of L. krameri occur in mountainous locations away from the BDDG (which has been suggested to be one of the main East Asia refugia for plants) (Chung et al., 2017), the species would have endured the glacial periods in favorable areas for the survival of the species. In addition to the BDDG, some of its ‘branches’ or associated mountain systems, such as ‘Nakdongjeongmaek’ (where Gajisan Mt. is located, home of populations LMR-1 to LMR-3) (Fig. 1), as well as Jeju Island (where LMR-5 is located) (Fig. 1), would have also acted as suitable refugia for many temperate and boreal plants; populations within refugia would have maintained relatively large N_e and connectivity among them (Dolezal et al., 2012; Chung et al., 2013, 2017). Orchids produce tiny seeds that have the potential for long distance dispersal by wind or storms (Arditti and Ghani, 2000; Trapnell and Hamrick, 2004; Yukawa et al., 2012; Takashima et al., 2016). As L. krameri occurs in southern Japan, one may consider that current Korean populations might have originated by postglacial expansion from outside of the peninsula. This scenario is, however, less likely than in situ survival in glacial refugia because L. krameri harbors high levels of H_eP. Analysis of ecological niche modeling would provide us insights into Pleistocene demographic history of L. krameri on the Korean Peninsula. However, given the limited number of presence records in South Korea and little information on the distribution of L. krameri in North Korea, it is unlikely that niche modeling will generate accurate and useful predictions.

According to population genetics theory, the “genetic health” of populations is characterized by large N_e, sufficient gene flow to counteract the effects of random genetic drift, and a low degree of inbreeding (Ottewell et al., 2016). These factors are closely related to the three major population genetic parameters G_ST, H_eP, and F_IS: in healthy populations or species, one would expect low degrees of G_ST and low levels of F_IS but high levels of H_eP. Ottewell et al. (2016) developed a plain language genetic assessment approach by setting eight possible combinations based on the rating of the three parameters as ‘low’ (L) or ‘high’ (H).

The F_ST value found in L. krameri (0.149) can be considered as ‘low’ because it is lower than the threshold value (F_ST = 0.150) proposed by Ottewell et al. (2016), but also because it is lower than the average values reported for plant species occurring mainly in the BDDG, orchid species, plant species with outcrossing-animal breeding system, and perennial herbs (Table 1). The value of H_eP (0.145) can be considered as ‘high’ because, despite being lower than the value found in L. makinoana (H_eP = 0.319) (Table 1), it is higher than the average values reported for orchids, plant species with outcrossing-animal breeding system, and herbaceous perennials (Table 1). Finally, we rate the third parameter F_IS (0.059) as ‘low’ because it does not deviate from the null hypothesis of F_IS = 0 and it is lower than its congener L. makinoana (F_IS = 0.199) (Table 1) and than what is expected for herbaceous perennials (F_IS = 0.174) (Duminil et al., 2009). Thus, in terms of all (pooled) populations examined, L. krameri in southern Korean Peninsula can be represented as L/H/L (which can be regarded as a “healthy” species). According to table 1 in Ottewell et al. (2016), L/H/L belongs to the first genetic management strategy which is characterized by the following species’ ecological and demographical scenarios: (1) historically connected populations maintaining gene flow contemporarily; (2) reduced effects of outbreeding depression and local adaptation because of somewhat homogeneous habitats (the species grows well on half-shaded slopes under trees forests rich in humus with a good drainage); (3) populations acting as “source/sink” metapopulations; and (4) avoidance of loss of genetic diversity and inbreeding due to large N_e. On the basis of these scenarios, Ottewell et al. (2016) recommend focusing on maintaining and facilitating metapopulation dynamics, although translocations between populations to boost population sizes is also a viable option.

In terms of individual populations, the only way to know whether a population is low (L) or high (H) regarding F_ST is by getting a matrix of pairwise F_ST values and check whether a given population has low or high levels of F_ST with the rest of populations. There were no significant differences among the five means (from 0.107 in LMR-2 to 0.154 in LMR-5) for F_ST pairwise values of each population with the rest of the populations (Kruskal-Wallis test, H = 5.287, p = 0.259). Thus, we considered that the allele frequencies are rather homogenous among the five populations and rated each F_ST mean as ‘L’ because Wright’s F_ST (0.149) estimated from five populations was considered as ‘L’. Except LMR-4, all populations of L. krameri could be expressed as L/H/L, and, thus, could be considered as “healthy” ones. The population LMR-4 could be represented as L/H/H because F_IS (0.186) is similar to that of L. makinoana (F_IS = 0.199) and slightly higher than the mean reported for herbaceous perennials (F_IS = 0.174) (Duminil et al., 2009). According to the ecological/demographical scenarios proposed by Ottewell et al. (2016), recent fragmentation or reduction in population size would have occurred in LMR-4, causing inbreeding. As management strategies, approaches to reduce inbreeding (e.g., facilitation of pollen/seed immigration and active translocations) should be initiated (Ottewell et al., 2016).

Finally, as L. makinoana could be represented as L/H/H, for the long-term survival of this species management strategies should be conducted as suggested above for LMR-4. As L. kumokiri is a genetically depauperate species at neutral loci and mainly occurs under forests showing a narrow ecological amplitude, it might lack adaptability to environmental changes and, thus, is prone to extinction. The species blooms during the rainy season and is rain-triggered self-pollinated (Suetsugu, 2019). From a conservation and restoration perspective, L. kumokiri is of higher priority regarding conservation and restoration actions, which should be focused on increasing the number of individuals and populations.