Inflorescences of two H. brasiliensis genotypes—Yunyan73477 (high-embryogenic, HE) and Reken628 (low-embryogenic, LE)—were collected during the peak flowering period in March. Immature male flowers (1.3–1.6 mm in diameter) were surface-sterilized with 0.1% HgCl₂ (BBI Life Sciences, Shanghai, China) for 10 min, followed by thorough rinsing with sterile distilled water. Under aseptic conditions, floral sepals were removed and anthers were excised and placed onto callus induction medium consisting of 1.0 mg L^-1 2,4-dichlorophenoxyacetic acid (2,4-D), 1.0 mg L^-1 kinetin (KT), 1.0 mg L^-1 naphthaleneacetic acid (NAA), 70 g L^-1 sucrose, and 50 mL L^-1 coconut water. Explants were incubated in the dark at 26–28°C for 7–8 weeks to induce loose, pale-yellow callus.

Induced calli were subsequently transferred to differentiation medium containing 1.0 g L^-1 activated charcoal, 2.0 mg L^-1 KT, 0.1 mg L^-1 NAA, 0.5 mg L^-1 gibberellic acid (GA₃), 0.2 mg L^-1 abscisic acid (ABA), 70 g L^-1 sucrose, and 50 mL L^-1 coconut water (all reagents from BBI Life Sciences). Cultures were maintained in the dark at 23–25°C. In the HE genotype, globular embryos formed after 3–4 weeks, and cotyledonary embryos developed after 8–9 weeks of continued culture. In contrast, the LE genotype produced only non-embryogenic callus under the same conditions, and no somatic embryos were observed throughout the culture period.

Samples representing different stages of SE were collected for metabolomic and transcriptomic analyses (transcriptome data reported previously in Li et al., 2023). For the HE genotype, four developmental stages were sampled: explants (HE-EX), induced callus (HE-IC), embryogenic callus (HE-EC), and cotyledonary embryos (HE-CE). For the LE genotype, three stages were collected: explants (LE-EX), induced callus (LE-IC), and non-embryogenic callus (LE-C). Three biological replicates were collected for each stage. All samples were immediately frozen in liquid nitrogen and stored at –80°C until further analysis.

Freeze-dried samples were prepared using a vacuum freeze dryer (Scientz-100F). The dried tissues were ground into a fine powder using an MM 400 mixer mill (Retsch) with zirconia beads at 30 Hz for 1.5 min. A 50 mg aliquot of powdered tissue was extracted with 1.2 mL of 70% methanol. Samples were vortexed for 30 s every 30 min, for a total of six cycles. After centrifugation at 12,000 rpm for 3 min, the supernatant was filtered through a 0.22 μm membrane (SCAA-104; ANPEL, Shanghai) and used for UPLC-MS/MS analysis.

Metabolomic profiling was performed on a UPLC-ESI-MS/MS system consisting of an ExionLC™ AD UPLC system coupled to an Applied Biosystems 6500 QTRAP mass spectrometer. Chromatographic separation was achieved using an Agilent SB-C18 column (1.8 μm, 2.1 × 100 mm) maintained at 40°C. The mobile phases were: solvent A, water containing 0.1% formic acid; and solvent B, acetonitrile containing 0.1% formic acid. The gradient program was: 0–9 min, 95–5% A; 9–10 min, 5% A; 10–11.1 min, 5–95% A; 11.1–14 min, 95% A. The flow rate was 0.35 mL min^-1 and the injection volume was 2 μL.

Electrospray ionization (ESI) was conducted in both positive and negative ion modes. Parameters were as follows: ion spray voltage, +5500 V (positive mode) or −4500 V (negative mode); source temperature, 500°C; ion source gas I (GSI), 50 psi; ion source gas II (GSII), 60 psi; curtain gas (CUR), 25 psi; and collision gas set to “high.” Metabolite detection was performed in multiple reaction monitoring (MRM) mode with nitrogen as the collision gas. Declustering potential (DP) and collision energy (CE) were optimized for each MRM transition.

