Genetic diversity and population structure of the endangered orchid Pelatantheria scolopendrifolia (Orchidaceae) in Korea

Due to substantial population decline, the Korean orchid P. scolopendrifolia is considered endangered and highly threatened. Like many endangered species, it is vulnerable to biological and anthropogenic threats that can lead to the loss of genetic diversity and, ultimately, extinction. Therefore, the assessment of genetic diversity and population genetic structure is imperative for conservation. In this study, we newly developed 15 polymorphic microsatellite markers. Analyses of genetic diversity and population genetic structure that included 182 samples from 11 populations were conducted using microsatellite markers and four noncoding regions of chloroplast DNA. Our study revealed a relatively low level of genetic diversity (Ho = 0.529, He = 0.356), albeit harboring with private alleles based on microsatellite genotyping data, and high haplotype diversities based on chloroplast DNA sequences data. The results of STRUCTURE and PCoA based on microsatellite genotyping data showed population differentiations. An AMOVA based on chloroplast DNA sequence data further corroborated these conclusions, indicating about 70% of variations found among populations. Low genetic diversity and divergence among the population might have been caused by factors, such as asexual reproduction, demographic events (bottleneck and population expansion), geographic isolation, and low gene flow. The development and implication of conservation strategies and management of P. scolopendrifolia are proposed based on these results.


Introduction
Species with small population size or small geographic ranges subsequently tend to have reduced gene flow and increased genetic divergence, due to genetic drift and inbreeding [1]. In particular, many endangered species have low genetic diversity and undergone population differentiation as a result of overexploitation [2,3] and habitat disturbance [4][5][6]. Because vulnerabilities to biological and anthropogenic threats are directly linked to loss of genetic diversity and extinction, the assessment of genetic diversity and population genetic structure is imperative for development and implication of conservation strategies. The family Orchidaceae is one of the most species-rich families, and this high diversity along with species-specific pollination has attracted numerous biologists [7][8][9]. Pollination success may affect the long-term survival of population or species, but it strongly depends on the presence of pollinator. Therefore, it may be vulnerable to anthropogenic threats, such as habitat disturbance, which affect both pollinator and orchids. The overexploitation of natural populations by gathering for ornamental and medicinal uses has been a major subject of conservation issues, and at least recently, awareness about the conservation of biodiversity has grown. In addition, terrestrial orchids are known to have small seeds, and germination is difficult due to requirements of specific nutrients and environmental conditions [10][11][12]. These factors may make restoration and conservation of plant at risk difficult and challenging. Given the conservation priority of highly threatened orchid species, several studies have been conducted, evaluating genetic diversity and structure based on various molecular markers, such as isozymes [13], allozymes [14,15], AFLPs [amplified fragment length polymorphisms; [16][17][18]], ISSRs [inter simple sequence repeats; [19]], and SSRs [simple sequence repeats or microsatellites; [20,21]]. In particular, studies of microsatellite markers remain popular in population genetics and conservation genetics because these markers are co-dominant, highly polymorphic, selectively neutral, and have high mutation rates [22,23]. Chloroplast (cp) DNA markers are also useful for evaluating population genetic diversity and structure [24][25][26] as well as phylogenetic analyses [27][28][29] because of non-recombinant and uniparental inheritance [30].
Pelatantheria scolopendrifolia (Makino) Aver. [= Cleisostoma scolopendrifolium (Makino) Garay] is an epiphytic or lithophytic evergreen perennial orchid, belonging to subtribe Aeridinae (Orchidaceae), and occurs across Korea, China, and Japan. In Korea, P. scolopendrifolia is highly restricted to southwestern corner of the Korean Peninsula, including coastal (i.e., Mokpo and Haenam) and inland (i.e., Naju) regions as well as islands of the southwest (i.e., Jeju, Bogil, Jindo, and Dongrae). With a narrow extent of occurrence (EOO, 2,000 km 2 ) and substantial population decline due to anthropogenic activities, such as illegal collection and habitat disturbances, the species is assessed as EN B2b(iii,iv)c(iii,iv,v) on Korean Red list and protected by law as Endangered Wildlife [31]. The rarity of this species and its need for conservation have influenced studies of P. scolopendrifolia. For example, given the fact that genetic diversity is particularly important to the conservation of rare or endemic species [32], allozyme-based and microsatellite-based studies were carried out, providing some baseline genetic information for conservation [2,33]. Their studies, however, yielded somewhat conflicting results; low genetic diversity within populations (He = 0.002) based on allozymes [2] versus relatively high genetic diversity (He = 0.4162) based on microsatellite markers [33]. In addition to population genetic studies, Son et al. [34] reported the male megachilid bee (Megachile yasumatsui) as a pollinator and self-incompatibility of this species. Furthermore, Son et al. [35] also presented a positive relationship between the frequencies of pollinator visits and fruiting rates. Given that preliminary studies have shown contradictory results, further genetic studies of P. scolopendrifolia using various markers and extensive sampling are necessary to provide more comprehensive genetic information and also to develop conservation guidances based on this genetic perspective.
Thus, the present study aimed to develop additional microsatellite markers of P. scolopendrifolia and to determine the genetic diversity, population genetic structure, and the factors that contribute to genetic patterns using the newly developed markers and four non-coding chloroplast DNA regions. Based on the findings of this study, we proposed several baseline guidelines or suggestions for conservation and management strategies of highly threatened species of P. scolopendrifolia in Korea.

