Integrated metabolomic and transcriptomic analyses reveal genetic regulation of anthocyanin accumulation in the peel of three jaboticaba cultivars

Jaboticaba (Myrciaria cauliflora Berg), an evergreen shrub in the family Myrtaceae native to Brazil, is commonly known as the 'tree grape'^[1]. Currently cultivated across various Chinese regions including Fujian, Guangdong, Hainan, and Chongqing, it has emerged as a novel tropical fruit in the domestic market. This dark-colored berry is characterized by its high concentration of phenolic compounds, most notably anthocyanins and ellagitannins^[2−6]. The exceptional abundance of anthocyanins (such as cyanidin and delphinidin) in the peel imparts a near-black hue, leading to its long-standing use as a natural coloring and flavoring agent in candies, jellies, and desserts^[7]. Research indicates that anthocyanins are water-soluble pigments and a subset of flavonoids; the six most common types found in higher plants are malvidin (Mv), pelargonidin (Pg), cyanidin (Cy), delphinidin (Dp), petunidin (Pt), and peonidin (Pn). Beyond their coloring properties, anthocyanins possess potent antioxidant characteristics that offer therapeutic value in preventing neurological diseases, cardiovascular conditions, cancer, and diabetes^[8−16]. Furthermore, these flavonoids can inhibit and reduce cellular damage caused by free radicals generated during metabolic processes.

Anthocyanin biosynthesis is a multi-step, highly regulated metabolic process typically categorized into early biosynthesis genes (EBGs) and late biosynthesis genes (LBGs) based on the functional order of enzymes^[17]. EBGs primarily include phenylalanine ammonia-lyase (PAL), cinnamate-4-hydroxylase (C4H), 4-coumarate: CoA ligase (4CL), chalcone synthase (CHS), and chalcone isomerase (CHI), which catalyze the initial steps of the phenylpropanoid and flavonoid pathways to form the basic flavanone skeleton. Subsequently, enzymes such as flavanone 3-hydroxylase (F3H), flavonoid 3',5'-hydroxylase (F3'5'H), and flavonoid 3'-hydroxylase (F3'H) perform hydroxylation modifications on this core structure to generate dihydrokaempferol and dihydroquercetin, thereby determining the branching direction of flavonoids and the diversity of downstream products. LBGs, such as dihydroflavonol 4-reductase (DFR), anthocyanidin synthase (ANS/LDOX), and UDP-glucose: flavonoid 3-O-glucosyltransferase (UFGT), are responsible for catalyzing the final steps, completing the reduction and glycosylation processes from dihydroflavonols to stable, colored anthocyanins^[18,19]. This 'EBGs-LBGs' stratification represents a classic model in anthocyanin biosynthesis research^[20]. Generally, anthocyanin synthesis is governed by two categories of genes: structural genes that directly participate in the biosynthetic pathway, and regulatory genes—predominantly the MBW complex—that exert transcriptional control over these structural genes^[21,22].

Integrating targeted metabolomics and transcriptomics analysis helps elucidate the coloring mechanisms in the peels of different jaboticaba cultivars. Unstable free anthocyanins are stabilized through a series of modifications, including UFGT-catalyzed glycosylation followed by acylation and methylation. These modifications significantly improve the stability and solubility of anthocyanins, ultimately forming functional pigments^[23]. The peel of jaboticaba features cyanidin-3-O-glucoside and delphinidin-3-O-glucoside as the major anthocyanin components. These pigments determine the typical dark purple to black-purple appearance of the peel and are key anthocyanins for peel coloration and antioxidant activity. Due to the high content and stable biological activity of anthocyanins in jaboticaba peel, it serves as a high-quality functional raw material rich in cyanidin-based anthocyanins^[24−26]. Cyanidin-type and delphinidin-type anthocyanins function not only in jaboticaba peel but also in various other fruit trees, where peel coloration is closely related to the accumulation of cyanidin, delphinidin, and their glycosides. For example, the peel of black-skinned pomegranate continuously accumulates high levels of anthocyanins throughout its development; its total peel anthocyanin content can be two orders of magnitude higher than that of common pomegranate varieties, and it is dominated by cyanidin-type and delphinidin-type anthocyanins, which serve as the direct pigment basis for its purple-black appearance^[27]. In dark-colored grape varieties such as black grapes, the peel is enriched with anthocyanins like delphinidin-3-O-glucoside, while the outer tissue of blackberry fruit features cyanidin-3-O-glucoside as the main anthocyanin component, with both exhibiting high antioxidant capacities^[28]. In fruit vegetables like eggplant, the deep purple to nearly black peel is similarly colored by epicarp anthocyanins, mainly consisting of delphinidin glycosides^[29]. In the Osmanthus fragrans cultivar 'Ziyingui', the degreening of the pericarp is synchronous with the temporal accumulation of cyanidin-type anthocyanins, and the color transition from green to purple is significantly and positively correlated with cyanidin accumulation levels^[30]. Similarly, deep purple berries such as blackcurrants possess peels rich in anthocyanins dominated by cyanidin and delphinidin^[31].

