Temporal dynamics of bacterial and fungal communities during the infection of Brassica rapa roots by the protist Plasmodiophora brassicae

The temporal dynamics of rhizosphere and root microbiota composition was compared between healthy and infected Chinese cabbage plants by the pathogen Plasmodiophora brassicae. When inoculated with P. brassicae, disease was measured at five sampling dates from early root hair infection to late gall development. The first symptoms of clubroot disease appeared 14 days after inoculation (DAI) and increased drastically between 14 and 35 DAI. The structure of microbial communities associated to rhizosphere soil and root from healthy and inoculated plants was characterized through high-throughput DNA sequencing of bacterial (16S) and fungal (18S) molecular markers and compared at each sampling date. In healthy plants, Proteobacteria and Bacteroidetes bacterial phyla dominated the rhizosphere and root microbiota of Chinese cabbage. Rhizosphere bacterial communities contained higher abundances of Actinobacteria and Firmicutes compared to the roots. Moreover, a drastic shift of fungal communities of healthy plants occurred between the two last sampling dates, especially in plant roots, where most of Ascomycota fungi dominated until they were replaced by a fungus assigned to the Chytridiomycota phylum. Parasitic invasion by P. brassicae disrupted the rhizosphere and root-associated community assembly at a late step during the root secondary cortical infection stage of clubroot disease. At this stage, Flavisolibacter and Streptomyces in the rhizosphere, and Bacillus in the roots, were drastically less abundant upon parasite invasion. Rhizosphere of plants colonized by P. brassicae was significantly more invaded by the Chytridiomycota fungus, which could reflect a mutualistic relationship in this compartment between these two microorganisms.


Introduction
All plant tissues including roots [1,2], leaves [3,4] and seeds [5,6] are surrounded by a large diversity of microorganisms assembled in microbial communities or microbiota. These microbial assemblies represent a continuum of symbiosis with the plant ranging from parasitic to PLOS  the first one involving a rapid recruitment of microorganisms in the vicinity of the root and a second step being their entry inside the root [29]. However, the second step is more complex than the first one, with each root niche playing a selective role in microbiota assembly [15]. Even though the composition and recruitment mechanisms of these communities are being extensively investigated, only few studies dealt with the stability of these assemblies during the plant growth and under the effect of biotic stresses. Among biotic stresses, soilborne plant pathogens cause major economic losses in agricultural crops. Most of them are adapted to grow and survive in bulk soil but can also invade the root tissues to establish parasitic relationships with the plant. Since soilborne pathogens are already present in the soil before sowing, infections are usually early and occur during the vegetative stages of plant growth. To infect root tissues, pathogens have to compete with other microorganisms of the rhizosphere microbiota for available nutrients and microsites. One of the major roles of the rhizosphere microbial communities could to provide a frontline defense for plant roots against infection by soilborne pathogens [19]. Some of the mechanisms involved in the activity of these beneficial rhizosphere microorganisms are well studied and include several direct interactions with plant pathogens as well as indirect interactions via the plant by stimulating the plant immune system [30][31]. These mechanisms are well documented, using specific strains, for some rhizobacteria like Pseudomonas sp. and Bacillus sp., and for some fungi like Trichoderma sp. and non-pathogenic Fusarium oxysporum. However, most of the responsible microbial networks underlying these defense mechanisms are currently largely unknown. Recently, some metagenomic approaches provided us new opportunities to enrich our knowledge about the strong interactions between telluric pathogens and their living environment [32].
Among soilborne pathogens, Plasmodiophora brassicae is responsible for clubroot disease, a serious disease for many members of the Brassicacae family. P. brassicae Woronin is an obligate protist within the class Phytomyxea (plasmodiophorids) of the protist supergroup Rhizaria [33]. The pathogen life involves three stages: survival in the soil as resting spores, root hair primary infection and finally secondary cortical infection [34]. This process is accompanied by the hyperplasia and the hypertrophy of infected roots, resulting in formation of clubshape galls on the root. The tissue disruption associated with large clubs reduces nutrient and water transport within the plant, and consequently reduces plant growth and yield. In our experiments, Brassica rapa subsp. Pekinensis (Chinese cabbage) was chosen as the plant model of Brassicaceae because the full P. brassicae life cycle was easily achieved under controlled conditions in this species. We specifically addressed the following questions: (i) what is the dynamics of root and rhizosphere communities of Chinese cabbage during the vegetative stages of plant growth? (ii) How does P. brassicae affect the composition of bacterial and fungal rhizosphere and root communities at each of its life cycle stages? and (iii) which microbial species are selected following the infection by P. brassicae? To address these questions, a timeseries experiment was conducted under controlled conditions. We described the dynamics of the root and rhizosphere (including rhizoplane) bacterial and fungal communities of noninoculated (also called "healthy") plants at successive time points during its vegetative stage. We also analysed the effects of a parasitic infection by Plasmodiophora brassicae, on the composition and dynamics of these microbial communities.

Materials
Soil. The experimental soil used for this study was collected at the INRA experimental site of La Gruche in Western Brittany (N: 48˚08.44', W: 01˚47.98'). The topsoil (0-5 cm) was removed and the layer between -5 and -30 cm was harvested, homogenized, sieved at 4 mm and subsequently stored in 500 L containers at ambient temperature in the dark until further used. Physical and chemical properties of the soil were determined at the Arras soil analysis laboratory (F-62000, Arras, France). These properties were determined as: 13.3% sand, 70.9% silt, 15.8% clay, pH 6.2, 12.0 g.kg -1 of organic carbon, 1.2 g.kg -1 of mineral N and 20.8 g.kg -1 of organic matter.
Pathogen. The eH selection isolate [35], used in this study, belongs to the most virulent P. brassicae pathotype P1 [36]. This isolate was kindly provided by J Siemens (University of Dresden, Germany). It was propagated on Chinese cabbage as root galls, harvested, washed and stored at -20˚C.

