Plant materials
Ziziphus jujuba spinosa were cultivated in the germplasm resources center of jujube at the Baijiazhuang Village, Taizijing Township, Xingtai City, Hebei Province, China. It is a common species that belongs to the plant family Rhamnaceae. The seed of Ziziphus jujuba spinosa (SZS) is oblate or oblong, and the surface is reddish brown or purplish red. We collected SZS from the same tree in the germplasm resources center of jujube every two weeks, on 15 August (T1), 1 September (T2), and 15 September (T3) 2020. We confirm that the collection of plant material is compliant with relevant institutional, national, and international guidelines and legislation. All SZS samples were frozen using liquid nitrogen and stored at -80 °C for RNA and metabolite extraction.
Metabolite extraction
Freeze-dried SZS samples were crushed using zirconia beads and a mixer mill (MM 400, VERDER RETSCH, Shanghai, China). The powder was weighed and 100 mg was resuspended in 1.0 mL 70% aqueous methanol and incubated overnight at 4 °C. The solution (i.e., extract) was centrifuged at 10,000 g for 10 min. The supernatant was added to the CNWBOND Carbon-GCB SPE Cartridge (250 mg, 3 mL; ANPEL, Shanghai, China; www.anpel.com.cn/) and then filtered (SCAA-104, 0.22 μm pore size; ANPEL) before UPLC-MS/MS analysis.
Ultra performance liquid chromatography (UPLC) conditions
The SZS sample extracts were analyzed using an UPLC-ESI-MS/MS system (UPLC, Shim-pack UFLC SHIMADZU CBM30A system, Shanghai, China). The analytical conditions were as follow: UPLC: column, Waters ACQUITY UPLC HSS T3 C18 (1.8 µm, 2.1 mm × 100 mm); solvent system, water (0.04% acetic acid): acetonitrile (0.04% acetic acid); gradient program, 100:0 (v/v) at 0 min, 5:95 (v/v) at 11.0 min, 5:95 (v/v) at 12.0 min, 95:5 (v/v) at 12.1 min, 95:5 (v/v) at 15.0 min; flow rate, 0.40 mL/min; temperature, 40°C; and injection volume: 2 μL. After UPLC, the effluent was alternatively connected to an ESI-triple quadrupole-linear ion trap (Q TRAP)-MS.
ESI-Q TRAP-MS/MS
The LIT and triple quadrupole (QQQ) scans were acquired using the API 6500 Q TRAP LC-MS/MS system equipped with an ESI Turbo Ion-Spray interface. The system was operated in the positive ion mode and was controlled using the Analyst 1.6 software (AB Sciex). The ESI source operation parameters were as follows: ion source, turbo spray; source temperature, 500 °C; ion spray voltage, 5,500 V; ion source gas I, gas II, and curtain gas were set at 55, 60, and 25.0 psi, respectively; and the collision gas was high. Instrument tuning and mass calibration were performed using 10 and 100 μmol/L polypropylene glycol solutions in the QQQ and LIT modes, respectively. The QQQ scans were acquired in MRM experiments, with the collision gas (nitrogen) set at 5 psi. The declustering potential and collision energy were optimized for individual MRM transitions. A specific set of MRM transitions was monitored for each period according to the metabolites eluted within the period.
Identification and quantitative analysis of metabolites
A publicly available standard metabolite database (Metware Biotechnology Co., Ltd., Wuhan, China), which was built on the basis of multiple ion monitoring–enhanced product ions, was used for the qualitative analysis of metabolites. The quantitative analysis of metabolites was completed using a published multiple reaction monitoring method37,38. An unsupervised principal component analysis, hierarchical cluster analysis, and partial least-squares discriminant analysis were performed using the prcomp function within R (www.r-project.org). Metabolites that were significantly differentially abundant between groups were identified according to the following criteria: variable importance in projection value ≥ 1 and an absolute log2(fold change) ≥ 1.
RNA extraction and Illumina sequencing
Total RNA was extracted from frozen leaves using the RNAprep Pure Plant Kit (TIANGEN Biotech, Beijing, China), after which RNA degradation and contamination were monitored by 1.2% agarose gel electrophoresis. The purified RNA concentrations were quantified using the NanoDropTM 1000 spectrophotometer (ThermoFisher Scientific, Shanghai, China). Poly (A) mRNA was enriched from the total RNA using Oligo (dT) magnetic beads. The poly (A) mRNA was subsequently fragmented by an RNA fragmentation kit (Ambion, Austin, TX, USA). The fragmented RNA was transcribed into first-strand cDNA using reverse transcriptase and random hexamer primers. Second-strand cDNA was synthesized using DNA polymerase I and RNase H (Invitrogen, Carlsbad, CA, USA). After end repair and the addition of a poly (A) tail, suitable length fragments were isolated and connected to the sequencing adaptors. The fragments were sequenced on the Illumina HiSeqTM 2500 platform39.
RNA sequencing (RNA-seq) data analysis and annotation
To acquire high-quality reads, the raw reads in fastq format were processed through in-house Perl scripts. Clean reads were obtained from raw data by removing adaptor sequences, low-quality reads, and reads containing ploy-N. All downstream analyses were based on clean data with high quality. Transcriptome assembly was accomplished using Trinity software (version2.5.1)40. Genes were functionally annotated according to the Kyoto Encyclopedia of Gene and Genome (KEGG) pathway database (http://www.genome.jp/kegg), the NCBI non-redundant (Nr) database (version 2018.4), the Swiss-Prot protein database (http://www.uniprot.org), the euKaryotic Clusters of Orthologous Groups (KOG) database (version 1.0), the Gene Ontology (GO) database (http://www.geneontology.org), and the Pfam database (version 33.0).
Gene expression levels were estimated by RSEM (version 1.2.26)41. Analysis of DEGs of the two groups was performed with the DESeq R package (1.10.1). DESeq provides statistical routines for determining DEGs using a model based on the negative binomial distribution. The results of all statistical tests underwent a Benjamini and Hochberg false discovery rate multiple testing correction. Genes were determined to be significantly differentially expressed at an adjusted P-value < 0.05. The GO enrichment analysis of the DEGs was performed using the topGO R package based on the Kolmogorov-Smirnov test. The pathways significantly associated with the DEGs were identified using the KEGG database (http://www.genome.jp/kegg)42,43,44 and KOBAS software45.
QRT-PCR expression analysis of genes involved in flavonoid biosynthesis
Total RNA of SZS was reverse-transcribed using the Quantscript Reverse Transcriptase Kit (TIANGEN Biotech, Beijing, China). The resulting cDNA was used as the template for analyzing gene expression by qRT-PCR. The specific primers involved in the flavonoid biosynthesis genes and actin gene (internal control) are listed in Table S3. The qRT-PCR analysis was performed using SYBR® Premix Ex Taq™ II (TaKaRa; http://www.takarabiomed.com.cn) and the ABI Prism 7500 system, including the software for the 7500 and 7500 Fast Real-Time PCR systems (version 2.0.1) (Applied Biosystems, Foster City, CA, USA). The 2−ΔΔCT method was used to quantify gene expression levels46.
Statistical analysis
Each experiment was three biological replicates. Data were analyzed using Microsoft Office Excel 2013 software. Data were presented as the means ± standard deviations. The significance of any differences was analyzed using p < 0.05 as the threshold.