Influence of Belowground Herbivory on the Dynamics of Root and Rhizosphere Microbial Communities

Ourry, Morgane; Lebreton, Lionel; Chaminade, Valérie; Guillerm-Erckelboudt, Anne-Yvonne; Hervé, Maxime; Linglin, Juliette; Marnet, Nathalie; Ourry, Alain; Paty, Chrystelle; Poinsot, Denis; Cortesero, Anne-Marie; Mougel, Christophe

doi:10.3389/fevo.2018.00091

ORIGINAL RESEARCH article

Front. Ecol. Evol., 26 June 2018

Sec. Behavioral and Evolutionary Ecology

Volume 6 - 2018 | https://doi.org/10.3389/fevo.2018.00091

This article is part of the Research Topic As Above so Below? Below-ground Interactions in Ecological Processes View all 19 articles

Influence of Belowground Herbivory on the Dynamics of Root and Rhizosphere Microbial Communities

$\r\nMorgane Ourry$ Morgane Ourry¹

Lionel Lebreton¹

Valérie Chaminade²

Anne-Yvonne Guillerm-Erckelboudt¹

Maxime Hervé²

Juliette Linglin¹

Nathalie Marnet^1,3

Alain Ourry⁴

Chrystelle Paty²

Denis Poinsot²

Anne-Marie Cortesero^2†

Christophe Mougel¹^*^†

¹IGEPP, INRA, Agrocampus Ouest, Université de Rennes 1, Le Rheu, France
²IGEPP, INRA, Agrocampus Ouest, Université de Rennes 1, Rennes, France
³BIA, Team “Polyphenols, Reactivity & Processes”, Le Rheu, France
⁴Normandie Université, UNICAEN, INRA, UMR 950 Ecophysiologie Végétale, Agronomie et nutritions N, C, S, Esplanade de la Paix, Caen, France

Recent studies are unraveling the impact of microorganisms from the roots and rhizosphere on interactions between plants and herbivorous insects and are gradually changing our perception of the microorganisms' capacity to affect plant defenses, but the reverse effect has seldom been investigated. Our study aimed at determining how plant herbivory influences the dynamics of root and rhizosphere microbial community assemblages and whether potential changes in root metabolites and chemical elements produced during herbivory can be related to microbial community diversity. We conducted our study on oilseed rape (Brassica napus) and its major belowground herbivore, the cabbage root fly (Delia radicum). We further assessed the influence of initial soil microbial diversity on these interactions. Different microbial diversities based on a common soil matrix were obtained through a removal-recolonization method. Root and rhizosphere sampling targeted different stages of the herbivore development corresponding to different perturbation intensities. Root bacterial communities were more affected by herbivory than some rhizosphere bacterial phyla and fungal communities, which seemed more resistant to this perturbation. Root herbivory enhanced the phylum of γ-Proteobacteria in the roots and rhizosphere, as well as the phylum of Firmicutes in the rhizosphere. Herbivory tended to decrease most root amino acids and sugars, and it increased trehalose, indolyl glucosinolates, and sulfur. Higher abundances of four bacterial genera (Bacillus, Paenibacillus, Pseudomonas, and Stenotrophomonas) were associated following herbivory to the increase of trehalose and some sulfur-containing compounds. Further research would help to identify the biological functions of the microbial genera impacted by plant infestation and their potential implications in plant defense.

Introduction

Interactions between plants and herbivorous insects are known to influence the evolution of plant defense mechanisms such as defensive and toxic metabolites or volatiles attracting natural enemies of the pest (Fürstenberg-Hägg et al., 2013; Nishida, 2014). Plants constitute an important crossroad in biological interactions since they are also involved in interactions with microorganisms that can improve plant growth, nutrient acquisition, and protection from different plant bioagressors (Richardson et al., 2009; Schnitzer et al., 2011). Recent studies are unraveling the impact of microorganisms from the rhizosphere on plant-insect interactions and are gradually changing our perception of the microorganisms' capacity to affect plant defenses (Pineda et al., 2017). A single strain or a whole microbial community has been shown to increase plant defenses, which in turn can negatively impact insect fitness (Hol et al., 2010; Pangesti et al., 2015a,b; Lachaise et al., 2017). However, the reverse effect has seldom been investigated while herbivore attacks on plants could represent important perturbations for belowground microbial communities.

Perturbations are seen as major drivers of ecosystems stability (Milchunas et al., 1988). Since the 90's, the positive relationship between an ecosystem's stability and its diversity was addressed in different systems and particularly in grasslands facing different abiotic stresses (Ives and Carpenter, 2007). It seemed that the recovery of an ecosystem following a stress could be conditioned by the ecosystem's initial diversity prior to this stress. For example, an ecosystem with high plant diversity that suffers a drought will better recover than an ecosystem with a low plant diversity (Tilman and Downing, 1994; van Ruijven and Berendse, 2010). Effects of such stresses were investigated on resistance (the ability to sustain a perturbation or stress) and resilience (the capacity to come back to a stable state) of different soil processes such as nitrogen and carbon cycles, organic matter decomposition, or respiration (Griffiths and Philippot, 2013). Resistance and resilience of grassland communities toward biotic stresses such as herbivory has also been investigated (see review by Cingolani et al., 2005). However, these studies concerned either grazing mammals or insect herbivory in grasslands, with the latter having variable effects on mycorrhiza (Johnson and Rasmann, 2015). To our knowledge however, this question has never been addressed for microbial communities of agroecosystem plants confronted to insect herbivory.

Insect herbivory can modify plant physiology (Bolton, 2009) and alter metabolite concentrations in plant tissues (Griffiths et al., 1994; Hopkins et al., 1995; Ponzio et al., 2017) as well as nutrient uptake (Katayama et al., 2014; Aziz et al., 2016). It can also modify the release of organic matter by plant roots (rhizodeposition), a major driver of soil microbial communities (Singh et al., 2004; Paterson et al., 2007). Some studies showed that root herbivory can influence plant carbon sources and rhizosphere chemistry, which in turn modify the abundance of bacteria and fungi, and community physiological profiles in the rhizosphere (Grayston et al., 2001; Dawson et al., 2004; Treonis et al., 2005). Katayama et al. (2014) showed that herbivory could also alter chemical element uptakes by the roots as well as nitrogen concentration of microbial origin. So far, when the effect of herbivory on plant-associated microorganisms was covered, the experiments either neglected to take into account the plant metabolites and chemical elements or did not integrate the notion of perturbation dynamics, or cultivable methods and fluorimetry were used to study microbial communities.

Our study aimed at determining the impact of insect herbivory on the plant-microorganism interactions. We hypothesized that (i) the dynamics of root and rhizosphere microbial communities would be influenced by belowground herbivory, (ii) these dynamics would depend on root metabolites and chemical elements induced by herbivory, (iii) the initial soil microbial diversity would influence plant chemistry and hence the dynamics of microbial communities. We conducted our experiment on oilseed rape (“OSR,” Brassica napus) and its belowground pest, the cabbage root fly (“CRF,” Delia radicum) (Ahuja et al., 2010). This fly is a specialist of brassicaceous plants. Females lay their eggs aboveground, at the base of plant stems and, upon hatching, which occurs a few days later, larvae feed by tunneling into the roots for 2–3 weeks, before pupating in the nearby soil and emerging. Herbivory of the CRF on brassicaceous plants is known to change the concentration of primary (e.g., sugars such as glucose and sucrose) and secondary metabolites (e.g., indolyl glucosinolates) in the roots (Hopkins et al., 1999; van Dam and Raaijmakers, 2006; Pierre et al., 2012; van Geem et al., 2015). However, nothing is known about how CRF herbivory and resulting biochemical changes in the roots influence the belowground microbial communities which interact with the plant.