Plant materials
We sampled P. scolopendrifolia from 11 populations, including two in Japan, to analyze genetic diversity and structure (Table 1). In Korea, we collected a total of 166 individuals, including inland population from Naju (NJ; 32 individuals), two coastal populations from Mokpo (MP-1, 5 individuals and MP-2, 13 individuals), and one in-land (HN-1; 14 individuals) and one coastal (HN-2; 7 individuals) population from Haenam region. Samples were also collected from the island of Jindo (JD-1, 58 individuals and JD-2, 11 individuals), Gwanmaedo (GM; 6 individuals), and Wando (WD; 20 individuals). Lastly, a total of 16 individuals from two populations (eight individuals per population) in Japan were also included as island populations to understand the origin of P. scolopendrifolia in Korea. These have been transplanted from natural population in Kyushu to the Batanic Garden for conservation. All necessary permits were obtained for the described study, which complied with all relevant regulations. The number of individuals randomly sampled per population was proportional to population size (raning from five in small populations to 58 in large populations), and a minimum collecting distance of 20 cm between sampled individuals was used to prevent resampling the same individual. All leaf samples were dried with silica gel, and voucher specimens for each population were deposited in the Chonnam National University Herbarium (CNU).

DNA extraction and microsatellite markers development
Total DNA was extracted using the DNeasy Plant Mini Kit (Qiagen, Carlsbad, California, USA), following the instructions of the manufacturer. The leaves sampled from two representative    [36]. The PrimerPipeline program which incorporates MIcroSAtellite identification software (MISA) [37] and Primer3 [38] was then used for mining the microsatellites. In order to search SSR loci, the following minimum thresholds were applied: di-, tri-, tetranucleotide with eight repeats, pentanucleotide with seven repeats, and hexanucleotide with six repeats. The primer set was designed according to the following parameters: annealing temperature of 57-63˚C, product size of 100-300bp, and length of 18-27bp. By comparing the results obtained from MP-2 and NJ individuals, 32 primer pairs with the same motif and different repeat number were preferentially selected, and then 22 primer pairs were randomly chosen from the SSR candidate list of MP-2 individual. Then, 54 primer pairs labeled with either 6-FAM or HEX fluorescent dye were used for initial polymorphism tests for eight individuals from six population (NJ, MP-1, MP-2, WD, HN-1, and JD-1). The 15 primers pairs which had confirmed polymorphism were finally applied to genotype the remaining 174 individuals. PCR amplifications were performed in a total volume of 15μL containing 0.5μL of genomic template DNA, 0.2μL of 10 pmol forward primer with fluorescence dye, 0.2μL of 10 pmol reverse primer, 1.0μL (with 2.5 mM MgCl 2 ) of PCR buffer, 0.2μL (each 10 mM) of dNTPs, and 0.05μL (5 U/μL) of Taq DNA polymerase (Inclone Biotech, Gyeonggi-do, Korea). PCR conditions were an initial denaturation at 95˚C for 3 min; followed by 35 cycles at 95˚C for 1 min, annealing at 60˚C for 1 min, extension at 72˚C for 2 min, and a final extension at 72˚C for 5 min. PCR amplicons were then separated with a GeneScan 500 LIZ Size Standard (Applied Biosystems, USA) at the Macrogen Company (Seoul, Korea) utilizing ABI 3730XL DNA Sequencer (Applied Biosystems, USA). Allele size was scored using GeneMapper v5.0 (Applied Biosystems).

