Genetic and reproductive characterization of distylous Primula reinii in the Hakone volcano, Japan: implications for conservation of the rare and endangered plant

Genetic and ecological evaluation are crucial in effective management of rare and endangered species, including those exhibiting complex breeding systems such as distyly. We studied a threatened distylous herb Primula reinii in the Hakone volcano, central Japan, to obtain baseline information of reproductive and genetic status towards conservation. In two representative populations inhabiting a central cone and somma of the volcano, population size, floral morph ratio, stigmatic pollen deposition, and fruit-set were measured. Using microsatellite markers, we evaluated genetic diversity, structure and differentiation of populations. Population bottlenecks and historical changes in population size were also estimated from genotype data. We found significant deviation from equal morph ratios in the central cone population, which also exhibited skewed mating success together with a high frequency of pollination within the same morph. These trends were not detected in the somma population. From genetic insights, the central cone population showed slightly lower genetic diversity, whereas no significant deviation from Hardy-Weinberg equilibrium was found in either population. The estimated moderate genetic differentiation and admixed genetic structure suggest recent lineage divergence and/or gene flow between populations. While robust evidence for a recent bottleneck was not obtained in our analyses, a clear signature of historical population contraction was detected in the central cone population. Our findings suggest that the skewed morph ratio strongly influenced the reproduction of small and isolated populations in the short-term, highlighting the vulnerability of distylous plant populations under ongoing anthropogenic pressure.


INTRODUCTION
Worldwide, numerous plants are already threatened by human-caused stress (e.g., habitat destruction) and climate changes (Jackson & Kennedy 2009). Among these plants, species having a sophisticated entomophilous breeding system such as distyly (heterostyly) are likely to be the most vulnerable to the detrimental effects of isolation and unreliable pollination service due to anthropogenic environmental alteration .
Distyly is a floral polymorphism, where populations have two floral morphs (a long-and short-styled morph; hereafter, referred to as the L-and S-morphs) that differ reciprocally in the heights of stigmas and anthers in flowers. Besides the morphological differences, distylous plants usually have a heteromorphic incompatibility system that prevents selfing and intramorph mating (Barrett 2002); only cross-pollination (i.e., legitimate pollination) between L-and S-morphs results in seed setting. In theory, such morphologically and physiologically disassortative mating between floral morphs generally leads to an equilibrium with equal morph ratios, as a result of negative frequencydependent selection and simple inheritance of distyly (Heuch 1979;Barret & Shore 2008). Accumulated evidence in distylous plants, however, has provided numerous examples of variation in population morph ratios (e.g., Kéry et al. 2003;Brys et al. 2008;Meeus et al. 2012). It has been advocated that floral morph bias can be governed by several factors, such as stochastic and deterministic events (Matsumura & Washitani 2000;Kery et al. 2003), maternal fitness differences between morphs (Hodgins & Barret 2006), and a combination of weak heteromorphic incompatibility and pollen limitation (Barret 1989;Brys et al. 2008).
Skewed morph ratios have often been found in small isolated distylous plant populations (e.g., Matsumura & Washitani 2000;Kery et al. 2003). Furthermore, it is well known that the deviation of morph frequencies from a 1:1 ratio can have negative reproductive and genetic consequences for populations. Indeed, skewed morph ratios result in the limited availability of compatible mates, which can contribute to reduced reproductive success (Kery et al. 2003;Wang et al. 2005;Pedersen et al. 2016) and increase the effects of genetic drift (Byers & Meagher 1992). Moreover, a combination of the loss of effective pollinators and the absence of a strict heteromorphic incompatibility system can increase inbreeding (selfing and biparental inbreeding) in morph-biased populations (Barret 1989;Wilcock & Neiland 2002;. Another consequence of skewed morph ratios is an increase in genetic drift and inbreeding, which can lead to a loss of genetic diversity Meeus et al. 2012) and may eventually result in the breakdown of distyly (Washitani 1996). Therefore, conservation studies on threatened distylous species which integrate ecological and genetic approaches are indispensable for assessing current status and predicting future extinction probability, as well as for planning effective conservation and/or restoration strategies .
The present paper investigated Primula reinii Franch. et Sav. (Primulaceae), an endemic primrose that inhabits mountainous regions of Japan. As in most primroses, the plant is a self-incompatible distylous species (Richards 2003). Because of its attractive and relatively large flowers with dwarf foliage, the primrose has suffered from anthropogenic activity (e.g. horticultural exploitation) in the wild. Based on the rarity and serious reductions in numbers and populations, the species is listed in the 'Vulnerable' category of the latest Japanese Red List (Ministry of the Environment 2019). Despite required effective management, especially for small populations exposed to ongoing human activities, little is known about their population status, reproductive success, or remnant genetic diversity. Simultaneously, this might provide an opportunity to study the immediate responses of a distylous plant population to demographic changes. In this study, we assessed the genetic and reproductive status of two neighborhood populations of P. reinii to provide baseline information pertinent to the conservation and preservation of this rare and endangered primrose. To discuss factors affecting reproductive success in natural habitats, we measured population size, morph frequency, stigmatic pollen deposition, and fruit-set within a population, and evaluated the genetic diversity and structure within and among populations.

