Greenhouse experiments and sampling
Wheat (Triticum turgidum cv. Negev) was cultivated in sandy loam soil (19% clay, 6% silt, 75% sand) classified as Calcic Haploxerept. The soil was obtained from intensive agriculture field located in Eshkol region, Israel (31.248,949, 34.379,872). Potatoes, wheat and peanuts were previously grown in this field. Initial soil parameters were: pH 8.78 ± 0.04, electrical conductivity 99 ± 1 (µS/m), NO3-N 0.22 ± 0.02 (mg/kg), NH4 0.30 ± 0.01 (mg/kg), P-PO4 0.09 ± 0.01(mg/kg), total soluble organic carbon 4.0 ± 0.04 (mg/kg) and total soluble nitrogen 0.70 ± 0.02 (mg/kg).
The plants were grown for 6 weeks (from December 2016 to February 2017) as described previously [25]. Briefly, 750 g of soil was distributed in a 700-mL plastic pot, with four seeds per pot. Those pots were able to sustain up to four wheat plants for six weeks under the experimental conditions. The wheat was grown in a greenhouse with two closed-system chambers at day/night temperatures of 25 °C/18 °C ± 1 °C, and with an automatically adjusted CO2-supply system (Emproco Ltd., Ashkelon, Israel). The photoperiod was 9 h and the daily light integral was 12.5 MJ/day. Wheat plants were grown in a sequence of three independent experimental cycles of 6 weeks each (five pots per treatment per cycle), with a 1-week shift between cycles. Plants were grown under either ambient (400 ppm) or elevated (850 ppm) atmospheric CO2 levels. Nutrient solution was prepared with 90% nitrogen supplied as nitrate and 10% supplied as ammonium using KNO3 and NH4NO3 to provide final concentrations of 30, 70 and 100 ppm nitrate [26]. Other macronutrients were supplied in each treatment at the following rate: P-15 ppm, K-150 ppm, Mg-24 ppm, Ca-120 ppm and S-40 ppm provide by NH4NO3, KNO3, CaCl2, KCl, MgCl2 and KH2PO4 salts. 40 ppm S and Ca were present in the tap water. Micronutrients were supplied at a rate of 1.3 ppm Fe, 0.7 ppm Mn, 0.3 ppm Zn, 0.05 ppm Cu, and 0.0375 ppm Mo using Korotin (Haifa Chemicals, Israel), a commercial micronutrient mix. Each pot was irrigated with 50 mL of the nutrient solution four times a week. The total amount of nitrogen in the 30 ppm nitrate treatment was 36 mg/pot (equivalent of ca. 73 kg N/ha), 70 ppm nitrate treatment was 84 mg/pot (equivalent of ca. 170 kg N/ha) and in the 100 ppm treatment, 120 mg/pot (equivalent of ca. 250 kg N/ha).
Soil and plant analyses
At the end of the 6th week of growth, 15 pots (5 pots per cycle) from each treatment were sampled for soil, shoots and roots, and the following parameters were measured: soil nitrate and ammonia content, soil EC and soil pH, shoot and root dry biomass, nitrogen concentration and content in shoot and roots. Soil properties and relevant methods were as described previously [25]. Briefly, soil EC and pH were determined in a solution of 1:5 air dry sieved soil:distilled water (w/v). Nitrate and ammonium concentrations were determined using an autoanalyzer (Lachat Instruments, Milwaukee, WI or Gallery Plus, Thermo Fisher Scientific, Waltham, MA, USA). Sampled shoots and roots were dried at 60 °C for 48 h, ground and weighed to obtain dry biomass. Total nitrogen concentration was determined using an autoanalyzer (Lachat Instruments or Gallery Plus) following digestion with sulfuric acid and peroxide [27].
Root DNA extraction for sequencing and qPCR
At the end of the 6th week of wheat growth, pots were randomly selected for DNA extraction. To obtain the root-surface-associated microbiome, wheat roots were collected in triplicate from each of the three cycles and were vortexed three time with 85% saline solution, until no visible soil particles were attached to the roots. Total DNA was extracted from 0.4 g of complete root system, using the Exgene Soil DNA mini isolation kit (GeneAll, Seoul, Korea) according to the manufacturer’s instructions.
Generation of qPCR plasmid standards
Plasmids containing the 16S rRNA gene were generated as described previously [28, 29]. Each PCR amplification product was ligated into pGEM-T Easy Vector (Promega, Madison, WI, USA) and plasmids were transformed into BioSuper Escherichia coli DH5α competent cells (Bio-Lab, Jerusalem, Israel). Circular plasmid DNAs were used as the standards to create calibration curves at 10-fold dilutions for gene quantification by real-time qPCR.
Assessment of gene copy numbers by qPCR
Copy numbers of the total bacterial community (16S rRNA gene) and translation elongation factor 1 (TEF, a plant housekeeping gene) were assessed using selected primers (Table S1) in roots of 6-week-old wheat plants with the StepOnePlus Real-Time PCR System (Applied Biosystems, Foster City, CA, USA). Triplicates from whole genomic DNA were diluted to 6 ng/µL and 1 µL was used in a 20-µL final reaction volume together with 50 µM forward and reverse primers and 10 µL 1X FAST MasterMix (Thermo Fisher Scientific). Three biological and three technical replicates were analyzed for each root DNA sample. Reaction efficiency was monitored in each run by means of an internal standard curve (constructed plasmids) using duplicates of 10-fold dilutions of standards ranging from 108–102 copies per reaction. Efficiency was 89–98% for all target genes and runs, and R2 values were greater than 0.99. Copy numbers of the target genes were calculated based on the relative calibration curve of the plasmid copy numbers. All data analyses were conducted using StepOne software v2.3 (Applied Biosystems).
Shotgun sequencing
Root DNA was extracted from each of the biological triplicates, in each of the three cycles. For sequencing, the DNA of the triplicates was combined, resulting in three biological replicates per treatment (one from each batch) and 18 samples altogether. Shotgun metagenome libraries were prepared using the Celero DNA-Seq library preparation kit (NuGen, Takara Bio, USA) with enzymatic shearing, according to the manufacturer’s instructions. All libraries were then pooled in equal volumes and size selection (350–400 bp fragments) was performed using a Blue PippinPrep instrument (Sage Scientific). The libraries were then sequenced using an Illumina MiniSeq instrument employing a mid-output kit. Based on the number of reads per sample, the samples were repooled with varying volumes, and size selection was performed again using the same size range. The final size-selected pool was sequenced on an Illumina NovaSeq instrument with an S4 flow cell, employing 2 × 150 base reads. Library preparation and pooling were performed at the University of Illinois at Chicago Sequencing Core (UICSQC), and sequencing was performed by Novogene Corporation (Chula Vista, CA, USA).
In total, we obtained 310 Gb of information, with 30–44 million sequences per root sample. These sequence data were submitted to the Sequence Read Archive (SRA) of the NCBI databases under accession numbers SUB6631533 and SUB8385777, BioProject: PRJNA592741.
All reads were subjected to quality control using FastQC v0.11.3 [30] and barcode trimming using Trimmomatics v0.32 [31]. Reads were mapped to the whole wheat metagenome using Bowtie2 v2.3.5.1 [32], and mapped reads were filtered out from each sample. Then, short Illumina reads from triplicates of each nitrate treatment (30, 70 and 100 ppm) were assembled using SPADES v3.13.0 [33] into longer contigs, to create three wheat root microbiome catalogs for each treatment separately. The 30 ppm nitrate catalog had 677,271 contigs with N50 of 964 bp, 70 ppm nitrate catalog had 644,394 contigs with N50 of 971 bp, and the 100 ppm catalog had 677,271 contigs with N50 of 964 bp. Those three catalogs were combined and Prodigal v2.6.2 [34] was used for protein-coding gene prediction. To create a non-redundant set of genes, we used CD-HIT-EST software v4.8.1 [35] with a similarity threshold of 95%. Those genes were used as the root gene catalog, which included 35 million partial genes. This gene catalog was searched against the non-redundant NCBI protein database using DIAMOND sensitive algorithm v0.9.24.125 [36] to assign taxonomic and functional annotations. Results were then uploaded to MEGAN Ultimate edition software v6.15.2 [37]. The LCA (lowest common ancestor) algorithm was applied (parameters used with minimum bit-score of 70, minimum support of 5% and 30% top threshold) to compute the assignment of genes to specific taxa. For functional annotation, the Kyoto Encyclopedia of Genes and Genomes (KEGG) database [38] was used. Following annotation, to generate taxonomic and functional count tables, each library was mapped to the gene catalog with Trinity mapping software v2.8.4 [39], with Bowtie2-modified parameters (–no-unal –gbar 99999999 -k 250 –dpad 0 –mp 1,1 –np 1 –score-min L,0,−0.9 -L 20 -i S,1,0.50).
Data analyses
Significance of interactions between CO2 and nitrate levels on soil and plant parameters was calculated using two-way ANOVA the least-squares method, in JMP 14 Pro software (SAS Institute Inc., Cary, NC, USA). Differences between soil and plant parameters as influenced by interactions between CO2 and nitrate levels was calculated using Student’s t test in JMP 14 Pro software and statistical significance was set at P < 0.05 The abundance of genes measured using qPCR, was calculated using Student’s t test and statistical significance was set at P < 0.05. All sequencing data analyses were performed using R statistical software. An nMDS (non-metric multidimensional scaling) ordination plot was constructed using R package VEGAN v.2.5-5 [40]. The data matrix was transformed using normalized count transformation, and the batch effect of three experimental cycles was removed using the DESeq2 v.1.22.2 package [41] limma::removeBatchEffect function. For community structure, ordination was generated using the Bray–Curtis dissimilarity matrix, and for community function of KEGG orthologous (KO) groups, ordination was generated using the Euclidean dissimilarity matrix. Permutational multivariate analysis was used to compute the variance between bacterial community structure or function and experimental parameters (CO2, nitrate and cycle), using the Adonis function in the Vegan R package [42]. For comparison of taxonomic changes and functional traits, differential abundance was estimated using DESeq2 [41] and was considered significant when the difference in abundance between genes had a FDR-adjusted P value < 0.05. For comparison of changes in taxonomic and functional traits, the bacterial read counts table was binned into KO groups, based on DIAMOND-MEGAN annotation. The links between gene functions (denitrification and type VI secretion system (T6SS)) and wheat root microbial community structure at the order level, represented as networks, were constructed using Cytoscape v.3.7.2 software [43].
Source: Ecology - nature.com