Microsatellite genotyping and chloroplast DNA sequencing
To analyze the genetic diversity and population structure of P. scolopendrifolia, 15 newly developed microsatellite markers in this study and additional three microsatellite markers developed by Han [33] were used ( Table 2). PCR amplification and detection of allele size were conducted in the same conditions as described previously.

Data analysis
Microchecker 2.2.3 [40] was used to detect evidence of null alleles and potential genotyping errors. Null allele frequencies for each locus were estimated with the expectation maximization (EM) algorithm described in Dempster et al. [41] using the program FreeNA [42]. Genetic diversity indices, such as the number of different alleles (Na), number of effective alleles (Ne), observed heterozygosity (Ho), expected heterozygosity (He), proportion of polymorphic loci (P), and private alleles (Pa), were estimated based on the microsatellite allele frequencies over 18 loci using GenAlEx v. 6.5 [43]. Allelic richness (Ar) was calculated based on minimum sample size of 5 individuals in each of the 11 populations using FSTAT 2.9.4. software [44]. Deviations from Hardy-Weinberg Equilibrium (HWE) and Linkage Disequilibrium (LD) were estimated using Genepop version 4.2 [45]. The genetic differentiation coefficient (F IS and Fst) was estimated for all populations, and an analysis of molecular variance (AMOVA) was conducted with 999 permutations at two hierarchical levels (among populations and within populations) using GenAlEx v. 6.5 [43]. Principal coordinate analysis (PCoA) was carried out based on a codominant genetic distance matrix to examine genetic similarities and relationships between individuals and among populations using GenAlEx v. 6.5 [43]. A Mantel test [46] was conducted to determine correlations between genetic and geographic distance using GenAlEx v. 6.5 [43] with 999 permutations. For this analysis, we utilized the two datasets: one included all populations while the other consisted of nine Korean populations only. The genetic structure of populations was analyzed based on a Bayesian algorithm by STRUCTURE v.2.3.4 [47]. The admixture model was used with correlated allele frequencies applying burn-in of 2 x 10 5 and runs of 2 x 10 5 repetitions. Twenty runs were conducted for each value of potential genetic clusters (K), with K ranging from 1 to 15. The best K was identified based on the approach of Evanno et al. [48] using the program STRUCTURE HARVESTER v.0.6.94 [49]. The graphical results were displayed using CLUMPAK [50]. Bottleneck analysis was performed using BOTTLENECK v.1.2.02 software [51]. Heterozygosity excess caused by recent population bottlenecks was tested using a onetailed Wilcoxon signed-rank test under a two-phase model (TPM) with a variance of 12% and a 70% proportion of a stepwise mutation model in the TPM. Mode-shift tests were also used to detect genetic bottlenecks. The four combined cp DNA data were used to construct an unrooted haplotype network using the statistical parsimony approach [52] as implemented in TCS 1.21 [53]. The analysis was run with gaps coded as missing data, and the connection limit excluding homoplastic changes set to 95% following Hart and Sunday [54]. Furthermore, consecutive point-mutations were considered as a single-step mutation event, and messy poly(A)/(T) regions were deleted in this analysis. Genetic diversity parameters such as nucleotide diversity and gene diversity were calculated using Arlequin version 3.5 [55]. Selective neutrality indices, Tajima's D [56] and Fu's Fs test [57], were conducted to assess population expansion. In addition, an analysis of molecular variance (AMOVA) implemented using Arlequin version 3.5 [55] was employed to estimate the genetic variation among and within populations. Haplotype sequences of P. scolopendrifolia were deposited in GenBank: MT333568-MT333587 for rpl20-rps12, MT333588-MT333607 for petA-psbJ, MT333608-MT333627 for accD-psaI, and MT333628-MT333647 for petL-psbE.

