2.1 Reagents and plant samples

LC-MS grade reagents such as acetonitrile, formic acid and water were purchased from Merck (Darmstadt, Germany), Sigma-Aldrich (St. Louis, MO, USA), and Watsons (Hong Kong, China), respectively. The other analytical reagents were purchased from Chengdu Kelong Chemical Reagent Factory.

The stems of 1-year and 3-year grown D. nobile (10 samples, respectively) were randomly collected from the Dendrobium Germplasm Garden of Chishui Xintian Company located in Chishui City, China in October of 2019 and authenticated by Associate Professor Daopeng Tan (Pharmacognosy, Zunyi Medical University). All fresh stems for secondary metabolite analysis were dried, ground into powder, passed through a sieve with 300 mesh, and stored at −80 ℃. The stems for transcriptome analysis were processed as described previously (Li et al., 2022). Briefly, all samples were divided into two groups according to growth years and each group had three biological replicates for sequencing.

2.2 UPLC-Q/TOF-MS conditions for secondary metabolite analysis

Agilent 1290 Infinity II UPLC liquid chromatograph system coupled with Agilent 6545 Q/TOF-MS system (Agilent, MA, USA) was employed for the present secondary metabolites analysis. The separation was performed by the Waters CORTECS UPLC C18 (100 mm × 2.1 mm, 1.6 μm) chromatographic column with 0.1 % formic acid water (A) − 0.1 % formic acid acetonitrile (B) as the mobile phase, and the gradient elution conditions was set as: 0–0.5 min, 5 % B; 0.5–4 min, 5 % − 40 % B; 4–5 min, 40 % − 75 % B; 5–5.1 min, 75 % − 95 % B; 5.1–6.5 min, 95 % B; 6.5–6.6 min, 95 % − 5 % B; 6.6–10 min, 5 % B. Flow rate: 0.4 mL/min, column temperature: 40 ℃, injection volume: 5 μL.

Mass spectrometry conditions: Agilent 6545 Q-TOF high-resolution mass spectrometry was used, with separate acquisition in positive and negative ion modes. Ion source parameters: capillary voltage (CV) of 4000 V, nozzle voltage (Liu et al.) of 1000 V, atomization gas pressure of 45 psi, dry gas flow rate of 10 L/min, sheath gas temperature of 350 ℃, sheath gas flow rate of 11 L/min; ion source temperature (TEMP) of 350 ℃, collision energy 10 V, mass scan range of m/z 50 ∼ 1200. The secondary mass spectrometry information was acquired in Auto MS/MS mode with CEs of 20, 30, and 40 V, respectively.

The pulverized powder samples (75 mg) were accurately weighed, added with 70 % methanol 1 mL, and then extracted by ultrasonication (400 W, 50 kHz) for 30 min, respectively. After cooling, the solution was centrifuged at 12000 rpm for 5 min, and the supernatant was loaded into the injection vial for UPLC-Q-TOF/MS analysis. An additional 0.1 mL aliquot of the multi-batch sample solution was added into a 5 mL centrifuge tube and mixed. The mixture was centrifuged at 12,000 rpm for 5 min. The supernatant was used as a quality control (QC) sample.

2.3 Identification of secondary metabolites of D. nobile

The mass spectrometry data were imported into Agilent MassHunter Profinder 10.0 software for peak matching, peak alignment, ion fusion and deconvolution processing. Based on the peak area, retention time and molecular weight, the fragmented peaks with false positives were excluded. The corresponding molecular formulae were obtained by fitting and calculating with Qualitative Analysis B.07.00 software, and matched with the local database for preliminary structure inference. The secondary metabolites were further identified based on the secondary mass spectrometry data with the information provided by references, Scifinder database, or in-house databases.

2.4 RNA extraction, cDNA library and Illumina sequencing

Total RNA of each sample was extracted by using Trizol (Invitrogen, CA, USA) on the basis of manufacturer’s instructions for the creation of individual cDNA libraries and Illumina sequencing. After purification and characterization, RNA with a RIN number greater than 7.0 was next purified from total RNA (5 μg) using poly-T oligo-attached magnetic beads. And then, the mRNA was split into small fragments of approximately 200 bp using divalent cations at high temperature. The cleaved RNA fragments were then reverse transcribed to create the final cDNA library following the protocol of the mRNASeq sample preparation kit (Illumina, San Diego, USA), with an average insert size of 300 bp (±50 bp) for the paired-end library. Paired-end sequencing was then performed on an Illumina Nova seq™ 6000.

2.5 De novo assembly, unigene annotation and functional classification