This study systematically elucidates the regulatory networks governing key anthocyanin metabolism in the fruit peels of three jaboticaba cultivars during ripening. An integrated analysis of targeted metabolomics and transcriptomics revealed the expression profiles of key anthocyanin biosynthesis-related genes and their potential transcription factors^[32,33]. The results indicate that McMYB4, McMYB24, McMYB73, and McbHLH130 likely participate in the regulation of anthocyanin biosynthesis by modulating the expression of relevant structural genes.

In this study, three jaboticaba cultivars—'White', 'Evergreen', and 'Sabara'—were selected due to the distinct differences in their pericarp coloration: 'White' exhibits a green peel, 'Evergreen' a reddish-black peel, and 'Sabara' a deep black peel. All samples were harvested at the peak of maturity from three vigorously growing trees cultivated by Qionghai Shengda Agricultural Science and Technology Development Co., Ltd. Following collection, the samples were immediately flash-frozen in liquid nitrogen and stored at −80 °C for subsequent analysis.

Targeted metabolomics analysis based on LC-MS/MS

Targeted anthocyanin detection was performed by Sanshu Biotech Co., Ltd. First, the sample solution was weighed into a centrifuge tube and mixed with 1 mL of extraction solvent (methanol : water : formic acid = 70:30:1, [v/v/v]) via vortexing. Extraction was conducted using ultrasonication for 20 min, followed by centrifugation at 12,000 rpm for 10 min; this procedure was repeated at least twice. The combined supernatants were purified using HLB solid-phase extraction (SPE). The HLB-SPE cartridges were pre-conditioned with 1 mL of 100% methanol and equilibrated with 1 mL of deionized water. After loading the supernatant, the column was washed with 1 mL of deionized water and eluted with 1 mL of methanol containing 5% formic acid. The eluate was evaporated under nitrogen and lyophilized, and the resulting residue was reconstituted in 0.2 mL of methanol for extraction. Instrumental analysis was conducted using a Thermo Scientific Vanquish UPLC system coupled with a Q Exactive high-resolution mass spectrometer. Chromatographic separation was achieved on a Waters HSS T3 column (50 mm × 2.1 mm, 1.8 μm) maintained at 40 °C. The mobile phase consisted of 0.1% formic acid in ultrapure water (A) and 0.1% formic acid in acetonitrile (B) at a flow rate of 0.3 mL/min with a 2 μL injection volume. The 10-min gradient program was set as follows: 0–1 min, 5% B; 6 min, 30% B; 7–8 min, 95% B; and 8.1–10 min, 5% B^[34]. Mass spectrometry was operated in ESI positive ion mode using Selected Ion Monitoring (SIM) with the following parameters: sheath gas 40 arb, auxiliary gas 10 arb, spray voltage 3,000 V, and capillary temperature 350 °C^[35]. Identification was performed using TraceFinder software by comparing retention times and key fragments with standards and databases. Quantitative analysis utilized the external standard method, with calibration curves (range 50–2,000 ng/mL) established via linear regression forced through the origin (Y = kX) using equal weighting^[36].

RNA-seq analysis

Utilizing the BGI (Sanya) sequencing platform, mRNA was isolated from total RNA via the Oligo(dT) magnetic bead enrichment method. The resulting RNA was fragmented into approximately 300 bp templates using ions, followed by random primer-mediated reverse transcription to synthesize double-stranded cDNA. The library underwent PCR amplification and size selection for 450 bp fragments, with quality control monitored using the Agilent 2100 Bioanalyzer system. Qualified libraries were pooled in equimolar ratios (2 nM) based on their respective indices. After alkaline denaturation, paired-end sequencing was performed on the Illumina platform. For the generated raw data, stringent quality filtering—including the removal of adapters and low-quality reads—was conducted using Fastp. High-quality reads were mapped to the target genome (jaboticaba) using HISAT2 (v2.2.1)^[37−39]. Quantitative assessment was performed with featureCounts based on the resulting BAM files. Finally, data normalization was executed using the Transcripts Per Million (TPM) algorithm to establish the gene expression profiles.

Differential expression analysis and functional annotation