Development of novel polymorphic microsatellite markers
Among the 54 paired primers designed and synthesized, 15 polymorphic and 3 monomorphic microsatellite markers were newly developed for P. scolopendrifolia (Table 2). Repeat motifs were detected 6 dinucleotides, 7 trinucleotides, 1 tetranucleotides, and 1 complex. The allele size ranged from 143 to 284 bp and the number of alleles per locus ranged from 2 to 9, with an average of 5.13 alleles per locus. Eighteen newly developed microsatellite markers were deposited at GenBank (MN592878-MN592895) and 15 polymorphic markers were used in this study.

Genetic diversity based on microsatellite markers
The missing data included in the final data set showed a low frequency of 0.15%. The presence of null alleles was checked using Microchecker 2.2.3 [40], and there were no regular patterns in the existence of null alleles in genotyping. While populations MP-1, MP-2, HN-1, GM, JPN-1, and JPN-2 detected no null alleles, the rest of the populations had null allele(s) at one or two loci; locus PS25 in HN-2, loci PS6 and CW3164 in JD-1, loci PS6 and CW3189 in JD-2, loci PS40 and PS43 in NJ, and loci PS25 and PS29 in WD ( Table 4). The program FreeNa [42] was used to calculate null allele frequencies and created global Fst values across all loci. The Fst values (Fst = 0.29697, 95% CI (confidence interval) = 0.24707-0.34803) after excluding null alleles (ENA) were similar to the Fst values (Fst = 0.30185, 95% CI = 0.25362-0.35281) derived from raw microsatellite data, suggesting that the effect of null alleles on the genetic structure of populations is likely negligible. Twelve loci deviated significantly from Hardy-Weinberg equilibrium (Table 3). Tests for the significance of linkage disequilibrium over all loci and populations showed linkage disequilibrium between all but 33 pairs of loci (S1 Table).
Genetic diversity was evaluated for 11 populations of P. scolopendrifolia (Table 4). At the population level, the number of alleles (Na) across loci varied from 1.278 (JPN-2) to 4.056 (JD-1). The mean of expected alleles (Ne) was 1.810. While the values of observed heterozygosity (Ho) varied from 0.278 (JPN-2) to 0.778 (HN-2), the values of expected heterozygosity (He) ranged from 0.139 (JPN-2) to 0.479 (JD-1). Ho was greater than He in all populations. Regarding the percentage of polymorphism (P), the smallest value was found in JPN-2 (28%) and the largest value (100%) was confirmed in the population JD-1. All populations other than HN-2 had one or more private alleles. Allele richness across loci showed that JD-1 had the highest mean value of 2.808. Negative F IS values were identified in all populations.

