Whole-Genome DNA Methylation Associated With Di ﬀ erentially Expressed Genes Regulated Anthocyanin Biosynthesis Within Flower Color Chimera of Ornamental Tree Prunus mume

: DNA methylation is one of the best-studied epigenetic modiﬁcations involved in many biological processes. However, little is known about the epigenetic mechanism for ﬂower color chimera of Prunus mume (Japanese apricot, mei). Using bisulfate sequencing and RNA sequencing, we analyzed the white (FBW) and red (FBR) petals collected from an individual tree of Japanese apricot cv. ‘Fuban Tiaozhi’ mei to reveal the di ﬀ erent changes in methylation patterns associated with gene expression leading to signiﬁcant di ﬀ erence in anthocyanins accumulation of FBW (0.012 ± 0.005 mg / g) and FBR (0.078 ± 0.013 mg / g). It was found that gene expression levels were positively correlated with DNA methylation levels within gene-bodies of FBW and FBR genomes; however, negative correlations between gene expression and DNA methylation levels were detected within promoter domains. In general, the methylation level within methylome of FBW was higher; and in total, 4,618 di ﬀ erentially methylated regions (DMRs) and 1,212 di ﬀ erentially expressed genes (DEGs) were detected from FBW vs. FBR. We also identiﬁed 82 DMR-associated DEGs, and 13 of them, including PmBAHD , PmCYP450 , and PmABC, were playing critical roles in phenylalanine metabolism pathway, glycosyltransferase activity, and ABC transporter. The evidence exhibited DNA methylation may regulate gene expression resulting in ﬂower color chimera of Japanese apricot.


Introduction
Prunus mume Sieb. et Zucc. (Japanese apricot, mei), being one of the traditional ornamental plants, has been domesticated for more than 3000 years in China [1,2]. The woody plant with great cultural and economic value is well known for its colorful flowers with their unique fragrance, the traits of which are the important breeding objectives for genetic improvement [2]. In the 1940s, breeders discovered varieties of flower color chimera with high ornamental and commercial value. During the past 70 years, at least five new cultivars of Japanese apricot with chimeric flowers for modern landscaping application were released through selective breeding; however, there has been little exploration into the formation mechanism of the flower color chimera.
At the cytological level, plant chimera organism is believed to be the result of the heterogeneous primordial cells proliferating mechanically in the shoot apical meristem [3]. In other word, the cell lines

Plant Materials
In the present study, both white petals (FBW) and red petals (FBR) were collected from one tree of 'Fuban Tiaozhi' mei, planted in Jiufeng Forest Farm of Beijing Forestry University, Haidian, Beijing, when the flowers were ready to bloom (Figure 1a). Then, the FBW, and FBR samples were frozen in liquid nitrogen immediately and stored at −80 • C.

Anthocyanin Determination
Part of each sample was finely ground, and 0.1 g powdered material was transferred to a calibration test tube (20 mL) containing 10 mL methanol with 1% HCl. The mixtures were vortexed for 1 min, and then the pigments were ultrasonically extracted by using a KQ 2200B ultrasonic cleaner (Kunshan, Jiangsu, China) at 4 • C for 15 min under dark condition. Then, the mixtures were centrifuged for 10 min at 12,000 rpm and 4 • C, and the supernatants were collected and stored at the refrigerator overnight at −20 • C. After filtration through 0.22 µm nylon membranes, the extracting solution was analyzed using a UV-visible spectrophotometer (UV-2450, Shimadzu Corporation, Tokyo, Japan). The anthocyanins concentrations were estimated according to the absorbance value and the coefficient of regression acquired by standard scale measurements [33]. The experiment was repeated biologically three times.

DNA Extraction and BS-Seq
Genomic DNA of FBW and FBR petals was extracted using a DNeasy Plant Mini kit (QiaGen, Shanghai, China), respectively. The degradation and contamination of the genomic DNA was monitored on 0.1% agarose gels, using a NanoPhotometer spectrophotometer (Implen, Calabasas, CA, USA) and a Qubit DNA Assay kit on a Qubit 2.0 fluorometer (Life Technologies, Carlsbad, CA, USA). Then, the DNA was mechanically fragmented in a Covaris S220 ultrasonicator, and the fragments were subjected to end reparation, adenylation, and methyl-treated adapter ligation. Next, the nonmethylated cytosines were transformed into uracil using an EZ DNA Methylation-Gold kit (Zymo Research). Utilizing the DNA fragments as templates, PCR amplification was carried out, and the products forming a sequencing library were quantified, and the insert size was assayed. The sequencing process was performed on the Hiseq 2500 (Illumina) platform to generate raw reads that were preprocessed with Trimmomatic v0.36 software and FastQC analysis [52]. Finally, the clean reads were obtained for subsequent analysis.