Materials and Methods

Soil Preparation and Inoculation

A batch of soil was collected in November 2014 from the layer −10 to −30 cm deep of a field in Brittany (La Gruche, Le Rheu, France, 48°08′44″N, 01°47′97″ W) where wheat was cultivated for 20 years, and it was stored in containers at ambient temperature in the dark. After a year of storage in containers, the soil was homogenized, ground, sieved at 4 mm to remove the macrofauna and mixed with 1/3 silica. As described in Lachaise et al. (2017), this mixture was sterilized at 35 kGy and left 2 months to stabilize before inoculation while the remaining unsterilized soil with no silica was ground and sieved at 2 mm before being suspended in water. Following the detailed protocol of Lachaise et al. (2017), this suspension was then undiluted (10⁰) or diluted at 10⁶ before inoculating the sterilized soil, hence creating the two levels of soil microbial diversity used in our experiment: respectively “high” and “low” initial soil microbial diversities, also referred to as soil microbial modalities (Figure 1). This dilution-inoculation method was performed three times in order to obtain three soil biological replicates per soil microbial modality. After inoculation, the soil was incubated in the dark for seven weeks at 18°C and 50% humidity. During this period, the bags containing the soil were opened weekly under sterile conditions using a laminar flow cabinet to homogenize the soil and facilitate microbial respiration and recolonization. This allowed the soil to reach optimal bacterial and fungal densities and similar abundances of Colony Forming Units between the two soil microbial modalities at the end of the recolonization period. A part of this soil was collected before sowing to evaluate initial soil microbial diversities (N = 18, 9 samples per soil microbial modality).

FIGURE 1

Figure 1. Experimental design. Collected soil was suspended in water and diluted or not, to obtain two solutions of respectively low and high microbial diversities. Sterilized soil was inoculated with one of these two solutions and incubated for 49 days. After sowing, plants were cultivated for 42 days and then two treatments were applied, obtaining healthy and infested plants. The experiment started after the infestation (i.e., D. radicum eggs deposited on the plant crown) and sampling was carried out according to the fly biological cycle and important development stages: egg hatching, 3rd instar larva, end of adult emergence, occurring around 1 day after infestation (DAI), 14 and 42 DAI respectively. These sampling times also corresponded to the beginning, the peak and the end of herbivory. Metagenomics, metabolomics and elemental analyses were performed on the obtained samples.

Insect Rearing

The population of cabbage root fly (“CRF,” Delia radicum) used in our experiment was collected in the field in 2015 (Le Rheu, Bretagne, France). In the laboratory, flies were fed on sugar, milk powder and yeast (ratio 1:1:1) and they were reared on rutabaga roots (Brassica napus subsp. rapifera) in a climatic room (16:8 LD, 21 ± 2°C; 60% ± 10% RH) as described in Neveu Bernard-Griffiths (1998). Adult flies were left to oviposit on rutabaga roots for 24h and black-headed eggs (i.e., ready to hatch) were collected 3 days later for our experiment.

Plant Growth, Infestation, and Sampling

Seeds of Brassica napus L. (subsp. oleifera cv. Tenor) were sown in individual pots using a layer of pozzolan at the bottom and the soil previously obtained (characterized either with a “high” or “low” initial soil microbial diversity). Plants were watered twice a week by sub-irrigation with a nutritive solution based on Hoagland and Arnon (1950) during the whole experiment. This solution was obtained by blending three separate solutions: a macronutrient solution (3 mM of KNO₃, 0.5 mM of KH₂PO₄, 1 mM of MgSO₄ 7H₂O, 2.5 mM of Ca(NO₃)₂ 4H₂O), a micronutrient solution (10 μM of MnSO₄ H₂O, 1 μM of ZnSO₄ 7H₂O, 0.5 μM of CuSO₄, 30 μM of H₃BO₃, 1 μM of Na₂MoO₄ 2H₂O, 0.5 μM of Co(NO₃)₂ 6H₂O), and a Fe EDTA solution (27 μM of Fe EDTA Na). Plants were cultivated during 6 weeks in a greenhouse under natural late 2016 winter photoperiod (mean temp. 15.3°C, min/max 5.5/25.6°C). Two treatments were applied on the plants: half of the plants were infested by depositing eight black-headed eggs on the crown per plant, the other half remained untreated. These two batches of plants were then referred to as “infested plants” and “healthy plants” respectively. All plants were then moved to a climatic chamber (photoperiod 16:8 LD and thermoperiod 20:18°C LD) for the rest of the experiment. Afterwards, healthy and infested plants were sampled at 1, 14, and 42 days after infestation (“DAI”), which respectively represented (i) the hatching stage which corresponds to the initiation of herbivory, (ii) the third larval instar which corresponds to the peak of herbivory, and (iii) the end of the fly emergence which corresponds to the end of herbivory (Figure 1). A total of 108 plants were harvested, with 9 plants per condition (2 soil microbial modalities, 2 treatments being healthy and infested, and 3 sampling times). Roots and rhizosphere were sampled as follows: the root fraction was collected, corresponding to the area from the crown to the tip roots, and was washed twice in 20 mL of sterile permuted water before being transferred to a clean Falcon tube; the whole root bath (i.e., 40 mL in a Falcon tube) corresponded to the rhizosphere fraction; both tubes were immersed in liquid nitrogen and stored at −80°C before being freeze-dried; only roots were ground using glass beads. In order to have sufficient material for molecular, metabolomic and elemental analyses, roots from three different plants of the same treatment, cultivated on the same soil microbial diversity and biological replicates, were pooled to make one sample, hence obtaining a total of 36 root samples (N = 3 samples per condition). The rhizosphere fraction (N = 36 rhizosphere samples) was treated similarly.

Molecular and Bioinformatic Analysis

Initial bulk soil (i.e., before sowing), as well as root and rhizosphere samples collected during the experiment, were analyzed.

