Integrative analysis of metabolome and transcriptome reveals the mechanism of color formation in cassava (Manihot esculenta Crantz) leaves

Cassava (Manihot esculenta Crantz) leaves are often used as vegetables in Africa. Anthocyanins possess antioxidant, anti-inflammatory, anti-cancer, and other biological activities. They are poor in green leaves but rich in the purple leaves of cassava. The mechanism of anthocyanin’s accumulation in cassava is poorly understood. In this study, two cassava varieties, SC9 with green leaves and Ziyehuangxin with purple leaves (PL), were selected to perform an integrative analysis using metabolomics and transcriptomics. The metabolomic analysis indicated that the most significantly differential metabolites (SDMs) belong to anthocyanins and are highly accumulated in PL. The transcriptomic analysis revealed that differentially expressed genes (DEGs) are enriched in secondary metabolites biosynthesis. The analysis of the combination of metabolomics and transcriptomics showed that metabolite changes are associated with the gene expressions in the anthocyanin biosynthesis pathway. In addition, some transcription factors (TFs) may be involved in anthocyanin biosynthesis. To further investigate the correlation between anthocyanin accumulation and color formation in cassava leaves, the virus-induced gene silencing (VIGS) system was used. VIGS-MeANR silenced plant showed the altered phenotypes of cassava leaves, partially from green to purple color, resulting in a significant increase of the total anthocyanin content and reduction in the expression of MeANR. These results provide a theoretical basis for breeding cassava varieties with anthocyanin-rich leaves.

Introduction

Cassava (Manihot esculenta Crantz) is the sixth staple food crop in the world (An et al., 2019; Cai et al., 2021), which feeding more than a billion people, particularly in sub-Saharan Africa (Kumba, 2012). Cassava is an energy-producing tuberous crop with starchy root which is not only used for industrial materials but also acts as food in southern China. The leaves are also a kind of forage full of protein, flavonoids, carotenoids, and so on. In Africa, cassava leaves are usually consumed as a vegetable in daily life. Most cassava leaves are green in color; however, different cassava germplasms have different characteristics in the leaves. Some specific cassava germplasms have purple leaves, which are rich in flavonoids, especially anthocyanin. Recent studies indicated that flavonoids are detected in the leaves or storage roots of cassava (Omar et al., 2012; Xiao et al., 2021; Fu et al., 2022). However, metabolomic and genetic differences in anthocyanin biosynthesis of cassava leaves in those specific germplasms are poorly understood.

Anthocyanin is an important component of flavonoids (Chen et al., 2014), which plays several functions in plant development. Anthocyanins are involved in plant physiological processes such as attracting pollinators and seed dispersers, combatting abiotic stress, and promoting plant defense (Forkmann, 2010; Zaynab et al., 2018; Sharma et al., 2019). As water-soluble pigments contribute to the majority of the purple, blue, and red colors, anthocyanins are involved in the color formation of leaves, flowers, and fruits (Bosse et al., 2021). As a strong antioxidant, anthocyanins also have a wide range of medicinal values such as antidiabetic, anticancer, antineurodegeneration, and antineuroinflammation, as well as treating cardiovascular diseases (Bai et al., 2019; Calis et al., 2020; Park et al., 2022; Ciumărnean et al., 2020). Therefore, understanding the biosynthesis and regulation mechanism of anthocyanins is very important for breeding anthocyanin-rich cassava varieties, especially for the nutritional intake of Africans.

The biosynthesis of anthocyanins is composed of three branches and is catalyzed by a series of enzymes in the phenylpropanoid and isoflavonoid biosynthetic pathways (Quattrocchio et al., 2006; Ferreyra et al., 2012). Anthocyanin is synthesized from malonyl CoA and coumaroyl-CoA. Naringenin is produced from 1 coumaroyl CoA and 3 malonyl-CoAs, desaturated by chalcone synthase (CHS), and then isomerized by chalcone isomerase (CHI). As a core intermediate for anthocyanin biosynthesis, naringenin is catalyzed by a series of enzymes to form different types of flavonoids, such as flavonoid 3′-hydroxylase (F3′H), flavonoid 3′5′-hydroxylase (F3′5′H), and naringenin 3-dioxygenase (N3D), then catalyzed by dihydroflavonol 4-reductase (DFR), anthocyanidin synthase (ANS), and anthocyanin reductase (ANR) to form a series of derivatives. However, the anthocyanins produced in the cytoplasm are unstable and require further glycosylation modification by UDP glucosyltransferase (UGT; named BZ1 in cassava), which is ultimately transferred to vacuoles and stored together with glutathione S-transferase (GST). (Winkel-Shirley, 2001; Jaakola, 2013). Nevertheless, the biosynthesis mechanism of anthocyanins in cassava is poorly understood. There are unknown reasons why purple leaves have more anthocyanins than green leaves. Therefore, it is necessary to explore the specific mechanism of anthocyanin biosynthesis in cassava.

Ultra-Performance Liquid Chromatography-tandem/Electron Spray Ionization-Quadrupole TRAP-Mass Spectrometry/Mass Spectrometry (UPLC/ESI-Q TRAP-MS/MS) is widely used in the identification and analysis of plant metabolite and has the advantages of high sensitivity and throughput, fast separation, and wide coverage. This technology has been widely applied to analyze the metabolites in tomatoes, strawberries, and asparagus (Paolo et al., 2018; Zhu et al., 2018). In recent years, an integrative analysis of metabolomics and transcriptomics was used to analyze the correlation between metabolites and genes expression in plants (Jiang et al., 2020; Liu et al., 2020; Xu et al., 2021; Zou et al., 2021).

