The seeds of A. manihot, sourced from a noncontaminated area in Sichuan Province, China, were subjected to surface sterilization by soaking in 0.05% sodium hypochlorite (NaClO) for 30 min followed by rinsing with deionized water. The sterilized seeds were directly planted in baskets filled with sterilized ceramsite stones and placed within the planting basket of an automatic hydroponic culture device ( Supplementary Figure S2 ). Germination occurred under controlled conditions maintaining a constant temperature of 25°C and a photoperiod of 14/10 h (day/night). After 2 weeks, the seedlings were precultured in a greenhouse using Hoagland nutrient solution, with a temperature setting of 28°C/24°C (day/night), a photoperiod of 14/10 h (day/night), and humidity maintained between 60% and 80%. The components of the nutrient solution were as follows: 4 mM Ca(NO₃)₂·4H₂O, 4 mM (NH₄)₂SO₄, 4 mM K₂SO₄, 4 mM KNO₃, 1.3 mM KH₂PO₄, 1 mM MgSO₄·7H₂O, 50 μM Fe-EDTA, 10 μM H₃BO₃, 5 μM MnSO₄·H₂O, 5 μM ZnSO₄·7H₂O, 1 μM CuSO₄·5H₂O, and 0.5 μM Na₂MoO₄·2H₂O (Mwamba et al., 2016). After 30 days of growth, the plants were transferred to Hoagland nutrient solution containing CdCl₂, while the control plants were not treated with CdCl₂ (CK). An additional preliminary experiment assessed the physiological parameters of A. manihot roots under a range of Cd concentrations, including 0, 50, 100, 200, and 400 μM (data not shown), which allowed the level of Cd exposure to be set at 100 μM in the present study. Plant samples were collected at four different time points: 12 h, 36 h, 72 h, and 7 days. At each time point, a total of 12 individuals were subdivided into three biological replicates. The harvested plants were washed with tap water, and the roots were treated with 20 mM ethylenediaminetetraacetic acid disodium salt (Na₂-EDTA) for 15 min to remove the Cd adsorbed on the root surface.

Dried samples of roots, stems, and leaves were ground, passed through a 0.15-mm mesh sieve, and digested using an acid mixture of nitric acid (HNO₃) and perchloric acid (HClO₄) (v:v, 4:1) (Wu et al., 2018). The Cd concentrations in the different plant tissues were analyzed by atomic absorption spectrometry (PinAAcle 900 H, Perkin Elmer, USA). Electrodeless discharge lamp at 228.8 nm (with a slit width of 0.7 nm) was used as the radiation source for the Cd.

Fresh root tissues were utilized directly for the detection of biochemical parameters. The malondialdehyde (MDA) content was measured using the thiobarbituric acid (TBA) reaction method (Peever and Higgins, 1989). The activities of superoxide dismutase (SOD), peroxidase (POD), and catalase (CAT) were determined following the procedure described by Wu et al. (2018). SOD activity was determined by nitroblue tetrazolium (NBT) reduction, POD activity through guaiacol oxidation, and CAT activity by measuring the decrease in H₂O₂ concentration.

The transcriptome and metabolome data were normalized and subjected to statistical analysis to elucidate the relationships between the genes and metabolites implicated in the Cd stress response. Pearson correlation analysis between DEGs and DEMs was performed using the COR function in R with normalized data. Subsequently, DEGs and DEMs with Pearson correlation coefficient (PCC) thresholds of ≥0.8 or ≤−0.8 were subjected to KEGG analysis.

Statistical analysis was performed using SPSS 20.0 (SPSS Inc., Chicago, IL, USA). Statistical significance was calculated using two- or one-way ANOVA followed by Tukey’s multiple comparison test with differences deemed significant at p < 0.05. Graphs were generated using GraphPad Prism 9 software (GraphPad Software Inc., La Jolla, CA, USA).