Transcriptome data for all samples used in this study were obtained from our previously published work (Li et al., 2023). Briefly, RNA extraction and library construction were performed using the same samples, followed by sequencing on the Illumina platform. After quality control of the raw sequencing data, reads were aligned to the rubber tree reference genome (Tang et al., 2016) using HISAT2. DEGs were identified using the DESeq2 package with thresholds of |log₂FC| ≥ 1 and false discovery rate (FDR) < 0.05. Genes were annotated through BLAST searches against the NR, Swiss-Prot, and KEGG databases. Gene Ontology (GO) and KEGG pathway enrichment analyses were performed using the clusterProfiler package in R.

Raw LC-MS/MS data were processed for peak detection, alignment, and quantification using a proprietary data processing pipeline. Metabolites were annotated by comparing retention times, mass spectra, and MS/MS fragmentation patterns with the MetWare Database (MWDB, Metware Biotechnology Co., Ltd., Wuhan, China) and public metabolite libraries.

Principal component analysis (PCA) was performed on the data that were unit variance standardized using the prcomp function in R software. DAMs were identified using orthogonal partial least squares discriminant analysis (OPLS-DA). Variable importance in projection (VIP ≥ 1.0) and fold change (FC ≥ 2 or ≤ 0.5) were used as selection thresholds. OPLS-DA models were constructed after log₂ transformation and mean centering of metabolite intensities, and model robustness was evaluated by 200 permutation tests (p < 0.05).

DAMs were mapped to KEGG pathways using KEGG Compound IDs (Kanehisa and Goto, 2000), and pathway enrichment analysis was conducted using a hypergeometric test with significance set at p < 0.05. The same criteria were applied for transcriptome-based pathway enrichment. For multi-omics integration, DEGs and DAMs were jointly mapped to KEGG pathways to identify co-regulated metabolic nodes.

Using the UPLC-MS/MS platform combined with the MWDB metabolite database, a total of 1,383 metabolites were identified across all developmental stages of H. brasiliensis SE. These metabolites were classified into 11 major categories, with flavonoids accounting for the largest proportion (23.50%), followed by phenolic acids (17.64%) and amino acids and derivatives (10.41%). Other abundant groups included lipids (9.26%), alkaloids (6.58%), organic acids (6.29%), and nucleotides and derivatives (5.42%) (Figure 1A).

To ensure data reliability, pooled quality control (QC) samples were included throughout the analytical sequence. The empirical cumulative distribution function (ECDF) analysis showed that the coefficient of variation (CV) distributions of biological replicates from all developmental stages largely overlapped. Across all sample groups, more than 75% of metabolites exhibited CV < 0.3, and over 85% displayed CV < 0.5, indicating excellent instrumental stability and technical reproducibility (Supplementary Figure S1). Pearson correlation analysis further confirmed the high consistency among biological replicates (Figure 1C). In the HE genotype, correlations ranged from r = 0.96–0.99 across the HE-EX, HE-IC, HE-EC, and HE-CE stages. Similarly, samples from the LE genotype showed correlations of r = 0.94–0.99.

PCA revealed clear metabolic transitions during the SE process. The first two principal components (PC1 and PC2) explained 47.68% and 15.91% of the total variance, respectively, accounting for 63.59% of overall metabolic variation. Samples clustered primarily according to developmental stage rather than genotype, with HE and LE samples at the same stage grouping more closely than samples at different stages within the same genotype (Figure 1B). Hierarchical clustering analysis supported this trend, showing tight clustering of biological replicates and stage-dependent separation of samples (Figure 1D). Notably, the cotyledonary embryo stage (HE-CE), which is unique to the HE genotype, formed an independent branch in both PCA and hierarchical clustering analyses, indicating distinct metabolic signatures associated with somatic embryo formation.

To further elucidate the metabolic mechanisms underlying the contrasting SE capacities of the two rubber tree genotypes, KEGG enrichment analysis was performed for all DAMs identified across the comparison groups.