Population genetic structure based on microsatellite markers
PCoA results for P. scolopendrifolia populations based on 18 microsatellite markers can be seen in Fig 1. The first two principal coordinates of individual genotypes explained 28.05% of the total genetic variance (15.92% and 12.13%, respectively). Although some individuals belonging to different populations had overlapped positions (e.g., HN-1 and MP-2, and MP-2 and NJ), most of the population showed generally distinct distributions from each other. The distribution of Korean and Japanese populations were particularly different. In AMOVA analysis, 38% of the total molecular variation was attributed to among regions while 30% and 32% of the molecular variation resided among and within populations, respectively (Table 5). We conducted an additional AMOVA analysis by setting three population types (costal, inland, and island) as regions. The result detected no variance among regions while 55% and 45% of molecular variation was respectively partitioned among and within population (S4 Table). The Mantel test showed a positive correlation between genetic and geographic distances in both datasets (R 2 = 0.5286, P = 0.001 using all populations; R 2 = 0.1867, P = 0.001 using nine Korean

PLOS ONE
populations) (S1 Fig). According to the STRUCTURE analysis, an optimal number of genetic clusters determined as ΔK reached a maximum value at K = 5. However, ΔK was also high for K = 2 and K = 8 (Fig 2). At K = 5, HN-1 consists of a population-specific genetic cluster (magenta). Although HN-1 and HN-2 are in the same region (Haenam), they showed different genetic structures. HN-2 was included in the same genetic cluster as the Wando and Naju populations. JD-1 had three genetic clusters shared with JD-2 and GM (green), and HN-2, NJ and WD (orange) as well as MP-1 and MP-2 (sky-blue). The two Japanese populations shared the same genetic cluster with each other, which was not found in any Korean populations. According to the Wilcoxon test, eight populations had a significant heterozygote excess under the TPM and SMM (p < 0.05), whereas nine populations showed shift mode, suggesting the occurrence of a recent bottleneck (Table 6).

Chloroplast haplotype frequency and distribution
The sequences of accD-psbI, petA-psbJ, petL-psbE and rps12-rpl20 with polymorphic sites were successfully amplified and sequenced in 182 individuals from the 11 populations. The total

PLOS ONE
Genetic patterns of P. scolopendrifolia length of the combined alignment was 2,370 bp. A total of 20 chloroplast haplotypes (A-T) were identified from this study (Table 1 and Fig 3). and ΔK [48] using STRUCTURE HARVESTER [49]. Best K = 5. (b) Bayesian clustering results obtained using STRUCTURE [47]. CLUMPAK [50] was used for visualization. Vertical black bar separate individuals of different populations and colors correspond to specific clusters. https://doi.org/10.1371/journal.pone.0237546.g002