The present study analyzed the Cd concentrations in the roots, stems, and leaves of A. manihot at four different time points following exposure to 100 μM Cd. The results revealed a significant increase in the Cd content in A. manihot with increasing duration of Cd exposure ( Figure 1 ). The average Cd concentrations in the roots were 33.21, 73.95, 119.00, and 137.38 mg·kg⁻¹ at 12 h, 36 h, 72 h, and 7 days after Cd treatment, respectively. Initially, the Cd concentrations in the shoots remained lower than those in the roots at 12–36 h. However, the Cd concentration in the shoots significantly increased to 203.81 mg·kg⁻¹ at 7 days indicating a 1.48-fold increase over that in the roots. These results verified the high Cd accumulation capability of A. manihot, as indicated by a previous study (Wu et al., 2018).

A. manihot demonstrated strong tolerance to Cd treatment showing no visible phytotoxicity symptoms. To evaluate the response of A. manihot roots to Cd stress, several biochemical indices associated with the oxidative stress response were monitored. MDA, one of the final products of cell membrane lipid peroxidation, is recognized as a well-known biomarker for oxidative lipid damage caused by increased ROS under abiotic stress (Wu et al., 2018; Yu et al., 2021). Compared with that in the control plants, the MDA content in the Cd-treated plants significantly increased by 64.0%–104.12% indicating that Cd induced cellular oxidative stress ( Figure 2A ). Plants respond to oxidative stress by enhancing their antioxidant defense systems (Wu et al., 2018; Zhang et al., 2021). In the roots of A. manihot, compared with those in the CK treatment, Cd stress significantly increased the activity of SOD, POD, and CAT by 19.73%–50%, 22.87%–38.89%, and 32.31%–45.40%, respectively, at 12 h, 36 h, 72 h, and 7 days. SOD activity significantly increased under Cd stress, peaking at 36 h, and then gradually decreased ( Figure 2B ). A similar pattern was observed for POD activity, which peaked at 72 h ( Figure 2C ). Additionally, the CAT activity in roots increased continuously with the duration of Cd exposure ( Figure 2D ).

After data filtering, approximately 40.61–66.83 million 150-bp high-quality clean reads were obtained from the RNA-seq data of 24 samples. The clean read Q20 and Q30 values for all test samples were greater than 97.51% and 92.74%, respectively ( Supplementary Table S2 ) demonstrating the high quality of sequencing and the feasibility of subsequent analyses. De novo assembly of the clean reads using the Trinity program generated a reference transcriptome comprising 73,122 unigenes, with an average length of 873 bp, a maximum length of 13,204 bp, a minimum length of 201 bp, and an N50 of 1,268 bp. Using the BLASTX tool with a cutoff E-value of 10⁻⁵, a total of 56,752 unigenes were functionally annotated ( Supplementary Table S3 ). Based on extensive database searches, Gossypium hirsutum (22.47%), Gossypium arboreum (22.14%), and Gossypium raimondii (21.91%) exhibited high similarity to A. manihot ( Supplementary Figure S3 ).

To further evaluate the validity of the RNA-seq data, 12 unigenes were randomly selected for qPCR analysis. The RT-qPCR results for these genes showed similar expression trends to the RNA-seq data ( Supplementary Figure S4 ) indicating the accuracy of the RNA-seq data. Furthermore, a strong positive correlation was observed between the qPCR and RNA-seq data (r = 0.7581, p < 0.001) ( Supplementary Figure S5 ). These results confirm the validity of RNA-seq and suggest that the RNA-seq data can be reliably used for further analyses.

The DEGs of A. manihot were identified by comparing RNA-seq data from Cd-treated root samples with CK samples at the same time points. Comparative DEG analysis revealed a total of 245 DEGs (122 upregulated and 123 downregulated) at 12 h, 5,708 DEGs (629 upregulated and 5,079 downregulated) at 36 h, 9,834 DEGs (3,606 upregulated and 6,228 downregulated) at 72 h, and 2,323 DEGs (647 upregulated and 1,676 downregulated) at 7 days ( Figure 3A ).