Genotype-based comparisons revealed that the flavonoid biosynthesis pathway was significantly enriched (p < 0.01) at all developmental stages (LE-EX vs. HE-EX, LE-IC vs. HE-IC, and LE-C vs. HE-EC) (Figures 5A–C). This consistent enrichment suggests that flavonoid metabolism is a key distinguishing feature between the two genotypes. Notably, during the LE-C vs. HE-EC stage, linoleic acid metabolism and glutathione metabolism exhibited the strongest enrichment, followed by flavonoid biosynthesis (Figure 5C).

Developmental stage-oriented comparisons showed that both genotypes exhibited strong enrichment of flavonoid-related pathways—particularly flavonoid biosynthesis and flavone and flavonol biosynthesis—during the explant-to-callus transition (LE-EX vs. LE-IC; HE-EX vs. HE-IC) (Figures 5D, F). In contrast, the subsequent callus-to-differentiation transition (LE-IC vs. LE-C; HE-IC vs. HE-EC) showed pronounced metabolic divergence between the two genotypes (Figures 5E, G). The LE genotype displayed significant enrichment in amino sugar and nucleotide sugar metabolism and nucleotide sugars biosynthesis (p < 0.01) (Figure 5E). Conversely, the HE genotype exhibited strong enrichment in valine, leucine and isoleucine biosynthesis and pantothenate and CoA biosynthesis (Figure 5G), pathways essential for amino acid metabolism, lipid metabolism, and cellular energy production.

In HE-EC vs. HE-CE, DAMs associated with the flavonoid biosynthesis pathway were predominantly downregulated (Figure 5H; Supplementary Figure S2), suggesting that cotyledonary embryo development may involve further suppression of secondary metabolism.

KEGG enrichment analysis revealed consistent and significant differences in the flavonoid biosynthesis pathway between the two genotypes across multiple developmental stages (Figure 5; Supplementary Figure S3). To clarify the regulatory roles of this pathway during SE, transcriptomic and metabolomic datasets were integrated, and gene-metabolite co-regulation patterns were compared during two key developmental stages: explant-to-callus induction and callus-to-differentiation (Supplementary Figure S3).

By integrating transcriptomic and metabolomic datasets, this study provides a comprehensive view of the metabolic dynamics underlying SE in H. brasiliensis, and highlights distinct metabolic behaviors between high- and low-embryogenic genotypes across developmental transitions. A total of 1,383 metabolites were identified, spanning 11 major metabolic categories. Among these, flavonoids, phenolic acids, and amino acids and their derivatives represented the dominant metabolite classes. PCA revealed that samples clustered primarily according to developmental transition rather than genotype, with the first two components explaining 63.59% of the total variance. This pattern indicates that developmental progression is the major driver of metabolic shifts during SE.

Among genotype-based comparisons, 240, 234, and 349 DAMs were detected at the LE-EX vs. HE-EX, LE-IC vs. HE-IC, and LE-C vs. HE-EC transitions, respectively, with the LE-C vs. HE-EC transition showing the highest metabolic divergence. More importantly, the functional categories of these metabolic differences underwent a marked transition. Early SE stages were dominated by secondary metabolic pathways, such as flavonoid and phenylpropanoid biosynthesis, whereas the LE-C vs. HE-EC transition was characterized by the activation of primary metabolic pathways, including linoleic acid metabolism, glutathione metabolism, and amino sugar and nucleotide sugar metabolism. This shift from secondary to primary metabolic dominance reflects a profound metabolic reprogramming process, which appears essential for successful SE progression.

Such metabolic reprogramming is consistent with observations from several other plant species. In multiple SE processes, early dedifferentiation is often associated with stress-related pathways (Fehér, 2015; Nic-Can and Loyola-Vargas, 2016; Ge et al., 2025). In Gossypium hirsutum, for example, high-embryogenic genotypes show DEGs enriched in fatty acid, tryptophan, and pyruvate metabolism, whereas low-embryogenic genotypes exhibit enrichment in pathways related to DNA conformational changes, resulting in divergent metabolic trajectories (Guo et al., 2020). A similar pattern has been reported in Coffea arabica, where SE involves a coordinated transition between stress-associated metabolism and growth-related primary metabolic pathways (Awada et al., 2019, Awada et al., 2023).

