Effects of the mycorrhizal fungus Ceratobasidium sp. AR2 on growth and flavonoid accumulation in Anoectochilus roxburghii

Isolation of mycorrhizal fungi

The roots of wild adult A. roxburghii were collected from tree hoods under the mountain streams in Dehua County, Fujian Province, southeast China in November 2017. Mycorrhizal fungi were isolated following the methods of Hoang et al. (2017) and Richardson, Currah & Hambleton (1993), with some modifications. Roots were rinsed with tap water to remove soil and impurities, and then treated with 1% sodium hypochlorite disinfectant for 1 min. Roots were disinfected for 30 s with double antibiotics (penicillin–streptomycin, both at a final concentration of 50 μg/mL) and then rinsed three times for 1 min each in sterile distilled water. The root surface was scraped using an anatomical needle, and mycelium pelotons were then released into the sterile distilled water. The mycelium pelotons were transferred to a small piece of double-resistant potato dextrose agar (PDA) medium with a pipette gun under a GL6545 stereomicroscope (Guilin Guiguang Instrument Co. Ltd., Guangxi, China), and the plates were incubated in a dark incubator at 22 °C. Samples were then incubated until hyphae were visible on the medium. Pure culture was obtained by transferring the hyphae onto PDA. The isolated strains were then taken to the Institute of Medical Plant Development, Chinese Academy of Medical Sciences/Peking Union Medical College (maintained at 4 °C and −80 °C).

Coculture of plant and mycorrhizal fungi

Tissue culture plantlets of A. roxburghii (approximately 3 months old and 2–3 cm in height) were used for coculture. Plantlets were transplanted into tissue culture bottles (9 cm in diameter, 12.5 cm in height) containing 125 mL H1 oat medium agar (H1OMA; 200 mg Ca(NO₃)₂·4H₂O; 100 mg KCl; 200 mg KH₂PO₄; 100 mg MgSO₄·7H₂O; 100 mg yeast extract; 2 g sucrose; 3 g finely ground oats; 5 g agar; and 1 L H₂O) (Clements, Muir & Cribb, 1986). Twelve bottles were used for each group (control group with no fungus and treatment group with fungus). One mycelial plug (approximately 1 cm × 1 cm) from 7-day-old fungi grown on PDA medium was placed into the center of each culture bottle, and four plantlets were then planted. The culture bottles were then kept in a conventional greenhouse with a 12-h light/12-h dark photoperiod at 24 °C ± 1 °C and an illumination intensity of 1,500 Lx. For contrastive evaluation of leaf color between the control and treatment groups, the RGB color pattern in Photoshop software was used to analyze the R, G, and B values of leaf color from both groups (Afshari-Jouybari & Farahnaky, 2011).

After 4 months, the root number, fresh weight, node number and tiller number of each plantlet were recorded. Section I tissues (the upper plantlets of the medium from the control group) and Section II tissues (the upper plantlets of the medium from the treatment group) were frozen with liquid nitrogen and stored at –80 °C until RNA-seq, qRT-PCR, and target metabolism analyses. Each treatment contained three replicates with six plantlets per replicate.

Histological study

One week after inoculation, fresh root segments were fixed in a solution of 2.5% glutaraldehyde and 1.6% paraformaldehyde in 50 mM phosphate buffer (pH 6.8) for 4 h at room temperature. After fixation, the samples were washed three times with phosphoric acid buffer for about 15 min each time. The samples were then dehydrated through a graded ethanol series (15% ethanol for 30 min; 30% ethanol for 30 min; 50% ethanol for 30 min; 70% ethanol for 1 h; 85% ethanol for 1 h; 95% ethanol for 1 h; absolute ethanol for 1 h). After dehydration, the samples were embedded in LRwhite gradient mixture (25% LRwhite for 24 h; 50% LRwhite for 24 h; 75% LRwhite for 24 h; 100% LRwhite for 24 h) and polymerized for 48 h at 60 °C. Next, 3-µm-thick sections were cut using glass knives with a Reichert-Jung 2040 Autocut rotary microtome. The sections were collected on slides and stained with 0.05% (w/v) toluidine blue O in benzoate buffer for general histology examinations. The sections were examined, and images were captured digitally using a digital camera attached to the microscope (Axio Imager A1; Carl Zeiss, Oberkochen, Germany).

