Integrated transcriptomic and metabolomic analysis reveals the molecular profiles of dynamic variation in Lilium brownii var. viridulum suffering from bulb rot

Lilium brownii var. viridulum, known as Longya lily, is a well-known medicinal and edible plant in China. Bulb rot is a common disease in Longya lily cultivation that severely affects the yield and quality of lilies. According field investigations, we found that different Longya lily plants in the same field had different degrees of resistance to root rot. To find the reasons leading to the difference, we performed metabolomic and transcriptomic analyses of Longya lily with different degrees of disease. The transcriptomic analyses showed that the number of differentially expressed genes increased in early and mid-stage infections (LYBH2 and LYBH3), while decreased in late-stage infection (LYBH4). A total of 2309 DEGs showed the same expression trend in diseased bulb compared healthy bulb (LYBH1). The transcription factors (TFs) analysis of DEGs showed that several common TFs, like WRKY, bHLH, AP2/ERF-ERF and MYB, were significantly activated in bulbs after decay. The metabolomic analyses showed that there were 794 differentially accumulated metabolites, and metabolites with significant changes in relative content largely were phenolic acids, followed by flavonoids and amino acids and derivatives. The combined analysis of transcriptome and metabolome indicated that phenylpropanoid biosynthesis pathway was crucial in Longya lily resistance to bulb rot. Therefore, we speculated that the different degree of resistance to bulb rot in Longya lily may be related to the transcript levels of gene and contents of metabolites in the phenylpropanoid biosynthesis pathway. Overall, these results elucidate the molecular responses of Longya lily to bulb rot and lay a theoretical foundation for breeding resistant varieties.

1 Introduction

Lilium brownii var. viridulum Baker, a perennial bulbous herb, is a mutant variety of L. brownii F. E. Brown ex Miellez within the Lilium genus of the Liliaceae family. Commonly known as the Longya lily, meaning dragon’s teeth in Chinese, this plant is characterized by bulbs made up of many thick scales, with wide bases and pointed, curved tips. Its bulbs contain phenols, saponins, sterols, and alkaloids, which known for their sedative, anti-inflammatory and antitussive properties (Yang et al., 2019). Nevertheless, with the expansion of cultivation area and the increase of continuous cropping years, the occurrence of Longya lilies disease has gradually increased, especially bulb rot. The bulb rot, is commonly caused by soilborne fungal pathogens such as F. oxysporum f. sp. lilii (F. oxysporum), Fusarium solani, F. proliferatum, and curvularia pseudobrachyspora that become an important factor limiting the large-scale planting of Longya lilies (Pan et al., 2003; Li et al., 2018; Ren et al., 2020; Zeng et al., 2020). During the early stage of bulb rot disease, the scales appear brown slightly concave spots, which gradually expand and necrosis. In the later stage, the bulb browns and decays, the scales rot and peel off from the bulb disk, the leaves wilt and eventually the whole plant dies (Jiao et al., 2021). According to the previous statistics, this disease is highly prevalent in Wanzai, Taihe and Gao’an in Jiangxi, with a loss rate of 20%–30% in production, and can reach more than 60% in severe cases, which causes huge economic losses. The disease is difficult to prevent and control due to it is affected by complex environmental factors, such as temperature for pathogen growth, moderate to high soil moisture content, soil compaction, poor drainage, and so on. To alleviate losses caused by bulb rot, chemical pesticides are used many times and in large quantities, thus making pathogenic bacteria resistant and polluting the environment (Liao et al., 2013; Peng et al., 2015). Therefore, understanding the molecular mechanism of disease resistance, and breeding resistant varieties are the most economical and effective measures to prevent bulb rot.

Plant responses to pathogen infection trigger significant changes in gene expression and metabolism (Lu et al., 2020; Shen et al., 2023; Zhang et al., 2023). Transcriptome analysis, facilitated by RNA-Seq technology, allows for gene expression profiling in pathogen-infected plants, particularly those lacking a reference genome (Sun et al., 2019; Zhang et al., 2021a). This approach serves as a valuable tool to investigate how plants respond to various stresses, both abiotic and biotic, shedding light on the molecular mechanisms underlying their ability to withstand such conditions (Sun et al., 2019; Lu et al., 2020; Zhang et al., 2021a; Shen et al., 2023; Zhang et al., 2023). Additionally, metabolomics, which focuses on small molecular components, the combined transcriptome and metabolome analysis is widely used to uncover molecular mechanisms of plant response to biological or abiotic stresses. For instance, in resistant eggplant varieties, the metabolisms of phytohormones, phenylpropanoids, and flavonoids were found to be differentially enriched, highlighting their crucial role in defense against bacterial wilt (Xiao et al., 2023). Similarly, resistant broccoli showed significant increases in three indolic glucosinolates (I3G, 4MI3G, and 1MI3G) following infection with Alternaria brassicicola, underscoring the importance of these compounds in defense mechanisms (Shen et al., 2023). Flavonoid metabolism was identified as a key player in conferring resistance to stem canker in Zanthoxylum bungeanum (Li et al., 2021). Previous research on bulb rot disease primarily focused on pathogen isolation and identification, and the pathogenic mechanism of the pathogen, and so on (Pan et al., 2003; Ai et al., 2017; Li et al., 2018; Ren et al., 2020; Zeng et al., 2020). However, the impact of bulb rot on Longya lily at the transcriptional and metabolic levels, as well as the underlying response mechanisms, remains unclear.

In this study, we collected Longya lily bulbs with different degrees of disease from the same field, and found the key genes and metabolites associated with disease response through the combined transcriptome and metabolome analysis. Furthermore, the potential resistance mechanisms was speculated through the function of key genes and metabolites. These identified genes and metabolites can not only serve as valuable references for early disease diagnosis of Longya lily bulbs, but also provide theoretical basis for breeding disease-resistant varieties.

2 Materials and methods

2.1 Sample collection of Longya lily bulbs

The Longya lily used in this study was harvested in 2023 from Baishui township, located in the northwest of Wanzai County, Jiangxi Province, China, a widely cultivated species in the region. This field has been naturally infected with bulb rot disease and has since spread to large areas. In present research, we selected underground bulbs of Longya lily as the experimental materials. Four disease levels were divided in the experiment: LYBH1 (healthy: no obvious brown spots), LYBH2 (mild: brown spots occupied about 1/4 of the entire bulb), LYBH3 (moderate: brown spots occupied about 1/2 of the entire bulb), LYBH4 (severe: approximately the entire bulb turned brown and rotted). The textural and sensory were evaluated immediately after harvest, thereby defining the degree of disease. The surface of each scale was washed with clean water and sterilized with 75% ethanol, and then rinsed three times with sterile water. All fresh scales were cut into small pieces and immediately frozen in liquid nitrogen. Three biological replicates were collected for each level. Each replicate contained mixed sampled bulbs from at least five sampled sites. All bulb samples were sealed and sent for transcriptomic and metabolomic analyses.