KEGG enrichment analysis revealed that the flavonoid biosynthesis pathway was significantly enriched in both developmental and genotype-based comparisons. Integrated transcriptome-metabolome analysis further demonstrated that this pathway is tightly co-regulated during early dedifferentiation, followed by pronounced divergence during the callus-to-differentiation transition. During the explant-to-callus transition, both genotypes exhibited highly consistent transcriptional and metabolic responses in the flavonoid biosynthesis pathway, with 67.5% of DEGs and 85.7% of DAMs shared between genotypes, all showing identical regulatory directions. In contrast, during the callus-to-differentiation transition, only 37.5% of DEGs and 6.7% of DAMs were shared, and 11 key genes displayed entirely opposite regulatory trends.

The early-stage conservation was primarily reflected in the coordinated activation of the anthocyanin precursor biosynthetic branch. From CHS (carbon backbone synthesis) through CHI, F3H, DFR, and LDOX (progressive structural modifications), multiple core enzymes were significantly upregulated in both genotypes. This activation was accompanied by pronounced depletion of flavanone precursors (naringenin and butin), together with strong accumulation of downstream dihydroflavonols (aromadendrin, taxifolin, and dihydromyricetin), consistent with enhanced metabolic flux through the anthocyanin precursor branch. Meanwhile, FLS genes were strongly downregulated (log₂FC –9.69 to –12.75), together with a marked reduction in flavonol metabolites such as galangin and quercetin. Flavone compounds (apigenin, luteolin) also decreased. This metabolic partitioning suggests that carbon flux during callus formation is preferentially directed toward anthocyanin and proanthocyanidin biosynthesis rather than flavonol production. Consistent with this trend, ANR and LAR showed coordinated upregulation, and the corresponding proanthocyanidin monomers epigallocatechin and gallocatechin accumulated robustly.

The conserved activation of the flavonoid pathway during early dedifferentiation likely reflects a shared stress response triggered by tissue culture conditions. Explant excision and callus induction generate substantial oxidative stress, necessitating enhanced antioxidant defenses (Cassells and Curry, 2001; Fehér, 2015; Zhou et al., 2016). Dihydroflavonols and proanthocyanidins possess strong antioxidant properties (Brunetti et al., 2013), and their accumulation may provide redox buffering capacity, protecting cells from oxidative damage associated with dedifferentiation (Zhou et al., 2016; Huang et al., 2019).

The strong upregulation of F3’5’H and the preferential accumulation of its products—highly hydroxylated flavonoids such as dihydromyricetin—further support this interpretation. High levels of B-ring hydroxylation are typically associated with enhanced antioxidant capacity (Agati et al., 2012), indicating that both genotypes synthesize more potent antioxidant metabolites as part of their dedifferentiation-associated stress response.

Similar patterns of flavonoid accumulation during SE have been reported across diverse plant species, including Korean pine (Peng et al., 2022), Gossypium hirsutum (Guo et al., 2019), Paeonia ostii (Zhang et al., 2023), and Silybum marianum (Khan et al., 2015), where flavonoids were suggested to promote differentiation and create a more favorable cellular environment for embryo formation. Collectively, these findings support a conserved role of flavonoid-mediated antioxidant activity during the initial phases of somatic embryogenesis.

In contrast to the conserved metabolic responses observed during early dedifferentiation, the callus-to-differentiation transition was characterized by striking metabolic divergence between the two genotypes. This divergence corresponded closely with their contrasting embryogenic capacities.