Molecular identification of fungus

After 7 days in PDA medium, an E.Z.N.A. Fungal DNA Mini Kit (Omegabiotek, Norcross, GA, USA) was used to extract genomic DNA from mycelia. The internal transcribed spacer (ITS) regions were amplified using primers ITS1 (5′-TCCGTAGGTGAACCTGCGG-3′) and ITS4 (5′-TCCTCCGCTTATTGATATGC-3′) (Zhou, Tang & Guo, 2018). PCR was performed in a 50 μL reaction volume containing 2 μL template DNA (200 ng/μL), 1 μL each primer (10 μM), 25 μL Taq PCR mix, and 21 μL sterile water. The reactions were performed in a Bio-Rad T100 Thermal Cycler (Bio-Rad, Hercules, CA, USA) using the following program: initial denaturation at 95 °C for 3 min; 35 cycles of denaturation at 94 °C for 1 min, annealing at 55 °C for 30 s, and extension at 72 °C for 1 min; and a final extension at 72 °C for 7 min. Amplification products were analyzed in 1% agarose gels (mixed with golden view) by electrophoresis and directly sequenced using a GENEWIZ sequencer (Suzhou, China). The sequence was submitted to the GenBank database (Accession Number: MN068847). For phylogenetic analysis, ITS sequences of the AR2 isolate and 18 Ceratobasidium sequences of known species from GenBank were analyzed, and the sequence of Waitea circinata was used as the outgroup taxon (Paduano et al., 2011). DNA sequences were aligned using Clustal W algorithm in MEGA 7.01, followed by manual adjustment. Distance trees were obtained using the neighbor-joining (NJ) method (Saitou & Nei, 1987) with the Tajima-Nei method. For assessing the relative robustness of branches, the bootstrap method was used with 1,000 replicates (Felsenstein, 1985).

Total RNA extraction

Total RNA from each sample was extracted with an RNeasy Plant Mini Kit (Qiagen, Hilden, Germany) according to the manufacturer’s instructions, followed by treatment with an RNase-free DNase I digestion kit (Aidlab, Beijing, China) to remove genomic DNA contamination. RNA quality and concentration were verified by 1% agarose gel electrophoresis and NanoDrop 2000 spectrophotometry (Thermo Fisher Scientific, Waltham, MA, USA), respectively. RNA integrity was assessed on an Agilent 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA). Library construction and sequencing were performed on an Illumina Hiseq Xten platform. Sequence data were deposited in the NCBI database (Accession Number: PRJNA562071).

cDNA library construction, sequencing, assembly and annotation

mRNA was isolated from total RNA with oligo (dT)-attached magnetic beads and randomly fragmented into short fragments using fragmentation buffer. With cleaved RNA fragments as template, the first cDNA strand was synthesized by random hexamers, and the second cDNA strand was synthesized by adding buffer, dNTPs, RNaseH, and DNA polymerase I. Double-stranded cDNA was then purified using a QIAEX II Gel Extraction kit (Qiagen). The ends were trimmed, and a poly (A) tail was added at the sequencing joint. PCR amplification was then performed, and PCR products were purified using an AMPure XP Kit (Beckman Coulter, Brea, CA, USA). Library quality was assessed on an Agilent 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA).