According to the protocol developed by the GenoSol platform (Dijon, France) and described in Plassart et al. (2012), soil DNA was extracted from 2 g of wet bulk soil (i.e., initial soil before sowing) or 1 g of freeze-dried rhizosphere respectively in 8 or 5 mL of lysis buffer containing 100 mM of Tris-HCl (pH 8), 100 mM of EDTA (pH 8), 100 mM of NaCl, 2% SDS, and ultrapure water. The following modifications were performed: tubes were vortexed at mid and at the end of incubation (i.e., bath at 70°C), then centrifuged at 3,500 rpm at 15°C for 10 min; during deproteinization, samples were centrifuged at 14,000 g at 4°C during 10 min; during DNA precipitation, tubes were placed at −20°C for 30 min. To obtain DNA pellets, tubes were centrifuged at 13,000 rpm at 4°C for 30 min and supernatant was discarded. DNA pellets were washed as follows: 400 μL of 70% ice-cold ethanol were added to samples, which were centrifuged at 13,000 rpm at 4°C for 5 min and supernatants were removed. Remaining traces of ethanol were removed by heating open tubes at 60°C for at least 15 min or more if needed. Pellets of DNA were resuspended with 100 μL of ultrapure water and the duplicated samples were finally pooled. Bulk soil and rhizosphere samples were purified twice. The first purification required Microbiospin (Biorad, Hercules, California, USA) columns of PVPP (PolyVinyl PolyPirrolidone, Sigma-Aldrich). To prepare the columns, their tips were removed and columns were placed in 2 mL Eppendorf tubes. Then, columns were washed with 400 μL of ultrapure water and centrifuged at 1,000 g and at 10°C during 2 min. After emptying the tubes, the previous step was carried out a second time, before centrifuging empty tubes at 1,000 g at 10°C for 4 min. Following Plassart et al. (2012), 100 μL of DNA were injected into the columns, previously transferred to a clean tube, however our samples were incubated on ice for 5 min, before a 4 min centrifugation at 1,000 g at 4°C. The obtained DNA (~95 μL) was used for the second purification, performed using the Geneclean^® Turbo kit (MP Biomedicals) with the following modifications: samples were centrifuged at 10,000 g for 10 s; after adding the GTE (GeneClean Turbo Elution Solution), samples were incubated on ice for 5 min and centrifuged at 10,000 g for 1 min; the GTE, incubation and centrifugation steps were repeated a second time, to finally obtain ~60 μL of clean DNA.

Root DNA was extracted using the NucleoSpin^® Plant II kit and protocol (Macherey-Nagel, Düren, Germany) with the following modifications: 30 g of freeze-dried root powder were used; (i) cell lysis was done with buffer PL1; (ii) incubation after adding RNase A lasted 30 min; (iii) the crude lysate was centrifuged before its filtration.

DNA quantification was performed using a Quantus fluorometer (Promega, Madison, Wisconsin, USA) and the Quantifluor kit (dsDNA: E2670).

PCR amplification and sequencing were performed at the GenoScreen platform (Lille, France) using the Illumina Miseq platform to a 2 × 250 bases paired-end version. We used PCR primer pairs 799F (5′-AACMGGATTAGATACCCKG-3′) and 1223R (5′-CCATTGTAGTACGTGTGTA-3′), and NS22B (5′-AATTAAGCAGACAAATCACT-3′) and SSU0817 (5′-TTAGCATGGAATAATRRAATAGGA-3′) (Lê Van et al., 2017) to amplify 16S and 18S rDNA genes, respectively.

To manage mismatch between reads 1 and 2 in the overlap region, we used bases trimming at Q30 with PRINSeq, trimming of specific primers with Cutadapt, assembly with FLASH starting with trimmed reads 1 and 2, a minimum of 30 bases overlapping was required with a 97% homology between reads. Regarding 18S rDNA, only read 1 was analyzed because i) amplicon size did not allow a sufficient overlapping area between read 1 and read 2 and ii) read 1 was of better quality. The GnS-PIPE bioinformatical pipeline was used for the bioinformatical analyses of 16S rDNA and 18S rDNA (Terrat et al., 2012). Raw data sets were deposited on the European Nucleotide Archive database system under the project number (PRJEB25217). Bulk soil samples accession numbers range from ERS2281263 to ERS2281298, those of root and rhizosphere samples range from ERS2255945 to ERS2256016 for 16S rDNA and from ERS2256770 to ERS2256841 for 18S rDNA.

Metabolites and Chemical Elements Analysis

Primary Metabolites: Amino Acids, Carbohydrates, Polyols, and Organic Acids

Quantification of free amino acids (AAs), non-structural carbohydrates, polyols, and organic acids (CPOs) was based on the method described by Gravot et al. (2010) and performed using 9 to 12 mg of freeze-dried root powder, with the same methanol-chloroform-water-based extraction. Minor adjustments were made to this protocol: after adding 100% chloroform, tubes were rapidly vortexed and then agitated for 10 min at room temperature; after adding water, samples were vortexed for 20 s and centrifuged for 5 min at 12,000 g and 15°C.

For AA derivatization and profiling, 50 μL of methanol-water extract were vacuum-dried. The dry residue was resuspended in 50 μL of ultrapure water and the tubes were rapidly vortexed, put in an ultrasonic bath for 5 min and centrifuged for 5 s at room temperature. Derivatization of AAs was performed using the AccQTag Ultra Derivatization kit (Waters) with the following modified volumes: 5 μL of the resuspended sample, 35 μL of AccQTag Ultra Borate Buffer and 10 μL of AccQTag Reagent were placed in a new tube, which was vortexed and placed in a water bath for 10 min at 55°C. The whole volume was transferred in a vial and derivatizated AAs were analyzed using liquid chromatography (Acquity UPLC-DAD system, Waters, Milford, MA, USA) according to Jubault et al. (2008). However, the column used for analyses was heated at 53°C and AAs were detected at 265 nm using a photodiode array detector. Identification of AAs was realized by comparison with a standard solution and quantification was made thanks to the internal standard BABA (3-aminobutyric acid).

Analysis of CPOs was performed by gas chromatography coupled to a flame ionization detector (GC-FID, Thermo Fisher Scientific, Waltham, MA, USA) and based on the protocol described by Lugan et al. (2009), which however required several modifications. The online derivatization for CPOs was performed with a Trace 1300 GC-FID (Thermo Scientific) equipped with a Tri Plus RSH (Thermo Scientific). Fifty microliters of the methanol-water extract were sampled in injection vials and vacuum-dried. This online derivatization was performed as follows: the dried extract was resuspended in 50 μL of pyridine containing 20 mg.mL⁻¹ of methoxyaminehydrochloride, under orbital shaking at 40°C for 90 min. Fifty microliters of MSTFA (N-methyl-trimethylsilyltrifluoroacetamide) were added before incubation at 40°C for 30 min. One microliter of the mixture was injected into the GC-FID with a split/splitless injector (split mode set to 1:20) at 260°C, on a TG-5MS column (30 m × 0.32 × 0.25 mm, Agilent Technologies) connected to a flame ionization detector at 300°C. The temperature gradient of the GC oven was: 4 min at 100°C followed by an increase of 10°C.min⁻¹ up to 198°C and maintained at this temperature for 2 min; an increase of 1°C.min⁻¹ up to 202°C; then an increase of 15°C.min⁻¹ ramp up to 268°C and held for 3 min followed by an increase of 1°C.min⁻¹ up to 272°C and raised to 210°C at 10°C.min⁻¹ maintained for 7 min. Identification of CPOs was realized by comparison with a standard solution, while quantification was achieved with adonitol as the internal standard.

Secondary Metabolites: Glucosinolates

Extraction and analysis of glucosinolates (GSLs) were performed based on the protocol from Hervé et al. (2014) with the following modifications: GSLs were extracted from 12 to 15 mg of freeze-dried root powder and tubes were centrifuged for sedimentation at 12,000 g and 15°C for 5 min. Analysis of GSLs was performed using liquid chromatography coupled with mass spectrometry (Acquity UPLC-TQD, Waters) with electrospray ionization in a negative mode. Chromatographic conditions, A and B-eluents and the gradient used are described in Hervé et al. (2014). Quantification of GSLs was realized using an external calibration with a standard solution containing glucoerucin, gluconasturtiin, and glucobrassicin in the range of 4 to 80 μmol.L⁻¹. These compounds were respectively used to quantify aliphatic, aromatic and indolyl GSLs.

Chemical Elements

