Genetic Diversity of Paeonia rockii (Flare Tree Peony) Germplasm Accessions Revealed by Phenotypic Traits, EST-SSR Markers and Chloroplast DNA Sequences

: Research Highlights: This study, based on the first collection of cultivated Paeonia rockii (flare tree peony, FTP) germplasm across the main distribution area by our breeding desires, comprehensively evaluates these accessions by using phenotypic traits, expressed sequence tag (EST)-simple sequence repeat (SSR) markers and chloroplast DNA sequences (cpDNA). The results show that these accessions collected selectively by us can represent the genetic background information of FTP as a germplasm of tree crops. Background and Objectives: FTP has high cultural, ornamental and medicinal value traditionally, as well as recently presenting a significance as an emerging edible oil with high α -linolenic acid contents in the seeds. The objectives of this study are to reveal the characteristics of the genetic diversity of FTP, as well as to provide scientific suggestions for the utilization of tree peony breeding and the conservation of germplasm resource. Materials and Methods: Based on the phenotypic traits, EST-SSR markers and chloroplast DNA sequence variation, we studied the diversity of a newly established population of 282 FTP accessions that were collected and propagated by ourselves in our breeding project in recent years. Results: (1) There was an abundant variation in phenotype of the accessions, and the phenotypic variation was evenly distributed within the population, without significant hierarchical structure, (2) the EST-SSR data showed that these 282 accessions had relatively high genetic diversity, in which a total of 185 alleles were detected in 34 pairs of primers. The 282 accessions were divided into three distinct groups, and (3) the chloroplast DNA sequences (cpDNA) data indicated that these accessions had a higher genetic diversity than the population level and a lower genetic diversity than the species level of wild P. rockii , and the existing spatial genetic structure of these accessions can be divided into two branches. Conclusions: From the results of the three analyses, we found that these accessions can fully reflect the genetic background information of FTP germplasm resources, so their protection and utilization will be of great significance for genetic improvement of woody peonies.