Deep transcriptome sequencing was performed using the Illumina Hiseq Xten. Raw data were filtered by removing the connector sequence and low-quality reads to yield high-quality clean data. Sequence alignment was conducted between the clean data for each sample and a transcript or unigene library was assembled. The obtained transcript and unigene reads were categorized as mapped reads and were used for subsequent analysis. Unigene function was annotated based on the following six databases: NCBI non-redundant protein sequences (Nr), NCBI non-redundant nucleotide sequences (Nt), Clusters of Orthologous Groups of proteins (COG), a manually annotated and reviewed protein sequence database (Swiss-Prot), Kyoto Encyclopedia of Genes and Genomes (KEGG), and Gene Ontology (GO).

Confirmation of flavonoid biosynthesis-related genes

Unigenes related to flavonoid biosynthesis were validated with qRT-PCR. RNA from the two groups (I and II) was reverse-transcribed into cDNA with a PrimeScript RT reagent kit (TaKaRa Bio, Shiga, Japan). qRT-PCR was performed with SYBR Premix Ex Taq (TaKaRa Bio, Shiga, Japan) on a LightCycler 480 II Real-Time PCR System (Roche Diagnostics, Basel, Switzerland). Three biological replicates and three technical replicates were analyzed, and all primer names and corresponding sequences are shown in Table S1. The elongation factor 1 alpha (EF-1α) gene in A. roxburghii was used as a reference control. The PCR conditions were as follows: denaturation at 95 °C for 30 s, followed by 40 cycles of amplification (95 °C for 5 s, 60 °C for 30 s). Melting curves were prepared by measurement at 95 °C for 5 s and 60 °C for 1 min. For each sample, three technical replicates were used, with a minimum of three biological replicates. Gene expression was evaluated using the 2^−ΔΔCt method (Livak & Schmittgen, 2001).

Extraction and analysis of total and main flavonoids

Flavonoids were extracted using an ultrasonic-assisted method. Each sample (approximately 20 mg) was soaked in 400 μL of 90% methanol for approximately 2 h according to the water-to-raw material ratio (20 mL/g) and then extracted three times for 20 min each time in an ice-water bath using the KH5200DE ultrasonic instrument (Kunshan Hechuang Ultrasonic Machinery Co., Jiangsu, China). Samples were then centrifuged at 10,000×g for 5 min, and each extract was collected and pooled.

The contents of total flavonoids were determined using colorimetry. Briefly, 150 μL of 5% NaNO₂ was added to one mL extract in a 10 mL volumetric tube, and the mixture was incubated for 6 min at room temperature. Next, 150 μL of 10% AlCl₃·6H₂O was added to the mixture, and the mixture was incubated for an additional 6 min, followed by the addition of 2 mL of 4% NaOH. After 15 min of incubation at room temperature, the absorbance at 510 nm was measured using an EnSpire Multimode Plate Reader (PerkinElmer, Waltham, MA, USA). Total flavonoid contents were calculated using a standard rutin curve. The linear regression equation was y = 1.9322x + 0.0501 (correlation coefficient (R²) = 0.9997).

To determine changes in monomer flavonoids, several available flavonoids known to be present in A. roxburghii were tested. Rutin (153-18-4; 98% purity), quercetin (117-39-5; 98% purity), isorhamnetin (480-19-3; 98% purity), isorhamnetin-3-o-rutinoside (604-80-8; 98% purity), kaempferol (520-18-3; 98% purity), and cyanidin-3-glucoside chloride (7084-24-4; 98% purity) were purchased from National Institutes for Food and Drug Control (Beijing, China). Measurement of the concentrations of these flavonoids was performed using a 1290 Infinity II series UHPLC System connected to the 6460 Triple Quadrupole Mass Spectrometer (Agilent Technologies, Santa Clara, CA, USA). UHPLC separation was completed using an Agilent 1290 Infinity II series UHPLC System equipped with a Waters ACQUITY UPLC HSS T3 column (100 × 2.1 mm, 1.8 μm). The mobile phase A was 0.1% formic acid in water, and the mobile phase B was methanol. The elution gradient was as follows: 0–1 min, 80–80% A; 1–9.5 min, 80–5% A; 9.5–11.2 min, 5–5%; 11.2–11.5 min, 5–80%; 11.5–15 min, 80–80%. The flow rate was 0.3 mL/min, and the injection volume was 1 μL. The column and autosampler temperatures were 35 °C and 4 °C, respectively. An Agilent 6460 triple quadrupole mass spectrometer (Agilent Technologies, Santa Clara, CA, USA) equipped with an AJS electrospray ionization interface was used for assay development. Typical ion source parameters were as follows: capillary voltage = +4,000/−3,500 V, nozzle voltage = +500/−500 V, gas (N₂) temperature = 250 °C, gas (N₂) flow = five L/min, sheath gas (N₂) temperature = 250 °C, sheath gas flow = 11 L/min, nebulizer = 45 psi. The six target flavonoids were used to obtain an external calibration curve, which was used for calculating the content of each flavonoid in each sample.