The HE genotype adopted a conversion-type metabolic pattern during the HE-IC vs. HE-EC transition. From the downregulation of the upstream rate-limiting enzyme CHS to the broad suppression of downstream glycosylation machinery, the flavonoid biosynthesis pathway was systematically repressed at the transcriptional level, accompanied by an overall attenuation of metabolite accumulation, despite the persistence of selected intermediates. This reconfiguration is consistent with KEGG enrichment results, which showed that the HE genotype preferentially activated primary metabolic pathways essential for growth and development, including valine, leucine and isoleucine biosynthesis and pantothenate and CoA biosynthesis.

By contrast, the LE genotype exhibited a maintenance-type metabolic pattern. Although the biosynthesis rate of anthocyanin precursors declined during differentiation, the glycosylation system became further activated, maintaining a relatively elevated level of secondary metabolic activity compared with the HE genotype, despite partial attenuation during differentiation. This pattern aligned with the enrichment of amino sugar and nucleotide sugar metabolism and biosynthesis of nucleotide sugars in the LE genotype, suggesting a metabolic configuration that favors secondary metabolism even as the tissue fails to undergo somatic embryo formation.

These contrasting patterns likely reflect distinct metabolic resource allocation strategies. The condensation reaction catalyzed by CHS requires malonyl-CoA (Winkel-Shirley, 2001; Austin and Noel, 2003), which is also a critical precursor for fatty acid biosynthesis (Rawsthorne, 2002). Developing somatic embryos demand substantial lipid production to support rapid cell division and membrane biogenesis (Baud and Lepiniec, 2010; Troncoso-Ponce et al., 2011). Thus, persistent high CHS expression in the LE genotype may create competition with lipid metabolism for shared precursors, potentially constraining embryo development. In contrast, the downregulation of CHS in the HE genotype during differentiation may alleviate this competition, thereby redirecting carbon resources toward primary metabolic processes necessary for embryogenesis.

Similarly, the widespread suppression of the glycosylation system in the HE genotype likely increases the availability of nucleotide sugar precursors. Glycosylation reactions consume substantial amounts of UDP-glucose, a central metabolic hub involved in a range of biosynthetic pathways (Bar-Peled and O’Neill, 2011). UDP-glucose serves not only as a glycosyl donor for secondary metabolites (Vogt and Jones, 2000), but also as the direct substrate for cellulose synthesis (Delmer, 1999). Through UDP-glucose dehydrogenase, it can be converted to UDP-glucuronic acid and subsequently to UDP-xylose and UDP-arabinose for hemicellulose and pectin biosynthesis (Lim and Bowles, 2004; Lerouxel et al., 2006). Allocation of UDP-glucose to flavonoid glycosylation may therefore compete with cell wall biosynthesis, potentially restricting cell division and tissue organization required for embryo development. The dual repression of precursor synthesis and glycosylation in the HE genotype could thus increase the pool of UDP-glucose and other nucleotide sugars available for primary cell wall formation during embryogenesis.

The concept of a metabolic trade-off between secondary metabolism and embryogenic capacity has been reported in several conifer and angiosperm SE systems. In Pseudotsuga menziesii, 48 transcripts associated with flavonoid metabolism were strongly upregulated (up to 34-fold) in non-embryogenic callus (NEC), accompanied by tissue browning and a complete loss of embryogenic competence (Gautier et al., 2018; Gautier et al., 2019). The authors suggested that the characteristic features of NEC included the upregulation of stress response metabolites and a shift in carbohydrate metabolism toward starch storage, indicating that excessive secondary metabolic activity may interfere with the coordinated metabolic patterns required for efficient cell differentiation. A similar conclusion was drawn in Larix kaempferi, where 39% of DAMs were flavonoids and phenolic acids predominantly accumulated in NEC, suggesting that excessive secondary metabolism may interfere with primary metabolism and delay differentiation (Wang et al., 2022). In Paeonia ostii, flavonoid-related genes—including F3’5’H and ANS—were downregulated in EC and SE relative to NEC (Zhang et al., 2023). Moreover, addition of activated charcoal, which adsorbs phenolics, improved embryogenic callus induction and somatic embryo development in Gossypium hirsutum (Kouakou et al., 2007; Thomas, 2008), suggesting that the excessive accumulation of phenolic compounds (including flavonoids) may have an inhibitory effect.