In this study, an integrative analysis of metabolomics and transcriptomics was performed on two cassava varieties, SC9 with green leaves and PL with purple leaves, from the National Cassava Germplasm Repository (NCGR) in China to investigate the mechanism of anthocyanin biosynthesis in cassava leaves. Significantly different metabolites (SDMs) and differentially expressed genes (DEGs) involved in anthocyanin biosynthesis were analyzed. The virus-induced gene silencing (VIGS)-induced phenotypes were also characterized. Our study provides the candidate genes and further explains the anthocyanin biosynthetic pathway in cassava. It also laid the foundation to develop functional substances and promote the cassava value chain.

Materials and methods

Plant material and sampling

Two cassava varieties, SC9 (Green leaves) and PL (Purple leaves), provided by the National Cassava Germplasm Repository were planted in the field (19°30’33.13”N, 109°30’19.34”E) located at Danzhou City, Hainan Province, China, in February. The young leaves of SC9 and PL were collected 4 months after planting. The leaves of SC9 and PL were frozen in liquid nitrogen and stored at -80°C for RNA and metabolite extraction. The cassava leaf sample was ground into powder (30 Hz, 1.5 min) in liquid nitrogen. For RNA extraction using the RNA extraction kit (Tiangen), 0.1g powder was weighted, and 0.05 g powder was weighted for metabolites extraction with 0.5 mL methanol/water/hydrochloric acid (500:500:1, V/V/V). The metabolite extract was vortexed for 5 min, ultrasound was performed for 5 min and was centrifuged for 3 min at 12,000 rpm. All the steps were done at a temperature of 4°C. The residue was re-extracted by repeating the above steps with the same conditions. The supernatants were collected and filtrated through a 0.22 μm membrane filter before LC-MS/MS analysis. Three biological replicates were performed for all experiments.

Library preparation and transcriptome sequencing

RNA extraction and qualification were performed according to the instruction by An et al. (An et al., 2022). NEBNext® UltraTM RNA Library Prep Kit (New England Biolabs, USA) was used to prepare sequencing libraries by Illumina® following the instruction and index codes. The library fragments were purified with the AMPure XP system (Beckman Coulter, Beverly, USA) to select cDNA fragments ranging from 250 to 300 bp in length. Purified (AMPure XP system) PCR products and library quality were assessed on the Agilent Bioanalyzer 2100 system.