Analysis of root chemical elements composition was performed according to Maillard et al. (2015). About 4 mg of freeze-dried root powder were used to measure total N and S contents, using a continuous flow isotope mass spectrometer (Nu Instrument, Wrexham, United Kingdom) linked to a C/N/S analyzer (EA3000, Euro Vector, Milan, Italy). For other elements such as K, Ca, P, Mg, Fe, Mn, Zn, Cu, Mo, Ni, and B, roots samples were submitted to microwave acid sample digestion (Multiwave ECO, Anton Paar, les Ulis, France) using 800 μL of concentrated HNO₃, 200 μL of H₂O₂ and 1 mL of Milli-Q water for 40 mg of root powder. All samples were previously spiked with two internal-standard solutions of Gallium and Rhodium, respectively, for a final concentration of 10 and 2 μg.L⁻¹. Mineralized samples were then diluted to 50 mL with Milli-Q water to obtain solutions containing 2.0% (v/v) of nitric acid, then filtered at 0.45 μm using a teflon filtration system (Filtermate, Courtage Analyses Services, Mont-Saint-Aignan, France). Samples were then analyzed by High Resolution Inductively Coupled Plasma Mass Spectrometry (HR ICP-MS, Thermo Scientific, Element 2^TM) and quantification of each element was performed using external standard calibration curves. Additional information about mineralization and the utilization of HR ICP-MS can be found in the annexes of Maillard (2016).

Statistical Analyses

Analyses were performed using the R software (R Core Team, 2016) and a 5% threshold for statistical significance.

Microbial Communities

When analyses were performed on microbial data, some samples had to be removed from the dataset for the following reasons: one root sample was lost before sequencing (from the low diversity-healthy plants-14 DAI condition), one rhizosphere bacterial sample (from the high diversity-healthy plants-1 DAI condition) showed an abnormal smaller total read count after normalization while one root fungal sample (from the high diversity-healthy plants-1 DAI condition) had abnormal phyla abundances compared to other samples.

Bacterial and fungal richnesses and diversities, represented by the number of observed Operational Taxonomic Units (OTUs) and the Shannon index (obtained with the diversity function in the “vegan” package) respectively, were obtained from non-rarefied OTU data and analyzed in roots and rhizosphere separately using a linear model. Models took into account the following factors: sampling time (1, 14, and 42 DAI), plant treatment (healthy vs. infested plants), soil microbial diversity (high vs. low diversity), and soil replicate, as well as paired interactions between the three first factors. Linear models were adjusted depending on the significance of interactions: interactions were all removed if none was significant while they were kept if at least one was significant or close, which was respectively the case of analyses on the rhizosphere and root compartments. A type II analysis of variance table was performed on the models, followed by comparisons based on least-squares means when possible.

In order to analyze the bacterial and fungal community structure, distance-based redundancy discriminant analysis (dbRDA, dbrda function in the “vegan” package) was performed on a Bray-Curtis dissimilarity matrix, obtained from OTU data, which were filtered using a 1%0 threshold and log2-transformed. A type II permutation test for constrained multivariate analyses was performed on the dbRDA using the “RVAideMemoire” package (Hervé, 2016a,b) to evaluate the contribution of each factor (i.e., compartment, time, treatment, soil microbial diversity, and soil biological replicates) to microbial community structure. Quantified contributions expressed as r-squared, were obtained from the varpart function (“vegan” package).

Phyla analyses were performed according to the scripts provided by Bulgarelli et al. (2015) but P-values were corrected using the False Discovery Rate (FDR) method. Plotting was done with the “ggplot2” package.

Genewise negative binomial generalized linear models were conducted using the “edgeR” package (Anders et al., 2013) on filtered data to recover genera differences between plant treatments within sampling time. The obtained P-values were corrected with the FDR method. Using “ggplot2” package, genus count differences between infested and healthy (Xi – Xh) were plotted when the treatment was significant on one hand and when these differences exceeded 250 or −250 counts.

Metabolites and Chemical Elements

For each sampling time, two different analyses were performed on metabolites and chemical elements, taking into account the plant treatment (healthy vs. infested plants), the soil microbial diversity (high vs. low diversity), and the interaction between both factors, as well as soil biological replicates. Both analyses required the data to be transformed using the fourth root. First, a redundancy analysis (RDA), associated to a permutation significance test based on cross-validation, was performed to determine the influence of the previously mentioned factors on the root metabolomic and elemental profiles. Second, linear models were used to assess precisely which metabolites and elements were affected by the above factors. Finally, a type II analysis of variance was performed on the linear models and the obtained P-values were corrected with the FDR method. The content differences between infested and healthy (Xi – Xh) was plotted.

Relationships Between Root Microbial Communities and Root Chemistry

To assess the relationship between root microbial communities and root chemistry under herbivory, we applied the DIABLO (Data Integration Analysis for Biomarker discovery using Latent cOmponents) method on our data with the plant treatment (i.e., infested vs. healthy plants) defined as the grouping factor, using the scripts provided by Hervé et al. (2018) and two R packages (“mixOmics and “RVAideMemoire”). Only metabolites or chemical elements significantly impacted by the treatment at a given time were kept for this analysis at this same given time. Bacterial and fungal genera used in this analysis were also significantly impacted by the plant treatment and bacterial count differences respected the 250/−250 threshold mentioned above. Hence, a three block-DIABLO was performed at 14 DAI on one hand, using three datasets (i.e., absolute abundance (expressed in counts of OTUs) of microbial genera, metabolites, and chemical elements) while a two block-DIABLO was performed at 42 DAI on the other hand with two datasets: microbial genera and combined metabolites and elements. In the latter case, an analysis on three blocks was not possible using a single column dataset for chemical elements, which was why it was included in the metabolite dataset. Discrimination of the grouping factor was assessed with a permutation significance test based on cross-validation.

Results

Our experiment was based on two soils of different bacterial and fungal alpha-diversities (Table S1) and phyla abundances (Figure S1).

Diversity and Structure of Microbial Communities

After 6 weeks of OSR culture in soils of different microbial diversities and 1 DAI, the bacterial and fungal alpha-diversities were similar in both root and rhizosphere compartments.

Bacteria

Root bacterial richness and diversity were significantly impacted by sampling time, initial soil microbial diversity and the interaction between time and treatment (Table 1). Diversity was also influenced by the interaction between time and initial soil microbial diversity.

TABLE 1

Table 1. Statistical output of bacterial and fungal alpha-diversity indices in the root and rhizosphere compartments, associated to Figure 2.

At high initial soil microbial diversity, the bacterial richness (observed OTUs) of healthy plants increased with time, reaching a peak of 900 observed OTUs at 14 DAI, and then decreased at 42 DAI (Figure 2A). Herbivory reversed this trend and richness was significantly increased at 42 DAI. Richness was significantly different between healthy and infested plants both at 14 and 42 DAI. Root bacterial diversity (Shannon index) had a similar profile at high and low initial diversities whether or not plants were infested. The Shannon index was significantly lower in infested plants than in healthy plants at 14 DAI but it was higher at 42 DAI (Figure 2A).

FIGURE 2

Figure 2. Bacterial and fungal alpha-diversity indices in the root and rhizosphere compartments. Mean value (± se) of bacterial (A) and fungal (B) richness (i.e., Observed OTUs) and diversity (i.e., Shannon index) for the root (left column) and rhizosphere (right column) compartments at 1, 14, and 42 days after infestation (DAI). In the root compartment, lowercase letters represent significant differences between modalities at one given time while asterisks represent significant differences between times for one given modality. Regarding the rhizosphere compartment, capital letters represent significant differences between times for all modalities.