It should be noted that the metabolic differentiation and morphological differentiation in our time-course analysis showed relatively synchronized patterns, with no significant metabolic changes preceding morphological changes. This suggests that changes in flavonoid metabolism may be coordinated with other developmental processes.

A central question emerging from our data is whether the differential regulation of the flavonoid biosynthesis pathway is a cause, a consequence, or a modulator of embryogenic capacity. While our analyses revealed strong associations between metabolic divergence and developmental outcomes, the underlying causal relationships remain unresolved.

Several studies suggest that upstream developmental regulators may coordinate metabolic reprogramming during SE. Key transcription factors such as BBM, WUS, and LEC1 are well-established determinants of embryogenic competence (Braybrook and Harada, 2008; Horstman et al., 2017). In Larix kaempferi, transcriptome profiling and predicted gene regulatory networks identified specific ERF, MYB, and DOF transcription factors as potential regulators of phenylpropanoid and flavonoid biosynthesis (Wang et al., 2022). MYB proteins are known to form MBW complexes with bHLH and WD40 partners to modulate flavonoid pathway gene expression (Li, 2014; Khan et al., 2015). These findings suggest that differences in upstream transcriptional regulatory networks associated with embryogenic competence may indirectly influence SE outcomes by coordinating the expression of multiple developmental and metabolic pathways, including flavonoid biosynthesis.

Conversely, metabolic states may exert feedback effects on developmental capacity. Several mechanisms have been proposed: (1) Resource competition, as discussed above, where flavonoid biosynthesis competes with primary metabolic pathways for shared precursors. (2) Modulation of hormone signaling and polar auxin transport, mediated by flavonoid effects on auxin efflux carriers and related regulatory proteins (Peer and Murphy, 2007). (3) Changes in cellular redox homeostasis, which influence signal transduction. In Picea species, a reduced redox state promotes somatic embryo induction, whereas a more oxidized environment favors embryo differentiation and development (Stasolla et al., 2003b). Furthermore, disturbances in ROS homeostasis interact with auxin signaling to regulate SE (Zhou et al., 2016; Kudełko and Gaj, 2019).

Although these mechanisms suggest plausible links between flavonoid accumulation and SE competence, direct experimental evidence connecting specific flavonoid species to embryogenesis remains limited.

Based on available findings, we propose that the relationship between flavonoid metabolism and embryogenic capacity represents a dynamic, bidirectional regulatory system rather than a linear causal pathway. Initial embryogenic competence—shaped by genetic and epigenetic factors—may trigger specific metabolic reprogramming events, including activation or suppression of flavonoid biosynthesis. In turn, these metabolic changes influence cell physiology through multiple layers, such as resource allocation, hormone transport, and redox balance, thereby positively or negatively affecting the progression of embryogenic development. Such reciprocal interactions may amplify small initial differences in competence into pronounced developmental outcomes.

In the context of our experimental system, one possible interpretation is that timely suppression of flavonoid biosynthesis in the HE genotype during the callus-to-differentiation transition reduces carbon and energy diversion into secondary metabolism, thereby favoring the biosynthetic demands of primary metabolism required for embryo formation. This could contribute to a self-reinforcing positive feedback loop that supports successful somatic embryogenesis. In contrast, the LE genotype maintains high levels of secondary metabolism, potentially causing persistent carbon and precursor diversion that interferes with the metabolic restructuring necessary for embryogenic development.

However, we emphasize that these interpretations remain hypothetical. Definitive resolution of the causal connections between flavonoid metabolism and embryogenic capacity will require targeted experiments, such as metabolic flux analysis, exogenous application or chemical inhibition of specific flavonoids, and isotope tracing to directly test how metabolic pathways influence SE outcomes. Future studies combining targeted manipulation of key DAMs (e.g., polyamine conjugates, flavonoid aglycones) with transcription factor binding assays (such as ChIP-qPCR for BBM, LEC1, or MYB regulators) will be essential to establish direct links between transcriptional regulation, metabolic flux, and embryogenic competence.

