The assessment of epigenetic diversity, differentiation, and structure in the ‘Fuji’ mutation line implicates roles of epigenetic modification in the occurrence of different mutant groups as well as spontaneous mutants

The ‘Fuji’ line includes many varieties with a similar genetic background and consistent inducement factors with epigenetic occurrence, thus it may be considered an ideal candidate for epigenetic research. In this study, 91 bud mutations of ‘Fuji’ apple were used as the test materials. Using the genetic variation within ‘Fuji’ as the control, the characteristics of epigenetic variation at different levels in both varieties and mutant groups were examined. The results showed that: (1) the global genomic DNA methylation level of the 91 bud mutants of ‘Fuji’ ranged from 29.120%-45.084%, with an average of 35.910%. Internal cytosine methylation was the main DNA methylation pattern. Regarding the variation of methylation patterns of ‘Fuji’ mutants, the vast majority of loci maintained the original methylation pattern existed in ‘Fuji’. CHG methylation variation was the main type of variation; (2) the variation in methylation patterns between the mutant groups was greater than that of methylation levels. Among these patterns, the variation in CHG methylation patterns (including CHG hypermethylation and CHG demethylation) was expected to be dominant. The observed variation in methylation levels was more important in the Color mutant group; however, the variation in methylation patterns was more obvious in both the early maturation and Spur mutant groups. Moreover, the range of variation in the Early-maturation group was much wider than that in the Spur mutant group; (3) epigenetic diversity and genetic diversity were both low between the mutant groups. In the ‘Fuji’ mutant groups, there was few correlation between genetic and epigenetic variation, and epigenetic differentiation resulted in more loci with moderate or greater differentiation; (4) the purifying selection seemed to play a major role in the differentiation of different groups of ‘Fuji’ mutants (65.618%), but epigenetic diversity selection still occurred at nearly 35% of loci. Sixteen epigenetic outlier loci were detected.

Introduction Malus domestica Borkh. cv. 'Fuji' is one of the most economically important apple varieties in China and the world. 'Fuji' apple is prone to mutation, through which abundant bud mutants have been derived, the majority of which have been adapted as superior varieties in production [1]. Studies have shown that epigenetics is an important cause of bud mutations in fruit crops, but it has not yet to be systematically studied in a large sample size. The study of the characteristics of epigenetic effects in bud mutants is the premise for further studying the epigenetic mechanism of mutations and carrying out epigenetic breeding.
Epigenetics studies the heritable changes in gene function that eventually lead to phenotypic variation with no changes in the underlying DNA sequence [2]. Epigenetics is involved in the regulation of gene expression, and changes are dynamic with respect to the endogenous and/or external environmental stimuli, thus affecting the phenotypic plasticity and environmental adaptability of organisms [3]. Epigenetic variation enhances biodiversity and complexity, especially in asexual organisms without gene recombination, which is helpful to the promotion of functional phenotypic diversity and has a greater impact on variation and evolution and more survival significance [4].
DNA methylation, representing the most important form of epigenetic modification, is ubiquitous in higher plants [5]. Different level of DNA methylation exists in normal organisms (5%-40%, [6]), while changes in DNA methylation can be caused by external environmental factors or various biotic and abiotic stresses [7]. Changes in DNA methylation level and methylation pattern may cause significant phenotypic variation in plant genome [8]. It can be expressed by the dynamics of hypermethylation and demethylation [9]. DNA methylation occurs on cytosines in different sequence contexts, and CG and CHG are the two main types of genome methylation. CG and CHG methylation are regulated by distinct enzymes and pathways [10]. Both CG and CHG methylation levels are different in a plant species, and CG methylation was reported predominately in many plant species [11].
The current research on DNA methylation is mainly focused on local and global aspects. The former mainly targets specific traits, starting with the changes directly related gene methylation in promoter or gene core sequence, and has more advantages in revealing the mechanism of occurrence of specific phenotypic variation; The latter studies the changes in the overall methylation of the genome, which is more conducive to the comprehensive analysis of the role of epigenetics in the variation and diversity. Many techniques have been developed to analyze global DNA methylation and its alterations [12]. The method of methylation-sensitive amplified polymorphism (MSAP) is one of the most efficient, economical, and widely one used for the detection of DNA methylation events [13][14]. This approach is a modified version of the amplified fragment length polymorphism (AFLP) method based on the differential sensitivity of isoschizomeric restriction endonucleases to site-specific cytosine methylation [15]. MSAP analysis employs two types of cleaved enzymes EcoRI (rare cutter) and HpaII/ MspI (frequent cutters) which recognize similar tetranucleotide 5'-CCGC representing a different sensitivity to methylation at the inner or outer cytosine. If one or both cytosines are methylated at both DNA strands, HpaII is inactive. If one or both cytosines are methylated in only one strand they are cleaved by HpaII. In contrast, MspI reacts when only the internal cytosine is hemi-or fully-(double strand) methylated [16]. So, on the basis of variant band patterns resulting from differential digestion of the genome by HpaII/MspI isozymes, variation in the DNA methylation level and pattern in the whole genome is detected [15,17]. MSAP has been successfully applied for the analysis of the variation in methylation levels and patterns in a variety of plant species [13].
In general, bud mutations are an important source of new varieties of fruit crops. Characteristics such as a perennial nature, long juvenile phase, heterozygosity, and sexual incompatibilities in fruit crops hamper their improvement through conventional breeding [18]. Compared with conventional hybridization, the selection of bud mutants gives the advantages of shortening the breeding cycle and reducing the workload and costs. Such an approach can be used to obtain excellent varieties by modifying individual traits without changing the desirable qualities of the parent plant [19]. A variety of perennial fruit trees of economic importance originated from bud changes [20][21]. However, Spontaneous bud mutation occurs at a very low frequency, moreover, many of these mutations may be deleterious, making the organism less adapted to its environment, and in some cases may even be lethal [18]. Therefore, mutants that survive in adverse environments and even present excellent phenotypes are considered to have good characteristics for adaption. Epigenetic regulation mediated by DNA methylation is considered as one of the important mechanisms in plant adaptive procedure [14,22]. Studying the molecular mechanism of bud mutation at the DNA level is thus of great significance.
'Fuji' is the most representative apple cultivar in which plenty of new bud mutations arise [23]. Multiple bud mutants generated from the standard cultivar 'Fuji'. Those represent highly similar genetic backgrounds and have abundant types of variation. Additionally, 'Fuji' bud mutants were frequently reported to be available in orchards at unusual locations, such as under dense high-voltage lines, at high altitudes, subjected to an abnormal climate, severe drought, waterlogging, frost, or sudden diseases or insect pest infestations [24]. This means that the conditions for induction of bud mutation are consistent with the factors inducing epigenetic effects. Therefore, in a summary of the above mentioned, 'Fuji' mutation line could be regarded as ideal sets for research on epigenetic regulation. However, to our knowledge, global DNA methylation in this line has not been reported yet.
A general understanding of the mechanisms of genome-wide DNA methylation in 'Fuji' mutants is a prerequisite for their utilization in epigenetic mechanism studies or breeding. Therefore, in the present study, nearly one hundred 'Fuji' bud mutants were used as study materials for the first time. Then, through genetic variation as control, we focused on the analysis of variations in DNA methylation levels and patterns within groups and between groups as well as epigenetic diversity.
We were interested in the following questions: (1) is DNA methylation involved in the occurrence of bud mutations in the 'Fuji' line? (2) what are the characteristics of the changes in both DNA methylation levels and patterns? (3) what are the effects of DNA methylation on not only the clustering status of the 'Fuji' line at the variety level but also genetic structure and differentiation at the group level? (4) how is natural selection affecting the occurrence of different mutant groups in the Fuji line?

