Genetic diversity and population structure of Euscaphis japonica, a monotypic species

Plant material

Based on the samples collected from the previous phenotypic diversity of E. japonica (Sun et al., 2019), we have expanded the scope of field investigation and collected more samples from August 2016 to March 2017. During field observation, we found that E. japonica had two major categories, deciduous and evergreen, and we sampled these two types separately. All observed individuals were sampled, and the latitude, longitude, altitude, and specific habitats were recorded. We collected 83 samples across nine populations in Wuyi Mountain (WYC1, WYC2, WYC3, WYC4, and WYL), Taimu Mountain (TML), Daiyun Mountain (DYC and DYL), and Xishui National Nature Reserve (ZYL) (Fig. 1). The Wuyi Mountain range is about 550 kilometers long from the northeast to southwest, so five populations were collected. We collected and dried 10 to 15 fresh leaves from each sample using silica gel in sealed bags (Chase & Hills, 1991). The dry leaves were stored at 4 °C in a refrigerator to be used for ISSR analysis. The characteristics of the different population localities and morphologies are shown in Table 1.

DNA isolation and purification

We extracted genomic DNA from the E. japonica leaves using a modified acetyl trimethyl ammonium bromide method (Uddin et al., 2014). Fresh leaves (50 mg) without veins were ground into a fine powder using liquid nitrogen and a mortar. We then quickly transferred the powder to 2-mL microcentrifuge tubes with 700 µL of extraction buffer (1.4 M NaCl, 100 mM Tris HCl (pH 8.0), 20 mM EDTA, 1% PVP, and 1% β-mercaptoethanol) that had been preheated at 65 °C. This mixture was incubated at 65 °C for 40 min with gentle shaking. After cooling for 10 min, we added chloroform: isopentyl alcohol mixture (24:1, equal volume: 700 µL) and centrifuged it at 12,000 rpm for 10 min. The supernatant was transferred into a clean 2-mL microcentrifuge tube, and two volumes of absolute ethanol were added and gently mixed using inversion. The DNA pellet was washed twice with 70% (v/v) ethanol and dissolved in TE buffer (pH 8.0) after air drying. We removed RNA contamination by adding 2 µL of RNase (10 mg/mL) to each tube, incubated the tubes for 30 min at 37 °C, and then added 400 µL of chloroform: isopentyl alcohol (24:1). Samples were centrifuged at 12,000 rpm for 10 min at 4 °C and the aqueous layer was transferred into a new tube. We added two volumes of chilled ethanol and 1% NaAc (3 mol/L, pH 5.2) to precipitate the DNA. The DNA was washed with 70% ethanol and dissolved in 100 µL of TE buffer after air drying. Each sample was diluted to 80 ng/mL using TE buffer and was stored at 4 °C for further ISSR analysis.

We measured the concentration of recovered DNA using NanoDrop™ spectrophotometry (Thermo Fisher Scientific, Waltham, MA, USA) at 260 nm and 280 nm, and verified the purity of the DNA using 0.8% (w/v) agarose gel. SYNGENE (BioRad, Hercules, CA, USA) was used to visualize and photograph the gel under a UV transilluminator.

ISSR amplification

PCRs contained 1 µL of 80 ng template DNA, 1 µL of 0.5 µmol/L primer, 0.4 µL of 0.2 µmol/L dNTP, 2 µL of 2.5 mmol/L MgCl₂, 0.3 µL of 1.5 U Taq enzyme, 2 µL of 10× PCR buffer, and sterile distilled water, for a final volume of 20 µL. ISSR primers were synthesized by Sangon Bioengineering (Shanghai, China) and all reagents were purchased from Solarbio (Beijing, China). We conducted PCR amplification using a Veriti™ 96-Well Thermal Cycler (Applied Biosystems, Foster City, CA, USA) as follows: initial denaturation at 95 °C for 5 min, 35 cycles at 94 °C for 1 min, at 55 °C for 30 s, and 72 °C for 1 min; final extension at 72 °C for 7 min, and the tubes were subsequently maintained at 4 °C before analysis. All ISSR primers had been initially tested, and 12 primers amplified DNA with polymorphic bands (Table 2). PCR products were analyzed on a 1.2% (w/v) agarose gel in 1 × TBE buffer and were visualized using the Gel Doc system (BioRad).

