Metabolic profiling during pear flesh development

To survey the changes in metabolites during pear flesh development, cv. Dangshansuli fruits were firstly collected at different time points, specifically at 4 (S1), 6 (S2), 8 (S3), 10 (S4), 12 (S5), 14 (S6), 16 (S7), 18 (S8), 20 (S9), 22 (S10), and 24 (S11) weeks after flower blooming (Fig. 1a). The typical characteristics of fruit development were observed during this period (Additional file 1: Fig. S1). These included a gradual increase in single fruit weight, longitudinal and transverse diameters, and soluble sugar content. Ethylene content showed a substantial increase at S10 and S11. Stone cell content increased from S1 to S3, but then gradually decreased. Flesh firmness decreased gradually from S8 to S11, while soluble solids increased gradually from S4 to S11. Secondly, metabolic measurements identified a total of 492 metabolites in pear flesh. On average, approximately 456 metabolites were detected at each stage, with 371 metabolites being commonly present across all stages (Additional file 2: Table S1). Principal component analysis confirmed distinct differences among the fruit flesh at the 11 stages (Fig. 1b). By comparing the relative content between each pair of stages (Additional file 2: Table S2), a total of 449 differentially accumulated metabolites (DAMs) were identified. These DAMs belonged to 32 different classes, including phytohormones, anthocyanins, amino and their derivatives, carbohydrates, flavone, flavonol, hydroxycinnamoyl derivatives, lipids, nucleotide and its derivatives, and organic acids, and presented eight clusters (I → VIII) based on their accumulation trends (Fig. 1c, Additional file 2: Table S3). Clusters I, II, and IV showed higher levels of metabolites in the fruit flesh at the early stages compared to later stages (Additional file 1: Fig. S2), indicating their association with early flesh development. An example is the sinapyl alcohol in cluster II, which was associated with the lignification of stone cells that are enriched in pear flesh from 15 to 55 days after flowering [42]. Clusters VII and VIII exhibited higher levels of metabolites in the fruit flesh at S7 and S6, respectively (Additional file 1: Fig. S2), indicating that these metabolites, including L-Alanine, L-Phenylalanine, L-Proline, and L-Valine in cluster VII, may be associated with fruit enlargement [43, 44]. It is reported that amino acids are significantly accumulated during pear fruit enlargement [45]. Clusters III, V, and VI showed higher levels of metabolites in the flesh at S10 and S11 compared to other stages (Additional file 1: Fig. S2), indicating their association with fruit ripening. For example, abscisic acid (ABA) in cluster V can increase the content of soluble sugars (including sucrose and glucose) and promote ethylene release, thereby promoting pear fruit ripening [46].

Fig. 1

Dynamics of metabolite during pear flesh development. a The pear fruits collected from 11 stages. b Principal component analysis of metabolome data in the fruit flesh at 11 stages. c Clustering analysis of metabolic profiling grouped the differentially accumulated metabolites into eight clusters, I → VIII. Z-score standardized values of each metabolite across all 11 stages were used for clustering analysis. The color bar denotes an increase in the contents of each metabolite, transitioning from blue to red color. d The number and cluster of the metabolites that were significantly increased in the next stage. S1 to S11 indicate the pear fruit at 4, 6, 8, 10, 12, 14, 16, 18, 20, 22, and 24 weeks after flower blooming, respectively

To mine the metabolites that promote flesh development, we further analyzed the DAMs between two adjacent stages (Additional file 1: Fig. S3a) and found that transition from one type of metabolite to another was more frequent in the early stages of flesh development compared to the later stages (Additional file 1: Fig. S3b). A total of 105 metabolites belonging to 27 classes were significantly accumulated in the flesh at the next stage compared to the previous stage, suggesting their potential importance in pear flesh development. Interestingly, the metabolites between any two adjacent stages belonged to at least two clusters (Fig. 1d). Moreover, the transition from S1 to S2 and from S2 to S3 involved 14 and 11 metabolite classes, respectively, while the transition from S9 to S10 and from S10 to S11 involved seven and four metabolite classes, respectively (Table 1). This result indicates that early flesh development is more dependent on the metabolite biosynthesis compared to fruit ripening. During fruit enlargement (from S3 to S9) [47], the transition from S4 to S5 and from S6 to S7 involved 13 and 12 metabolite classes, respectively, while the transitions from S3 to S4, from S5 to S6, from S7 to S8, and from S8 to S9 involved seven, five, five, and eight metabolite classes, respectively (Table 1). These results indicate that the biosynthesis of metabolites from S4 to S5 and from S6 to S7 may play important roles in fruit enlargement.