Sequence Mapping and Differentially Methylated Regions (DMRs) Detection
The clean reads were aligned to the reference genome of P. mume [53][54][55]. And the sequencing depth and coverage were estimated according to duplicates that aligned to a unique genomic region [54,55]. The non-conversion rate (r) of bisulfate treatment was defined as a rate of the number of sequenced cytosines at all cytosine reference positions divided by the number in the lambda genome. We also defined m C and C as numbers of methylcytosine and cytosine, respectively. Then, the methylation level was calculated by an equation described as ML( m C) = reads( m C)/(reads( m C) + reads(C)). Meanwhile, the parameter ML( m C) was corrected to ML(corrected) as ML(corrected) = (ML( m C) − r)/(1 − r) [56]. The mCG, mCHG, and mCHH contexts (methylcytosine occurs at CG, CHG, and CHH regions, respectively) and their densities, and distributions in each chromosome were analyzed using previous methods, respectively [39,45,57,58]. Between samples, the global methylation levels and distributions were also tested [42]. And DMRs and differentially methylated loci were identified by using a DSS software [59][60][61].

Functional Annotation of DMR-Related Genes
The GOseq R soft package was used for Gene Ontology (GO) enrichment analysis of DMR-related genes that were aligned with GO Ontology (p-value < 0.05) [62]. While the Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis was used to explore pathways of DMR-associated genes [63]. We also used KOBAS software to test DMR-related genes' statistical enrichment in the KEGG pathways [64].

RNA Isolation and Sequencing
The FBW and FBR materials were processed for RNA extraction using an RNeasy Plant Mini kit (QiaGen), respectively. RNA quality and quantity were determined using 1% agarose gels and a NanoPhotometer spectrophotometer (Implen, Calabasas, CA, USA), respectively. Also, an RNA Nano 6000 Assay kit of the Bioanalyzer 2100 system (Agilent, Carpinteria, CA, USA) was used to assess the RNA integrity. A total amount of 3 µg RNA was used as the input material for sequencing library construction. By using a NEBNext Ultra RNA Library Prep Kit for Illumina (NEB, Ipswich, MA, USA), first-and second-strand cDNA synthesis was generated, and the fragments were purified using an AgencourtAMPure XP system (Beckman Coulter, Brea, CA, USA). Then, the purified fragments were subjected to 3 -end adenylation and the addition of a poly-A tail and a NEBNext adapter. The libraries and 150-200 bp cDNA fragments were selected and assessed by the Bioanalyzer 2100 system. The libraries were sequenced on an IlluminaHiSeq system to generate 125-bp/150-bp paired-end reads (clean bases > 6 Gb).

Methylation Modification of Gene Expression
In this section, we explore the relationship between cytosine methylation and gene expression level according to methylomes and transcriptomes. The correlation between methylation level and gene expression was progressively analyzed on the four steps below: firstly, we mapped the gene methylation level and gene expression densities onto each of the chromosomes of P. mume and investigated the methylation levels at different regions (promoter, gene-body region, 2 kb upstream and downstream regions, transcription start site, and transcription end site) that affected genes expression [47,55,68]; secondly, we identified DEGs (including promoters, and 2 kb upstream and downstream regions of gene-body) that modified with CG, CHG, and CHH methylation patterns according to the methods described by Ng et al. [69] and Zhou et al. [70]; thirdly, the expression levels for DMR-related genes (including CG, CHG, CHH patterns) were detected [71,72]; and at last, we identified the set of DMR-related genes that overlapped with DEGs that had been subjected to GO and KEGG enrichment.

Content of Anthocyanins
As shown in Figure 1b, the anthocyanins content in white petals of 'Fuban Tiaozhi' mei (FBW; 0.012 ± 0.005 mg/g, fresh weight) was less than anthocyanins content in red petals of 'Fuban Tiaozhi' mei (FBR; 0.078 ± 0.013 mg/g, fresh weight). And it indicated that the anthocyanins were the main pigment observed directly by the eyes in petals of 'Fuban Tiaozhi' mei.