2.2 RNA extraction and library construction

Total RNA was extracted using the Spin Column Plant Total RNA Purifcation Kit (Sangon Biotech, Shanghai, China), following the manufacturer’s instructions. The RNA quality was assessed by analyzing integrity with 1% agarose gel electrophoresis and the Agilent 2100 bioanalyzer, and concentration by the Qubit 2.0 Fluorometer. For cDNA library construction, total RNA ≥1 µg was used as a template with the Illumina NEBNext^® UltraTM RNA Library Prep Kit (NEB, Ipswich, MA, United States). Library quality was evaluated on the Agilent Bioanalyzer 2100 (Agilent, Santa Clara, CA, United States). Once qualified, all cDNA libraries were sequenced on the Illumina Hiseq 2000 platform (San Diego, CA, United States).

2.3 Reads mapping and functional annotation

Raw reads underwent quality control to obtain clean data by removing adaptors, reads with N content >10%, and low-quality bases (Q ≤ 20) exceeding 50%. High-quality reads were further filtered based on sequencing error rate and GC content distribution, then mapped to the assembled transcriptome for gene expression quantification in the absence of reference genomes. Clean reads were assembled into a reference sequence using Trinity (version 2.13.2) (Friedman et al., 2011).

Gene expression levels were quantified as fragments per kilobase per transcript per million mapped reads (FPKM) after normalization. Gene functions were annotated using databases such as the Kyoto Encyclopedia of Genes and Genomes (KEGG), Gene Ontology (GO), and NCBI non-redundant protein sequences (Nr). Differentially expressed genes (DEGs) were identified between groups using DESeq2 software (version 1.22.2) with criteria of |Log₂fold change (FC)| ≥ 1 and false discovery rate (FDR) < 0.05 (Love et al., 2014; Varet et al., 2016). Subsequently, DEGs were further analyzed through clustering using the R software package (pheatmap: version 1.0.12). Gene functions were annotated using databases such as Kyoto Encyclopedia of Genes and Genomes (KEGG), Gene Ontology (GO) and NCBI non-redundant protein sequences (Nr). PlnTFDB (Perez-Rodriguez et al., 2010) and PlantTFDB (Jin et al., 2014) database was integrated to predict plant transcription factors (TFs) using iTAK software (version 1.7a) (Zheng et al., 2016).

2.4 Metabolites extraction

The freeze-dried samples were crushed using a mixer mill (MM 400, Retsch Technology, Haan, Ger-many) with a zirconia bead for 1.5 min at 30 Hz. The 50 mg of lyophilized powder was dissolved in 1.2 mL of 70% methanol solution and vortexed for 30 s at 30 min intervals for a total of 6 times, then centrifuged (12,000 rpm for 3 min). The supernatant was aspirated, and the sample was filtered through a microporous membrane (0.2 µm pore size) and stored in the sample vial for UPLC-MS/MS analysis.

2.5 Quantitative analysis of metabolomics

Sample extracts were analyzed using ultra-performance liquid chromatography electrospray ionization tandem mass spectrometry (UPLC-ESI-MS/MS) system (UPLC, Shim-pack UFLC SHIMADZU CBM30A system; MS/MS, Applied Biosystems 6500 QTRAP). Chromatographic separation was conducted on an Agilent SB-C18 column (1.8 µm, 2.1 mm × 100 mm, Agilent, New York, NY, United States) at 40°C. The UPLC analysis utilized a solvent system of ultrapure water (0.1% formic acid) and acetonitrile (0.1% formic acid) with a gradient program of 95:5 v/v at 0.0 min, transitioning to 5:95 v/v from 0.0 to 9.0 min, maintaining 5:95 v/v from 9.0 to 10.0 min, and returning to 95:5 v/v from 10.0 to 11.1 min, continuing until 14.0 min. The flow rate was set at 0.35 mL min⁻¹, temperature at 40°C, and injection volume at 2 µL. The effluent was directed to an ESI-triple quadrupole-linear ion trap (QTRAP)-MS for further analysis. The ESI source operated at a source temperature of 550°C, ion spray voltage of 5,500 V (positive mode)/- 4,500 V (negative mode), and gas pressures of 50, 60, and 25 psi for ion source gas I, gas II, and curtain gas, respectively. Mass spectrometry data processing was performed using Analyst software (Version 1.6.3, AB Sciex, Fram-ingham, MA, United States). Metabolite analysis involved qualitative assessment based on secondary spectral information from the self-built database MWDB (Metware database) and quantitative analysis using multiple reaction monitoring (MRM) modes of triple quadrupole mass spectrometry.

2.6 Statistical analysis

Multivariate statistical analysis, including principal component analysis (PCA) and orthogonal partial least squares discriminant analysis (OPLS-DA), was performed using the R software package. The OPLS-DA models were evaluated based on their R2X, R2Y and Q2 values (R2: interpretation rate of the model, Q2: predictability of the model). The PCA plots were used to analyze the variability between groups and within groups. The OPLS-DA was performed to maximize the metabolome differences between sample pairs. The relative importance of each metabolite to the OPLS-DA model was tested using the variable importance in projection (VIP) as a parameter. For the VIP (VIP >1) value from OPLS-DA model, combined with absolute Log₂FC (|Log₂FC| ≥ 1.0) were used as the thresholds for screening differentially accumulated metabolites (DAMs). The functions of DAMs were further annotated using the KEGG compound database.

3 Results

3.1 Morphology of Longya lily bulbs with different degrees of disease

According to the field observations, we discovered that the Longya lily plants and bulbs with different disease levels had large phenotypic differences. The healthy plants were morphologically greener in color and exhibited normal growth, with no disease spots on the bulbs. In the early stages of the occurrence of the bulb rot disease, the edges of some leaves began to yellow, but the plant still grew normally, and a small number of brown spots appeared on the outer scales of the bulbs. With the development of diseases, a large number of leaves yellowing expanded and eventually covered the entire leaf, the surface of the bulbs had several sunken dark brown lesions and multiple lesions could merge with each other to appear water-stained. In the late stage of the disease, the leaves of the whole plant turned yellow or purple from bottom to up, and then wilted and dried up, and the base of the stem shrank, which was easy to break; the color of the lesions on the bulbs gradually darkened and deepened into the bulb tissue, creating a cavity, and then the bulbs completely deteriorated and formed mushy rot, eventually lead to death of the whole plants (Figures 1A, B). These performances were consistent with previous reports showing that lily infected bulb rot disease caused by F. oxysporum.

Figure 1