Subsequently, the functions of the DEGs were classified according to GO classifications. The downregulated DEGs were predominantly associated with “structural molecule activity” in the molecular function category, “intracellular ribonucleoprotein complex” and “ribonucleoprotein complex” in the cellular component category, and “phenylpropanoid metabolic process,” “lignin metabolic process,” and “cell wall organization or biogenesis” in the biological process category ( Figure 3B , Supplementary Table S4 ). These results suggest that Cd stress impacts the cell wall organization or biogenesis of A. manihot. In addition, the upregulated DEGs were significantly enriched in molecular function categories, such as “tetrapyrrole binding,” “catalytic activity,” and “oxidoreductase activity,” while the cellular component category was dominated by “intracellular ribonucleoprotein complex,” “ribonucleoprotein complex,” and “ribosome.” In addition, the biological process category was the most enriched in “defense response,” “metabolic process,” and “salicylic acid-mediated signaling pathway” ( Figure 3B , Supplementary Table S4 ). These results indicate that A. manihot can enhance catalytic activity and tetrapyrrole binding and activate defense mechanisms and salicylic acid-mediated signaling pathways to achieve stable Cd detoxification or improve Cd tolerance (Yan et al., 2020).

To analyze the temporal expression patterns of Cd-responsive genes, STEM analysis was employed to cluster DEGs with similar expression patterns ( Supplementary Figure S6 ). In the STEM analysis, profiles 1, 5, 7, 8, and 11 with p < 0.05 were considered significant profiles ( Figure 4A ). Subsequently, GO enrichment analysis of the DEGs within these significant profiles was performed to elucidate the functional significance of the transcriptional changes induced by Cd stress ( Figure 4B ). In profiles 1 and 7, genes associated with cell wall biosynthetic and metabolic pathways were particularly enriched, including “glucuronoxylan metabolic process,” “xylan metabolic process,” “lignin metabolic process,” “phenylpropanoid metabolic process,” and “polysaccharide catabolic process.” Gene expression in profile 7 decreased after 36 h, whereas in profile 1, there was a decrease in gene expression from 12 to 72 h, which subsequently increased after 72 h. Profiles 5 and 8 also shared a similar expression pattern to profile 1, where genes were primarily associated with energy and substance metabolism. Profile 5, enriched in GO terms related to “glycerol-3-phosphate metabolic process,” showed a decrease in gene expression levels from 12 to 36 h followed by an increase thereafter. Profile 8, enriched in GO terms, such as “ATP biosynthetic process,” “purine nucleoside triphosphate biosynthetic process,” and “purine ribonucleoside triphosphate biosynthetic process,” demonstrated a decrease in gene expression from 36 to 72 h but an increase afterward. In contrast, profile 11 exhibited genes that were upregulated between 36 and 72 h and subsequently downregulated after 72 h. The enriched pathways represented in this profile included “lipid homeostasis,” “acylglycerol homeostasis,” “response to toxic substances,” and “response to zinc ion” ( Figure 4B ).

WGCNA was used to investigate various physiological parameters and key genes involved in Cd detoxification and accumulation in A. manihot. This analysis identified 17 modules with an unsigned TOM through the dynamic tree cutting ( Supplementary Figure S7 ). Notably, the purple module exhibited significant positive correlations with Cd accumulation in both roots (r = 0.77, p < 0.05) and shoots (r = 0.74, p < 0.05), as did MDA (r = 0.67, p < 0.05), POD (r = 0.50, p < 0.05), and CAT (r = 0.71, p < 0.05) ( Figure 5A ). We speculated that genes in the purple module were associated with enhanced resistance and Cd accumulation in A. manihot. Therefore, the purple module was selected for further analysis.

There were 577 genes in the purple module. Based on the gene expression heatmaps and eigengene histograms, the genes within the purple module were upregulated at 72 h and 7 days following Cd treatment ( Figure 5B ). Among them, 321 and 85 genes in the purple module were differentially expressed at 72 h and 7 days, respectively ( Figure 5C ). Subsequently, KEGG pathway analysis was conducted on the overlapping DEGs (344 genes in total). The overlapping DEGs predominantly mapped to pathways, such as “phenylpropanoid biosynthesis,” “glutathione metabolism,” “vitamin B6 metabolism,” and various lipid metabolism-related pathways, including “glycerolipid metabolism,” “glycerophospholipid metabolism,” “fatty acid elongation,” and “linoleic acid metabolism” ( Figure 5D ).