Table 1 The classes of metabolites involved in the transition from one stage to the nextGene regulation of pear flesh metabolism

To establish a relationship between gene activity and pear flesh metabolism, mass spectrometry was conducted to investigate the proteomes at each stage. A total of 5106 proteins were detected in all samples combined. On average, approximately 4074 proteins were detected at each stage, with 2860 proteins being consistently present across all stages (Additional file 2: Table S1). By comparing protein abundance between stages, a total of 3469 differentially expressed proteins (DEPs) were identified in the fruit flesh (Fig. 2a). Integrative analysis of the proteome and metabolome showed that 1038, 1012, 167, 843, 927, 287, 426, and 162 DEPs were positively (Pearson coefficient > 0.85 and false discovery rate < 0.05) were positively correlated with metabolites in clusters I → VIII, respectively, while 834, 883, 68, 509, 922, 60, 458, and 58 DEPs were negatively (Pearson coefficient <  − 0.85 and false discovery rate < 0.05) correlated with metabolites in the same clusters (Additional file 2: Fig. S4a). Overall, 2513 DEPs were correlated with 423 DAMs (Additional file 2: Table S4). Among these, 1581 DEPs were correlated with the metabolites in a single cluster, while 932 DEPs were correlated with metabolites in at least two clusters (Additional file 1: Fig. S4b).

Fig. 2

Integrative analysis of metabolome, proteome, and transcriptome data of pear flesh. a Clustering analysis of proteome data shows the three reliable replicates of each sample at any stage. Z-score standardized value of each protein across all 11 stages (from S1 to S11) were used for clustering analysis. b Clustering analysis of transcriptome data shows the three reliable replicates of each sample at any stage. Z-score standardized value of each gene across all 11 stages was used for clustering analysis. The color bar indicates the increasing expression levels of both protein and gene from blue to red. c Peonidin-3-galactoside biosynthesis pathway. d Sugar biosynthesis and metabolic process. e Abscisic acid biosynthesis pathway. The pathway genes with the background of red color were positively correlated to the target metabolites (highlighted by blue box). f The metabolic regulatory network for each metabolite cluster. I → VIII indicate the different clusters

To compensate for the limited number of proteins detected through mass spectrometry, transcriptome sequencing was conducted to investigate the expression patterns of all predicted genes in the pear fruit flesh. A total of 3.06 billion raw reads were generated from all samples (Additional file 2: Table S5), with 931.09 million clean reads successfully assembled into 33,051 genes. On average, approximately 28,748 genes were detected at each stage, and among them, 24,806 genes were consistently present across all 11 stages (Additional file 2: Table S1). By comparing the expression levels between each pair of stages, a total of 22,055 differentially expressed genes (DEGs) were identified in fruit flesh (Fig. 2b). Integrative analysis of the transcriptome and metabolome showed that 9524, 7120, 382, 6257, 2757, 505, 1818, and 433 DEGs were positively (Pearson coefficient > 0.85 and false discovery rate < 0.05) correlated with metabolites in clusters I → VIII, respectively, while 2138, 2564, 80, 1144, 5372, 202, 4050, and 481 DEGs were negatively (Pearson coefficient <  − 0.85 and false discovery rate < 0.05) correlated with metabolites in the same clusters (Additional file 2: Fig. S5a). Overall, 14,399 DEGs, including 2218 DEP-coding genes, were correlated with 439 DAMs (Additional file 2: Table S6). Among these DEGs, 9139 DEGs were correlated with metabolites in a single cluster, while 5260 DEGs were correlated with metabolites in at least two clusters (Additional file 1: Fig. S5b).

