- Research article
- Open Access
Transcriptomic responses of a simplified soil microcosm to a plant pathogen and its biocontrol agent reveal a complex reaction to harsh habitat
BMC Genomics volume 17, Article number: 838 (2016)
Soil microorganisms are key determinants of soil fertility and plant health. Soil phytopathogenic fungi are one of the most important causes of crop losses worldwide. Microbial biocontrol agents have been extensively studied as alternatives for controlling phytopathogenic soil microorganisms, but molecular interactions between them have mainly been characterised in dual cultures, without taking into account the soil microbial community. We used an RNA sequencing approach to elucidate the molecular interplay of a soil microbial community in response to a plant pathogen and its biocontrol agent, in order to examine the molecular patterns activated by the microorganisms.
A simplified soil microcosm containing 11 soil microorganisms was incubated with a plant root pathogen (Armillaria mellea) and its biocontrol agent (Trichoderma atroviride) for 24 h under controlled conditions. More than 46 million paired-end reads were obtained for each replicate and 28,309 differentially expressed genes were identified in total. Pathway analysis revealed complex adaptations of soil microorganisms to the harsh conditions of the soil matrix and to reciprocal microbial competition/cooperation relationships. Both the phytopathogen and its biocontrol agent were specifically recognised by the simplified soil microcosm: defence reaction mechanisms and neutral adaptation processes were activated in response to competitive (T. atroviride) or non-competitive (A. mellea) microorganisms, respectively. Moreover, activation of resistance mechanisms dominated in the simplified soil microcosm in the presence of both A. mellea and T. atroviride. Biocontrol processes of T. atroviride were already activated during incubation in the simplified soil microcosm, possibly to occupy niches in a competitive ecosystem, and they were not further enhanced by the introduction of A. mellea.
This work represents an additional step towards understanding molecular interactions between plant pathogens and biocontrol agents within a soil ecosystem. Global transcriptional analysis of the simplified soil microcosm revealed complex metabolic adaptation in the soil environment and specific responses to antagonistic or neutral intruders.
The soil microbial community is probably the most important component of the soil biota and it is responsible for a wide range of ecologically significant functions . Soil microorganisms are key determinants of biological soil quality, supporting soil fertility and indirectly plant health . Soil bacteria and fungi have undergone evolution and niche differentiation, interacting in different ways with various outcomes in soil ecosystems . Microbe-microbe interactions have an important impact on microbial fitness in the soil  and they can positively or negatively affect other participants, including organisms of higher taxa . These interactions can be synergistic, neutral or antagonistic, and can affect the growth, metabolism and differentiation of members of the microbial community . For example, direct and indirect competition, like antagonism and nutrient consumption, have an adverse effect on interacting members of the population . Plant pathogens have to compete with members of the rhizosphere microbiota for available nutrients and microsites in order to infect root tissue, and they are significantly restricted in growth by the antagonistic activities of biocontrol microorganisms in suppressive soils . However, microbes have evolved strategies not only to fight each other, but in some cases relationships to adapt or support each other, increasing the overall fitness of the community . Despite the widespread occurrence of interspecific microbial interactions in nature and their crucial relevance, there is a lack of understanding about how each species perceives and responds to interactions with other microbes within a complex community .
Molecular interplay between soil microorganisms has mainly been studied in dual cultures of plant pathogens and biocontrol agents [2, 7–10]. Soil phytopathogenic fungi cause severe diseases in several crops, and Armillaria spp. are among the most serious plant root pathogens, and cause significant economic losses to perennials, such as fruit crops, timber and ornamental trees in boreal, temperate and tropical regions of the world . In particular, one of the most aggressive species (Armillaria mellea) colonises living woody roots, causing damping off and root rot . Infected plants show a decline in vigour and production and generally die within a few years of infection . Armillaria mellea can persist saprotrophically in roots after the death of the host and it can also spread below ground through specialised and robust structures known as rhizomorphs . Fumigants and fungicides cannot completely eradicate the Armillaria spp. inoculum in soil, and post-infection control is commonly unsuccessful, as in the case of several soil phytopathogens . Therefore, biocontrol microorganisms have been increasingly screened for antagonism against this soil phytopathogenic fungus. Trichoderma spp. are one of the most widely applied biocontrol fungi and the Trichoderma atroviride SC1 strain has been identified as an effective biocontrol agent of A. mellea . The biocontrol properties of this strain have been studied in vitro, demonstrating that its mechanism of action relies mainly on hyperparasitization of A. mellea hyphae . The biocontrol properties of Trichoderma spp. have been extensively studied using several approaches, including transcriptomics . Most transcriptomic studies of microbial biocontrol mechanisms have focused on molecular responses of the biocontrol agent and the target plant pathogen in dual cultures [2, 7, 8, 10, 18]. However, microbial properties in complex consortia often differ from those of their single components . In particular, biocontrol agents should compete with indigenous communities to occupy niches and show appropriate biocontrol properties against the target pathogens . Introduction of T. atroviride changes the microbial structure of the soil microflora during the first two weeks following inoculation , suggesting strong interaction between this biocontrol strain and indigenous microbial communities. The objective of this study was to clarify the transcriptomic response of a simplified soil microcosm composed of 11 soil microorganisms to the introduction of a plant pathogen (A. mellea) and its biocontrol agent (T. atroviride) in a soil matrix. Greater insight into these mechanisms can help to understand the dynamics underlying interactions among soil microorganisms following the application of biocontrol agents against soil phytopathogenic fungi, and the impact of biocontrol agents on indigenous non-target microorganisms.
Growth conditions for soil microorganisms
Azorhizobium caulinodans ORS 571 (DSM 5975, ), Bacillus subtilis subsp. subtilis 168 (ATCC 23857, ), Cupriavidus metallidurans CH34 (DSM 2839, ) and Pseudomonas protegens Pf-5 (ATCC BAA477, ) were grown in 20 ml of Luria-Bertani broth (LB; Oxoid) at 25 °C overnight, with orbital shaking at 250 rpm. To prepare the inoculum suspension of each bacterial strain, 15 ml of the culture grown overnight were centrifuged at 5000 × g for 10 min at room temperature. The supernatant was discarded and the bacterial cells were washed three times by suspending in 20 ml of sterile isotonic solution (0.9 % NaCl) and centrifuging at 5000 × g for 10 min. The concentration of each bacterial suspension was assessed by measuring the optical density (OD) at 600 nm with a spectrophotometer (Ultrospec 3100, GE Healthcare Life Sciences) and adjusted to 5 × 108 cells ml−1.
Debaryomyces hansenii CBS 767 , Pichia stipitis CBS 6054 , Schizosaccharomyces pombe 972 h (DSM 70576, ), Saccharomyces cerevisiae S288c (ATCC 204508, ) and T. atroviride SC1 (available in our strain collection and in the commercial product Vintec, Bi-PA, Londerzeel, Belgium) were grown on malt extract agar (MEA; Oxoid) at 25 °C. Aspergillus niger CBS 513.88 , Fusarium oxysporum f. sp. lycopersici 4287 (CBS 123668, ) and Penicillium chrysogenum Wisconsin 54–1255 (DSM 1075, ), were grown on potato dextrose agar (PDA; Oxoid) at 25 °C. To prepare the inoculum suspension of filamentous fungi and yeasts, conidia and cells were collected from 21-day-old cultures by washing each plate with 3 ml of sterile isotonic solution, using a glass rod under sterile conditions. The suspension of each strain was filtered with sterile cloth, and cells were washed three times by suspending in 6 ml of sterile isotonic solution and centrifuging at 5000 × g for 10 min. The concentration of each conidia and cell suspension was adjusted to 5 × 108 cells ml−1 by counting with a hemocytometer under a light microscope. For the inoculum suspension of A. mellea M6132 (kindly provided by Dr. Simone Prospero, Swiss Federal Research Institute, Birmensdorf), 1 g of mycelium was collected from 21-day-old culture grown on MEA at 25 °C and suspended in 4 ml of sterile isotonic solution. Sterile steel beads were added and the A. mellea mycelium was ground in a mixer-mill disruptor (MM 400, Retsch) at 25 Hz for 3 min. The ground mycelium was washed three times by suspending in 4 ml of sterile isotonic solution and centrifuging at 5000 × g for 10 min, and it was then suspended in 2 ml of sterile isotonic solution.
The viability of counted cells in each inoculum suspension was validated through a dilution plating method (Additional file 1). Basically, 0.1 ml of each inoculum suspension was subjected to 10-fold serial dilution, 0.2 ml of each dilution was plated in triplicate on plates containing the appropriate media for each strain, bacterial and fungal colony forming units (CFUs) were then assessed after incubation at 25 °C for 24 and 48 h respectively.
Set-up of the simplified soil microcosm and sample collection
The soil matrix was obtained as described by Ellis , by mixing 35 g of sand (Sigma-Aldrich), 10 g of kaolinite (Sigma-Aldrich), 5 of bentonite (Sigma-Aldrich) and 0.1 g of CaCO3 (Sigma-Aldrich) in a 200 mL glass vessel for plant tissue culture (Sigma-Aldrich). The soil matrix was sterilised at 121 °C for 30 min, dried in an oven at 80 °C for 24 h and cooled down to room temperature. Subsequently, 1 g of humic acid (Sigma-Aldrich) was added under sterile conditions and the soil matrix was mixed by vigorous shaking.
To prepare the simplified soil microcosm, the soil matrix was inoculated with a suspension obtained by mixing 1 ml of the suspension of each microbial strain to a final volume of 15 ml in sterile isotonic solution, and it was incubated at 25 °C with orbital shaking at 150 rpm. Five conditions were analysed in triplicate: the simplified soil microcosm incubated for 24 h either without exogenous fungi (SSM), with the biocontrol agent T. atroviride (SSM+T), with the plant pathogen A. mellea (SSM+A) or with both (SSM+T+A), and the simplified soil microcosm collected at the beginning of the experiment just after inoculation of the 11 soil microorganisms (SSM0). For each condition, 20 aliquots of 2 g of soil were collected and placed into 15 ml sterile tubes, immediately frozen in liquid N2 and stored at −80 °C. In an preliminary experiment the viability of each microbial strain after 24 h of incubation in the soil matrix was validated (Additional file 1). Specifically, each microbial strain was incubated in the soil matrix at 25 °C, after 24 h the CFUs were assessed by a dilution plating method on plates containing the above mentioned media and no significant changes (t-test, P > 0.05) on microbial CFUs were found for each soil microorganism before and after incubation in the soil matrix.
RNA isolation and mRNA enrichment
Total RNA was extracted from 2 g of soil using the RNA PowerSoil Total RNA Isolation Kit (MoBio Laboratories) with slight modifications. Briefly, beads were directly added to the frozen soil sample, and the bead solution was added simultaneously with Solution SR1. After vortexing, the sample was centrifuged at 2500 × g for 20 min at room temperature and the following steps of the extraction protocol were carried out, according to the manufacturer’s instructions. Total RNA was quantified using a Nanodrop 8000 (Thermo Fisher Scientific) and its quality was checked using a Bioanalyzer 2100 (Agilent Technologies). In a preliminary experiment, each soil microorganism was separately inoculated in the soil matrix and the equivalent quality and yield of total RNA were assessed 24 h after incubation at 25 °C.
For each of the five conditions, twelve independent extractions of total RNA were obtained from the simplified soil microcosm, and three replicates deriving from the pool of four RNA extractions were used in further analysis. Total RNA was purified and concentrated using the RNeasy Plus Micro kit (Qiagen), which allowed depletion of genomic DNA, and 3 μg total RNA was subjected to mRNA enrichment. The RiboMinus Eukaryote Kit for RNA-Seq (Invitrogen, Life Technologies) was used, and probes of the Ribominus Transcriptome Isolation Kit (Yeast and Bacteria) (Invitrogen, Life Technologies) were added in the hybridization step, in order to simultaneously remove the ribosomal RNA of fungi, bacteria and yeasts. Purified mRNA was quantified using a Nanodrop 8000 (Thermo Fisher Scientific) and its quality was checked using a Bioanalyzer 2100 (Agilent Technologies) to verify the absence of ribosomal RNA.
Library construction and Illumina sequencing
Three replicates for each condition were subjected to RNA-Seq library construction, using the TruSeq SBS v3 protocol (Illumina), and paired-end reads of 100 nucleotides were obtained using an Illumina HiSeq 2000 at Fasteris (Plan-les-Ouates, Switzerland). Briefly, 200 ng of the purified mRNA was fragmented using Zn-catalysed hydrolysis and converted into double-stranded cDNA with random priming. Following end repair, indexed adapters were ligated and cDNA fragments of 200 ± 50 bp were purified. Purified cDNA was amplified by PCR and validated by Illumina sequencing, after which mRNA-Seq libraries were multiplexed (four libraries per lane) and sequenced according to the manufacturer’s instructions.
Mapping of sequenced reads and assessment of gene expression
Quality trimming and sequence filtering were carried out using Trimmomatic (version 0.32)  to remove sequencing adapters, low-quality bases at the end of the reads (Phred quality score lower than 15 in a sliding window of 4 bases) and short reads (shorter than 50 bases or with an average Phred quality score below 25). The FastQC (version 0.10.1) was used for sequence quality control .
The reference genomes of the 13 microorganisms were downloaded from public databases (Table 1), unique identifiers of contigs and genes were generated by adding a two letter code for each microorganism and the microcosm genome was built by merging the genomic data. Filtered read pairs were mapped to the microcosm genome using TopHat2 (version 2.0.10)  with the default settings, except for the inner distance between mate pairs of 200 ± 50, maximum intron length of 4000, and minimum intron length of 10 . The results of the alignments were analysed with Samtools (version 0.1.19)  and visualised using the Integrative Genomics Viewer . Unambiguously mapped read pairs were used to improve the specificity of gene expression estimation, and read pairs mapped to gene regions were counted using HTSeq . As for transcriptomic studies of microbial consortia [33, 34], global normalisation of read mapping to the microcosm genome was used and gene expression levels were calculated as fragments per kilobase of each gene per million of unambiguously mapped read pairs (FPKM; ) of the microcosm genome. Using this approach, changes in gene expression levels are related to transcriptional regulation and the biomass abundance of each microorganism, compared with other members of the microcosm .
Identification of differentially expressed genes
Differentially expressed genes (DEGs) were identified with the DESeq2 package . Briefly, counts of unambiguously mapped reads were normalised for sequencing depth following negative-binomial distribution and an empirical Bayes shrinkage approach was used for estimation of data dispersion and fold changes . A false discovery rate (FDR) of Benjamini-Hochberg multiple tests lower than 5 % (p < 0.05) with a minimum Log2 fold change of one and a minimum expression level of 0.6 FPKM (corresponding to a coverage of 1× of the gene sequence in our samples) in at least one condition were imposed to identify DEGs through pairwise comparison. Four pairwise comparisons of the simplified soil microcosm were analysed: SSM vs. SSM0, SSM+A vs. SSM0, SSM+T vs. SSM0, and SSM+T+A vs. SSM0. DEGs were then grouped according to their expression profiles and clusters were visualised using MultiExperiment Viewer software (version 4.9) . Two additional pairwise comparisons (SSM+A vs. SSM+T+A and SSM+T vs. SSM+T+A) were considered to identify the DEGs of A. mellea and T. atroviride. Fold changes of FPKM values were calculated, with the minimum expression value of 0.002 FPKM imposed on non-expressed genes.
Functional annotation of differentially expressed genes
DEGs were annotated using the ARGOT2 function prediction tool , and grouped into 17 functional categories of Gene Ontology (GO) biological process terms , based on the BlastX results against the UniProtKB database . Genes not associated with any biological process category were assigned to the GO root (GO:0008152, namely biological process). Gene descriptions were obtained from the annotation files available in public databases for the 13 soil microorganisms (Table 1). To improve gene descriptions, sequences of DEGs were aligned using a BlastX search against the Swiss Prot database , and the best hits (E-value lower than 1E−5) were manually checked.
GO terms significantly overrepresented in the DEGs of the simplified soil microcosm were identified using the Biological Networks Gene Ontology (BiNGO) tool , and biological networks were visualised with Cytoscape version 3.2.1 . Basically, sequences of DEGs were aligned against the protein database of S. cerevisiae S288c [44, 45] using a BlastX search and homologous S. cerevisiae proteins were identified (E-value lower than 1E−5). Redundancies were removed and significantly overrepresented (P < 0.05) GO biological process terms were identified for the up- and down-regulated genes of each cluster in comparison with the whole S. cerevisiae annotation as a reference set.
DEGs were also assigned to KEGG (Kyoto Encyclopedia of Genes and Genomes) ortholog groups, redundancies were removed and they were mapped to KEGG pathways using the KEGG Automatic Annotation Server (KAAS) [46, 47] with the BBH method and default gene data set [48, 49].
Gene expression analysis by quantitative real-time RT-PCR
First-strand cDNA was synthesized from 600 ng of purified RNA using SuperScript VILO with random hexamers (Invitrogen, Life Technologies). cDNA was tested for contamination by genomic DNA through amplification of a translocon alpha subunit Sec61 of S. pombe with intron-spanning primers (Additional file 2). PCR reaction was carried out with a T-Professional Thermocycler (Biometra) using the DreamTaq DNA Polymerase (Fermentas, Life Technologies), and amplification of the expected cDNA fragment of Sec61 was verified with gel electrophoresis.
Quantitative real-time PCR (RT-qPCR) was carried out on ten times-diluted cDNA with Platinum SYBR Green qPCR SuperMix-UDG (Invitrogen, Life Technologies) and specific primers (Additional file 2), using a LightCycler 480 instrument (Roche, Branford, CT, USA). Reactions were set with two initial steps at 55 °C for 10 min and 95 °C for 2 min, followed by 50 cycles at 95 °C for 15 s and 60 °C for 1 min. Each sample was examined in three technical replicates and dissociation curves were analysed to verify the specificity of each amplification reaction. To validate the amplification of the target genes, real-time PCR products were treated with Exosap-IT (Usb Corporation) and sequenced on both strands using an AB3730xl instrument (Applied Biosystems) at the Sequencing Platform facility (Fondazione Edmund Mach). Sequences were examined using FinchTV Version 1.4.0 (Geospiza Inc.) and aligned with the expected target gene using T-Coffee Multiple Sequence Alignment . The specificity of designed primers and amplified sequences was further validated with a BlastN search against all genes of the microcosm composed by 13 soil microorganisms, and specific alignments to each target gene were verified.
Cycle threshold Ct) values were extracted with LightCycler 480 SV1.5.0 software (Roche) using the second derivative calculation. For each gene, the reaction efficiency (Eff) was calculated with the LinRegPCR 11.1 software  and used to calculate relative quantities (RQ), according to the formula: RQ = Eff(Ct–Ct′), where Ct is the threshold cycle and Ct’ is the average Ct of all the conditions analysed . Four housekeeping genes were selected, as their expression was not significantly affected in the five conditions (Additional file 2), and their stability was validated using the ∆Ct method described by Silver et al. . Normalised relative quantities (NRQ) were then calculated by dividing the RQ by a normalisation factor, based on the RQ value of the two most stable reference genes : Actin 1 of S. pombe and Elongation factor 1α of F. oxysporum, for global normalisation of the microcosm transcriptome.
Experimental design and RNA-Seq analysis of the simplified soil microcosm
The simplified soil microcosm was based on a soil matrix simultaneously inoculated with 11 soil microorganisms whose genomes were already sequenced and available that were selected to mimic a community of bacteria (Gram-positive and negative), yeasts (fission and budding yeasts) and filamentous fungi commonly found in soil [44, 54–66], having beneficial (biocontrol agent or nitrogen fixing microorganisms), pathogenic, saprophytic or neutral properties for crops (Table 1). The plant pathogen (A. mellea) and its biocontrol agent (T. atroviride) were introduced into the simplified soil microcosm and samples were collected 24 h after incubation under controlled conditions, to analyse the early transcriptional response of the microbial community to these exogenous fungi. Five conditions were analysed in triplicate: the simplified soil microcosm incubated for 24 h either without exogenous fungi (SSM), with the biocontrol agent T. atroviride (SSM+T), with the plant pathogen A. mellea (SSM+A) or with both (SSM+T+A), and the simplified soil microcosm collected at the beginning of the experiment just after inoculation of the 11 soil microorganisms (SSM0). From 46 to 80 million paired-end reads of 100 nucleotides were obtained for each replicate and more than 84 % of them passed the quality control test (filtered read pairs; Additional file 3). More than 96 % of filtered read pairs were aligned (mapped read pairs) to the microcosm genome obtained by merging the genomes of the 11 soil microorganisms, the plant pathogen and the biocontrol agent. From 18 to 27 % of mapped read pairs were unambiguously aligned to unique locations (unique read pairs, Additional file 3) and reads mapping to multiple locations were not considered for gene expression analysis, as they possibly matched orthologous genes with high sequence similarity among the 13 soil microorganisms analysed.
The majority (more than 77 %) of unique read pairs were mapped to gene regions of the microcosm genome (Additional file 3) and the proportion of read pairs mapping to genes revealed the different transcriptional activity of soil microorganisms (Additional file 4A–D). Specifically, high (more than 10 %) and low (less than 0.1 %) proportions of read pairs were mapped to the P. protegens and B. subtilis genes respectively, in agreement with the high and low number of CFUs inoculated in the soil matrix (Additional file 1). However, high proportions of expressed genes (more than 76 %) were obtained for underrepresented transcriptomes (Additional file 4E, F and Additional file 5), such as those of A. caulinodans and D. hansenii (less than 2 % of unique pairs mapping to genes; Additional file 4D). Proportion of expressed genes of each soil microorganism was consistent among the five conditions, indicating that transcriptional activity and viability were not affected by incubation in the artificial soil in the presence or absence of A. mellea and T. atroviride. As expected, read pairs of A. mellea and T. atroviride were found either in the SSM+A or SSM+T and in the SSM+A+T condition, confirming the specificity of the alignment results (Additional file 4B, D and F).
Identification and functional annotation of differentially expressed genes
Major transcriptional changes occurred during incubation of the soil microorganisms in the soil matrix, independently of A. mellea and T. atroviride introduction (Additional file 6). Pearson’s correlation coefficients of gene expression levels among the simplified soil microcosm collected at the beginning of the experiment (SSM0) and the four conditions at 24 h ranged from 0.85 to 0.87, and they were lower than those found among the four conditions at 24 h (ranging from 0.98 to 0.99). A total of 28,309 DEGs were found using DESeq statistical analysis, with a FDR of 5 %, a minimum fold-change of two and an expression level greater than 0.6 FPKM (Additional file 7). About 40 % of DEGs were expressed in all conditions, and a total of 3023 DEGs were specifically modulated by A. mellea and T. atroviride introduction (SSM+A+T; Fig. 1a). Moreover, eight and 88 DEGs of A. mellea and T. atroviride were identified in comparison of SSM+A and SSM+T vs. SSM+T+A respectively, and they represent transcriptional changes during pathogen-biocontrol interaction within the simplified soil microcosm. The proportion of DEGs ranged from 12 to 82 % of the expressed genes for each strain and it was not correlated with the proportion of expressed genes (Pearson’s correlation coefficient of 0.53, P = 0.06, Additional file 8). For example, C. metallidurans showed high transcriptional activity (91 % of expressed genes) and reduced transcriptional changes (38 % of DEGs), while P. protegens and P. stipitis were both highly transcriptionally active (more than 98 % of expressed genes) and responsive to the tested conditions (more than 64 % of DEGs).
Functional annotation of DEGs revealed regulation of microbial processes related to gene expression, nucleobase-containing compounds (DNA and RNA metabolism) and macromolecule metabolic processes during incubation in the soil matrix (Fig. 1b). Transcriptional regulations of the simplified soil microcosm included primary metabolism processes (carbohydrate, protein and lipid metabolic process), transport and cellular responses (response to stimulus and signal transduction), and they were stimulated by adaptation to the soil matrix environment and by A. mellea or T. atroviride introduction. However, a large fraction of DEGs (17 and 12 % respectively) were associated with poorly defined processes (single-organism cellular process) or were not associated with known functions (assigned to the Gene Ontology root, biological process).
Expression profiles of the differentially expressed genes
The DEGs of the simplified soil microcosm were grouped into 16 clusters based on the expression profiles (clusters 1–16; Fig. 2), and two additional clusters included the DEGs of T. atroviride (cluster 17) and A. mellea (cluster 18; Additional file 9). The majority (53 %) of DEGs of the simplified soil microcosm were modulated by incubation in the soil matrix, with similar expression profiles in the presence or absence of A. mellea and T. atroviride, alone and combined (cluster 1), and reflected adaptation processes to the soil matrix environment not influenced by the presence of the plant pathogen and the biocontrol agent. A large fraction (17 and 4 % respectively) of DEGs were modulated by incubation in the soil matrix, with similar expression profiles in the presence or absence of A. mellea and T. atroviride, with reinforced modulation (cluster 15) or no modulation (cluster 16) in the presence of A. mellea and T. atroviride combined. Of the DEGs of the simplified soil microcosm, 322 and 327 genes were modulated by T. atroviride alone (cluster 2) or combined with A. mellea (cluster 3), respectively. Moreover, 346 and 1202 genes were modulated by A. mellea alone (cluster 4) or combined with T. atroviride (cluster 5), respectively. Cluster 6 comprised 1407 genes whose expression was affected by the presence of A. mellea and T. atroviride combined. The expression of 1993 genes was modulated by the presence of A. mellea and T. atroviride alone and combined (cluster 7), while that of 478 genes was affected by the presence of A. mellea or T. atroviride alone (cluster 8). In addition, six minor clusters comprised DEGs affected by 24 h incubation in the soil matrix, with specific profiles in the presence of A. mellea and/or T. atroviride (clusters from 9 to 14).
To validate the RNA-Seq results, the expression levels of 15 DEGs were analysed by real-time RT-PCR (Fig. 3a–q). Genes with different expression profiles belonging to different soil microorganisms were analysed, including genes associated with different functional categories (Additional file 2). Although the extent of modulation revealed by real-time RT-PCR and RNA-Seq may differ [67, 68], the expression profiles obtained using the two methods agreed totally for ten genes. The expression profiles generated by real-time RT-PCR and RNA-Seq differed in one or all of the conditions for three (Fig. 3e, g and i) and two (Fig. 3d and c) genes, respectively. These differences could be due to differences in sensitivity, particularly in distinguishing members of multigene families and orthologous genes.
Metabolic pathways differentially affected by Armillaria mellea and Trichoderma atroviride
The metabolic pathways of the simplified soil microcosm were affected overall by 24 h incubation in the soil matrix, with similar expression profiles in the presence or absence of the plant pathogen and biocontrol agent (clusters 1 and 15; Additional file 10A and B). Complex reprogramming of primary metabolic processes (organic acid, carbohydrate and lipid metabolism) indicated adaptation of the microbial communities, from the nutrient conditions of the rich microbiological media to the harsh habitat of the soil matrix (Additional file 11). For example, genes involved in the metabolism of aromatic compounds (quinate, vanillate, ferulate, and benzene metabolism) and in energy production under anaerobic conditions (formate and acetate metabolism) were upregulated by incubation in the soil matrix as possible adaptation of the simplified soil microcosm to humic acid metabolism under semi-anaerobic conditions. Likewise, growth regulators and conidiation-related genes were modulated by incubation of the simplified soil microcosm as a consequence of adaptation processes from in vitro to soil conditions (Additional file 7).
The introduction of T. atroviride caused modulation of genes related to secondary metabolic pathways (metabolism of folate and toxic compounds) and regulators of cell growth and differentiation (controllers of meiotic processes and regulators of hyphal development) in the simplified soil microcosm (cluster 3; Additional file 10C). Oxidation-reduction processes (Fig. 4a), genes involved in signal transduction (kinases and phosphatases) and defence response (multidrug- and heavy metal- efflux proteins) were activated by the simplified soil microcosm as a reaction to T. atroviride (Additional file 7). Conversely, nucleic acid and RNA metabolic processes were mainly activated by the simplified soil microcosm during incubation with A. mellea (cluster 5), as well as processes related to response to stimuli (receptors and protein kinases) and to oxidative stresses (catalases and oxidoreductases; Fig. 4b, Additional files 7 and 10D).
The presence of both A. mellea and T. atroviride (cluster 7) upregulated expression regulators (transcriptional regulators NRG1 and NahR, a transcriptional activator of hydrogenase-4), genes implicated in defence processes (drug resistance proteins, protein KRE1, and ABC transporters for drug resistance) and secondary metabolism (Additional file 10E). Metabolisms of amino acids and amines were upregulated, as well as the regulation processes of signal transduction and protein localization (Fig. 4c). Conversely, genes related to regulation of cell division (spindle pole body-associated proteins, telomere length regulators, actin and myosin genes), cell cycle (mitotic cyclins, meiotic activators and cell division control proteins) and microbial growth (negative regulator of sporulation PMD1, morphogenesis-related gene MSB1, genes WHI4 and DSE2) were downregulated in the simplified soil microcosm (Additional file 7), as a possible consequence of competition with A. mellea and T. atroviride.
Molecular interactions between biocontrol agents and target pathogens have mainly been studied in dual culture assays [2, 7, 8, 10, 18] and little is known about the more intricate molecular interplay of a complex microbial community. Microorganisms can sense and respond to the presence of neighbouring strains in the environment and alter their behaviour and performance accordingly . In this study, we took advantage of RNA-Seq technology to study the interaction of a biocontrol agent (T. atroviride) with a ubiquitous plant pathogen (A. mellea) in a simplified soil microcosm containing 11 microorganisms. Complex microbial responses were found in the simplified soil microcosm. In particular, major transcriptional regulation occurred when single colonies grown on rich media in Petri dishes were transferred to the soil matrix, reflecting microbial acclimatisation to a different environment and the presence of other members of the microbial community. However, specific molecular responses were activated by the simplified soil microcosm in response to the introduction of the plant pathogen or the biocontrol agent, upregulating neutral adaptation processes and active defence mechanisms, respectively (the overview of main processes is reported in Fig. 5 and key DEGs are summarised in Additional file 12).
The results confirmed that RNA-Seq was an appropriate method for specifically estimating gene expression levels within a complex microbial consortia [33, 34]. However, the use of unambiguously mapped read pairs for gene expression quantification maximises specificity in the sequence count assignment, but could underestimate the expression levels of homologous genes within the microbial community. Moreover, gene expression levels in metatranscriptomic analysis are determined not only by transcriptional regulation, but also by the biomass of each microorganism, compared with other members of the community . For example, high and low proportions of P. protegens and B. subtilis reads were obtained, respectively, in agreement with their respective biomass in the soil matrix. S. pombe had the same concentration as the mean of other soil microorganisms and accounted for the greatest proportion of sequenced read pairs, resulting the most transcriptionally active microorganisms of the simplified soil microcosm.
Metabolic adaptation of the simplified soil microcosm in the soil matrix was not modified by the introduction of Armillaria mellea and Trichoderma atroviride
Major transcriptional changes of genes related to complex metabolic pathways occurred through incubation in the soil matrix and showed similar expression profiles in the presence or absence of A. mellea and T. atroviride (cluster 1; Additional file 12), suggesting more significant adaptation to the soil environment than to the microbial intruders. Metabolic shifts from the nutrient conditions of the rich microbiological media to the soil matrix were highlighted by the transcriptional activation of genes implicated in the metabolism of aromatic compounds, formate and acetate (Fig. 5) possibly involved in the catabolism of humic acid  and energy production under anaerobic conditions .
The expression of defence genes was activated by the soil microorganisms to compete with each other, such as six hydrolases of A. niger, C. metallidurans and F. oxysporum, responsible for the detoxification of epoxides , and three polyketide synthases (PKS) of A. niger, P. stipitis and D. hansenii, responsible for the biosynthesis of toxic molecules [5, 7]. The upregulation of genes encoding 23 flagellar proteins, two flagellar regulators, three chemotaxis proteins and two pilus assembly proteins suggested the activation of bacterial cell motility, which could contribute to the initial phases of colonisation in the soil , as well as for exploreing soil niches and getting away from competitors . Based on these results, metabolic adaptation and niche colonisation processes seem to be based on the metabolic plasticity and competition/defence capacity of each microorganism against neighbouring competitors in the soil matrix and these transcriptional changes dramatically increase the complexity of microbial community reprogramming.
The simplified soil microcosm responded to the introduction of Trichoderma atroviride by activating defence-related processes
At transcriptional level, the overall response of the simplified soil microcosm to the introduction of T. atroviride (328 DEGs of cluster 3) revealed a reaction to an antagonistic intruder, such as the activation of signal transduction processes, metabolism of toxic compounds, oxidation-reduction and defence processes, as well as shifts in cell growth and differentiation (Fig. 5, Additional file 12). Specifically, the activation of genes implicated in signal transduction processes (five protein kinases, a protein phosphatase and a calcium-binding protein) and transcriptional regulation (the transcriptional regulatory proteins HoxA and PHO23, a RNA polymerase and a zinc finger domain protein) indicated recognition and reaction of the simplified soil microcosm to the biocontrol agent. The microbial response included the upregulation of genes encoding enzymes implicated in the metabolism of toxic compounds, such as a nonribosomal peptide synthetase (NRPS) of A. niger, a methenyltetrahydrofolate cyclohydrolase and a riboflavin biosynthesis protein of A. caulinodans, a formylbenzenesulfonate dehydrogenase and a phytoene synthase of C. metallidurans, and an isovaleryl dehydrogenase of F. oxysporum. Furthermore, other defence-related genes were upregulated, possibly to protect the soil microbial community against the biocontrol agent, such as three resistance proteins and a gallate transporter of C. metallidurans, a peptidase of P. protegens and a protease of P. stipitis. Proteases and enzymes responsible for the detoxification of toxic molecules have already been found to be transcribed in the plant pathogenic fungi Rhizoctonia solani and A. niger in response to antagonistic strains belonging to Bacillus and Serratia genera [2, 5], suggesting that biocontrol agents activate typical defence processes in neighbour microorganisms. Oxidoreductases of A. niger and F. oxysporum were upregulated in the simplified soil microcosm by T. atroviride to remove active oxidants, as already observed for R. solani in response to Serratia spp. . Conversely, stress-related genes were downregulated in A. niger in response to B. subtilis attachment , suggesting that microbial responses could differ in mutualistic and antagonistic interactions.
Biocontrol agents actively compete with indigenous microorganisms through the production of antimicrobial compounds  and T. atroviride can antagonize and/or parasitize target fungi , suggesting that this fungal strain may have inhibited the growth of the soil microbial community. This hypothesis is in line with the inhibition of growth and differentiation processes of the simplified soil microcosm by T. atroviride introduction that was suggested by the downregulation of six regulators of cell division and four genes implicated in hyphal growth of A. niger, D. hansenii, F. oxysporum, S. cerevisiae and S. pombe. Likewise, a regulator of sporulation and two genes related to cell wall biogenesis were downregulated in D. hansenii and S. cerevisiae. Transcriptional changes of the growth- and morphology-related genes of plant pathogens in response to antagonistic bacteria are known to be linked to stress responses , and they could be proposed as markers of competitive interactions among the microbiota. Two flagellar transcriptional regulators (Rmet_1690 and Rmet_1691) were downregulated in response to T. atroviride introduction, and chemotactic motility genes were modulated in Collimonas fungivorans and Serratia plymuthica incubated with A. niger  and R. solani , respectively. Taken together, these expression profiles reflected the reaction of the simplified soil microcosm to the competition displayed by the biocontrol agent T. atroviride through antagonism and production of toxic molecules.
The simplified soil microcosm responded to the introduction of Armillaria mellea by activating neutral adaptation processes
The transcriptional response of the simplified soil microcosm to A. mellea introduction (1203 DEGs of cluster 5) activated more neutral adaptation processes as compared to T. atroviride introduction (Fig. 5, Additional file 12). Six protein kinases, two receptors and 21 transcription regulators were upregulated in the simplified soil microcosm and are candidate activators of the microcosm response to A. mellea. Similar expression changes were found for 16, 14 and 13 genes responsible for RNA metabolism, ribosome biogenesis and regulation of oxidative stress, respectively. Signal transduction  and oxidation-reduction processes  are commonly activated during microbe-microbe communications and expression profiles reflect the response of the simplified soil microcosm to A. mellea. However, genes implicated in defence-related processes were mainly down- and up-regulated in response to A. mellea and to T. atroviride, respectively. These opposite regulations indicate that the simplified soil microcosm may recognise the antagonistic properties of T. atroviride and react with early activation of defence processes, which are not required in the presence of the plant pathogen. Specifically, the A. mellea introduction downregulated the expression of genes responsible for cell wall reinforcement, stress-related processes and the metabolism of toxic and aromatic compounds in A. niger, A. caulinodans, C. metallidurans, P. chrysogenum, S. cerevisiae and S. pombe. In nature microbes have evolved mechanisms not only to fight, but also to adapt, and in some cases to support each other . The microbial transcriptional response differs strongly during antagonistic, neutral and beneficial interactions  and defence-related processes are generally downregulated during mutualistic interactions . Therefore, the simplified soil microcosm was probably able to recognise A. mellea as a non-competitive intruder and it did not stimulate early antagonistic reactions, but rather neutral adaptation processes to this exogenous fungus. However, microbial communities of suppressive soils could activate antagonistic traits against soil phytopathogens after long incubation periods .
The response mechanisms to Trichoderma atroviride dominated over those to Armillaria mellea following simultaneous introduction to the simplified soil microcosm
When A. mellea and T. atroviride were simultaneously introduced to the simplified soil microcosm (cluster 7), the activation of defence-related genes indicated an attempted defence response by the microbial community, mediated by the upregulation of genes implicated in cell wall reinforcement, lytic processes (proteases and peptidases) and drug resistance (Fig. 5, Additional file 12). Secondary metabolic processes play a crucial role in the biocontrol activities of T. atroviride , and genes for the metabolism of aromatic compounds and biofilm formation were upregulated, possibly to defend the microbial community against this fungal intruder. Likewise, two ABC transporter genes were upregulated in A. caulinodans and C. metallidurans, and their activation has been associated with detoxification pathways during biocontrol processes . Increased protein synthesis during biocontrol interactions was linked to the defence mechanisms of the prey microorganisms , and genes implicated in the metabolism of amino acids were upregulated by the simplified soil microcosm. However, microbial distress was indicated by the downregulation of growth-related genes, such as 32 regulators of division processes, 31 regulators of the cell cycle, ten genes implicated in DNA replication processes and nine in chromatin remodelling processes.
In addition to defence processes, signal transduction and transcription regulations were also activated. In particular, genes encoding protein kinase, protein phosphatases and signalling proteins were activated in C. metallidurans, D. hansenii, F. oxysporum, P. chrysogenum, P. protegens and P. stipitis. Likewise, the upregulation of 23 transcription factors and 22 genes implicated in the RNA metabolism reflected the response of the microbial community. Specifically, the upregulated transcription factors prtT (Ps_47947), YjiE (Rmet_0397), NRG1 (Ps_32504) and NahR (Rmet_0894) are possibly implicated in the regulation of defence processes, such as protease biosynthesis, response to environmental stress , regulation of biofilm formation  and secondary metabolism , respectively. Although different genes were implicated, key functional processes modulated by the simplified soil microcosm in response to the simultaneous introduction of A. mellea and T. atroviride revealed expression profiles similar to those for the condition in which only T. atroviride was introduced. The response to the antagonistic fungus may therefore dominate neutral reactions to the less competitive plant pathogen A. mellea. However, 180 genes with an unknown function were upregulated, indicating that several yet-to-be identified proteins may have been involved in the response of the simplified soil microcosm.
Trichoderma atroviride biocontrol processes were already activated by the simplified soil microcosm and Armillaria mellea did not expresses antagonistic traits
T. atroviride DEGs included 30 and 58 genes up and down-regulated during incubation in the presence of A. mellea, respectively (cluster 17). In particular, upregulated genes of calcium and cAMP-mediated signalling were possibly implicated in the response of T. atroviride to A. mellea (Fig. 5, Additional file 12). The activation of two and five genes related to growth regulation and metabolic processes, respectively suggest stimulation of T. atroviride growth in the simplified soil microcosm containing A. mellea, in agreement with the transcriptional response during mycoparasitic interaction with R. solani . Genes implicated in the metabolism of aromatic compounds were upregulated (a 6-hydroxy-D-nicotine oxidase, a 3-hydroxybenzoate 6-hydroxylase and a 4-hydroxy-2-oxovalerate aldolase) and they are possibly involved in the catabolism of microbial compounds. A T. atroviride peptidase was upregulated and is potentially implicated in biocontrol activity against A. mellea. Key lytic enzymes in biocontrol activity against fungal plant pathogens [76, 82, 83] were already expressed by T. atroviride during incubation in the simplified soil microcosm (SSM+T sample), and they were not further modulated by A. mellea introduction (SSM+T+A sample), as in the case of 23 proteases, 13 ß-glucanases, ten chitinases and two glucosaminidases. Likewise, two PKSs, eight NRPSs, 17 enzymes of toxin metabolism and four of terpenoid metabolism (terpene, trichothecene and trichodiene), six siderophore transporters and 23 ABC transporters possibly implicated in the production of toxic biocontrol molecules and antifungal components [84, 85] were highly expressed by T. atroviride in the simplified soil microcosm. These expression profiles suggested that biocontrol processes were already activated by T. atroviride incubated in the soil matrix to compete with soil microorganisms and the simplified soil microcosm consequently reacted by activating defence processes and detoxification mechanisms. On the other hand, the expression of genes encoding lytic enzymes and toxic secondary metabolites was underrepresented in A. mellea as compared with T. atroviride. Specifically, only 12 proteases, seven ß-glucanases, eight chitinases were expressed by A. mellea during incubation in the simplified soil microcosm in the presence (SSM+T+A) or absence (SSM+A sample) of T. atroviride and no PKSs, NRPSs and genes of siderophore transport and toxin efflux were expressed (Additional file 12), supporting that the activation of defence mechanisms by the simplified soil microcosm was not needed in response to this non-competitive intruder. Moreover, A. mellea genes responsible for sugar transport, energy and carbohydrate metabolism were downregulated in response to T. atroviride introduction, possibly reflecting the inhibition of A. mellea growth in the presence of the biocontrol agent. Likewise, genes implicated in signal transduction processes (a calmodulin and a Ras-related protein) and cell wall reinforcement were downregulated, suggesting that no appropriate cellular reactions were activated by A. mellea in response to T. atroviride in the simplified soil microcosm.
This transcriptome analysis of a simplified soil microcosm is a useful approach to model the complex microbial responses of a microbial community. The change in the growing environment, including both the nutrient-limiting conditions of the soil matrix and the presence of individuals of other species, profoundly impact the transcriptome and the metabolic plasticity of the microbial community. The transcriptional responses of the microbial community differed according to the habitus of microbial intruders (plant pathogen or biocontrol agent). The two intruders were specifically recognised by the other members of the simplified soil microcosm: the biocontrol agent stimulated early defence reactions in the soil microbiota, while the plant pathogen seemed to be perceived as a non-competitive intruder, resulting in the downregulation of defence mechanisms. This opposite response of the simplified soil microcosm agreed with the stronger expression of lytic enzymes and antimicrobial pathways by the biocontrol agent as compared with the plant pathogen. Key biocontrol processes of T. atroviride were already activated, most probably by the other microorganisms of the simplified soil microcosm, and they were not further enhanced by A. mellea introduction, in agreement with the generic antagonistic properties of T. atroviride rather than specific activity against A. mellea. These results represent a step towards understanding the molecular interplay underlying interactions among soil microorganisms and the early transcriptional impact of biocontrol agents on non-target microorganisms. However, further functional and metabolomic studies with a broad range of different soil microorganisms and sampling time points are required in order to deeply characterize the recognition mechanisms of competitive and non-competitive intruders.
Colony forming units
Differentially expressed genes
Efficiency of the quantitative real-time PCR
False discovery rate
Fragments per kilobase per million of unambiguously mapped read pairs
KEGG pathways using the KEGG Automatic Annotation Server
Malt extract agar
Nonribosomal peptide synthetase
Potato dextrose agar
Quantitative real-time PCR
Simplified soil microcosm incubated for 24 h either without exogenous fungi
Simplified soil microcosm incubated for 24 h with the plant pathogen Armillaria mellea
Simplified soil microcosm incubated for 24 h with the biocontrol agent Trichoderma atroviride
Simplified soil microcosm incubated for 24 h with the biocontrol agent T. atroviride and the plant pathogen A. mellea
Simplified soil microcosm incubated collected at the beginning of the experiment just after inoculation of the 11 soil microorganisms
Hofman J, Bezchlebová J, Dušek L, Doležal L, Holoubek I, Anděl P, et al. Novel approach to monitoring of the soil biological quality. Environ Int. 2003;28:771–8.
Gkarmiri K, Finlay R, Alstrom S, Thomas E, Cubeta M, Hogberg N. Transcriptomic changes in the plant pathogenic fungus Rhizoctonia solani AG-3 in response to the antagonistic bacteria Serratia proteamaculans and Serratia plymuthica. BMC Genomics. 2015;16:630.
Nazir R, Warmink JA, Boersma H, Van Elsas JD. Mechanisms that promote bacterial fitness in fungal-affected soil microhabitats. Fems Microbiol Ecol. 2010;71:169–85.
Wargo MJ, Hogan DA. Fungal-bacterial interactions: a mixed bag of mingling microbes. Curr Opin Microbiol. 2006;9:359–64.
Benoit I, van den Esker MH, Patyshakuliyeva A, Mattern DJ, Blei F, Zhou M, et al. Bacillus subtilis attachment to Aspergillus niger hyphae results in mutually altered metabolism. Environ Microbiol. 2015;17:2099–113.
Chapelle E, Mendes R, Bakker PAHM, Raaijmakers JM. Fungal invasion of the rhizosphere microbiome. ISME J. 2015: doi:10.1038/ismej.2015.82.
Mathioni SM, Patel N, Riddick B, Sweigard JA, Czymmek KJ, Caplan JL, et al. Transcriptomics of the rice blast fungus Magnaporthe oryzae in response to the bacterial antagonist Lysobacter enzymogenes reveals candidate fungal defense response genes. PLoS One. 2013;8:e76487.
Mela F, Fritsche K, de Boer W, van Veen JA, de Graaff LH, van den Berg M, et al. Dual transcriptional profiling of a bacterial/fungal confrontation: Collimonas fungivorans versus Aspergillus niger. ISME J. 2011;5:1494–504.
Massart S, Perazzolli M, Höfte M, Pertot I, Jijakli MH. Impact of the omic technologies for understanding the modes of action of biological control agents against plant pathogens. BioControl. 2015;60:725–46.
Neupane S, Finlay RD, Alström S, Elfstrand M, Högberg N. Transcriptional responses of the bacterial antagonist Serratia plymuthica to the fungal phytopathogen Rhizoctonia solani. Environ Microbiol Rep. 2015;7:123–7.
Baumgartner K, Coetzee MP, Hoffmeister D. Secrets of the subterranean pathosystem of Armillaria. Mol Plant Pathol. 2011;12:515–34.
Travadon R, Smith ME, Fujiyoshi P, Douhan GW, Rizzo DM, Baumgartner K. Inferring dispersal patterns of the generalist root fungus Armillaria mellea. New Phytol. 2012;193:959–69.
Perazzolli M, Bampi F, Faccin S, Moser M, De Luca F, Ciccotti AM, et al. Armillaria mellea induces a set of defense genes in grapevine roots and one of them codifies a protein with antifungal activity. Mol Plant Microbe Interact. 2010;23:485–96.
Aguin O, Mansilla JP, Sainz MJ. In vitro selection of an effective fungicide against Armillaria mellea and control of white root rot of grapevine in the field. Pest Manag Sci. 2006;62:223–8.
Longa CMO, Pertot I, Tosi S. Ecophysiological requirements and survival of a Trichoderma atroviride isolate with biocontrol potential. J Basic Microbiol. 2008;48:269–77.
Pellegrini A, Corneo PE, Camin F, Ziller L, Tosi S, Pertot I. Studying trophic interactions between a plant pathogen and two different antagonistic microorganisms using a 13C-labeled compound and isotope ratio mass spectrometry. Rapid Commun Mass Spectrom. 2012;26:510–6.
Lorito M, Woo SL, Harman GE, Monte E. Translational research on Trichoderma: from ‘omics to the field. Annu Rev Phytopathol. 2010;48:395–417.
Barret M, Frey-Klett P, Boutin M, Guillerm-Erckelboudt AY, Martin F, Guillot L, et al. The plant pathogenic fungus Gaeumannomyces graminis var. tritici improves bacterial growth and triggers early gene regulations in the biocontrol strain Pseudomonas fluorescens Pf29Arp. New Phytol. 2009;181:435–47.
Tarkka MT, Sarniguet A, Frey-Klett P. Inter-kingdom encounters: recent advances in molecular bacterium-fungus interactions. Curr Genet. 2009;55:233–43.
Perazzolli M, Antonielli L, Storari M, Puopolo G, Pancher M, Giovannini O, et al. Resilience of the natural phyllosphere microbiota of the grapevine to chemical and biological pesticides. Appl Environ Microbiol. 2014;80:3585–96.
Savazzini F, Longa CMO, Pertot I. Impact of the biocontrol agent Trichoderma atroviride SC1 on soil microbial communities of a vineyard in northern Italy. Soil Biol Biochem. 2009;41:1457–65.
Leibniz Institute DSMZ - German Collection of Microorganisms and cell Culture. www.dsmz.de. Accessed 01 Sept 2013.
ATTC Collection. www.lgcstandards-atcc.org. Accessed 01 Sept 2013.
CBS Collection. www.cbs.knaw.nl. Accessed 01 Sept 2013.
Ellis RJ. Artificial soil microcosms: a tool for studying microbial autecology under controlled conditions. J Microbiol Methods. 2004;56:287–90.
Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30:2114–20.
FASTQC: a quality control tool for high throughput sequence data. http://www.bioinformatics.bbsrc.ac.uk/projects/fastqc. Accessed 1 Jan 2014.
Kim D, Pertea G, Trapnell C, Pimentel H, Kelley R, Salzberg SL. TopHat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions. Genome Biol. 2013;14:R36.
Sibthorp C, Wu H, Cowley G, Wong PW, Palaima P, Morozov I, et al. Transcriptome analysis of the filamentous fungus Aspergillus nidulans directed to the global identification of promoters. BMC Genomics. 2013;14:847.
Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, et al. The sequence alignment/map format and SAMtools. Bioinformatics. 2009;25:2078–9.
Robinson JT, Thorvaldsdottir H, Winckler W, Guttman M, Lander ES, Getz G, et al. Integrative genomics viewer. Nat Biotechnol. 2011;29:24–6.
Anders S, Pyl PT, Huber W. HTSeq - a Python framework to work with high-throughput sequencing data. Bioinformatics. 2015;31:166–9.
Jung JY, Lee SH, Jin HM, Hahn Y, Madsen EL, Jeon CO. Metatranscriptomic analysis of lactic acid bacterial gene expression during kimchi fermentation. Int J Food Microbiol. 2013;163:171–9.
Johansson H, Dhaygude K, Lindstrom S, Helantera H, Sundstrom L, Trontti K. A metatranscriptomic approach to the identification of microbiota associated with the ant Formica exsecta. PLoS One. 2013;8:e79777.
Mortazavi A, Williams BA, McCue K, Schaeffer L, Wold B. Mapping and quantifying mammalian transcriptomes by RNA-Seq. Nat Methods. 2008;5:621–8.
Di Bella JM, Bao Y, Gloor GB, Burton JP, Reid G. High throughput sequencing methods and analysis for microbiome research. J Microbiol Methods. 2013;95:401–14.
Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15:550.
Saeed AI, Bhagabati NK, Braisted JC, Liang W, Sharov V, Howe EA, et al. TM4 microarray software suite. Methods Enzymol. 2006;411:134–93.
Falda M, Toppo S, Pescarolo A, Lavezzo E, Di Camillo B, Facchinetti A, et al. Argot2: a large scale function prediction tool relying on semantic similarity of weighted Gene Ontology terms. BMC Bioinformatics. 2012;13:S14.
Ashburner M, Ball CA, Blake JA, Botstein D, Butler H, Cherry JM, et al. Gene ontology: tool for the unification of biology. The gene ontology consortium. Nat Genet. 2000;25:25–9.
UniProt. http://www.uniprot.org/. Accessed 01 Jan 2016.
Maere S, Heymans K, Kuiper M. BiNGO: a cytoscape plugin to assess overrepresentation of gene ontology categories in biological networks. Bioinformatics. 2005;21:3448–9.
Shannon P, Markiel A, Ozier O, Baliga NS, Wang JT, Ramage D, et al. Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res. 2003;13:2498–504.
Cherry JM, Adler C, Ball C, Chervitz SA, Dwight SS, Hester ET, et al. SGD: saccharomyces genome database. Nucleic Acids Res. 1998;26:73–9.
Saccharomyces cerevisiae S288c. http://fungi.ensembl.org/Saccharomyces_cerevisiae/Info/Index. Accessed 01 Sept 2014.
Moriya Y, Itoh M, Okuda S, Yoshizawa AC, Kanehisa M. KAAS: an automatic genome annotation and pathway reconstruction server. Nucleic Acids Res. 2007;35:W182–5.
KAAS: KEGG Automatic Annotation Server. http://www.genome.jp/tools/kaas. Accessed 01 Jan 2016.
Yamada T, Letunic I, Okuda S, Kanehisa M, Bork P. iPath2.0: interactive pathway explorer. Nucleic Acids Res. 2011;39:W412–5.
iPath2.0: interactive pathway explorer. http://pathways.embl.de/. Accessed 01 Jan 2016.
T-Coffee Multiple Sequence Alignment. http://www.ebi.ac.uk/Tools/msa/tcoffee/. Accessed 01 Mar 2016.
Ruijter JM, Ramakers C, Hoogaars WM, Karlen Y, Bakker O, van den Hoff MJ, et al. Amplification efficiency: linking baseline and bias in the analysis of quantitative PCR data. Nucleic Acids Res. 2009;37:e45.
Hellemans J, Mortier G, De Paepe A, Speleman F, Vandesompele J. qBase relative quantification framework and software for management and automated analysis of real-time quantitative PCR data. Genome Biol. 2007;8:R19.
Silver N, Best S, Jiang J, Thein S. Selection of housekeeping genes for gene expression studies in human reticulocytes using real-time PCR. BMC Mol Biol. 2006;7:33.
Lee K-B, De Backer P, Aono T, Liu C-T, Suzuki S, Suzuki T, et al. The genome of the versatile nitrogen fixer Azorhizobium caulinodans ORS571. BMC Genomics. 2008;9:271.
Barbe V, Cruveiller S, Kunst F, Lenoble P, Meurice G, Sekowska A, et al. From a consortium sequence to a unified sequence: the Bacillus subtilis 168 reference genome a decade later. Microbiology. 2009;155:1758–75.
Janssen PJ, Van Houdt R, Moors H, Monsieurs P, Morin N, Michaux A, et al. The complete genome sequence of Cupriavidus metallidurans strain CH34, a master survivalist in harsh and anthropogenic environments. PLoS One. 2010;5:e10433.
Paulsen IT, Press CM, Ravel J, Kobayashi DY, Myers GSA, Mavrodi DV, et al. Complete genome sequence of the plant commensal Pseudomonas fluorescens Pf-5. Nat Biotech. 2005;23:873–8.
Dujon B, Sherman D, Fischer G, Durrens P, Casaregola S, Lafontaine I, et al. Genome evolution in yeasts. Nature. 2004;430:35–44.
Kumar S, Randhawa A, Ganesan K, Raghava GP, Mondal AK. Draft genome sequence of salt-tolerant yeast Debaryomyces hansenii var. hansenii MTCC 234. Eukaryot Cell. 2012;11:961–2.
Jeffries TW, Grigoriev IV, Grimwood J, Laplaza JM, Aerts A, Salamov A, et al. Genome sequence of the lignocellulose-bioconverting and xylose-fermenting yeast Pichia stipitis. Nat Biotech. 2007;25:319–26.
Wood V, Gwilliam R, Rajandream MA, Lyne M, Lyne R, Stewart A, et al. The genome sequence of Schizosaccharomyces pombe. Nature. 2002;415:871–80.
Pel HJ, de Winde JH, Archer DB, Dyer PS, Hofmann G, Schaap PJ, et al. Genome sequencing and analysis of the versatile cell factory Aspergillus niger CBS 513.88. Nat Biotech. 2007;25:221–31.
van den Berg MA, Albang R, Albermann K, Badger JH, Daran J-M, M Driessen AJ, et al. Genome sequencing and analysis of the filamentous fungus Penicillium chrysogenum. Nat Biotech. 2008;26:1161–8.
Ma L-J, van der Does HC, Borkovich KA, Coleman JJ, Daboussi M-J, Di Pietro A, et al. Comparative genomics reveals mobile pathogenicity chromosomes in Fusarium. Nature. 2010;464:367–73.
Kubicek C, Herrera-Estrella A, Seidl-Seiboth V, Martinez D, Druzhinina I, Thon M, et al. Comparative genome sequence analysis underscores mycoparasitism as the ancestral life style of Trichoderma. Genome Biol. 2011;12:R40.
Collins C, Keane TM, Turner DJ, O’Keeffe G, Fitzpatrick DA, Doyle S. Genomic and proteomic dissection of the ubiquitous plant pathogen, Armillaria mellea: toward a new infection model system. J Proteome Res. 2013;12:2552–70.
Zenoni S, Ferrarini A, Giacomelli E, Xumerle L, Fasoli M, Malerba G, et al. Characterization of transcriptional complexity during berry development in Vitis vinifera using RNA-Seq. Plant Physiol. 2010;152:1787–95.
Perazzolli M, Moretto M, Fontana P, Ferrarini A, Velasco R, Moser C, et al. Downy mildew resistance induced by Trichoderma harzianum T39 in susceptible grapevines partially mimics transcriptional changes of resistant genotypes. BMC Genomics. 2012;13:660.
Tyc O, Wolf AB, Garbeva P. The effect of phylogenetically different bacteria on the fitness of Pseudomonas fluorescens in sand microcosms. PLoS One. 2015;10:e0119838.
Stout JD. The role of protozoa in nutrient cycling and energy flow. In: Alexander M, editor. Advances in microbial ecology, vol. 5. Boston: Springer US; 1980. p. 1–50.
Leonhartsberger S, Korsa I, Bock A. The molecular biology of formate metabolism in enterobacteria. J Mol Microbiol Biotechnol. 2002;4:269–76.
Turnbull GA, Morgan JA, Whipps JM, Saunders JR. The role of bacterial motility in the survival and spread of Pseudomonas fluorescens in soil and in the attachment and colonisation of wheat roots. Fems Microbiol Ecol. 2001;36:21–31.
Hibbing ME, Fuqua C, Parsek MR, Peterson SB. Bacterial competition: surviving and thriving in the microbial jungle. Nat Rev Microbiol. 2010;8:15–25.
Pal KK, McSpadden GB. Biological control of plant pathogens. Plant Health Instr. 2006;10:1094–117.
Deveau A, Barret M, Diedhiou A, Leveau J, de Boer W, Martin F, et al. Pairwise transcriptomic analysis of the interactions between the ectomycorrhizal fungus Laccaria bicolor s238n and three beneficial, neutral and antagonistic soil bacteria. Microb Ecol. 2015;69:146–59.
Atanasova L, Crom SL, Gruber S, Coulpier F, Seidl-Seiboth V, Kubicek CP, et al. Comparative transcriptomics reveals different strategies of Trichoderma mycoparasitism. BMC Genomics. 2013;14:121.
Seidl V, Song L, Lindquist E, Gruber S, Koptchinskiy A, Zeilinger S, et al. Transcriptomic response of the mycoparasitic fungus Trichoderma atroviride to the presence of a fungal prey. BMC Genomics. 2009;10:567.
Gebendorfer KM, Drazic A, Le Y, Gundlach J, Bepperling A, Kastenmuller A, et al. Identification of a hypochlorite-specific transcription factor from Escherichia coli. J Biol Chem. 2012;287:6892–903.
Braun BR, Kadosh D, Johnson AD. NRG1, a repressor of filamentous growth in C. albicans, is down‐regulated during filament induction. EMBO J. 2001;20:4753–61.
Park W, Padmanabhan P, Padmanabhan S, Zylstra GJ, Madsen EL. nahR, encoding a LysR-type transcriptional regulator, is highly conserved among naphthalene-degrading bacteria isolated from a coal tar waste-contaminated site and in extracted community DNAb. Microbiology. 2002;148:2319–29.
Reithner B, Ibarra-Laclette E, Mach RL, Herrera-Estrella A. Identification of mycoparasitism-related genes in Trichoderma atroviride. Appl Environ Microbiol. 2011;77:4361–70.
Viterbo A, Ramot O, Chernin L, Chet I. Significance of lytic enzymes from Trichoderma spp. in the biocontrol of fungal plant pathogens. Antonie Van Leeuwenhoek. 2002;81:549–56.
Steindorff AS, Ramada MHS, Coelho ASG, Miller RNG, Pappas GJ, Ulhoa CJ, et al. Identification of mycoparasitism-related genes against the phytopathogen Sclerotinia sclerotiorum through transcriptome and expression profile analysis in Trichoderma harzianum. BMC Genomics. 2014;15:1–14.
Mukherjee PK, Horwitz BA, Kenerley CM. Secondary metabolism in Trichoderma – a genomic perspective. Microbiology. 2012;158:35–45.
Zeilinger S, Gruber S, Bansal R, Mukherjee PK. Secondary metabolism in Trichoderma – Chemistry meets genomics. Fungal Biol Rev. 2016;30:74–90.
Azorhizobium caulinodans ORS 571. http://bacteria.ensembl.org/azorhizobium_caulinodans_ors_571/Info/Index. Accessed 01 Sept 2014.
Bacillus subtilis subtilis 168 http://bacteria.ensembl.org/Bacillus_subtilis_subsp_subtilis_str_168/Info/Index. Accessed 01 Sept 2014.
Cupriavidus metallidurans CH34. http://bacteria.ensembl.org/cupriavidus_metallidurans_ch34/Info/Index. Accessed 01 Sept 2014.
Pseudomonas protegens Pf-5. http://bacteria.ensembl.org/pseudomonas_protegens_pf_5/Info/Index. Accessed 01 Sept 2014.
Debaryomyces hansenii CBS767. http://genome.jgi.doe.gov/Debha1/Debha1.home.html. Accessed 01 Sept 2014.
Pichia stipitis CBS6054. http://genome.jgi-psf.org/Picst3/Picst3.home.html. Accessed 01 Sept 2014.
Schizosaccharomyces pombe 972h. http://fungi.ensembl.org/Schizosaccharomyces_pombe/Info/Index. Accessed 01 Sept 2014.
Aspergillus niger CBS 513.88. http://fungi.ensembl.org/Aspergillus_niger/Info/Index. Accessed 01 Sept 2014.
Fusarium oxysporum lycopersici 4287. http://fungi.ensembl.org/Fusarium_oxysporum/Info/Index. Accessed 01 Sept 2014.
Penicillium chrysogenum Wisconsin54-125. http://genome.jgi.doe.gov/PenchWisc1_1/PenchWisc1_1.home.html. Accessed 01 Sept 2014.
Trichoderma atroviride IMI 206040. http://genome.jgi.doe.gov/Triat2/Triat2.home.html. Accessed 01 Sept 2014.
Armillaria mellea DSM 3731. http://genome.jgi.doe.gov/Armme1_1/Armme1_1.home.html. Accessed 01 Sept 2014.
The authors thank Chiara Pellegrini, Carmela Sicher and Denise Ress for technical support in preparing the simplified soil microcosm and Simone Prospero (Swiss Federal Research Institute WSL, Birmensdorf) for providing the haploid Armillaria mellea strain M6132.
This project has received funding from the European Union’s Seventh Framework Programme under grant agreement: 324416 (project INNOVA, subprogramme: FP7-PEOPLE-2012-IAPP).
Availability of data and materials
The sequences were deposited at the Sequence Read Archive of the National Center for Biotechnology (https://www.ncbi.nlm.nih.gov/sra/?term=SRP065901) under accession number SRP065901 and BioProject number PRJNA281302.
NH and AP carried out the simplified soil microcosm experiments, sample collection and RNA extraction. MP, LS and YVdP analysed the RNA-Seq data through mapping, gene expression estimation, annotation and statistical analysis. LL and MP carried out the RT-qPCR experiments. MP, GP and IP contributed to data interpretation and manuscript writing. IP conceived the study, designed the experiment and coordinated all research activities. All the authors approved the final manuscript.
The authors declare that they have no competing interests.
Consent to publish
Ethics approval and consent to participate
Concentration of the soil microorganism in the simplified soil microcosm. (DOCX 18 kb)
Primer sequences of the microcosm genes analysed by real-time RT-PCR. (XLSX 18 kb)
RNA-Seq sequencing and mapping results for each replicate. (PDF 14 kb)
Distribution of read pair alignments to the genomes of the soil microorganisms. (A, B) Distribution of read pair alignments to the 13 soil microorganisms calculated with the Samtools software  and expressed as a percentage (%) of total alignments to the microcosm genome. (C, D) Distribution of unique read pairs mapping to genes of the 13 soil microorganisms, counted using HTSeq  and expressed as a percentage (%) of total unique read pairs mapping to genes in the microcosm genome. (E, F) Percentage (%) of expressed genes (more than one read pair) calculated as compared to the total predicted genes for each soil microorganism. Mean and standard error values of three replicates are reported for each condition: the simplified soil microcosm collected at the beginning of the experiment (SSM0) and 24 h after incubation either without exogenous fungi (SSM), with the biocontrol agent Trichoderma atroviride (SSM+T), with the plant pathogen Armillaria mellea (SSM+A) or with both (SSM+T+A). (PDF 23 kb)
Expression levels of genes of the simplified soil microcosm. (XLSX 39274 kb)
Pearson’s correlation coefficients among replicates and conditions for RNA-Seq analysis. (PDF 6 kb)
Clustering and functional annotation results of differentially expressed genes. (XLSX 15048 kb)
Proportion of expressed and differentially expressed genes for each soil microorganism. (PDF 8 kb)
Distribution of differentially expressed genes of each soil microorganism in 18 clusters, based on the expression profiles. (XLSX 12 kb)
Metabolic pathways of the simplified soil microcosm, modulated by incubation in the soil matrix. Metabolic pathways deactivated (left panels) and activated (right panels) by 24 h incubation in the soil matrix (A) without reinforced modulation (cluster 1) and (B) with reinforced modulation (cluster 15) in the presence of Armillaria mellea and Trichoderma atroviride combined. Metabolic pathways modulated by the introduction of (C) T. atroviride (cluster 3), (D) A. mellea (cluster 5), or (E) both (cluster 7). KEGG pathways were visualised using the iPath2 tool , the pathways of upregulated (green) and downregulated (red) genes were highlighted, and a section of the most relevant pathways is reported for each panel. (PDF 2815 kb)
Biological networks of Gene Ontology (GO) terms. GO biological process terms of the simplified soil microcosm, upregulated by incubation in the soil matrix with similar expression profiles in the presence or absence of Armillaria mellea and Trichoderma atroviride (cluster 1). Significantly enriched GO terms (P < 0.001) were identified using the BiNGO tool  and visualised with Cytoscape software . The colour scale legend indicates the level of significance for enriched GO terms. White nodes indicate not significantly overrepresented categories. (PDF 111 kb)
Key differentially expressed genes discussed in the manuscript, based on their functional categories and expression profiles. Each sheet contains the genes in each cluster discussed. (XLSX 23094 kb)
About this article
Cite this article
Perazzolli, M., Herrero, N., Sterck, L. et al. Transcriptomic responses of a simplified soil microcosm to a plant pathogen and its biocontrol agent reveal a complex reaction to harsh habitat. BMC Genomics 17, 838 (2016) doi:10.1186/s12864-016-3174-4
- Soil microbial community
- Soil transcriptome
- Biological control
- Plant pathogen
- Microbial interaction
- Gene expression