Figure 1. (A) Growth of healthy and diseased Longya lilies in the field. (B) Each disease level of Longya lily bulbs. LYBH1: Healthy; LYBH2: Mild; LYBH3: Moderate; LYBH4: Severe. Scale bars = 1.5 cm.

3.2 RNA sequencing and quality control

Based on quality assessment of the RNA extraction, the extracted concentrations were between 10 and 32 ng μL⁻¹, and the integrity were from 6.3 to 8.6, indicating good RNA quality, which met the requirements for sequencing and library construction (Supplementary Table S1). A total of 99,711,225 raw reads were obtained by transcriptome sequencing, and 96,627,540 clean reads were derived after filtering the original data; altogether, 173.92 Gb of clean data was acquired. Each sample had over 12 Gb of clean data and the GC content was about 50.38%. The percentage of Q20 base was about 97.92%, and the Q30 base was about 94.09%, indicating the relatively high quality of transcriptome data (Supplementary Table S2). Ultimately, the sequence and expression information of 127,621 unigenes was obtained for further data analysis.

3.3 Differentially expressed genes analysis

To investigate the gene expression changes induced by bulb rot disease in Longya lily, all samples were divided into three comparison groups: LYBH2 vs. LYBH1, LYBH3 vs. LYBH1 and LYBH4 vs. LYBH1. The differentially expressed genes (DEGs) for different comparisons were identified by DESeq based on FPKM values. A total of 9,976 DEGs were identified in three comparison groups. In LYBH2 vs. LYBH1, there were 3,922 DEGs (1,861 upregulated and 2,061 downregulated). Similarly, there were 7,595 DEGs (4,068 upregulated and 3,527 downregulated) were identified in LYBH3 vs. LYBH1, whereas the LYBH4 vs. LYBH1 group contained 6,590 DEGs (3,622 upregulated and 2968 downregulated) (Figure 2A). Venn diagram analysis showed that 2,309 (986 upregulated and 1,323 downregulated) DEGs were shared in three comparison groups, 996, 1,955 and 1,203 DEGs were unique to the LYBH2 vs. LYBH1, LYBH3 vs. LYBH1 and LYBH4 vs. LYBH1, respectively (Figure 2B). This result indicated that there was a significant change in gene expression levels between the healthy and infected Longya lily bulbs.

Figure 2

Figure 2. (A) Number of up/downregulated DEGs in the three comparison groups. (B) Venn diagram showing the DEGs numbers in each comparison group.