After removing the debatable genes presenting an inconsistent correlation with a metabolite in transcript and translation levels, a comprehensive gene-metabolite database was constructed (Additional file 3: Table S7). In this database, 220 metabolites could be detected in the Kyoto Encyclopedia of Genes and Genomes (KEGG) database. However, the transcriptome data of fruit flesh only included survey information for the pathway genes of only 101 metabolites. Out of these 101 metabolites, 63 positively or negatively correlated with their pathway genes (Additional file 3: Table S8). To further evaluate the effectiveness of this database, we conducted tests using previously reported metabolic pathways. For instance, peonidin-3-galactoside, a component of anthocyanin, was positively correlated with various anthocyanin biosynthetic genes, including Cinnamate 4-hydroxylase, 4-Coumarate:CoA ligase, Chalcone synthase, Chalcone isomerase, Flavonoid 3-hydroxylase, Leucoanthocyanidin dioxygenase, UDP-glucuronosyl/UDP-glucosyltransferase, O-methyltransferase, and components of the MYB-bHLH-WD40 complex (Fig. 2c, Additional file 3: Table S7). Additionally, sucrose and glucose are positively correlated with sugar metabolism-related genes including Sorbitol dehydrogenase, Fructokinase, Hexokinase, Phosphoglucomutase, UTP-glucose-1-phosphate uridylyltransferase, sucrose-phosphate synthase, and tonoplast monosaccharide transporter (Fig. 2d, Additional file 3: Table S7). Moreover, ABA was positively correlated with its biosynthetic genes, including Phytoene synthase (PSY), Zeta-carotene desaturase (ZDS), beta-ring hydroxylase (HYB), Short-chain dehydrogenase/reductase (SDR), and Aldehyde oxidase (AO) (Fig. 2e, Additional file 3: Table S7). These results suggest that the gene-metabolite database accurately identifies crucial genes involved in previously reported metabolic pathways and thus provides valuable insights into the metabolic regulatory network within fruit flesh.

To survey the potential regulation of various metabolism pathways, the metabolic regulatory networks for different classes of metabolites within each cluster were constructed (Fig. 2f). The result showed that 4764 genes were correlated with at least two classes of metabolites, and among them, 2390 genes were involved in at least two clusters (Additional file 3: Table S7). For instance, out of the 10 ABA biosynthetic genes detected, 8 of them were also correlated with sucrose and glucose (Fig. 2f; Additional file 3: Table S7). This result is consistent with the previous report that ABA can increase soluble sugar content [46]. In addition, by mapping the genes in the database onto the chromosomes of the cv. Dangshansuli genome, visualization of this database provided evidence that the accumulation of each metabolite class is governed by the cooperatively regulation of genes located on multiple chromosomes (Additional file 1: Fig. S6).

Effect of DNA methylation on gene expression

To investigate the dynamics of DNA methylation in the fruit flesh, genome-wide bisulfite sequencing was performed on flesh samples at each stage and generated 866.75 million raw reads across all samples. Approximately 97.80% of these raw reads were used to examine DNA methylation sites in the whole genome, the number of clean bases obtained for each stage was more than 30 times larger than the size of the pear genome (Additional file 3: Table S9). As a result, a total of 39,016,459 sites were cytosine-methylated in the fruit flesh (Additional file 2: Table S1). On average, each stage displayed cytosine methylation at 27,345,777 sites, with 18,634,952 sites being consistently methylated across all stages. During the fruit flesh development, the DNA methylation level gradually increased in the whole genome and at the CHH region, while it gradually decreased at the CG region, and hardly changed at the CHG region (Additional file 3: Table S9). These trends were also observed in the gene body and its flanking regions, with the exception of a gradual decrease in cytosine methylation level within the gene body (Fig. 3a). Moreover, the CG methylation level gradually decreased in the two flanking regions of transposable elements (TE), while the CHH methylation level gradually increased in the TE body as well as its flanking regions (Fig. 3a). These results suggest that the increase in DNA methylation mainly is primarily driven by an elevation in CHH methylation. This finding is consistent with a previous report on sweet oranges [26].