Methods for de novo assembly, unigene annotation and functional classification were processed as previous studies (Li et al., 2022). Briefly, de novo assembly of transcriptomes was performed by Trinity 2.4.0 (Grabherr et al., 2011). Trinity groups transcripts into clusters based on shared sequence content. Such transcript clusters are very loosely referred to as “genes”. The longest transcript in the cluster was selected as the “gene” sequence (aka Unigene). All assembled Unigenes were aligned using DIAMOND with the Non-Redundant (Nr) Protein Database (https://www.ncbi.nlm.nih.gov/), Kyoto Encyclopedia of Genes and Genomes (KEGG) (https://www.genome.jp/kegg/), Gene Ontology (GO) (https://www.geneontology.org), the SwissProt (https://www.expasy.ch/sprot/) and eggNOG (https://eggnogdb.embl.de/) databases with a threshold E value < 0.00001.

2.6 qRT-PCR analysis

To verify the accuracy of the transcriptome data, we selected nine differential genes (DEGs) for qRT-PCR validation. Total RNA was extracted by BIO FIT polysaccharide polyphenol bioreagent kit and cDNA was synthesized by TIANGEN reverse transcription kit. For each sample, three biological replicates were used for qRT-PCR assay with three technical replicates. Gene expression was calculated by the 2^-ΔΔCt method. The primers used for qRT-PCR were listed in Table S2.

2.7 Molecular docking

The crystal structures of the candidate glycosyltransferase-related proteins from Dendrobium were downloaded from the AlphaFold Protein Structure Database (AlphaFold DB, https://alphafold.ebi.ac.uk) and modified using the Autodock tools 1.5.6 software (Jumper et al., 2021; Varadi et al., 2022). AlphaFold DB provides programmatic access to and interactive visualization of predicted atomic coordinates, per-residue and pairwise model-confidence estimates and predicted aligned errors. The 3D structures of sesquiterpene glycosides and their aglycones were built by ChemBioDraw Ultra14.0 and then converted to PDBQT coordinateds using AutoDockTools. The rotatable bonds in the ligand were assigned with AutoDock Tools, and the ligand docking was performed with the AutoDock Vina. Construct 2D maps of protein–ligand interactions using LigPlot software to analyze the interaction forces between protein and ligand binding. Construct 3D maps of protein–ligand interactions using PyMOL software to view the binding sites between proteins and ligands.

2.8 Data analysis

All visualizations analysis such as PCA, OPLS-DA and Volcano plot were performed in the R program (version 4.1.1) using the packages of ggplot2, pheatmap, and others or SIMCA-P 14.0 software (Umetrics AB, Umea, Sweden). Differences between two groups were calculated by t-test in R, p-values < 0.05 were identified as statistically significant differences, and all data were expressed as “mean ± standard error (SEM)”.

3.1 Profile analysis of secondary metabolites of D. nobile

In the present study, UPLC-Q/TOF-MS coupled with multivariate analysis were employed for comparative analysis of secondary metabolites in D. nobile with 1 or 3 grown years. A total of 511 ion signals were detected by both positive and negative modes, excluding some that could be false signals. First, principal component analysis (PCA) was used to identify differences in the profiles of the secondary metabolites in 1-year-grown and 3-year-grown D. nobile. It was obvious that the 1-year-grown and 3-year- grown samples were located in two different regions (Fig. 1A), which indicated that the composition of secondary metabolites in 1-year-grown and 3-year-grown D. nobile was significantly different. Next, Volcano plot were employed to screen for significantly different secondary metabolites in 1-year-grown and 3-year-grown D. nobile (Fig. 1B). The results showed that 345 of the detected secondary metabolites were higher in the 3-year-grown D. nobile, only 97 secondary metabolites were higher in the 1-year-grown D. nobile, and another 69 were not significantly different (Fig. 1C).

Close

Fig. 1

Effect of growth years on the profile of secondary metabolites in D. nobile. (A) PCA plot of 1-year-grown and 3-year-grown D. nobile. Each spot represented one sample. (B) Volcano plot of 1-year-grown and 3-year-grown D. nobile. Each spot represented one secondary metabolite. (C) Amounts of secondary metabolites with different accumulation changes.

Export to PPT

Generally, secondary metabolites in plants accumulate positively with increasing growth years, as in the case of these 345 secondary metabolites in D. nobile. Instead, these 97 secondary metabolites that showed the opposite accumulation pattern aroused our interest. Based on the structural identification by mass spectrometry, 29 amino acids, 7 phenolic acids, 6 alkaloids and 2 sesquiterpenes were identified among these 97 secondary metabolites that decreased with the growth year (Fig. 1C).

3.2 Identification of alkaloids and sesquiterpenes in D. nobile

The alkaloids represented by Dendrobine are the main active components in D. nobile (Huang et al., 2019; Lv et al., 2020). It has been reported that sesquiterpenes have a common biosynthetic pathway with alkaloids in D. nobile (Gong et al., 2021; Li et al., 2017; Li et al., 2022). In order to determine the structure features of differential alkaloids and sesquiterpenes with growth-dependent expression styles in D. nobile, annotation of secondary metabolites was performed by comparing accurate masses and fragments with literature data or in-house databases. Finally, the 6 alkaloids were identified as dendrobine, N-methyl-dendrobine, 6-hydroxyl-dendrobine, nobilonine, N-isopentenyl-dendrobine, and dendrobine-N-oxide, respectively. And the 2 sesquiterpenes were elucidated as confertifolin and elemol (Fig. 2). Their detailed information including chemical name, adducts, experimental and theoretical m/z, molecular formula, and MS/MS data were shown in Table S1. Furtherly, t-test was performed on the peak area data of the 6 alkaloids and 2 sesquiterpenes to verify that these structurally confirmed secondary metabolites had a significant trend and statistical significance. The results indicated that those secondary metabolites were indeed significantly higher in the 1-year-grown D. nobile (P < 0.05) as shown in Fig. 3.

Close

Fig. 2

Chemical structure of alkaloids and sesquiterpenes in D. nobile that content decreased with growth year.

Export to PPT

Close

Fig. 3

Comparing intensity of 6 alkaloids and 2 sesquiterpenes in D. nobile that content decreased with growth years. The asterisk indicated significant differences (*P < 0.05).

Export to PPT

3.3 Identification of sesquiterpene glycosides in D. nobile

Although it has been suggested in the literature that the alkaloids in D. nobile were furtherly biosynthesized from sesquiterpenes as precursor compounds (Gong et al., 2021; Li et al., 2017; Li et al., 2022), however, sesquiterpenes found in previous chemical composition investigations were rare, instead a large number of sesquiterpene glycosides have been isolated from D. nobile (Tan et al., 2023a; Tan et al., 2023b). In the present study, only 2 sesquiterpenes were found its content decreased with growth years (Fig. 3). However, no sesquiterpene glycosides were identified from the secondary metabolites that decreased with growth years. Not surprisingly, 14 sesquiterpene glycoside components were furtherly found in secondary metabolites that increased with growth years through the same identification method as the alkaloids by comparing accurate masses and fragments with literature data or in-house databases, including dendronobiloside A, dendronobiloside B, dendronobiloside C, dendronobiloside D, dendronobiloside E, dendromoniliside A, dendromoniliside B, dendromoniliside C, dendromoniliside D, dendroside A, dendroside C, dendroside E, dendroside F, dendroside G (Fig. 4). Their detailed information was shown in Table S1. In the same way, t-test was performed on the peak area data of the 14 sesquiterpene glycosides to verify that these structurally confirmed secondary metabolites had a significant trend and statistical significance. The results indicated that those secondary metabolites were indeed significantly higher in the 3-year-grown D. nobile (P < 0.05) as shown in Fig. 5.

Close

Fig. 4

Chemical structure of sesquiterpene glycosides in D. nobile that content increased with growth years.

Export to PPT

Close

Fig. 5

Comparing intensity of 14 sesquiterpene glycosides in D. nobile that content increased with growth years. The asterisk indicated significant differences (*P < 0.05).

Export to PPT

3.4 Hypothesis on the biosynthetic pathways of alkaloids and sesquiterpene glycosides in D. nobile

Since Dendrobine is the main active component in D. nobile and it has a special sesquiterpene backbone, the biosynthesis of Dendrobine has been a research hotspot focused by scientists. As early as 1966, the biosynthesis of Dendrobine was confirmed in the literature based on isotope tracing by the sesquiterpenoid biosynthetic pathway (Gong et al., 2021). Furtherly, the complete biosynthetic pathway of Dendrobine was proposed by infection with Mycena sp. (MF23) (Li et al., 2017a). However, sesquiterpenes were isolated from D. nobile are relatively rare, on the contrary, the sesquiterpene glycosides, which are present in large quantities, were not included in this biosynthetic pathway. Therefore, we proposed a new hypothesis on the biosynthetic pathway of alkaloids and sesquiterpene glycosides in D. nobile based on the previous biosynthetic pathway (Fig. 6). In the present biosynthetic pathway, a portion of the sesquiterpenes were glycosylated to form sesquiterpene glycosides that were stored in D. nobile, which could reasonably explain why the accumulation of alkaloids decreases with the growth year, while the sesquiterpene glycosides show the opposite accumulation pattern.

Close

Fig. 6

Hypothesized biosynthetic pathways of alkaloids and sesquiterpene glycosides in D. nobile.

Export to PPT

3.5 UDP-glycosyltransferase gene analysis

Although it has been reported in the literature that the variation of Dendrobine content in D. nobile was related to the expression of genes such as Cytochrome P450, Aminotransferase, N-methyltransferase, etc. (Gong et al., 2021; Li et al., 2017; Li et al., 2022). However, it could not explain why sesquiterpene glycosides showed opposite accumulation patterns compared to alkaloids in D. nobile. The conversion of sesquiterpene glycoside from sesquiterpenes in organisms must be biocatalyzed by UDP-glycosyltransferases. Therefore, we hypothesized that the increase of sesquiterpene glycosides in D. nobile is due to the high expression of UDP-glycosyltransferases.

In the present study, 184 UDP-glycosyltransferase genes in transcriptome data of D. nobile were collected by KEGG (Kyoto Encyclopedia of Genes and Genomes; https://www.kegg.jp/ or https://www.genome.jp/kegg/). First, OPLS-DA was used to show the difference in the expression of UDP-glycosyltransferase genes in different growth years. The results showed that 1-year and 3-year-grown samples were located in two different regions, indicating that expression profiles involved in UDP-glycosyltransferase genes were significantly affected by growth years (Fig. 7A). The VIP value distribution map in the OPLS-DA model was used to screen the points with VIP values greater than 2 as the differential genes between the two groups (Fig. 7B). At last, nine potential key UDP-glycosyltransferases were screened, including 5GT, AS, LOC110100769, PF3R4, RhGT1, SB.20O07.14, UGT73C3, UGT74F1, UGT86A1.

Close

Fig. 7

Alteration of UDP-glycosyltransferase genes induced by growth years. (A) OPLS-DA analysis of UDP-glycosyltransferase genes. (B) VIP plot showing UDP-glycosyltransferase genes altered by growth years. (C) Phylogenetic analysis of 9 UDP-glycosyltransferase genes discovered from D. nobile.

Export to PPT

The UDP-glycosyltransferase genes of Arabidopsis thaliana, which is a model plant, have been well studied for their classification and functional annotation. Therefore, 20 UDP-glycosyltransferase genes from different families of Arabidopsis thaliana were selected for multiple sequence alignment with the above-mentioned nine UDP-glycosyltransferase genes from D. nobile, and a phylogenetic tree analysis was constructed for functional annotation (Fig. 7C). The results showed that all the nine Dendrobium UDP-glycosyltransferase genes were located in the branch of UDP-glycosyltransferase family I. The donor molecule of glycosyltransferase family I is uridine diphosphate glucose, which plays the role of glucose glycosylation of ligand molecules, which is consistent with the fact that sesquiterpene glycosides in D. nobile are mainly glucose glycosides.

3.6 qRT-PCR validation of DEGs related to growth years

Based on the results of the exploration of UDP-glycosyltransferase DEGs related to growth years, to confirm the dependability of the UDP-glycosyltransferase genes expression profiles for DEGs, 9 genes with the largest upregulation fold were validated by real-time quantitative RT-PCR using gene-specific primers (Table S2). Based on the qRT-PCR results, 6 of 9 candidate genes, including 5GT, AS, LOC110100769, PF3R4, SB.20O07.14, and UGT73C3, consistent with the gene expression levels obtained from RNA-seq data (Fig. 8).

Close

Fig. 8

Qrt-pcr validation of differential gene expression under growth years. The asterisk indicated significant differences (*P < 0.05), NS indicated no significant differences.

Export to PPT

3.7 Molecular docking analysis

To further predict the glycosylation catalytic activity of candidate UDP-glycosyltransferases on the sesquiterpenes of D. nobile, molecular docking was used to simulate the binding modes between 14 differential sesquiterpene glycosides in D. nobile and 6 candidate UDP-glycosyltransferase including 5GT, AS, LOC110100769, PF3R, SB.20O07.14, and UGT73C3. After docking, the conformational cluster with the lowest binding free energy is shown in the complex of the docked compounds with UDP-glycosyltransferases (Table S3, Fig. 9). Based on the binding ability of the candidate UDP-glycosyltransferases to the sesquiterpene glycosides, the results showed that these 6 candidate enzymes had good affinity with all 14 sesquiterpene glycosides. Some recent studies have shown that many natural product biosynthetic genes exist in clusters (plant biosynthetic gene clusters, BGCs) in plant genomes (Liu et al., 2020). Glycosyltransferases usually have broad-spectrum glycosylation catalytic effects (Wang et al., 2020). Therefore, the sesquiterpene glycosides in D. nobile may be formed by the joint participation of multiple UDP-glycosyltransferases.

Close

Fig. 9

Sesquiterpene glycosides in D. nobile were docking with UDP-glycosyltransferases (UGT73C3). (Dendronobiloside A, Dendronobiloside B, Dendronobiloside C, Dendronobiloside D, Dendronobiloside E, Dendromoniliside A, Dendromoniliside B, Dendromoniliside C, Dendromoniliside D, Dendroside A, Dendroside C, Dendroside E, Dendroside F, Dendroside G, in order).

Export to PPT