Plant material and MSAP, AFLP analyses
This study was conducted on 'Fuji' and it's 91 elite varieties arising from bud mutations. Three major types: Color mutant group, Early-maturation mutant group, and Spur mutant group could be classified (Table 1). We sampled plant materials with the same age on the same dates (12 June 2016) and at their identical phenological stage (bearing fully expanded leaves), trying to avoid that possible developmental variation in DNA methylation would confound variety or group differences in methylation patterns. Fully expanded fresh leaves from pooled individuals (ca. 5 plants being clonally propagated from a single mother variety) of each variety were collected. Young leaves immediately frozen in liquid nitrogen and then stored at -80˚C prior to DNA isolation. Total genomic DNA from leaf tissue was extracted by using DNeasy Plant Mini Kit (Qiagen). MSAP molecular marker was used for the analysis of epigenetic variation in plants examined. For comparison, genetic variation was also conducted by AFLP molecular marker. Both methods performed on the same set of plants, and EcoRI primers were labeled with or green (JOE), blue (FAM), or yellow (NED) fluorescent dyes.
MSAP and AFLP were performed as described by Xiong et al. [17] and Vos et al. [25], respectively. MSAP was essentially the same as the AFLP protocol. The difference between MSAP and AFLP procedure was replacing the MseI enzyme with the enzyme either HpaII or MspI in MSAP. Thus, differences in the PCR products detected with EcoRI/HpaII and EcoRI/ MspI would reflect different methylation states. The primer sequences of MSAP and AFLP were listed in the S1 Table. Based on previous pilot tests, we selected 23 optimal primer combinations for MSAP analysis (S2 Table) and 19 optimal primer combinations for the AFLP analysis (S3 Table). HM-, E-, and M-corresponded to the sequence of H/M00, E00, and M00, respectively (S1 Table). The repeatability of banding patterns assessed by conducting two sets of independent MSAP and AFLP analyses and only the consistent bands were included.
Using an automated sequencer (ABI PRISM 1 3730 Genetic Analyzer, Applied Biosystems) to separated and detected MSAP and AFLP PCR products. GeneScan Rox-500 labeled with a red (ROX) dye was an internal size marker. A peak size between 80 and 500 bp was selected to study the polymorphic DNA fragments. MSAP products were scored as present '1' or absent '0' on the chromatogram to create a binary matrix.