Fig. 3

The involvement of DNA methylation in gene transcript during pear flesh development. a DNA methylation levels of genes and TEs in the fruit flesh at 11 stages (from S1 to S11). b Analysis of the differentially expressed genes (DEGs) correlated with DNA methylation in the promoter. The percentage represents the proportion of the DEGs that were correlated or uncorrelated with DNA methylation in the promoter. c Flesh callus was treated by the DNA methylation inhibitor 5-azacytidine (5’-Aza). Control is the untreated callus. d Genome-wide bisulfite sequencing of 5’-Aza-treated and control calli revealed the changes in DNA methylation. Standard error bars were calculated based on three replicates. Analysis of variance were calculated by Student’s t-test. Single asterisk stand for the level of significance at P-value < 0.05. e Chromosome locations of DNA methylation sites and all predicted genes revealed the genes modified by DNA methylation of all three context (CG, CHG, and CHH). Chr 1ꟷ17 represent the 17 chromosomes in pear genome. f A diagram showed the ratio of upregulated and downregulated genes between S1 and later stages. These genes have been identified to be repressed by 5’-Aza treatment in fruit callus and to be modified by DNA methylation in developing fruits. g A venn diagram showed the number of the DEGs modified by CG, CHG, and/or CHH methylations

To survey the effect of DNA methylation on gene expression, we investigated and found that cytosine methylation was present in the promoters of 38,159 genes including 20,322 DEGs (Additional file 2: Fig. S7a). Among these DEGs, the expression levels of 6209 DEGs were correlated with the average levels of cytosine methylation in their corresponding promoters, either positively or negatively (false discovery rate < 0.05; Fig. 3b, Additional file 3: Table S10). This suggests that DNA methylation may have both positive and negative effect on gene expression in pear fruit flesh. To verify this, pear flesh callus was treated with 50 mM 5-azacytidine (5’-AZa), a DNA methylation inhibitor (Fig. 3c). Genome-wide bisulfite sequencing of the 5’-Aza-treated and control (untreated) calli generated a total of 634.06 million raw reads (Additional file 4: Table S11) and showed that cytosine methylation level was reduced upon treatment with 5’-Aza (p-value = 0.0071 < 0.05; Fig. 3d). Subsequently, RNA sequencing of the 5’-Aza-treated and control calli generated a total of 267.34 million raw reads (Additional file 4: Table S12) and identified 1426 DEGs (Additional file 4: Table S13). A total of 58,943 differential methylation regions (DMRs) were detected from the promoters of 1328 DEGs, which included 705 DEGs present in the gene-metabolite database (Additional file 4: Table S13). These results suggest that the expression levels of these 1328 genes result from the change in DNA methylation in pear flesh callus.

To identify the methylation regions within the promoters of DEGs in fruit flesh, the methylation levels in a region enriched with cytosine methylation sites were analyzed between stages. A total of 178,499 DMRs were identified within the promoters of 35,176 genes in the fruit flesh (Fig. 3e, Additional file 2: Table S2). Among these DMRs, 92,912 were located within the promoters of 18,672 DEGs (Additional file 4: Table S14), including 12,396 DEGs that were correlated with 437 DAMs (Additional file 1: Fig. S7b, Additional file 4: Table S15). Out of the 12,396 DEGs, 647 were differentially expressed between the 5’-Aza-treated and control calli (Additional file 4: Table S13). After excluding the genes that were not differentially expressed between S1 and other stages, we found that 56.67 to 78.17% of genes, which were repressed by 5’-Aza treatment in fruit callus (Additional file 4: Table S13), were higher expressed in the developing fruits at later stages (from S3 to S11) compared to S1 (Fig. 3f). These results indicate that the expression of these 18,672 DEGs is likely modified by DNA methylation in the fruit flesh. Notably, 95.58% of these DEGs were modified by CHH DMR (Fig. 3g), suggesting that CHH DMR may play a more important role in modifying genome-wide gene expression in the fruit flesh compared to CG and CHG DMRs.