Differentially expressed genes (DEGs) were identified using DESeq2 (v1.48.1), with screening criteria set at |log2FoldChange| > 1 and p-adj < 0.05. For the metabolomic data, differentially accumulated metabolites (DAMs) were determined based on |log2FoldChange| > 1, p < 0.05, and FDR < 0.05. Protein sequences were identified using eggNOG-mapper (https://eggnog-mapper.embl.de/) and TBtools, while transcription factors were predicted via the iTAK database (https://itak.feilab.net/cgi-bin/itak/index.cgi). To further uncover biological pathways, we utilized the KEGGREST (v1.48.0) package. Subsequently, KEGG pathway enrichment analysis for both DEGs and DAMs was executed using clusterProfiler (v4.16.0) to identify significantly overrepresented biological processes.

Combined transcriptomic and metabolomic analysis and network construction

Soft clustering analysis was performed on the differentially expressed genes among jaboticaba varieties using the Mfuzz (v2.68.0) package to characterize dynamic expression trends. Subsequently, WGCNA (v1.73) was employed to identify functional modules closely associated with anthocyanin content, with the soft-thresholding power set to nine to achieve a scale-free topological network. To pinpoint core regulatory nodes, components were screened based on a weight threshold of 0.5 and high gene significance (GS) for anthocyanin levels. Finally, an integrated regulatory network was constructed by incorporating transcription factor (TF) data and visualized using Cytoscape (v3.8.0).

Analysis of cis-acting regulatory elements

Cis-acting elements within the promoter regions of target genes were identified using the FIMO (Find Individual Motif Occurrences) tool from the MEME Suite (v5.5.5). Transcription factor binding sites (TFBS) were scanned based on motif reference data retrieved from the JASPAR CORE (Plantae) database (https://jaspar.genereg.net/). This analysis enabled the high-sensitivity detection and precise mapping of the spatial distribution of candidate motifs along the promoter sequences.

Dual-LUC reporter assay

To verify the regulatory effects of McMYB4 on the promoter activities of McPAL_2 and McC4H, the promoter sequences of McPAL_2 and McC4H were cloned into the pGreenII 0800-LUC reporter vector to construct the reporter vectors pGreenII 0800-LUC-McPAL_2pro and pGreenII 0800-LUC-McC4Hpro, respectively, while the coding sequence of McMYB4 was cloned into the pGreenII-62SK effector vector to yield the effector vector pGreenII-62SK-McMYB4. For the experimental setup, the co-transformation of each promoter with McMYB4 served as the experimental group, whereas promoter + empty vector, empty vector + McMYB4, and empty vector + empty vector served as negative controls. The constructed vectors were introduced into Agrobacterium tumefaciens strain GV3101 (pSoup-p19) via electroporation, and single colonies were inoculated into YEB liquid medium containing the appropriate antibiotics and cultured at 28 °C. After harvesting the bacterial cells by centrifugation at 4,000 rpm for 5 min, the pellet was resuspended in a buffer containing 10 mM MgCl₂ and 120 μM AS to adjust the OD₆₀₀ to approximately 1.0, and the suspension was allowed to stand for over 1 h at room temperature (28 °C) or over 3 h at 4 °C prior to use. Equal volumes of the effector and reporter bacterial cultures were mixed and infiltrated into the abaxial surface of tobacco (Nicotiana benthamiana) leaves. The infiltrated plants were then cultured under low-light conditions for 48 h. Subsequently, the leaves were harvested, uniformly coated on the abaxial side with a 0.2 mg/mL potassium luciferin solution, and incubated in the dark for 5–10 min before fluorescence signals were captured using a Tanon 5200 Automated Chemiluminescence Imaging System. Dual-luciferase activity was measured using the Dual-Luciferase Reporter Assay System (Vazyme), where the firefly luciferase (LUC) activity served as the reporter signal and the Renilla luciferase (REN) activity served as the internal control; the regulatory effect of McMYB4 on the transcriptional activities of the McPAL_2 and McC4H promoters was evaluated by calculating the LUC/REN ratio^[40−43].

In this study, targeted anthocyanin detection was employed to analyze the absolute content of anthocyanins in the pericarp of three jabuticaba cultivars: 'White' (green pericarp), 'Evergreen' (red-black pericarp), and 'Sabara' (dark black pericarp) (Fig. 1a). The Total Ion Chromatogram (TIC) of different samples exhibited complete peak shapes and stable retention times, confirming the reliability of the method and providing a solid foundation for subsequent data analysis. During instrumental analysis, quality control (QC) samples were uniformly inserted into the sequence, with all QC relative standard deviations (RSD) below 30%, indicating a stable detection system (Supplementary Fig. S1; Supplementary Table S1). PCA showed that PC1 and PC2 collectively explained 99.1% of the total variation (66.9% and 32.2%, respectively), reflecting significant metabolic differences between samples and high data reliability due to the tight clustering of biological replicates. Cluster analysis further demonstrated distinct separation between 'White', 'Evergreen', and 'Sabara', with minimal intra-group variation (Fig. 1b). Among the 15 metabolites detected, cyanidin-3,5-diglucoside, cyanidin-3-galactoside, and procyanidin B4 had an absolute content of zero across all cultivars. The remaining 12 substances primarily included six anthocyanins (delphinidin, cyanidin, petunidin, pelargonidin, peonidin, and malvidin), three flavonols (quercetin, kaempferol, and isorhamnetin), and one each of procyanidins (procyanidin B2), flavones (luteolin), and flavonoid glycosides (rutin) (Fig. 1c). Significant differences were observed between cultivars: except for quercetin and rutin, the concentrations of the other 10 substances were higher in 'Sabara' than in 'White'; meanwhile, 'Evergreen' showed higher levels of only rutin, luteolin, malvidin, delphinidin, cyanidin, and isorhamnetin compared to 'White'. Notably, 'Sabara' contained the highest levels of all anthocyanins except for pelargonidin and malvidin (Fig. 1d; Supplementary Table S2). Based on the screening criteria of |log2FoldChange| > 1 and p-value < 0.05, various differentially accumulated metabolites (DAMs) were identified. In the Sabara_vs_White comparison, nine DAMs were identified, all of which were upregulated. In Evergreen_vs_White, eight DAMs were identified, including five upregulated (cyanidin, malvidin, rutin, luteolin, and isorhamnetin) and three downregulated (petunidin, pelargonidin, and kaempferol). In Sabara_vs_Evergreen, seven DAMs were found; while pelargonidin and peonidin were absent in 'Evergreen', six metabolites were upregulated in 'Sabara', and only rutin was downregulated (Fig. 1e–g; Supplementary Table S3).

Figure 1. 

Targeted anthocyanin data analysis of the pericarp from three jaboticaba cultivars. (a) Phenotypic differences in pericarp coloration among 'White', 'Evergreen', and 'Sabara' cultivars. (b) Principal component analysis on the anthocyanin targeted metabolomics data of the fruit peels from these three cultivars using the "prcomp" function in the R stats package. (c) Classification histograms of 12 detected metabolites plotted using the ggplot2 package in R. (d) Line graph representing the targeted anthocyanin data across the three cultivars (log-transformed and z-standardized). (e)–(g) Radar plots plotted using the ggplot2 package in R to show the differentially accumulated metabolites (DAMs) compared to Evergreen and White, Sabara and White, and Sabara and Evergreen.

Transcriptome analysis of three cultivars of jaboticaba

Each sample yielded at least 6 GB of high-quality clean reads, with Q20 and Q30 values exceeding 97% and 92%, respectively, and GC content ranging from 50.16% to 51.28%, indicating high sequencing data quality. A total of 31,235 genes were detected across all samples, and after filtering out null/missing values and invariant genes, 26,090 genes were retained for downstream analysis (Supplementary Table S4). Principal component analysis (PCA) demonstrated that the three biological replicates for each jaboticaba cultivar clustered tightly, reflecting excellent experimental reproducibility and minimal intra-group variation. Furthermore, the 'White', 'Evergreen', and 'Sabara' cultivars were clearly separated without overlap; PC1 and PC2 explained 47.4% and 27.3% of the total variance, respectively (74.7% cumulative), indicating significant transcriptomic divergence among the cultivars (Fig. 2a). Differential expression gene (DEG) screening was performed across the comparisons. The 'Evergreen'_vs_'White' comparison yielded 3,913 DEGs (2,435 upregulated, 1,478 downregulated), while 'Sabara'_vs_'White' exhibited 6,309 DEGs (3,240 upregulated, 3,069 downregulated), significantly higher than the other groups. In the 'Sabara'_vs_'Evergreen' comparison, 4,055 DEGs were identified (1,741 upregulated, 2,314 downregulated) (Fig. 2b), suggesting that 'Sabara' possesses the most distinct gene expression profile. A core set of 823 DEGs was common to all three comparisons (Fig. 2c; Supplementary Table S5). To further characterize expression profiles, soft clustering categorized the DEGs into six clusters (Fig. 2d; Supplementary Table S6), including continuous decline (Clusters 1 and 5), increase-then-decline (Clusters 2 and 6), and continuous increase (Cluster 4). Notably, the expression trend of Cluster 4 aligned with the anthocyanin content across cultivars. KEGG enrichment analysis of Cluster 4 (Fig. 2e; Supplementary Table S7) revealed significant enrichment in pathways closely related to anthocyanin synthesis, particularly 'phenylpropanoid biosynthesis' and 'flavonoid biosynthesis', which form the core metabolic framework for anthocyanin production. Additionally, hormone-related and upstream signaling pathways (e.g., 'plant hormone signal transduction', 'MAPK signaling pathway–plant') were enriched, suggesting that hormonal regulation and signaling networks may indirectly modulate anthocyanin accumulation by influencing structural genes and transcription factors. Multiple carbohydrate metabolism pathways (e.g., 'starch and sucrose metabolism', 'glycolysis/gluconeogenesis') were also enriched, likely providing the necessary energy and carbon sources. These KEGG results suggest that anthocyanin biosynthesis may be governed by a comprehensive network involving direct biosynthetic pathways, hormonal signaling, stress responses, and primary metabolism.

Figure 2. 

Transcriptome analysis of the pericarp from three jabuticaba cultivars. (a) Principal component analysis (PCA) of the transcriptomes for 'White', 'Evergreen', and 'Sabara'. (b) Differential expression analyses for the three comparison groups—Evergreen_vs_White, Sabara_vs_White, and Sabara_vs_Evergreen—using DESeq2, and the identity and number of differentially expressed genes (DEGs) for each. (c) Venn analysis of the differentially expressed genes (DEGs) across the Evergreen_vs_White, Sabara_vs_White, and Sabara_vs_Evergreen comparison groups using the VennDiagram package in R. (d) Mfuzz soft clustering analysis on all differentially expressed genes (DEGs) using the Mfuzz package in R. (e) KEGG enrichment analysis on the differentially expressed genes in the fourth cluster using the clusterProfiler and KEGGREST packages in R, with red-highlighted entries representing pathways related to anthocyanin biosynthesis.

WGCNA analysis

To identify genes related to pericarp pigmentation in three jaboticaba cultivars, we implemented a weighted gene co-expression network analysis (WGCNA) incorporating the TPM values of all genes and the absolute contents of DAMs (cyanidin, delphinidin, pelargonidin, luteolin, and kaempferol) from the flavonoid biosynthesis pathway. The genes were categorized into 33 distinct co-expression modules, and the turquoise module was found to be highly correlated with anthocyanin synthesis, with correlation coefficients greater than 0.9 and reaching extreme significance (Fig. 3a). Through KEGG enrichment of the genes in the turquoise module (Fig. 3b; Supplementary Table S8), a total of 127 pathways were significantly enriched. Among them, 16 structural genes (three PAL, one F3'H, four 4CL, one FLS, one C4H, two CHI, one CHS, one F3H, one F3'5'H, and one ANS) in the phenylpropanoid biosynthesis and flavonoid biosynthesis pathways were found to be involved in anthocyanin biosynthesis. Except for FLS and 4CL_2, all these structural genes were significantly positively correlated with anthocyanins, especially delphinidin and cyanidin, which showed extremely high correlations with the structural genes (r > 0.9, p < 0.05). Furthermore, the absolute module membership (MM) values of these structural genes in the turquoise module were all greater than 0.8 (Supplementary Table S9). 'Sabara' exhibited significantly higher gene expression in the turquoise module, showing stronger activity or a more positive regulatory state in the biological processes related to this module, while 'White' and 'Evergreen' showed lower gene expression, though 'Evergreen' was higher than 'White' (Fig. 3c). The module eigengene expression values of 'Sabara' were significantly higher than other samples and the baseline, all characterized by high expression; meanwhile, the eigengene expression values of 'Evergreen' and 'White' were significantly below the baseline, with 'White' being lower than 'Evergreen', implying that some related biological processes in this module might not be effectively activated in 'Evergreen' and 'White' (Fig. 3d).

Figure 3. 

Identification of regulatory modules associated with anthocyanin accumulation using weighted gene co-expression network analysis (WGCNA). (a) Heatmap of WGCNA module-trait correlations utilizing the WGCNA, pheatmap, and Hmisc R packages. (b) KEGG pathway enrichment analysis for genes in the MEturquoise module using the clusterProfiler and KEGGREST packages in R. (c) Heatmap of gene expression patterns within the MEturquoise module using the pheatmap package in R. (d) Bar chart showing the expression levels of the module eigengene for the MEturquoise module.

Integrated transcriptomic and metabolomic analysis

By integrating targeted metabolomics and transcriptomics analysis, this study provides an in-depth resolution of the key regulatory mechanisms underlying anthocyanin biosynthesis in three jaboticaba cultivars. Within the turquoise module, 16 structural genes involved in anthocyanin synthesis were identified. Further correlation with Cluster 4 from the Mfuzz analysis highlighted nine key structural genes (including PAL_2, 4CL_4, C4H, CHS, CHI_1, CHI_2, F3H, F3'H, and ANS) that exhibited significant differential expression across the cultivars (Fig. 4). These genes were most highly expressed in the 'Sabara' cultivar, which possesses the highest anthocyanin content, followed by 'Evergreen', while the lowest expression levels were observed in the 'White' cultivar. Metabolomic profiling further confirmed that the accumulation patterns of major anthocyanins, such as delphinidin, cyanidin, and pelargonidin, as well as related flavonoids like luteolin, are highly consistent with the expression trends of these structural genes.

Figure 4. 

Combined analysis. Biosynthetic pathways of five DAMs (indicated in red font), illustrating the metabolite accumulation levels and expression profiles of key enzyme genes across the 'White', 'Evergreen', and 'Sabara' cultivars (normalized using log-transformed and z-standardized methods). Key structural genes identified in Cluster 4 are highlighted in red. Enzyme abbreviations are as follows: CHS, chalcone synthase; CHI, chalcone isomerase; F3H, flavanone 3-hydroxylase; F3'H, flavonoid 3'-hydroxylase; F3'5'H, flavonoid 3',5'-hydroxylase; FLS, flavonol synthase; and ANS, anthocyanidin synthase.

Regulatory network of key structural genes, transcription factors, anthocyanin accumulation, and analysis of cis-acting elements

By employing iTAK to predict transcription factors (TFs) within the turquoise module and conducting rigorous co-expression analysis (weight > 0.5), a high-confidence regulatory network for anthocyanin biosynthesis was constructed. This network is centered around 9 key structural genes and 19 differentially expressed TFs, comprising 7 McMYBs and 12 McbHLHs, with a total of 143 high-confidence interactions identified between these TFs and structural genes (Supplementary Table S10). The expression levels of these structural genes exhibited a significant positive correlation with the accumulation of delphinidin, pelargonidin, and cyanidin (Fig. 5a). Notably, these TFs displayed contrasting expression patterns in 'White' and 'Sabara', while the expression levels in 'Evergreen' remained intermediate (Fig. 5b). Furthermore, FIMO motif scanning of the promoter regions of the nine key structural genes revealed a widespread distribution of MYB and bHLH binding sites; specifically, potential binding sites for McMYB4, McMYB24, and McMYB73 were identified in the promoters of genes such as 4CL_4, C4H, CHI, PAL_2, F3'H, F3H, and ANS, whereas McbHLH130 binding sites were primarily located within the promoters of CHI_1, F3H, and 4CL_4 (Fig. 5c; Supplementary Table S11).

Figure 5. 

Co-expression network and cis-acting element analysis. (a) A co-expression network of McMYB and McbHLH transcription factors, structural genes, and metabolites within the MEturquoise module was constructed using the WGCNA package in R and then beautified using Cytoscape. (b) TPM expression heatmaps of the McMYB and McbHLH transcription factor families were generated using the pheatmap package in R. (c) Motif scans of promoter regions of nine key structural genes were performed using the JASPAR CORE plant database to predict their cis-regulatory elements (plotted with TBtools).

McMYB4 represses the promoter activities of McPAL_2 and McC4H

To verify the regulatory effect of McMYB4 on the promoter activities of McPAL_2 and McC4H, a dual-luciferase reporter assay was performed using a tobacco transient expression system. The results showed that, compared with the negative control, the LUC luminescence signals of both McPAL_2 and McC4H promoters were markedly weakened upon co-expression with McMYB4 (Fig. 6a, c). Further determination of the dual-luciferase activity revealed that McMYB4 significantly decreased the LUC/REN ratios of the McPAL_2 and McC4H promoters (p < 0.001) (Fig. 6b, d), indicating that McMYB4 can inhibit the transcriptional activities of these two promoters (Supplementary Tables S12−S13). Combined with the gene expression analysis, the expression of McMYB4 was negatively correlated with that of McPAL_2 and McC4H. Taken together, these results demonstrate that McMYB4 may function as a negative transcriptional regulator to modulate the expression of McPAL_2 and McC4H, thereby influencing the expression of key structural genes in the phenylpropanoid metabolic pathway, and potentially further participating in the biosynthesis of secondary metabolites such as flavonoids and anthocyanins.

Figure 6. 

Dual-luciferase reporter assay demonstrating the regulatory effects of McMYB4 on the promoters of McPAL_2 and McC4H. (a) Transient dual-luciferase assay of McMYB4 on the McPAL_2 promoter in Nicotiana benthamiana leaves. Region 1 represents the experimental group co-infiltrated with pGreenII 0800-LUC-McPAL_2pro and pGreenII-62SK-McMYB4; Region 2 represents the negative control co-infiltrated with pGreenII 0800-LUC-McPAL_2pro and the empty pGreenII-62SK vector; Region 3 represents the negative control co-infiltrated with the empty pGreenII 0800-LUC vector and pGreenII-62SK-McMYB4; Region 4 represents the negative control co-infiltrated with the empty pGreenII 0800-LUC and empty pGreenII-62SK vectors. (b) Quantification of dual-luciferase activity (LUC/REN) driven by the McPAL_2 promoter. (c) Transient dual-luciferase assay of McMYB4 on the McC4H promoter in N. benthamiana leaves. The four infiltration regions correspond to those described in panel (a), with McPAL_2pro replaced by McC4Hpro. (d) Quantification of dual-luciferase activity (LUC/REN) driven by the McC4H promoter. Data are presented as the mean ± SE of three biological replicates. Asterisks indicate significant differences between treatments (*** p < 0.001).

Anthocyanins are key compounds determining the diverse hues of plant organs. Specifically, the high accumulation of cyanidin-type anthocyanins provides the biochemical foundation for the black or dark skin observed in various plants^[44]. The intensity of the color is directly proportional to the concentration of anthocyanin^[45], as seen in black pomegranate, blackberry, blackcurrant, and black mulberry^[27,46,47]. In this study, cyanidin-type derivatives (cyanidin-3,5-diglucoside and cyanidin-3-galactoside) were not detected in the three jaboticaba cultivars; however, the levels of cyanidin were exceptionally high. This suggests that cyanidin exists in its free form within the pericarp of jaboticaba. Although free cyanidin is highly unstable in aqueous solutions and easily influenced by pH, metal ions, or co-pigments^[48], high concentrations of free cyanidin under specific cellular conditions—such as co-pigmentation, metal complexation, high pigment aggregation, and optimal vacuolar pH—can significantly enhance light absorption and shift the color toward deep purple or near-black^[49−51]. Furthermore, delphinidin and its O-glycoside derivatives are also vital pigment bases for deep purple or near-black phenotypes in many plants, as they increase molecular absorbance and cause a bathochromic shift, deepening the tissue color^[52]. The variation in the relative content of delphinidin-type and cyanidin-type anthocyanins, along with their synergistic effects in enhancing light absorption efficiency and promoting metal complexation or co-pigmentation, ultimately dictates the presentation of dark or black colors^[53]. In contrast, pelargonidin and its derivatives typically impart red or orange tones and show a lower correlation with dark or black phenotypes than cyanidin and delphinidin^[54,55]. Among the three jaboticaba cultivars, pelargonidin content was relatively low compared to cyanidin and delphinidin—it was even undetectable in 'Evergreen', while being relatively higher in 'Sabara'.

Targeted metabolomics analysis revealed that the absolute contents of free cyanidin and delphinidin were high across all three jaboticaba cultivars, whereas free pelargonidin and flavonoids such as luteolin and kaempferol were present in lower concentrations. Previous studies have indicated that while free anthocyanins possess poor intrinsic stability, they can undergo copigmentation with luteolin and kaempferol to enhance both color stability and intensity^[56,57]. In the 'White', 'Evergreen', and 'Sabara' cultivars, the primary anthocyanins (cyanidin and delphinidin) exhibited a gradient distribution, with their concentration trends aligning with the observed peel coloration. Notably, in 'Sabara', the cyanidin content reached as high as 82,477 ng/100 mg, with log2FC values of 10.71 and 5.06 in the 'Sabara_vs_White' and 'Sabara_vs_Evergreen' comparisons, respectively. Conversely, pelargonidin levels remained low in all three cultivars and were undetectable in 'Evergreen'. Our research demonstrates that cyanidin and delphinidin are the predominant anthocyanins influencing peel color in these three jaboticaba cultivars; specifically, the exceptionally high concentration of cyanidin is a critical factor underlying the deep black coloration of the 'Sabara' peel.

Combined analysis revealed that 16 structural genes participate in anthocyanin biosynthesis. Among these, nine structural genes (one PAL, one 4CL, one C4H, one F3H, CHS, two CHI, one F3'H, and one ANS) were consistently identified in both Mfuzz Cluster 4 and the WGCNA turquoise module. Genes in Cluster 4 mirrored the accumulation trends of major anthocyanins across the three cultivars, while the turquoise module was identified as the primary module regulating anthocyanin and flavonoid biosynthesis. These nine structural genes exhibited significant positive correlations with anthocyanin content (p < 0.05), with Gene Significance (GS) values exceeding 0.9 for cyanidin and delphinidin and module membership (MM) values surpassing 0.8, although GS values for pelargonidin were lower. These genes function at different stages of the biosynthetic pathway: upstream genes such as PAL_2, C4H, and 4CL_4 convert phenylalanine into p-coumaroyl-CoA, providing precursors for downstream synthesis; CHS, CHI_1, CHI_2, and F3H drive the conversion from chalcones to dihydroflavonols, serving as critical branch points; F3'H influences the hydroxylation of dihydroflavonols to determine anthocyanin types and promote cyanidin synthesis^[58], consistent with the finding that cyanidin levels far exceed delphinidin; and ANS acts as the core enzyme catalyzing the formation of anthocyanin monomers. In the three jaboticaba cultivars, these nine genes showed significant differential expression: they were highly expressed in 'Sabara' (highest anthocyanin content), lowest in 'White' (lowest content), and intermediate in 'Evergreen'. Metabolomic data confirmed that the accumulation patterns of major anthocyanins (delphinidin, cyanidin, pelargonidin) and flavonoids (Luteolin) were highly consistent with these expression trends. Based on a weight > 0.5, seven differentially expressed MYBs and 12 bHLH transcription factors (TFs) co-expressed with these structural genes were identified in the turquoise module. FIMO analysis of the 2,000 bp upstream promoter regions revealed potential binding sites for McMYB4, McMYB24, McMYB73, and McbHLH130, suggesting direct transcriptional regulation. Research has shown that MYB and bHLH transcription factors often regulate downstream gene expression through the formation of protein complexes rather than acting independently^[59]. In plants, MYB–bHLH complexes are widely involved in secondary metabolism, developmental regulation, and stress responses^[60]. A typical regulatory mechanism involves MYB proteins recognizing specific DNA cis-acting elements via their R2R3 domains, while bHLH proteins bind to E-box or G-box elements. Meanwhile, the two proteins form stable complexes through protein–protein interaction domains, thereby enhancing transcriptional activation efficiency or improving regulatory specificity^[61]. For example, MYB4 belongs to the R2R3-MYB transcription factor family and generally functions as a transcriptional repressor in Arabidopsis thaliana, where it participates in the regulation of phenylpropanoid and flavonoid metabolism. MYB4 can interact with bHLH and other proteins to interfere with the activation function of the MBW complex and has been demonstrated to directly repress the expression of key structural genes, thereby reducing the accumulation of flavonoid and anthocyanin metabolites^[62−65]. In our study, the expression of McMYB4 exhibited a significant negative correlation (r < −0.9, p < 0.05) with several key structural genes, including CHI_1, ANS, CHI_2, F3H, C4H, and PAL_2. To validate the regulatory role of McMYB4, dual-luciferase reporter assays were performed using the promoters of the upstream structural genes McPAL_2 and McC4H. The results demonstrated that the co-expression of McMYB4 significantly suppressed the activities of both promoters, as evidenced by a marked decrease in the LUC signal and the LUC/REN ratio. These findings are highly consistent with the negative correlations observed in the transcriptomic data, suggesting that McMYB4 functions as a negative regulator in the phenylpropanoid biosynthetic pathway. By inhibiting the expression of these pivotal upstream structural genes, McMYB4 likely restricts the metabolic flux toward flavonoid and anthocyanin biosynthesis, thereby limiting the accumulation of these secondary metabolites. Consequently, while the coordinated upregulation of structural genes drives anthocyanin accumulation, McMYB4 functions as a crucial negative regulatory component that fine-tunes anthocyanin biosynthesis, which may maintain the homeostasis of flavonoid and anthocyanin accumulation in jaboticaba pericarp.

In this study, a comprehensive analysis of the transcriptome data and targeted anthocyanin content was conducted across three jaboticaba cultivars: 'White', 'Sabara', and 'Evergreen'. The results indicate that cyanidin and delphinidin are the primary anthocyanins governing fruit peel coloration, with cyanidin playing the most significant role; furthermore, 'Sabara' exhibited significantly higher anthocyanin levels than both 'White' and 'Evergreen'. Based on Mfuzz clustering and WGCNA analysis, nine structural genes involved in the entire anthocyanin biosynthetic pathway were identified, showing expression levels highly consistent with anthocyanin accumulation—peaking in 'Sabara', remaining lowest in 'White', and showing intermediate levels in 'Evergreen'. Simultaneously, the co-expression network identified seven McMYB and 12 McbHLH transcription factors (TFs) highly co-expressed with these structural genes, displaying contrasting expression patterns between 'White' and 'Sabara'. Notably, potential binding sites for McMYB4, McMYB24, McMYB73, and McbHLH130 were discovered within the promoter regions of the structural genes. Furthermore, dual-luciferase reporter assays demonstrated that McMYB4 significantly suppressed the activities of the McPAL_2 and McC4H promoters, and the results were consistent with those obtained from the WGCNA co-expression network analysis. These findings regarding key structural genes and potential regulatory factors provide a critical foundation for further elucidating the molecular regulatory network of peel coloration in jaboticaba.