Statistical analysis
The resulting data of MSAP and AFLP were processed using Excel 2016. According to the scoring approach [26], MSAP raw data were transformed into a binary data matrix before running statistical analyses and computation. Because the enzymes used for MSAP analysis recognize the same restriction site (5'-CCGG) but have different sensitivities to methylation modifications, the final products of MSAP present four types of DNA methylation at the 5'-CCGG sites, namely, no methylation (H1M1, condition I), CHG methylation (H1M0, condition II), CG methylation (H0M1, condition III), and CG/CHG methylation (H0M0, condition IV). DNA methylation level% = (condition II+condition III)/(condition I+condition II+condition III+condition IV).
Methylation-Sensitive Polymorphism (MSP) matrix data, which were converted using R program msap [27], were used for all analyses of DNA methylation. The change of patterns of cytosine methylation at CCGG sites in the 'Fuji' mutation line was listed in Table 2. Accordingly, contrasting to the standard cultivar 'Fuji', the changes of DNA methylation pattern in its series mutations could be summarized into four categories: (a) CG hyper: H1M1 to H0M1, H1M0 to H0M0, H1M0 to H0M1, H1M1 to H0M0; (b) CHG hyper: H1M1 to H1M0, H0M1 to H0M0, H1M1 to H0M0, H0M1 to H1M0; (c) CG hypo: H0M1 to H1M1, H0M0 to H1M0, H0M1 to H1M0; and (d) CHG hypo: H1M0 to H1M1, H0M0 to H0M1, H1M0 to H0M1, H0M0 to H1M1.

DNA methylation analysis and variation coefficient calculation of different mutant groups in 'Fuji'
The epigenetic relationship of different mutant groups in 'Fuji' was analyzed with 12 indexes related to DNA methylation levels and patterns. Indicators include CHG hyper-methylation frequency (V5), CG hyper-methylation frequency (V6), CG hypo-methylation frequency (V7), CHG hypo-methylation frequency (V8), the total hyper-methylation frequency (V9), the total hypo-methylation frequency (V10), total genome methylation frequency (V11), the frequency of condition I (Non-methylation frequency, V12), the frequency of condition II methylation (CHG methylation frequency, V13), the frequency of condition III methylation (CG methylation ratio, V14), the frequency of condition IV methylation (V15), and the total amplified loci of varieties (V16). The Coefficient of Variation (CV) of each of the above indexes was calculated in the three mutation groups in order to investigate whether varieties in different mutant groups behave in different ways. CV = (standard deviation SD/Mean)×100%. Correlation analysis of 12 major epigenetic parameters (V5-V15) was performed by software IBM SPSS Statistics 22 [28].

Epigenetic and genetic similarity, clustering and principal component analysis
The similarity coefficient between varieties was calculated using the SM coefficient through the SimQual procedure of NTSYSpc2.11 software package [29] and unweighted pair group method average (UPGMA) method was used for cluster analysis. Circle diagrams were drawn using TBtools [30]. MSAP-PCA and AFLP-PCA analysis were carried out by R program package msap 27 and Adegenet 2.1.1 [31], respectively.

Epigenetic and genetic diversity and molecular variation analysis
Various diversity measurement parameters including polymorphism site proportion (PPL), effective allele variance (Ne), Shannon diversity index (I), expected heterozygosity (He), and hierarchical analysis of molecular variance (AMOVA) analysis between and within mutant groups were estimated in GenAlex 6.51 [32], using 999 random permutations.  Epigenetic and genetic Structural analysis STRUCTURE 2.3.4 [33] was utilized to analyze genetic and epigenetic structure. Admixture and correlated allele frequencies model were chosen. Ten independent runs were made with values of K set from 2 to 4, with three iterations for each value of K. The length of the burn-in period was set at 10,000, and the number of Markov chain Monte Carlo (MCMC) repeats after burn-in was set at 100,000. The result from STRUCTURE output file was performed online by STRUCTRE Harvester 0.6.8 [34]. The results were averaged for a particular K using CLUMPP 1.1.2 [35] and visualized by DISTRUCT [36].  [38] scale of evidence, a log 10 PO>2.0 is interpreted as 'strong evidence' of selection. For our analysis, the estimation of model parameters was set as 10 pilot runs of 5,000 interaction each, followed by 100,000 interactions [37]. Outliers were calculated using a burn-in of 50,000 interactions, a thinning interval of 20, and a sample size of 5000. FDR = 0.05 was used.

Correlation analysis between epigenetic diversity and genetic diversity
Using NTSYS2.0 software, the correlations between epigenetic distance and genetic distance of three 'Fuji' mutant groups as well as varieties calculated by MSAP and AFLP were evaluated using the Mantel test implemented through NTSYSpc 2.11 software package [28].
All of the statistical significance in this study was determined by IBM SPSS Statistics 22 [32]. A value of P < 0.05 was considered significant.

MSAP and AFLP amplification
A total of 2954 CCGG loci were detected via the genome-wide methylation analysis of 92 'Fuji' varieties with 23 pairs of MSAP primer combinations, and 129 CCGG loci were amplified on average with each pair of primers. Among these loci, 2752 CCGG loci showed polymorphism, for a polymorphic ratio of 93.162%. Different MSAP amplification patterns were obtained from different varieties ( Table 1). The number of amplified methylated loci ranged from 1793 to 2136, and that of polymorphic loci ranged from 1613 to 1956, with an average polymorphic ratio of 63.566% (58.612% to 71.076%). Among the 2954 CCGG loci, 1748 were methylationsensitive loci (MSL), accounting for 59.174% of the total amplified loci. A total of 1627 polymorphic MSL loci were obtained, accounting for 93.078% of the total amplified loci (S2 Table).
Nineteen pairs of AFLP primer combinations amplified 1745 total loci in the same test varieties used for MSAP, among which 1620 were polymorphic loci, accounting for 92.837% of the total amplified loci. Different primer combinations produced different amplification results. The total number of amplified loci obtained with a single pair of primers ranged from 36 to 166, the number of polymorphic loci ranged from 28 to 150, and the polymorphic ratio ranged from 57.143% to 100%. On average, the total numbers of loci and polymorphic loci generated by amplification with each pair of primers were 92 and 86, respectively. The four primer combinations (E-AGG+M-CAG, E-AGG+M-CTA, E-ACG+M-CAC, E-ACG+M-CTG) produced 100% polymorphic loci. Primer amplification details are shown in S3 Table.

Analysis of DNA methylation levels and variation patterns at the variety level
As shown in Table 1, the global genomic DNA methylation level of the 91 varieties was 29.120%-45.084%, with an average of 35.929%. Among these modifications, the internal cytosine methylation level was 22.311% (16.528%-29.448%), and the external cytosine methylation level was 13.599% (10.242%-22.870%). The former was significantly higher than the latter (P<0.01), implying that internal cytosine methylation was the main DMA methylation way in the 'Fuji' line.
According to the cutting profiles of the HpaII and MspI methylation-sensitive endonucleases in the original 'Fuji' variety, the banding patterns could be divided into four types: A, B, C, and D ( Table 2). In comparison with the original 'Fuji', there were many types of possible locus variation in the mutant varieties, so the variation pattern of methylation loci could be subdivided into several subcategories: A1, A2, A3, A4, B1, B2, B3, B4, C1, C2, C3, C4, D1, D2, and D3. As shown in S1 Fig., in DNA methylation variation patterns A and C, loci identical to the original methylation pattern presented the highest proportion, indicating that during the occurrence of bud mutations in 'Fuji' line, the majority of loci maintained the original methylation pattern, while only a few loci exhibited methylation variation (S5 Table). All types of methylation variation patterns (15 subclasses) were detected in the test varieties.
The above 15 subclasses of methylation variation patterns were further categorized into four types: CG hypermethylation (CG-hyper), CHG hypermethylation (CHG-hyper), CG demethylation (CG-hypo) and CHG demethylation (CHG-hypo). As shown in S4 Table, CHG-hypo, CHG-hyper, CG-hyper, and CG-hypo displayed differences among different tested varieties, but the general trend basically showed the following correlations: CHG-hypo>CHG-hyper>CG-hyper>CG-hypo (S2 Fig). (CHG-hypo+CHG-hyper) was significantly higher than (CG-hypo+CG-hyper) (P<0.01). The relative trend between the total demethylation frequency and the hypermethylation frequency in different varieties also exhibited diversity, including the following findings: (1) the demethylation frequency was higher than the frequency of methylation, (2) the hypermethylation frequency was higher than that of demethylation, and (3) hypermethylation and demethylation frequencies were approximately the same. These results indicated that methylation pattern variations were not fixed during the occurrence of the 'Fuji' mutation.

Analysis of the DNA methylation level and variation pattern at the mutant group level
The epigenetic relationships of the genomes of different mutant groups in the 'Fuji' line were analyzed using 12 parameters related to DNA methylation levels and patterns (V5-V16). Table 3 showed that the variation coefficients (CVs) of V5, V6, V7, V8, V9, and V10 were generally greater than those of V11, V12, V13, V14, V15, and V16. Among these 12 parameters, V5-V10 and V11-V16 reflected variations in the DNA methylation pattern and DNA methylation level, respectively. Therefore, it could be deduced that the variation of methylation patterns among varieties was greater than that of methylation levels. V11-V16, V11, V12, V15, and V16 exhibited the same degree of variation among groups, while V13 and V14 showed relatively higher CV than these groups, indicating that variation in CHG and CG methylation levels is abundant among varieties of different mutant groups in the 'Fuji' line.
The CVs of each parameter were compared between the three mutant groups and analyzed (Table 3). It was shown that the emphases of epigenetic variation were different in different groups; for example, in the Early-maturation group, several parameters, including V5, V6, V7, V9, and V10, presented higher values than in the other groups; in the Spur group, V8 was much more prominent; in the Color group, V13 showed a remarkably high value. Based on the above results, it seemed that the variation in the methylation level was more important in the Color group. The variation in the methylation pattern was obvious in the Early-maturation and Spur groups; moreover, the extent of the variation in the Early-maturation group was much wider.
The correlations of 12 DNA methylation parameters were calculated and plotted using the R packages psych and corrplot. As presented in Fig 1, a total of 43 pairs (P<0.01) were significantly correlated with each other, among which V6 and V9 exhibited the highest significant positive correlation (0.93), followed by V5 and V9 (0.91). The lowest correlation was found for V10 and V13 (0.25), followed by V9 and V12 (0.31). The significant negative correlation was highest between V5 and V11 (-0.24), followed by V10 and V14 (-0.27). The lowest negative correlation was found for V15 and V16 (-1.00), followed by V10 and V15 (-0.74). The results showed that there was a direct positive correlation between the total hypermethylation frequency and the frequency of CG hypermethylation (that is, the higher the frequency of CG hypermethylation, the higher the frequency of total hypermethylation). The total demethylation frequency presented a small correlation with the proportion of CHG methylation but a significant negative correlation with the proportion of CG methylation (that is, the higher the proportion of CG methylation, the lower the total demethylation frequency). Furthermore, the above 12 DNA methylation parameters of the three mutant groups were compared and analyzed. The data distribution of all 12 variables conformed to a normal distribution and the assumptions of the statistical analysis. The variance homogeneity test showed that, with the exception of V9 (the total hypermethylation frequency), the other 11 parameters all met the requirements of homogeneity of variance. Therefore, the Tamhane's T2 and Games-Howell tests were chosen for V9, whereas the Duncan and LSD tests were chosen for the other parameters for multiple comparison analysis with SPSS software. Multiple comparative analysis among groups showed that among the three mutant groups, V5, V6, V7, V8, V9, V10, V15, and V16 presented significant differences (P<0.01), but no significant difference was found between V11, V12, V13, and V14. The average multiple comparison results revealed that the differences between the Spur and Early-maturation groups were not significant, showing similar epigenetic characteristics, but different degrees of significant differences existed between the Spur and Color mutant groups as well as between the Early-maturation and Color groups. For example, V5, V6, V7, V8, V9, V10, V15, and V16 were significantly different between the Spur and Color mutant groups; the difference in V6, V7, V8, V9, V10, V15, and V16 between the Early-maturation and Color mutant groups was very significant (Fig 2). Taken together, the above results showed that among the three 'Fuji' mutant groups, the Spur and Early-maturation groups showed similar epigenetic patterns; however, between the Color groups and either the Spur group or the Early-maturation group, epigenetic differences were apparent.

Epigenetic similarity calculation, clustering and principal component analysis
The results of the genetic similarity analysis based on the MSAP-MSL data sets (S6 Table) showed that the epigenetic similarity coefficient between varieties ranged from 0.515 to 0.807; the lowest similarity was found between 'Yishuizhongqiu' and 'Wengao No.3', and the highest . The results showed that there was apparently epigenetic variation between mutant varieties and their 'Fuji' mother, and different degrees of epigenetic differentiation also occurred between different mutant varieties.
The UPGMA results showed (Fig 3(A)) that nearly 84% (43/51) of the varieties in the Color mutant group clustered closely together to form an independent cluster (Cluster 1). The Earlymaturation mutant group, the Spur-maturation mutant group, 'Fuji', and the remaining eight

PLOS ONE
The epigenetic assessment of Fuji mutation line implicates roles of epigenetic modification varieties of the Color mutant group formed another group, designated Cluster 2. In comparison with the position of 'Fuji', other varieties showed obvious epigenetic variation, and the epigenetic variation of the Color mutant group was most prominent. We also found that within Cluster 2, the Early-maturation mutant group and Spur mutant group were mixed and distributed without distinct grouping, reflecting the close epigenetic relationship between them. In addition, although 'Chang Fujia', 'Zhaofuwang', 'Yishuihong', and 'Jihong' in the Color mutant group were concentrated in Cluster 2, they were relatively closer to Cluster 1, but 'Akifu No.5', 'Line I Fuji', 'Iwate line I', and 'Iwafu No.10' were located far from Cluster 1. These four varieties were all selected in Japan and belonged to the Color mutant group. Moreover, they were the earliest mutant cultivars obtained from 'Fuji', and few clonal descendants are selected from these cultivars at present. Therefore, their clustering results might be related to their original region of selection or breeding generation. In conclusion, the UPGMA results implied that the greatest epigenetic variation existed in the Color mutant group, which might present unique epigenetic and evolutionary mechanisms compared with the other two variant groups.
The results of principal component analysis further supported the UPGMA clustering results. Principal components PC1 and PC2 accounted for 13.6% and 4.5% of the total variation, respectively. As shown in Fig 4(A), the Spur mutant and Early-maturation groups were mixed, showing similar epigenetic consistency. 'Fuji' was included in this group, displaying a close epigenetic relationship with its members. The Color mutant group was basically independent and was located far from 'Fuji' and the other two mutant groups, indicating its unique epigenetic variation.  Table). The mean similarity coefficient between the standard 'Fuji' variety and its mutants was 0.80, ranging from 0.684 (between 'Fuji' and 'Su Fuji') to 0.877 (between 'Fuji' and 'Yanfu No.7'), indicating that the genetic variation of the mutants differed from that of 'Fuji".