Plant species and study site
Primula reinii is a diploid (2n = 24) perennial herb occurring as a chasmophyte on wet shaded rocky cliffs in the mountains (Image 1). The natural populations are rare and usually small and isolated from each other because of their narrow edaphic niche and low dispersal ability arising from the nature of the species. A single ramet produces one to two pink flowers from J TT mid-April to early May. Although flower visitors are not well-known in this species, its narrow corolla tube and recent pollinator observations in their related species (Yamamoto et al. 2018) imply that long-tongued flying insects, such as bumblebees and bee flies can be effective pollinators for the species. Under cultivation conditions, the generation time (between seed germination and first flowering) of the species is estimated to be 2-3 years (Yamamoto et al. 2017).
Fieldwork was conducted from April 2013 to October 2014. We selected P. reinii populations from two sites in the Hakone volcano within a special protected zone of the Fuji-Hakone National Park in central Japan (Fig.  1a). These sites were severely isolated by a volcanic landform, i.e., caldera (distance: ca. 7 km) (Fig. 1b). One was on Mt. Kintoki (KIN population,35.289'N,139.004'E,1,212m) and the other was on Mt. Komagatake (KOM population,35.228'N,139.021'E;1,275m). In both sites, the primrose occurs on rock outcrops near the mountain top in which populations are composed of scattered patches within an area of approximately 20 × 20 m 2 . Although a few smaller patches also located alongside the studied populations, these rim populations were sparse and separated geographically from the studied center population (>100 m). Thus, each studied population was considered to be one reproductive unit.
Currently, the Hakone volcano is an attractive destination to tourists. There is more than 10,000 people living within the caldera, and an approximate average of 60,000 tourists visit the area every day. Including Primula reinii, approximately 80% of endangered plants

J TT
found in the volcano are presumed to have decreased due to habitat destruction and horticultural exploitation (Osawa & Inohara 2008).

Measurements of basic population traits, pollen deposition, and fruit-set
During the flowering season, the number of all flowering individuals and flowers within each population, along with flower morph (L-or S-morph), were recorded. Whether the two morphs were equally frequent within a population (i.e., deviations from a 1:1 morph ratio) was investigated with Chi-square goodness-of-fit tests.
To elucidate the role of pollen limitation and morphratio variation on female reproductive success, we measured the stigmatic pollen load. In the natural fields, 20 fully-opened flowers were collected from each population in the afternoon and carefully transported to the laboratory the same day. In the laboratory, flowers were dissected, and stigma removed and mounted on a microscope slide in aniline blue staining solution (0.1% aniline blue in 0.1 M K 3 PO 4 ). Under a compound microscope (Olympus), we directly counted the number of legitimate (pollen from the opposite morph) and illegitimate (pollen from the same morph) pollen grains on stigmas of each floral morph based on pollen size differences (L-morph: 18.0 ± 0.1 μm; S-morph: 28.6 ± 0.1 μm; Y. Ojima, unpublished data).
At the beginning of the fruiting season (August), the population mean for fruit-set per flower was measured with the exception of some flowers that were used for the measurements of stigmatic pollen loads. Even if the fruit was set, the reproductive success will strongly depend on fruit predation (e.g., Matsumura & Washitani 2000;Yamamoto et al. 2013). Thus, we continued observations until October that was immediately before avulsion of the capsules.
All field surveys described above were conducted in 2013 and 2014 for KOM and KIN population, respectively.