Statistical analysis

SPSS 12.0 software (IBM, Chicago, IL, USA) was used for data analysis, and values are presented as the mean ± standard deviation. An independent Student’s t-test was used to compare the effects of mycorrhizal fungi on the response variables considered in this study. A value of P ≤ 0.05 was considered statistically significant.

Fungal isolation and identification

AR2 cultures showed rapid growth (9–10 mm/day) on PDA medium. The colonies were yellowish-white in color, and the mycelium grew densely (Fig. 1A). Colony microstructure observations showed that the mycelia had right-angled branches and hyphae constriction at branching points (Fig. 1B). There was mycelial fusion between hyphae: short lateral branches were produced between two hyphae, and short lateral branches grew in opposite directions to complete mycelial fusion (Fig. 1C). Ceratobasidium and Thanatephorus are both members of genus Ceratorhiza (Dearnaley, Martos & Selosse, 2012). The main difference between these organisms is that the mycelial cells of the former are binucleate, whereas those of the latter are multinucleated. Our observation of the AR2 nuclear chamber showed that its hyphal cells were binucleate (Fig. 1D). A BLAST search of the ITS sequences based on the extracted rDNA showed 98% identity with the ITS sequences of Ceratobasidium sp. Moreover, the results of the NJ phylogenetic tree analysis showed that AR2 was closely related to Ceratobasidium sp. and aggregated into the same branch with Ceratobasidium sp. (KJ495974) isolated from A. formosanus (Jiang et al., 2015) (Fig. 1E). Thus, AR2 isolated from the roots of wild A. roxburghii was identified to be a member of the genus Ceratobasidium based on morphological features and phylogenetic tree analysis of the ITS sequences.

Variations in biomass parameters in A. roxburghii

After 1 week of coculture, the roots in the AR2 treatment group showed a typical structure (intracellular pelotons) of orchid mycorrhizae (Fig. 2A). After 4 months, the treatment group showed better growth than the control group (Fig. 2B), and new buds appeared in the treatment group (Fig. 2C). Additionally, the leaf color differed significantly between the two groups; the leaves of mycorrhizal A. roxburghii were redder than those of nonmycorrhizal plants (Fig. 2D; Table S2). The growth indices of A. roxburghii infected with AR2 are shown in Table 1. The fresh weight, dry weight, height, root number, and bud number were higher in the mycorrhizal plants than those in the control plants (Table 1).

Transcriptome sequencing, de novo assembly and unigene annotation