DNA Methylation Landscape in Japanese Apricot cv. 'Fuban Tiaozhi' Mei
In total, 12.48 and 12.40 Gb clean bases were produced within sequenced genomes of FBW and FBR, respectively. It was estimated that approximately 25 million clean reads (~59%) of each methylome could be mapped to the reference genome (Supplementary Table S1). It also revealed that the maximum coverage was displayed at a sequencing depth of approximately 25× to 30× of the reference genome ( Figure 2a). Within methylomes of FBW and FBR, 12.88-13.36% of cytosine sites (C sites) were methylated, including three patterns revealed by mCG (35.43-35.61%), mCHG (19.25-19.26%), and mCHH (8.65-9.30%) contexts (Table 1). Meanwhile, it was shown that chromosome Pm3 was highly methylated (Figure 2b). Among the total mC sites, the relative proportions of mCG, mCHG, and mCHH contexts were 26.83-27.97%, 20.41-21.06%, and 50.86-52.74%, respectively ( Figure 2c). Simultaneously, there existed a tendency toward mCG and mCHG contexts with high methylation levels (70-100%); however, mCHH context performed the opposite trend, with low methylation level (0-30%) taking up a large proportion ( Figure 2d). As shown in Figure 3a,b, the DNA methylation (divided into mCG, mCHG, and mCHH contexts) density and level of 'Fuban Tiaozhi' mei within each chromosome were consistent with the distribution of low-gene-density areas.

Difference in DNA Methylation Level between the Two Methylomes
The status and level can vary among different genomes, although the distributions of DNA methylation density and level showed similarly tendencies within corresponding chromosomes of FBW and FBR. We first calculated that the methylation levels of CG (2.43-86.26%), CHG (1.42-60.32%), and CHH (0.48-4.81%) sites within genome of FBW were different from methylation levels of CG (1.95-86.21%), CHG (1.15-59.57%), and CHH (0.39-4.46%) sites within genome of FBR. Additionally, variation in the majority of mCG, mCHG, and mCHH contexts with strongly methylated regions occurred similarly at paired chromosomal positions (Figure 3c). Then, we analyzed that gene functional regions were differentially methylated and TEs domains presented high methylation levels. The mCHG and mCHH levels showed low methylation levels within exon and intron domains; however, CG sites were highly methylated within exon and intron regions. In the promoter region, mCG and mCHG displayed a declining trend (Figure 3d). The mCG context displayed an upward tendency within gene-body domains compared with 2 kb up-and down-streams within transcribed regions, but mCHG and mCHH contexts showed the opposite trend (Figure 3e). In general, the methylation levels were higher within methylome of FBW (Figure 3d,e).

DMR-Related Genes Participate in Color Regulation Pathway
We detected 4,618 DMRs within the two methylomes in order to find DMR-related genes that take part in flower color regulation. The CG-DMRs were methylated with levels between 1.59% and 97.92%, and averaged at 49.49% and 48.22% within FBW and FBR samples, respectively. CHG-DMRs were averaged at 39.31% and 39.13% within genomes of FBW and FBR, respectively. DNA methylation levels of both of CG-and CHG-DMRs were well-distributed approximately. However, the majority of CHH-DMRs were distributed with low methylation levels ( Figure 4a). Heat maps and circos plots (Figure 4b,c) showed wide variations between methylomes of FBW and FBR, respectively. Simultaneously, we annotated the DMRs within gene functional regions that presented the minority of DMRs were distributed in the TSS (transcriptional start site) and TES (transcriptional end site) areas. For the CG-DMRs, the hypermethylated ones were less than the hypomethylated (FBW vs. FBR) ones within promoter and exon domains; while the hypermethylated CHG-DMRs were larger than the hypomethylated CHG-DMRs. Obviously, the DMRs displayed more hypermethylated CHH contexts than that of hypomethylated ones all through the gene function regions (Figure 4d).