Genetic similarity calculation, clustering and principal component analysis
The UPGMA results showed (Fig 3(B)) that 92 samples could be clustered into four groups (Clusters 1-4), and 'Aki-fu No.5' was clustered in the outermost area and separated into  AFLP (B). Blue, green, and red indicate the materials from the Spur group, the Early-maturation group, and the Color group, respectively. https://doi.org/10.1371/journal.pone.0235073.g004 different groups, indicating distant genetic similarity. In addition, 'Hongwangjiang', 'Gai Fuji', and 'Su Fuji' were grouped together to form Cluster 3. Cluster 1 and Cluster 2 were closely related to each other. 'Fuji' was clustered in Cluster 1, which was closely related to 'Shanfu No.6', 'Qiufuhong', 'Lele Fuji', and 'Qingnonghe No.1'. In Cluster 1 and Cluster 2, local clustering and accumulation phenomena of the same 'Fuji' mutation type were found, as observed in the regions between 'Yanfu No.7' and 'Duanzhi Fuji' (Spur mutant group), 'Aki Fuji' and 'Nagafu No.36' (Color mutant group), and 'Miyazaki' and 'Chengji No.1' (Spur mutant group) in the former group and the region between 'Qianxuan No.1' and 'Shoufu No.3' and that between 'Yanfu No.3' and 'Nagafu No.4' in the latter group. However, from the overall perspective of the cluster diagram, the three major 'Fuji' mutant groups were mixed and arranged, and there was no obvious boundary between different groups. Fig 4(B) showed the PCA results of the 92 test varieties. The three dimensions of PC1 and PC2, PC1 and PC3, and PC2 and PC3 all showed results that were consistent with UPGMA clustering; that is, the genetic similarity among the three major mutation types of 'Fuji' was very high, without obvious boundaries.