Statistical analysis

We scored amplification based on the amplicon bands from the gel photographs. The presence of bands at an amplicon level was scored as 1 and the absence was scored as 0, which we used to calculate the raw data of the 0–1 matrix (Uddin et al., 2014). We used POPGENE 1.3.1 to analyze various genetic parameters, such as the number of alleles (Na), the number of effective alleles (Ne), Nei’s genetic diversity (H), Shannon’s information index (I), the percentage of polymorphic loci (PPB), the level of gene flow (Nm), and the Dice’s similarity coefficient (GS) (Yeh et al., 1997). We calculated polymorphism information content (PIC) using the following formula: PIC = 2 Pi (1–Pi), where Pi is the frequency of polymorphic band occurrence in different primers (Khaleghi et al., 2017). The unweighted pair group method with the arithmetic average (UPGMA) and NTSYS 2.1 were used for dendrogram construction. To verify the UPGMA clustering results, we performed a cophenetic correlation analysis using NTSYS 2.1. Principal component analysis (PCA) was performed to show the E. japonica multiple dimensional distributions in a scatter plot (NTSYS 2.1) (Wang et al., 2009). We used GenAlex 6.5 to calculate the variance components within and between populations (Peakall & Smouse, 2012). The genetic structure of the E. japonica population was determined using STRUCTURE 2.3.1 based on Bayesian model-based clustering (Pritchard, Stephens & Donnelly, 2000). The procedure was carried out by selecting the correlated allele frequencies among populations and the admixture model. The presumed populations (K) denoted from 1 to 13 and estimated 20 independent runs for each K. The operating parameter was a burn-in period of 100,000 and 100,000 Markov Chain Monte Carlo replicates after burn-in. The optimal number of clusters was identified using Structure Harvester (http://taylor0.biology.ucla.edu/struct_harvest/) (Earl & vonHoldt, 2012).

ISSR polymorphism and genetic diversity

We amplified a total of 122 identifiable bands, 121 of which were polymorphic bands. The polymorphic ratio (PPB) was 99.18%. The size of the PCR products ranged from 250 to 2,000 bp. The largest number of bands (14) was obtained using the UBC 808 primer, while only six bands were obtained using the UBC 890 primer. The highest PIC value was 0.62 (UBC 809), the lowest PIC value was 0.48 (UBC 856), and the mean PIC per primer was 0.56 (Table 2). The high PPB and PIC values suggest that E. japonica has abounding genetic information, and that ISSR markers are appropriate to use in the genetic diversity analysis of E. japonica populations

Table 3 shows the statistical analysis of the E. japonica genetic parameters when using GenAIEx 6.5. At the species level, the average Ne was 1.4974, H was 0.3016, and I was 0.4630. At the population level, the ranges of H and I across different populations were 0.0994–0.3688 and 0.1594–0.3971, respectively. Five populations had H values greater than 0.2: WYC2 (0.2607), WYC3 (0.2512), WYL (0.2360), DYL (0.2301), and WYC4 (0.2141). Four populations had I values greater than 0.3: WYC2 (0.3971), WYC3 (0.3751), ZYL (0.3548), and DYL (0.3471). The five populations with the largest Ne values were WYC2 (1.4090), WYC3 (1.4276), ZYL (1.4034), DYL (1.3857), and WYC4 (1.3610). The five populations with the largest population polymorphic rates (PPB) were WYC (86.07%), WYC4 (72.13%), WYC3 (70.49%), ZYL (68.85%), and DYL (67.21%). Considering each genetic diversity index collectively, including the H, I, and Ne values, the WYC2 population had the highest genetic diversity, followed by WYC3 and WYC4.

Genetic structure

The UPGMA tree, based on individual samples, consisted of four major clusters. Most ZYL population samples and all WYL samples were grouped into cluster A. The remaining ZYL samples, all TML samples, and DYL samples were grouped into cluster B. Cluster C contained WYC3 and WYC4 samples. Cluster D included WYC1, WYC2, and DYC samples. Clusters A and B consisted of deciduous E. japonica, and clusters C and D consisted of evergreen E. japonica (Fig. 2A). We based the PCA on genetic similarities to determine the genetic relationships between E. japonica populations and to create a two-dimensional display of these relationships. The PCA of the cumulative data grouped the nine populations into clusters I, II, and III, where cluster I included deciduous E. japonica and clusters II and III included evergreen E. japonica. Cluster II contained WYC1, WYC2, and DYC samples, and cluster III included WYC3 and WYC4 samples (Fig. 2B). The STRUCTURE analysis results indicated that the number of optimal clusters was four (Figs. 2C and 2D) based on the maximum delta K = 4. Therefore, we divided the populations into four clusters (Fig. 2E): ZYL, TML, WYL, and DYL (purple); WYC2 (green); WYC1, DYC, and some WYC4 (blue); and the remaining WYC4 and WYC3 (yellow). The samples in the purple groups are deciduous E. japonica, and the samples in the green, blue, and yellow groups are evergreen E. japonica. Interestingly, the great isolation of population ZYL is not reflected in an unusual level of genetic distance (Figs. 2A and 2E). However, this population forms a distinct entity in the dendrogram. These results showed that the genetic differentiation between different geographic populations is not obvious. Furthermore, from Fig. 2E, there was little genetic infiltration between the evergreen and deciduous populations. The lower average gene flow (Nm) of E. japonica (Nm = 0.7804) also confirms the previous results (Table 2). In addition, we found that the deciduous (Ne = 1.450, H = 0.271, I = 0.413, PPB = 84.43%) and evergreen E. japonica had rich genetic variation. (Ne = 1.474, H = 0.290, I = 0.446, PPB = 97.54%) (Table 3).

Genetic differentiation

The Dice’s genetic similarity coefficient ranged from 0.729 to 0.938. Two populations had relatively high Dice’s similarity coefficients, indicating a close genetic relationship and small genetic differences between the two populations. The Dice’s similarity coefficients were higher between WYC1and WYC2 (0.954), WYC2 and WYC4 (0.921), WYC3 and WYC4 (0.936), and DYL and ZYL (0.938), while the Dice’s similarity coefficients between DYC and TML (0.729), DYC and ZYL (0.793), WYL and TML (0.733), WYL and DHC (0.764), and WYL and WYC1 (0.766) were lower (Table 4). Genetic differentiation coefficient (Fst) values greater than 0.25 indicates greater genetic differentiation across populations (Wright, 1978). The Fst values for WYC1 and WYC2 (0.0883), WYC1 and WYC4 (0.2249), WYC2 and WYC4 (0.2070), WYC2 and WYC3 (0.2350), WYC3 and WYC4 (0.1509), and WYC3 and DYL (0.2172) populations were less than 0.25, and the Fst values for other populations were all greater than 0.25 (Table 4). Overall, there was a high similarity coefficient and a low Fst value across deciduous or evergreen populations in the same mountain, and there was a lower similarity coefficient and high Fst value between deciduous and evergreen populations.

The AMOVA analysis showed that 30% of total genetic variation existed across populations and 70% existed within populations (Table 5). The genetic variation within deciduous and evergreen populations accounted for 64% and 77% of the total variation, respectively. The genetic variation across different deciduous and evergreen populations accounted for 36% and 23% of the total variation, respectively. These results showed that for both deciduous E. japonica and evergreen E. japonica, most genetic variation was from within populations rather than from among populations.

Genetic diversity and differentiation

We determined polymorphisms in nine E. japonica populations from China by estimating the genetic variability across different populations. Compared to previous reports, E. japonica’ s average Ne, H, and I values (Ne = 1.4974, H = 0.3016, I = 0.4630) were higher than those observed in Staphylea bumalda (Staphylea, Staphyleaceae) (Ne = 1.341, H = 0.227, I = 0.370) (Chen, Wang & Wang, 2014). These results indicate that E. japonica has high genetic variation at the species level. There are two types of E. japonica, deciduous and evergreen E. japonica. We evaluated the genetic diversity of E. japonica, deciduous E. japonica, and evergreen E. japonica, and found that they all have abundant genetic variation, and the genetic variation within the population was much greater than that between populations. Numerous previous studies have shown that most woody species have more variation within individual populations than across different populations (Hamrick, Godt & Sherman-Broyles, 1992; Wei et al., 2012). Woody species with long-lived, outcrossing breeding systems and animal seed dispersal have higher variations within populations than across populations (Albert, Raspe & Jacquemart, 2005; Schaal et al., 2010; Ramírez-Valiente, Valladares & Aranda, 2014). E. japonica is a long-lived plant with a mixed mating system. Bee is the primary pollinator for E. japonica, meaning that pollen transfer can occur between adjacent populations (Sun et al., 2017), such as the WYC1 and WYC2 populations, WYC3 and WYC4 populations located in the Wuyi Mountain. The red pericarp of E. japonica attracts birds to eat the fruit that helps with seed dispersal. The lower interpopulation genetic diversity of E. japonica, deciduous E. japonica, and evergreen E. japonica and higher intrapopulation genetic diversity may be the results of its long life, mixed mating, and seed dispersal systems.

The lower Dice’s similarity coefficient and higher Fst value between the deciduous and evergreen populations confirmed that the deciduous and evergreen E. japonica experienced differentiation. Frequent gene flow can prevent genetic drift and reduce genetic differentiation (Tremblay & Ackerman, 2001). The STRUCTURE analysis and the lower average gene flow indicated that there was less gene flow between the evergreen and deciduous populations. According to our long-term field observations, the florescence of deciduous E. japonica was observed from April to May, and the florescence of evergreen E. japonica was observed from May to June (Sun et al., 2017; Sun et al., 2019). Although bees are the primary pollinators, there is a small possibility of genetic interactions between the deciduous and evergreen populations in the same Mountain, such as DLY and DYC populations in Daiyun Mountain, and WYC1 and WYL populations in Wuyi Mountain. We speculate that little gene exchange between deciduous and evergreen E. japonica may be the result of the differentiation between these two type populations.

Implications for E. japonica conservation

The deciduous E. japonica is mainly distributed in southern China, Japan, Korea, and the evergreen E. japonica is distributed in Fujian, Jiangxi, Guangxi provinces in southern China (Li et al., 2018; Sun et al., 2019). According to our field observations, the influence of human activities has been the main challenge affecting the survival of E. japonica populations (Sun et al., 2019). Therefore, based on the genetic diversity information of the population, an effective E. japonica conservation strategy is proposed to create favorable conditions for the innovation and effective utilization of E. japonica germplasm resources. We detected the highest genetic diversity in the WYC2 population (PPB = 86.07%, H = 0.2607, I = 0.3971), followed by WYC3 (PPB = 70.49%, H = 0.2512, I = 0.3751). The WYC2 and ZYL populations were located in the Jiangshi Nature Reserve and the Xishui Natural National Reserve, respectively, and the DYL population was located in sparsely populated natural villages. These habitats are well-preserved without frequent human interventions. Based on our UPGMA clustering, genetic diversity index, and population location, we suggest that the ZYL, DYL, WCY3, and WCY2 populations in clusters A, B, C, and D should be conserved in-situ. In-situ management of genetic resources can ensure that the majority of extant variation is preserved (Asddisalem et al., 2016; Negri & Tiranti, 2010).

Furthermore, we should work on recovering populations with low diversity, destroyed habitats, and small population sizes. TML and DYC habitats were shrinking due to human interventions. The WYC1 population is located near a city, and the WYC3 population was transplanted from local forests to forest farms because of habitat destruction. The populations with the lowest genetic diversity were WYL (PPB = 20.49%, H = 0.0803, I = 0.1174) and DYC (PPB = 21.31%, H = 0.0738, I = 0.1115) population. Therefore, we recommend an ex-situ conservation strategy of these populations. In ex-situ conservation, genetic resources are preserved outside their natural habitats in facilities such as seed banks and botanical gardens. These strategies aim to preserve genetic material in collections (Brush, 2000; Negri & Tiranti, 2010). Li & Pritchard (2009) and Richards et al. (2010) confirmed that seed collections are a useful and effective way to maintain the size of most ex-situ populations, and that vegetative propagation can rapidly and effectively expand ex-situ population sizes (Li et al., 2018). In natural E. japonica populations, collecting and storing seeds are of significantly important, especially for populations with low genetic diversity and severely damaged habitats, such as TML, WYL, and WYC1. Other populations with very small individuals, such as DYC, can benefit from vegetative propagation and in situ conservation of genetic resources.