Increased DNA methylation during flesh development is associated with decreased expression of the genes involved in DNA demethylation

To elucidate the underlying mechanism of increased DNA methylation during pear flesh development, we initially employed phylogenetic analyses to identify genes involved in the RdDM pathway (Additional file 1: Fig. S8a-d). Among these genes, PbDCL2.2 was not expressed in fruit flesh (Additional file 4: Table S16). Overall, the expression levels of PbAGO4.1, PbAGO4.2, PbAGO4.3, PbAGO6, PbDCL2.1, PbDCL3.1, PbDCL3.2, PbDCL3.3, PbDCL4, PbNRPD1, PbNRPE1, PbRDR2.1, PbRDR2.2, and PbRDR2.3 gradually decreased during flesh development (Additional file 1: Fig. S8e), and of which, most were negatively correlated with C and CHH methylation levels (Additional file 1: Fig. S8f). This result indicates that the decreased expression of these genes is not consistent with the increased DNA methylation during pear flesh development. In contrast, the expression levels of PbDCL2.3, PbDCL2.4, PbDCL2.5, and PbDCL2.6 gradually increased during flesh development (Additional file 1: Fig. S8e), and notably, the expression level of PbDCL2.3 was positively correlated with C and CHH methylation levels (Additional file 1: Fig. S8f). Theoretically, the increased expression of these DCL genes raises an abundance of RdDM-dependent 24-nt siRNA. To test this hypothesis, small RNAome sequencing was performed on fruit flesh at 11 stages. A total of 620.56 million raw reads were generated from the small RNA-Seq dataset, with 95.63% (clean reads) used for the analysis of 24-nt siRNA clusters (Additional file 4: Table S17). In total, 2,279,387 siRNAs were detected in the fruit flesh with an average of approximately 173,425 siRNAs detected at each stage, and of the 139,462 siRNAs were consistently present across all stages (Additional file 2: Table S1). As a whole, the levels of 24-nt siRNAs in the whole genome, as well as in CHH DMRs, CG DMRs, and CHG DMRs, were reduced in fruit flesh from S5 to S11 compared to those from S1 to S4 (Additional file 1: Fig. S9). This result shows that the expression of the four DCL genes does not result in an abundance of 24-nt siRNAs in the fruit flesh. Therefore, the changes in gene expression observed in RdDM pathway are not associated with the increase in DNA methylation during pear flesh development.

The maintenance of plant DNA methylation is governed by DNA methytransferases [27]. To examine the potential association between DNA methytransferases and the increased DNA methylation during pear flesh development, we isolated the 12 pear orthologs of Arabidopsis DNA methyltransferase genes from the pear genome (Additional file 1: Fig. S10a). The expression levels of the pear DNA methyltransferase genes, PbDRM1.1, PbDRM1.2, PbDRM3, PbCMT2.1, PbCMT2.2, PbCMT2.3, PbCMT2.4, PbCMT3.1, and PbCMT3.2, were higher in the fruit flesh at S1 and S2 compared to other stages (Additional file 1: Fig. S10b), and of which, most were negatively correlated with C and CHH methylation levels (Additional file 1: Fig. S10c). Moreover, the expression level of PbMET1.2 was higher in the fruit flesh at S3 and S5 compared to S10 and S11 (Additional file 1: Fig. S10b), and was not correlated with C and CHH methylation levels (Additional file 1: Fig. S10c). These results indicate that the differential expression of these genes does not consistent with the increased DNA methylation during pear flesh development. In contrast, the expression level of PbMET1.1 gradually increased during flesh development (Additional file 1: Fig. S10b) and was positively correlated with C and CHH methylation levels (Additional file 1: Fig. S10c). Notably, MET1.1 is responsible for maintaining CG cytosine methylation [27]. However, CG methylation level gradually decreased during pear flesh development (Fig. 4a) and was negatively correlated with the expression level of PbMET1.1 (false discovery rate = 0.002 < 0.05). Therefore, the changes in the expression of the DNA methyltransferase genes are also not associated with the increased DNA methylation during pear flesh development.