In addition to pathway-level differences, several individual DAMs showed pronounced stage- or genotype-specific accumulation patterns, suggesting potential functional roles during SE. Although their precise mechanisms remain unclear, these metabolites represent important clues to the metabolic framework underlying embryogenic competence.

N1,N10-bis(p-coumaroyl)spermidine, which accumulated across multiple stages in the HE genotype, is a derivative of polyamines such as spermidine and spermine—molecules known to be essential for plant growth, development, and in particular for cell division and differentiation (Hasegawa et al., 1984; Kusano et al., 2008). Acylation of polyamines can alter their biological activity and cellular localization (Bassard et al., 2010). In Arabidopsis, polyamine acyltransferases participate in floral organ development (Grienenberger et al., 2009). The specific accumulation of this compound in the HE genotype may therefore reflect enhanced polyamine-related regulatory activity supporting its strong proliferative and differentiation capacity.

Vomifoliol, a degradation product of ABA and carotenoid oxidative cleavage (Hasegawa et al., 1984), accumulated specifically in LE-IC vs. LE-C. ABA plays a critical role in SE, particularly in embryo maturation. Exogenous ABA promotes SE maturation in Quercus suber and Vitis vinifera, increasing dry mass accumulation (Garcia-Martin et al., 2005; Acanda et al., 2020). In Arabidopsis, ABA has been reported to play complex and sometimes contrasting roles in SE. Some studies indicate that modulation of ABA signaling affects SE efficiency, with ABA receptor overexpression repressing SE and ABA-deficient mutants showing increased SE capacity (Chen et al., 2021). Vomifoliol accumulation may therefore indicate altered ABA metabolism or signaling in the LE genotype during differentiation.

The LE genotype also exhibited significant accumulation of purine metabolism intermediates, including xanthine and hypoxanthine. These compounds arise from nucleotide catabolism but can also be recycled through salvage pathways for nucleotide synthesis (Ashihara et al., 2018). Purine nucleotides are essential for DNA/RNA synthesis, ATP-driven energy metabolism, and signaling (Stasolla et al., 2003a). Studies in Picea glauca SE showed that proliferating cells rely on both de novo and salvage pathways to meet nucleotide demand (Stasolla and Thorpe, 2004). The accumulation of these catabolic intermediates in the LE genotype may reflect an imbalance in nucleotide metabolism.

L-ascorbic acid (vitamin C) was upregulated during the HE-EC vs. HE-CE stage. As one of the most abundant and functionally important soluble antioxidants in plants, ascorbic acid participates in reactive oxygen species (ROS) scavenging, regeneration of glutathione and other antioxidants, and serves as a cofactor for various dioxygenases, thereby playing essential roles in plant development and stress responses (Ishikawa et al., 2018; Celi et al., 2023; Conklin et al., 2024). In plant tissue culture systems, numerous studies have demonstrated that ascorbic acid profoundly influences somatic embryogenesis and regeneration. It effectively suppresses the oxidation of phenolic compounds released from wounded explants, thereby reducing tissue browning and cell death, which enhances explant survival and the viability of embryogenic cells (Amente and Chimdessa, 2021; Wang et al., 2024). Moreover, modulation of the glutathione/ascorbate redox state has been shown to significantly affect the efficiency of embryogenic callus induction, the establishment of embryo polarity, and the formation of meristematic primordia (Stasolla, 2010; Ishikawa et al., 2018; Conklin et al., 2024). Therefore, we hypothesize that the upregulation of L-ascorbic acid during somatic embryo formation in the HE genotype may help protect rapidly dividing embryogenic cells from ROS-induced damage and facilitate cell-wall remodeling and other growth processes associated with embryo patterning, ultimately supporting normal somatic embryo development and meristem formation.