Introduction
Tree peony or Mudan, known as "the king of flowers" in China, has high cultural, ornamental and medicinal value traditionally, as well as recently presenting a significance as an emerging edible oil with high unsaturated fatty acid (>90%) and α-linolenic acid (>40%) contents in the seeds [1,2].
Tree peony belonging to the section Moutan in the genus Paeonia and all wild species are endemic to China. It has been cultivated since the Tang Dynasty (618-906 AD) and is now widespread in temperate regions of China and other countries. At present, there are more than 1500 germplasm accessions of tree peony in the world, included in about 17 cultivar groups. In China, long-term and repeated domestication and selection breeding have formed ten groups, about 1000 germplasm accessions, with diverse genetic background [3].
As the representative of the most-widely cultivated tree peonies, Paeonia suffruticosa Andrews or the common tree peony (CTP), includes the germplasm accessions traditionally cultivated in China and Japan, while the germplasm accessions originated from P. rockii or the flare tree peony (FTP) have established a distinct group from CTP and are known well, with a colored flare at the base of petals [3]. Compared to CTP, FTP is more fertile in most cases as well as more resistant to stress conditions like coldness, drought and poor soil. FTP is distinct from CTP in many other features, like more fragrant blossoms, and taller and stronger plants with more longevity (Figure 1) [1,4]. Since many studies [1,[5][6][7][8][9][10][11] have confirmed that FTP is a valuable genetic resource of tree peonies, greatly promising in the breeding and industry of peonies, the further in-depth study of FTP has become obviously required for the utilization and protection of FTP as a germplasm resource of tree crops. The analysis of genetic diversity can provide beneficial data for collection of germplasm resources, breeding application and researches on origin and evolution [12]. Phenotypic variation analysis is one of the most important, basic, intuitive but indispensable methods for the development and application of plant germplasm resources. In tree peonies (Paeonia Sect. Moutan), researchers have focused on the genetic diversity of phenotypic traits in different species and cultivated germplasm, including P. delavayi [13], P. ostii [14,15], P. delavayi var lutea [16], P. suffruticosa [17][18][19] and so on. In recent years, the emergence of molecular markers has greatly accelerated the development of genetics and genomics and is expected to be a powerful tool to accelerate the research of tree peony breeding. Researchers had used random amplification polymorphic DNA (RAPD), inter-simple sequence repeat (ISSR), sequence-related amplified polymorphism (SRAP), conserved DNA-derived polymorphism (CDDP), inter-primer binding site (iPBS), target region amplified polymorphism (TRAP), expressed sequence tags microsatellite markers (EST-SSR), chloroplast DNA sequences (cpDNA) and EST-SSR with cpDNA to analyze the genetic diversity of P. suffruticosa [20][21][22][23][24][25][26], P. x yananensis [27], P. rockii [28,29], P. delavayi [30], P. ludlowii [31], Paeonia subsect. Delavayanae [32,33], P. decomposita [34], P. jishanensis [35], P. ostii [36] and P. quii [37]. In P. rockii, Yuan et al. [28,29] used cpDNA and SSR to analyze the genetic diversity in 20 wild populations, showing that wild P. rockii had high genetic diversity at the level of species and can be divided into three genetically distinct clusters corresponding to three geographically distributed regions. Xu et al. [37] used cpDNA and EST-SSR to analyze the genetic diversity of P. rockii, P. jishanensis and P. qiui to discover that the genetic diversity of P. rockii among the three tree peony species was relatively low. In the population of 462 cultivated P. rockii seedlings, Wu et al. [38] used 40 EST-SSR markers to reveal a medium genetic diversity and the genetic structure composed of three subgroups. These studies provided us with very useful references to carry out further research on the genetic diversity of tree peony germplasm.
With high polymorphism, codominance and ease of operation, EST-SSR molecular markers have been applied more and more in the study of genetic diversity [39]. Chloroplast genome is maternal lineage and retains ancestral patterns of genetic diversity longer than nuclear DNA [40]. Moreover, in the cross-breeding system, the effective population size of cpDNA is only half of that of nuclear DNA, which may lead to higher genetic differentiation and more significant genetic structure [41]. Chloroplast non-coding region sequences have higher mutation rates than coding region sequences, which can be used for population genetics research [42,43].
Based on the phenotypic traits, EST-SSR markers and chloroplast DNA sequence variation, this paper reports the research of genetic diversity in the newly established population of 282 FTP accessions that were collected, propagated and mostly named as cultivars by us in our breeding project in recent years. Our objectives are to reveal the characteristics of genetic diversity of these accessions, as well as to provide scientific suggestions for the utilization in peony breeding and the conservation of FTP germplasm resource.

Plant Materials and DNA Extraction
Germplasm resources of FTP had been selectively collected by our breeding desires across the main cultivation area of the FTP, Lanzhou city, Zhang county, Longxi county, Lintao county and Linxia county of Gansu province in Northwest China since 1990, propagated by grafting in Guose Peony Garden of Yan Qing district, Beijing, China (40°46′ N, 116°07′ E) since 2011 ( Figure 2). The plants of the 282 accessions used in this study are randomly planted in the same nursery field and maintained under the same condition of water and fertilizer throughout the year. In March 2018, the young leaves of all accessions were individually collected and dried in silica gel, from which total genomic DNA was extracted by using a DNAsecure plant kit (Tiangen Biotech, Beijing, China). The extracted DNA was detected by electrophoresis using 2% agarose gels and a UnicoUV-visible Spectrophotometer (Agilent, Palo Alto, CA, USA), respectively. The samples were diluted to 20-30 ng/μL with deionized water and then stored in the refrigerator at −20 °C.
A total of 282 accessions, of those collected above, were used for EST-SSR analysis, and the phenotypic traits were investigated only in the completely matured 15-year-old plants of 180 accessions. The color of flowers is genetically stable for every specific accession and in the studied accessions can be classified into five groups: 89 white, 33 pink, 32 red, 116 purplish red and 12 purplish black. By this color group, we selected the accessions for cpDNA sequencing, including 24 white, 10 pink, 8 red, 35 purplish red and 3 purplish black accessions, with the proportions of 26.97%, 30.30%, 25.00%, 30.17% and 25.00%, respectively (Supplementary Table S1).

Measurement of Phenotypic Traits
We recorded 28 phenotypic traits according to the measurement standard (Table 1) during the flowering time from April to May and the seed ripening time from August to September in 2019. In April to May, three plants were randomly selected for each accession as three replicates and three branches with flowers were randomly selected for each plant. The number of flowers, plant height, crown breadth and tiller number of three plants were measured first, then nine flower diameters were measured, and nine carpel numbers were counted. A random outer petal was selected from each flower to measure the length and width of each petal and flare. Then, the number of petals of each flower was counted. Finally, the second or third compound leaf from bottom to top of each branch was selected for leaf trait investigation. From August to September, we chose the three plants and three branches of each plant that we had measured in the blooming period. First, we counted the fruit number of three plants and fruit weight per plant. Then, we measured individual fruit weight, number of seeds per fruit, seed weight per fruit and effective carpel number. Lastly, we randomly selected a single carpel of each follicle to measure fruit length, fruit width, fruit height and pericarp thickness. Petal width The width of the outer petal (cm) 3 Flare length The length of flare (cm) 4 Flare width The width of flare (cm) 5 Flower diameter The maximum width of flower at full bloom (cm) 6 Petal number The number of petals in the whole flower (score) 7 Carpel number The number of carpels in the whole flower (score) Branch and leaf trait 8

Plant height
The aboveground height of the whole plant (cm) 9 East-west crown breadth The width of a plant in the east-west direction (cm) 10 North-south crown breadth The width of a plant in the north-south direction (cm) 11 Fruit number The number of follicles per plant (score) 12 Flower number The number of flowers per plant (score) 13 Tiller number The number of branches at the base of a plant (score) 14 Compound leaf length The length from tip of compound leaf to petiole (cm) 15 Compound leaf width The length at widest part of compound leaf (cm) 16 Petiole length The length of petiole (cm) 17 Leaflet number The number of leaflets in a compound leaf (score) 18 Terminal leaflet length The length of terminal leaflet (cm) 19 Terminal leaflet width The width of terminal leaflet (cm) Fruit trait 20 Fruit weight per plant The weight of fruits per plant (g) 21 Individual fruit weight The weight of one fruit (g) 22 Number of seeds per fruit The number of seeds per fruit (score) 23 Seed weight per fruit The weight of seeds per fruit (g) 24 Effective carpel number The number of carpels that produce seeds (score) 25 Fruit length The length of a single carpel of mature follicles (mm) 26 Fruit width The width of a single carpel of mature follicles (mm) 27 Fruit height The height of a single carpel of mature follicles (mm) 28 Pericarp thickness The thickness of pericarp of mature follicles (mm)

Chloroplast DNA Sequences
Since the same chloroplast gene fragment has different evolutionary rates in different plants, the optimal chloroplast fragment should be selected for different species [46][47][48]. From the reported chloroplast fragments, three pairs of chloroplast gene spacer fragments, petB-petD, accD-psaI and psbE-petL, which can be amplified stably and have high polymorphism, were screened out on the basis of existing common primers [26,49]. Primer information is shown in Table 3. Three chloroplast DNA primers were selected to access the genetic diversity of 80 accessions that can represent the typical flower colors from 282 accessions. The PCR amplification reaction was conducted in a 50 μL solution, including: 25 μL 2× Power Taq PCR Master MIX (Aidlab Biotechnologies, Beijing, China), 2.5 μL and 10 μmol/L each of forward and reverse primer, 2 μL and 20-25 ng/μL genomic DNA and 18 μL ddH2O. The PCR procedure was as follows: 5 min at 94 °C, followed by 35 cycles of 30 s at 94 °C, 30 s at 52-54 °C and 50 s at 72 °C, and lastly, 7 min at 72 °C. The amplified fragment results were detected by capillary electrophoresis using an ABI3730xl DNA Analyzer (Applied Biosystems, Carlsbad, CA, USA).

Analysis of Phenotypic Data
The average values per plant were used in the statistical analysis and the data were analyzed using the statistics software SPSS version 18.0 (IBM Inc., Chicago, IL, USA) [50]. One-way analysis of variance (ANOVA; with post hoc Duncan' s multiple range test with a probability p < 0.05 and p < 0.01) was used to analyze the 28 traits' variation and Pearson's correlations between traits were calculated. Prior to the ANOVA and Pearson's correlation analysis, all data were tested for normality with the Shapiro-Wilk W test and for homogeneity of variance with the Levene's test, and the nonnormal data were logarithmic transformed. We adjusted the p-value for Pearson's correlation using the Benjamini-Hochberg (BH) false discovery rate (FDR) correction. Coefficient of variation (%) = (standard deviation/mean) × 100. Principal component analysis (PCA) was performed for 28 traits. Maximum deviation method was used for factor rotation, and the principal components were extracted according to Kaiser criterion (characteristic root > 1). Then, with the first principal component (PC-1) and the second principal component (PC-2) as the main coordinates, the distribution diagrams of 28 quantitative traits and 180 accessions were plotted, respectively.

Analysis of EST-SSR Data
GeneMarker V2.2.0 was used to analyze the capillary electrophoresis data. GenAIEx version 6.5 [51] was used to calculate the following indicators: Number of Different Alleles (Na), Number of Effective Alleles (Ne), Shannon's Information Index (I), Observed Heterozygosity (Ho), Expected Heterozygosity (He), Inbreeding coefficients (FIS) and Nei's genetic diversity (GD). Polymorphism information content (PIC) was calculated for each locus by using a Microsatellite Toolkit. GENEPOP version 4.2 [52] was used to detect whether microsatellite loci deviated from the Hardy-Weinberg equilibrium (HWE). Additionally, a Neighbor-Joining phylogenetic tree based on Nei's unbiased genetic distance was constructed with MEGA-X [53]. Finally, Principal Coordinates Analysis (PCoA) was carried out using GenAIEx version 6.5 based on Nei's genetic distance of all accessions.

Analysis of Chloroplast DNA Sequences
MAFFT version 7 [54] was used to conduct multi-sequence alignment and adjustment. SequenceMatrix 1.7.8 [55] was used to splice three chloroplast non-coding regions, accD-psaI, psbE-petL and petB-petD, into a sequence. Values of haplotype diversity (Hd), haplotype number (Hap) and nucleotide diversity (Pi) were calculated by DnaSP 5.0 [56]. PopART [57] and the Median-Joining network (MJ) algorithm were used to construct the association map of the chloroplast DNA haplotypes. MEGA-X [53] was used to construct the phylogenetic tree of 3 haplotypes based on the maximum likelihood method (ML). Meanwhile, phylogenetic trees were constructed for the sequences of 80 accessions based on the Maximum Parsimony method (MP) with PAUP 4.0 [58]. The consistency index was 1.000000, the retention index was 1.000000, and the composite index was 1.000000 for all sites and parsimony-informative sites (in parentheses).

Phenotypic Traits Variation
ANOVA analysis on phenotypic variation showed that in addition to the number of flowers per plant, the other 27 quantitative traits showed very significant differences in 180 accessions (p < 0.01) ( Table 4). Since all investigated plants were adult with stable characteristics growing in a consistent cultivation condition, it could be considered that the phenotypic variation was mainly due to genotypic differences. Descriptive statistical analysis was conducted at the whole population level, and it was found that the coefficient of variation of flower traits, stem and leaf traits and fruit traits ranged from 12.33% to 108.73%, with a very high degree of variation. Among the flower traits, the coefficient of variation ranged from 12.33-108.73%. The variation of petal number (108.73%) was very obvious, but the other flower traits varied weakly in various accessions, among which the variation coefficients of flower diameter (12.33%) and petal length (12.33%) were the smallest. As for the branch and leaf traits, the variation coefficients of tillers' number (90.53%), fruits (66.86%) and flowers (57.43%) per plant were relatively higher. The variation coefficient of the fruit traits ranged from 17.84% to 99.93% (average 48.50%). Among them, the variation coefficient of seed weight per fruit was the largest (99.93%), followed by seed number per fruit (94.30%). So, significant phenotypic differences of phenotypic variation in the FTP accessions provide valuable germplasm resources for selective breeding of tree peonies.

Distribution of Phenotypic Traits Variation in the Population
PCA on 28 quantitative traits showed that the cumulative contribution rate of 8 principal components reached 70.516% (Table 5). The explanatory contribution rate of the first principal component was 17.722%. The traits that determined the first principal component were, in order, seed weight per fruit, individual fruit weight, number of seeds per fruit, fruit width, fruit height, fruit length and fruit weight per plant, all of which were fruit traits. Therefore, the first principal component represented the variation of fruit traits, which can be defined as the fruit character factor. The traits that determined the second principal component were fruit number, flower number, eastwest crown breadth, north-south crown breadth and fruit weight per plant. The third principal component mainly consisted of petal length, petal width, flower diameter, flare length and flare width, all of which were floral traits, which can be defined as the flower character factor. The fourth principal component consisted of compound leaf length, petiole length, compound leaf width, terminal leaflet length and terminal leaflet width, all of which were stem and leaf traits. Length and width of the terminal leaflet majorly affected the character of the fifth principal component. The fourth and fifth principal components can be defined as the branch and leaf character factor. The PCA of quantitative traits explained the variation distribution of phenotypic traits in the whole population. In addition, among the principal coordinate distributions of all quantitative traits, we found that fruit traits were distributed along the horizontal axis, stem and leaf traits along the vertical axis, and flower traits were distributed near the origin (Figure 3a). The 180 accessions were evenly distributed on the coordinate axes composed of the first principal component and the second principal component (Figure 3b), indicating that the variation of these phenotypic traits was evenly distributed in the population and no obvious hierarchical structure was formed.

Correlation Analysis of Quantitative Traits
Correlation analysis on 28 quantitative traits showed that there were different degrees of correlation among all traits (Table 6). Among the 378 pairs of combinations, 205 pairs of all traits showed significant correlation (p < 0.05), among which 191 pairs showed very significant correlation (p < 0.01). In tree peony breeding for ornamental uses, we mainly focused on the improvement of the flower diameter and the petal number. There were 17 traits that were very significantly positively correlated with flower diameter, including all floral traits except petal numbers and all fruit traits except effective carpel number. Petal number was negatively correlated with flower diameter, and only positively correlated with effective carpel number. This indicated that the traits of flower and fruit were closely related and the two ornamental traits we want to improve are opposite. For the breeding of oil tree peony, we need accessions that have more fruit number per plant and more seed weight per fruit. In Table 6, we can find that the trait that was very significantly positively correlated with these two traits was fruit weight per plant. For conventional breeders, it is not necessary to peel off the pericarp to weigh the seeds, but simply to weigh the fruit weight per plant to make a preliminary selection of high seed yield accessions.  Note: ns = not significant, * p < 0.05, ** p < 0.01. The traits corresponding to the serial number are shown in Table 1.

EST-SSR Polymorphism
Amplification of 282 accessions using 34 SSR pairs generated 185 alleles and the number of alleles detected per locus was in a range of two to thirteen, with an average of 5.441 alleles ( Table 2). Observed heterozygosity (Ho) and expected heterozygosity (He) ranged from 0.004 to 0.993 and 0.004 to 0.815, with means of 0.537 and 0.489, respectively. Among all 34 SSR loci, PS004, PS026, PS047, PS068, PS073, PS074, PS095, PS119, PS139, PS157, PS166, PS187, PS221, PS260, PS265, PS296, PS309, PS311, PS337, PS345, PS356 and PS367 had lower He than Ho. Overall, the mean of Ho was higher than that of He, indicating that there was a heterozygote excess in these accessions. The effective number of alleles (Ne) was 2.291 ± 0.166. Shannon diversity index (I) was 0.908 ± 0.072, and the range of PIC at 34 SSR pairs of 282 accessions was 0.004~0.792. According to the Shannon's Information Index, the genetic diversity of PS095 (I = 1.857) was the highest, followed by primers PS074 (I = 1.457), PS166 (I = 1.455) and PS356 (I = 1.476). Primer PS323 has the lowest genetic diversity. The mean value of SSR marker polymorphism information content (PIC) in this study was 0.611. In addition, 24 SSR sites deviated significantly from the HWE (p < 0.001), and 2 SSR sites deviated significantly from the HWE (p < 0.01).

Genetic Relationships among 282 Accessions
It was found that 282 accessions were firstly divided into 3 major branches (Figure 4a), comprising 175, 87 and 20 accessions. Branch Ⅰ included 61 white, 25 pink, 21 red, 59 purplish red and 9 purplish black accessions, branch Ⅱ included 23 white, 8 pink, 9 red, 44 purplish red and 3 purplish black accessions, and branch Ⅲ included 5 white, 2 red and 13 purplish red accessions. The accessions of each color were relatively evenly distributed on each cluster, which showed the same results as the PCoA analysis ( Figure 4b) and showed no regular distribution. Some accessions with similar traits have a relatively long genetic distance, while some with great phenotypic differences have relatively short genetic distances. For example, accessions '73' and '74', both of which had single white flowers, were very distantly related, while accessions '40' and '75', one purplish red double and one white single, were particularly closely related. This also indicated that there was no direct relationship between the relative distance and the color of the accessions.

Genetic Diversity Based on the cpDNA
The accD-psaI, psbE-petL and petB-petD fragments of three chloroplast non-coding regions were used to detect 80 accessions. After combination, the total length of the three fragments was 1580 bp. A total of 4 mutation sites were detected in the combined fragment, all of which were single base mutations and corresponding to three haplotypes. The petB-petD fragment had 2 variation loci, the other two fragments were each one, respectively (Table 7). According to the calculation results of the three combined fragments by DnaSP 5.0, the total Hd was 0.164 and Pi was 0.28 × 10 −3 . Table 7. The variations in accD-psaI, psbE-petL and petB-petD regions among 80 accessions.

Domain
Variation Site Hap 1 Hap 2 Hap 3 accD-psaI Note: A, T, C and G are bases, respectively.
In the haplotype spectrum diagram (Figure 5a), Haplotype 1 was located in the center of the network, which was an obvious widespread haplotype. It appeared in 73 accessions, which may be the original haplotype. Haplotype 2 and 3 contained 6 and 1 accessions, which were '98', '120', '140', '146', '150', '187' and '128', respectively. According to Figure 5a, Haplotype 1 was closer to Haplotype 3 than Haplotype 2, showing the same result with phylogenetic relationship of 3 haplotypes (Figure 5b). Three haplotypes were clustered into two branches in the tree, and Haplotype 2 formed a single branch. Meanwhile, phylogenetic trees were constructed for the sequences of 80 accessions. The accessions were divided into two distinct branches ( Figure 6). The first branch consisted of 6 germplasms of Haplotype 2, including '98', '120', '140', '146', '150' and '187', and the second branch consisted of 74 accessions of Haplotype 1 and Haplotype 3, which also indicated that the phylogenetic relationship of Haplotype 1 and 3 were closer than that of Haplotype 2.

Diversity of Phenotypic Traits
Because the phenotype is the result of the interaction of genotype and external environment and the diversity of phenotypic traits is the manifestation of genotype difference in morphology, the phenotypic variation in the population mainly depends on the genotype difference when all accessions used in this study were grown in the same condition. In FTP, Pang et al. [59] reported phenotypic variation coefficients ranging from 10% to 30% based on 32 traits of 150 accessions, but Wu et al. [60], on the diversity of 29 quantitative traits of 462 accessions, found that the coefficient of genetic variation was 9.52-112.1%. By the genetic diversity of 28 quantitative traits of 180 accessions, this study revealed a genetic variation coefficient of 12.33-108.73%, which was similar to Wu et al. [60]. As the plants used both by Wu et al. and by this study are widely selected and collected from the main cultivation area of the FTP, Gansu province ( Figure 2), we felt that the accession population of this study should be large enough to represent the overall situation of FTP as a germplasm resource. Moreover, each accession collected in this study has been clonally propagated and most of them were named as the cultivars (Supplementary Table S1) and can be released into the industry as crop germplasm resources.
Studies on phenotypic variation can not only improve our understanding of the biological basis of plants, but also have very important breeding value [61]. In FTP used in this study, the fact that the number of petals varied from 9 petals to 238.33 petals (Table 4) indicates that these accessions have diverse flower types, from single through to semi-double to double [1]. Such various flower types are really valuable for breeding new cultivars for ornamental plants like tree peonies. Moreover, the tree peonies have been recently developed rapidly as an emerging oil crop and the selection of high-yield accessions has become very crucial for increasing cultivation. Among the accessions used in this study, some of them with single flower, in most cases, are fertile for the formation of fruits (seeds) and can be screened as the high-yield oil germplasms, in which the maximum fruit weight per plant was up to 2661.64 g (Table 4).

Genetic Diversity Based on SSR Markers
The primary task of evaluating all kinds of samples in germplasm resources is to determine their genetic constitution in the repository. We can also compare newly collected cultivated or wild plants from various locations with existing germplasm accessions to facilitate identification. With the advent of molecular technologies and methods, it has become more efficient and scientific to use genetic information at the genome level to analyze accession differences in genetic resources and compare these samples with germplasm collections elsewhere [62].
This study detected a total of 185 alleles in 34 pairs of EST-SSR primers in 282 FTP accessions, and the average number of different alleles (Na) was 5.441 (Table 2), which was smaller than in 20 populations (Na = 9.15) of 335 P. rockii wild individuals [29], but larger than in 40 pairs of primers in 462 P. rockii seedlings (Na = 4.5) [38]. Such differences might be caused by the different primers used in the different studies, but it was very clear that the wild germplasm resources of P. rockii do have more abundant allele variation than that of cultivated FTP. The average inbreeding coefficients (FIS) was negative, and the expected heterozygosity (He) of 22 pairs of primers was lower than the observed heterozygosity (Ho) ( Table 2), indicating that there was a surplus of heterozygotes in this FTP population [29]. Yuan et al. [63] believed that the higher observed heterozygote meant higher diversity of the population, as their results showed that the wild P. rockii (Ho = 0.475) has significant diversity compared to CTP (Ho = 0.61) and FTP (Ho = 0.58). In FTP used in this study, the observed heterozygote was 0.537, which was larger in wild P. rockii (Ho = 0.475) but smaller in CTP (Ho = 0.61), and closer to FTP of Yuan et al. (Ho = 0.58). This study by Neighbor-joining tree based on microsatellite analysis identified the 282 FTP accessions into three distinct groups and there was no direct relationship between the genetic distance and the flower color. The grouping was consistent with the studies on 335 wild P. rockii individuals of 20 populations and 462 P. rockii seedlings, which were also divided into three groups [29,38]. Guo et al. [64] used SRAP to classify 16 accessions of P. suffruticosa, which also indicated that the genetic distance had no correlation with the color: white, green, yellow and black accessions formed a group, while blue, pink and/or multi-colored flowers formed another group.
To sum up, SSR data showed that the number of alleles for FTP used in this study was less than that of wild P. rockii, but more than that of cultivated P. rockii seedlings. The heterozygosity of FTP is less than that of CTP, larger than that of wild P. rockii, and close to that of cultivated P. rockii seedlings. This also showed that wild plants have abundant allele variation but low population heterozygosity, and cultivated germplasm have less allele variation but high population heterozygosity, which just indicates that gene exchange is more frequent and more heterozygosity is retained artificially under cultivation conditions [65]. Most woody perennial species are obligately outcrossing, resulting in high heterozygosity [66]. The longer the culture history of tree peony, the higher the heterozygous, and these individuals maintained high heterozygosity by asexual reproduction or clonal propagation under cultivation conditions. The FTP accessions were divided into three groups, which were consistent with the wild P. rockii and cultivated P. rockii seedlings. So, we speculate that the germplasm resources of FTP in this study, to a large extent, can represent the genetic background information of species tree peony P. rockii germplasm.

Genetic Diversity Based on the cpDNA
According to Yuan et al. [28], the mean nucleotide haplotype diversity (Hd) and the mean nucleotide diversity (Pi) of the chloroplast (cpDNA) in 335 P. rockii wild individuals at the population level was 0.0686 and 0.765 × 10 −4 respectively, while Hd and Pi at the species level was 0.887 and 0.185 × 10 −2 , which was significantly higher than population level. Based on the analysis of three cpDNA sequences, Xu et al. [37] measured the diversity of 214 individuals in 24 populations of P. rockii and obtained the same results as Hd and Pi at the species level, 0.874 and 0.294 × 10 −2 , which were significantly higher than Hd (0.1097) and Pi (0.383 × 10 −4 ) at the population level. This shows that wild P. rockii has a low genetic diversity at the population level and a higher genetic diversity at the species level. This situation was consistent with other tree peony species, P. delavayi [33], P. ludlowii [33], P. qiui [37] and P. jishanensis [37]. However, Hd and Pi in FTP accessions in this study were 0.164 and 0.28 × 10 −3 , which were higher than the population level and lower than the species level in populations of wild P. rockii. Considering these accessions as a population to compare with wild P. rockii, we found that the genetic diversity of the FTP under cultivation conditions is really different from its ancestor species in the habitats in forming genetic diversity or evolution. This is probably because wild plants of P. rockii reproduce by seeds in the same population, while gene exchange does not occur due to geographical isolation between populations. The chloroplast diversity of all individuals in the same population was lower, while it was higher at the species level because each population has distinct haplotypes. In this study, the accessions of FTP were selected and collected from various places of Gansu province without geographic isolation, and gene exchange may have occurred during the formation of germplasms through repeated hybridization, which finally resulted in a higher chloroplast diversity at the population level after long cultivation.
Meanwhile, based on the results of phylogenetic tree of chloroplast gene fragments ( Figure 6) and haplotype spectrum (Figure 5b) analysis, the existing spatial genetic structure of the FTP can be divided into two branches, which was consistent with Yuan et al. [27] using chloroplast fragment to divide the genetic structure of wild P. rockii into two branches.

Conversation and Utilization of FTP Germplasm Accessions
Research on genetic diversity and structure is the theoretical basis for the protection and utilization of species. One of the goals of species protection is to maximize the conservation of genetic diversity of the species [67,68]. The accessions we investigated have abundant phenotypic variations and represent FTP germplasm, from which we have selected to name some cultivars for ornamental oil uses. The next step is to establish multi-site repositories and a multi-year, comprehensive phenotypic database of FTP germplasm resources. Germplasm accessions are crucial resources for studying crucial characteristics like plant resistance and flowering time through multi-year and multi-site observations. Similarly, using these germplasms to study phenotypic variation may be significant to researching phenotypic limberness as they react to changing climatic and different living environments.
With the results from SSR and cpDNA analysis, the FTP accessions can represent the genetic background information of P. rockii germplasm resources and have relatively high genetic diversity. In order to protect this diversity, a germplasm resource bank should be established to conserve as many accessions as possible. Analysis of chloroplast sequence variation showed that three haplotypes were found in 80 accessions. One haplotype consists of 73 accessions, but the other two haplotypes only consist of six and one, respectively. This suggests that the presence of seven accessions is closely related to the genetic diversity of P. rockii germplasms, and their disappearance will bring about an irreversible loss of diversity. Therefore, special protection should be given to these unique accessions.
Due to the rapid change of global climate caused by human activities, the selection pressure of woody ornamental plant commercialization is increasing, and the diversity of existing germplasm resources has more important breeding value [69]. As a kind of high-quality germplasm resource with strong resistance, high ornamental value and great potential for oil use, P. rockii accessions play an important role in the improvement of tree peony. At the level of DNA, this cultivated P. rockii population we have established can represent the basic level of P. rockii germplasm resources. In terms of phenotypic diversity, the variation of each trait is also extremely rich. Therefore, the value of these accessions is especially obvious for breed improvement and development of tree peonies. Meanwhile, the germplasm resources contain various plant genetic resources, which is of great significance to the basic research of plant biology. In order to explore the potential excellent characters within these accessions, we are required to continue to explore the phenotypes, genotypes and preservation of these resources. In this study, we planted these accessions from the traditional cultivation area of Gansu province to Beijing. It is necessary to continuously observe and measure these characters to explore the adaptability of the FTP to the different environment. Ultimately, germplasm collections should not simply be viewed as tools for conservation, but as resources with unique benefits that hold immense value for sustaining the future of woody perennials and their wild relatives. The comprehensive evaluation and utilization of these cultivated accessions is to protect the germplasm of P. rockii tree peony and its derivatives.

Conclusions
This study is the first to analyze the genetic diversity of FTP accessions newly selected by us using phenotypic traits, EST-SSR markers and chloroplast DNA sequences. Phenotypic data indicated that these accessions have abundant phenotypic variations and are representative of woody peony P. rockii germplasm. SSR data showed that the FTP accessions used in this study had relatively high genetic diversity and were composed of three distinct groups, like their wild relatives. The cpDNA data indicated that the genetic diversity of chloroplast genome in the FTP accessions had a higher genetic diversity than the population level and a lower genetic diversity than the species level of wild P. rockii, and they can be divided into two branches based on this. The FTP accessions fully reflected the genetic diversity and current situation of species tree peony P. rockii germplasm to a large extent, so they can be regarded as the core germplasm resources and their protection and utilization will be of great significance both for genetic improvement and biological studies of tree peonies and for sustaining development of the peony industry.