Fig. 4

The decreased expression of demethylase genes is consistent with the increased DNA methylation during pear flesh development. a Expression patterns of demethylase genes in the fruit flesh at all 11 stages (S1ꟷS11). Z-score standardized values of each DEG across all 11 stages were used for clustering analysis. The color bar indicates the increasing expression levels of gene from blue to red. b Correlation analysis of demethylase genes with C (left panel) and CHH (right panel) methylations. The dotted line with red color represents the false discovery rate at 0.05. c A diagram showing the establishment, maintenance, and elimination of DNA methylation in fruit flesh

The regulation of plant methylome involves the dynamic interplay between DNA methylation and DNA demethylation processes [27]. To determine the expression patterns of DNA demethylase genes in pear fruit flesh, we identified the eight pear orthologs of the Arabidopsis ROS1 gene from the pear genome (Additional file 1: Fig. S11). Among these genes, PbDME1.1, PbDME1.2, and PbDME1.3 were not expressed in the fruit flesh (Additional file 4: Table S16). Overall, the expression levels of PbDME1.4, PbDME1.5, and PbROS1.3 gradually decreased during flesh development (Fig. 4a) and were negatively correlated with C and CHH methylation levels (Fig. 4b). Additionally, the expression levels of PbROS1.1 and PbROS1.2 gradually decreased in the fruit flesh from S5 to S11 (Fig. 4a). These results suggest that the decreased expression of these five genes is consistent with the increased DNA methylation during pear flesh development (Fig. 4c). Moreover, considering the role of Arabidopsis IDM1 in catalyzing histone H3K18 acetylation to create a chromatin environment that facilitates ROS1 function [35], we identified the three pear orthologs of Arabidopsis IDM1 gene from the pear genome (Additional file 1: Fig. S11). The expression levels of PbIDM1.1, PbIDM1.2, and PbIDM1.3 also gradually decreased during flesh development (Fig. 4a) and were negatively correlated with C and CHH methylation levels (Fig. 4b). Consequently, the decreased expression of all three IDM genes is consistent with the increased DNA methylation during flesh development. In conclusion, these findings suggest that the increased DNA methylation during pear flesh development may be contributed by the decreased expression of genes involved in DNA demethylation.

DNA methylation is involved in pear fruit metabolism

It is determined that 12,396 DMR-modified DEGs correlated with 437 DAMs (Additional file 4: Table S15). To clarify the relationship between these DAMs and DNA methylation in the promoter of the DEGs, a correlation analysis was performed, revealing that DNA methylation in the promoters of 3987 DEGs was correlated (False discovery rate < 0.05) with 316 DAMs present in all eight clusters (Additional file 4: Table S18). Moreover, it was observed that DNA methylation in the whole genome was negatively correlated with DAMs in clusters I, II, and IV, but was positively correlated with DAMs in cluster V (Additional file 4: Table S19). These results indicate that DNA methylation likely plays a role in regulating fruit flesh metabolites by impacting expression of genes involved in metabolite production.

Furthermore, a majority of the DMR-modified DEGs associated with the synthesis and metabolic processes of ABA, sucrose, and glucose were identified (Additional file 4: Table S15), and the DNA methylation in the promoters of partial DEGs, even in the whole genome, correlated with ABA, sucrose, and glucose (Additional file 4: Tables S18 and S19). ABA is an important phytohormone that promotes pear fruit ripening by increasing soluble sugar (including sucrose and glucose) content and ethylene release [46]. Therefore, DNA methylation may be involved in pear flesh development by influencing ABA biosynthesis. To test this hypothesis, pear fruit flesh was treated with 5’-Aza beneath the pericarp at S9 (Fig. 5a). Genome-wide bisulfite sequencing of the 5’-Aza-treated and control (H2O) pericarp generated a total of 674.03 raw reads (Additional file 4: Table S11) and showed a reduction in cytosine methylation level due to the treatment (p-value = 0.038 < 0.05; Fig. 5b). Compared to the control, the 5’-Aza treatment promoted chlorophyll accumulation in the pericarp (Fig. 5c), leading to a green pericarp (Fig. 5a). Moreover, genome-wide bisulfite sequencing of the 5’-Aza-treated and control (H2O) flesh generated a total of 679.69 raw reads (Additional file 4: Table S11) and showed a reduction in cytosine methylation level caused by the treatment (p-value = 0.012 < 0.05; Fig. 5b). Compared to the control, the 5’-Aza treatment promoted ABA, β-carotene, and xanthophyll biosynthesis in the flesh (Fig. 5c). This result suggests that DNA methylation is involved in carotenoid metabolism and perhaps accelerates fruit ripening by increasing ABA production.

