Experimental design and seedling culture
We used wild-type (Col-0) and cca1-1 and toc1-101 mutants (both deficient in central transcriptional-translational feedback loops of the circadian clock) of Arabidopsisthaliana to determine the role of plant circadian clock in determining the nature of the rhizospheric microbiome. The surfaces of Arabidopsis seeds were sterilized to prevent bacterial contamination of solid media and were vernalized after that. The seeds were then cultured in Petri dishes containing Murashige and Skoog (MS) medium under sterile conditions at 25 °C at a light intensity of 300 μmol photons/m2/s and a 12:12-h light:dark photoperiod. The MS medium contained 3% sugar and 0.5% agar and was autoclaved at 115 °C for 30 min before use. Two-week-old aseptic seedlings were transplanted into polycarbonate pots (400 mL) containing autoclaved potting soil (Sun Gro Horticulture, Avagam, Massachusetts, USA) and soil slurry. Approximately 30 g of grassland soil collected near the Zhejiang University of Technology, China (30°17′45.11″N, 120°09′50.07″E) was mixed with 200 mL of sterile water by vigorous shaking for 60 s. Eighteen milliliters of soil slurry were added to each pot before transplantation. Plants were grown in an artificial greenhouse at 25 ± 0.5 °C and 80% relative humidity under cool-white fluorescent light (300 μmol photons/m2/s) with a 12:12-h light:dark cycle. Samples for transcriptomic and metabolomic analyses and 16S rRNA gene sequencing were collected 14 days after transplantation. The timepoints for the analysis of circadian rhythmicity were 07:00 (one hour before the start of the light exposure, Zeitgeber time (ZT) 0), 13:00 (ZT6), 19:00 (one hour before the dark cycle, ZT12), 01:00 (next day, ZT18) and 07:00 (next day, ZT24). The details of the sample collection are provided below.
Sampling of plant shoots and RNA isolation
About 0.1 g of aboveground plant tissue was weighed at each timepoint and placed into a 1.5-mL centrifuge tube, frozen in liquid nitrogen and stored at -80 ℃. Three replicate samples were collected for each genotype for transcriptomic analysis.
TRIzol® RNA reagent was used to extract total RNA from the samples following the manufacturer’s instructions (Invitrogen, Carlsbad, California, USA). Genomic DNA was removed using DNase I (TaKaRa), RNA quality was assessed using a 2100 Bioanalyzer (Agilent) and RNA was quantified using an ND-2000 spectrophotometer (NanoDrop Technologies). Only samples with high-quality RNA (OD260/280 = 1.8-2.2, OD260/230≥2.0, RIN≥6.5, 28S:18S≥1.0, total RNA>2μg) were used to construct a sequencing library.
Library preparation and RNAsequencing
RNA purification, reverse transcription, library construction and sequencing were all performed at Shanghai Majorbio Bio-pharm Biotechnology Co., Ltd. (Shanghai, China) following the manufacturer’s instructions (Illumina, San Diego, USA). One microgram of total RNA for each genotype was used to prepare the RNA-seq transcriptomic library using the TruSeqTM RNA kit (Illumina, San Diego, USA) for sample preparation. Messenger RNA was isolated by using the polyA-selection method and oligo(dT) beads and then fragmented in a fragmentation buffer. A SuperScript double-stranded cDNA synthesis kit (Invitrogen, Carlsbad, California, USA) with random hexamer primers (Illumina) was used to synthesize double-stranded cDNA. The synthesized cDNA was end-repaired and phosphorylated, and adenines were added to the 3′ end following Illumina’s protocol for constructing libraries. Libraries were size-selected for cDNA target fragments of 200-300 bp on 2% Low Range Ultra Agarose followed by PCR amplification using Phusion DNA polymerase (New England Biolabs, Boston, Massachusetts, USA) for 15 cycles with the default parameters mentioned in the instruction of Phusion DNA polymerase. A paired-end RNA-seq library was sequenced using the Illumina Novaseq 6000 system (2×150 bp read length) after quantification using a TBS380 mini-fluorometer.
Analysis of the RNA-seq data
The sequences were trimmed and quality controlled using SeqPrep (https://github.com/jstjohn/SeqPrep) and Sickle (https://github.com/najoshi/sickle) using default parameters, and clean reads were separately aligned to the A. thaliana reference genome using TopHat (http://tophat.cbcb.umd.edu/, version2.1.1) in orientation mode. Only reads uniquely matched to the genome, allowing ≤2 mismatches without insertions or deletions, were mapped. The numbers of fragments per kilobase of exon per million mapped reads (FPKM) were calculated to represent the expression level. The ARSER algorithm [18] was used to detect circadian rhythmic genes using the parameters: minimal period, 20; maximal period, 28; default period, 24 and p = 0.05. A functional-enrichment analysis of the selected gene sets was performed to detect significantly enriched genes in the Kyoto Encyclopedia of Genes and Genomes (KEGG) relative to the whole-transcriptome background. The KEGG pathway analysis was carried out using KOBAS 2.1.1 (http://kobas.cbi.pku.edu.cn/download.php).
Preparation of exudate samples
Six replicates were taken from each genotype, and six seedlings were taken from each parallel. The plant was removed by cutting the polycarbonate pots which Arabidopsis planted in and its rhizosphere was gently nondestructively removed by shaking. The remaining soil attached to the roots was removed by washing with ddH2O, and the roots were then placed in a prepared petri dish (with an orifice plate) containing 40 mL of sterile water. Each root was allowed to secrete for 6 h each replicate. The secretions were then collected in a 50-ml centrifugal tube using a 0.45-μm aqueous filter head, quickly frozen in liquid nitrogen and stored at -80 ℃. The weight of the plant and the volume of the secretions were recorded. The samples were thawed on ice, and the exudates were extracted in a 50% methanol buffer. The exudates were centrifuged at 4000 g for 20 min, and the supernatants were transferred to new 96-well plates. Samples for quality control (QC) were prepared by pooling equal volumes (10 μL) from all samples. The samples were stored at -80 °C prior to analysis by liquid chromatography and mass spectroscopy (LC-MS).
Measurement and identification of root exudates
All chromatographic separations were performed using an ultraperformance liquid chromatograph (UPLC) (SCIEX, UK). An ACQUITY UPLC T3 column (100 mm×2.1 mm, 1.8 µm, Waters, UK) was used for the reversed-phase separation. The column oven was maintained at 35 °C, and the flow rate was 0.4 ml/min. The mobile phase consisted of solvents A (0.1% formic acid) and B (acetonitrile, 0.1% formic acid) with the following conditions for gradient elution: 0-0.5 min, 5% B; 0.5-7 min, 5 to 100% B; 7-8 min, 100% B; 8-8.1 min, 100 to 5% B and 8.1-10 min, 5% B. Exudates eluted from the column were detected using a TripleTOF5600plus high-resolution tandem mass spectrometer (SCIEX, Boston, USA) in both positive- and negative-ion modes. The Ionspray voltage floating was set at 5000 and -4500 V for the positive- and negative-ion modes, respectively. The MS data were acquired in IDA mode, and the TOF mass ranged from 60 to 1200 Da. The survey scans were acquired in 150 ms, and as many as 12 product ion scans were collected if exceeding a threshold of 100 counts per second and with a +1 charge. Total cycle time was fixed to 0.56 s. Four-time bins were summed for each scan at a pulse frequency of 11 kHz by monitoring the 40-GHz multichannel TDC detector with four-anode/channel detection. Dynamic exclusion was set for 4 s. Mass accuracy was calibrated during the acquisition every 20 samples. A QC sample was collected after every 10 samples to evaluate the stability of the LC-MS throughout the acquisition. Raw LC-MS data files were converted to mzXML format and processed using the XCMS, CAMERA and metaX toolboxes in R. Each ion was identified by combining the retention times and m/z data. Exudates with mass differences between the observed and database values of <10 ppm were annotated using the online KEGG and HMDB databases. The molecular formulae of the exudates were determined and validated by measurements of isotopic distribution.
Analysis of metabolomic data
The intensity of the peaks of the exudates was calculated by the area of the peaks and preprocessed using metaX, and peaks were removed when detected in <50% of the QC samples or <80% of the biological samples. The remaining peaks with missing values were calculated with the k-nearest neighbor algorithm to further improve the quality of the data. The preprocessed data set was used for detecting outliers and evaluating batch effects using a principal component analysis (PCA). A robust QC-based LOESS signal correction was fitted to the QC data in the order of injection to minimize the drift of signal intensity over time. Metabolic features with relative standard deviations >30% across all QC samples were removed. A principal coordinate analysis (PCoA) combined with a permutational analysis of variance (PERMANOVA) were performed to evaluate the dissimilarity between genotypes using Bray-Curtis distances. The ARSER algorithm was used to detect circadian rhythmic metabolites using the parameters: minimal period, 20; maximal period, 28; default period, 24 and p = 0.05. Spearman’s rank correlations (r) between genus were calculated using the psych package in R. Only robust (r > 0.8 or r <−0.8) and statistically significant (P-value < 0.05) correlations were incorporated into co-occurrence network. Network visualization and modular analysis were made with Gephi (v 0.9.2). Node-level topological properties were also calculated using Gephi.
Collection of rhizospheres and bulk soils
Exposed topsoil (about 1 cm thick) was discarded, and the soil around the roots was gently removed to minimize mechanical damage to the roots. The roots were then placed in a 50-mL centrifuge tube filled with 20 ml of PBS buffer (137 mM NaCl, 10 mM phosphate buffer, 2.7 mM KCl, pH 7.3-7.5). The tube was then placed on a shaker at 180 rpm for 20 min, and the roots were discarded using sterile tweezers. The washing solution was centrifuged at 3900g for 20 min. The supernatant was discarded, retaining the rhizosphere [19]. The bulk soil (soil without roots) from the same pot as each sample of rhizosphere was also collected. Samples were stored at -80 °C prior to 16S rRNA gene sequencing. Four parallel samples were collected for each genotype for sequencing.
16S rRNA gene sequencing
A FastDNA SPIN Kit for Soil (MP Biomedicals, California, USA) was used to extract total DNA from the microorganisms in the rhizosphere and bulk soil following the manufacturer’s instructions. The V3-V4 region of the 16S rRNA gene was amplified using the primers 338F (5′-ACTCCTACGGGAGGCAGCAG-3′) and 806R (5′-GGACTACHVGGGTWTCTAAT-3′) and a PCR thermocycler (GeneAmp 9700, ABI, USA). Purified amplicons were sequenced on an Illumina MiSeq platform (Illumina, San Diego, California, USA). Sequences with ≥97% similarity were assigned to the same operational taxonomic unit (OTU). Taxa were assigned using the SILVA database (https://www.arb-silva.de/) based on the RDP classifier (Version 2.11, https://sourceforge.net/projects/rdp-classifier/) algorithm. The microbial compositions were also assessed using the UniFrac distance, which is a phylogenetically based measure of the degree of similarity between microbial communities [20]. Pairwise distances were determined for all pairs of samples based on either the taxonomic representation (unweighted) or the types and relative abundances (weighted). A PCoA of the UniFrac distances was used to visualize the differences among samples. A PERMANOVA was performed to evaluate the variance between genotypes. The ARSER algorithm was used to detect genera with circadian rhythms using the parameters: minimal period, 20; maximal period, 28; default period, 24 and p = 0.05.
Analysis of transcriptomes, root exudates and 16S rRNA gene sequences
The Procrustes test was used to analyze the correlations between rhythmic genes in the wild type and all exudates using the vegan package in R (3.6.2). Multiple linear regression was used to analyze the regulation of all exudates from the 10 most abundant phyla (Proteobacteria were divided into three classes) using SPSS (v 20.0.0).