Orchidaceae mycorrhizal infection is often confined to underground roots and rhizomes (Singh & Varma, 2000). In order to avoid the interference of genes and compounds from mycorrhizal fungi, we chose uninfected portions of A. roxburghii above the culture medium for analysis. cDNA libraries were constructed using Section II of A. roxburghii inoculated with AR2 and Section I of normal A. roxburghii as a control (Fig. 2A). We sequenced and pooled equivalent quantities of RNA extracts isolated from the two sections, and all samples were sequenced in three biological replicates. After removing adaptors, ambiguous nucleotides, and low-quality sequences, 247.76 Mb clean reads (total data: 37.17 Gb) were obtained. Q20 and Q30 base percentages were greater than or equal to 95.68% and 89.28%, respectively. The GC contents of Section I and Section II were in the ranges of 47.62–47.74% and 48.76–49.00%, respectively (Table S3). The assembled sequences showed that the unigenes were in the range of 61,217–67,640 bp with an N50 length of 1,242–1,440 bp and 42.81–44.02% GC content, whereas transcripts were in the range of 92,355–105,866 bp with an N50 length of 1,140–1,330 bp and 42.79–44.07% GC content (Table S4). These data indicated that the constructed transcriptome library had high-quality sequencing data and could be used for subsequent analyses, providing abundant data resources for subsequent gene research.

Among the 133,321 identified unigenes, 63,136 (47.36%), 31,497 (23.62%), 46,082 (34.56%), 45,505 (34.13%), 25,349 (19.01%) and 26,447 (19.84%) were aligned to the Nr, Nt, SwissProt, KEGG, COG, and GO databases, respectively (Table S5). According to GO analysis, 130,618 non-redundant unigenes were classified into three major functional ontologies: biological process, cellular component, and molecular function (Fig. 3). For biological process, the dominant subcategories were cellular process (13,243) and metabolic process (13,595). In the category of cellular component, cell (10,254), cell part (10,160), and membrane (8,516) were the most highly represented. Among molecular function terms, binding (11,457) and catalytic activity (13,083) were the most highly represented. Moreover, within each of these three categories, a few genes were assigned the subcategories of “cell killing,” “extracellular matrix component,” and “protein tag.”

In addition, all unigenes were searched against the COG database for functional prediction and classification. In total, 36,198 nonredundant unigenes were assigned to COG classifications and clustered into 24 functional categories (Fig. 4). Among these, the cluster of general function prediction only (6,366) was the largest group, followed by replication, recombination, and repair (3,441); transcription (3,249); post-translational modification, protein turnover, and chaperones (2,713); and signal transduction mechanisms (2,646). Only a few unigenes were assigned to the subcategories of cell motility (13) and nuclear structure (10). Furthermore, 869 genes were classified as having unknown functions. These GO and COG annotations showed that the expressed unigenes in A. roxburghii encoded many metabolism-related proteins.

Analysis of the metabolic pathways modulated by unigenes was performed using KEGG annotations. According to KEGG analysis, 68,046 unigenes were assigned to different pathways (Fig. 5). The results showed that metabolic pathways were the most enriched, followed by “genetic information processing.” The subcategory of human diseases was the least enriched. Flavonoid (ko00941), flavone and flavonol (ko00944), and anthocyanin metabolism (ko00942) were among the detected metabolic pathways, and 109 unigenes encoding proteins involved in flavonoid biosynthesis were identified (Table S6), indicating that the expression levels of these genes could be responsible for flavonoid accumulation.

Analysis of differentially expressed genes

Genes with |log₂(fold change)| > 1 and q < 0.001 were defined as differentially expressed genes (DEGs). In this study, the six samples were divided into two groups, that is, groups I (I1, I2, I3) and II (II1, II2, II3). The different expression patterns of genes were then analyzed between groups I and II. In total, 9,645 DEGs were identified between groups I and II, with 6,223 up-regulated genes and 3,422 down-regulated genes. GO annotations showed that 14,669 genes were classified into 50 GO terms (Fig. 6), with metabolic process and catalytic activity as the dominant terms. KEGG pathway analysis of the DEGs from groups I and II revealed 136 enriched pathways. These pathways were involved in six processes, including cellular processes, environmental information processing, genetic information processing, human diseases, metabolism, and organismal systems (Fig. 7). In the top 20 pathways, the number of unigenes involved in metabolic pathways was higher (Table S7). Moreover, 23.90% and 14.95% of the top 20 KEGG pathways were involved in metabolic pathways and biosynthesis of secondary metabolites, respectively. These results suggested that AR2 had important effects on the metabolism of secondary metabolites in A. roxburghii.