The clustering of the index-coded samples was performed on a cBot Cluster Generation System using TruSeq PE Cluster Kit v3-cBot-HS following the instructions. Then, the library preparations were sequenced on an Illumina Hiseq platform and 125 bp/150 bp paired-end reads were generated. In this study, the coding sequences (CDSs) predicted from the Phytozome (https://phytozome-next.jgi.doe.gov/) were annotated through KEGG (https://wwwkegg.jp/kegg), GO (https://geneontology.org), NR (https://www.ncbi.nlm.nih.gov/), Swiss-Prot (https://www.uniprot.org/), trEMBL (https://www.uniprot.org/), and KOG (https://www.ncbi.nlm.nih.gov/research/cog-project/) databases. PL was used as an indicator to measure transcript or gene expression levels. DESeq2 was used to identify DEGs with absolute log₂fold change of equal or greater than 1 and false discovery rate (FDR) less than 0.05 (Varet et al., 2016). GO and KEGG enrichment analyses of DEGs were further implemented by employing the cluster Profiler R package (version 4.1.3). After 4 months of planting, cassava leaves were collected and three biological replicates were used for all experiments.

The identification and quantification of anthocyanin content in cassava leaves

The identification and quantification of anthocyanin metabolites were performed at the Wuhan MetWare Biotechnology Co., Ltd. (Wuhan, China). The extraction of samples was analyzed using a UPLC-ESI-MS/MS system (UPLC, ExionLC™ AD, https://sciex.com.cn/; MS, Applied Biosystems 6500 Triple Quadrupole, https://sciex.com.cn/). The analytical conditions of UPLC were as follows: column is ACQUITY BEH C18 (Waters, 1.7 µm, 2.1 mm*100 mm); the solvent system is water (with 0.1% formic acid) and methanol (with 0.1% formic acid) with four gradient programs 95:5 V/V, 50:50 V/V, 5:95 V/V, and 95:5 V/V hold for 10 min, 6 min, 14 min, 14 min, and 2 min, respectively; flow rate is 0.35 mL/min; temperature is 40°C; injection volume is 2 μL. Analyst 1.6.3 software (Sciex) was used for data acquisitions. The MetWare MWDB database was used to identify the metabolites through different parameters such as retention time, m/z value, and fragmentation mode. Differentially accumulated metabolites of filtering conditions were as follows: p-value < 0.05, |log₂ (FC)| is equal or greater than 1, and variable importance in projection (VIP) is equal or greater than 1. To study the accumulation of specific metabolites, the R package (version 4.1.3) was used to perform principal component analysis (PCA) and orthogonal partial least squares-discriminant analysis (OPLS-DA). After 4 months of planting, cassava leaves were collected and three biological replicates were performed for each experiment.

qRT-PCR involved in anthocyanin biosynthetic pathway

The total cDNA of cassava leaves were synthesized using the HiScriptIII 1^st Strand cDNA Synthesis Kit (Vazyme, Code No. Q312-01). The cDNA diluted to 10-fold by MillQ was used as a template to measure the genes’ relative expression level. The specific primers for genes involved in the anthocyanin biosynthetic pathway and the cassava internal control gene (actin) are listed in Table S1. qRT-PCR was conducted with a real-time PCR ABI 6800 system using the ChamQ SYBR qPCR Master Mix (Vazyme, Code No. Q311-02). The comparative CT method (2^-ΔΔCT) was used to quantify the genes’ relative expression (Schmittgen and Livak, 2008). Three biological replicates were performed for each experiment.

Combined analysis of the transcriptome and metabolome

According to the results of the analysis of the combination of metabolites and transcriptome differential gene detected in this study, the differential genes and all metabolites were mapped into the KEGG pathway map to better understand the relationship between genes and metabolites. Pearson correlation coefficients were calculated to integrate transcriptome and metabolome data. For the joint analysis between the transcriptome and metabolome, the screening criterion was a Pearson correlation coefficient great than 0.8.

VIGS in cassava

For vector construction, 300 bp MeANR was cloned into pCsCMV-NC (negative control) as described by Tuo et al., 2021. The experiment materials used in VIGS were the stems of SC9 which were planted in nutrient-rich soil in an air-conditioned growth room (16-h photoperiod; 2000 lux light intensity) for 20 days. The inoculation method was followed as reported by An et al., 2022. Transformed Agrobacterium tumefaciens was precipitated, resuspended, and injected into the cassava leaves and allowed to grow in the greenhouse for 28 dpi. After that, the expression level of the targeted genes in VIGS-MeANR silenced lines was detected. Three MeANR silencing plants and SC9 were tested. All the primers involved are shown in Table S1.

Statistical analysis

Statistical analysis was performed using Excel 2010 software (Microsoft Office, USA). Data are presented as means ± standard deviations (SD). The levels of statistical significance were analyzed by the least significant difference (Kruskal-Wallis test, p < 0.05).

Results

The total anthocyanin content in the leaves of different cassava varieties

An important pigment, anthocyanin, plays a vital role in the formation of color in cassava leaves. In this study, two different varieties, SC9 and PL, were selected as materials. The content of the total anthocyanins in SC9 and PL was extracted using the methanol-HCl method. The results indicated that the extraction color of SC9 and PL leaves displayed yellow-green and purple-red colors, respectively. The content of total anthocyanins of PL leaves was 413.83 ng/g of fresh weight (FW) which was significantly higher than the total anthocyanins content of SC9 (12.39 ng/g FW) (Figure 1) (Kruskal-Wallis test; p < 0.05).

FIGURE 1

Figure 1 The phenotype and the total anthocyanin content between SC9 and PL are characterized. (A) The color of the cassava leaf and extraction between SC9 and PL are shown. Bar: 1 cm. (B) The total content of anthocyanin between SC9 and PL was measured. Error bars indicate the standard deviation (SD). Significant differences are indicated by an asterisk (*p < 0.05).

The identification of flavonoids and anthocyanin metabolites in cassava leaves

To investigate the differential flavonoids and anthocyanin metabolites between SC9 and PL, three replications of each sample were analyzed by UPLC/ESI-Q TRAP-MS/MS analysis. A total of 107 flavonoid metabolites were identified in different cassava leaves. PCA analysis indicated that the metabolites of SC9 and PL could be clearly distinguished (Figure S1). Among these metabolites of flavonoids, rutin is the most abundant in cassava leaves, followed by kaempferol-3-O-rutinoside and quercetin-3-O-glucoside, which account for 97% of flavonoids in cassava leaves; however, this was not significantly different between SC9 and PL.

SDMs were screened using the VIP value > 1 in the OPLS-DA model and FC ≥ 2 or FC ≤ 0.5 as the identification criterion. A total of 25 SDMs in PL were upregulated compared to SC9. These SDMs were generally classified into eight categories, including cyanidins, delphinidins, flavonoids, malvidins, peonidins, petunidins, procyanidins, and pelargonidins. Interestingly, six categories of the SDMs belonging to the anthocyanins were highly accumulated in PL compared with SC9 (Figure 2; Table 1).

FIGURE 2

Figure 2 The heatmap of different classes of anthocyanin. The horizontal axis stands for sample processing information, and the vertical axis stands for metabolite information. The cluster tree on the left is the metabolite clustering tree. The scale is the expression quantity obtained after standardized processing; the red color indicates the higher expression quantity and the grey indicates that the expression quantity is N/A. The group includes SC9 and PL.

TABLE 1

Table 1 The contents of differential anthocyanidins between SC9 and PL.

Differentially expressed genes in cassava leaves

To further understand the anthocyanin biosynthesis in different cassava leaves, transcriptome sequencing was performed. A total of 277,823,272 clean reads were produced by the cassava leaves. The average value of Q20 was 97.46% and the Q30 value was over 92%; the average content of GC was 44.83%, and the sequencing error rates were 0.03%, which was lower than 0.1%. These data indicated that the quality of the sequencing data was good enough for further analysis (Table S2).

With the filter criteria |Log₂FC| ≥1 and FDR < 0.05, 4,846 DEGs were identified in the cassava leaves, including 2,349 up-regulated genes and 2,497 down-regulated genes in PL compared to SC9 (Table S4). The volcano plots were performed to show the overall distribution of gene expression levels (Figure 3A). In order to characterize the function of DEGs, GO and KEGG enrichment analyses were performed. The results showed that these DEGs were associated with biological processes, cellular components, and molecular functions (Figure 3B). Combined with the KEGG annotation database, a number of DEGs were enriched in many metabolic pathways, including secondary metabolites biosynthesis and flavonoid biosynthetic pathways (Figure 3C). The homologous genes associated with anthocyanin biosynthesis in PL were calculated and compared with SC9 (Table S3).

FIGURE 3

Figure 3 Screening and functional enrichment analysis of DEGs between SC9 and PL. (A) The overdistribution of gene expression levels and fold differences between SC9 and PL. The x-axis shows the change in gene multiples, and the y-axis shows the significance level of differentially expressed genes. The red dots stand for up-regulated differential genes, the green dots stand for down-regulated differential genes, and the blue dots stand for non-differentially expressed genes. (B) KEGG pathway enrichment of DEGs. The x-axis shows the Rich factor, and the y-axis shows the KEGG pathway. The larger the Rich factor the greater the degree of enrichment. The larger the dot the more the number of enriched differential genes in those pathways. The redder the dot the more significant the enrichment. (C) GO enrichment analysis of the 4,846 DEGs in three functional groups (Biological process, Cellular component, and Molecular function) between SC9 and PL. The x-axis shows the secondary entries, and the y-axis shows the number of entry differential genes.

The combined analysis of the transcriptome and metabolome

To further elucidate the relationship between metabolite accumulation and gene expression in cassava leaves, the analysis in combination with transcriptomic and metabolomic data was performed and illustrated. In anthocyanin biosynthesis, 25 SDMs and 10 DEGs were screened and analyzed using correlation analysis after the criteria |Correlation| ≥ 0.8. These results demonstrated that the anthocyanin metabolite changes were significantly associated with the gene expression in the anthocyanin biosynthetic pathway (Figure 4). The SDMs such as delphinidin-3-O-sophoroside, malvidin-3-O-sambubioside, and peonidin-3,5-O-diglucoside had a higher correlation with Manes.16G016400 (MeANR), while malvidin-3-O-glucoside, delphinidin, and peonidin-3-O-rutinoside had a higher correlation with Manes.01G070200 (MeANS) (Figure 5). These results indicated that these SDMs and DEGs may participate in anthocyanin biosynthesis.

FIGURE 4

Figure 4 SDMs and DEGs were altered in the anthocyanin biosynthetic pathway between SC9 and PL. Green and red colors indicate the down-regulation and up-regulation of DEGs, respectively.

FIGURE 5

Figure 5 Network analysis of SDMs and DEGs in the anthocyanin biosynthetic pathway between SC9 and PL. Green squares indicate metabolites in the anthocyanin biosynthetic pathway, and red circles indicate genes expressed in the anthocyanin biosynthetic pathway. These SDMs and DEGs were analyzed by the criteria |Pearson’s correlation coefficient| ≥ 0.8.

Alterations of anthocyanin biosynthesis in PL compared to those in SC9

The relative expression levels in 10 genes of the anthocyanin biosynthesis including MeCHS, MeCHI, MeN3D, MeF3’H, MeF3’5’H, MeDFR, MeANS, and MeBZ1, except MeANR, were significantly up-regulated in PL than in SC9 leading to high anthocyanin accumulation. To verify the credibility of the transcriptome information, we further selected ten DEGs to validate the sequencing results. The qRT-PCR results were consistent with the RNA-Seq results, suggesting that those DEGs may be involved in anthocyanin biosynthesis (Figure 6). The homologous genes associated with anthocyanin biosynthesis were also confirmed by qRT-PCR analysis (Figure S2).

FIGURE 6

Figure 6 qRT-PCR analysis of the expression of anthocyanin biosynthetic pathway-associated genes. Leaves from cassava plants were used for each RNA extraction (three RNA extractions per plant; n = 3). Data indicate means ± SD. Asterisks indicate a significant reduction or increase compared to the control plants. Significant differences are indicated by an asterisk (*p < 0.05).

Effects of transcription factors on metabolite accumulation between SC9 and PL leaves

Transcription factors (TFs) regulate the key enzyme expression and participate in secondary metabolism biosynthesis processes in plants. More than 15 TF families including 210 differentially expressed TFs were found in both SC9 and PL leaves in the study (Table S5). The most abundant difference expressed transcription factors (DETFs) were MYB-dominant protein (23.24%, 34 were up-regulated and 9 were down-regulated), basic helix-loop-helix protein (bHLH) (15.68%, 15 were upregulated and 14 were downregulated), NAC domain-containing proteins (12.97%, 16 were up-regulated and 8 were down-regulated), followed by AP2/ERF (9.19%, 5 were up-regulated and 12 were down-regulated), WRKY (6.49%, 10 were up-regulated and 2 were down-regulated), Zinc finger (6.49%, 6 were up-regulated and 6 were down-regulated), GARS (5.95%, 5 were up-regulated and 6 were down-regulated), and other TFs (Table 2). These TFs may participate in anthocyanin and flavonoid biosynthesis or other secondary metabolism biosynthesis in cassava leaves.

TABLE 2

Table 2 Transcription factors found between SC9 and PL leaves.

MeANR is necessary for anthocyanin biosynthesis and color formation in cassava leaves

To investigate whether MeANR was a key enzyme negatively regulating anthocyanin biosynthesis in cassava, VIGS-MeANR silenced plants were synthesized and characterized in SC9 after 28 dpi. VIGS-MeANR silenced plant presented the altered phenotype of cassava leaves, partially from green to purple color (Figure 7A). The transcript level of MeANR and the total anthocyanin content were also measured in VIGS-MeANR silenced plant. The results showed that the expression level of MeANR was reduced by approximately 60% compared with the control plant (Figure 7B). The content of the total anthocyanin in the VIGS-MeANR silenced plant was 52.46 ng/g, which was higher than that in the control plant (Figure 7C). Furthermore, the metabolites in anthocyanin biosynthesis were analyzed, and it was found that these metabolites were highly accumulated in VIGS-MeANR silenced plants (Table 3). These results indicated that MeANR is a key factor that negatively regulates anthocyanin biosynthesis and color formation in cassava leaves.

FIGURE 7

Figure 7 Silencing of MeANR in SC9 and its effects on anthocyanin biosynthesis. (A) Phenotypes of cassava SC9 leaves and plants displayed by VIGS of MeANR. pCsCMV-NC was used as control. The upper photos show a young leaf while the lower photo shows the plant. Bars: 1 cm. (B) qRT-PCR analysis on the expression level of MeANR in VIGS-MeANR plants. Leaves from cassava plants were used for each RNA extraction (three RNA extractions per silencing plant; n = 3). Error bars indicate the standard deviation and significant differences are indicated by asterisk (*p < 0.05). (C) The total anthocyanin content was measured in VIGS-MeANR plants. Error bars indicate the standard deviation and significant differences are indicated by an asterisk (*p< 0.05).

TABLE 3

Table 3 Significantly different metabolites in the anthocyanin biosynthetic pathway of the VIGS-MeANR plant in comparison with the control plant.

Discussion

Anthocyanins are secondary metabolites in plants and have important biological functions in antioxidation, defense regulation, and protection against environmental stresses (Shen et al., 2018; Jiang et al., 2020). In the study, we found that the anthocyanin contents in PL were higher than that in SC9 (Figure 1) which was consistent with the finding of the high anthocyanins content in purple leaves of tea cultivar “Zijuan” (Li et al., 2017) and “Baitang” tea (Rothenberg et al., 2019). Interestingly, a total of 20 significantly different anthocyanin metabolites were identified in cassava leaves by UPLC/ESI-QTRAP-MS/MS, which were highly accumulated. Among these anthocyanins, the content of three anthocyanins, i.e., cyanidin-3-(6-O-p-caffeoyl)-glucoside, cyanidin-3-O-rutinoside, and delphinidin-3-O-rutinoside was the most abundant in PL. The cyanidin-3-O-sophoroside content was 668.55-fold higher in PL compared with SC9 and the content of delphinidin-3-O-rutinoside was highest in PL with 295.935 ng/g of fresh weight (Table 1). These results suggest that delphinidin-3-O-rutinoside is the major anthocyanin in PL, which may result in the purple pigmentation of cassava leaves.

Anthocyanin biosynthesis was regulated by a series of enzymes. In the study, differentially expressed genes involved in anthocyanin biosynthesis were identified by transcriptome analysis in cassava leaves (Table S3). In the anthocyanin biosynthesis, the expression levels of MeCHS, MeCHI, MeF3’H, MeF3’5’H, MeDFR, MeANS, and MeBZ1 were up-regulated in PL than in SC9 (Figure 6), except MeANR. ANS is a key enzyme of the anthocyanin pathway that plays an important role in anthocyanidin formation (Ben-Simhon et al., 2015). In our study, the expression level of MeANS was the highest among these genes in PL vs. SC9. Some studies also showed that overexpression of SmANS or StANS increased anthocyanin content in Salvia miltiorrhiza and Solanum tuberosum, respectively (Li et al., 2019a; Zhang et al., 2020). ANR genes (MrANR1, 2 and MiANR1-1, 1-2, and 1-3) are the critical factors in proanthocyanidins (PA) biosynthesis, and there is competition between the anthocyanin and PA biosynthetic pathways (Tan et al., 2018; Li et al., 2019b). Simultaneously, the expression level of MeANR was significantly down-regulated in PL. We also used the VIGS system to rapidly analyze the function of MeANR. VIGS-MeANR silenced plant regulated the altered phenotype of cassava leaves, partially from green to purple color (Figure 7A), resulting in a significant increase of the total anthocyanin content (Figures 7B, C). The metabolites of anthocyanin were also highly accumulated in VIGS-MeANR silenced plants (Table 3). The results demonstrated that MeANR is a negative regulator of anthocyanin biosynthesis in cassava. In grape berries, VIGS-VvANR silencing increases the expression of VvANS and promotes anthocyanin accumulation (Yang et al., 2022). In short, the high expression level of MeANS and the low expression level of MeANR may promote the accumulation of anthocyanin, leading to the formation of purple leaves. At present, most of the cassava varieties have green leaves. In order to develop the varieties with higher anthocyanin content in leaves, silencing or knock-out MeANR could be done in cassava varieties with green leaves so that the varieties with higher anthocyanin leaves might be available and, thereby, fulfill the requirements for food and feeding in the future.

TFs play crucial roles as regulators in anthocyanin biosynthesis. MYB, bHLH, WRKY, WD40, zinc finger, and MADs proteins have been reported for regulating anthocyanin biosynthesis (Terrier et al., 2008; Lloyd et al., 2017; Trainin et al., 2021). To date, the regulation of anthocyanin biosynthesis was known to be not only controlled by a single TF but also regulated by TF/TFs combined with some proteins via an anthocyanin biosynthesis network. A specific example is the MYB-bHLH-WD40 (MBW) ternary complex of transcription factors involving anthocyanin accumulation and inhibition between Arabidopsis, Maize, and Cauliflower (Hichri et al., 2011; Petroni and Tonelli, 2011; Chiu and Li, 2012; Shi et al., 2018). In kiwifruit, an MYB/bHLH complex regulates tissue-specific anthocyanin biosynthesis in the inner pericarp of red-centered kiwifruit Actinidia chinensis cv. Hongyang (Wang et al., 2019). Our study found 210 TFs including MYB, bHLH, NACs, AP2/ERFs, WRKY, GARS, Zinc-finger, MADs, HSF, B3, Trihelix, bZIP, HD-ZIP, SBP, and GATA by analyzing the transcriptome data, which showed the significantly different expression levels in PL vs. SC9 (Table S5). We speculated that the anthocyanin biosynthesis may be regulated by these transcription factors in cassava and their biological functions need to be further studied.

Data availability statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://www.ncbi.nlm.nih.gov/PRJNA897109.

Author contributions

This study was designed by JC and SC. The experiment was performed by XL and FA. Data analysis was finished by JX, WZ, ZW, WO, and KL. Author XL wrote the paper. JC and SC modified the manuscript. All authors reviewed and approved the manuscript.

Funding

This work was financially supported by the National Natural Science Foundation of China Projects (No.32101811 and No. 31871687), the National Key Research and Development Program of China (No. 2019YFD1000500), and the Central Public-interest Scientific Institution Basal Research Fund for Chinese Academy of Tropical Agricultural Science (1630032022007).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpls.2023.1181257/full#supplementary-material

References

An, F. F., Chen, T., Li, Q. X., Qiao, J. J., Zhang, Z. W., Carvalho, L. J., et al. (2019). Protein cross-interactions for efficient photosynthesis in the cassava cultivar SC205 relative to its wild species. J. Agric. Food Chem. 67. doi: 10.3389/fpls.2022.901128

CrossRef Full Text | Google Scholar

An, F. F., Xiao, X. H., Chen, T., Xue, J. J., Luo, X. Q., OU, W. J., et al. (2022). Systematic analysis of bHLH transcription factors in cassava uncovers their roles in postharvest physiological deterioration and cyanogenic glycosides biosynthesis. Front. Plant Sci. 13. doi: 10.3389/fpls.2022.901128

CrossRef Full Text | Google Scholar

Bai, L., Li, X., He, L., Zheng, Y., Lu, H., Li, J., et al. (2019). Antidiabetic potential of flavonoids from traditional chinese medicine: a review. Am. J. Chin. Med. 47 (5), 933–957. doi: 10.1142/S0192415X19500496

PubMed Abstract | CrossRef Full Text | Google Scholar

Ben-Simhon, Z., Judeinstein, S., Trainin, T., Harel-Beja, R., Bar-Ya’akov, I., Borochov-Neori, H., et al. (2015). A”White” anthocyanin-less pomegranate (Punica granatum l.) caused by an insertion in the coding region of the leucoanthocyanidin dioxygenase (LDOX, ANS) gene. PloS One 10 (11), e0142777. doi: 10.1371/journal.pone.0142777

PubMed Abstract | CrossRef Full Text | Google Scholar

Bosse, M. A., Silva, M. B., Oliveira, N. G. M. O., Araujo, M. A. A., Rodrigues, C., Azevedo, J. P., et al. (2021). Physiological impact of flavonoids on nodulation and ureide metabolism in legume plants. Plant Physiol. Biochem. 166, 512–521. doi: 10.1016/j.plaphy.2021.06.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Cai, J., Zhang, J., Ding, Y., Yu, S., Lin, H. X., Yuan, Z. Q., et al. (2021). Different fertilizers applied alter fungal community structure in rhizospheric soil of cassava (Manihot esculenta crantz) and increase crop yield. Front. Microbiol. 12. doi: 10.3389/fmicb.2021.663781

CrossRef Full Text | Google Scholar

Calis, Z., Mogulkoc, R., Baltaci, A. K. (2020). The roles of flavonols/flavonoids in neurodegeneration and neuroinflammation. Mini Rev. Med. Chem. 20 (15), 1475–1488. doi: 10.2174/1389557519666190617150051

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, W., Gao, Y. Y., Xie, W. B., Gong, L., Lu, K., Wang, W. S., et al. (2014). Genome-wide association analyses provide genetic and biochemical insights into natural variation in rice metabolism. Nat. Genet. 46 (7), 714–721. doi: 10.1038/ng.3007

PubMed Abstract | CrossRef Full Text | Google Scholar

Chiu, L. W., Li, L. (2012). Characterization of the regulatory network of BoMYB2 in controlling anthocyanin biosynthesis in purple cauliflower. Planta 236 (4), 1153–1164. doi: 10.1007/s00425-012-1665-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Ciumărnean, L., Milaciu, M., Runcan, O., Vesa, S., Răchișan, A., Negrean, V., et al. (2020). The effects of flavonoids in cardiovascular diseases. Molecules 25 (18), 4320. doi: 10.3390/molecules25184320

PubMed Abstract | CrossRef Full Text | Google Scholar

Ferreyra, M. F., RiuS, S. P., Casati, P. (2012). Flavonoids: biosynthesis, biological functions, and biotechnological applications. Front. Plant Sci. 3. doi: 10.3389/fpls.2012.00222

CrossRef Full Text | Google Scholar

Forkmann, G. (2010). Flavonoids as flower pigments: the formation of the natural spectrum and its extension by genetic engineering. Plant Breeding. 106, 1–26. doi: 10.1111/j.1439-0523.1991.tb00474.x

CrossRef Full Text | Google Scholar

Fu, L. L., Ding, Z. H., Tie, W. W., Yang, J. H., Yan, Y., Hu, W. (2022). Integrated metabolomic and transcriptomic analyses reveal novel insights of anthocyanin biosynthesis on color formation in cassava tuberous roots. Front. Nutr. 9. doi: 10.3389/fnut.2022.842693

CrossRef Full Text | Google Scholar

Hichri, I., Barrieu, F., Bogs, J., Kappel, C., Delrot, S., Lauvergeat, V. (2011). Recent advances in the transcriptional regulation of the flavonoid biosynthetic pathway. J. Exp. Bot. 62 (8), 2465–2483. doi: 10.1093/jxb/erq442

PubMed Abstract | CrossRef Full Text | Google Scholar

Jaakola, L. (2013). New insights into the regulation of anthocyanin biosynthesis in fruits. Trends Plant Sci. 18, 477–483. doi: 10.1016/j.tplants.2013.06.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Jiang, T., Zhang, M. D., Wen, C. X., Xie, X. L., Tian, W., Wen, S. Q., et al. (2020). Integrated metabolomic and transcriptomic analysis of the anthocyanin regulatory networks in Salvia miltiorrhiza bge. flowers. BMC Plant Biol. 20 (1), 349. doi: 10.1186/s12870-020-02553-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Kumba, Y. K. (2012). Genetic characterization of exotic and landraces of cassava in Ghana (Peking: University of Science and Technology).

Google Scholar

Li, H. Y., Liu, J. L., Pei, T. L., Bai, Z. Q., Han, R. L., Liang, Z. S. (2019a). Overexpression of SmANS enhances anthocyanin accumulation and alters phenolic acids content in Salvia miltiorrhiza and Salvia miltiorrhiza bge f. alba plantlets. Int. J. Mol. Sci. 20 (9), 2225. doi: 10.3390/ijms20092225

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, J., Lv, X. J., Wang, L. X., Qiu, Z. M., Song, X. M., Lin, J. K., et al. (2017). Transcriptome analysis reveals the accumulation mechanism of anthocyanins in ‘Zijuan’ tea (Camellia sinensis var. asssamica (Masters) kitamura) leaves. Plant Growth Regul. 81, 51–61. doi: 10.1007/s10725-016-0183-x

CrossRef Full Text | Google Scholar

Li, H., Tian, J., Yao, Y. Y., Zhang, J., Song, T. T., Li, K. T., et al. (2019b). Identification of leucoanthocyanidin reductase and anthocyanidin reductase genes involved in proanthocyanidin biosynthesis in Malus crabapple plants. Plant Physiol. Biochem. 139, 141–151. doi: 10.1016/j.plaphy.2019.03.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, Y., Lv, J., Liu, Z., Wang, J., Yang, B., Chen, W., et al. (2020). Integrative analysis of metabolome and transcriptome reveals the mechanism of color formation in pepper fruit (Capsicum annuum l.). Food Chem. 306, 125629. doi: 10.1016/j.foodchem.2019.125629

PubMed Abstract | CrossRef Full Text | Google Scholar

Lloyd, A., Brockman, A., Aguirre, L., Campbell, A., Bean, A., Cantero, A., et al. (2017). Advances in the MYB-bHLH-WD repeat (MBW) pigment regulatory model: addition of a WRKY factor and co-option of an anthocyanin MYB for betalain regulation. Plant Cell Physiol. 58 (9), 1431–1441. doi: 10.1093/pcp/pcx075

PubMed Abstract | CrossRef Full Text | Google Scholar

Omar, N. F., Hassan, S. A., Yusoff, U. K., Abdullah, N. A., Wahab, P. E., Sinniah, U. (2012). Phenolics, flavonoids, antioxidant activity and cyanogenic glycosides of organic and mineral-base fertilized cassava tubers. Molecules 17 (3), 2378–2387. doi: 10.3390/molecules17032378

PubMed Abstract | CrossRef Full Text | Google Scholar

Paolo, B., Saverio, O., Mirko, M., Matteo, B., Lara, G., Azeddine, S. A. (2018). Gene expression and metabolite accumulation during strawberry (Fragaria×ananassa) fruit development and ripening. Planta 248 (5), 1143–1157. doi: 10.1007/s00425-018-2962-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Park, M. Y., Kim, Y., Ha, S. E., Kim, H. H., Bhosale, P. B., Abusaliya, A., et al. (2022). Function and application of flavonoids in the breast cancer. Int. J. Mol. Sci. 23 (14), 7732. doi: 10.3390/ijms23147732

PubMed Abstract | CrossRef Full Text | Google Scholar

Petroni, K., Tonelli, C. (2011). Recent advances on the regulation of anthocyanin synthesis in reproductive organs. Plant Sci. 181 (3), 219–229. doi: 10.1016/j

PubMed Abstract | CrossRef Full Text | Google Scholar

Quattrocchio, F., Baudry, A., Lepiniec, L., Grotewold, E. (2006). “The regulation of flavonoid biosynthesis,” in The science of flavonoids. Ed. Grotewold, E. (New York, NY: Press), 97–122.

Google Scholar

Rothenberg, D. O. N., Yang, H. J., Chen, M. B., Zhang, W. T., Zhang, L. Y. (2019). Metabolome and transcriptome sequencing analysis reveals anthocyanin metabolism in pink flowers of anthocyanin-rich tea (Camellia sinensis). Molecules 24 (6), 1064. doi: 10.3390/molecules24061064

PubMed Abstract | CrossRef Full Text | Google Scholar

Schmittgen, T. D., Livak, K. J. (2008). Analyzing real-time PCR data by the comparative CT method. Nat. Protoc. 3 (6), 1101–1108. doi: 10.1038/nprot.2008.73

PubMed Abstract | CrossRef Full Text | Google Scholar

Sharma, A., Shahzad, B., Rehman, A., Bhardwaj, R., Landi, M., Zheng, B. (2019). Response of phenylpropanoid pathway and the role of polyphenols in plants under abiotic stress. Molecules 24 (13), 2452. doi: 10.3390/molecules24132452

PubMed Abstract | CrossRef Full Text | Google Scholar

Shen, J. Z., Zou, Z. W., Zhang, X. Z., Zhou, L., Wang, Y. H., Fang, W. P., et al. (2018). Metabolic analyses reveal different mechanisms of leaf color change in two purple-leaf tea plant (Camellia sinensis l.) cultivars. Hortic. Res. 5, 7. doi: 10.1038/s41438-017-0010-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Shi, H., Liu, G., Wei, Y., Chan, Z. (2018). The zinc-finger transcription factor ZAT6 is essential for hydrogen peroxide induction of anthocyanin synthesis in Arabidopsis. Plant Mol. Biol. 97 (1-2), 165–176. doi: 10.1007/s11103-018-0730-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Tan, L., Wang, M., Kang, Y. F., Azeem, F., Zhou, Z. X., Tuo, D. C., et al. (2018). Biochemical and functional characterization of anthocyanidin reductase (ANR) from Mangifera indica l. Molecules 23 (11), 2876. doi: 10.3390/molecules23112876

PubMed Abstract | CrossRef Full Text | Google Scholar

Terrier, N., Torregrosa, L., Ageorges, A., Vialet, S., Verries, C., Cheynier, V., et al. (2008). Ectopic expression of VvMybPA2 promotes proanthocyanidin biosynthesis in grapevine and suggests additional targets in the pathway. Plant Physiol. 149 (2), 1028–1041. doi: 10.1104/pp.108.131862

PubMed Abstract | CrossRef Full Text | Google Scholar

Trainin, T., Harel-Beja, R., Bar-Ya’akov, I., Ben-Simhon, Z., Yahalomi, R., Borochov-Neori, H., et al. (2021). Fine mapping of the “black” peel color in pomegranate (Punica granatum l.) strongly suggests that a mutation in the anthocyanidin reductase (ANR) gene is responsible for the trait. Front. Plant Sci. 12. doi: 10.3389/fpls.2021.642019

PubMed Abstract | CrossRef Full Text | Google Scholar

Tuo, D. C., Zhou, P., Yan, P., Cui, H. G., Liu, Y., Yang, X. K., et al. (2021). A cassava common mosaic virus vector for virus-induced gene silencing in cassava. Plant Methods 17 (1), 74. doi: 10.1186/s13007-021-00775-w

PubMed Abstract | CrossRef Full Text | Google Scholar

Varet, H., Brillet-Gueguen, L., Coppee, J. Y., Dilles, M. A. (2016). SARTools: a DESeq2- and EdgeR-based r pipeline for comprehensive differential analysis of RNA-seq data. PloS One 11 (6), e0157022. doi: 10.1371/journal.pone.0157022

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, L. H., Tang, W., Hu, Y. H., Zhang, Y. B., Sun, J. Q., Guo, X. H., et al. (2019). A MYB/bHLH complex regulates tissue specific anthocyanin biosynthesis in the inner pericarp of red-centered kiwifruit Actinidia chinensis cv. hongyang. Plant J. 99 (2), 359–378. doi: 10.1111/tpj.14330

PubMed Abstract | CrossRef Full Text | Google Scholar

Winkel-Shirley, B. (2001). Flavonoid biosynthesis. a colorful model for genetics, biochemistry, cell biology, and biotechnology. Plant Physiol. 126 (2), 485–493. doi: 10.1104/pp.126.2.485

PubMed Abstract | CrossRef Full Text | Google Scholar

Xiao, L., Cao, S., Shang, X. H., Xie, X. Y., Zeng, W. D., Lu, L. Y., et al. (2021). Metabolomic and transcriptomic profiling reveals distinct nutritional properties of cassavas with different flesh colors. Food Chem. (Oxf). 2, 100016. doi: 10.1016/j.fochms.2021.100016

PubMed Abstract | CrossRef Full Text | Google Scholar

Xu, X. Y., Lu, X. N., Tang, Z. L., Zhang, X. N., Lei, F. J., Hou, L. P., et al. (2021). Combined analysis of carotenoid metabolites and the transcriptome to reveal the molecular mechanism underlying fruit colouration in zucchini (Cucurbita pepo l.). Food Chem. (Oxf). 2, 100021. doi: 10.1016/j.fochms.2021.100021

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, B., Wei, Y., Liang, C., Niu, T., Zhang, P. F., Wen, P. F. (2022). VvANR silencing promotes expression of VvANS and accumulation of anthocyanin in grape berries. Protoplasma 259 (3), 743–753. doi: 10.1007/s00709-021-01698-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Zaynab, M., Fatima, M., Abbas, S., Sharif, Y., Umair, M., Zafar, M. H., et al. (2018). Role of secondary metabolites in plant defense against pathogens. Microb. Pathog. 124, 198–202. doi: 10.1016/j.micpath.2018.08.034

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, H. L., Zhao, X. J., Zhang, J. P., Yang, B., Yu, Y. H., Liu, T. F., et al. (2020). Functional analysis of an anthocyanin synthase gene StANS in potato. Sci. Hortic. 272, 109569. doi: 10.1016/j.scienta.2020.109569

CrossRef Full Text | Google Scholar

Zhu, G. T., Wang, S. C., Huang, Z. J., Zhang, S. B., Lia, Q. G., Zhang, C. Z., et al. (2018). Rewiring of the fruit metabolome in tomato breeding. Cell 172 (1-2), 249–261. doi: 10.1016/j.cell.2017.12.019

PubMed Abstract | CrossRef Full Text | Google Scholar

Zou, Q. J., Wang, T., Guo, Q. S., Yang, F., Chen, J. M., Zhang, W. Y. (2021). Combined metabolomic and transcriptomic analysis reveals redirection of the phenylpropanoid metabolic flux in different colored medicinal Chrysanthemum morifolium. Ind. Crop Prod. 164, 113343. doi: 10.1016/j.indcrop.2021.113343

CrossRef Full Text | Google Scholar