PLOS ONE
Genetic patterns of P. scolopendrifolia geographically close populations (HN-2 and WD). Although HN-1 and HN-2 were in the same region on the administrative district, they had different haplotypes. Notably, haplotypes J and M of population NJ were shared with JD-2 and GM, yet these haplotypes were not identified in other populations between NJ and Jindo region (Fig 4).  Table 7). AMOVA analysis showed significant genetic variation in among population within regions (Korea and Japan) (

PLOS ONE
using three population types as regions, most genetic variation occurred among populations within regions (63.06%) (S4 Table).

Development of microsatellite markers for population genetic study
This study reported microsatellite markers newly developed in Pelatantheria scolopendrifolia, one of highly threatened and endangered orchid species in Korea. A total of 15 polymorphic microsatellite markers were successfully amplified in all individuals from 11 populations, showing observed heterozygosity values from 0.156 to 1.000 (mean Ho = 0.481) and expected heterozygosity from 0.119 to 0.527 (mean He = 0.325). Although the values of heterozygosity were similar to or lower than those of other orchid species, for example, Caladenia huegelii (mean Ho = 0.551 and He = 0.690) [58] and Cypripedium kentuckiense (mean Ho = 0.436 and He = 0.448) [20], the presence of private alleles and the average number of alleles per locus of 5.13 suggest that the markers were appropriate for evaluating genetic diversity in this study. We could not verify the transferability of microsatellite markers into closely related species due to there being only one genus one species in Korea. Even though microsatellite markers were known to have species-specific characterization, there are many studies that successfully apply microsatellite markers into a closely related species. For example, this has been shown in Pinus resinosa (Pinaceae) [59], Allium sativum (Lilliaceae) [60], Artocarpus altilis (Moraceae)

PLOS ONE
Genetic patterns of P. scolopendrifolia [61], and Prunus mongolica (Rosaceae) [62]. Therefore, these sets of microsatellite markers will likely prove useful in future population genetic studies related to the genus Pelatantheria as well as P. scolopendrifolia specially.

Patterns of genetic diversity
Since species with small population size or limited geographic range are vulnerable to reduced genetic diversity, studies on the genetic diversity of these species are of particular importance for conservation [63][64][65]. To assess genetic diversity of the endangered orchid Pelatantheria scolopendrifolia, we used two types of molecular markers, 18 microsatellite markers and four non-coding cpDNA markers. In this study, we found that the values of observed heterozygosity based on microsatellite genotyping dataset ranged from 0.278 to 0.778 with an average value of 0.529. The values of expected heterozygosity (mean Ho = 0.356) were all lower than those of observed heterozygosity, suggesting heterozygosity excess. The mean He and Ho values observed in P. scolopendrifolia are lower when compared to other endangered or endemic orchids. For example, Caladenia huegelii had a mean He = 0.69 and mean Ho = 0.551 [58], and Cymbidium tortisepalum had a mean He = 0.653 and mean Ho = 0.619 [21]. On the other hand, values of expected heterozygosity of P. scolopendrifolia were similar to Gastrodia elata which was evaluated as having low genetic variation [66]. Rare and endemic species as well as geographically restricted species are expected to have lower levels of genetic diversity than widespread species [67][68][69][70]. Compared to these results, genetic diversity of P. scolopendrifolia is lower to the mean He and Ho values reported for narrow and regional ranged species [69]. Given that genetic diversity is deeply influenced by the life cycle, distribution range, the pattern of reproduction, and gene flow [71][72][73], life cycle can be ruled out as a cause of low genetic diversity due to the longevity of this species. In accordance with the vulnerability of small populations to genetic drift [1] and subsequent increases of homozygosity [74], these factors likely play a large role in the low genetic diversity of both Korean and Japanese P. scolopendrifolia populations. Historical demographic events, such as founder effects, bottlenecking, population expansion and contraction, and habitat fragmentation, can all impact genetic diversity and structure [63,75,76]. The previous genetic study of P. scolopendrifolia using allozymes suggested that bottleneck and founder effects are possible causes for extremely low levels of allozyme variation [2]. In this study, evidence of bottleneck effects were found via BOTTLENECK analysis (Table 6). Nine out of 11 populations showed a distribution of alleles with a shifted mode, suggesting the presence of a recent bottleneck. In contrast, the large populations NJ and JD-1 indicated a normal L-shaped mode, as expected in the absence of a recent bottleneck. Also, in spite of high haplotype diversity, gene and nucleotide diversity values are low (Fig 3,  Table 7), indicating rapid population expansion [76]. In Tajima's D and Fu's Fs test to determine the presence of the population expansion, six populations showed negative values on Tajima's D and Fu's Fs, but only NJ was significant (p < 0.05), corroborating the aforementioned result. The star-like shape of the haplotype network also supported this assumption. Therefore, this study suggests that demographic events including bottlenecks and population expansion have played important roles in reducing genetic diversity of P. scolopendrifolia. Reproductive pattern may also contribute to this low genetic diversity. Mating system can be very influential, as reproduction through outcrossing is known to lead to a high level of genetic diversity [77]. Chung et al. [2] reported that P. scolopendrifolia is self-compatible, while Son et al. [35] came to the opposite conclusion based on the correlation between the pollinator visits (male Megachile bee, Megachile yasumatsui) and fruit set. The diversity of the populations MP-2, HN-2, and WD confirmed pollinator visits were not significantly higher than the populations without a pollinator. Furthermore, gene diversity, nucleotide diversity, and haplotype diversity based on cpDNA sequences were unable to reveal any relationship with the presence or absence of a pollinator. These results are contrary to the population genetic theory which generally holds that the genetic diversity of outcrossing species tends to be higher [69]. This may be because P. scolopendrifolia can be reproduced through clonal ramets in addition to outcrossing. Vegetative reproduction is one of the causes low genetically variation [78,79]. Conversely, the patterns of allele richness and private allele abundance in P. scolopendrifolia demonstrate that genetic diversity tends to be higher in larger population sizes due to the higher likelihood of maintaining a robust genetic pool. The gene flow is another factor to explain low genetic diversity. In general, pollen and seed dispersal are important for maintaining gene flow. Gene flow between populations prevent from differentiation and increase genetic diversity between populations. Despite the presence of pollinators, the reason for low genetic diversity may be related to the behavior of the pollinator itself. It is not known how far Megachile bees fly, but male bees typically carry pollen farther distances than female bees [80], and Euglossine bees can return to their nest from as far 23 km in tropical rain forest [81]. The distance between the populations of P. scolopendrifolia with and without pollinators is approximately 30 km. Therefore, it is probably that little to no pollinator activity occurs between populations of P. scolopendrifolia. Based on these facts, it is inferred that gene flow by pollen dispersal is weak. This is consistent with PCoA and STRUCTURE results indicating the population differentiation (Figs 1 and 2). However, the haplotypes shared between or among populations may provide evidence of gene flow through seed dispersal by wind (Fig 4). The fruits of P. scolopendrifolia are matured during late summer and autumn. Because Korea has a monsoon climate, winds from the southeast or south blow in summer and winds from the northwest blow in winter. In the fall, there is no distinct wind direction as it corresponds to the transitional period of the two monsoons. Therefore, the direction of seed dispersal can be random and haplotype diversity might not distribute in one direction. Various factors, including small population size, demographic history, asexual reproduction, and low gene flow, likely contribute to the low levels genetic diversity seen in P. scolopendrifolia, as well as causing significant LD values (S1 Table).

Population differentiation and genetic structure
The mean Fst value (0.299, S2 Table) among populations measured by microsatellite markers was lower than the Fst value (0.899) based on allozymes [2]. Comparison of population differentiation generated via different types of markers might be not appropriate; however, the mean Fst of P. scolopendrifolia is higher than that of long-lived perennial, narrow-ranged, and outcrossing species (0.19, 0.23, and 0.22, respectively) [69]. This result was supported by a PCoA based on microsatellite data (Fig 1) and AMOVA analysis based on cpDNA data, with both showing that the majority of genetic variation existed among populations within regions ( Table 8, S4 Table). This paints a picture of differentiation among populations, that is consistent with the previous study [2]. These results might be explained by limited gene flow among populations and their isolation. The results of the Mantel test also indicate positive correlations between geographic and genetic distances, supporting population differentiation by isolation (S1 Fig). Significant pairwise genetic differentiation was detected between NJ and JD-1/JD-2 (S2 Table). The two populations are approximately 70 km apart; however, they had the smallest different in Fst values, further supported by genetically similarity based on their haplotype distributions (Fig 4). Although it is not supported by STRCUTURE analysis, haplotypes J and M were shared with each other. This may additionally raise the possibility of seed dispersal from JD-1 or JD-2 to NJ that has recently undergone population expansion. Nonetheless, a more thorough and clear understanding of these genetic pattern will require future study involving more rigorous comparison of NJ and the Jindo populations. The Bayesian analysis grouped the populations into 5 clusters (K = 5; Fig 2). The admixture structure existed only weakly, and most populations showed a simple structure with geographically close populations (e.g., GM and JD-1/JD-2) having the same structure. This is also somewhat consistent with the patterns of shared haplotypes between geographically close populations other than JPN-1 and JPN-2 (Fig 4). Although identical genetic structure and shared haplotypes indicate pollen or seed dispersal between or among populations, bottleneck or founder event and geological barriers to prevent gene flow could be factors to affect simple genetic structures (Fig 2) and populationspecific haplotypes (Fig 4). Notably, JD-1 is composed of three genetic clusters, indicating the presence of three genetically distinct subpopulations. Therefore, the founders of this population either had substantially different genetic compositions, or each subpopulation had a different origin.

Haplotype distribution and dispesal of P. scolopendrifolia
According to Qian and Ricklefs [82], during periods of glacial maxima (about 18,000 years ago) in the Pleistocene China, Korea, and Japan were connected by continental shelf. Simulated vegetation maps of East Asia during the last glacial maximum [82,83] showed vegetation distributions among East China, the Korean Peninsula and Japan that coincide with the current distribution pattern of P. scolopendrifolia. This raises the possibility that the influx of P. scolopendrifolia to Korea was not only from China, but also from Japan. Three major haplotypes (A, J, and O) in the haplotype network support this notion. In addition, founder events may also have influenced the haplotype distribution. After the influx of P. scolopendrifolia from China or Japan, genetic drift related to founder effects or bottlenecks in the process of dispersal from southern to northern populations would have contributed to the low level of genetic diversity and simple population structure. In this process, the Jindo populations (JD-1 and JD-2) with relatively high genetic variation might be situated as a geographic center of seed dispersal. However, after colonization, factors such as excessive collection, asexual reproduction, geographic isolation, and the loss of genetic diversity due to low gene flow provide strong reasons for the absence of a specific pattern in the haplotype distribution map (for example, the relevance of distance and genetic structure).

Implications for conservation
P. scolopendrifolia is categorized as an endangered species according to IUCN Red List criteria: B2b(iii,iv)c(iii,iv,v). The species is also listed as vulnerable in Japan following 'Japanese Red Lists' (https://www.nationalredlist.org). Therefore, the development of policies for the conservation of P. scolopendrifolia is essential. Fortunately, since 1998, it has been legally protecting P. scolopendrifolia in Korea. Information on genetic diversity and population genetic structure, as well as ecological knowledge, should be considered when developing and designing conservation strategies for species or populations. Subsequently, this study provides several suggestions. First, population size must be stabilized. Most of the populations were not large, and MP-2, HN-1, JD-1, and JD-2 are easily accessible along the trail, leaving them vulnerable to illegal collection. Therefore, periodic monitoring and strict prohibition of illegal collection seem necessary. Second, the establishment of more efficient conservation strategies will require more extensive research on P. scolopendrifolia ecology. Han [33] reported that the lichen Ramalina intestiniformis was commonly found in association with all populations of P. scolopendrifolia. Though this relationship is not well understood, it appears related to strategies for preserving or obtaining moisture. Factors such as these can be used to create more preferable environments for P. scolopendrifolia reintroduction, thus increasing the success rate of settlement. Lastly, genetic diversity and structure should be carefully considered. Genetic diversity was evaluated as relatively low in this study as a result of low diversity in allozyme data [2]. However, the existence of population-specific private alleles and haplotypes, and the unique genetic structures of population can guide the selection of the most ideal population from which to source individuals for reintroduction to a specific location. Realizing this will require the combination of seed collection from each population and ex-situ conservation. In addition, Son et al. [34,35] confirmed the existence of pollinators and revealed the relationship between pollinator visits and fruit set. Given that self-incompatibility and need of biotic pollinator for pollination, human-mediated pollination can subsequently be used as a strategy to increase the fruit set rate and genetic diversity. Pollinator protection strategies such as pollinator habitat enhancement and pest management are also necessary due to the local distribution of the pollinator.

Conclusions
In this study, we newly developed 15 polymorphic microsatellite markers of P. scolopendrifolia and provided information related to genetic diversity and population structure based on microsatellite and cpDNA data. Low levels of genetic diversity and population differentiation were observed. In addition, the identification of the private alleles and population-specific haplotypes, further enhances the need for conservation. However, it was not possible to find a connection between genetic structure or haplotype distribution and specific factors such as the presence of pollinator. Various causes such as small population size, demographic events, clonal reproduction, and low gene flow all likely have role in shaping the observed genetic patterns of P. scolopendrifolia. Lastly, taking these patterns into account informed several suggestions for the development of future conservation strategies.
Supporting information S1