In the purple module, the eight genes with the most connections in the network were considered hub genes ( Figure 5E , Supplementary Table S5 ). Among these genes, four hub genes were involved in the phenylpropanoid biosynthesis pathway, including three POD homologs (Unigene0046043, Unigene0046044, and Unigene0010159) and a coniferyl-aldehyde dehydrogenase (REF1) homolog (Unigene0006924). Additionally, Unigene0053629 is a homolog of the ethylene-responsive transcription factor ERF98, which belongs to the AP2/ERF transcription factor family. Other hub genes included an ABC transporter homolog (Unigene0055296), a cytochrome P450 CYP82D47 homolog (Unigene0041804), and a RING-H2 finger protein ATL3 homolog (Unigene0035257). Taken together, these core physiological processes and the hub genes that regulate them may be involved in the regulation of A. manihot’s response to Cd in terms of uptake, translocation, and detoxification.

To investigate the simultaneous changes in the metabolome and identify the metabolic adaptations of A. manihot to Cd stress, a broadly targeted metabolomic analysis of its roots was performed. Using an LC-MS/MS-based widely targeted metabolomics approach, we detected a total of 722 metabolites ( Supplementary Table S6 ). Multivariate PCA and OPLS-DA were subsequently applied to identify differentially abundant metabolites across various time points. The PCA results showed that biological repeat samples within each group clustered together demonstrating the stability and reliability of the sequencing results and confirming good reproducibility among samples in each group ( Figure 6A ). Notably, as the duration of Cd exposure increased, distinct differences emerged among the samples. According to predefined criteria, a total of 66 (20 upregulated, 46 downregulated), 62 (16 upregulated, 46 downregulated), 156 (89 upregulated, 67 downregulated), and 90 (55 upregulated, 35 downregulated) DEMs were identified at 12 h, 36 h, 72 h, and 7 days, respectively ( Figure 6B ). These findings highlight the differences in metabolic responses at various stages of Cd stress, suggesting the occurrence of adaptive metabolic adjustments in A. manihot as the stress progresses.

Among the identified DEMs, six metabolites consistently appeared across all time points, including D-pipecolinic acid, adenosine O-ribose, melatonin, N,N-didesvenlafaxine, LysoPE 18:3, and LysoPE 18:1 ( Figure 6C ). The upregulated DEMs were categorized into 12 major groups, predominantly including alkaloids, lipids, flavonoids, phenolic acids, organic acids, amino acids, and nucleotide and its derivatives ( Supplementary Table S7 ). Furthermore, the analysis revealed the participation of several phytohormones, including abscisic acid and gibberellin, in the plant’s adaptive response to Cd stress ( Supplementary Table S7 , Supplementary Figure S8 ). Further analysis through KEGG pathway analysis of these upregulated DEMs highlighted their involvement in fatty acid-related pathways such as “linoleic acid metabolism,” “α-linolenic acid metabolism,” “fatty acid biosynthesis,” and “fatty acid degradation” ( Supplementary Figure S8 ).

To gain further insight into the potential regulatory mechanisms influenced by Cd stress, we conducted an integrated analysis of the metabolome and transcriptome data from the roots of A. manihot. KEGG enrichment analysis of correlated DEGs and DEMs revealed that phenylpropanoid biosynthesis was significantly enriched in both the metabolome and transcriptome datasets across the four time points ( Figure 7 , Supplementary Table S8 ). Furthermore, lipid metabolism-related pathways, which were enriched at all four time points, exhibited temporal differences. In summary, the analysis of metabolome and transcriptome data identified lipid metabolism, especially α-linolenic acid metabolism, and phenylpropanoid biosynthesis as key metabolic pathways involved in the response of A. manihot to Cd stress. Hence, a detailed analysis of the lipid metabolism and phenylpropanoid biosynthesis pathways was conducted ( Figures 8 , 9 ).

In the lipid metabolism pathway, significant variations in the levels of lysophospholipids, including lysophosphatidylcholines (LysoPCs), lysophosphatidylethanolamines (LysoPEs), and lysophosphatidylserines (LysoPSs), were observed. At 12 h, there were significant increases in the levels of lysoPC C16:0 (2) and C18:3 (2). Interestingly, most lysoPEs showed significant decreases relative to the controls, except for lysoPE 18:0 (2), which increased at 72 h. Furthermore, phospholipase A (PLA), a crucial gene involved in the metabolism of fatty acids and lysophospholipids, was significantly downregulated at 36 h and subsequently upregulated at 72 h. At common time points, significant upregulation of metabolites, such as 13-hydroperoxyoctadecatrienoic acid (13-HPOTrE), 9-hydroxy-10E,12Z,15Z-octadecatrienoic acid (9-HOTrE), and stearidonic acid in α-linolenic acid metabolism, was also observed under Cd stress. Consistent with these findings, the reaction products of upstream genes, including fatty acid desaturase (FAD) and lipoxygenase (LOX), were significantly upregulated at 72 h in response to Cd stress.