To evaluate the function of the shared DEGs between the three comparison groups, functional enrichment analyses using GO and KEGG database were carried out to identify possible biological processes or pathways involved in bulb rot disease infection response. All DEGs annotated in the GO database were divided into three categories: biological process (the top three GO terms were cellular process, metabolic process, and response to stimulus), molecular function (catalytic activity and binding, transcription regulator activity), and cellular component (cellular anatomical entity and protein-containing complex (Figure 3A). Furthermore, all DEGs were mapped to the KEGG pathways. According to the annotation results, we discovered that most DEGs were annotated to metabolic pathway, biosynthesis of secondary metabolites, MAPK signaling pathway-plant, plant hormone signal transduction, phenylpropanoid biosynthesis, starch and sucrose metabolism and plant-pathogen interaction pathways. The top 20 pathways with the lowest p-values showed that flavonoid biosynthesis, MAPK signaling pathway-plant, flavone and flavonol biosynthesis, phenylpropanoid biosynthesis, isoquinoline alkaloid biosynthesis, plant hormone signal transduction, biosynthesis of various alkaloids, starch and sucrose metabolism, tyrosine metabolism, mismatch repair and circadian rhythm-plant pathways were significantly enriched (Figure 3B), some of which have been previously reported to be related to plant disease resistance, like flavonoid biosynthesis, phenylpropanoid biosynthesis, plant hormone signal transduction, starch and sucrose metabolism, and so on. In addition, we performed KEGG enrichment analysis for DEGs unique to each group (Supplementary Table S3). In LYBH2 vs. LYBH1, a total of 113 pathways were enriched by specific DEGs, of which the glutathione metabolism and monoterpenoid biosynthesis pathways were significantly enriched. In LYBH3 vs. LYBH1, the number of pathways for specific DEGs enrichment was 124, and the phenylpropanoid biosynthesis, cyanoamino acid metabolism and plant hormone signal transduction pathways were significantly enriched. Correspondingly, a total of 109 pathways were enriched, and the protein processing in endoplasmic reticulum and carotenoid biosynthesis were significantly enriched among them. This indicated that the Longya lily in the different stage of infection defended against pathogen invasion by activating different metabolic pathways.

Figure 3

Figure 3. (A) GO enrichment analysis of DEGs from the three comparisons. (B) Top 20 enriched KEGG pathways of DEGs shared in the three comparison groups.

3.4 Transcription factors analysis of DEGs

Transcription factors (TFs) and transcriptional regulators (TRs) play an important regulatory role in plant growth, development and stress defence response. Therefore, we also analysed the TFs (TRs) in DEGs. In total, 52 TFs (TRs) from 270 DEGs were identified in the comparisons of LYBH2 vs. LYBH1; 59 TFs (TRs) from 524 DEGs were observed in the comparisons of LYBH3 vs. LYBH1, while 56 TFs (TRs) from 425 DEGs were found in the comparisons of LYBH4 vs. LYBH1. 46 common differentially expressed TFs (TRs) were identified in the three comparisons, such as WRKY, NAC, MYB, AP2/ERF-ERF, bHLH, bZIP, Tify, etc. These TFs (TRs) families have previously been reported to be involved in the regulation of various disease resistance pathways. Among them, the first two families with more TFs were bHLH and AP2/ERF-ERF family, which were all largely upregulated in the three comparisons. The C2C2-CO-like, HB-other, MADS-M-type and zf-HD TFs encoded by downregulated DEGs were unique in the LYBH2 vs. LYBH1 and LYBH3 vs. LYBH1 comparison groups. Some TFs (TRs) were only in the LYBH3 vs. LYBH1 and LYBH4 vs. LYBH1 com-parison groups, like AP2/ERF-RAV, BES1, C2C2-YABBY, E2F-DP, HB-KNOX and NF-YC. Also, CAMTA and TUB, were downregulated. Similarly, the RAF, LFY and DBP TFs (TRs) were downregulated only in LYBH3 compared to LYBH1, while the BBR-BPC and HMG, both upregulated and downregulated, were just discovered in the LYBH4 vs. LYBH1 comparison group. All data for the TFs can be found in Supplementary Table S4.

3.5 Overview of the metabolome

To gain insight into the metabolic changes and the possible defense mechanisms in the Longya lily bulbs after decay, widely targeted metabolomics by UPLC-MS/MS was performed to broadly evaluate metabolite profiles in healthy and diseased Longya lily bulbs. Based on UPLC-MS/MS analysis and a self-building database, a total of 1,603 metabolites were detected in all samples and divided into 11 categories. Among them, five categories account for more than 10% of the metabolites, including phenolic acids (13.91%), alkaloids (12.04%), flavonoid (11.04%), amino acids and derivatives (10.92%) and steroids (10.04%) (Figure. 4A). Principal component analysis (PCA) analysis of metabolomic profiles was conducted to provide information about distinctly separate groups among these samples of the four disease levels of Longya lily bulb and the three replicates of each disease level together. The two PC together explained 70.81% of the chemical variables. Based on PC1, 62.44% of the information could clearly distinguish LYBH1, LYBH2 and LYBH3, LYBH4 samples. The two-dimensional plot showed obvious variations in the metabolites from the bulbs of Longya lily after decay. However, the LYBH3 and LYBH4 samples had similar metabolic profiles perhaps due to the relatively similar severity of onset. Furthermore, three duplicate samples of each disease level grouped together, indicating high reproducibility among replicates (Figure 4B).

Figure 4

Figure 4. (A) Ring diagram of metabolites classification. (B) Principal component analysis (PCA) of all metabolites.

3.5.1 Identification of differentially accumulated metabolites

Based on VIP (VIP >1) and absolute Log₂FC (|Log₂FC | ≥ 1.0), there were 794 differentially accumulated metabolites (DAMs) were identified in three comparisons, of which 37 DAMs only accumulated in infected lily bulbs (LYBH2, LYBH3 and LYBH4), including 10 flavonoids, 5 lignans and coumarins, 4 steroids, 4 alkaloids, 3 phenolic acids, 2 terpenoids, 2 quinones, 1 organic acid, 1 nucleotides derivative and 5 other compounds. In LYBH2 vs. LYBH1, there were 440 DAMs (406 upregulated and 34 downregulated). Correspondingly, there were 646 DAMs (603 upregulated and 43 downregulated) in LYBH3 vs. LYBH1, and 668 DAMs (621 upregulated and 47 downregulated) in LYBH4 vs. LYBH1 (Figure 5A). And the most upregulated metabolites were phenolic acids, followed by lipids (Table 1). This indicated that phenolic compounds played key roles in in lily bulb rot resistance. In addition, a total of 363 DAMs were shared in three comparison groups, the number of DAMs that were unique to the comparison group of LYBH2 vs. LYBH1, LYBH3 vs. LYBH1 and LYBH4 vs. LYBH1 were 38, 74 and 85, respectively (Figure 5B).

Figure 5

Figure 5. (A) Number of differentially accumulated metabolites detected in the three comparison groups. (B) Venn diagram showing DAMs numbers shared in three comparisons.

Table 1

Table 1. The detailed amount of difference metabolites in the three comparisons of lily.

In this study, the top 10 DAMs with the largest absolute value of Log₂FC were screened in three comparison groups, respectively (Supplementary Table S5). Among them, the upregulated phenolic acids (p-Coumaroylmalic acid and 1-O-p-Coumaroyl-β-D-glucose), amino acids derivatives (N-[malonyl] phenylalanine), and other substances (juglanoside D) were shared in the three comparison groups. Apart from this, the top 10 DAMs in LYBH2 vs. LYBH1 group included three upregulated phenolic acids (1-Feruloyl-3-p-coumaroyl-2′,6′-diacetylsucrose, vanilloloside and 6-O-p-Coumaroyl-β-D-glucose*), one upregulated flavonoid (kaempferol-3-O-rutinoside-7-O-rhamnoside), one upregulated quinone (7-hydroxy-4-methoxyphenanthrene-2-O-glucoside) and one upregulated alkaloid (N′,N″-Diferuloylspermidine); the top 10 DAMs in LYBH3 vs. LYBH1 had two upregulated amino acids derivatives (jasmonoyl-L-isoleucine and 2′-deoxyadenosine), one upregulated flavonoid (choerospondin), one upregulated lipid (linoleoyl ethanolamine) and one upregulated alkaloid (N′, N″-Diferuloylspermidine), and one downregulated flavonoid (quercetin-3-O-glucoside (isoquercitrin)); the top 10 DAMs in LYBH4 vs. LYBH1 included three upregulated flavonoids (galangin-7-O-glucoside, apigenin-5-O-glucoside and choerospondin), one upregulated phenolic acid (ethyl ferulate) and one upregulated lipid (linoleoyl ethanolamine), and one downregulated flavonoid (quercetin-3-O-galactoside (hyperin)). Most of these compounds metioned above belong to phenolic acids, followed by flavonoids and amino acids and derivatives, and a few belong to lipids and alkaloids.

3.5.2 KEGG enrichment analysis of DAMs

Based on KEGG annotations of the shared DAMs from the three comparisons, we discovered that a total of 63 pathways were annotated by all shared metabolites, the pathways with the highest number of metabolites were metabolic pathways, biosynthesis of secondary metabolites, biosynthesis of cofactors, linoleic acid metabolism, ABC transporters, biosynthesis of amino acids, alpha-linolenic acid metabolism, aminoacyl-tRNA biosynthesis, phenylpropanoid biosynthesis and purine metabolism. The KEGG functional enrichment analysis showed that linoleic acid metabolism, alpha-linolenic acid metabolism, phenylpropanoid biosynthesis, plant hormone signal transduction, glutathione metabolism, cyanoamino acid metabolism and ABC transporters were significantly enriched, of which linoleic acid metabolism, alpha-linolenic acid metabolism and phenylpropanoid biosynthesis were extremely significantly enriched (p < 0.01) (Figure 6). Thus, we speculate DEMs involved in these metabolic pathways were essential in the Longya lily response to bulb rot.

Figure 6

Figure 6. Top 20 enriched KEGG pathways of DAMs shared in the three comparison groups.

3.6 Integrated transcriptomics and metabolomics analysis

To further explore the relationship between DEGs and DAMs in Longya lily responsive to bulb rot disease, the co-joint KEGG pathway enrichment analysis of the shared 2309 DEGs and 363 DAMs in three comparison groups was conducted. Results showed that these DEGs and DAMs enriched to 115 and 62 pathways, respectively, whereas 51 pathways were shared by DEGs and DAMs. Then, we showed the 20 pathways most significantly enriched in terms of DEGs and DEMs of the pathways shared in three comparison groups with the lowest p-value (Figure 7). Among them, we found that linoleic acid metabolism, alpha-linolenic acid metabolism and phenylpropanoid biosynthesis were significantly enriched by DAMs, but the pathways significantly enriched by DEGs were phenylpropanoid biosynthesis, plant hormone signal transduction, ABC transporters, stilbenoid, diarylheptanoid and gingerol biosynthesis and tyrosine metabolism. Interestingly, only phenylpropanoid biosynthesis was significantly enriched in both the transcriptomic and metabolomic data. This result suggested that the phenylpropanoid biosynthesis pathway may be the main regulatory pathway involved in the Longya lily defenses against bulb rot disease. A network of the diagram was constructed to further analyze the relationship between DAMs and DEGs in the phenylpropanoid biosynthesis pathways (Figure 8). We found that there were five upregulated DAMs (L-tyrosine, trans-5-O-[p-coumaroyl] shikimate, 1-O-sinapoyl-β-D-glucose, sinapinaldehyde and coniferin), and 55DEGs involved in this pathway. Furthermore, we calculated Pearson correlations between DEGs and DAMs involved in phenylpropanoid biosynthesis pathway, retaining only strong (|r| (correlation coefficient) ≥ 0.8) and statistically significant (p < 0.05) associations (Supplementary Table S6). The results showed that the cinnamate 4-Hydroxylase (CYP73A) encoded by the cluster-65569.4 and cluster-65569.10 was negatively correlated with L-tyrosine, trans-5-O-(p-coumaroyl) shikimate, 1-O-sinapoyl-β-D-glucose, sinapinaldehyde and coniferin, while the CYP73A encoded by the cluster-76247.11 and cluster-76247.15 was positively correlated with trans-5-O-(p-coumaroyl) shikimate, 1-O-sinapoyl-β-D-glucose, sinapinaldehyde and coniferin. The shikimate O-hydroxycinnamoyl transferase (HCT) encoded by the cluster-34314.1 and cluster-57950.0 was negatively correlated with the sinapinaldehyde, trans-5-O-(p Coumaroyl) shikimate and 1-O-Sinapoyl-β-D-glucose, and the cinnamyl-alcohol dehydrogenase (CAD) encoded by the cluster-47087.0 was negatively correlated with the 1-O-Sinapoyl-β-D-glucose. The phenylalanine ammonia-lyase (PAL), 4-coumarate-CoA ligase (4CL), caffeoyl-CoA O-methyltransferase (CCoAOMT), cinnamoyl-CoA reductase (CCR), caffeic acid 3-O-methyltransferase (COMT), caffeoylshikimate esterase (CSE) and peroxidase (POD) encoded by the remaining DEGs were positively the five DAMs involved in the pathway. In addition, transcriptome analysis showed that the expression of some genes encoding the CYP73A, 4CL, HCT, CCR and CAD decreased with the increase of disease severity, like the Cluster-65569.1, Cluster-65569.10, Cluster-65569.4, Cluster-71471.0, Cluster-34314.1, Cluster-57950.0, Cluster-47087.0 and Cluster-9648.19, while the expression of the genes encoding the PAL, CCoAOMT, CSE, COMT and POD increased.

Figure 7

Figure 7. Co-enrichment analysis of DAMs and DEGs shared in the three comparison groups.

Figure 8

Figure 8. The schematic diagram of key genes and metabolites in the phenylpropanoid biosynthesis involved in Longya lily the response to bulb rot. The words in bold black indicated the DAMs in each sample. The words in red indicated the DEGs in each sample.

4 Discussion

It has been reported that the bulb rot, mainly caused by F. oxysporum, is a common issue in lily cultivation and production, reducing harvest yield and quality, and subsequent death of severely infected plants (Li et al., 2018). This is a big challenge considering the multi-pathogenic cause of this disease and poor technologies for its management. Currently, this disease is controlled by a variety of methods, such as chemical pesticides, disease management in the field, plant breeding for resistant varieties/genetically modified plants, as well as biological control agents (Liao et al., 2013; Peng et al., 2015; Wu et al., 2015). With the rapid development of science and technology, the breeding of resistant varieties is considered the most economical and environmentally friendly long-term strategy for plant disease management. Therefore, it will become increasingly important to explore internal resistance genes to improve lily resistance. In this study, we performed transcriptome and metabolome analysis of Longya lily with different degrees of disease. Transcriptome analysis showed that there were 3922 DEGs in LYBH2 vs. LYBH1, 7595 DEGs in LYBH3 vs. LYBH1, 6590 DEGs in LYBH4 vs. LYBH1, and 2309 DEGs shared in three comparison groups, indicating complex and contrasting transcriptional changes in Longya lily induced by bulb rot disease. At the same time, we could see from the above data that the numbers of DEGs increased during early and mid-stage infections (LYBH2 and LYBH3), while decreased during late-stage infection (LYBH4). The reason for the results may be that plants required more transcriptional regulatory genes expression in response to pathogen invasion and spread during the early and mid-stages, while plants already began to resist the pathogens attack or enter the stage of host death in the late-stage infection, leading to reduced gene expression. TFs are key regulators of plant responses to environmental signals, directly influencing target genes by binding to specific cis-regulatory elements in their promoters (Riechmann, 2018). Therefore, we further performed TFs analysis of all DEGs. The results showed that a total of 52, 59 and 56 TF families annotated by DEGs identified in LYBH2 vs. LYBH1, LYBH3 vs. LYBH1 and LYBH4 vs. LYBH1, respectively. The most common TFs detected included AP2/ERF-ERF, bHLH, NAC, MYB, Tify, and WRKY. The WRKY family is one of the largest families of specific TFs that regulate various defense processes and play crucial roles in regulating defense-related genes by binding to W-Box cis-elements (Rushton et al., 1996; Rushton et al., 2010). Certain WRKY members have been shown to be important in lily plants defense against F. oxysporum. For instance, Deng et al. reported that the LrWRKY11 positively regulated Lilium regale resistance to F. oxysporum by interacting with salicylic acid/jasmonic acid signaling pathways and modulating LrDIR1 expression to increase lignin/lignans accumulation (Deng et al., 2023). Additionally, Karre. et al. identified that silencing the HvWRKY23 gene increased the proportion of diseased spikelets and fungal biomass (Karre et al., 2019). Our study identified bHLH and AP2/ERF-ERF as first two families with more TFs, with significant upregulation in all three comparisons. bHLH, the second largest TF family, is known for its role in transcription and protein-interacting regulation in plants (Jin et al., 2021). Similarly, AP2/ERF TFs are crucial in regulating various disease resistance pathways (Licausi et al., 2013). The AP2/ERF family can be categorized into ERF, DREB, AP2, and RAV subfamilies based on conserved DNA binding domains, with ERF specifically involved in plant defence responses to pathogens via SA and ethylene (ET)/JA-dependent signal transduction pathways (Ma R. et al., 2017). In F. oxysporum-inoculated lilies, genes encoding WRKY and ERF families showed increased expression levels (He et al., 2019). Consequently, it was suggested that the differential regulation of these TFs may trigger the production of antimicrobial compounds as a defence mechanism against pathogens.

The changes of this transcriptome were related to the metabolome during the infection and spread of pathogenic bacteria. To prevent the invasion of pathogens, plants will mobilize a large number of defense factors, especially the second metabolites, such as reactive oxygen species and phytoalexins (Kumudini et al., 2018). Phytoalexins are primarily synthesized through the phenylpropanoid metabolism pathway, including flavonoids, terpenoids, alkaloids, and phenolic acids (Ahuja et al., 2012). Metabolomics results indicated that the relative content of two phenolic acids (p-Coumaroylmalic acid and 1-O-p-Coumaroyl-β-D-glucose), N-(malonyl) phenylalanine and juglanoside D continuously increased with the severity of the disease. In addition, some metabolites were identified only in one of the stages of infection, such as the kaempferol-3-O-rutinoside-7-O-rhamnoside was especially discovered in LYBH2. Ma et al. found that the kaempferol 3-O-rutinoside, a flavonoid and derivative of quercetin, had certain anti-bacterial effect (Ma Y. et al., 2017). The quercetin-3-O-glucoside (isoquercitrin) was only in LYBH3 and downregulated, which was consistent with the results of previous studies (Wei et al., 2021). The quercetin-3-O-galactoside (hyperin) and apigenin-5-O-glucoside were unique to LYBH4. Previous studies have shown that flavonoids including quercetin, apigenin, rutin and naringenin are known for their antiviral properties, their mechanism of action involves inhibiting viral neuraminidase, proteases, and DNA/RNA polymerases that interact with viral proteins (Ninfali et al., 2020). Additionally, they can enhance plant tissue responses to pathogenic microorganisms, thereby improving disease resistance (Gatto et al., 2011). Therefore, the changes in the relative content of these metabolites may lead to enhance plant disease resistance.

The combined transcriptomic and metabolomic analyses revealed that phenylpropanoid biosynthesis was significantly enriched in both the transcriptomic and metabolomic data, indicating this pathway was crucial in Longya lily response to bulb rot. In response to pathogen infection, there were 55 DEGs involved in this pathway, and found these DEGs were involved in the synthesis of various important proteins, such as PAL, CSE, CCoAOMT, COMT and POD, which were upregulated in LYBH2, LYBH3 and LYBH4. However, some DEGs encoding the cinnamate 4-Hydroxylase (CYP73A), 4-coumarate-CoA ligase (4CL), shikimate O-hydroxycinnamoyl transferase (HCT), cinnamoyl-CoA reductase (CCR) and cinnamyl-alcohol dehydrogenase (CAD) were downregulated in LYBH2, LYBH3 and LYBH4. Phenylpropanoid metabolism is a critical secondary pathway in plants, giving rise to flavonoids, lignans, cinnamic acids, amides, and other metabolites. This pathway initiates with key enzymes PAL, CYP73A and 4CL, leading to the production of cinnamic acid and p-coumaroyl-CoA, which further branches into the phenylalanine metabolic pathway and the flavonoid metabolic pathway (Dixon et al., 2002). The PAL, CYP73A and 4CL enzymes play a significant role in the synthesis of phenolic compounds involved in plant defense mechanisms. And the activities of these enzymes aligned closely with the expression of their respective genes, potentially leading to consistent accumulation of the catalytic products. The POD is involved in lignin synthesis, which provides a structural barrier for the plant against the pathogen (Sun et al., 2024). And the increase of its activity can promote the oxidation of phenol to quinone, which is harmful to pathogenic bacteria (Lu et al., 2020). Phenylpropanoid metabolites play a key role in both plant development and interactions with the environment (Yadav et al., 2020). In this study, we discovered five upregulated metabolites associated with bulb rot in Longya lily, including L-tyrosine, trans-5-O-(p-Coumaroyl) shikimate, sinapaldehyde, 1-O-sinapoyl-β-D-glucose and coniferin. L-tyrosine is the precursor of lignin, which can prevent the penetration of bacteria through strengthening the cell wall, affecting the formation of the haustorium within the cell (Kutschera et al., 2019). Coniferyl alcohol is converted to coniferin by UGT72E. Coniferyl alcohol is one of the major monomers of lignin, and coniferin is a plausible monolignol precursor in the lignification of plant cell walls (Guan et al., 2022), and one study was found to inhibit fungal growth and melanization (Kumar et al., 2020). Interestingly, despite significant changes in gene expression profiles, there were no notable differences in metabolite levels between healthy and infected bulbs, suggesting that environmental factors, structural genes, and TFs may play a crucial role in the biosynthesis and accumulation of metabolites.

5 Conclusion

We studied the variations of gene expression and metabolite in Longya lily bulbs with different degrees of disease using transcriptome and metabolomics analyses. The results indicated differential expression of genes during disease progression, influencing the composition and content of certain DAMs like quercetin, apigenin, kaempferol, and coniferin. The combined KEGG enrichment results indicated that phenylpropanoid biosynthesis pathways may play essential roles in fighting against bulb rot disease. This study lays the foundation for understanding the molecular mechanisms of Longya lily’s resistance to bulb rot disease and offers insights for breeding resistant lily varieties. The potential roles of the identified key regulatory genes and metabolites in disease resistance outlined in this study may need further exploration and validation.

Data availability statement

The datasets presented in the study are publicly available. This data can be found here: NCBI SRA repository under the project number PRJNA1143388 (https://www.ncbi.nlm.nih.gov/sra/PRJNA1143388).

Author contributions

NC: Formal Analysis, Funding acquisition, Investigation, Writing–original draft. LZ: Formal Analysis, Investigation, Writing–original draft. YX: Data curation, Investigation, Writing–original draft. CW: Data curation, Investigation, Writing–original draft. HL: Project administration, Supervision, Writing–review and editing. YW: Supervision, Writing–review and editing.

Funding

The author(s) declare that financial support was received for the research, authorship, and/or publication of this article. This research was supported by the General project of Jiangxi Provincial Traditional Chinese Medicine Science and Technology Plan (2023B1364), the Fundamental Research Funds for the Central Public Welfare Research Institutes (ZZ16-ND-12, ZZ17-ND-12).

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/fgene.2024.1432997/full#supplementary-material

References

Ahuja, I., Kissen, R., and Bones, A. M. (2012). Phytoalexins in defense against pathogens. Trends Plant Sci. 17 (2), 73–90. doi:10.1016/j.tplants.2011.11.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Ai, C., Hui, J., Wang, J., and Ma, D. (2017). “Study on the activity of cell wall degrading enzymes in Fusarium oxysporum,” in The 2017 annual meeting of the chinese society of plant pathology.

Google Scholar

Deng, J., Wang, Z., Li, W., Chen, X., and Liu, D. (2023). WRKY11 upregulated dirigent expression to enhance lignin/lignans accumulation in Lilium regale Wilson during response to Fusarium wilt. J. Integr. Agric. doi:10.1016/j.jia.2023.07.032

CrossRef Full Text | Google Scholar

Dixon, R. A., Achnine, L., Kota, P., Liu, C. J., Reddy, M. S., and Wang, L. (2002). The phenylpropanoid pathway and plant defence-a genomics perspective. Mol. Plant Pathol. 3 (5), 371–390. doi:10.1046/j.1364-3703.2002.00131.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Friedman, N., Regev, A., Grabherr, M. G., Haas, B. J., Yassour, M., Levin, J. Z., et al. (2011). Full-length transcriptome assembly from RNA-Seq data without a reference genome. Nat. Biotechnol. 29 (7), 644–652. doi:10.1038/nbt.1883

PubMed Abstract | CrossRef Full Text | Google Scholar

Gatto, M. A., Ippolito, A., Linsalata, V., Cascarano, N. A., Nigro, F., Vanadia, S., et al. (2011). Activity of extracts from wild edible herbs against postharvest fungal diseases of fruit and vegetables. Postharvest Biol. Technol. 61 (1), 72–82. doi:10.1016/j.postharvbio.2011.02.005

CrossRef Full Text | Google Scholar

Guan, M., Li, C., Shan, X., Chen, F., Wang, S., Dixon, R. A., et al. (2022). Dual mechanisms of coniferyl alcohol in phenylpropanoid pathway regulation. Front. Plant Sci. 13, 896540. doi:10.3389/fpls.2022.896540

PubMed Abstract | CrossRef Full Text | Google Scholar

He, X., Li, W., Zhang, W., Jin, X., Shenkute, A. G., Aynalem, T., et al. (2019). Transcriptome sequencing analysis provides insights into the response to Fusarium oxysporum in Lilium pumilum. Evol. Bioinforma. 15, 1176934319838818. doi:10.1177/1176934319838818

PubMed Abstract | CrossRef Full Text | Google Scholar

Jiao, X., Zhang, X., and Zhou, Y. (2021). Advances in the study of medicinal Lilium bulb diseases. China Plant Prot. 41 (10), 30–37.

Google Scholar

Jin, J., Zhang, H., Kong, L., Gao, G., and Luo, J. (2014). PlantTFDB 3.0: a portal for the functional and evolutionary study of plant transcription factors. Nucleic Acids Res. 42 (Database issue), D1182–D1187. doi:10.1093/nar/gkt1016

PubMed Abstract | CrossRef Full Text | Google Scholar

Jin, R., Kim, H. S., Yu, T., Zhang, A., Yang, Y., Liu, M., et al. (2021). Identification and function analysis of bHLH genes in response to cold stress in sweet potato. Plant Physiology Biochem. 169, 224–235. doi:10.1016/j.plaphy.2021.11.027

CrossRef Full Text | Google Scholar

Karre, S., Kumar, A., Yogendra, K., Kage, U., Kushalappa, A., and Charron, J. B. (2019). HvWRKY23 regulates flavonoid glycoside and hydroxycinnamic acid amide biosynthetic genes in barley to combat Fusarium head blight. Plant Mol. Biol. 100 (6), 591–605. doi:10.1007/s11103-019-00882-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Kumar, V., Hatan, E., Bar, E., Davidovich-Rikanati, R., Doron-Faigenboim, A., Spitzer-Rimon, B., et al. (2020). Phenylalanine increases chrysanthemum flower immunity against Botrytis cinerea attack. Plant J. 104, 226–240. doi:10.1111/tpj.14919

PubMed Abstract | CrossRef Full Text | Google Scholar

Kumudini, B. S., Jayamohan, N., Patil, S., and Govardhana, M. (2018). “Chapter 11-primary plant metabolism during plant-pathogen interactions and its role in defense,” in In plant metabolites and regulation under environmental stress. Editor P. Ahmad (Cambridge, MA, USA: Academic Press), 215–229. doi:10.1016/B978-0-12-812689-9.00011-X

CrossRef Full Text | Google Scholar

Kutschera, A., Dawid, C., Gisch, N., Schmid, C., Raasch, L., Gerster, T., et al. (2019). Bacterial medium-chain 3-hydroxy fatty acid metabolites trigger immunity in Arabidopsis plants. Sci. (New York, N.Y.) 364 (6436), 178–181. doi:10.1126/science.aau1279

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, P., Ruan, Z., Fei, Z., Yan, J., and Tang, G. (2021). Integrated transcriptome and metabolome analysis revealed that flavonoid biosynthesis may dominate the resistance of Zanthoxylum bungeanum against stem canker. J. Agric. Food Chem. 69 (22), 6360–6378. doi:10.1021/acs.jafc.1c00357

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, R., Cheng, H., Zhang, P., and Ai, F. (2018). Identification, biological characteristics and screening of bacteriostatic agents of Fusarium wilt. Genomics Appl. Biol. 37 (9), 3844–3851. doi:10.13417/j.gab.037.003844

CrossRef Full Text | Google Scholar

Liao, H., Ding, J., Dong, L., Jiang, Q., and Li, W. (2013). Preliminary report of lily anthracnose and epidemic chemical agent prevention and control test. Chin. melon Veg. 26 (1), 3. doi:10.3969/j.issn.1673-2871.2013.01.010

CrossRef Full Text | Google Scholar

Licausi, F., Ohme-Takagi, M., and Perata, P. (2013). APETALA2/Ethylene responsive factor (AP2/ERF) transcription factors: mediators of stress responses and developmental programs. New Phytol. 199 (3), 639–649. doi:10.1111/nph.12291

PubMed Abstract | CrossRef Full Text | Google Scholar

Love, M. I., Huber, W., and Anders, S. (2014). Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 15 (12), 550. doi:10.1186/s13059-014-0550-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Lu, C., Guo, N., Yang, C., Sun, H., and Cai, B. (2020). Transcriptome and metabolite profiling reveals the effects of Funneliformis mosseae on the roots of continuously cropped soybeans. Bmc Plant Biol. 20 (1), 479. doi:10.1186/s12870-020-02647-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Ma, R., Xiao, Y., Lv, Z., Tan, H., Chen, R., Li, Q., et al. (2017a). AP2/ERF transcription factor, Ii049, positively regulates lignan biosynthesis in isatis indigotica through activating salicylic acid signaling and lignan/lignin pathway genes. Front. Plant Sci. 8, 1361. doi:10.3389/fpls.2017.01361

PubMed Abstract | CrossRef Full Text | Google Scholar

Ma, Y., Liu, Y., Sun, A., Du, Y., Ye, M., Pu, X., et al. (2017b). Intestinal absorption and neuroprotective effects of kaempferol-3-O-rutinoside. RSC Adv. 7, 31408–31416. doi:10.1039/C7RA05415G

CrossRef Full Text | Google Scholar

Ninfali, P., Antonelli, A., Magnani, M., and Scarpa, E. S. (2020). Antiviral properties of flavonoids and delivery strategies. Nutrients 12 (9), 2534. doi:10.3390/nu12092534

PubMed Abstract | CrossRef Full Text | Google Scholar

Pan, Q., Zhu, M., Deng, J., and Tang, X. (2003). “Study on the occurrence and control of Fusarium lily wilt,” in East China plant pathology academic symposium and the 10th member congress of Jiangsu plant pathology society.

Google Scholar

Peng, J., Li, Z., and Zhong, Y. (2015). Study on the control effect of different chemical agents on lily blight and gray mold. Mod. Hortic. 10, 3. doi:10.3969/j.issn.1006-4958.2015.10.007

CrossRef Full Text | Google Scholar

Perez-Rodriguez, P., Riano-Pachon, D. M., Correa, L. G., Rensing, S. A., Kersten, B., and Mueller-Roeber, B. (2010). PlnTFDB: updated content and new features of the plant transcription factor database. Nucleic Acids Res. 38 (Database issue), D822–D827. doi:10.1093/nar/gkp805

PubMed Abstract | CrossRef Full Text | Google Scholar

Ren, S., Mei, H., Li, Y., Du, P., Wen, T., Qiao, Y., et al. (2020). Isolation and identification of lily Fusarium wilt pathogen and the effect of spore suspension concentration on the extent of disease. J. Nanjing Agric. Univ. 43 (6), 1097–1106. doi:10.7685/jnau.202001038

CrossRef Full Text | Google Scholar

Riechmann, J. L. (2018). “Transcription factors of arabidopsis and rice: a genomic perspective,” in Annual plant reviews online, 28–53. doi:10.1002/9781119312994.apr0297

CrossRef Full Text | Google Scholar

Rushton, P. J., Somssich, I. E., Ringler, P., and Shen, Q. J. (2010). WRKY transcription factors. Trends Plant Sci. 15 (5), 247–258. doi:10.1016/j.tplants.2010.02.006

PubMed Abstract | CrossRef Full Text | Google Scholar

Rushton, P. J., Torres, J. T., Parniske, M., Wernert, P., Hahlbrock, K., and Somssich, I. E. (1996). Interaction of elicitor-induced DNA-binding proteins with elicitor response elements in the promoters of parsley PR1 genes. Embo J. 15 (20), 5690–5700. doi:10.1002/j.1460-2075.1996.tb00953.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Shen, Y., Wang, J., Shaw, R. K., Sheng, X., Yu, H., Branca, F., et al. (2023). Comparative transcriptome and targeted metabolome profiling unravel the key role of phenylpropanoid and glucosinolate pathways in defense against Alternaria brassicicola in broccoli. J. Agric. Food Chem. 71 (16), 6499–6510. doi:10.1021/acs.jafc.2c08486

PubMed Abstract | CrossRef Full Text | Google Scholar

Sun, J., Wang, Y., Zhang, X., Cheng, Z., Song, Y., Li, H., et al. (2024). Transcriptomic and metabolomic analyses reveal the role of phenylalanine metabolism in the maize response to stalk rot caused by Fusarium proliferatum. Int. J. Mol. Sci. 25 (3), 1492. doi:10.3390/ijms25031492

PubMed Abstract | CrossRef Full Text | Google Scholar

Sun, L., Yu, D., Wu, Z., Wang, C., Yu, L., Wei, A., et al. (2019). Comparative transcriptome analysis and expression of genes reveal the biosynthesis and accumulation patterns of key flavonoids in different varieties of Zanthoxylum bungeanum leaves. J. Agric. Food Chem. 67 (48), 13258–13268. doi:10.1021/acs.jafc.9b05732

PubMed Abstract | CrossRef Full Text | Google Scholar

Varet, H., Brillet-Gueguen, L., Coppee, J. Y., and Dillies, 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

Wei, Y., Zhang, Y., Meng, J., Wang, Y., Zhong, C., and Ma, H. (2021). Transcriptome and metabolome profiling in naturally infested Casuarina equisetifolia clones by Ralstonia solanacearum. Genomics 113 (4), 1906–1918. doi:10.1016/j.ygeno.2021.03.022

PubMed Abstract | CrossRef Full Text | Google Scholar

Wu, Z., Li, Y., Wang, Y., Yang, L., Pan, H., Xie, Z., et al. (2015). Isolation and identification of an antagonistic bacterium against Fusarium on lily and its antifungal substances. Microbiol. China 42 (7), 1307–1320.

Google Scholar

Xiao, X. O., Lin, W., Feng, E., and Ou, X. (2023). Transcriptome and metabolome response of eggplant against Ralstonia solanacearum infection. Peerj 11, e14658. doi:10.7717/peerj.14658

PubMed Abstract | CrossRef Full Text | Google Scholar

Yadav, V., Wang, Z., Wei, C., Amo, A., Ahmed, B., Yang, X., et al. (2020). Phenylpropanoid pathway engineering: an emerging approach towards plant defense. Pathogens 9 (4), 312. doi:10.3390/pathogens9040312

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, L., Fu, Y., Fan, J., and Huang, Y. (2019). Morphological characteristics and karyotype analysis of three local varieties of Lilium brownii var.viridulum. Plant Sci. J. 37 (5), 559–568. doi:10.11913/PSJ.2095-0837.2019.50559

CrossRef Full Text | Google Scholar

Zeng, H., Li, R., Lu, Q., and Yan, D. (2020). First report of bulb rot on lily caused by Curvularia pseudobrachyspora in China. Plant Dis. 104 (1), 291. doi:10.1094/PDIS-05-19-0943-PDN

CrossRef Full Text | Google Scholar

Zhang, C., Luo, Q., Tang, W., Ma, J., Yang, D., Chen, J., et al. (2023). Transcriptome characterization and gene changes induced by Fusarium solani in sweet potato roots. Genes 14 (5), 969. doi:10.3390/genes14050969

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, Z., Chen, Q., Tan, Y., Shuang, S., Dai, R., Jiang, X., et al. (2021a). Combined transcriptome and metabolome analysis of alfalfa response to thrips infection. Genes 12 (12), 1967. doi:10.3390/genes12121967

PubMed Abstract | CrossRef Full Text | Google Scholar

Zheng, Y., Jiao, C., Sun, H., Rosli, H. G., Pombo, M. A., Zhang, P., et al. (2016). iTAK: a program for genome-wide prediction and classification of plant transcription factors, transcriptional regulators, and protein kinases. Mol. Plant 9 (12), 1667–1670. doi:10.1016/j.molp.2016.09.014

PubMed Abstract | CrossRef Full Text | Google Scholar