At low initial soil microbial diversity, the bacterial richness (observed OTUs) of healthy plants did not change with time but increased significantly after herbivory (i.e., 42 DAI) (Figure 2A). At this low initial soil diversity, root bacterial diversity (Shannon index) increased in healthy plants at 14 DAI and remained stable at 42 DAI while infested plants displayed an opposite profile (Shannon index decreased at 14 DAI and increased at 42 DAI).

Rhizosphere bacterial richness and diversity were both influenced by sampling time and initial soil microbial diversity (Table 1, Figure 2A). Richness and Shannon index decreased over time in both healthy and infested plants, the latter being similar between 14 and 42 DAI. These indices showed higher values in plants grown on high soil microbial diversity than on low soil diversity.

A dbRDA confirmed that bacterial community structure was driven by compartment, time, treatment and initial soil microbial diversity (Figure 3A). Differences in bacterial community structures between healthy and infested plants were highly significant.

FIGURE 3

Figure 3. Driving factors of bacterial (A) and fungal (B) community structure. Presented for each microbial kingdom is an ordination plot and the output of the type II permutation test performed on the dbRDA, using the Bray-Curtis dissimilarity index and OTU data set. It includes the variation (adjusted r-squared) explained by each factor, the F-value, the degrees of freedom (df), and the P-value (in bold when significant). Compartment refers to roots and rhizosphere, time refers to 1, 14, and 42 days after infestation (DAI), treatment refers to healthy and infested plants and soil microbial diversity refers to high and low soil microbial diversities.

Fungi

Root fungal richness was four times lower than root bacterial richness and ranged from 150 to 300 observed OTUs. It was only influenced by sampling time (Table 1). The richness of healthy and infested plants at both soil diversities similarly increased over time and reached 250 OTUs at 42 DAI (Figure 2B).

Root fungal diversity (Shannon index) was impacted by time and by the interaction between time and treatment but not by soil microbial diversity (Table 1). Both richness of healthy and infested plants increased over time but the Shannon index was higher in infested plants at 14 DAI, hence during the peak of herbivory (Figure 2B).

Rhizosphere fungal richness was similar to root richness, with 250 to 350 observed OTUs, and it did not vary as much as bacterial richness. Rhizosphere fungal richness was impacted by initial soil microbial diversity and by the soil biological replicates (Table 1). Richnesses of healthy and infested plants remained similarly stable over time (Figure 2B). Richness was greater in plants grown on high soil microbial diversity.

Rhizosphere fungal diversity was influenced by both sampling time and initial soil microbial diversity (Table 1). The Shannon index of healthy and infested plants increased similarly over time and diversity was greater in plants grown on high soil microbial diversity (Figure 2B).

Fungal community structure drivers were compartment, time, and treatment, which were highly significant, as well as initial soil microbial diversity (dbRDA, Figure 3B).

Identification of Microbial Phyla and Genera Associated With Herbivory and Initial Soil Microbial Diversity

Relative Abundance of Microbial Phyla

The roots (Figure 4A) and rhizosphere (Figure 4B) contained four major bacterial phyla, of which Bacteroidetes, Firmicutes, and Proteobacteria were the most abundant. The root bacterial phyla were mainly influenced by infestation but not at 1 DAI. Significant differences occurred during the peak of herbivory at 14 DAI: the phylum of γ-Proteobacteria increased under herbivory, while the phyla of α-, δ-Proteobacteria, Actinobacteria, Bacteroidetes, and “Others” decreased. At 42 DAI, infestation still influenced the phyla of γ-Proteobacteria, Actinobacteria and “Others.” Initial soil microbial diversity impacted β-Proteobacteria at 14 and 42 DAI with a larger abundance in high diversity-infested plants while on the contrary Actinobacteria was more abundant in low diversity-infested plants at 42 DAI.

FIGURE 4

Figure 4. Dominant bacterial and fungal phyla detected in the roots and rhizosphere. Mean relative abundance are presented for bacterial (A,B) and fungal (C,D) phyla from the root (A,C), and rhizosphere (B,D) compartments. All modalities are represented, taking into account the sampling time (1, 14, and 42 days after infestation), the treatment (healthy and infested plants), and initial soil microbial diversity (high and low levels). Black asterisks indicate significant differences between healthy and infested plants at a given time and soil microbial diversity, while white asterisks indicate significant differences between soil microbial diversities at a given time and for a given treatment.

The bacterial phyla of the rhizosphere (Figure 4B) were also influenced by infestation as well as the initial soil microbial diversity and their abundances were in the same range as in the roots. Overall, the phyla of Actinobacteria, β-, δ-, γ-Proteobacteria, and “Others” differed between high and low microbial diversities at the different times. Infestation started to impact phyla at 14 DAI, with more γ-Proteobacteria and less α-Proteobacteria in infested plants. At 42 DAI, the phyla of α-Proteobacteria and Firmicutes were still affected by infestation and more abundant in infested plants.

In the roots and the rhizosphere, there were four major fungal phyla among which Ascomycota (divided into three sub-phyla), Basidiomycota and Chytridiomycota were the most abundant. Relative abundances of root fungal phyla (Figure 4C) were very variable. At 1 DAI, there was no difference between modalities. Pezizomycotina and Basidiomycota abundances increased in infested plants at 14 and 42 DAI respectively while Taphrinomycotina abundance decreased at 42 DAI. The phyla of Basidiomycota and Chytridiomycota were influenced by the soil microbial diversity at 42 DAI in infested plants.

Fungal phyla (Figure 4D) seemed to vary less in the rhizosphere than in the roots but proportions of Chytridiomycota seemed to decrease in favor of Basidiomycota, Pezizomycotina, and Saccharomycotina. There was no effect of herbivory on fungal phyla at 1 and 14 DAI, but Chytridiomycota abundance decreased in infested plants at 42 DAI. The phyla of Blastocladiomycota and “Others” were lightly influenced by the initial soil microbial diversity, with larger abundances in plants of high diversity.

Absolute Abundance of Microbial Genera

A total of 2,031 bacterial genera were detected in the roots and in the rhizosphere.

The root compartment presented 74 bacterial genera from nine different phyla, which varied significantly following herbivory. Most of these genera showed a different abundance between healthy and infested plants at 14 and/or 42 DAI but not at 1 DAI (Table S2). Out of these 74 genera, Bacillus, Clostridium, Paenibacillus, Pseudomonas, and Stenotrophomonas were the most abundant genera influenced by herbivory (Figure 5A). Bacillus decreased in infested plants while the other four genera increased and none of them returned to an abundance similar to the one detected in healthy plants at 42 DAI. Abundance varied with the initial soil microbial diversity. At low soil diversity, the increase of Paenibacillus abundance following herbivory was greater than at high soil diversity while that of the other genera were similar between both diversities.

FIGURE 5

Figure 5. Absolute abundance variation of root and rhizosphere bacterial genera between healthy and infested plants. Mean differential absolute abundance (± se) are presented for bacterial genera from the root (A) and rhizosphere (B) compartment, associated to either high or low soil microbial diversity, at 14 and 42 days after infestation (DAI). Genera represented in the figure were significantly impacted by the plant treatment (i.e., healthy vs. infested) at a given time and exceeded the 250/–250 threshold. Genus abundance is expressed as the difference between absolute abundance of infested plants (Xi) and absolute abundance of healthy plants (Xh).