13-HPOTrE serves as the key intermediate in the biosynthesis pathway of jasmonic acid (JA) from the α-linolenic acid metabolic pathway (Feussner and Wasternack, 2002). The sustained upregulation of 13-HPOTrE suggested that the synthesis of JA from α-linolenic acid is strongly induced by Cd stress. Moreover, the expression levels of several JA pathway-related structural genes, including allene oxide synthase (AOS) and 12-oxophytodienoate reductase (OPR), and β-oxidation-related genes, such as acyl-CoA oxidase (ACX), enoyl-CoA hydratase/3-hydroxyacyl-CoA dehydrogenase (MFP2), and acetyl-CoA acyltransferase (ACAA1), were upregulated at 72 h. These results provide further evidence that the synthesis of JA was strongly induced by Cd stress.

In the phenylpropanoid biosynthesis pathway, Cd stress triggered the upregulation of phenylalanine ammonia lyase (PAL), cinnamate 4-hydroxylase (C4H), caffeic acid O-methyltransferase (COMT), ferulate 5-hydroxylase (F5H), 4-coumarate-CoA ligase (4CL), and POD expression at 72 h or 7 days. This resulted in marked increases in the production of p-coumaric acid, trans-cinnamic acid, ferulic acid, and sinapic acid. In addition, REF1 exhibited upregulated expression at 72 h and 7 days. Concurrently, there were significant decreases in the concentrations of three primary monolignols: sinapyl alcohol, p-coumaryl alcohol, and coniferyl alcohol. Given the critical role of POD in the final steps of phenylpropanoid biosynthesis, particularly in lignin polymerization through the oxidation of monolignols, the reductions in these precursors coupled with the upregulation of POD may lead to an increase in lignin content.

Exposure to Cd disturbs the redox homeostasis of plants, leading to pronounced increases in ROS production and lipid peroxidation, consequently triggering oxidative stress (Asgari Lajayer et al., 2017; Haider et al., 2021). In response to Cd stress, physiological traits, which are pivotal for plant survival and adaptation under challenging conditions, play significant roles in orchestrating defense strategies (Pan et al., 2021; Xu et al., 2022). A noteworthy example of such physiological responses is the activation and regulation of antioxidant enzymes under Cd stress, which play essential roles in restoring cellular functions and maintaining homeostasis (Yu et al., 2021; Liu et al., 2023). Our findings indicate that Cd stress induced a significant increase in the content of the lipid peroxidation product MDA suggesting increased lipid oxidative damage in plant cells ( Figure 2A ). Moreover, exposure to Cd stress significantly elevated the activity of antioxidant enzymes (SOD, POD, and CAT) in A. manihot roots in this study ( Figures 2B–D ), consistent with findings in Cd hyperaccumulators such as Youngia japonica (Yu et al., 2021), Pterocypsela laciniata (Zhong et al., 2019), and Arabis paniculate (Liu et al., 2023). Among these enzymes, POD plays pivotal roles in regulating diverse plant physiological processes under Cd stress (Rui et al., 2016; Loix et al., 2017). In the roots of A. manihot, there was no significant difference in POD activity compared to that in the CK group at 12 h, but a significant increase was observed thereafter ( Figure 2C ). This finding aligns with previous research suggesting that POD becomes particularly important in the later stages of exposure to Cd toxicity (Yang et al., 2007; Rui et al., 2016).