Epigenetic and genetic diversity and molecular variation analysis
As shown in Table 4, the mean Shannon index, the number of effective alleles, and the mean expected heterozygosity based on MSL data among the three mutant groups were 0.415, 1.424, and 0.266, respectively. In general, the three genetic diversity indexes showed a consistent variation trend. For mutants with a high Shannon index (I), the effective allelic variance (Ne) and mean expected heterozygosity (He) were also higher. The I and He indices presented significant differences among the three mutant groups (P<0.01). The analysis of molecular variance showed that most of the variation occurred within the mutant groups (85%), and only a very small proportion (15%) occurred among the mutant groups (P< 0.01).
The average Shannon index, number of effective alleles and average expected degree of heterozygosity in the three 'Fuji' mutant groups reflected by AFLP analysis were 0.290, 1.303 and 0.186, respectively. In general, there was no significant difference in genetic diversity indexes between the three groups basing on AFLP data. AMOVA showed that most of the genetic variation existed within the mutant groups (95%); a very small portion (5%) existed between the populations (P< 0.01).

Epigenetic and genetic structure analysis
Structure analysis was carried out based on MSL epigenetic data and AFLP-based genetic data, and 10,000 simulations were run with population numbers K = 2 to 4. The results showed that LnP (D) values were generated for K = 3 for both types of data, and no values were generated for K = 2 or K = 4. Therefore, K = 3 was selected as the optimal number of clusters. Fig 5(A) showed that in the 'Fuji' mutation line, the Spur mutant group and the Early-maturation mutant group were highly consistent, while the Color mutant group showed a large proportion of heteromorphic lineages, indicating the unique epigenetic structure of the Color mutant group. As shown in Fig 5(B), there was no clear division between any of the mutant groups, which reflected the similarity of their genetic compositions.