Experimental design
Plant growth assay and inoculation. ECD5 plants were cultivated in pots filled with 400 g of the experimental soil mixed with sterilized sand in the ratio 2:1. The experiment was conducted under a randomized complete block design using three blocks consisting of three replicates of four plants each. In each block, replicates were randomly distributed and placed in a greenhouse under the following conditions: 16 hours light (day) at 22˚C and 8 hours dark (night) at 19˚C. A mean photosynthetically active photon flux density of 150 μmol.m -2 .s -1 at plant level during the 16 hours daylight was maintained. Some pots without plants were designated "bulk soil".
Inoculum was prepared from three galls stored at -20˚C as described previously [37]. In brief, spores were extracted by thawing the frozen galls at room temperature, and then homogenizing in 100 mL of sterilized water at high speed for 2 min. The resulting spore suspension was filtered through two sieves (250 and 100 μm pore diameters). The spore concentration was determined with a Malassez cell and adjusted to 1 × 10 7 spores.mL -1 . Ten-day-old seedlings were inoculated by pipetting 1 mL of spore suspension containing 1 × 10 7 spores.mL -1 onto the soil surface at the base of each seedling. Non-inoculated plants and bulk soil were poured with sterile water. All pots, including bulk soil controls, were watered periodically every three days from the top with 8 mM Hoagland solution to maintain a water retention capacity between 70-100%.
Symptom development and clubroot severity measurement. Disease severity was assessed in inoculated plants during the vegetative stage of plant growth at 0 (T1), 7 (T2), 14 (T3), 23 (T4) and 35 (T5) days after inoculation (DAI) with P. brassicae, corresponding to 10, 17, 24, 33 and 45 DAS, respectively. Clubroot severity was recorded using the scale: 0, no visible swelling; 1, very slight swelling usually confined to lateral roots; 2, moderate swelling on lateral roots and taproot; 2+, severe clubs on all roots, but some roots remain present; and 3, no root left, only one big gall. A disease index (DI) was calculated as described by [42]: DI = (n1 � 25 + n2 � 50 + n2+ � 75 + n3 � 100)/N, where "ni" is the number of plants in the symptom class "i" and N is the total number of plants tested. Disease data were analyzed using a likelihood ratio test on a cumulative link model [CLMM; function "clmm", package "RVAideMemoire", [43]]. Pairwise comparisons of LSMeans were then computed. To measure the hypertrophy of infected root, taproot width was also assessed at each date of sampling at 1 cm under the soil surface. Taproot width data were compared between non-inoculated (or healthy) and inoculated (or diseased) plants using a linear model [LMM; function "lmer", package "lme4"]. Pairwise comparisons of LSMeans [function "lsmeans", package "lsmeans"] and FDR correction for p-values were then performed.
Sampling of "rhizosphere", "root" and "bulk soil" compartments Rhizosphere and root compartments from healthy and diseased plants were sampled at 10 (T1), 17 (T2), 24 (T3), 33 (T4) and 45 (T5) DAS. The "rhizosphere compartment" defined as the soil particles firmly attached to roots was collected by centrifugation of root washings. The "root compartment" was defined as the root tissues depleted of soil particles and epiphytic bacteria by sequential washing and sonication treatments and was therefore enriched in rootinhabiting bacteria.
Rhizosphere and root samples were collected from planted pots in a soil depth of -1 to -6 cm from the surface. Roots were separated from non-adhering soil particles, collected in 15 mL Falcon containing 20 mL sterile water and vortexed for 1 min. Seminal and nodal roots were included in the analysis. After vortexing, roots were transferred into a sterile Petri dish and subjected to a second washing treatment with 5 mL sterile water. Double washed roots were transferred in 5 mL sterile water and sonicated twice for 3 s at 40 Hz to detach microbes living in close association with root tissues. Roots were transferred into a Petri dish, cut into fragments smaller than 5 mm, ground to a powder with a pestle in liquid nitrogen-chilled mortar with Fontainebleau sand and stored at -80˚C until further analysis. The soil suspensions collected in Falcon tubes or in the Petri dishes after the first, the second washing treatments and the sonicated solution were pooled, centrifuged at 4,000 g for 20 min and the pellet, referred to as the rhizosphere, was frozen in liquid nitrogen and stored at -80˚C until further analysis.
Soil samples were collected from unplanted pots at T1, T3 and T5 in a soil depth of -1 to -6 cm from the surface. The soils from four pots were pooled, transferred in 10 mL sterile water and vortexed for 1 min. The soil suspension was centrifuged at 4,000 g for 20 min and the pellet, referred as the bulk soil, was frozen in liquid nitrogen and stored at -80˚C until further analysis.