In this study, transcriptome analysis revealed significant changes in the transcripts of A. manihot roots from 12 h to 7 days under Cd stress. At the initial stage of Cd stress, only a limited number of genes exhibited differential regulation. Notably, the transcript-level responses of A. manihot to Cd stress became more pronounced during the medium term (36–72 h) ( Figure 3A ). Moreover, more downregulated DEGs were identified in A. manihot roots than upregulated DEGs under Cd stress at different time points ( Figures 3A and 6B ). These downregulated DEGs were mainly enriched in metabolic pathways related to cell wall organization or biogenesis ( Figure 3B , Supplementary Table S4 ). Interestingly, DEGs involved in secondary cell wall synthesis exhibited diverse temporal expression patterns in response to Cd stress. For example, DEGs related to hemicelluloses, such as xylan and glucuronoxylan, decreased after 36 h, while those linked to lignin metabolic processes declined continuously from 12 to 72 h before increasing thereafter ( Figure 4 ). The cell wall is a dynamic structure whose composition rapidly changes in response to Cd stress (Loix et al., 2017). In general, elevated levels of polysaccharide components result in greater accumulation of metal ions within the cell wall, thereby impeding their entry into the protoplast (Guo et al., 2021; Wang et al., 2022). This serves as a significant defense mechanism against trace metal stress across many plant species. However, studies on hyperaccumulators have revealed that the fixation of metal ions by root cell walls is less stable than that in non-hyperaccumulator plants potentially leading to the transportation of metal ions to the shoots (Loix et al., 2017; Guo et al., 2020). For example, the Cd non-hyperaccumulating ecotype of S. alfredii exhibits abundant functional groups and xyloglucan in hemicellulose resulting in tighter Cd binding within the root cell wall. Conversely, the Cd-hyperaccumulating ecotype of S. alfredii shows lower levels of hemicellulose-bound Cd, which allows easier transport into the shoot (Guo et al., 2020). Therefore, during a period of Cd exposure, the dynamic inhibition of the secondary cell wall in A. manihot may facilitate the upward transport of Cd ions.

Similarly, we found that other metabolic processes were also dynamically regulated under Cd stress in A. manihot. For example, during the middle of Cd exposure, the STEM analysis revealed initial increases followed by decreases in the expression levels of DEGs related to lipid and acylglycerol homeostasis, while DEGs in pathways involved in ATP synthesis exhibited decreased expression followed by increased expression ( Figure 4 ). Lipid homeostasis involves a balance between the accumulation of membrane lipids and the accumulation of storage lipids such as triacylglycerol. Previous studies have demonstrated that lipid homeostasis is known to exert systemic effects capable of influencing plant survival, growth, and development during stress (Boutté and Jaillais, 2020). In A. manihot, the dynamic control of lipid homeostasis responds to disruptions in lipid metabolism induced by Cd stress. Taken together, these findings demonstrate that Cd stress dynamically affects the gene expression of A. manihot in a temporally specific manner.

We utilized physiological indicators as trait files to conduct further analysis of these DEGs. By performing WGCNA, we identified eight hub genes that likely play key roles in Cd stress responses ( Figure 5 , Supplementary Table S5 ). Interestingly, among these hub genes, four were homologous genes related to lignin biosynthesis, including three POD genes and one REF1 homolog. In terms of the antioxidant system, POD assists plants in coping with the excess H₂O₂ induced by Cd exposure. In addition to its role in redox homeostasis, this enzyme family participates in various cellular processes, including cell wall loosening, cross-linking, and lignification (Rui et al., 2016; Loix et al., 2017). POD is involved in a complex redox network in A. manihot playing a regulatory role in response to Cd in terms of accumulation and detoxification. Furthermore, a crucial transcription factor has been identified. ERF98, an ethylene response factor, serves as a central regulatory hub in plant responses to abiotic stresses (Wu et al., 2022). Among the identified hub genes, we also identified a homolog of an ABC transporter, a homolog of cytochrome P450 CYP82D47, and a homolog of the RING-H2 finger protein ATL3. ABC transporters, a significant portion of the membrane protein family, play vital roles in the transport and detoxification of heavy metals. In Arabidopsis, AtABCC3 acts as a transporter for PC–Cd complexes influencing Cd accumulation and tolerance (Brunetti et al., 2015). Similarly, OsABCG36 is activated under heavy metal stress serving as one of the root hub genes responsible for exporting Cd or Cd conjugates, thereby increasing Cd tolerance in Oryza sativa (Fu et al., 2019).