Relative expression of flavonoid synthesis-related genes

Next, we evaluated the mechanisms of flavonoid accumulation in mycorrhizal A. roxburghii. The metabolic pathways of six common flavonoids in A. roxburghii are shown in Fig. 8, and 14 structural genes were found to be involved in its synthesis. We obtained nine phenylalanine ammonia lyase (PAL) genes, two cinnamate-4-hydroxylase (C4H) genes, twelve 4-coumarate CoA ligase (4CL) genes, seven chalcone synthase (CHS) genes, four chalcone isomerase (CHI) genes, fifteen flavanone 3-hydroxylase (F3H)/flavonol synthase (FLS) genes, ten flavonoid 3′-hydroxylase (F3′H) genes, seven dihydroflavonol 4-reductase (DFR) genes, seven leucoanthocyanidin dioxygenase (LDOX) genes, six O-methyltransferase (OMT) genes, twelve anthocyanin 3-O-glucosyltransferase (GT1) genes, seven 3-O-rutinoside-glucosyltransferase (GT2) genes, nine UDP-glucosyltransferase (UFGT) genes, and two 3-O-glucoside rhamnosyltransferase (UFRT) genes from transcriptome databases (Table S6). Among these genes, fourteen were analyzed using qRT-PCR to validate the accuracy of the transcriptome data. The chosen genes were all expressed in both sections to a greater or lesser extent (Fig. 9). PAL (CL13945. contig6), 4CL (Unigene 36812), CHS (CL3547. contig2), F3′H (CL5534. contig1), DFR (CL4355. contig1), LDOX (CL2219. contig2), OMT (CL13313. contig2), and UFGT (CL13667. contig1) were significantly upregulated, while the genes CHI (CL6389. contig4), F3H/FLS (CL13017. contig3), and GT2 (Unigene3073) were significantly downregulated in Section II compared with Section I. Furthermore, C4H (CL7400. contig2), GT1 (CL449. contig1), and UFRT (Unigene1776) showed no significant changes between Section I and Section II. These validation results were consistent with the obtained transcriptome data, indicating that RNA-seq provided reliable data for the identification of key genes related to flavonoid metabolism in A. roxburghii.

Contents of total and selected flavonoids of A. roxburghii

Using A. roxburghii with uninoculated AR2 as the control, the contents of total flavonoids and six known flavonoid monomers (quercetin, rutin, isorhamnetin, isorhamnetin-3-O-rutinoside, kaempferol, and cyanidin-3-glucoside chloride) in A. roxburghii were measured by colorimetric analysis and UHPLC-MS-MS, respectively. The recovery rates of quercetin (94.9%; relative standard deviation% (RSD%) = 2.7%), rutin (98.9%; RSD% = 11.2%), isorhamnetin (108%; RSD% = 4.3%), isorhamnetin-3-O-rutinoside (107.1%; RSD% = 4.4%), kaempferol (106.3%; RSD% = 2.6%), and cyanidin-3-glucoside chloride (105.0%; RSD% = 4.7%) were determined by comparing them with an external standard method. The multiple reaction monitoring parameters and calibrations for these standards are shown in Table 2. The R² was higher than 0.9986. The lower limits of detection ranged from 1.53 nM to 24.41 nM, and the lower limits of quantitation were from 3.05 nM to 48.83 nM. As shown in Table 3, mycorrhizal fungi AR2 could significantly increase the accumulation of total flavonoids of A. roxburghii. Rutin, isorhamnetin, and cyanidin-3-glucoside chloride contents were significantly higher after inoculation with AR2 than in the control, showing increases of 62.59%, 260.00%, and 214.29%, respectively. In contrast, the yields of quercetin and isorhamnetin-3-O-rutinoside in Section II were lower than those in Section I. The content of kaempferol in Section I was 0.79 ± 0.01 nmol/g, although this compound was undetectable in Section II.