DNA extraction and pathogen quantification
Root and soil DNA extraction. The GnS-GII protocol was used for root and DNA extraction [44]. For root samples, five to 150 mg of each root sample were homogenized for 3 x 30 s at 4 m.sec -1 in a FastPrep-24 (MP-Biomedicals, Illkirch, France) in 2 mL of the "lysing matrix E" solution from MP-Biomedicals (Illkirch, France) containing 100 mM Tris (pH 8.0), 100 mM EDTA (pH 8.0), 100 mM NaCl, and 2% (wt/vol) sodium dodecyl sulphate. The samples were incubated for 30 min at 70˚C, and then centrifuged at 7,000 g for 1 min at 20˚C. To remove proteins from the extracts, 1 mL of the collected supernatant was incubated for 10 min on ice with 1/10 volume of 3 M potassium acetate (pH 5.5) and centrifuged at 14,000 g for 5 min at 4˚C. Finally, after precipitation with 900 μL of ice-cold isopropanol, the nucleic acids were washed with 70% ice-cold ethanol and DNA was resuspended in 200 μL ultrapure sterile water. DNA was separated from the residual impurities, particularly humic substances, by centrifuging through two types of minicolumns. Firstly, aliquots (100 μL) of crude DNA extract were first loaded onto Microbiospin (Biorad, Hercules, California, USA) columns of PVPP (PolyVinyl PolyPyrrolydone) and centrifuged at 1,000 g for 2 min at 10˚C. Secondly, the eluate was purified with the Geneclean turbo kit (Q-Biogene, Illkirch, France). DNA concentration and purity were determined with a Nanodrop (Agilent Technologies, Les Ulis, France).
The same protocol was used to extract DNA from soil samples except that, before homogenization in the FastPrep-24, 2 g of each soil sample were mixed with 5 mL of a solution containing 100 mM Tris (pH 8.0), 100 mM EDTA (pH 8.0), 100 mM NaCl, and 2% (w/v) sodium dodecyl sulphate in a 15 mL "lysing matrix E" Falcon tube from MpBio.
Measurement of pathogen DNA amount in roots by real-time qPCR. Plant root colonization by P. brassicae was also monitored by quantitative PCR. The predicted 18S gene was used to estimate P. brassicae DNA amount per ng of total extracted DNA. A portion of this gene sizing 164 bp was amplified with the primers PbK1F/PbK1R (5'-TTGGGTAATTTGCG CGCCTG-3'/5'-CAGCGGCAGGTCATTCAACA-3'). All reactions were performed in 20 μL qPCR reaction using 10 μL of SYBR Green Master Mix (Roche), 1 μL of DNA (2.5 ng) and 0.08 μL of each primer (100 μM). Amplification conditions were as follows: 5 min at 95˚C, followed by 45 two-step cycles at 95˚C (10s) and 60˚C (40s). Standard curves were constructed using serial dilutions of P. brassicae DNA extracted from resting spores. A linear model [LMM; function "lmer", package "lme4"] was used to analyze the pathogen DNA quantification data. Pairwise comparisons of LSMeans [function "lsmeans", package "lsmeans"] and FDR correction for p-values were then performed.

Bacterial and fungal community composition and diversity
Sequencing of 16S and 18S rDNA genes. The structure of microbial communities associated to soil and root samples collected during the experiments was assessed though amplification and subsequent sequencing of bacterial (16S) and fungal (18S) rDNA genes. PCR amplification and sequencing were performed at GenoScreen (Lille, France) using the Illumina Miseq platform to a 2 × 300 bases paired-end version with an adequate read assembly method. For soil and root DNA extracts, a 420 bp fragment of the V5-V7 region of the bacterial 16S rDNA gene was amplified using the universal bacterial primers 799F_16S (5'-AACMGGATTA GATACCCKG-3') and 1223R_16S (5'-CCATTGTAGTACGTGTGTA-3') [45,46]. Before sequencing, PCR products were purified to eliminate a 760 bp fragment corresponding to plant mitochondrial DNA amplification. A 530bp fragment of the fungal 18S rDNA that includes the variable regions V4 (partial) and V5 was also amplified using the primer pair NS22B (5'-AATTAAGCAGACAAATCACT-3') and SSU0817 (5'-TTAGCATGGAATAATR RAATAGGA-3') [47,48].
Analysis of MiSeq sequencing data. After reads assembly, sequences were processed with GnS-PIPE bioinformatics platform developed by GenoSol platform and optimized for amplicons analysis [49,50]. The reads were filtered and eliminated if they harbored one or more ambiguities (Ns) or an average quality score below 30. A PERL program was applied to obtain strict dereplication (i.e., clustering of strictly identical sequences). After this initial quality filtering step, the reads were aligned with INFERNAL alignments [51] and clustered at 97% sequence similarity into operational taxonomic units (OTU) using another PERL program. All single-singletons (reads detected only once and not clustered) were then deleted to eliminate PCR chimeras and large sequencing errors. These final sequences were used to produce rarefaction curves. The retained high-quality reads were used for taxonomy-based analysis of each OTU using similarity approaches against dedicated reference databases from SILVA [52]. The raw data sets are available on the European Nucleotide Archive database system under the project accession number PRJEB26948. Root and soil samples accession numbers range from ERS2513216 to ERS2513353 for 16S and 18S rDNA.
Alpha diversity. To compare bacterial or fungal composition among bulk soil, rhizosphere soil and root from healthy and diseased plants, the richness was characterized by the number of OTUs found in each sample. As metric of taxonomy diversity, Shannon diversity was also determined using the "vegan" package in R, version 2.2-1 [53]. Since values were conformed to normality assumptions, two-way Anova and post-hoc Tukey's HSD test were used to examine pairwise differences between samples for these measures.
Beta diversity. After normalization by sample size, OTU counts without at least a mean of one read per sample were removed from the analysis. The genera OTU counts were also rarefied to 1,000 counts per sample and Log2-transformed rarefied values were used to calculate a Bray-Curtis distance dissimilarity matrix using the function "vegdist" of the R package "Vegan". The beta diversity distance matrices were plotted using a bi-dimensional Principal Coordinates Analysis (PCoA) using the function "plot" of the R package "Vegan". To quantify the influence of each factor on the beta diversity, a canonical analysis of principal coordinates (CAP, [54]) followed by a permutation-based ANOVA (PERMANOVA) was performed using the R package "vegan" according to the method described by [55].
Statistical analysis on phyla counts. To identify phyla enriched in rhizosphere and root microhabitats compared to unplanted soil and to compare phyla composition between samples collected from healthy and diseased plants, we employed linear statistics on Relative Abundances (RA) values (log2 > 5‰ threshold) using a script developed from the R package "Limma". Differentially abundant phyla between two samples were calculated using moderated t-tests. The resulting p-values were adjusted for multiple hypotheses testing using the Benjamini-Hochberg (BH) correction.
Detection of differentially enriched OTUs. EdgeR is a workflow largely based on the free open-source R language and Bioconductor software [56]. This workflow was originally used to analyze count-based differential expression of RNA sequencing as part of transcriptome studies [57] and was recently adapted to metagenomic data analysis [27]. OTU counts without at least a mean of one read per sample were removed from the analysis. To normalize the data for each sample OTU count, the trimmed mean of M values normalization method (TMM) was used according to the method described by [58]. A Log2-transformation was performed on the normalized data for statistical comparisons. Threshold, normalization and transformation steps were performed using a custom R script. To identify differentially abundant genus in bacterial and fungal communities between sampling dates and treatments (non-inoculated or inoculated) in root or soil samples, EdgeR was used to fit a model with treatment (non-inoculated or inoculated) � sampling date (T1 to T5) terms to the count data in each compartment by using glmFit and glmLRT with tagwise dispersion and to test for significant effects of each term. EdgeR employs statistical methods supported on negative binomial distribution as a model for count variability. Data from root and rhizosphere soil were not analyzed together because composition biases between samples from these two compartments were not eliminated by TMM normalization. To examine whether having a diverged or conserved communities composition was associated with treatment � time effect, the model was fitted to subsets of the normalized counts data and used "contrasts" to identify genera with significant differential abundances in pairwise comparisons. A likelihood ratio test (LRT) was performed to specify the difference of interest and the resulting p-values were adjusted for multiple hypotheses testing using the Benjamini-Hochberg (BH) correction.