Fig. 5

The involvement of DNA methylation in fruit flesh metabolism. a Phenotype change was detected in pear fruit treated with 5’-Aza compared to H2O (Control). White arrow indicates the dark green pericarp. b Genome-wide bisulfite sequencing of 5’-Aza-treated and control fruit revealed the changes in DNA methylation in pericarp and flesh of pear fruits. c Measurement of chlorophyll, ABA, β-carotene, and xanthophyll in the flesh and pericarp of the pear fruits treated with 5’-Aza and H2O (Control). d Expression levels of ABA biosynthetic genes and the positive correlated TFs were determined in the fruit flesh treated with 5’-Aza and H2O (Control). The standard error bars were calculated based on three replicates. Analysis of variance was calculated by Student’s t test. Single and double asterisks stand for the level of significance at P-value < 0.05 and < 0.01, respectively

To clarify the regulation of DNA methylation on carotenoid metabolism, quantitative real-time polymerase chain reaction (qRT-PCR) analysis was performed on the fruit flesh treated with 5’-Aza and H2O. Among the DMR-modified DEGs, 10 pathway genes and 35 TFs in expression level were positively correlated with ABA content (Additional file 4: Table S15). Out of these 10 pathway genes, PbPSY, PbZDS1, PbHYB1, PbSDR2, and PbAAO were expressed at higher levels in pear flesh from S5 to S11 compared to those from S1 to S4 (Additional file 1: Fig. S12). The expression levels of PbPSY, PbHYB1, and PbAAO were upregulated in the fruit flesh with the 5’-Aza treatment compared to the control, while no significant difference was observed for PbZDS1 and PbSDR2 expressions (Fig. 5d). These results suggest that DNA methylation represses the expressions of PbPSY, PbHYB1, and PbAAO to limit ABA production in fruit flesh. Of the 35 TFs, 18 were relatively highly expressed in the fruit flesh from S5 to S11 compared to other stages within their respective families (Additional file 1: Fig. S12). The expression levels of PbZFP1, PbHB3, PbERF, PbB3.1, PbZFP4, PbbZIP1, PbGRAS3, PbbHLH3, and PbNAC1 were upregulated in the fruit flesh with the 5’-Aza treatment compared to the control, while no significant difference was observed for PbZFP2, PbHB1, PbHB2, PbIAA2, PbIAA3, PbbHLH2, and PbbHLH1 expressions (Fig. 5d). These results suggest that the decreased DNA methylation may promote the expressions of several TFs to induce the expression of PbPSY, PbHYB1, and PbAAO. In addition, PbIAA1 and PbMYB1 were downregulated in the fruit flesh with the 5’-Aza treatment compared to the control (Fig. 5d), suggesting that the decreased DNA methylation can also inhibit gene expression.

Verification of a novel transcription factor with roles in regulating ABA biosynthesis

To clarify the potential involvement of the tested TFs in ABA biosynthesis, a dual-luciferase assay was conducted using the promoters of ABA biosynthetic genes. A total of 12 TFs were selected based on expression and DNA methylation analyses (Fig. 5d, Additional file 4: Table S15). The result showed that the LUC gene, under the control of the PbPSY, PbZDS1, PbHYB1, PbSDR2, and PbAAO promoters, was significantly activated by at least three of the twelve TFs tested (Fig. 6a). This suggests that these 12 TFs are likely to activate a subset of genes involved in ABA biosynthesis.

Fig. 6