In the rhizosphere compartment, fewer bacterial genera (46 genera from 9 phyla) were impacted by the infestation at 14 and/or 42 DAI (Table S3). Compared to the roots, the range of bacterial absolute abundance in the rhizosphere was lower. This time, only Bacillus, Clostridium, Paenibacillus, and Pseudomonas were the most abundant genera influenced and increased by herbivory (Figure 5B). As previously shown in the roots, the rhizosphere was also characterized by a change of range in the differential abundance associated to the initial soil microbial diversity, where the range was lower at low diversity.

A total of 2,593 fungal genera were detected in the roots and in the rhizosphere. The fungal genera varied as much as their phyla (Table S4). In the roots, four out of the five herbivory-influenced genera were more abundant in infested plants, including Ajellomyces and Filobasidiella from the subphylum of Pezizomycotina and the phylum of Basidiomycota respectively. No herbivory-influenced genus was significantly impacted at low microbial diversity in the roots. In the rhizosphere, eight genera were influenced by infestation, Torulaspora and Tuber being more abundant in infested plants at 1 and 42 DAI respectively (Table S4).

Root Compounds Modulated by Herbivory

Primary and Secondary Metabolites

Metabolomic profile of OSR roots was significantly influenced by the treatment at 1, 14, and 42 DAI (Table 2) where infested and healthy plants showed two different profiles, associated to metabolite variations (Figure S2). Their profiles were also impacted by the initial soil microbial diversity but only at 14 DAI. The constrained variance (i.e., part of the variance explained by our variables) had the highest value (73%) at 14 DAI.

TABLE 2

Table 2. Statistical output of multivariate analysis performed on root metabolites and elements.

Overall, 33 metabolites (17 AAs, 9 CPOs, and 7 GSLs presented in Figures 6A–C respectively) were significantly affected by infestation at one or several moments of the plant-insect interaction and 2 metabolites were impacted by the soil microbial diversity (Table S5). Many of them differed at 1 and/or 14 DAI and most of these 33 metabolites were less abundant in infested plants except 2 AAs (β-alanine, SMCSO), 2 CPOs (trehalose, glycerate), and 2 GSLs (glucobrassicin, neoglucobrassicin), which production increased. Twenty-three of these metabolites reached back a stable state, close to the profile of healthy plants at 42 DAI. All GSLs varied due to infestation but six were still influenced at 42 DAI.

FIGURE 6

Figure 6. Variation of root metabolite contents between healthy and infested plants. Mean differential content (± se) of the root amino acids (A), carbohydrates, polyols and organic acids (B), and glucosinolates (C) at 1, 14, and 42 days after infestation (DAI) that were significantly impacted by the plant treatment (i.e., healthy vs. infested), are represented in this figure. Metabolite content is expressed as the difference between content of infested plants (Xi) and content of healthy plants (Xh). Asterisks show significant differences (P < 0.05) between healthy and infested plants at a given time. High and low initial soil microbial diversities are represented in black and gray bars respectively. Metabolites are sorted out in the same order as in Table S5, with only the significant metabolites remaining.

Chemical Elements

Elemental profiles of OSR roots were only impacted by the treatment at 14 and 42 DAI (Table 2, Figure S3). At 14 DAI, the profile of infested plants was highly modified compared to healthy plants. As for metabolites, the constrained variance was the highest at 14 DAI (61%).

Only 8 chemical elements (4 macroelements, 3 microelements, and 1 heavy metal) were significantly impacted by the infestation, while soil microbial diversity had no effect on chemical element content (Figure 7, Table S6). Six elements (Mg, K, Fe, Na, V, Cd) decreased in infested plants at 14 DAI while N content increased. All of these elements reached a state similar to healthy plants at 42 DAI. Sulfur also increased in infested plants, but only at 42 DAI.

FIGURE 7

Figure 7. Variation of root chemical element content between healthy and infested plants. Mean differential content (± se) of the macro- (Mg, N, K, S), microelements (Fe, Na, V), and heavy metal (Cd) at 1, 14, and 42 days after infestation (DAI) that were significantly impacted by the plant treatment (i.e., healthy vs. infested), are represented in this figure. Chemical element content is expressed as the difference between content of infested plants (Xi) and content of healthy plants (Xh). Asterisks show significant differences (P < 0.05) between healthy and infested plants at a given time. High and low initial soil microbial diversities are represented in black and gray bars respectively. Chemical elements are sorted out in the same order as in Table S6, with only the significant chemical elements remaining.

Relationships Between Root Compounds and Microbial Communities

At the peak of herbivory (i.e., 14 DAI), the DIABLO method confirmed the significant differences between infested and healthy plants (CER = 0, P = 0.002) with the three score plots showing similar patterns (Figure 8). The axis 1 of microbiome data was positively correlated to axes 1 of metabolite and element data (r = 0.86 and 0.95, respectively) while metabolite data and element data were also positively correlated (r = 0.84). Infested plants seemed to be associated to (i) four bacterial genera (Clostridium, Paenibacillus, Pseudomonas, and Stenotrophomonas) and one fungal genus (Torulaspora), that were overexpressed, (ii) trehalose (CPO) and neoglucobrassicin (GSL), also overexpressed, (iii) nitrogen but to a lesser extent. The remaining AAs, CPOs, GSLs and seven elements were underexpressed in infested plants.

FIGURE 8

Figure 8. DIABLO graphical outputs at 14 days after infestation (DAI). Score plots of microbiome, metabolite and chemical element datasets are positioned on the left side while the correlation circles associated to each of these score plots are positioned on the right side of the figure. Correlations between axes 1 (axes 2) of microbiome and metabolite datasets: 0.86 (0.86); of microbiome and chemical element datasets: 0.95 (0.60); of metabolite and element datasets: 0.84 (0.77).

At 42 DAI, the discrimination between infested and healthy plants was also significant (CER = 0, P = 0.002) and axis 1 of microbiome data was highly positively correlated to axis 1 of combined metabolite and element data (r = 0.95), confirmed by similar patterns on the two score plots (Figure 9). At the end of herbivory, infested plants seemed to be associated to (i) four bacterial genera (Clostridium, Paenibacillus, Pseudomonas, Stenotrophomonas) and three fungal genera (Ajellomyces, Filobasidiella, Torulaspora), which were overexpressed and (ii) SMCSO (AA), glycerate (CPO), glucobrassicin and neoglucobrassicin (GSLs), and sulfur. Conversely, one bacterial and one fungal genus, as well as four metabolites were underexpressed in these infested plants.

FIGURE 9

Figure 9. DIABLO graphical outputs at 42 days after infestation (DAI). Score plots of microbiome dataset and combined metabolite and chemical element datasets are positioned on the left side while the correlation circles associated to each of these score plots are positioned on the right side of the figure. Correlations between axes 1 (axes 2) of microbiome and metabolite/chemical element datasets: 0.95 (0.90).

Discussion

