{'fraction_non_alphanumeric': 0.03893949617633828, 'fraction_numerical': 0.029183535762483132, 'mean_word_length': 4.977983193277311, 'pattern_counts': {'":': 0, '<': 13, '<?xml version=': 0, '>': 2, 'https://': 1, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Integrative analysis of transcriptome and metabolome reveals the effect of DNA methylation of chalcone isomerase gene in promoter region on Lithocarpus polystachyus Rehd flavonoids

Metabolite biosynthesis is regulated by gene expression, which is altered by DNA methylation in the promoter region. Chalcone isomerase (CHI) gene encodes a key enzyme in the Lithocarpus polystachyus Rehd flavonoid pathway, and the expression of L. polystachyus CHI (LpCHI) is closely related to the synthesis of flavonoid metabolites. In this study, we analyzed the DNA methylation site of the LpCHI promoter and its effect on gene expression and metabolite accumulation. The proportions of three types of LpCHI promoter DNA methylation are 7.5%, 68.75%, 18.75%, determined by bisulfite sequencing. Transcriptome sequencing shows that LpCHI is strongly up-regulated in LpCHI promoter methylation Type A but down-regulated in LpCHI promoter methylation Type B and Type C. The expression of LpCHI shows no significant difference between Type B and Type C. Moreover, nine kinds of differentially expressed transcription factors (DETFs) bind to seven CpG-sites of the LpCHI promoter region to regulate LpCHI expression. The results of metabolomics show that differentially accumulated flavonoids are higher in LpCHI promoter methylation Type A than in LpCHI promoter methylation Type B and Type C. Additionally, a positive correlation was found between the LpCHI expression and flavonoids accumulation. These results show that the effect of CpG site-specificity on gene transcription is great than that of overall promoter DNA methylation on gene transcription. The mechanisms of flavonoid genes regulating metabolite accumulation are further revealed.


Introduction
As a traditional Chinese herbal medicine, Lithocarpus polystachyus Rehd belongs to Fagaceae and is evergreen. Due to its young leaves can be used as tea leaves, L. polystachyus is called "sweet tea" in folk [1,2]. In addition, L. polystachyus is rich in flavonoid metabolites, and trilobatin and phlorizin of dihydrochalcone are its main active components [3,4]. Leaf aqueous extract of L. polystachyus can inhibit proliferation, migration, and invasion of cancer cells and has anticancer activity [5]. The flavonoid fraction of L. polystachyus could reduce fasting blood glucose, improve serum lipid levels, and exerts therapeutic and preventive effects on hyperlipidemia and hyperglycemia in diabetes rats [6]. Flavonoids are common secondary metabolites in plants and are synthesized by phenylalanine derivatives through the phenylpropane synthesis pathway [7,8]. Chalcone synthase (CHS), the first rate-limiting enzyme in the biosynthesis pathway of flavonoids, can catalyze three malonyl CoA and one coumarin CoA to produce naringin chalcone (6 ′ -hydroxychalcone). Chalcone isomerase (CHI) is the second key enzyme in the biosynthesis pathway of flavonoids. Naringin chalcone (6 ′ -hydroxychalcone) isomerizes to 5,7,4 ′ -trihydroxyflavanone (naringenin) via CHI. The CHI superfamily is divided into four families: Type I CHI widely exists in vascular plants and has activity; Type II CHI with activity is mainly present in legumes; Type III CHI (FAP) has fatty acid catalytic activity; Type Ⅳ CHI (CHIL) has no catalytic activity [9,10]. Since the cDNA of CHI was first cloned from pea in 1987, scholars have successively analyzed the CHI gene sequence of Mirabilis himalaica [11], Peer review under responsibility of KeAi Communications Co., Ltd.
DNA methylation is an important epigenetic mechanism. DNA methylation of plant genes is involved in X-chromosome inactivation, transposable element silencing, genomic imprinting, and genome stability. Additionally, it regulates embryogenesis, flowering time, and leaf morphology [18]. DNA methylation within the promoter region of plant functional genes prevents transcription factors (TFs) from binding to methylated promoters, thereby inhibiting gene expression [19]. The promoter methylation level of the osmotic stress Wdreb2 TF is negatively correlated with its gene expression level in wheat [20]. The different methylation sites and levels of anthocyanidin synthase (ANS) gene promoter in red and white cultivars may result in different expression levels of ANS, leading to different flower colors of two lout cultivars [21]. At the developmental stages of Malus halliana, different methylation levels of R2 and R8 in the MtMYB10 promoter could cause differentially expressed MtMYB10, thus leading to different anthocyanin contents in flowers [22]. The evaluation of the interconnection of the DNA methylome and transcriptome shows melatonin effects on promoter methylation of flavonoid biosynthesis genes and gene expression and controls flavonoid accumulation [23]. Additionally, whole-genome DNA methylation patterns of apple fruit skin show that hypomethylation of m CG and m CHG contexts in flavonoid biosynthesis pathway genes leads to transcriptional activation and promotes anthocyanin accumulation [24]. There have been reports about the effect of DNA methylation of L. polystachyus CHI (LpCHI) gene promoter on its gene expression. However, the molecular mechanism of LpCHI promoter DNA methylation regulation has not been thoroughly studied. In this study, we analyzed the DNA methylation in the LpCHI promoter region by bisulfite sequencing. Further, we revealed the influence of the CHI gene on flavonoid synthesis based on transcriptomics and metabolomics of L. polystachyus. Technical and theoretical supports for promoting the industrial synthesis of L. polystachyus flavonoids can be provided by artificially regulating the DNA methylation state of the key enzyme gene promoter in engineered bacteria.

Bisulfite sequencing and analysis of DNA methylation of the LpCHI promoter
According to Lin et al. [24], total DNA was extracted from the mature leaves of L. polystachyus. Afterward, the PCR amplification was performed using the methylation-specific kit according to the manufacturer's instructions. A 25 μL reaction mixture containing 1 μL of methylation DNA, 1 μL of CHIjhS1 (CHIjhS3), 1 μL of CHIjhX1 (CHIjhX3), 0.4 μL of MSP DNA polymerase, 2 μL of 10 × Msp PCR buffer, and 13 μL of ddH2O was prepared. The following cycling conditions were used for amplification: 94 • C for 5 min, (94 • C for 20 s, 55 • C for 30 s, 72 • C for 20 s (10 s)) × 35 cycles, 72 • C for 5 min. The sequences of the primers are shown in Supplementary Table S6. PCR products of bisulfite-treated DNA sequences were analyzed by ethidium bromide-stained 1% agarose gel, and the corresponding band was sequenced at the Nuosai Genome Research Center (Beijing, China). Sequencing results of bisulfite-treated DNA sequences were analyzed using BiQ Analyzer 2.0, and DNA sequencing [17] was performed to screen DNA methylation sites.

RNA-seq analysis
Total RNA was extracted from frozen mature L. polystachyus leaves, and the mRNA libraries for each sample were constructed and sequenced in the Illumina Nova Seq 6000 system. After filtering and quality control of the original data obtained by sequencing, high-quality clean reads were obtained. Sequences from each library were spliced using the Trinity (http://trinityrnaseq.github.io/) software. In addition, all isoforms with dad FPKM (Fragments Per Kilobase of transcript per Million mapped reads) < 1, TPM (Transcripts Per Million) < 1, or IsoPct (Isoform PerCen Tage) < 5.0% were removed. The TGICL program (version 2.0.6) (https://sourceforge.net/projects/tgicl/) was used to cluster transcripts into unigenes. The BLAST algorithm was used to compare the unigene sequences with sequences from the NCBI nonredundant (NR), Swiss-Prot, Gene Ontology (GO), Clusters of Orthologous Groups (COG) of proteins, and Kyoto Encyclopedia of Genes and Genomes (KEGG) databases. Moreover, KOBAS v2.0 (http://kobas.cbi.pku.edu.cn) was applied to identify enriched KEGG pathways associated with unigenes. After predicting the amino acid sequences encoded by the unigenes, the sequences were compared with those in the Pfam database using the HMMER (http://www.hmmer.org/) program for functional annotation of the unigenes.

Sample extraction and UPLC-MS/MS analysis
The freeze-dried mature leaves were pulverized using a mixer mill (MM 400, Retsch, Haan, German) at 30 Hz for 1.5 min with zirconia beads. A total of 100 mg powder was weighed and extracted with 1.0 mL 70% aqueous methanol at 4 • C overnight. After centrifugation at 10,000 g for 10 min, the supernatant was absorbed and filtered by 0.22 μm Pro size microporous membrane for UPLC-MS/MS analysis.
Based on MWDB (Metware Database), a qualitative analysis was conducted according to secondary spectrum information. Metabolites were quantitatively analyzed by multiple reaction monitoring (MRM) with QQQ mass spectrometry. The above analysis was completed by Metware Biotechnology Co., Ltd. (Wuhan, China).

RT-qPCR (real-time quantitative PCR)
Total RNA was extracted and purified from L. polystachyus leaves of three types of DNA methylation at different stages using the RNAprep Pure Plant kit (Tiangen, Beijing, China). Afterward, the purity, concentration, and integrity of the RNA were determined with the Nano-Drop spectrophotometer (Thermo Fisher, Waltham, MA, USA) and a 2100 Bioanalyzer (Agilent, Santa Clara, CA, USA). According to the manufacturer's instructions, 11 μL total RNA was reverse transcribed using the RevertAid Frist Strand cDNA synthesis kit (Beijing Bayerdi Biotechnology). Primer 5.0 (Premier Biosoft Interpairs, Palo Alto, CA, USA) was used to design primers specific for the selected genes (Supplementary Tables S6) and a dissolution curve analysis was used to determine specificity. The L. polystachyus GAPDH (glyceraldehyde-3phosphate dehydrogenase) gene was used as the internal reference gene. The Applied Biosystems 7900HT Fast Real-Time PCR System (Life Technologies, Waltham, MA, USA) with a 96-well plate was used for amplification reactions consisting of 3 min at 95 • C, followed by 40 cycles of 5 s at 95 • C, and 15 s at 60 • C in a volume of 10 μL. The qRT-PCR assay was performed in three biological replicates, and relative expression levels were calculated based on the 2 − ΔΔCt method.

Metabolite extraction and determination
The leaves were oven-dried to a constant weight and ground. A total of 10 mL extracting solution (60% ethanol) was added to 0.100 g of leaf powder for ultrasonic extraction. After adding 1 mL of extracting solution (70% ethanol), mixture A was extracted by ultrasonic method (ultrasonic power: 300 W, crushing time: 5 s, interval time: 8 s) at 60 • C for 30 min. The supernatant was collected after centrifugation at 12,000 rpm for 10 min. The supernatant was then processed according to the instructions of the kit for the determination of total flavonoids in plants. Total flavonoid content was calculated using the standard rutin curve.
A1-A3, B1-B3, and C1-C3 denote the three biological replicates of Type A, Type B, and Type C, respectively. Q20: percentage of bases with a Phred quality score greater than 20. Q30: percentage of bases with a Phred quality score greater than 30.

Transcriptome sequencing and unigene annotation
We collected samples of LpCHI promoter methylation Type A, LpCHI promoter methylation Type B, and LpCHI promoter methylation Type C for RNA-seq analysis. On average, 25.0 million clean reads were obtained for each sample. The GC content of sequences for the nine samples is about 44%, and more than 92% of the bases have a Phred quality score greater than 30 (Table 1). Transcripts were spliced by Trinity, and 254,492 transcripts and 216,116 unigenes were obtained. The total length of unigenes is 107,119,884 bp, with the longest being 15,721 bp and the shortest being 200 bp, and the average length is 497 bp. N50 is 641 bp, and N90 is 240 bp.

Expression of flavonoid synthesis-related genes
Based on significance level (FDR ≤0.01) and log2 fold change (FC), the screening conditions for differentially expressed genes (DEGs) serve as the threshold for up-regulation (FC≥ + 0.5) and down-regulation (FC ≤ -0.5). A total of 20,593 DEGs were obtained in the transcriptome sequencing (Fig. 3A). There are 14,121 DEGs from Type A_vs_Type B, including 6852 up-regulated genes and 7269 down-regulated genes.
Additionally, 14,137 DEGs, including 7547 up-regulated genes and 6590 down-regulated genes, were identified between LpCHI promoter methylation Type A and Type C. We identified 10,301 DEGs, including 7547 up-regulated and 5757 down-regulated genes in LpCHI promoter methylation Type B compared with LpCHI promoter methylation Type C. The Venn diagram shows the overlap between DEGs in the three comparisons (Type A_vs_Type B, Type A_vs_Type C, and Type B_vs Type_C). A total of 1745 DEGs are different in all three comparisons, with 1891 DEGs only present in Type A_vs_Type B and 905 DEGs only present in Type B_vs_Type C (Fig. 3A). In addition, there are 8570 DEGs between the group of Type A_vs_Type B and Type A_vs_Type C, suggesting that the transcript expressions of LpCHI promoter methylation Type A, Type B, and Type C differ significantly.
In plants, phenylalanine can be converted by phenylalanine ammonia-lyase (PAL) to produce cinnamic acid, catalyzed by cinnamate 4-hydroxylase (C4H) and 4-coumarate CoA ligase (4CL) to produce ρ-coumaroyl CoA, and then enters the flavonoid biosynthesis pathway [27]. According to the KO numbers of the three genes and the flavonoid biosynthesis pathway (ko00941), 53 DEGs were identified in all comparisons (Fig. 3B, Supplementary Table 1). Specifically, TRINI-TY_DN22537_c0_g2 and TRINITY_DN256813_c0_g1 annotate to CHI in KEGG, Pfam, and Nr; TRINITY_DN22537_c0_g2 is identical to the LpCHI (MZ733721) [28] sequence cloned in the previous study; TRINI-TY_DN256813_c0_g1 is the gene encoding Type Ⅳ CHI (CHIL). Their expressions are high in LpCHI promoter methylation Type A but are not significantly different in LpCHI promoter methylation Type B and Type C (Fig. 3B). Two unigenes encoding LpPAL exhibit high expression in Type A but do not differ significantly in Type B and Type C. In addition, LpC4H1 transcription is not significantly different in Type A and Type B, while LpC4H2 is up-regulated in Type A but not significantly different in Type B and Type C. Four transcriptions encoding 4CL show the same expression trend as LpPAL, except that Lp4CL4 has high expression in Type A but low expression in Type B (Fig. 3B). CHS catalyzes three malonyl CoA and one coumaroyl CoA to naringenin. The two LpCHS identified from transcriptome sequencing data are differentially expressed in the three types of LpCHI promoter methylation and show the same aggregation trend as LpCHI (Fig. 3B). Downstream genes of LpCHI in the flavonoid biosynthetic pathway, including LpDFR1 (dihydroflavonol 4-reductase), LpF3 ′ H (flavonoid 3 ′ -monooxygenase), LpANS1 (anthocyanidin synthase), LpANS2, LpANR2 (anthocyanidin reductase), and three LpLARs (leucoanthocyanidin reductase), are down-regulated in Type B and Type C. However, LpFLSs (flavonol synthase) are principally expressed in Type B (Fig. 3B).
In order to validate the reliability of the transcriptome data, realtime PCR (qRT-PCR) was conducted to DEGs of the flavonoid biosynthetic pathway. The results show that the relative expression trend of DGEs is consistent with transcriptome sequencing data, proving the reliability and accuracy of transcriptome sequencing data (Fig. 4).
The transcript expression of DEGs in the flavonoid biosynthetic pathway in young and mature leaves was assessed by qRT-PCR, and the changing transcript expression at different developmental stages was analyzed. Expressions of LpFLSs involved in the flavone and flavonol biosynthesis pathways and LpANS2 involved in the anthocyanin biosynthesis pathway are greater at the young stage compared with the mature stage. Expressions of LpPALs, LpCHSs, LpCHI, LpCHIL are lower at the mature stage than at the young stage (Fig. 4).

Analysis of TFs associated with LpCHI
Transcription factor binding sites in the LpCHI promoter region were predicted using Transcriptional Regulatory Map with Castanea mollissima as a reference. The results suggest that 33 kinds, 148 TFs of C. mollissima can bind to the LpCHI promoter region with 284 binding sites. There are nine kinds TFs that can bind DNA methylation sites in LpCHI promoter CpG-island. Each TF can bind 1 to 3 CpG-sites (Fig. 5,  Supplementary Table 2). Additionally, TF binding sites are absent in -1455, -480, -1483, and -1761 sites, and the remaining 12 CpG-sites have TF binding sites (Fig. 5). The -1469 and -380 CpG-sites occur in LpCHI promoter DNA methylation Type A and Type B; the -1469 site has C2H2, ERF, NAC, and MYB TF binding sites; the − 1380 site has ERF binding site. CpG-sites are only present in Type C, where LBD, TCP, HSF, GATA, bZIP, C2H2, and ERF binding sites are found.
Multiples sequences of C. mollissima 148 TF amino acid sequences and L. polystachyus amino acid sequences in the transcriptome database were aligned by BLAST, and 427 DETFs were identified. Correlation analysis was performed between 427 DETFs and LpCHI expressions by the R package. PCC (Pearson correlation coefficient) > 0.9 and P < 0.01 were set as the cutoff criterion. A total of 148 DETFs were filtered, which strongly correlated with LpCHI expression, including bHLH, bZIP, ERF, MYB, and NAC. Among them, 85 DETFs were positively correlated, and 63 DETFs were negatively correlated ( Fig. 6A; Supplementary Table 3). The protein-protein interaction network of 148 DETFs and LpCHI was reconstructed with STRING 11, and 91 proteins had protein-protein    A. thaliana [29] (Supplementary Fig. 1A). LqTCPs are grouped with CmTCP1, LqGATA1 displays closer relationships with CmGATA1, and LqGATA4 is closely related to CmGATA2. Five LpHSFs expressions are negatively correlated with LpCHI. Only LpHSFB clusters with CmHSF1.
LpC2H2s, LpLBDs, and predicted TFs of C. mollissima are not clustered ( Supplementary Figs. 1B-F). Eighteen LpMYBs were classified into eight groups. Among them, LpMYB4, LpMYB11, CmMBY1, and CmMBY2 are clustered in Group II; LpMYB8, LpMYB6, LpMYB12, and CmMBY3 are in Group Ⅷ. However, LpMYB1 and LpMYB10 in Group III do not cluster with A. thaliana TFs (Supplementary Fig. 2). Sixteen LpERFs and ten CmERFs are classified into nine branches. Seven branches cluster with LpERFs and CmERFs apart from Group I and Group II. LpERF15 shows 97.07% sequence similarity to CmERF9 in Group Ⅶ. The sequence similarity between LpERF6 and CmERF8 in group VIII is 90.64%. LpERF15 and LpERF6 positively correlate with LpCHI expression (Supplementary Fig. 3). Three in the identified six branches of LpNACs have negative relations to LpCHI and are clustered with CmNAC. LpNAC13 exhibits relatively high sequence similarity with CmNAC at 95.95%, and it clusters with ANAC042, a TF involved in the regulation of camalexin biosynthesis [30] (Supplementary Fig. 2).

Identified metabolites involved in flavonoid biosynthesis
To further evaluate the effects of different LpCHI promoter DNA methylation on flavonoid biosynthesis, metabolites in L. polystachyus leaves of three types of LpCHI promoter DNA methylation were analyzed using UPLC-MS/MS. A total of 441 metabolites containing 11 classes were detected and quantified. Among them, 107 (24.26%) are flavonoids, followed by 75 (17.01%) phenolic acids (Fig. 7A).
The results of principal component analysis (PCA) show that metabolites from different types of LpCHI promoter DNA methylation are clearly separated in the score plots, where the first principal component (PC1) is plotted against the second principal component (PC2). PC1 and PC2 represent 32.36% and 28.99% of the total variations, respectively (Fig. 7B). These results suggest that the samples between groups are scattered and clustered within groups, indicating good duplication within groups and significant differences between groups (Fig. 7B).
The screening criterion for differentially accumulated metabolites (DAMs) is VIP≥1, FC ≥ 2 or ≤0.5. All 184 kinds of DAMs are filtered out (Fig. 8A). In the comparison between Type A and Type B, 66 and 46 metabolites present up-regulation and down-regulation, respectively. Additionally, 105 kinds of DAMs are filtered out in Type A_vs_Type C, of     In addition, DAMs were mapped to the KEGG database and tested for enrichment. The results show that four pathways ("Flavonoid biosynthesis", "Anthocyanin biosynthesis", "Isoflavonoid biosynthesis", and "Flavone and flavonol biosynthesis") related to flavonoid metabolism are enriched in the three comparisons. The "Anthocyanin biosynthesis" pathway is significantly enriched in Type A_vs_Type B (Fig. 8B-D).
Cluster analysis for flavonoids of DAMs shows that 26 of 79 differentially accumulated flavonoids (DAFs) are significantly different among LpCHI promoter DNA methylation Typer A, Type B, and Type C; the accumulation of 8 flavonoids (DSK0207-DSK1158) is lower in Type A than in Type B and Type C; 8 flavonoids (DSK1116-DSK1417) have a higher accumulation in Type A than in Type B and Type C (Fig. 9A, Supplementary Table 4).

Correlation of metabolites and DEG of the flavonoid biosynthesis pathway
To investigate the effect of flavonoid synthesis gene expression on flavonoid content accumulation in L. polystachyus, we screened transcriptome data and qRT-PCR data for flavonoid synthesis genes consistent with LpCHI gene expression trends. The PCC (PCC>0.8, P < 0.01) between DEG and DAFs was calculated, and the correlation network diagram between genes and metabolites was constructed. Among

Flavonoids content, expression of key genes, and enzyme activity
The total flavonoid content is 47.5713 (±0.1653) g/100 g to 8.3803 (±0.1210) g/100 g for the leaves of L. polystachyus. The total flavonoid content of young leave is more than twice that of mature leaves. At different developing stages, the total flavonoids content of the three types of LpCHI promoter DNA methylation leave is the highest in Type B, followed by Type A and Type C (Fig. 11A). Naringenin is a CHI catalyzed product with 1.4568 (±0.0246) g/100 g for Type Ay and 0.2592 (±0.1273) g/100 g for Type Cy. The content of other samples is less than 0.1000 g/100 g. The naringenin content was significantly higher in the early stage than in the mature stage (Fig. 11B). In an integrative analysis of LpCHI expression and metabolome, kaempferol and apigenin are significantly correlated with LpCHI. These results are consistent with our metabolome results, indicating that Type A has higher kaempferol and apigenin content than Type B and C at the mature stage. The kaempferol content of Type A is significantly lower than that of Type B and Type C at early stages ( Fig. 11C and D). Trilobatin, phlorizin, and phloretin are the major chalcones and medicinal constituents in L. polystachyus. The content of trilobatin range from 0.0380 (±0.0015) g/100 g to 16.0065 (±0.1115) g/100 g, and the content of phloretin range from 0.0065 (±0.0023) g/100 g to 0.5345 (±0.0091) g/100 g. Both substances decrease during leaf maturation. Phlorizin content range from 1.4331 (±0.01113) g/100 g to 3.6023 (±0.0816) g/100 g and increase during leaf maturation (Fig. 11E-G).   Table 5). The relative expression of LpCHI is significantly correlated with the content of apigenin and kaempferol, which is consistent with our integrative analysis of transcriptome and metabolite (P < 0.01). The correlation between the relative expression of LpCHI and total flavonoid content is highly significant (P < 0.01). The relative expression of LpCHI is highly correlated with the content of trilobatin and phloretin, as well as the total content of flavonoid, trilobatin, and apigenin (P < 0.01). In addition, the PCC of LpCHIL is greater than that of LpCHI. The relative expression of Lp4CL1-3, LpCHS2, LpC4H2, and LpF3 ′ H are also significantly correlated with the total content of flavonoids, trilobatin, and apigenin (P < 0.01).
The activities of PAL, C4H, 4CL, CHS, and CHI enzymes involved in flavonoid biosynthesis are also determined. The results show little difference in the enzyme activities of PAL and CHI in different samples of L. polystachyus at different periods (P < 0.05). There is a large difference in enzyme activities for 4CL, followed by CHS and C4H (Fig. 11H to L). In addition, the results of SPSS correlation analysis suggested that PAL, C4H, 4CL, CHS, and CHI activities are not correlated with the total content of flavonoids, naringenin, kaempferol, apigenin, trilobatin, phlorizin, and phloretin (Supplementary Table 5).

Discussion
One of the regulatory mechanisms of gene expression is the binding of TFs to DNA. TFs recognize specific sequences in the promoter region and activate or suppress their expression [31]. MYB proto-oncogene protein (CMYB) binds to Sirtuin 4 (SIRT4) promoter region and suppresses gene transcriptional activity, and nuclear respiratory factor 1 (NRF1) binds to the SIRT4 promoter region and suppresses gene transcriptional repression [32]. DNA methylation of promoter CpG islands causes silencing of gene expression by blocking TF binding [31]. However, recent studies have shown that three models of DNA methylation affect TF binding: DNA methylation suppresses TF binding, DNA methylation promotes TF binding, and DNA methylation does not affect TF binding [33]. Methylation of 5-methylcytosine within the AP-2, MYC binding sites can directly suppress TF binding to nucleotide sequence [34]. DNA methylation of Kaiso and C/EBPα binding sites facilitates the binding to recognition sequences [35,36]. CpG methylation does not affect the binding of YY1 protein to the bi-directional promoters of the Surf-1 and Surf-2 genes [37]. The core promoter region of the CHI gene in Paeonia lactiflora Pall has 16 CpG sites, with methylation levels ranging from 59.8% to 65.7% at different growth stages. The DNA Fig. 11. Flavonoid content and enzymes activities of three types of LpCHI promoter DNA methylation at different growth stages. The "y" indicates young leaves, and "m" indicates mature leaves. Different lowercase letters represent significant differences (P < 0.05). methylation levels are negatively correlated with CHI expressions, and total flavonoid content in petals presents a similar trend to CHI expression at development stages, indicating that DNA methylation in the promoter region of CHI affects flavonoid content [38]. Previous studies have also shown that the methylation level of LpCHI gene promoter DNA is negatively correlated with the relative expression level of L. polystachyus [17]. In this study, the DNA methylation level of the LpCHI promoter is 37.5% for Type A, 68.75% for Type B, and 18.75% Type C. LpCHI expression is higher in Type A than in Type B, while there is no difference between Type B and Type C. Furthermore, no correlation is found between the DNA methylation level of LpCHI promoter and gene expression (P < 0.05). Therefore, we further analyzed the CpG-sites and found 3 different sites between the DNA methylation Type A and the Type C of the LpCHI gene promoter (-1469, -1480, -1380). The recognition sites for C2H2, MYB, and NAC TFs are at -1469 and ERF at -1469 and -1380, respectively (Fig. 4). We speculate that after DNA methylation, TFs are binding to the -1469 and -1380 sites that transcriptionally activate the expression of the LpCHI gene. Five different GpG-sites between the Type A and Type B, as well as -1797, -1783, -1735, -1723 sites, can bind to LBD, TCP, HSF, C2H2, GATA, ERF, and bZIP (Fig. 4). These 4 sites may contain TFs that bind DNA methylation to suppress and promote gene expression. The expression of the LpCHI gene is higher in Type A than in Type B. However, no significant difference is found between Type B and Type C, which can be reasonably explained when transcriptional inhibition is greater than transcriptional activation. The sensitivity of human TFs to DNA methylation by methyl-SELEX and bisulfite-SELEX show that most C2H2 and GATA family members preferentially recognize the binding sites of DNA methylation. Moreover, DNA methylation in the binding region decreases the binding of MYB, bZIP, and HSF families to the recognition sequences [39].
Secondary metabolites are important in the defense mechanisms of plants against abiotic and biotic stresses, improving flower colors and fruit quality and attracting pollinators and seed dispersers. They are defensive metabolites and toxins, which fight invasive pathogens and herbivores and protect plants from abiotic stresses, such as UV stress. In addition, secondary metabolites are a source of natural flavors and medicines [40]. TFs are DNA binding proteins with sequence specificity. It interacts with the regulatory regions of the target gene (usually promoters) to activate or suppress the expression patterns and levels of structural genes, regulating the synthesis of plant secondary metabolites. Several families of TFs do not directly bind to the target gene but interact with other cofactors to form complexes to regulate the expression of the target gene [40]. The MYB TF family is one of the largest families of transcription factors that regulate flavonoid biosynthesis in plants. The R2R3-MYB gene SbMYB8 of Scutellaria baicalensis can bind to the GmMYB92 BS3 sequence in the SbCHS promoter region, thereby regulating flavonoid biosynthesis [41]. TT2 (MYB), TT8 (bHLH), and TTG1 (WDR) form a stable ternary complex that can bind to the BANYULS promoter to directly regulate the expression of BANYULS, which in turn promotes the anthocyanin biosynthesis in A. thaliana [42]. Camalexin levels are much lower in mutants of ANAC042 plants than in the wild type, and NAC042 regulates the camalexin biosynthetic genes P450s then induce phytoalexin biosynthesis in A. thaliana [30]. TCP1 binds to the DWF4 promoter and regulates DWF4 expression, which modulates brassinosteroid biosynthesis in A. thaliana [43]. Three AP2/ERF TFs are identified as positive regulators of the CitCHIL1 expression. CitERF32 and CitERF33 activate the transcription by directly binding to the CitCHIL1 promoter, while CitRAV1 forms a transcriptional complex with CitERF33, strongly enhancing the activation efficiency and flavonoid accumulation in Ougan (Citrus reticulata cv. Suavissima) [44]. G-BOX BINDING FACTOR 1 (GBF1) of the bZIP family binds to CATALASE2 (CAT2) and inhibits CAT2 expression during senescence [45]. GBF1 and MYC2 (bHLH) physically interact in the nucleus, bind to HYH promoter, and mutually inhibit each other through non-DNA binding bHLH-bZIP heterodimers [46].
The CHI is one of the key enzymes required for flavonoid biosynthesis, regulating flavonoid content accumulation. Sun et al. reported that the expressing OjCHI of Ophiorrhiza japonica in Arabidopsis tt5 mutant restores the phenotype of anthocyanins and flavonols, indicating its function as a Type I CHI for flavonoid biosynthesis in vivo [47].
In vitro functional demonstrates that DcCHI1 of Dracaena cambodiana catalyzes the formation of naringenin from naringenin chalcone, while DcCHI4 lacks this catalytic activity. Moreover, overexpression of DcCHI1 or DcCHI4 enhances flavonoid production in D. cambodiana and tobacco [13]. CHIL does not possess the activity of intramolecular cyclization of naringenin chalcone. However, studies have shown that CHIL functions as a flavonoid enhancer in plants, such as Japanese morning glory [48], Humulus lupulus Linn. [49], and Arabidopsis [50]. In this study, the comprehensive analysis of the transcriptomic and metabolomic data showed that the expression of LpCHI and LpCHIL are significantly and positively correlated with the total content of flavonoids, kaempferol, apigenin, trilobatin, phlorizin, and phloretin (P < 0.01). Yang et al. measured the total content of flavonoids, trilobatin, and phlorizin of L. polystachyus leaves at different development stages. Total flavonoid content is highest in the early stage and decreases at the mature stage [4]. Comprehensive transcriptome and metabolome analysis of L. polystachyus leaves revealed that the expression of key genes in flavonoid biosynthesis pathways is consistent with the content accumulation pattern of the flavonoid, with high expression in young leaves and low expression in mature leaves [1]. In this study, the total flavonoids content is determined by UV-spectrophotometry, and the 6 flavonoids are determined by UPLC. The results show that the total flavonoid content is highly accumulated at the early stage and less accumulated at the mature stage. The expression of key genes in flavonoid biosynthesis pathways at different stages is examined by qRT-PCR, indicating a higher expression in young leaves than in mature leaves. These results of this study are in line with previous results. The environment where plants grow can influence enzyme activities in vivo and regulate metabolic activities. The CO 2 stress enhances phenolic accumulation in the root of Panax ginseng through induction of glucose-6-phosphate (G6PDH) activity [51]. Blue light-treated strawberries induce increased activities of four key enzymes in anthocyanin production: PAL, TAL (tyrosine ammonia-lyase), C4H, and 4CL [51]. In this study, the enzyme activities of PAL, C4H, 4CL, CHS, and CHI involved in flavonoid synthesis are measured. We found some differences in enzyme activities of different samples and different stages. However, there is no correlation between enzyme activity and flavonoid content. The differences in flavonoid metabolites in different types of LpCHI promoter DNA methylation samples are mainly caused by differences in the expression levels of genes involved in flavonoid metabolism in different samples.

Conclusion
DNA methylation at a specific site in the LpCHI promoter CpG island region affects the activity of TFs bind to the promoter, which then regulates LpCHI expression. In this study, DNA methylation at specific sites in the promoter region is found to have a more significant regulatory effect on gene expression than the DNA methylation ratio. Among 16 CpG-sites, 9 TFs binding to them are detected. DNA methylation of -1469 and -1380 cytosines lead to greater transcriptional activation than transcriptional inhibition. In addition, DNA methylation of -1797, -1783, -1735, and -1723 cytosines resulted in more significant transcriptional inhibition than transcriptional activation. Since different TFs differ in their sensitivities to DNA methylation, further experiments are needed to identify the TFs causing differences in LpCHI transcription. Under the same growth environment, the enzyme activities of PAL, C4H, 4CL, CHS, and CHI are different in various samples and leaves at different growth and development stages. However, these activities are not the key factors affecting the accumulation of flavonoids. The transcriptional differences of genes involved in flavonoid biosynthesis are the main factors regulating the accumulation of flavonoids.

Data availability
The raw RNA-seq data supporting the findings of this study have been deposited in the NCBI SRA database under accession number PRJNA771767.

CRediT authorship contribution statement
Limei Lin: designed the research, conducted most of the experiments, analyzed the transcriptome and metabolome, Formal analysis, Writingoriginal draft, analyzed the data and drafted the manuscript.

Declaration of competing interest
The authors report no declarations of interest.