DEGs Predicted by Comparative Transcriptomics Analysis
On the basis of the RNA sequencing, 52.79 and 56.32 million clean reads were generated for FBW and FBR (Table 2), respectively. In transcriptome of FBW, 90.5%, 0.9%, and 8.6% of the clean reads were mapped to the exon, intron, and intergenic areas, respectively. In comparison, 88.3%, 1.2%, and 10.4% of the clean reads were annotated within exon, intron, and intergenic domains of FBR transcriptome (Figure 6a). The numbers of clean reads aligned to "+" and "−" strands of the reference genome were approximately equal (Figure 6b, Table 2). It was calculated that the gene expression level of FBR was significant higher than the expression level of FBW genes (p = 0.027; Figure 6c). Selected by DEGseqv1.12.0 package, 1,212 DEGs (up-regulated, 704; down-regulated, 508) presenting differentially expressed patterns (Figure 6d,e) were detected from FBW vs. FBR. According to GO enrichment (p < 0.05), these up-and down-regulated genes were enriched into 162 and 133 GO terms, respectively (Supplementary Figure S3). We also used KEGG pathway analysis to reveal further biological functions of the DEGs, with those that presented 54 and 4 up-regulated DEGs taking part in biosynthesis of secondary metabolites (ko pmum01110) and ABC transporters (ko pmum02010) pathways, respectively. Interestingly, 4 and 38 down-regulated DEGs were found in flavonoid biosynthesis (ko pmum00941) and biosynthesis of secondary metabolites (ko pmum01110) pathways, respectively (Figure 7).

Correlation Between DNA Methylation and Gene Expression Levels
In order to explore the relationship between gene expression and DNA methylation levels, the Pearson correlation test (2-tailed) was performed. It was found that the gene expression levels were positively correlated with DNA methylation levels (p < 0.001) in FBW (r = 0.062) and FBR (r = 0.099) genomes within gene-body regions. However, negative correlations between gene expression and DNA methylation levels (p < 0.001) of FBW (r = −0.137) and FBR (r = −0.118) were detected within promoter domains (Figure 8a). In the purpose of understanding the relationship between gene expression and mC contexts in detail within gene-body and its 2 kb up-and down-stream areas, we detected the distribution of the gene expression levels that divided into four groups (high level, medium level, low level, and none expressed). The cytosine sites were lacking DNA methylation at TSS and TES sites. In gene-body regions, the mCG/CG levels were increasing, and the none-expressed genes were methylated with low mCG/CG levels. However, the methylation levels of mCHG and mCHH contexts were decreased in gene-body domains of all of the none-expressed and expressed genes. Simultaneously, these none-expressed genes along with highly mCHG/CHG and mCHH/CHH methylated levels were presented. In the 2 kb up-or down-stream of the gene-body areas, high methylation levels were detected at CG and CHG sites of the none-expressed genes. Therefore, gene expression levels were positively correlated with mCG/CG levels, but they were negatively correlated with mCHG/CHG and mCHH/CHH levels. And within mCHG context, none-expressed genes always displayed highly methylated. We could also observe that the distributions of gene expression levels within mC contexts presented the similar tendencies of FBW vs. FBR (Figure 8b). On the other hand, we also detected the frequency distribution of genes with different methylation levels in gene-body and promoter areas classified into six groups that were projected into a coordinate system of gene expression level (x-coordinate) and frequency (y-coordinate). At low expression levels estimated as log 2 (FPKM + 1) <1, more genes lacked DNA methylation within gene-bodies, but more genes rich in DNA methylation within promoters were detected. However, when the threshold value setting was 2 < log 2 (FPKM+1) < 8, more and minority genes were highly methylated within gene-bodies and promoters, respectively (Figure 8c). Thus, it suggested that genes rich in mCG sites within gene-bodies and genes lacking methylation within promoters were highly expressed in whole genomes of FBW and FBR. Gene expression levels were classified into four groups: none expression (purple; FPKM < 1); low expression (blue; 1 ≤ FPKM < FPKM_25%); medium expression (green; FPKM_25% ≤ FPKM < FPKM_75%); and high expression (red; FPKM ≥ FPKM_75%). FPKM stands for fragments per kilobase of exon per million fragments mapped, and FPKM_25% and FPKM_75% refer to values at the boundary of 25th and 75th percentiles of expression levels, respectively. Next, the different gene regions (gene-body and 2 kb up-and down-stream) were divided into 50 bins, and the methylation levels of each bin were averaged. The xand y-axes represent gene-body regions and DNA methylation levels by mCG, mCHG, and mCHH contexts, respectively; TSS, transcription start site; TES, transcription end site; (c) Frequencies of expressed genes that anchored by mC sites within gene-body and promoter domains of each sample. DNA methylation levels were grouped into six sections: group 1 (red; 0 < methylation level < level_20%); group 2 (yellow-green; level_20% ≤ methylation level < level_40%); group 3 (green; level_40% ≤ methylation level < level_60%); group 4 (turquoise; level_60% ≤ methylation level < level_80%); group 5 (blue; methylation level ≥ level_80%), and unmethylated group (purple; methylation level = 0). Level_20%, _40%, _60%, and _80% represent values at the boundaries of 20th, 40th, 60th, and 80th percentiles of methylation levels. The xand y-axes represent gene expression levels and gene frequencies, respectively.