Our study showed that infestation of OSR by a belowground herbivore increased root bacterial alpha-diversity, while it had no effect on rhizosphere bacteria and very little on fungi. Interestingly, herbivory was associated with an increase of γ-Proteobacteria and Firmicutes (and five of their most dominant affiliated genera: Pseudomonas, Stenotrophomonas and Bacillus, Clostridium, Paenibacillus respectively) as well as Ajellomyces and Filobasidiella from the fungal phyla of Ascomycota and Basidiomycota. Trehalose, sulfur-containing metabolites and sulfur contents increased in infested plants and could explain the variations observed in bacterial phyla and genera. Herbivory seemed to have only a short-term negative effect on the richness and diversity of the microbial communities, which were both restored when the perturbation ended. The chemical composition of roots matched this restoration. Initial soil microbial diversity itself had little impact on microbial communities of the root and its rhizosphere and did not significantly influence root chemistry.

Herbivory Influences Microbial Community Diversity and Composition

Our study showed that root herbivory decreased bacterial richness and diversity in roots at the peak of herbivory but globally increased both at the end of herbivory. Opposite trends were found in grassland plant communities, where plant diversity decreased after the end of herbivory (e.g., Collins et al., 1998; Pucheta et al., 1998) but microbial communities may be hard to compare with plant communities. Further from the roots, bacterial and fungal communities of the rhizosphere were respectively not affected and lightly affected by herbivory in our study. In other studies on the putative influence of herbivory on bacterial communities, herbivory marginally increased the abundance of soil nitrifying bacteria and archeae (Le Roux et al., 2008) and, (in contradiction with our own results) Kong et al. (2016) found that white fly herbivory significantly decreased bacterial richness in the rhizosphere of pepper. Differences between these two earlier studies may be due to differences in the biological models studied but also to the nature and length of the perturbation and type of experiments. Compared to bacteria, we observed that fungal communities were very variable and only lightly influenced by herbivory. These results are consistent with previous studies showing that fungi are quite resistant to herbivory, such as leaf mining in Ageratina altissima (Asteraceae), which did not influence the communities of endophytes colonizing leaves (Christian et al., 2016).

Roots and rhizosphere microbial communities of healthy OSR were mainly composed of bacterial taxa Actinobacteria, Bacteroidetes, Firmicutes, and Proteobacteria (more specifically of Bacillus, Paenibacillus, and Pseudomonas at the genera level), while fungal taxa were mainly Ascomycota, Chytridiomycota and Basidiomycota. Our results on bacteria are in line with those of Bulgarelli et al. (2012) and Bodenhausen et al. (2013) in Arabidopsis thaliana. Regarding fungi, Ascomycota and Basidiomycota accounted for 96% of the community in Microthlaspi spp. populations, another brassicaceous plant (Glynou et al., 2018). A recent study on OSR showed that Proteobacteria and Actinobacteria as well as Ascomycota and Basidiomycota were the most abundant bacterial and fungal phyla respectively in the roots and rhizosphere, however with different dominant genera (Gkarmiri et al., 2017). While the phylum of Chytridiomycota was only present in roots in the previous study, it was part of the rhizosphere composition in ours. These differences in phyla and genera might be due to the use of isotopic labeling by Gkarmiri et al. (2017), representing only the active portion (i.e., assimilating labeled photosynthetates) of the plant communities. Furthermore, many drivers have been shown to shape plant microbial communities: plant root exudates (Turner et al., 2013; Tsunoda and van Dam, 2017); plant development (de Campos et al., 2013; Chaparro et al., 2014; Wagner et al., 2016); accession (Micallef et al., 2009); genotype (Wagner et al., 2016); species (Miethling et al., 2000); soil (Vandenkoornhuyse et al., 2015); and agricultural management practices (Hartmann et al., 2015; Rathore et al., 2017). It is therefore not surprising to find some differences in our study.

Root herbivory increased the abundance of γ-Proteobacteria in the roots and rhizosphere, and that of Firmicutes in the rhizosphere but it decreased abundance of α-Proteobacteria in both compartments and δ-Proteobacteria and Bacteroidetes in the roots only. These changes were accompanied by an increased abundance of Bacillus, Paenibacillus, Pseudomonas, and Stenotrophomonas. However, some of these phyla in infested plants at 42 DAI reached back a stable state similar to that found in healthy plants at 42 DAI, hence highlighting the resilience of microbial communities after the end of insect herbivory. These results are consistent with those of Kong et al. (2016) on whitefly herbivory on pepper showing an increase in the abundance of γ-Proteobacteria and Stenotrophomonas in the rhizosphere and with those of Kim et al. (2016) where aphid herbivory increased the abundance of Paenibacillus. According to Card et al. (2015), Bacillus, Paenibacillus, Pseudomonas, and Stenotrophomas are either inhibitors or antagonists of plant-pathogenic bacteria or fungi, while Bacillus and Paenibacillus are known to be also entomopathogenous (Monnerat et al., 2009; Grady et al., 2016; Zhao et al., 2017).

Herbivory also affected fungal phyla and genera but to a lesser extent than bacterial taxa: it decreased the abundance of Chytridiomycota in the rhizosphere and increased the abundance of Ajellomyces (slightly) and Filobasidiella (strongly) in the roots. Tkacz et al. (2015) demonstrated an increase of a Chytridiomycota in the rhizosphere of Brassica rapa over time, while this phylum progressively colonized the roots during plant development, to finally dominate fungal communities (Lebreton, personal communication). Most fungi are oligotrophic and grow slowly, with limited carbon sources (Ho et al., 2017). The lesser effect of herbivory on fungal communities observed in our study might be due to the timeframe of our experiment, probably too short compared to the growth rate of fungi, which is based on their utilization of complex trophic resources.

Herbivory Influences Root Metabolite and Chemical Element Composition

Herbivory tended to decrease root AAs and CPOs, which is consistent with the results of Hopkins et al. (1999), van Leur et al. (2008), and Lachaise et al. (2017). Herbivory increased glycerate and trehalose, a sugar assumed to play a role in plant defenses against aphid infestation (Singh et al., 2011). Sucrose was not affected by herbivory, contrary to the study of Pierre et al. (2012). S-methyl cysteine sulfoxide (SMCSO), a toxic compound detectable in brassicaceous crops (Edmands et al., 2013) was increased by herbivory. We therefore hypothesize that both SMCSO and trehalose might play a role in OSR defense against root herbivory.

Our study showed that herbivory also modulated GSLs by increasing two indolyl (glucobrassicin and neoglucobrassicin) but decreasing four aliphatic and one aromatic GSLs consistently with van Dam and Raaijmakers (2006). This is also coherent with the same decrease of total GSL contents found by van Dam and Raaijmakers (2006) and Lachaise et al. (2017). A recent study demonstrated the detrimental effect of glucobrassicin on aboveground generalist and specialist herbivores: a lower larval weight and faster development time associated with an increased mortality (Santolamazza-Carbone et al., 2017). However, in another study, aliphatic GSLs affected larval weight in generalists but not in specialist insects (Arany et al., 2008). Glucosinolates are generally considered as defensive compounds against herbivores but their impact on insect life history traits appear very species specific and their influence on microbial communities is difficult to predict.

Herbivory increased sulfur and nitrogen contents. The increase of sulfur content in root tissues probably corresponded to higher amounts of indolyl GSLs (glucobrassicin and neoglucobrassicin) and SMCSO, which are all sulfur-containing metabolites. The increase in nitrogen could be linked to insect frass (Kagata and Ohgushi, 2011). These authors demonstrated that there was more nitrogen excreted by larva than ingested and they suggested that the excess of nitrogen might originate from the plant organic nitrogen (e.g., AAs or proteins), which could not be digested by the insect. These observations could explain our results, especially since the nitrogen increased occurred only at 14 DAI, the peak of herbivory when the CRF was still at the larval stage.