Outlier analysis
Bayesian population genomic analysis showed that out of 1748 methylation-sensitive MSAP loci, 601 exhibited positive alpha values, and 1147 exhibited negative alpha values. A positive alpha value indicates positive or diversified selection, while negative values indicate the negative or purifying selection, also known as neutral selection. Thus, purifying selection seemed to play a major role in the differentiation of the different groups of 'Fuji' mutants (65.618%). However, epigenetic diversity selection still occurred at nearly 35% of loci. Sixteen epigenetic outlier loci were detected, accounting for 0.915% of the total DNA methylation-sensitive loci (Fig 6(A)). The epigenetic genetic frequencies of each locus in each mutant group were calculated individually. If the frequency of a locus in a group was greater than 50%, it was set as a specific outlier locus of that population. Therefore, Fig 7 showed that there were 13, 13, and 3 methylated outlier loci in the Spur, Early-maturation, and Color mutant groups, respectively. According to the profile of the frequency of loci, the Spur and Early-maturation mutant groups presented the same pattern, which basically maintained the state of the existence of the majority of loci in 'Fuji'. However, the pattern of the Color mutant group was completely different from that of the other two mutant groups as well as 'Fuji'. In the Color mutant group, 13 locus deletions and 3 new locus mutations occurred, indicating different epigenetic patterns of the outlier loci among the different mutant groups and suggesting unique epigenetic characteristics of the Color mutant group. According to AFLP analysis, among the 1745 amplified loci, Bayesian population genomic analysis generated 568 loci with positive alpha values and 1177 loci with negative values. Purifying selection played a major role (67.450%) in the differentiation of the different 'Fuji' mutant groups. Genetic diversity selection occurred at only a few loci (32.550%). At the level of FDR = 0.05, the number of outlier loci detected was 0, indicating that there was no significant difference in genetic differentiation between the three 'Fuji' mutant groups (Fig 6(B)).