Metabolites represent the biochemical end products of gene activity reflecting the adaptability of organisms (Wu et al., 2023). In this study, we observed a trend of increasing changes in the number of DEMs with prolonged Cd stress peaking at 72 h before decreasing thereafter. ( Figure 6B ). The upregulated DEMs were predominantly composed of alkaloids, lipids, flavonoids, phenolic acids, organic acids, amino acids, nucleotides, and their derivatives ( Supplementary Table S7 ). Interestingly, at the early stage of Cd stress, there was a significant increase in the level of abscisic acid potentially triggering intricate stress-adaptive signaling cascades (Shen et al., 2022). Additionally, melatonin was significantly upregulated at all time points ( Figure 6C ) consistent with previous findings (Ahammed et al., 2024). Remarkably, exogenous melatonin application has been demonstrated to enhance Cd tolerance in Solanum lycopersicum through increased antioxidant potential and phytochelatin biosynthesis (Hasan et al., 2015). Therefore, the sustained increase in melatonin levels may have contributed to the enhanced tolerance of A. manihot to Cd stress. Subsequently, KEGG analysis of the metabolome revealed that the response of A. manihot to Cd stress involves various metabolic pathways at different stages ( Supplementary Figure S8 ). These upregulated DEMs were mainly enriched in fatty acid-related pathways indicating the importance of fatty acid metabolism in the response of A. manihot to Cd stress.

In this study, analysis of time-course metabolome and transcriptome data, along with their integrated results, revealed the crucial role of lipid metabolism-related pathways in the defense response of A. manihot to Cd stress ( Figures 4 , 5 , 7 , and S8 ). Similarly, a study on Tamarix hispida reported that lipid synthesis and metabolism play key roles in enhancing Cd tolerance (Xie et al., 2023). Lipids, which are fundamental components of plant cell biofilms, participate in a wide range of biological processes, including growth, development, and stress responses (Hou et al., 2016; Hu et al., 2023; Liang et al., 2023). Furthermore, accumulating evidence suggests a role for lipid molecules, such as lysophospholipids, fatty acids, and phosphatidic acid, in plant signaling processes (Liang et al., 2023).

In this study, the levels of various types of lysophospholipids exhibited diverse patterns of change ( Figure 8 ). LysoPCs [C16:0 (2) and C18:3 (2)] significantly increased under Cd exposure at 12 h. LysoPCs are typically present in trace quantities in plant tissues, but their levels substantially increase under abiotic or biotic stress conditions (Okazaki and Saito, 2014). Previously, it was reported that the amounts of lysoPCs (C14:0, C15:0, C16, C17:1, and C18:3) significantly increased under cold stress (Sun et al., 2021). At an early stage, the increased generation of lysoPCs in A. manihot roots may mediate signaling cascades leading to alterations in gene expression during Cd stress responses (Cappa and Pilon-Smits, 2014). Interestingly, the majority of the detected LysoPEs, specifically LysoPE 18:3 and LysoPE 18:1, decreased consistently at all treatment time points. This result is consistent with findings on T. hispida under Cd stress (Xie et al., 2023). However, there have been few reports regarding the functions of lysophospholipids in response to Cd stress in higher plants.

In addition, we observed that with continuous Cd stimulation, the gene expression of PLA was initially downregulated at 36 h followed by upregulation at 72 h ( Figure 8 ). PLA is critically important because it hydrolyzes the sn-1 and sn-2 positions of glycerophospholipids to produce lysophospholipids and free fatty acids (Hou et al., 2016). Previous studies have indicated a strong correlation between PLA and plant stress defense through its roles in the biosynthesis of JA and oxylipins, as well as the activation of downstream defense (Hou et al., 2016; Hu et al., 2023; Liang et al., 2023). Moreover, genes involved in α-linolenic acid metabolism, including FAD, LOX, AOS, OPR, and β-oxidation genes, exhibited significant upregulation at 72 h under Cd stress corresponding to significant increases in metabolites such as 13-HPOTrE, 9-HOTrE, and stearidonic acid ( Figure 8 ). These findings suggest that phospholipase is activated in the plasma membrane promoting α-linolenic acid metabolism in A. manihot roots. The upregulation of the LOX and AOS genes, along with the increased levels of the crucial intermediate metabolite 13-HPOTrE, implies activation of the downstream JA pathway. Previous studies have indicated that under heavy metal stress, JA regulates antioxidant activity to enhance plant tolerance to heavy metal stress (Hewedy et al., 2023).