Relationship Between Bacterial Genera and Root Chemistry in Response to Herbivory: Focus on Dominant Genera and Compounds Increased by Herbivory

Root herbivory—obviously—generates root degradation, which creates habitat spatial heterogeneity and as a consequence could modify microbial communities. Changes in microbial diversity and composition following herbivory is often hypothesized to be based on modifications of the plant chemistry, such as plant exudates (Kim et al., 2016; Kong et al., 2016). However, variations in root metabolite and rhizodeposits occur during the plant's life cycle. In addition, bacteria and fungi are known to differ in their ability to use organic compounds (Boer et al., 2005). In our study, we focused on the vegetative stage of OSR, a stage when fungal communities show only limited changes (Mougel et al., 2006).

Increased trehalose following herbivory was associated with an overexpression of Stenotrophomonas and Pseudomonas. Stenotrophomonas appears to produce trehalose as an osmoprotective substance (Wolf, 2002) while Pseudomonas is able to use this compound as a carbon source for its growth (De La Fuente et al., 2007) as well as other sugars found in root exudates, a trait that makes this genus quite competitive in microbial communities (Lugtenberg et al., 1999). A modification of plant microbial communities mediated by plant chemistry was also found in Kim et al. (2016) where bacteria such as Paenibacillus were able to grow on medium-based root exudates, which came from aphid-infested plants. Moreover, glycerate was also associated with microbial changes under herbivory. This could be explained by the fact that glycerate-derived compounds are usually accumulated by microorganisms under abiotic stress and this process could be similar for biotic stress (see review on abiotic stress by Empadinhas and da Costa, 2011). We suggest that the root chemistry disturbed by CRF herbivory might be beneficial to colonization of microorganisms such as Pseudomonas, exhibiting competitive traits to exploit trophic resources (Ho et al., 2017).

Following herbivory, enhancement of sulfur-containing compounds (i.e., SMCSO, glucobrassicin, neoglucobrassicin, sulfur) was associated to an overexpression of Bacillus, Paenibacillus, Pseudomonas, and Stenotrophomonas. Aziz et al. (2016) demonstrated that plant biomass loss caused by aboveground herbivory was lessened in presence of Bacillus. This decrease seemed mediated by an increase of indolyl GSLs. The fact that these bacteria are associated to sulfur-containing compounds might be linked to their potential role in the sulfur cycle. Pseudomonas for example dominates the arylsulfatase bacterial communities in the OSR rhizosphere (Cregut et al., 2009) and this enzyme (produced by microorganisms in order to mineralize organic sulfur) is present both in the roots and in the rhizosphere of this plant (Knauff et al., 2003). A plant undergoing herbivory could send external signals to recruit bacteria, which could stimulate the synthesis of plant defense compounds (e.g., GSLs). An increase of trehalose in the roots and rhizosphere via root exudates (Paterson et al., 2007) could represent an external signal, attracting beneficial symbionts and promoting bacterial growth. These bacteria could then enhance enzymatic activities (e.g., arylsulfatase) allowing sulfur mineralization and stimulation of GSL production.

Following herbivory, nitrogen was also associated with microorganisms but to a lesser extent. Nitrogen is a limiting factor for plant growth because of its low bioavailability in the soil, since eukaryotes are not able to fix atmospheric nitrogen (Grady et al., 2016). To have access to nitrogen sources, plants rely on nitrogen-fixing microorganisms (i.e., diazotroph). Clostridium and Paenibacillus are free-living diazotrophic species (Choudhary and Varma, 2016; Grady et al., 2016) while a nitrogen fixing trait was recently discovered in the genome of a Pseudomonas species (Yan et al., 2008). We suggest that plants stressed by herbivory might be able to select better microorganisms (i.e., in terms of fitness and cost in nutrients, such as carbohydrates, to the plant), which will provide more nitrogen and enable them to better resist this perturbation.

Soil Microbial Diversity Did Not Influence the Plant-Microorganism Interactions

Initial soil diversity had respectively little and no effect on root bacterial and fungal alpha-diversities. Only rhizosphere bacterial and fungal diversities were impacted by soil microbial diversity, but not by herbivory. When evaluating abiotic and biotic stresses (Kissoudis et al., 2014), some studies established a positive relationship between an ecosystem stability-resilience and its diversity, while others found opposite relationships or none at all (Griffiths and Philippot, 2013). In our study, the resilience of microbial communities did not seem to come from their initial diversity: microbial results were globally similar at high and low initial diversities and plant chemistry was not impacted by these levels of initial diversity either, as found by Lachaise et al. (2017). We hypothesize that our “low” microbial diversity condition might still have been too rich to properly assess the effect of soil diversity and that a stronger dilution may be necessary in future studies.

Conclusion

In summary, herbivory led to root chemical changes, involving carbohydrates and sulfur-containing compounds, which partly shaped belowground microbial communities and particularly the phyla of γ-Proteobacteria and Firmicutes and a couple of their affiliated genera. It indicates that a plant suffering from herbivory emits either defensive and/or nutritive compounds that influence the recruitment of soil microorganisms by the rhizosphere and the roots.

Our results encourage the determination of the precise identity and functions of microorganisms responding to herbivory, using more accurate primers and metatranscriptomic approaches, and the understanding of the feedback-loop existing between these identified microorganisms and the plant chemistry modified by herbivory. This future research could represent the next step toward the transition from correlation to causation in order to develop sustainable and innovative plant protection strategies.

Author Contributions

MO, CM, and A-MC conceived and designed research. LL and A-YG-E prepared the soil used in the experiment. MO, CP, and VC conducted plant harvest. JL and CM conducted respectively the molecular and bioinformatics analyses. MO performed metabolomics analysis, advised by NM. MO did the statistical analyses and interpreted data, advised by LL, MH, AO, and DP. MO wrote the manuscript, which was commented and approved by all authors.

Funding

This work was supported by grants from the Plant Health and Environment division of the French National Institute for Agricultural Research (INRA) (AAP2015) and a grant from Rennes Métropole for CM (300 01003).

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

The authors thank the GenoScreen platform for sequencing analyses and the Metabolic and Metabolomics Profiling Platform (P2M2) for metabolomics analyses. We are most grateful to PLATIN' (Plateau d'Isotopie de Normandie) core facility for all element analysis used in this study. We also thank Christophe Lunel for his help during the plant sampling.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fevo.2018.00091/full#supplementary-material

Abbreviations

AA, Amino Acid; CPO, Carbohydrate, Polyol, and Organic Acid; CRF, Cabbage Root Fly; DAI, Days After Infestation; dbRDA, distance-based Redundancy Discriminant Analysis; DIABLO, Data Integration Analysis for Biomarker discovery using Latent cOmponents; FDR, False Discovery Rate; GSL, Glucosinolate; OSR, Oilseed Rape; RDA, Redundancy Analysis; SMCSO, S-Methyl Cysteine Sulfoxide.

References

Ahuja, I., Rohloff, J., and Bones, A. M. (2010). Defence mechanisms of Brassicaceae: implications for plant-insect interactions and potential for integrated pest management. A review. Agron. Sustain. Dev. 30, 311–348. doi: 10.1051/agro/2009025