Mantel correlation test of epigenetic and genetic similarity between different mutant groups
The Mantel test was conducted for the genetic and epigenetic similarity coefficient matrix obtained based on MSAP and AFLP analysis. The results showed that the r values of the genetic and epigenetic correlation coefficients of the Spur mutant group, Early-maturation mutant group, and Color mutant groups were -0.044 (P = 0.368), 0.327 (P = 0.978), and 0.107 (P = 0.924), respectively. The overall correlation coefficient was 0.108 (P = 0.987). The above results showed that there was a slight correlation between genetic and epigenetic variation in the 'Fuji' mutation lines, suggesting that the epigenetic variation was independent of genetic variation.

Comparison of epigenetic and genetic diversity between different mutant groups
The average polymorphic loci ratio, the polymorphic loci ratio, as well as epigenetic/genetic diversity parameters of each 'Fuji' mutant group based on MSAP were significantly higher than those obtained based on AFLP analysis (P<0.01), indicating that the epigenetic diversity of each mutant group in the 'Fuji' mutation line was significantly higher than the genetic diversity. For AMOVA between as well as within the mutant group, a certain degree of differentiation under both markers was revealed, and among groups, the epigenetic differentiation ratio (15%) was greater than the genetic differentiation ratio (5%).

Comparison of epigenetic and genetic structure between different mutant groups
The results of the UPGMA clustering analysis, PCA, and structure analysis based on AFLP and MSAP data showed that the three 'Fuji' mutant groups were genetically similar and that, regardless of what algorithm the analysis was based on, they could not be clearly separated. In terms of epigenetics, there was a close relationship between the Spur and the Early-maturation mutant groups, which were always mixed, similar to their genetic relationship. However, for the varieties of the Color mutant group tended to cluster together, and this group was independent of the other two mutant groups, showing its unique epigenetic structure.

Comparison of outlier loci between different mutant groups
The numbers of loci with positive and negative alpha values were approximately the same when generated on the basis of AFLP and MSAP. For each marker type, the proportion of loci with a negative alpha value was greater than that with a positive alpha value. Regarding outlier loci, MSAP analysis detected 16 loci, whereas AFLP analysis detected 0, indicating that epigenetics played a large role in the differentiation of the 'Fuji' mutant groups; and from the viewpoint of the frequency of outlier loci between different groups, epigenetics in the Color mutant group contributed greatly to the differentiation of the 'Fuji' mutant groups.