DMR-Associated DEGs Took Part in Anthocyanin Formation
We also detected the methylation status of the DEGs, and it showed that the mCG, mCHG, and mCHH densities and mCG/CG, mCHG/CHG, and mCHH/CHH levels presented similar distribution tendencies. The decreasing mCG densities and mCG/CG levels within DEG's promoters were the least among the three functional regions (p < 0.001) at CG sites; however, the densities and levels were the richest of DNA methylation (p < 0.001) at both the CHG and CHH site in promoter areas. Unexpectedly, information was obtained that the methylation densities and levels of DEGs between FBW and FBR were not significantly different; however, the down-regulated DEGs were higher methylated than up-regulated DEGs (p ≤ 0.010) within each functional region, except for the CG context within exon region (Figure 9a,b). It suggested that the expression of DEGs was affected by DNA methylation. We then detected the DMR-associated DEGs to investigate the patterns of the DEGs regulated by cytosine methylation.  Table S2). Similarly, 14 promoter hypermethylated and 9 promoter hypomethylated DMR-related genes were associated with up-regulated DEGs, as well as 17 and 7 down-regulated DEGs were performed hypermethylation and hypomethylation within promoters, respectively (Figure 9c, Supplementary Table S3). After functional annotation, we identified that 13 putative hyper-or hypomethylated DEGs at CHH sites within gene-bodies (Supplementary Table S2) and promoters (Supplementary Table S3) were playing critical roles in phenylalanine metabolism pathway (Pm009966, Pm011003, Pm011258, Pm017164, Pm019289, Pm020893, and Pm025210), glycosyltransferase activity (Pm000414, Pm001802, Pm004453, Pm020721, and Pm027780), and ABC transporter (Pm013365).