Our findings have practical implications for improving somatic embryogenesis efficiency in H. brasiliensis and potentially other recalcitrant species within the experimental framework examined here. Identifying the callus-to-differentiation transition as a critical metabolic checkpoint suggests that interventions during this stage are likely to be most effective. Definitive resolution of causal relationships within this framework will require targeted functional validation. In particular, genetic manipulation of key flavonoid biosynthetic nodes such as CHS and flavonoid glycosyltransferases, together with in vitro supplementation or inhibition assays, would provide a direct means to test whether modulation of flavonoid flux influences embryogenic competence.

The integration of metabolomic and transcriptomic datasets highlights the importance of evaluating metabolic phenotypes, not only gene expression patterns, when optimizing SE systems. Medium composition may be refined to favor primary metabolism over secondary metabolism during the callus-to-differentiation transition—for example, through adjusting nitrogen sources, carbon-nitrogen ratios, or supplementing metabolic intermediates that redirect flux away from the phenylpropanoid pathway. It should be noted that the metabolic reprogramming captured in this study likely represents downstream outputs of broader transcriptional regulation during SE. Key embryogenesis-associated transcription factors, such as BBM and LEC, are known to orchestrate large-scale developmental transitions and may indirectly shape the observed accumulation patterns of core DAMs through regulation of hormone- and redox-related genes. Future integration of transcription factor binding information with transcriptomic and metabolomic datasets would facilitate the construction of a conceptual “transcription factor–gene–metabolite” regulatory network, providing a systems-level framework for linking developmental regulators with metabolic state transitions during somatic embryogenesis.

Several questions remain. First, the causal relationship between flavonoid metabolism and embryogenic capacity requires direct experimental validation, such as genetic manipulation or chemical inhibition of specific flavonoids. Second, the proposed resource competition model involving malonyl-CoA, UDP-glucose, and other precursors remains theoretical without supporting metabolic flux evidence. Third, the upstream regulators that coordinate multiple metabolic programs are still unknown, and may involve higher-order regulatory layers such as epigenetic mechanisms, including chromatin accessibility, DNA methylation, and histone modification. Finally, while flavonoid biosynthesis emerged as the pathway with the most consistent enrichment, other metabolic processes—including amino acid, nucleotide, and lipid metabolism—also differed between genotypes and may play important roles.

Importantly, these metabolic classes should not be viewed as independent of flavonoid metabolism but as interconnected components of a coordinated metabolic reprogramming underlying somatic embryogenesis. Amino acid, nucleotide, lipid, polyamine, and antioxidant pathways are metabolically linked to phenylpropanoid metabolism through shared carbon precursors, redox homeostasis, and cellular energy demands, and their genotype- and stage-specific alterations likely reflect shifts in resource allocation accompanying the acquisition of embryogenic competence.

In this context, flavonoid biosynthesis emerged as the most consistently enriched and clearly differentiated pathway between genotypes, justifying its detailed discussion in this study. Nevertheless, the accumulation patterns of non-flavonoid metabolites indicate that embryogenic potential is regulated by an integrated metabolic network rather than by a single compound class, with flavonoids representing a prominent but not exclusive node within this framework. Future studies should therefore focus on how multiple metabolic pathways are integrated and co-regulated during somatic embryogenesis, rather than on individual pathways in isolation.

Although the present study focused on two genotypes with contrasting embryogenic capacities using a single explant type, this design enabled clear resolution of genotype-dependent metabolic divergence under highly controlled developmental conditions. Embryogenic competence, however, is likely a quantitative trait rather than a binary state, and inclusion of genotypes with intermediate embryogenic capacity will be important for determining whether the observed suppression of flavonoid biosynthesis during the callus-to-differentiation transition follows a gradual or threshold-dependent pattern. Furthermore, because somatic embryogenesis can be initiated from diverse explant sources with distinct physiological and stress-response backgrounds, comparative metabolomic analyses across different explant types (e.g., shoot tips, leaves) will be valuable for assessing whether flavonoid-associated metabolic reprogramming reflects a general regulatory principle or an explant-specific response. Such expanded experimental frameworks will help define the broader applicability of the metabolic trade-offs proposed in this study.