Discussion
Extensive epigenetic research has been carried out in many plant species, whereas study on the fruit crops is still in progress. So far, the available reports mainly consist of several popular fruit crops such as apple [39][40][41][42][43][44][45][46][47][48], orange [49][50], banana [51], pear [52][53], strawberry [54], and grape [55][56]. Epigenetic studies in bud mutations are just in the initial stages since limited literature is available [40-41, 45-47, 52-53]. In 'Fuji' apple bud mutation occurrence, epigenetic mechanisms mediated by DNA methylation begun only last year [46][47][48], and all of them were done from the perspective of local DNA methylation. Few studies were performed on global DNA methylation in bud variation with a large sample size even on fruit crops. In the present study, we have incorporated a large-scale collection of 'Fuji' mutants from around the world as plant materials. Moreover, investigations of the genetic and epigenetic perspectives, at both the individual variety and mutant groups were conducted for the first time. Our studies revealed wide DNA methylation alterations between bud mutations in the 'Fuji' line. The results would be helpful for comprehensively understanding the epigenetic characteristics and differentiation of the 'Fuji' line. It also benefits for future further research on the epigenetic mechanism of bud mutations.
From an adaptive perspective, the modification of methylation status may allow trees to rapidly respond to abrupt changes in environmental conditions and contribute to their longterm responses to more general environmental scenarios [57]. Bud mutations often occur under stresses such as drought, pathogen infection, extreme weather and climate conditions, limited nutrient availability, human activities, and natural or artificial selection and are a result of adaptation to abnormal external environmental conditions [58,59]. The contribution of epigenetic modification to the ability of plants to adapt to various stresses has been well demonstrated [26,22]. Therefore, it can be theoretically deduced that a certain mutation is associated with epigenetic variation to some extent. This study detected nearly 32% epigenetic differentiation between the mutants and their original mother 'Fuji', showing many types of methylation pattern variants, and indicating that the variation in DNA methylation might be involved in the occurrence of 'Fuji' mutations. Our results explain the differentiation of epi-phenotypes, related to the similar genome structure, highlighted by Guarino et al. [60]. This may also verify the hypothesis of significant plant promotion to environmental adaptation, by phenotype diversification, increasing the probability of the plant to thrive when faced with a changing environment [3].
Different levels and degrees of methylation occur in higher plants to maintain normal plant development [59]. In general, the total genome DNA methylation level in different species detected by MSAP analysis is between 4.7% and 60.0% [61]. In this research, the DNA methylation level detected by MSAP was approximately 36%, similar to the findings of Li [62] in apple, in accordance with the relatively stable characteristics of DNA methylation levels in this species.
CG and CHG are two main forms of DNA methylation at CCGG sites. Their relative proportions in the organism vary among different species, and the majority of available reports present results show that the proportion of the former is greater than of the latter [11]. In this study, the levels of CG and CHG methylation were found to be 22.311% and 13.599%, respectively. The level of the former was significantly higher than that of the latter, indicating full methylation derived was the main form of DNA methylation in the 'Fuji' line, which conformed to the results in the majority of plant species.
Natural mutations in plants are usually accompanied by changes in the levels of genetic and epigenetic. AFLPs and MSAPs are the two marker techniques commonly used for the detection of plant genetic and epigenetic variation [13][14]. Expected heterozygosity (He), also known as gene diversity [63], is an important parameter for the measurement of the general genetic diversity among plant groups and shows wide applicability to any system of polyploidy, self-crossing or asexual reproduction in populations [64]. The lower the genetic diversity, the higher the degree of the genetic homogeneity in a population. It is generally known that a heterozygosity value higher than 0.5 indicates genetic diversity between population individuals [65]. In this study, the average genetic diversity values of the two types of molecular markers were 0.186 (AFLP) and 0.266 (MSAP), suggesting that the genetic and epigenetic diversity levels of the mutant groups in the 'Fuji' line were relatively low. Comparatively, in the three examined mutant groups, the epigenetic expected heterozygosity of the Color mutant group was higher than those of the other two groups (P<0.01), implying its relatively high epigenetic diversity. Additionally, in combination with the results from UPGMA cluster, STRUCTURE analysis, and outlier loci detection, the Color mutant group displayed unique epigenetic characteristics, which contributed largely to the epigenetic differentiation between the three mutant groups in the 'Fuji' line. This phenomenon might be related to the degree of difficulty in bud mutation of certain traits, the breeding and selection objectives in a given period, or the degree of artificial intervention.
Epigenetic variations may exhibit phenotypic differences in different environments or growth periods and even in different tissues and plant organs [61]. However, in some species such as walnut, the global genomic DNA methylation level in different tissues and organs does not differ significantly [65]. The leaf can provide information about epigenetic modification and adaptation in response to different environmental conditions [60]. In our study, the leaves were considered as the plant material, as such, the results and conclusions of the epigenetic investigation were limited to the analysis of DNA matrix purified from particular tissue collection. The differences between the Early-maturation mutant group and the Color mutant group were mainly reflected in the performance of fruits. Therefore, if fruit samples of these two mutant types were used as the examined material, the results would not be similar to those obtained for leaves and remain inconclusive since no similar research has previously been conducted in apple. The tissue specificity of epigenetic modifications is also unknown.
Additionally, we found that the mean hypermethylation and demethylation frequencies of CHG type were both significantly higher than those of CG type, indicating that variation in CHG methylation pattern played a key role in DNA methylation pattern variation in the 'Fuji' lines. In plants, cytosine methylation is a context-dependent process [9]. Our findings supported it. DNA methylation status reflects the outcome of the dynamic regulation of establishment, maintenance, and active-removal activities. CG and CHG methylation are catalyzed by various different enzymes and are under control by different pathways [9]. Hence, in view of our findings, the mechanism of epigenetic involvement in the occurrence of 'Fuji' bud mutation could be further explored, focusing on the act of enzymes and detailed pathways related to CHG methylation.
In conclusion, the present study uncovered abundant changes in methylation levels and patterns between not only bud mutants and their mother 'Fuji' but also bud mutants, indicating that it may be possible that epigenetics mediated by DNA methylation was involved in the occurrence of the 'Fuji' bud mutation line. The epigenetic mechanism of the Color bud mutant group was unique, which could be the focus of further research. The correlation between epigenetic variation and genetic variation was weak, showing their independence, which provided further support for using 'Fuji' line as ideal sets for epigenetic studies in the future.