Interestingly, stearidonic acid, a nontraditional fatty acid, significantly increased at multiple time points. This may be attributed to the upregulation of FAD genes, which increases the conversion of α-linolenic acid to γ-linolenic and stearidonic acids (Lee et al., 2019). In B. napus, the transcript levels of the FAD genes were greater in the Cd-tolerant cultivar than in the sensitive cultivar under 250 μM Cd stress (Xu et al., 2019). Moreover, FAD genes have been extensively identified in diverse plant species; these genes are activated under various abiotic stress conditions and thereby augment plant stress resistance (Zhang et al., 2005; Yu et al., 2009; Xu et al., 2019). For instance, in tomatoes, LeFAD3 overexpression significantly enhanced resistance to both salt and cold stress (Zhang et al., 2005; Yu et al., 2009). Therefore, we speculate that lipid metabolism maintains the stable cell membrane state in A. manihot under Cd stress and participates in the signaling processes involved in the plant stress response, thereby enhancing the Cd tolerance of plants.

Recent studies have highlighted the important role of phenylpropanoid biosynthesis in mitigating the adverse effects of Cd stress (Chen et al., 2020; Wang et al., 2022; Yu et al., 2023). In our studies, significant increases in the levels of trans-cinnamic acid and ferulic acid were detected at 72 h followed by pronounced increases in p-coumaric acid and sinapic acid levels at 7 days ( Figure 9 ). Similarly, the accumulation of phenolic acids was detected in the roots of Kandelia obovata under Cd stress (Chen et al., 2020). On the one hand, enhancing phenolic acid biosynthesis in response to Cd exposure plays a crucial role in scavenging ROS due to their potent antioxidant activity (Sharma et al., 2019; Chen et al., 2020). On the other hand, phenolic acids are capable of chelating Cd or rendering Cd biologically unavailable, thereby mitigating its toxicity (Chen et al., 2020). It should be noted that the composition and concentration of phenolic acids differ among plants exposed to diverse heavy metal stresses (Anjitha et al., 2021). Such variations may be determinants of the specific defense strategies plants employ to cope with these environmental stresses. Consistent with these findings, structural genes involved in the synthesis of these compounds, including PAL, C4H, COMT, F5H, and 4CL, were significantly upregulated at 72 h or 7 days ( Figure 9 ). PAL is the first enzyme in the general phenylpropanoid pathway, catalyzing the nonoxidative elimination of ammonia to yield trans-cinnamic acid, which is then transformed into p-coumaric acid through a process catalyzed by C4H (Dong and Lin, 2021). Both trans-cinnamic acid and p-coumaric acid serve as precursors for a wide range of organic compounds and can be influenced by the metabolic efficiency of the phenylpropanoid, lignin, and flavonoid pathways (Chen et al., 2020; Dong and Lin, 2021).

In addition, significant reductions in the levels of three primary monolignols, namely, sinapyl alcohol, p-coumaryl alcohol, and coniferyl alcohol, were observed in the later stages of Cd stress ( Figure 9 ). This phenomenon could be associated with the activation of REF1, which enzymatically converts coniferaldehyde and sinapaldehyde into sinapic acid and ferulic acid, respectively, thereby depleting the substrates necessary for monolignol synthesis (Nair et al., 2004). Interestingly, REF1 was identified as a hub gene in the present study ( Figure 5C ). Nevertheless, the precise mechanism by which REF1 influences the phenylpropanoid biosynthesis pathway in response to Cd stress requires further investigation. Additionally, an increase in the POD transcript level might also influence the monolignol content. The expression of the GhPER8 gene, a ligninolytic peroxidase in tobacco leaves, significantly reduced the levels of coniferyl alcohol and sinapic acid, which are substrates for G-lignin and S-lignin biosynthesis, respectively (Gao et al., 2019). In summary, phenylpropanoid biosynthesis, through the production of phenolic acids and the modulation of lignin, plays a vital role in detoxifying A. manihot against Cd.