Population sampling, DNA extraction, and microsatellite genotyping
In each population, 32 plants were sampled randomly (without distinction of floral morph types) from the entire area occupied by each population for genetic analysis. Leaf materials were collected and dried in silica gel. Before DNA extraction, leaves were homogenized with a disposable homogenizer (Biomasher 2; Nippi Co., Tokyo, Japan) to a fine powder. Total DNA was extracted from 40 to 80 mg silica-dried leaf tissue using the grassfiber filter method (Takakura 2011). The extracted DNA was dissolved in a TE solution and stored at 4 0 C until use.
After the preliminary marker screening, the genotypes of each individual were characterized following seven microsatellite markers that were originally developed for Primula sieboldii E. Morren: ga0161, ga0218, ga0580, ga0691, ga1140, Pri0141, and 2ca135 (Ueno et al. 2003(Ueno et al. , 2009Kitamoto et al. 2005). PCR amplifications and allele-size determination of fragment analysis were performed in accordance with the methods described by Yamamoto et al. (2017).

Population genetic analysis
For all seven microsatellite loci, the absence of linkage disequilibrium (LD) and the presence of null alleles were tested using Genepop v4.2 (Raymond & Rousset 1995). The LD test was verified using a Markov chain method with 1,000 dememorization steps, and 1,000 iterations per batch. Null allele frequencies were estimated by maximum-likelihood estimator based on the expectationmaximization algorithm (Dempster et al. 1977) with the default setting.
The following measures were calculated for each population: number of alleles (A), effective number of alleles (A E ), number of private alleles (A P ), expected heterozygosity (H E ), and inbreeding coefficient (F IS ). Deviations from the Hardy-Weinberg equilibrium were determined by the exact test and permutations. All measurements were calculated using GenoDive v2.0 (Meirmans & van Tienderen 2004). GenoDive was also used to compute the population's genetic differentiation pairwise F ST and G' ST indices (Hedrick 2005), and F ST was tested for significance using 10,000 permutations.
To estimate genetic structure of P. reinii populations in Hakone volcano, we used the model-based clustering method STRUCTURE 2.3.4 (Pritchard et al. 2000) and non-model-based method principal component analysis (PCA). STRUCTURE analysis was conducted for all samples across the two populations. Under an admixture model with correlated allele frequency, 20 independent simulations were run for each K (K = 1-5) with 5 × 10 5 Markov chain Monte Carlo (MCMC) steps and a burn-in period of 10 5 interactions. The most likely value of K was determined by the ΔK method (Evanno et al. 2005) with STRUCTURE HARVESTER 0.6.94 (Earl & vonHoldt 2012). CLUMPAK (Kopelman et al. 2015) was used to average the outputs from multiple STRUCTURE runs and produce the graphical results. The F value, the amount of genetic drift between each cluster and a common ancestral population, was also calculated for each cluster. The PCA analysis was performed using the package adegenet 2.0.1 (Jombart 2008) in R 3.5.2 (R core Team 2018).
To detect a genetic imprint of past population J TT bottlenecks, we first used the heterozygosity excess method (Cornuet & Luikart 1996) implemented within the program BOTTLENECK v1.2 (Piry et al. 1999). This method is suitable to detect very recent and less severe bottlenecks, and has low false positive error rates (Williamson-Natesan 2005). All simulations were done with mutation-drift equilibrium conditions (2,000 replicates) under the stepwise mutation model (SMM), infinite allele model (IAM), and two-phase mutation model (TPM: 70% SMM and 30% IAM). A two-tailed Wilcoxon signed-rank test was used to determine a significant excess of heterozygosity.
We also calculated the M-ratio (Garza & Williamson 2001) for each population using Arlequin (Excoffier et al. 2005). The M-ratio test is considered to have a greater detection power for ancient and moderate-to-severe population declines in comparison with the former method (Williamson-Natesan 2005). M-ratio represents the number of alleles relative to the range in allele sizes. After a severe bottleneck, the number of alleles should reduce faster than the allelic size range, which results in a reduced M-ratio (Garza & Williamson 2001). Thus, the magnitude of the decrease reflects the severity and duration of the reduction in population size, and generally an M-ratio <0.68 is indicative of the presence of a bottleneck (Garza & Williamson 2001).
Finally, we conducted a Bayesian demographic analysis using the R package, Vareff (Nikolic & Chevalet 2014). In contrast to the first two moment-based methods, this coalescent-based approach can examine temporal changes in the effective population size (Ne). The function Vareff simulates prior changes in the effective population size from microsatellite data by resolving coalescent theory and using an approximate likelihood MCMC (Nikolic & Chevalet 2014). After a series of preliminary runs, we used the prior parameter settings for each population (Table S1), following recommendations from Nikolic & Chevalet (2014). We set the estimated mutation rate to 5 × 10 -4 (Estoup et al. 2002) for all loci, and ran each analysis under a twophased mutation model with a proportion of 0.22 for multiple mutations (Peery et al. 2012), for 10 5 MCMC steps (NumberBatch = 1,000,000, LengthBatch = 10), sampling every 10 steps (SpaceBatch = 10) with an acceptance ratio of 0.25 (AccRate = 0.25), after burning of 10,000 steps. Estimations of sizes were searched for from sampling time to 5,000 and 500 generations ago.

Population traits and morph ratio
Each population trait is summarized in Table 1. A total of 72 flowering individuals and 99 flowers were found in the KIN population, whereas the KOM population had fewer (52 individuals and 69 flowers). The morph ratio in KIN did not deviate significantly from a 1:1 ratio (L-morph ratio = 0.54), even in the number of flowers (L-morph ratio = 0.52). In contrast, the number of flowering individuals of the L-morph was significantly higher than that of the S-morph in KOM (L-morph ratio = 0.65) and even higher for the number of flowers (L-morph ratio = 0.70).

Pollen deposition
Stigmas of both floral morphs received pollen grains in each population, but the numbers varied greatly between the individuals, ranging from zero to 321. In the KIN population, no differences in stigmatic pollen loads were detected between morphs (Fig. 2a). In addition, the proportion of deposited legitimate pollens was not significantly different between both morphs (Fig. 2c), while the proportions varied greatly among the L-morph stigmas in comparison with the S-morph. This is complemented by the result that the S-morph stigmas received significantly more legitimate pollen grains than the L-morph stigmas (Fig. 2b), implying that S-morph stigmas were pollinated more effectively than the opposite morph.

J TT
In contrast, in the KOM population, L-morph stigmas received a significantly greater number of pollen grains than the S-morph stigmas (Fig. 2a). After classifying pollen grains, however, we found no legitimate pollen grains loaded on the L-morph stigmas (Fig. 2b); that is, most L-morph stigmas were covered with a large quantity of illegitimate pollen.
Although several S-morph stigmas were legitimately pollinated, similar to other populations (Fig. 2c), there was no significant difference in the number of legitimate deposited pollen grains between the two morphs (Fig. 2b).

Fruit-set
At the population level, fruit-set ratio was much higher in the KIN population than in the KOM population (32.3 % and 14.3 %, respectively) (Fig. 2c). Within a population, both L-and S-morph scored comparable values in the KIN population (37.1 % and 26.7 %, respectively). In contrast, fruit-set of L-morph in the KOM population was less than half of the opposite morph (10.5 % and 27.3 %, respectively). We continued monitoring until October, but no evidence of fruit predation was found in either population, namely fruit-set was almost unchanged throughout the fruiting season.

Genetic diversity
LD between locus pairs was not significant. Although the frequencies of the majority of the null alleles were lower than 0.1, higher frequencies of null alleles were detected on 2ca135 and ga1140 loci in the KOM

J TT
population (0.227 and 0.114, respectively). As the presence of null alleles may affect the estimation of genetic diversity or differentiation, we excluded the two loci and repeated several analyses to compare results between seven and five microsatellites. This trial revealed no clear difference in the results based on 5 vs. 7 loci (Table S2). Thus, seven loci were used in all analyses described below.
Genetic diversity parameters for the two populations are presented in Table 2. In total, 63 alleles were amplified by seven microsatellite markers, with an average 9.0 alleles per locus. All diversity measurements were slightly higher in the KIN population (A = 6.3, A E = 3.3, A P = 2.9 and H E = 0.652) than in the KOM population  (b) (a) J TT (A = 6.1, A E = 2.6, A P = 1.9 and H E = 0.544). The inbreeding coefficient value (F IS ) was positive and comparable between the populations, but each value did not deviate significantly from zero.

Genetic structure and evidence a recent bottleneck
A moderate genetic differentiation was detected between populations (F ST = 0.115, P < 0.001; G' ST = 0.286). The STRUCTURE analysis based on the ΔK method indicated that ΔK was 462.5 for K = 2 and ΔK were <3 for other values of K. Therefore, the optimal ΔK for K = 2 showed that the best-fit model for the 64 sampled individuals of P. reinii revealed two clusters (Fig. 3a). Although several admixed individuals were found in each population, all samples formed a clear genetic structure between the two populations. The F values of clusters produced by STRUCTURE analysis were higher in the KOM population (F = 0.186) than in the KIN population (F = 0.086), indicating that the KOM population had undergone a larger genetic drift compared to that of the KIN population. In the PCA (Fig. 3b), the first two axes explained 17.0% and 10.7% of the variances in the experimental data, respectively. The results also distinguished the two populations, suggesting the existence of two genetic units corresponding to each population.
In BOTTLENECK analyses, the two-tailed Wilcoxon signed-rank test provided statistical support (P < 0.025) to the presence of a recent bottleneck in the KOM population under the SMM and TPM, whereas no  (Table 2). On the contrary, the M-ratio test indicated that both populations experienced a reduction in population size. The M-ratio values were 0.312 and 0.338 for the KIN and KOM populations, respectively (Table 2). A clear signature of historical population contraction was detected only in KOM via the third method, the Bayesian population demographic analysis. The bottleneck began approximately 1,000 generations ago (Fig. 4a), whereas a gradual decline was settled at least 100 generations ago (Fig. 4b). In contrast, the KIN population seems to have historically had a large constant population size (Fig. 4a); however, recent changes were unclear due to the broad confidence levels (Fig. 4b).

Reproductive status and genetic diversity
The observed low reproduction in KOM is congruent with reports that morph-biased populations experience reduced reproduction (Byers & Meagher 1992;Kery et al. 2003;Wang et al. 2005;Pedersen et al. 2016). Given that almost all stigmas were covered with L-morph pollen grains (Fig. 2), it is plausible that frequent self-or intra-morph (i.e., illegitimate) pollination had occurred among the KOM L-morphs. Therefore, our ecological data indicate that the low fruiting success in KOM L-morphs was caused by stigmatic clogging (Yeo 1975) as a consequence of the skewed morph ratio. Because L-morph flowers generally produce greater amounts of pollen grains than S-morph flowers (Richards 2003), it is apparent that the total pollen pool within KOM was occupied by a large amount of L-morph pollen. Similar to our results, previous studies in distylous plants showed higher female reproductive success in the relatively less abundant morph than the dominant morph (e.g., Wyatt & Hellwig 1979;Thompson et al. 2003;Wang et al. 2005;García-Robledo & Mora 2007). Thus, these results may demonstrate negative frequency-dependent patterns of reproductive success in the distylous primrose.
The indices of genetic diversity were relatively high and comparable between the two populations (Table 2), despite the skewed morph ratio observed. In addition, each population exhibited low F IS levels with no significant deviation from the Hardy-Weinberg equilibrium. These results allow for the conclusion that Primula reinii growing in the volcano had maintained sufficient genetic diversity as a result of outbreeding.
Overall, this study suggests the persistence of distylous self-incompatibility system in the P. reinii populations. Nevertheless, determination of the exact causes of floral morph bias in KOM was not possible based on the limited ecological and genetic data currently available. Because skewed morph ratios are often explained by several biotic and abiotic factors as discussed in the Introduction, there is a need for future studies investigating the ability of selfing and intra-morph mating, maternal fitness differences between morphs, pollinator assemblage, and population demography.

Genetic differentiation and structure
Our molecular analysis showed that genetic differentiation was moderate between the two populations (F ST = 0.115). Additionally, signs of genetic admixture between the populations were detected in PCA and STRUCTURE analyses (Fig. 3). There are at least two non-exclusive explanations for this: recent lineage divergence and gene flow. According to accumulated geographical surveys, Mt. Kintoki (locality of KIN pop.) and Mt. Komagatake (locality of KOM pop.) formed approximately ca 350-300 ka (Nagai & Takahashi 2008) and ca 27-20 ka (Kobayashi 1999;Nagai & Takahashi 2008), respectively. Formation of the central cone (i.e., Mt. Komagatake) clearly corresponded to the period of the last glacial maximum (LGM; ca 25-15 ka), suggesting that the KOM population was established at least after the last glacial period. The observed high F value (STRUCTURE analysis) and low private alleles in KOM may support a migration scenario that the population experienced a founder effect arising from a post-glacial refugial isolation and subsequent migration from the lowland of the caldera to the high-altitude areas of the central cone during the late Pleistocene and Holocene. Hence, it is plausible that the detected genetic admixture between populations suggests incomplete lineage sorting (i.e., sharing ancestral polymorphism between populations) due to recent lineage divergence.
Given the geographically close relationship between the populations (Fig. 1b), the presence of contemporary gene flow will also be taken into consideration. Because the two populations are severely isolated by a volcanic landform, gene flow mediated by pollen would be a plausible hypothesis. Moreover, in the flowering season we found claw marks, a useful indicator for the pollination services provided by bumblebees (Washitani et al. 1994), on the petals of each population. This may suggest that the bumblebees have a key role in pollination within the Primula reinii populations. Although bumblebees are known as strong-flying insects (e.g., Rao & Strange 2012), previous observations in other Primula species have demonstrated that pollen transfer by bumblebees J TT generally occurs within short distances (e.g., Ishihama et al. 2006). Therefore, we determined that the pollen flow between the populations might occur contemporarily but on very rare occasions. Nevertheless, deciding among the possible explanations for the genetic composition of the primrose in the Hakone volcano is difficult due to the weak evidence based on an insufficient number of loci.

Recent and historical demography
The two tests for a recent bottleneck yielded mixed results (Table 2). Based on the BOTTLENECK analysis, only the KOM population exhibited excess heterozygosity. In contrast, the M-ratio test supported a recent population size reduction in both populations. As mentioned above, however, because these inconsistent results might be attributed to the low statistical power of our sample size (e.g., number of loci or individuals), our results should be interpreted with caution. Nevertheless, such conflicting results often indicate the severity or timing of the reduction in population size (Williamson-Natesan 2005;Marshall et al. 2009;Padilla et al. 2015;Tóth et al. 2019), and were expected due to the differences in power detecting a bottleneck (Peery et al. 2012).
Considering the robust results in KOM, it is likely that the morph-biased population may have undergone more recent and severe bottlenecks in comparison with another population. In theory, the BOTTLENECK analysis can demonstrate population bottlenecks over a period of 0.2-4.0 Ne generations (Cornuet & Luikart 1996). Assuming for KOM population of Ne = 100 (Fig.  4) and a generation time of 2-3 years, it translates into approximately 50-1000 years before the present. On the other hand, a clear sign of recent (within 100 generations) population bottleneck was not found in the Bayesian demographic analysis (Fig. 4). Therefore, based on results from a series of demographic analyses, it is difficult to draw a definitive conclusion on whether recent bottlenecks occur or not, and thus, we defer a final conclusion until more genetic data are available in the future.
Contrary to this, the Bayesian demographic analysis provided strong evidence in support of a historical population bottleneck in KOM inhabiting the central volcanic cone. The first signs of population decline would have occurred 2-3 ka (assuming a generation time of 2-3 years). This timeframe post-dates a climatic warming, known as the Jomon optimum transgression, that occurred approximately 6ka, implying that historical population bottlenecks were likely due to volcanic activities as opposed to climatic events. According to geological records, the last major eruption of the volcano was from the central cone in 2.7-2.9 ka (Kobayashi 1997;Kobayashi et al. 2006), and intermittent phreatic eruptions continued until present-day.
Although speculative, these evidences may support the idea that the historical population declines experienced by the KOM could have been associated with repeated eruptive activities in the central cone. Perhaps, the detected recent bottlenecks in KOM are caused by eruptive activities rather than human activities.
On the other hand, the estimated effective population size in the KIN population inhabiting the somma mountains was large and constant in the long term, suggesting that the population has been maintained without suffering from volcanic eruptions occurring in the central cone. Further studies for the lineage divergence and demographic history of P. reinii in this region, using more informative datasets (e.g., single nucleotide polymorphisms), will be valuable because volcanism is one of the key abiotic factors in the plant's diversification and distribution in Japan (e.g., Yoichi et al. 2017;Nagasawa et al. 2020), located in the Pacific Ring of Fire.

Implication for conservation
Our study suggests that morph imbalances are striking effects on the reproduction of P. reinii population in the short-term. Accordingly, a measure of morph ratio should be given top priority in conservation management of the species, and enhancement of habitat monitoring should be considered as in situ managements to protect remnant individuals and to maintain optimum morph frequencies from horticultural exploitation. Considering the observed negative frequency-dependent patterns of reproductive success, if heteromorphic selfincompatibility is totally strict in P. reinii, the skewed morph ratio in KOM may be improved in the future when regeneration is successful. However, the exact breeding system of the species remains poorly understood. Therefore, in addition to other examinations (e.g., the germination requirements and the effect of storage time of seeds) towards a future ex situ conservation strategy, the levels of within morph fertility and selfing ability should be resolved immediately to evaluate the medium-to long-term risk of extinction in the remnant populations across species distribution ranges.
The two surveyed populations in the Hakone volcano were distinguished by two genetic clusters, suggesting that each population should be divided into a different management unit to maintain evolutionary distinctiveness and ecological viability (Moritz 1994;Frankham et al. 2002). The moderate genetic J TT differentiation and the presence of large amounts of private alleles between the populations highlight this suggestion; thus, artificial inter-population crossing should be avoided in this case. Nevertheless, the lack of samples from other parts of the volcano will influence the estimated genetic structure. Thus, an exhaustive population sampling, including other remnant small population, is required to elucidate the genetic structure and demographic history of P. reinii occurring in the Hakone volcano as is also needed for planning conservation strategies.
To our knowledge, this is the first conservation genetics study on threatened plants in the Hakone volcano, which harbors approximately 1,800 plant species (Tanaka 2008). Thus, the results discussed here will be useful for designing both in situ and ex situ conservation strategies for P. reinii as well as other plants inhabiting the volcano and shed light on the instability of plant populations due to the impacts of volcanism and human activities. Our study highlights the importance of studies in conservation, integrating ecological and genetic approaches to accurately assess the population status of endangered species and draw up effective conservation strategies.  Table S1. Prior parameter settings for each population using VarEff software. www.threatenedtaxa.org

J TT
The Journal of Threatened Taxa (JoTT) is dedicated to building evidence for conservation globally by publishing peer-reviewed articles online every month at a reasonably rapid rate at www.threatenedtaxa.org. All articles published in JoTT are registered under Creative Commons Attribution 4.0 International License unless otherwise mentioned. JoTT allows allows unrestricted use, reproduction, and distribution of articles in any medium by providing adequate credit to the author(s) and the source of publication.