Comparison of communities from root, rhizosphere and bulk soils in healthy plants: The rhizosphere effect
In the samples collected at T1, T3 and T5 from healthy plants and bulk soil, the greatest numbers of bacterial/fungal OTUs were detected in bulk and rhizosphere soils (2,240/1,242 and 2,280/1,220 OTUs on average, respectively). A significant reduction of richness was observed in root compartment (530/677 OTUs on average) compare to bulk and rhizosphere soils (S1 and S2 Figs). A significant reduction of bacterial and fungal diversities in the root samples compared to bulk and rhizosphere soils was also observed at each sampling date (S1 and S2 Figs). A temporal effect on bacterial richness and diversity was measured but only in the root compartment where the number of OTUs and the Shannon index were higher at T3. In each compartment, no temporal variations of fungal richness and diversity was measured.
When looking at the microbial composition, we found that root bacterial and fungal communities were clearly distinct from rhizosphere and bulk soil communities at each sampling date ( Consistently, we observed at T5 a clear separation between root microhabitat and soil samples followed by segregation of the rhizosphere and bulk soil samples. To explain the variance observed, the significant effect of the sampling date was weaker than the compartment.

Composition and dynamics of healthy Chinese cabbage rhizosphere and root microbiota
In the rhizosphere of healthy plants. In the rhizosphere of healthy plants, the most abundant bacterial phyla found were Proteobacteria, Firmicutes, Actinobacteria and Bacteroidetes, with 86% to 90% abundances at each sampling date between T1 and T5 (Fig 1). Within the rhizosphere-inhabiting Proteobacteria, the α-Proteobacteria were over-represented compared to the β-, γ-and δ-Proteobacteria (Fig 1). Between T1 and T5, a significant increase of α-Proteobacteria and a decrease of Firmicutes were measured in rhizosphere soil while no temporal variation of bulk soil composition at phylum level was observed (S5 Fig). At T5, the enrichment of members from the Proteobacteria and Bacteroidetes phyla significantly discriminated rhizosphere from bulk soil samples (S5 Fig). We tried to narrow down the bacterial community to those OTUs (� 97% sequence similarity), which showed a minimum relative abundance of 0.1% at least in one of the rhizosphere samples. A total of 429 OTUs were identified in the rhizosphere of healthy plants (S1 Table). At the genus level, OTU1 assigned as Bacillus (Firmicutes) dominated these rhizosphere communities at each sampling date with 12 to 18% relative abundances between T1 and T5. OTU4 (Sphingomonas, α-Proteobacteria), OTU7 (Pseudolabrys, α-Proteobacteria), OTU9 (Sporosarcina, Firmicutes), OTU6 (Bradyrhizobium, α-Proteobacteria), and OTU10 (Rhodopseudomonas, α-Proteobacteria) were also highly represented (S1 Table). No temporal variation of these dominant OTUs was observed between T1 and T5. However, several minor OTUs with significant relative abundance variations between two sampling dates were detected in these bacterial communities ( Table 1).
The fungal rhizosphere communities from healthy plants were largely dominated by Ascomycota, with 64% to 69% relative abundances between T1 and T4 (Fig 1). Three other phyla, Mucoromycotina (12% to 16%), Basidiomycota (3% to 8%) and Chytridiomycota (2% to 10%) were also detected but at lower frequencies. Until T4, the proportions of these four phyla compared to the whole fungal microbiota were relatively stable. From the 168 fungal OTUs identified with at least a relative abundance of 0.1% in one of these rhizosphere samples, OTU2 and OTU4, assigned as two Sordariomycetes, were detected at high frequencies between T1 and T4, but no temporal variations of these dominant OTUs were observed. At T5, fungi from the Chytridiomycota phylum were drastically more abundant (47.5%) than at the beginning of the kinetics and a decrease of Ascomycota (64.6% at T4 to 37.9% at T5) was measured (Fig 1). At this date of sampling, OTU1 and OTU16, assigned as two Chytridiomycota, were the most abundant OTUs with 18.3 and 16.4% relative abundances, respectively (S1 Table). Variations of other less dominant OTUs were also observed between T1 and T5 (Table 1). While the proportion of Chytridiomycota fungi strongly increased in the rhizosphere compartment at T5, their relative abundances remained low in bulk soil samples during all the kinetics.
Inside the roots of healthy plants. In the roots of healthy plants, bacterial communities were also dominated by Proteobacteria (41% to 51% RA between T1 and T5) and Bacteroidetes (21% to 33%) as in the rhizosphere (Fig 1). They also contained Actinobacteria (1% to 10%) and Firmicutes (3% to 7%), although to a lesser extent than in the rhizosphere samples. In the root samples, the α-and γ-Proteobacteria were over-represented compared to the β-and δ-Proteobacteria. Cyanobacteria were also detected. Between T1 and T5, more important variations of phylum frequencies occurred in the roots of healthy plants than in their rhizosphere (Fig 1). Actinobacteria increased significantly in frequencies while Proteobacteria decreased. A total of 202 genera were identified in these communities. The dominant OTUs were OTU2 assigned to a Flavisolibacter (Bacteroidetes) which relative abundances varied between 8% and 17%, OTU3 assigned to an unknown Cyanobacterium, OTU19 (Devosia, α-Proteobacteria), OTU12 (Pseudomonas, γ-Proteobacteria) and OTU5 (Flavobacterium, Bacteroidetes) (S1 Table). While the proportion of OTU2 strongly increased in the rhizosphere compartment of healthy plants between T1 and T5, the proportion of OTU19 decreased from 4.2% to 1.0% and no temporal variations was observed for the other main OTUs. Several other bacterial OTUs with significant variations in abundances from one date of sampling date to another were also detected in healthy plant roots ( Table 1).
The root fungal communities were dominated by Ascomycota (85.1%) at T1 and replaced progressively by fungi from the Chytridiomycota phylum during the kinetics of plant growth. At T5, OTU1 assigned to the Chytridiomycota phylum was detected in the roots of healthy plants at a very high frequency with a mean of 53% relative abundance (S1 Table). Variations of several minor OTUs were also observed in the fungal communities of diseased plants during the time-series experiment ( Table 1).
To conclude, weak fluctuations were measured in the composition of rhizosphere and root communities of healthy plants before T4, whereas important changes occurred in these communities between T4 and T5. According to the variations of OTU relative abundances, these changes were first observed in the bacterial and fungal communities of plant roots and then in the rhizosphere (Table 1).  The impact of a pathogen on the plant root and rhizosphere microbiota

Symptom development and clubroot severity
Differences of taproot width between healthy and diseased plants appeared at T3 and increased drastically between T3 and T5 (Fig 2). Disease index was low at T3 (DI = 16.7%), increased rapidly to 68.5% at T4 and reached a maximum of 86% at the end of the experiment (Fig 2). The amount of P. brassicae DNA followed a similar evolution profile (Fig 2). However, although DI increased between T4 and T5, there were no significant variations of P. brassicae DNA amount in roots at these time-points. Throughout the time of the experiment, no differences of leaf number and leaf area, plant height and shoot biomass were observed between healthy and diseased plants. In contrast, root length and biomass of inoculated plants decreased significantly but only between T4 and T5 (S2 Table). At the end of the experiment, some galls had become brownish and some mature resting spores were observed in gall tissues. According to these results, the duration of the life cycle of P. brassicae infection in Chinese cabbage was approximately 35 days in our experimental conditions. The root hair infection and the beginning of the cortical infection stages occurred before 14 DAI, and clubroots formed gradually between 14 and 35 DAI.

Effect of P. brassicae on the rhizosphere and root microbiota of Chinese cabbage
In the rhizosphere of diseased plants compared to healthy plants. In the rhizosphere compartment, no significant variation of bacterial richness and diversities was measured In bacterial communities, the sampling date explained 21.2% of the overall variance of the data (p = 0.001, 95% CI = 19.5%, 24.3%) and the inoculation condition (I vs NI) 4.4% of this variance. This proportion of the variation, albeit small, was found significant (p = 0.02, 95% CI = 3.3%, 5.3%). The microbial dynamics of healthy and diseased plant communities clearly diverged from T3 to T5 as visualized in PCoA (Fig 3). At T4, there is no variation of bacterial phyla_subphyla relative abundances between healthy and diseased plants (Fig 4), but two OTUs (OTU35 and OTU188) assigned to two genera from the α-Proteobacteria phylum (Sphingopyxis and Rhodobacter, respectively) and two non-assigned β-Proteobacteria (OTU54 and OTU151) became more abundant in the rhizosphere of diseased than healthy plants (Fig 5A). At T5, Proteobacteria (α, β and γ) and Bacteroidetes were consistently more abundant in the rhizosphere of diseased than healthy plants, while both Firmicutes and Acidobacteria were less abundant (Fig 4). At this sampling date, 20 OTUs belonging mainly to the Proteobacteria, Bacteroidetes and Firmicutes phyla were significantly more abundant in inoculated than in non-inoculated plant samples and 8 less abundant (Fig 5B). Among these 28 rhizospheric OTUs, the more frequent ones were OTU1 (Bacillus) that decreased between T1 and T5 in the rhizosphere of all plants but more drastically in diseased plants especially at T5, and OTU5 (Flavobacterium), OTU14 (Dokdonella), OTU17 (Pseudomonas), OTU35 (Sphingopyxis), OTU54 (unknown β-Proteobacteria) which were all significantly more abundant in inoculated than non-inoculated samples at T5 (Fig 6). In fungal communities, the sampling date explained a higher proportion of the variance than in bacterial communities (35%, p = 0.001, 95% CI = 26.1%, 49%), while the inoculation condition (inoculated vs non-inoculated) had no significant effect (3.9%, p = 0.077, 95% CI = 2.8%, 5.6%). Until T4, no variation of fungal phylum frequencies was observed (Fig 4). At T5, while Ascomycota, Basidiomycota and Mucoromycotina were less abundant in the rhizosphere of diseased than healthy plants, no significant variation of Chytridiomycota was observed (Fig 4). However, the major OTU (OTU1) assigned to the Chytridiomycota phylum significantly increased in diseased plant, while four minor OTUs also varied: higher relative abundances for OTU55 and OTU60 but lower for OTU11 and OTU20 in diseased than healthy plant samples at T5 (Fig 7). Higher changes of OTU relative abundances occurred in diseased than healthy plant rhizosphere communities during the time-series experiment ( Table 1).

Inside the roots of diseased plants compared to healthy plants.
In the root compartment, no clear significant differences of bacterial and fungal richness and diversity of communities from healthy and diseased plants were found at each date of sampling (S6 and S7 Figs). In bacterial communities, the sampling date explained 24.4% of the overall data variance (p = 0.001, 95% CI = 20.9%, 28.1%) and the inoculation condition (inoculated vs non-inoculated) 6.2% of this variance (p = 0.002, 95% CI = 4.4%, 8.2%). No significant differences in community composition between inoculated and non-inoculated root samples were observed until T4 when one bacterial OTU (OTU2) assigned to the Flavisolibacter genus decreased drastically in relative abundances, while six minor OTUs (OTU54, OTU151, OTU122, OTU150 and OTU422 and OTU440) were slightly more abundant in inoculated than non-inoculated samples (Fig 5A). At T5, Actinobacteria were less abundant in the roots of diseased than (T4) and 35 (T5) days after inoculation (DAI) by P. brassicae. Phyla with relative abundances below 1% were grouped as "others". At each sampling date, significant (p-values � 0.05) and non-significant differences between NI and I plants are indicated by stars and "ns", respectively. In bacterial communities, the Proteobacteria phylum was divided into four subphyla: α-, β-, γ-and δ-Proteobacteria. https://doi.org/10.1371/journal.pone.0204195.g004 The impact of a pathogen on the plant root and rhizosphere microbiota healthy plants but β-Proteobacteria were more abundant at the phyla-subphyla level (Fig 4). We observed significant differences in relative abundances of 28 OTUs between inoculated and non-inoculated root samples (Fig 5B). Among these 28 OTUs, OTU2 (Flavisolibacter), OTU21 (Streptomyces) and OTU44 (Pseudomonas), were the main OTUs which relative abundances had decreased in diseased plants on one hand (Fig 8). On the other hand, the main OTUs which frequencies increased in inoculated vs non-inoculated plants were OTU17 (Pseudomonas) but also the two non-assigned β-Proteobacteria OTU54 and OTU62 (Fig 8). Regarding fungal communities, the date of sampling accounted for a higher proportion of the variance than in bacterial communities (36.6%, p = 0.001, 95% CI = 26.2%, 52.9%), while the condition (inoculated vs non-inoculated) had no significant effect (2.7%, p = 0.55, 95% CI = 1.7%, 3.9%) as in the rhizosphere. At each date of sampling, there was no difference in fungal phylum (Fig 4) and OTUs frequencies between diseased and healthy root samples. However, changes of OTU relative abundances occurred in root communities of healthy and diseased plants during the time-series experiment ( Table 1).

Discussion
In our study, the stability of assembled root and rhizosphere communities of Chinese cabbage was investigated by a time-series experiment, during the plant growth and under the effect of the parasitic invasion by P. brassicae. During the plant growth, for healthy plants, most of Ascomycota fungi previously recruited by the plant were replaced, mainly in the root compartment, by a Chytridiomycota fungus. The root and rhizosphere-associated community assemblies were also strongly modified by P. brassicae infection during the secondary cortical infection stage of clubroot disease.
Framed numbers indicate the number of OTUs with significant different frequencies between I and NI plants. (B) At T5, bacterial OTUs that significantly differed in their relative abundances (expressed in %) in the root and rhizosphere samples between I and NI plants are represented. Each histogram bar represents the mean RA (± SEM) of three replicates. Only significant differences (p-values � 0.05) between I and NI plants (represented by white and black bars, respectively) are shown in both compartments. Framed numbers indicate the number of OTUs with significant different frequencies between I and NI plants. https://doi.org/10.1371/journal.pone.0204195.g005

A weak but significant rhizosphere effect
Clearly, the communities that assembled in the rhizosphere and bulk soils of healthy plants were very different from the communities found in the roots. These results are consistent with earlier findings on other plant species [2,10,15,59,60]. We found a significant "rhizosphere effect". Indeed, the alpha diversity analysis showed that the microbiota diversities of the bulk and rhizosphere soils were not distinct from each other. These observations are similar to the findings of several authors in Arabidopsis thaliana who reported the resemblance of bacterial communities between rhizosphere and bulk soil samples in multiple soil types [2,10]. At each sampling date, bulk soil and rhizosphere compartments shared a large proportion of OTUs. However, the enrichment of OTUs assigned to the Proteobacteria and Bacteroidetes bacterial phyla, but also to the Chytridiomycota fungal phylum, significantly discriminated rhizosphere from bulk soil samples at the end of the experiment.

The structure of microbial communities associated with the rhizosphere and roots of healthy plants evolved over time
Roots and rhizosphere of healthy plants were preferentially colonized by Proteobacteria and Bacteroidetes bacterial phyla. Rhizosphere bacterial communities also contained Actinobacteria and Firmicutes but to a higher extent than in the roots. This result was expected because Proteobacteria, Bacteroidetes and Actinobacteria phyla were also highly abundant in the rhizosphere soil of many Brassicaceae species like A. thaliana [2] and B. napus [61][62][63], with the exception of Bacteroidetes being present at low frequencies in the rhizosphere of B. napus cultivated in a Podzol [61] and in a soil collected from an organically managed field [63]. Furthermore, higher frequencies of Firmicutes were observed in rhizosphere communities of Chinese cabbage than in other Brassicaceae species. Actinobacteria were detected at lower frequencies in the roots of Chinese cabbage than in the roots of A. thaliana [2,10] and B. napus [61,62,64]. As in the roots of all Brassicaceae, Cyanobacteria were also abundant in the root of Chinese cabbage. The impact of a pathogen on the plant root and rhizosphere microbiota

A fungus belonging to the Chytridiomycota phylum became dominant in the roots and rhizosphere of non-inoculated plants
The variations of fungal OTU frequencies in the communities of healthy plants were observed mainly at the two last sampling dates. The main changes in fungal communities concerned the relative abundances of an unknown Chytridiomycota which increased drastically at the end of the experiment in the two plant compartments, but especially in roots. This fungus replaced Ascomycota fungi that were previously dominant. In contrast to bacterial 16S sequences, fewer fungal 18S The impact of a pathogen on the plant root and rhizosphere microbiota sequences were available to use in taxonomic assignment. However, [65] described the fungal rhizosphere microbiota succession of B. rapa plants in compost over three plant generations by sequencing of ITS regions. From the second generation, the Chytridiomycota fungus assigned as Olpidium brassicae became dominant in the rhizosphere fungal communities. The organism in our samples could be O. brassicae or a close relative, but ultimately this would require confirmation by culturing or more extensive sequencing. This fungus is considered as a soilborne obligate parasite that invades Brassica rhizosphere, infects roots and reduces production of pods and seeds [66,67]. Its resting spores can remain dormant in the soil for up to 20 years before infecting roots. However, no symptom was observed on non-inoculated plant roots in our study indicating an asymptomatic interaction or a non-pathogenic fungus able to colonize root. For a non-mycorrhizal plant, we hypothesized that such endophytic fungus could play a role in plant nutrition as it was observed for phosphorus nutrition with Colletotrichum tofieldiae in natural population of A. thaliana [68] and with a fungus belonging to the Helotiales order in Arabis alpina [69].

Clubroot disease altered microbial community structure from the Chinese cabbage roots, then from its rhizosphere
To analyse how the soilborne pathogen affects bacterial and fungal root communities, Chinese cabbage seedlings were inoculated by P. brassicae resting spores ten days after sowing. Noninoculated and inoculated plants were cultivated in controlled conditions for 45 days after sowing. The bacterial and fungal metagenomes from the roots and rhizosphere of healthy and diseased plants were compared at several sampling dates after inoculation. We demonstrated that the invasion by a soilborne parasite changed root and rhizosphere microbial communities already assembled from the soil. Such results about the impact of a soilborne pathogen on the indigenous plant-associated microbiome was also described for Rhizoctonia solani on the lettuce microbiome [70] and for Ralstonia solanacearum on the tomato rhizosphere microbiota [71].
After inoculation, resting spores of P. brassicae released zoospores, which invaded the plant rhizosphere, reached the surface of the root hair and penetrated through the cell wall inside root hairs to form primary plasmodia. After nuclear divisions, the primary plasmodia differentiated into zoosporangia and secondary zoospores were formed in each zoosporangium to be released into the rhizosphere soil [34]. During this primary infection stage, P. brassicae was not detected in roots by quantitative PCR indicating that the amount of protist was very low. The interactions between the primary and secondary zoospores and the plant microbiota by direct or indirect mechanisms did not result in detectable changes in bacterial and fungal communities, neither in the roots, nor in the rhizosphere. After being released, the secondary zoospores penetrate the taproot cortical tissues. Inside invaded taproot cells, the pathogen develops into secondary plasmodia which are associated to cellular hypertrophy, followed by gall formation in the tissues [34]. This secondary infection stage was localized inside the roots. During this cortical infection, the amount of P. brassicae increased drastically and multiple direct interactions between the protist and the endosphere communities could occur.
We demonstrated that when P. brassicae developed inside the roots during its secondary infection stage, it strongly modified the endophytic bacterial communities and lightly the fungal communities. Then, probably as a consequence of the disturbances caused by the interactions between P. brassicae and the endophytic communities inside the roots, shifts in rhizosphere communities of diseased plants occurred only at the last date of sampling. In our study, several taxa (including Bacillus sp., Flavisolibacter sp. and Steptomyces sp.) showed their occurrence reduced in infested plants. These taxa have been associated with plant growth promotion or biocontrol functions such as Bacillus strain QST713 which is known to suppress clubroot on canola probably via antibiotics production and systemic resistance induction [72]. Changes in plant microbiota probably occurred by direct microbe-microbe interactions, mainly in the root compartment and then by direct microbial exchanging between the two compartments. However, changes in microbial composition following plant-parasite interactions are often hypothesized to be based on some modifications of the plant chemistry. Salicylic Acid (SA) and Jasmonic Acid (JA) are important hormonal regulators of the plant immune signalling network in which it is commonly accepted that SA is effective against biotrophic and JA against necrotrophic pathogens. However, P. brassicae is a biotrophic parasite and both SA and JA signalling pathways could play a role in partial inhibition of clubroot development in compatible interactions between A. thaliana and P. brassicae [73]. The defense-related phytohormones SA and JA are also known to important modulators of microbiota assembly of A. thaliana [74,75]. The accumulation of SA and/or JA or both in plant roots in response to P. brassicae infection could lead to modify the composition of plant rhizodeposits and to stimulate specific microbiota in the roots and rhizosphere.
These direct microbe-microbe or indirect microbe-plant interactions could drive the selection of a plant protective microbiome. In few situations, the competitive interaction between soilborne pathogens and root microbiota for available nutrients and microsites could lead to a strong restriction of the pathogen by the activities of specific microorganisms. These situations were already described in suppressive soils for soilborne [76] and foliar parasites [77]. A sequence of events taking place in the rhizosphere of sugar beet seedlings growing in a disease suppressive soil infected by R. solani was proposed as a model [76]. The fungus may induce, directly or indirectly via the plant, stress responses in the rhizosphere microbiome by the production of oxalic and phenylacetic acid and lead to shifts in microbiome composition by the activation of Oxalobacteraceae, Burkholderiaceae, Sphingobacteriaceae and Sphingomonadaceae families present in the suppressive rhizosphere microbiome. This stress in turn could trigger a response in these bacterial families, leading to the activation of antagonistic traits that restrict pathogen infection [74].
Several results have focussed on (i) the interplay between the plant immune system and root microbiota [78], (ii) the rule of root exudate metabolites as a driver of the rhizosphere microbiota assembling in relation to plant pathogen [26,79] and (iii) the modification of root plant metabolites during root fly attack in relation with modification of root microbiota diversity [27]. In light of these elements and our results, we propose a theoretical model to link modifications of microbiota diversity during Plasmodiophora brassicae infections (Fig 9). In this theoretical framework, P. brassicae during the first step of its life cycle crosses the plant rhizosphere and infect the root hair without inducing changes in microbiota composition as a consequence of plant metabolism modification. Then, the parasite penetrates inside the roots during the second step of its life cycle and induces due to gall growth strong modifications of the root microbiota. This invasion leads to modification of plant metabolites and root exudates, but also to induction of the plant immune system. As consequences of these trophic and defense modifications, we observed the selection in the root microbiota of specific microorganisms that could (i) use new metabolites, (ii) produce a signal triggering defense responses in plants and (iii) activate, directly or indirectly, other microorganisms in the root and rhizosphere microbiota to control the protist. Future studies will focus on investigating these hypotheses by, among others, selecting other soils and plant genotypes to promote the mechanisms that lead to the restriction of parasitic infection. We also want to develop more functional analysis of the plant-microbiota interactions to identify the underlying mechanisms.
The importance of microbiome for the functioning of plant has been widely recognized. Understanding the complex interactions between the pathogen or more generally biotic stress, the plant and its rhizosphere microbiome network are also key elements in shaping a plant-protective microbiome to improve the efficacies of biocontrol agents and partially resistant plants in controlling soilborne plant diseases. By this, plant microbiome is expected to have an important impact in biotechnology and will be a key point for the next Green Revolution as a harbinger to draw a new model for sustainable agriculture.  The impact of a pathogen on the plant root and rhizosphere microbiota (RS) and bulk soil (BS) samples at different sampling dates are represented. Fungal diversities were estimated with OTUs count data normalized by sample size and rarefied to 1,000 counts. Sampling date refers to 10 (T1), 24 (T3) and 45 (T5) days after sowing (DAS). For each sample, the number of replicates was n = 3. At each sampling date, lowercase letters indicate significant differences (p-values � 0.05) between conditions, which were assessed by ANOVA followed by post hoc Tukey's HSD test. Fungal diversity was estimated with OTUs count data normalized by sample size and rarefied to 1,000 counts. Richness and diversity associated to NI and I plants at each sampling date (T1 to T5) were compared. Sampling date refers to 0 (T1), 7 (T2), 14 (T3), 23 (T4) and 35 (T5) days after inoculation (DAI) with P. brassicae. At each sampling date, lowercase letters indicate significant differences (p-values � 0.05) between conditions, which were assessed by ANOVA followed by post hoc Tukey's HSD test. (TIF) S1 Table. Comparison of bacterial (B) and fungal (F) OTUs relative abundances in the roots (R) and rhizosphere (RS) of healthy and diseased Chinese cabbage plants. This table is organized into four tabs corresponding to the description of bacterial OTUs i) from the roots, ii) from the rhizosphere, fungal OTUs iii) from the roots and iv) from the rhizosphere. Mean values of abundance (expressed in %) were obtained from three replicates per condition and sampling date. Condition refers to non-inoculated (NI) and inoculated (I) plants. Sampling date refers to 0 (T1), 7 (T2), 14 (T3), 23 (T4) and 35 (T5) days after inoculation (DAI) by P. brassicae. OTUs with relative abundances below 1% were not shown. Significant differences (p-values � 0.05) of OTU frequencies between two samples are indicated by crosses. For example, T2 NI/I refers to the comparison of each OTU frequencies between NI and I plants at