Discussion
Flower color is an important quality character for ornamental plants, and the uncommon color individuals are favored. Therefore, 'Fuban Tiaozhi' mei, an attractive cultivar, can produce one tree with both red and white petals and has been widely applied in landscape architecture. Recently, people considered that the chimeric color petals are induced by diverse anthocyanin contents and differentially expressed genes [9,10]. And cytosine methylation is also associated with the coloration of flower petal, fruit skin, and so forth, by regulating gene expression [34][35][36][37]. However, little is known about the epigenetic mechanism regulating flower color chimera of 'Fuban Tiaozhi' mei. The present research shed light on the relationship between methylcytosine and gene expression that takes part in the pathway of anthocyanins, the concentrations of which were significantly different between white and red petals.
Phenotype is determined by gene expression, and previous studies showed that structural and transcription factor genes control the plant coloration. Provably, anthocyanins start from the phenylpropanoid biosynthetic pathway when the phenylalanine ammonia lyase (PAL) catalyzes phenylalanine into cinnamic acid, and 4-coumaroyl CoA ligase (4CL) activates 4-coumaric acid into 4-coumaroyl-CoA [10,17]. Then, the flavonoids metabolites come around, and the branch pathway is regulated by a series of enzymes including CHS, CHI, F3H, and DFR [17,[21][22][23][24]. ANS, involved downstream of the flavonoids-branched pathway, is one key enzyme to transform leucoanthocyanidin to colored anthocyanidin before the final glycosylation steps controlled by UFGT [10,25,26]. For instance, in Petunia hybrida with bicolored flowers, enhancement of CHS expression induces blue or red coloration, respectively, in blue-white or red-white sectors of variegated flowers [8,73]. By using RNA sequencing method, it has demonstrated that the different expression of CHS, CHI, F3H, DFR, LDOX, ANS, and/or UF3GT genes taking part in the anthocyanin biosynthetic pathway was responsible for purple flesh coloration in a Dioscoreaalata cultivar, bicolored petals in lily, and variegated petals in 'Fuban Tiaozhi' mei [10,74,75]. Meanwhile, the expression of regulatory genes MYB, bHLH, and WD40 were significantly correlated with anthocyanin content in black carrot taproot [19]. In our work, we also found putative genes including PmCHS, PmANS, and PmUFGT, and so forth showed obvious consistency with the changing pattern of anthocyanin. These up-and down-regulated DEGs take part in flavonoid biosynthesis (ko pmum00941) and biosynthesis of secondary metabolites (ko pmum01110). The evidence suggested a possible mechanism that genes within FBW and FBR petals were differently expressed, resulting in color chimera flowers of Japanese apricot.
According to previous research, the gene expression is also regulated by epigenetic modifications results in variations; and DNA methylation is one of the best-studied modifications and closely interwoven with the process of endogenous gene transcription [12]. By using the single-base resolution technique, we studied the landscape of the methylomes of 'Fuban Tiaozhi' mei. In general, the high methyl-cytosine density and level within each chromosome had adistribution consistent with that of low-gene-density areas. And mCG context played an upward tendency within gene-body domains, but mCHG and mCHH contexts showed the opposite trend. Similarly, Zhang et al. [76] described that body-methylated genes were constitutively expressed at a higher level, whereas promoter-methylated genes tend to be expressed in a tissue-specific manner. We also determined that methylcytosine levels were positively correlated with gene expression within gene-body regions, but negative correlations were detected within promoter domains. The result is consistent with the finding that genes methylated in transcribed regions are highly expressed and constitutively active in A. thaliana [76] and the fact that mCG methylation is often linked to increased gene expression [77][78][79].
Further than that, we investigated the relationship between methylcytosine and gene expression level, and DMR-related genes including DEGs within FBW vs. FBR. It was found that the methylation levels were higher within methylome of FBW containing a lower level of anthocyanin. And 4,618 DMRs were associated with 1,029 and 1,383 genes within the functional regions of gene-bodies and promoters, respectively. In total, 82 DMR-associated DEGs (44 up-regulated, 38 down-regulated) were detected and annotated, such as down-regulated genes PmBAHD (Pm011003) and PmCYP450 (Pm025210), and PmABC (Pm013365) were enriched in phenylalanine metabolism pathway and ABC transporter, respectively. Cytochrome P450, which the flavonoid 3 -hydroxylase (F3 H) and flavonoid 3 ,5 -hydroxylase (F3 5 H) belong to, determining the hydroxylation pattern of the B-ring is considered to play a key role in flavonoid biosynthesis and anthocyanins compounds essentially [7,[80][81][82][83]. ABC transporter was determined to take part in membrane transport of flavonoid [84,85] and anthocyanin pigment [86,87]. Our research suggested these genes could be regulated by methylcytosine modification, though little knowledge was explored. Previously, the demethylation of the 5 -upstream promoter region of OgCHS generated mosaic red anthocyanin in lip crests, sepals, and petals of yellow flowers of Oncidium 'Gower Ramsey' [35]. And MYB, bHLH74 genes with demethylation are associated with the formation of red-skinned pears and apples by inhibiting anthocyanin synthesis [36,37,49]. Thus, the evidence suggested that DNA methylation levels regulated key genes participating in anthocyanins biosynthesis, leading to flower color chimera of 'Fuban Tiaozhi' mei.

Conclusions
In present work, we shed light on the regulatory mechanism for flower color chimera of ornamental tree Japanese apricot by using single-base resolution methylome detection (BS-seq) and transcriptome sequencing (RNA-seq) methods. In phenotype, the white and red petals from an individual tree resulted from the significant difference in anthocyanins accumulation. According to RNA-seq, we found DEGs taking part in flavonoid biosynthesis, secondary metabolites, and ABC transporters pathways. It was also revealed that the DEGs displayed various DNA methylation statuses, and down-regulated DEGs were higher methylated than up-regulated DEGs, in general. Furthermore, we detected DMR-related genes associated with DEGs, such as PmBAHD, PmCYP450, and PmABC, contributed to anthocyanins accumulation. The evidence suggested that genome DNA methylation associated with DEGs regulated bicolor flowers of 'Fuban Tiaozhi' mei. These findings will enhance our understanding of epigenetics in flower color regulation and be beneficial to molecular-assisted selection (MAS) for ornamental plants.
Sequencing Data: All of the BS-seq and RNA-seq data were submitted to NCBI (https://www.ncbi.nlm.nih.gov/) with an accession number PRJNA597717.