Abstract Background Platostoma palustre is a kind of plant resource with medicinal and food value, which has been differentiated into many different varieties after a long period of breeding. The cultivars of Taiwan(TW) and Pingyuan(PY) are widely grown in Guangdong, but a clear basis for species differentiation has not yet been established, resulting in the mixing of different species which limits their production and application. Results Regarding leaf surface morphology, the TW exhibited greater leaf area, non-glandular hairs, and the number of stomata than the PY. Regarding chemical activities, the TW exhibited higher total flavonoid content and antioxidant activity than the PY. In metabolomics, a total of 85 DAMs were detected, among which four flavonoid DAMs were identified, all of which were up-regulated in TW expression. Transcriptome analysis identified 2503 DEGs, which were classified according to their functional roles. The results demonstrated that the DEGs were primarily involved in amino acid metabolism, carbohydrate metabolism, sorting and degradation. A total 536 transcription factors (TFs) were identified, of which bHLH and MYB were the top two most abundant TFs families. Combined analysis of metabolome and transcriptome indicated that the phenylpropanoid pathway plays a significant role in flavonoid synthesis. Furthermore, real-time fluorescence qrt-PCR validation demonstrated that the expression trend of 10 DEGs was consistent with the transcriptomics data. Conclusion The phenylpropanoid pathway affects the synthesis of secondary metabolites, resulting in functional differences. In this study, metabolomic and transcriptomic analyses were performed to elucidate the regulatory mechanisms of flavonoid synthesis in P. palustre and to provide a theoretical basis for the identification, differentiation and breeding cultivation of different cultivars of P. palustre. Supplementary Information The online version contains supplementary material available at 10.1186/s12870-024-05909-5. Keywords: Platostoma palustre, Flavonoids, Metabolomics, Transcriptomics, Antioxidant activity Background Platostoma palustre (or Mesona chinensis Benth.), also known as “Xiancao”, is a herbaceous plant in the family Lamiaceae, which is distributed in the provinces of Guangdong, Guangxi, Zhejiang, Jiangxi and Taiwan, in China [[36]1], and is one of the important medicinal and edible cash crops [[37]2]. P. palustre has a lengthy history of serving as a resource for medicinal and food purposes, and it has been utilized in ancient folk medicine to treat hypertension, diabete and liver diseases [[38]3]. P. palustre is rich in many active ingredients, including polysaccharides, flavonoids, triterpenoids, phenolic acid and other compounds [[39]4]. It is widely used as an antioxidant drug [[40]5–[41]8] and has the effects of anti-oxidative stress damage [[42]9–[43]11], anti-inflammation [[44]12, [45]13] and lowering blood lipids [[46]8]. P. palustre is recognized for its refreshing flavor and high nutritional value as a food ingredient. It is widely used in the production of grass jelly, which has a significant consumer following and is gradually becoming a hot research topic. Due to the increasing demand both domestically and overseas, the cultivation area of P. palustre has expanded annually. After a prolonged period of natural selection and deliberate breeding, plant morphology diverges resulting in the formation of distinct cultivars. Among them, the two cultivars Taiwan(TW) and Pingyuan(PY) are more widely planted in Guangdong. Xiaoying Lu et al. utilized headspace gas-chromatography-mass spectrum(HS-GC–MS) to analyze the volatile constituents of TW and PY, which were distinctly categorized into two groups [[47]14]. In addition to variations in chemical composition, the two cultivars also exhibited differences in their agronomic traits and herb quality. TW plants were observed to be taller, with thicker leaves, denser trichome growth, higher yields, and fewer diseases. In contrast, PY plants were found to be shorter, possessed thinner leaves, shorter trichomes, lower yields, and more diseases, yet displayed an enhanced aroma after the herbs were stored. A clear basis for species differentiation has not been established, which has resulted in the two cultivars being often mixed in commercial production [[48]15]. This has led to variability in quality, which has seriously affected the stability of the quality of the herb and restricted its use in clinical applications and food production. Flavonoids represent the principal active ingredients in P. palustre, which can be employed as the primary components of anti-inflammatory, antibacterial, antiviral and stomachic and elimination drugs [[49]16, [50]17]. They are closely related to the quality of P. palustre herbs. However, the majority of current studies focus on the isolation, extraction and in vitro functional application of the active components of P. palustre, with fewer detailed studies on the mechanism of flavonoid synthesis and regulation in vivo. The correlation analysis of metabolomics and transcriptomics can analyze the co-expression of metabolites and genes related to metabolic pathways, to explore the relationship between functional characteristics of plant varieties and metabolites or genes [[51]18–[52]20]. This study will use non-targeted metabolomics and transcriptomics to analyze the TW and PY cultivars. The aim is to explore the differences in metabolites and genetic characteristics between these two cultivars, determine their respective total flavonoid content, and assess antioxidant activity to provide a strong foundation for the differentiation, breeding, and cultivation of P. palustre across various cultivars. Materials and methods Plant materials The Pingyuan grass (PY) and Taiwan grass (TW) seedlings used in this experiment were sourced from the planting base of P. palustre in Shizheng Town, Pingyuan County, Meizhou City, Guangdong Province, China. They were planted in Shi Zhen Mountain, Guangzhou University of Traditional Chinese Medicine (GZTCM) in early May 2022 and were identified as belonging to the family Labiatae, specifically Platostoma palustre, by Associate Professor Zhang Guifang of GZTCM. The leaves of the plant were harvested on 10th December 2022 and all tissues were promptly frozen in liquid nitrogen and transferred to -80 °C for subsequent use. Metabolome and transcriptome analyses included six and three biological replicates, respectively, and the sample numbers in the two omics were as follows: TW1 ~ TW6, PY1 ~ PY6; TWL1 ~ TWL3, PYL1 ~ PYL3. Screening micromorphological traits using scanning electron microscopy The sample leaves were collected in the same way as above, the leaf abaxial (AB) and adaxial (AD) surfaces of the species under this study were mounted onto stubs with double-sided adhesive tape, coated for 2 min with gold in a polaron JFC-1100E coating unit, and then examined and photographed with a JEOL JSM-IT200 in Guangzhou University of Chinese Medicine [[53]21]. The quantitative and qualitative features of the epidermal cells of both leaf surfaces (AB and AD) were recorded, including data on leaf surface structures, both closed and opened stomata, i.e., length and width, and measured by ImageJ analysis software [[54]22]. Extract preparation The resulting leaf samples were collected, dried in the shade to a constant weight, and then crushed through a 40-mesh sieve to obtain a fine sample powder. Then 0.6 g of sample fine powder was added to 100 mL of 60% ethanol solution and refluxated at 75 °C for 4 h, extraction of supernatant, sample liquid was obtained. Determination of total flavonoid content The content of total flavonoids was determined by sodium nitrite-aluminum nitrate method. Add 1.0 ml of sample liquid to 2.0 mL of 60% ethanol, followed by the addition of 0.5 mL of a 5% sodium nitrite solution. After allowing the reaction to proceed for 6 min, introduce 0.5 mL of a 10% aluminum nitrate solution and shake the mixture well. After a further 6 min, add 4.0 mL of a 4% sodium hydroxide solution and dilute the solution to a final volume of 10 mL with 60% ethanol and allow the mixture to stand for 15 min.The absorbance was determined at 510 nm using a spectrophotometer (A390, AOELab, CHINA) and the total flavonoid content was calculated using the rutin standard curve, and the results were expressed as the rutin equivalent mg/g of dry plant matter. Determination of antioxidant capacity According to the instructions, The DPPH reagent test kit (Solarbio, China, BC4755) was used to measure the DPPH radical scavenging ability of the sample solution of Mesona chinensis [[55]23]. 10 μL of sample liquid was mixed with 190 μL of working liquid and added to a 96 microplate. The mixture was incubated at room temperature under the dark for 30 min, then absorbance values at 517 nm were measured using a microplate reader (EPOCH2, BioTek, USA). The DPPH radical scavenging activities of each sample were calculated as the percent inhibition according to the following equation: DPPH radical scavenging activity (%) = [[A3-(A1—A2)] ÷ A3] × 100%. Where A1 is the absorbance of the test sample, A2 is the absorbance of the positive control, and A3 is the absorbance of the blank. The scavenging capacity for ABTS radical cations was measured using a total antioxidant capacity assay kit( Nanjing Jiancheng Bioengineering Institute, China, A015-2–1)following the manufacturer’s instructions [[56]24]. 10 μL of sample liquid was mixed with 190 μL of working liquid and added to a 96 microplate. After mixing, the solution was left at room temperature under dim light for 6 min, then the absorbance values at 405 nm were measured with using a microplate reader (EPOCH2, BioTek, USA). The antioxidant capacity of the sample is expressed as Trolox-Equivalent Antioxidant Capacity. FRAP scavenging capacities were determined using a total antioxidant capacity assay kit( Nanjing Jiancheng Bioengineering Institute, China, A015-3–1)following the manufacturer’s instructions [[57]25]. 5 μL of sample liquid, 180 μL of FRAP solution were added to the wells of a 96-well microplate. The mixture was incubated at 37℃ for 5 min, then the absorbance values at 593 nm were measured with using a microplate reader (EPOCH2, BioTek, USA). The antioxidant capacity of the sample is expressed as Trolox-Equivalent Antioxidant Capacity. Metabolite extraction and profiling analysis An appropriate amount of sample was accurately weighed into a 2 mL centrifuge tube, according to other research methods [[58]26]. 600 µL of MeOH (stored at -20℃), containing 2-Amino-3-(2-chloro-phenyl)-propionic acid (4 ppm), was added, and the mixture was vortexed for 30 s. 100 mg of glass beads were then added, and the sample was placed in a tissue grinder for 90 s at 60 Hz, followed by room temperature ultrasound for 15 min. The sample was centrifuged for 10 min at 12,000 rpm and 4℃ using a refrigerated centrifuge (H1850-R, Hunan Xiangyi Laboratory Instrument Development Co., Ltd., Hunan City, China). The supernatant was filtered through a 0.22 µm membrane and transferred into a detection bottle for LC–MS detection. LC–MS conditions The extracts were analyzed using an Ultra Performance Liquid Chromatography system(Thermo Vanquish, Thermo Fisher Scientific, USA)coupled with a tandem mass spectrometer(Thermo Orbitrap Exploris 120, Thermo Fisher Scientific, USA). The experiment conditions were as follows: The UPLC system was equipped with ACQUITY UPLC HSS T3 C18 column (Waters, 1.8 μm × 2.1 mm × 100 mm); solvent system, water. The column maintained at 40 ℃. The flow rate and injection volumewere set at 0.25 mL/min and 2 μL, respectively. For LC-ESI ( +)-MS analysis, the mobile phasesconsisted of (B2) 0.1% formic acid in acetonitrile (v/v) and (A2) 0.1% formic acid in water (v/v). Separation was conducted under the following gradient: 0 ~ 1 min, 2% B2; 1 ~ 9 min, 2% ~ 50% B2;9 ~ 12 min, 50% ~ 98% B2; 12 ~ 13.5 min, 98% B2; 13.5 ~ 14 min, 98% ~ 2% B2; 14 ~ 20 min, 2% B2. ForLC-ESI (-)-MS analysis, the analytes was carried out with (B3) acetonitrile and (A3) ammoniumformate (5 mM). Separation was conducted under the following gradient: 0 ~ 1 min, 2% B3; 1 ~ 9 min, 2% ~ 50% B3; 9 ~ 12 min, 50% ~ 98% B3; 12 ~ 13.5 min, 98% B3; 13.5 ~ 14 min, 98% ~ 2% B3;14 ~ 17 min, 2% B3 [[59]27]. Mass spectrometric detection of metabolites was performed on Orbitrap Exploris 120 (ThermoFisher Scientific, USA) with ESI ion source. Simultaneous MS1 and MS/MS (Full MS-ddMS2 mode, data-dependent MS/MS) acquisition was used. The parameters were as follows: sheath gaspressure, 30 arb; aux gas flow, 10 arb; spray voltage, 3.50 kV and -2.50 kV for ESI( +) and ESI(-), respectively; capillary temperature, 325 ℃; MS1 range, m/z 100–1000; MS1 resolving power, 60,000 FWHM; number of data dependant scans per cycle, 4; MS/MS resolving power, 15,000 FWHM; normalized collision energy, 30%; dynamic exclusion time, automatic [[60]28]. Qualitative and quantitative analysis of metabolites The raw data were firstly converted to mzXML format by MSConvert in ProteoWizard software package (v3.0.8789) [[61]29] and processed using XCMS [[62]30] for feature detection, retention time correction and alignment. The metabolites were identified by accuracy mass (< 30 ppm) and MS/MS data which were matched with HMDB [[63]31] ([64]http://www.hmdb.ca), massbank [[65]32] ([66]http://www.massbank.jp/), LipidMaps [[67]33] ([68]http://www.lipidmaps.org), mzcloud [[69]34] ([70]https://www.mzcloud.org) and KEGG [[71]35] ([72]http://www.genome.jp/kegg/). The robust LOESS signal correction (QC-RLSC) [[73]36] was applied for data normalization to correct for any systematic bias. Afternormalization, only ion peaks with relative standard deviations (RSDs) less than 30% in QC were kept to ensure proper metabolite identification. The Ropls [[74]37] software was used for all multivariate data analyses and modelings. After scaling data, models were built on principal component analysis (PCA) and partial least-square discriminant analysis (OPLS-DA). The metabolic profiles could be visualized as score plot, where each point represents a sample.The corresponding loading plot and S-plot were generated to provide information on the metabolites that influence clustering of the samples. All the models evaluated were tested forover fitting with methods of permutation tests. OPLS-DA allowed the determination of discriminating metabolites using the variable importance on projection (VIP). The P value and Variable importance projection (VIP) produced by OPLS-DA and applied to discover the contributable-variable for classification. Finally, P value < 0.05 and VIP values > 1were considered to be statistically significant metabolites. The identified metabolites in metabolomics were then mapped to the KEGG pathway for biological interpretation of higher-level systemic functions. The metabolites and corresponding pathways were visualized using KEGG Mapper tool. RNA extraction, library construction, and sequencing Total RNA was isolated using an RNAprep Pure Plant kit (Tiangen, Beijing, China) according to the instructions. Using the NanoPhotometer®spectrophotometer (IMPLEN, CA, USA) and the RNA Nano 6000 Assay Kit of the Agilent Bioanalyzer 2100 system (Agilent Technologies, CA, USA. RNA purity and integrity were assessed. cDNA libraries were created using the NEBNext® Ultra™ RNA Library Prep Kit for Illumina® (NEB, USA). Firstly, mRNA was purified from total RNA using poly-T oligo-attached magnetic beads. Fragmentation was carried out using divalent cations under elevated temperature in an Illumina proprietary fragmentation buffer. First strand cDNA was synthesized using random oligonucleotides and Super Script II. Second strand cDNA synthesis was subsequently performed using DNA Polymerase I and RNase H. Remaining overhangs were converted into blunt ends via exonuclease/polymerase activities and the enzymes were removed. After adenylation of the 3′ ends of the DNA fragments, Illumina PE adapter oligonucleotides were ligated to prepare for hybridization. To select cDNA fragments of the preferred 400–500 bp in length, the library fragments were purified using the AMPure XP system (Beckman Coulter,Beverly, CA, USA). DNA fragments with ligated adaptor molecules on both ends were selectively enriched using Illumina PCR Primer Cocktail in a 15 cycle PCR reaction. Productswere purified (AMPure XP system) and quantified using the Agilent high sensitivity DNA assayon a Bioanalyzer 2100 system (Agilent). The sequencing library was then sequenced on NovaSeq6000 platform (Illumina) by Suzhou PANOMIX Biomedical Tech Co., LTD. Transcriptome assembly and Gene functional annotation Samples are sequenced on the platform to get image files, which are transformed by the software of the sequencing platform, and the original data in FASTQ format (Raw Data) is generated. Then we use Cutadapt(v1.15) software filter the sequencing data to get high quality sequence (Clean Reads) for further analysis. For the transcriptome sequencing project without reference genome, we use Trinity(v2.5.1) software to montage Clean Reads for the transcripts for analyzing later. After the completion of the montaging, FASTA format Transcript sequence files can be obtained. The longest transcript of each gene was extracted as the representative sequence of the gene, called Unigene. Gene function for Unigene has been annotated. The database used in gene function annotation includes NR (NCBI non-redundant protein sequences), GO (Gene Ontology), KEGG (Kyoto Encyclopedia of Genes and Genome), eggNOG (evolutionary genealogy of genes: Non-supervised Orthologous Groups), Swiss-Prot, Pfam. Differential expression analysis The DEGs were identified with a threshold of |log2fold change|> 1 and p value < 0.05. ClusterProfiler (3.4.4) software was used to carry out the enrichment analysis of the KEGG pathway of DEGs, focusing on the significant enrichment pathway with p value < 0.05. Verification of candidate genes by quantitative real-time PCR (qRT-PCR) The DEGs identified in this study were verified by qRT-PCR (C1000 Touch PCR instrument, Bio-Rad, USA) using primers synthesized by Shanghai Sangon Bioengineering Co., Ltd (Tabel S1). The reaction mixture consisted of 10.0 µL of 2 × ChamQ SYBR qPCR Master Mix, 0.4 µL of forward primer (10 µM), 0.4 µL of reverse primer (10 µM), 0.4 µL of 50 × ROX Reference Dye 1, 2.0 µL of diluted cDNA, and 6.8 µL of ddH2O. qRT‒PCR was carried out using the following program: 95℃ for 30 s, followed by 40 amplification cycles at 95℃ for 5 s and 60℃ for 30 s. Finally, the relative expression levels of genes were analyzed using the 2^−ΔΔCT [[75]38]method with the reference gene GAPDH [[76]39]. Each sample was analyzed using three technical replicates. Statistical analysis All results were presented as the mean ± SD of 3 independent biological replications. SPSS 22.0 software was used to process the experimental data and analyze using analysis of variance (ANOVA) to identify significant differences between the data. Graphpad Prism 9 software and the free BIODEEP online platform were used to draw the experimental graph. TBtools were used to analyse the transcriptomic and metabolomic data representing the expression values and draw the heatmap. Results The leaf appearance and shape of TW and PY Fresh leaves from various cultivars were compared for measurement. The leaf blade traits of TW are ovate-lanceolate, while those of PY are long ovate (Fig. [77]1A, B). The leaf area and thickness were measured, and the resulting averages were 11.18 ± 3.16 cm^2 and 1.24 ± 0.11 mm for TW, 9.54 ± 1.34 cm^2 and 0.93 ± 0.16 mm for PY, respectively. Subsequently, the scanning electron microscope (SEM) of the dried leaves showed that (Fig. [78]1C-F) the non-glandular hairs of the leaves of the two cultivars were concentrated near the leaf veins, and the non-glandular hairs of TW were more and longer than PY, and the density of stomata TW was higher than that of PY. Fig. 1. Fig. 1 [79]Open in a new tab Leaf of P. palustre between TW and PY. The phenotypical characteristics of TW(A)and PY(B). Distribution characteristics of non-glandular hair in leaves of different varieties under 50 fold(C, D). Distribution characteristics of glandular hair and stomata in leaves under 200 fold(E, F) Flavonoid contents and antioxidant capacities between cultivars The total flavonoid content and antioxidant activity of TW and PY were determined to better understand the differences between the two cultivars (Fig. [80]2). The results showed that the total flavonoid content of TW was significantly higher than that of PY (120.97 and 99.76 mg/g of dry weight, respectively; p < 0.05). TW showed a statistically significant (p < 0.05) increase in the capacity to scavenge 2,2-diphenyl-1-picrylhydrazyl (DPPH) free radicals and ferric reducing antioxidant power (FRAP) (61.49% and 482.44 µmol Trolox/g in TW, respectively; 40.29% and 345.33 µmol Trolox/g in PY, respectively). In terms of scavenging free radical 2,2-azinobis [3-ethyl-benzothiazoline-6-sulfonic acid] (ABTS), the activities of TW and PY were 163.75 µmol Trolox/g and 131.31 µmol Trolox/g, respectively. However, there was no statistically significant difference between them. Fig. 2. Fig. 2 [81]Open in a new tab Total flavonoid contents and antioxidant activities between the TW and PY. A Total flavonoid contents. B DPPH radical scavenging activity. C ABTS scavenging capacity. D FRAP scavenging activity The total flavonoid content and antioxidant capacity of TW were higher than those of PY, and there was a positive correlation between them. Qualitative analysis of metabolites LC–MS/MS-based widely non-targeted metabolomics was carried out, and the mass spectrometry data were substituted into the metabolic component information database for detection and analysis. Further using ppm as a screening standard, compounds within ± 10 ppm were retained, and a total of 501 compounds were identified. The metabolites were clustered, in the heatmap (Fig. [82]3A), and all biological replicates were grouped, indicating a strong correlation among the replicates and high data reliability. The 501 metabolites can be classified in detail according to their properties (Table [83]S2): 82 lipids, 78 organic acids, 58 saccharides and alcohols, 53 amino acids, 26 phenols, 24 nucleotides and their derivatives, 22 esters, 21 flavonoids, 17 terpenoids, 16 steroids, 13 amines, 13 lignans and coumarins, 12 aldehydes, 11 alkaloids, 10 vitamins, 9 ketones and 36 others. Fig. 3. [84]Fig. 3 [85]Open in a new tab A Results of cluster analysis of metabolites detected in all samples. B, C The score and loading results of PCA Principal component analysis (PCA) To understand the differences in the composition of metabolites in P. palustre among cultivars, PCA was used to perform pattern recognition on all P. palustre samples, which can visualize grouping trends and outliers. The PCA scores of the mass spectrometry data for each treatment quality control sample (Fig. [86]3B) showed that two principal components (PC1 and PC2) accounted for 20.4% and 16.8%. From the figure (Fig. [87]3C), Metabolites contributing significantly to PC1 included saccharopine, mesaconate, geniposidic acid, 3,4-methylenedioxyamphetamine, etc. Metabolites contributing significantly to PC2 included (S)-2, 3-epoxysqualene, cellotetraose, 3-(2-hydroxyphenyl)propanoic acid, naringenin, etc. PCA results showed that the content of compounds identified in TW and PY samples was different to some extent. Discriminant analysis of orthogonal partial least squares(OPLS-DA) However, the detailed differences in each cluster remained unknown. OPLS-DA is a supervised pattern recognition method that can be used to analyze, classify and reduce the dimensionality of complex datasets. The OPLS-DA scores plot(Fig. [88]4A) shows the TW and the PY groups were clearly separated, indicating that the metabolites of TW and PY tended to separate. The OPLS-DA model was validated by the permutation test (Fig. [89]4B). In the figure, the R2 and Q2 values on the left were lower than those on the right, indicating good classification and predictability of the model. Fig. 4. [90]Fig. 4 [91]Open in a new tab A, B The score result and permutation test plot of OPLS-DA Differentially accumulated metabolites(DAMs) analysis According to the variable importance of the inter-group prediction (VIP) value and the significance test of the OPLS-DA model, the potential metabolites causing resource differences were screened. The variables, only under the condition of VIP > 1.0 and p < 0.05 (t-test), were considered to have a meaningful contribution to the OPLS-DA model. Using this model as a basis, we conducted a screening for distinct metabolites within the TW and PY. 85 DAMs were screened out, of which 39 were up-regulated and 46 were down-regulated. The DAMs can be classified in detail according to their properties (Fig. [92]5A): 13 organic acids, 12 amino acids, 10 lipids, 7 saccharides and alcohols, 6 esters, 4 alkaloids, 4 phenols, 4 flavonoids, 3 aldehydes, 3 terpenoids, 3 lignans and coumarins, 2 amines, 2 steroids, 1 nucleotides and their derivatives, 1 vitamin and 10 others. The DAMs between TW and PY are mainly organic acid compounds, accounting for 15.29%, most of which are involved in amino acid metabolism, followed by amino acids and lipids, accounting for 14.12% and 11.76%, respectively. It can be seen that the DAMs of flavonoids were all up-regulated in TW varieties, including Naringenin (14.78 times), Cyanidin 3-glucoside (4.91 times), 5, 7-dihydroxyflavone (3.95 times), Formononetin (1.39 times). Fig. 5. [93]Fig. 5 [94]Open in a new tab A Classification of significantly differentially regulated metabolites between TW and PY. B Scatterplot of the KEGG pathway enriched by different metabolites between TW and PY. The vertical axis represents the name of the pathway, and the horizontal axis represents the P value. The size and color of bubbles indicate the number and degree of enrichment of DAMs, respectively Significance analysis of the Kyoto Encyclopedia of Genes and Genomes (KEGG) can elucidate the primary biological processes associated with DAMs. The DAMs were involved in 56 pathways, and the enriched results also presented in a bubble chart (only the results of the top 20) are shown (Fig. [95]5B). The figure showed that the top five significantly enriched KEGG pathways were tyrosine metabolism, phenylalanine biosynthesis, Synthesis and degradation of ketone bodies, vitamin B6 metabolism, and phenylalanine metabolism. Naringenin and 5, 7-dihydroxyflavone were enriched in the KEGG pathway of Flavonoid biosynthesis. Transcriptomics analysis From Illumina sequencing, raw data and raw reads were obtained from the TW and PY cultivars. Q30 scores of all samples ranged from 93.91% to 94.46%(Table [96]1), indicating high sequencing data quality. After trimming of connector sequences and low-quality reads, 199,574 clean sequences were obtained and 71,220 Unigene was concatenated by Trinity (Table [97]2). The maximum and average lengths of the assembled genes were 16568 bp and 993.30 bp, respectively. The length of the annotated gene is shown in (Fig. [98]S1). Table 1. Raw data analysis Sample Reads No Bases (bp) Q30 (bp) N (%) Q20 (%) Q30 (%) PYL1 39378418 5946141118 5597142022 0.003467 98.02 94.13 PYL2 45207198 6826286898 6419085392 0.003431 97.98 94.03 PYL3 40596726 6130105626 5783013118 0.003458 98.1 94.33 TWL1 44939142 6785810442 6409878654 0.003443 98.13 94.46 TWL2 41389498 6249814198 5882005998 0.003477 98.03 94.11 TWL3 42570664 6428170264 6037231591 0.003429 97.93 93.91 [99]Open in a new tab Table 2. Sequence statistics Transcript Unigene Total Length (bp) 258023747 70742520 Sequence Number 199574 71220 Max. Length (bp) 16568 16568 Mean Length (bp) 1292.87 993.3 N50 (bp) 1879 1549 N50 Sequence No 44375 13270 N90 (bp) 574 415 N90 Sequence No 137542 51160 GC% 41.79 41.2 [100]Open in a new tab Differentially expressed genes(DEGs) analysis According to the conditions that |log2FoldChange|< 1 and p value < 0.05, 3384 DEGs were screened out, including 1454 up-regulated genes and 1930 down-regulated genes. Through matching with the GO, SwissProt, PFAM, KEGG and eggNOG database information, 2503 DEGs were identified (Fig. [101]6A). According to functional classification, 798 DEGs with known functions were further partitioned into 5 initial categories, including Metabolism, Genetic Information Processing, Cellular Processes, Environmental Information Processing, Organismal Systems. Then, we performed KEGG enrichment analysis of the second KEGG pathway category, which revealed the following pathways: Amino acid metabolism (92), Biosynthesis of other secondary metabolites (37), Carbohydrate metabolism (144), Energy metabolism (32), Environmental adaptation(47), Folding, sorting and degradation (61), Glycan biosynthesis and metabolism (14), Lipid metabolism (53), Membrane transport (21), Metabolism of cofactors and vitamins (32), Metabolism of other amino acids (27), Metabolism of terpenoids and polyketides (18), Nucleotide metabolism(13), Replication and repair (17), Signal transduction (64), Transcription(17), Translation (72), Transport and catabolism (37). KEGG enrichment annotation (Fig. [102]6B) showed that DEGs were mainly enriched in plant-pathogen interaction, ABC transporters, MAPK signaling pathway—plant, etc. Fig. 6. [103]Fig. 6 [104]Open in a new tab A Distribution and classification of DEGs between TW and PY. B Scatterplot of the KEGG pathway enriched by DEGs between TW and PY. The vertical axis represents the name of the pathway, and the horizontal axis represents the P value. The size and color of bubbles indicate the number and degree of enrichment of DEGs, respectively Transcription factor analysis of differentially expressed Genes Transcription factors (TFs) can bind specific cis-acting elements in the promoter region to regulate gene expression levels, thereby regulating the synthesis of metabolites in plants. A total of 536 DEGs were identified as TFs, and the most abundant TF family was bHLH, followed by MYB. 53 DEGs were annotated as bHLH TFs, of which 36 were up-regulated and 17 were down-regulated in TW (Fig. [105]S2). bHLH TFs can form a complex with other types of TFs such as MYB and WD40 complex proteins to co-regulate the expression of genes coding for key enzymes in the flavonoid biosynthetic expression of genes encoding key enzymes in the pathway. A total of 35 DEGs were identified to be annotated as MYB TFs, 19 were up-regulated and 16 were down-regulated in TW (Fig. S3). MYB TFs play important roles in the biosynthesis of plant secondary metabolites, especially flavonoids and anthocyanins. The results showed that two TF families, bHLH and MYB, showed significant differential expression, and it was hypothesised that these differentially expressed single genes in these two TF families might significantly affect the synthesis of flavonoids and other secondary metabolites in the two cultivated lines of P. palustre. The detailed gene name, predicted gene family and TFs description were listed in Table S3. There were also 42 DEGs annotated as MYB-related TFs, 33 DEGs as NAC TFs, 31 DEGs as ERF TFs, 26 DEGs as WRKY TFs, 25 DEGs as C2H2 TFs, and 25 DEGs as bZIP TFs. Correlation analysis between the Metabolomic and transcriptomic data We used the O2PLS (Two-way Orthogonal Partial Least Squares) method to assess the intrinsic correlation between transcriptomics and metabolomics by calculating the score for each sample and obtaining the joint score, and by calculating the loading values for each differential gene and differential metabolite to obtain the loading plot and the absolute values of the load values of the first 20 DAMs/DEGs are selected to construct the histogram (Fig. [106]7A,B). From the results it was seen that DAMs such as Chorismate, Naringenin, Cyanidin 3-glucoside are associated with the phenylpropanoid pathway and are involved in the synthesis of flavonoids by the phenylpropanoid pathway. Fig. 7. [107]Fig. 7 [108]Open in a new tab Metabolite and transcript correlation analysis. A The top 20 DEGs and DAMs loading histogram. Orange bars represent genes and red bars represent metabolites. The abscissa represents the combined loading value pq1, and the ordinate represents the DEGs/DAMs. B Assess the intrinsic correlation between the metabolites and transcripts by the two-way orthogonal partial least squares (O2PLS) method. The orange dots represent genes and the red dots represent metabolites. The abscissa and ordinate represent the combined loading value, each gene/metabolite has a relative coordinate point in pq1 and pq2, p represents the loading value of the gene, and q represents the loading value of the metabolite. C Histogram of the KEGG pathway enriched by the diferentially expressed genes and metabolites. The ordinate represents the pathway name, and the abscissa represents the pathway to impact Correlation analysis of metabolomics and transcriptomics of TW and PY cultivars of P. palustre was performed in combination with KEGG pathway enrichment analysis, and the same pathway was enriched in the transcriptomics and metabolomics(Fig. [109]7C). 5 KEGG pathways are involved in biological processes related to phenylpropanoid metabolism, including Phenylpropanoid biosynthesis (map00940), Phenylalanine metabolism (map00360), and Flavonoid biosynthesis (map00941), Tyrosine metabolism (map00350), Tryptophan metabolism(mapmap00380). In the Phenylpropanoid biosynthesis pathway, it was discovered that 10 DEGs have also been involved in the biological processes of Cyanoamino acid metabolism and Starch and sucrose metabolism. In transcriptomics, DEGs that were annotated in KEGG pathways related to flavonoid biosynthesis and phenylpropanoid biosynthesis were presented in heat maps to illustrate their expression levels (Fig. [110]8A). 7 DEGs were found to be enriched in KEGG pathways related to flavonoid biosynthesis. The annotated enzymes were ANS (TRINITY_DN30600_c0_g1), CHS (TRINITY_DN90504_c0_g1), and FLS (TRINITY_DN14690_c0_g1), UGAT (TRINITY_DN14126_c0_g1) and other enzymes were down-regulated in TW. In the flavonoid precursor synthesis pathway, 22 DEGs were found to be enriched in KEGG pathways related to phenylpropanoid biosynthesis. Between TW vs PY, we found that there were multiple DEGs that co-encode an enzyme, such as beta-glucosidase (E3.2.1.21) with 5 DEGs, and TRINITY_DN23872_c0_g1, TRINITY_DN15452_c0_g2, TRINITY_DN12046_c0_g1, and TRINITY_DN10064_c0_g1 expression were up-regulated, TRINITY_ DN75332_c0_g1 expression downgraded. Among them, TRINITY_DN15452_c0_g2, TRINITY_DN12046_c0_g1 and TRINITY_DN75332_c0_g1 were identified as C2H2 transcription factors. HCT had 2 DEGs, TRINITY_DN48336_c0_g1 expression was up-regulated and TRINITY_DN11381_c0_g1 expression was down-regulated. Fig. 8. [111]Fig. 8 [112]Open in a new tab A Heatmap of gene expression associated with flavonoid biosynthesis and phenylpropanoid biosynthesis. B The flavonoids biosynthesis pathway potentially used by the two cultivars based on gene expression levels and pathway annotations in the KEGG database. In the heatmap, the red block represents up-regulation and green, and blue represent downregulation of expression, respectively. The metabolic components and genes were mapped to be generated based on phenylpropanoid biosynthesis (k000940) and flavonoid biosynthesis (ko00941) To better understand the molecular mechanism of flavonoid accumulation in the two cultivars of P. palustre, we performed a predictive reconstruction of flavonoid biosynthesis based on the results of the screened flavonoid biosynthesis-related differential genes and KEGG pathway analyses (Fig. [113]8B). Although flavonoid DAMs, including naringenin, 5, 7-dihydroxyflavone, cyanidin 3-glycoside, and Formononetin, showed significant upregulation in TW, functional genes related to flavonoid synthesis such as PAL, CHS, and ANS showed higher expression levels in PY. Notably, no significant changes were observed in the expression of direct catalytic genes in these DAMs. Typically, CHS is a key enzyme in the process of flavonoid synthesis, while ANS is an important enzyme in the process of anthocyanin synthesis. However, the accumulation of 5, 7-dihydroxyflavanone and cyanidin-3-glucoside observed in TW was in contrast to the trend of expression of the enzymatic activities of CHS and ANS. This phenomenon may be due to the fact that the higher flavonoid content in TW triggers the feedback regulation mechanism, which keeps the genes of these related enzymes at a low expression level to maintain the internal equilibrium state of the system [[114]40–[115]42]. Confirmation of DEGs via qRT-PCR To validate the reliability and stability of the RNA sequencing (RNA-seq) data of the DEGS, we selected 9 DEGS for qRT-PCR validation (Fig. [116]9). Among them 3 DEGs (ANS, FLS, CHS) were associated with flavonoid biosynthesis, 2 DEGs (HPD, AOC) were associated with phenylalanine metabolism and tyrosine metabolism, and 4 DEGs (CAD, UGT72E, PAL, blgx) were associated with phenylpropanoid biosynthesis. The results are consistent with the transcriptomic gene expression trends, suggesting that the transcriptomic data were reliable. Fig. 9. Fig. 9 [117]Open in a new tab qRT-PCR validation of the relative expression levels of 9 selected genes from different cultivars of P. palustre Discussion In this experiment, 85 DAMs (39 up-regulated and 46 down-regulated) and 3384 DEGs (1454 up-regulated and 1930 down-regulated genes) were identified between TW and PY cultivars. In metabolomics, 21 flavonoids and 16 terpenoids, as well as 11 alkaloids, were identified among the secondary metabolites (Table S4). Plants can produce a variety of secondary metabolites with different chemical characteristics in response to the changing environment [[118]43]. Flavonoids are secondary metabolites of a phenolic nature that exist in nature and possess a broad spectrum of pharmacological actions [[119]44], such as antibacterial activity [[120]45], anti-inflammatory activity [[121]46–[122]48], antiviral activity [[123]49], antioxidant activity [[124]50–[125]52], and play a defense role in plants in response to various biological and abiotic stresses through dynamic changes [[126]53]. We have detected 21 flavonoids in total, comprising six flavonoids, five anthocyanins, four flavonols, three dihydroflavonols, two isoflavones, and one dihydroisoflavone. All four flavonoid DAMs were up-regulated in TW. Naringenin, the DAMs measured in the metabolome, is a significant intermediate in the flavonoid synthesis pathway [[127]54], as it is a flavanone produced through enzymatic catalysis of chalcone. Formononetin is the conversion of chalcone to isoflavones after a series of enzyme-catalyzed conversions, and liquintigenin generates Formononetin under the regulation of HIDH [[128]55]. Flavanones are catalysed by F3H to form dihydroflavonols, which are then regulated by DFR, ANS, etc. to form Cyanidin [[129]56, [130]57], and Cyanidin combines with glycosides to form Cyanidin 3-glucoside; 5,7-Dihydroxyflavone is formed by cinnamoyl-coenzyme A under the constant regulation of CHS, CHI, and FNS II to form salicin, which is formed by dehydrogenation of salicin [[131]58]. Subsequently, we demonstrated via sodium nitrite-aluminum nitrate spectrophotometry and antioxidant activity experiments that the TW possesses a greater total flavonoid content and antioxidant capacity than PY. Furthermore, we confirmed a positive correlation between antioxidant activity and total flavonoid content in P. palustre. We hypothesized that the differences in the ability of TW and PY cultivars to synthesize flavonoids were related to the phenylpropanoid pathway. It is widely accepted that flavonoid synthesis primarily occurs via the phenylpropanoid pathway [[132]59]. Coumaroyl coenzyme A, with phenylalanine and tyrosine as precursors, undergoes catalysis by CHS and CHI to form dihydroflavonoids. These are then processed by various enzymes to produce different types of flavonoid compounds [[133]60–[134]62]. As evidenced by the transcriptomics and metabolomics findings of this experiment, there was a significant enrichment of DAMs and DEGs in the phenylpropanoid pathway. In transcriptomics, 5, 6, and 16 DEGs were enriched in KEGG pathways for tyrosine, phenylalanine, and phenylpropanoid biosynthesis, respectively. In metabolomics, 5, 3, and 4 DAMs were enriched in KEGG pathways for tyrosine, phenylalanine, and phenylpropanoid biosynthesis, respectively. Chorismate, Naringenin, and Cyanidin 3-glucoside were also found to be associated with the phenylpropane pathway and involved in flavonoid synthesis in O2PLS analysis. Enrichment analyses aim to identify biological pathways that are significant in various biological processes, phenylpropanoid metabolism is one of the important secondary metabolic pathways in plants. Therefore, it is reasonable to hypothesize that the functional differences in flavonoid synthesis in P. palustre may be linked to the phenylpropanoid pathway. The two cultivars of P. palustre have different resistance traits that can be interpreted and explored from a variety of perspectives. As can be seen from the SEM images, TW has more non-glandular hairs than PY. Non-glandular hairs, as specialised accessory structures embedded in the plant epidermis, separating the external environment from the plant epidermis. They perform important functions, including defense [[135]63] and resistance [[136]64]. These functions are particularly important in the context of insect and pathogen resistance. Non-glandular hairs can affect insect development, feeding, movement, and delay or limit the access of phytophagous insects to the plant epidermis [[137]65]. Additionally, research studies have indicated that non-glandular hairs may possess disease resistance. For instance, the legume red clover has been observed to exhibit resistance to powdery mildew, which facilitates normal plant development [[138]66]. Non-glandular hairs safeguard the apical shoots of plant stems from external harm during the preliminary stages of growth. It has been observed that plants with longer non-glandular hairs exhibit greater resilience towards cold temperatures [[139]67]. Concurrently, flavonoids have significant antioxidant effects and reduce oxidative damage caused by the accumulation of reactive oxygen species, which play a role in plant growth, development, and defense [[140]68]. To a certain extent, the antioxidant activity of plants can reflect the strength of stress tolerance and disease resistance. It has been demonstrated that maize seedlings treated with exogenous abscisic acid (ABA) showed an increase in antioxidant enzyme activities CAT and SOD, resulting in improved resistance to drought [[141]69]. Additionally, papaya treated with methyl jasmonate (MEJA) displayed increased antioxidant activity, leading to enhanced tolerance to low temperatures [[142]70]. Transfecting tobacco with maize Cat2 resulted in higher catalase (CAT) activity in comparison to untransfected plants, thus inhibiting pathogen growth [[143]71]. In this experiment, the flavonoid content and antioxidant activity of the TW cultivar exceeded that of the PY cultivar in this investigation. These may explain why the TW cultivar experiences less disease. In addition to flavonoids, we have identified 11 types of alkaloids, including isoquinoline, pyridine, and indole alkaloids, in the study. 4 of the DAMs, Isocorypalmine, Isoquinoline, 17-O-Acetylnorajmaline, and Nicotine, were found to be up-regulated in TW expression. The composition of the 16 terpenoids includes 7 monoterpenoids, 2 sesquiterpenes, 3 diterpenes, and 4 triterpenoids. DAMs xanthotoxin was significantly up-regulated in TW, while that of Soyasapogenol A was down-regulated. The phenylpropanoid metabolic pathway is a significant pathway for the synthesis of plant secondary metabolites. Numerous studies have uncovered the enzymatic reaction process of phenylpropanoid metabolism in plants, along with the entire metabolic pathway’s regulatory mechanism. PAL, C4H, and 4CL, which are the key enzymes in the phenylpropanoid metabolic pathway [[144]72], successively reacted to produce cinnamic acid, p-hydroxycinnamic acid, and p-coumaroyl coenzyme A. These substrates are eventually converted into a variety of phenylpropanoid compounds, including flavonoids, lignans, terpenoids, alkaloids, and other secondary metabolites. Research has shown that the phenylpropanoid metabolic pathway is linked to plant resistance [[145]68]. When plants are under stress, enzymes in this pathway, including PAL, 4CL, and C4H, become more active. This increase in activity leads to the metabolic synthesis of lignin, which promotes the degree of cellular lignification. Additionally, the pathway produces a variety of metabolites such as phenols, flavonoids, and terpenes, which further synthesize phytophysical proteins. It has been shown that the phytopropane pathway plays a role in producing salicylic acid, which is essential for activating the plant’s immune pathway and enhancing its defense against diseases. As a result, this regulates the plant’s disease resistance and defense ability [[146]73]. At the same time, the study of transcriptome and metabolome data can provide targeted improvement strategies, such as genetic engineering to improve the synthesis of specific secondary metabolites, to optimize the efficacy and stability of swamp grass [[147]74–[148]76]. This study revealed the relationship between flavonoid gene regulatory network and metabolic network, provided insights into the coordinated metabolic pathways of key enzymes and metabolites, provided scientific basis for subsequent variety selection and breeding work [[149]77], and played a theoretical guidance for ensuring the stability of P. palustre in the process of picking, production, storage and use [[150]78, [151]79]. Conclusion Differences in phenotypic characteristics, total flavonoid content, and antioxidant activity of P. palustre samples from TW and PY cultivars were identified. Analyzing the transcriptome and metabolome, we observed the impact of the phenylpropanoid pathway on the synthesis of secondary metabolites in P. palustre, which may account for functional variations, and explored the correlation between antioxidant activity and plant resistance. Our investigation has provided a new understanding of the regulatory mechanisms behind flavonoid biosynthesis in P. palustre. This has led to critical theoretical data that significantly aids the process of identifying P. palustre from various cultivars, as well as facilitating its selection and cultivation. Supplementary Information [152]Supplementary Material 1.^ (51.6KB, xlsx) [153]Supplementary Material 2.^ (1.3MB, docx) Acknowledgements