Skip to main content

A novel plant-fungal association reveals fundamental sRNA and gene expression reprogramming at the onset of symbiosis

Abstract

Background

Beneficial associations between plants and microbes are widespread in nature and have been studied extensively in the microbial-dominant environment of the rhizosphere. Such associations are highly advantageous for the organisms involved, benefiting soil microbes by providing them access to plant metabolites, while plant growth and development are enhanced through the promotion of nutrient uptake and/or protection against (a)biotic stresses. While the establishment and maintenance of mutualistic associations have been shown to require genetic and epigenetic reprogramming, as well as an exchange of effector molecules between microbes and plants, whether short RNAs are able to effect such changes is currently unknown. Here, we established an interaction between the model grass species Brachypodium distachyon (Bd, Pooideae) and the beneficial fungal root endophyte Serendipita indica (Si, syn. Piriformospora indica, Sebacinales) to elucidate RNA interference-based regulatory changes in gene expression and small (s)RNA profiles that occurred during establishment of a Sebacinalean symbiosis.

Results

Colonization of Bd roots with Si resulted in higher grain yield, confirming the mutualistic character of this interaction. Resequencing of the Si genome using the Oxford Nanopore technique, followed by de novo assembly yielded in 57 contigs and 9441 predicted genes, including putative members of several families involved in sRNA production. Transcriptome analysis at an early stage of the mutualistic interaction identified 2963 differentially expressed genes (DEG) in Si and 317 in Bd line 21-3. The fungal DEGs were largely associated with carbohydrate metabolism, cell wall degradation, and nutrient uptake, while plant DEGs indicated modulation of (a)biotic stress responses and defense pathways. Additionally, 10% of the upregulated fungal DEGs encode candidate protein effectors, including six DELD proteins typical for Sebacinales. Analysis of the global changes in the sRNA profiles of both associated organisms revealed several putative endogenous plant sRNAs expressed during colonization belonging to known micro (mi)RNA families involved in growth and developmental regulation. Among Bd- and Si-generated sRNAs with putative functions in the interacting organism, we identified transcripts for proteins involved in circadian clock and flowering regulation as well as immunity as potential targets of fungal sRNAs, reflecting the beneficial activity of Si.

Conclusions

We detected beneficial effects of Si colonization on Bd growth and development, and established a novel plant-mutualist interaction model between these organisms. Together, the changes in gene expression and identification of interaction-induced sRNAs in both organisms support sRNA-based regulation of defense responses and plant development in Bd, as well as nutrient acquisition and cell growth in Si. Our data suggests that a Sebacinalean symbiosis involves reciprocal sRNA targeting of genes during the interaction.

Background

Mutualistic associations between beneficial microbes and plants are widespread and highly advantageous, especially in the microbial-dominant environment of the rhizosphere. This relationship benefits soil microbes by providing them access to plant metabolites; in return, they enhance plant growth and development by promoting nutrient uptake and/or protection against (a)biotic stresses [1, 2]. The beneficial or parasitic outcome of a plant-microbe interaction is governed by the genotype and physiological status of the host, identity of the microbe, and environmental factors such as soil type and nutrient availability [3, 4]. The establishment and maintenance of mutualistic associations (here called symbiosis) require genetic and epigenetic reprogramming and metabolome modulation by the exchange of effector molecules between the beneficial microbe and the plant [5, 6]. Beneficial microbes have a significant impact on crop production, due to their effects on plant health and yield. However, considerable gaps in knowledge prior to their establishment in agricultural practice remain, including systemic identification of microbial abundance and diversity in various ecosystems, understanding the influence of climate, soil conditions, management practices, and, lastly, elucidating the intricacies of molecular mechanisms governing establishment of colonization and nutrient acquisition [7].

Crucial for regulation of gene expression, RNA interference (RNAi) is a well-known eukaryotic gene silencing mechanism [8], mediated by small RNAs (sRNAs) of 20–24 nucleotides (nt) in size and RNAi-associated proteins, primarily from the Argonaute (AGO), Dicer-like (DCL) and RNA-dependent RNA polymerase (RdRP) families [9]. DCLs generate sRNAs from longer RNA molecules, whereas AGOs bind sRNAs within an RNA-induced silencing complex (RISC). In the context of plant-microbe interactions, microbial protein effectors are known to promote pathogen colonization by suppressing host immune responses [10] and have been described in mutualistic associations as well [5], including the Sebacinalean symbiosis [11]. Recent findings suggest that sRNAs, through RNAi-based regulatory mechanisms, also can serve as effectors of pathogenic microbes [12], whereby the sRNA is secreted to suppress translation of a host mRNA via RNAi. Conversely, plants can secrete sRNAs that target virulence-associated mRNAs in the microbe [13]. This transfer of sRNAs and subsequent gene silencing in the target organism is called cross-kingdom RNAi [12].

We studied the association of the beneficial fungus Serendipita indica (Si) with the model grass Brachypodium distachyon (Bd, purple false brome, Pooideae [14];). Si is an endophytic fungus belonging to the order Sebacinales that colonizes the rhizodermis and cortex of a large spectrum of plants [15]. Si serves as an excellent model for beneficial microbes as it (i) primes plants for disease resistance against biotrophic [16] and necrotrophic [17] fungi, oomycetes [18], and viruses [19]; (ii) enhances the tolerance of plants against abiotic stress [20]; (iii) promotes growth and yield [21]; (iv) has a sequenced 25 Mb genome [22]; and (v) is genetically transformable and culturable in axenic conditions [23]. Si initially undergoes a biotrophic growth phase during Arabidopsis (Arabidopsis thaliana) and barley (Hordeum vulgare) colonization, with suppression of innate immune responses [24, 25] and activation of induced systemic resistance [16]. Subsequently, Si colonization of barley enters a cell-death associated phase and switches to a saprophytic lifestyle [26, 27].

Bd is a temperate grass species belonging to the Pooideae subfamily and a model for genetic studies of stress resistance and yield parameters of cereals [28]. Bd is self-pollinating, genetically transformable, easy to cultivate, and has a sequenced genome of 272 Mb [29,30,31]. It shares evolutionary proximity and broad synteny with complex crop genomes, such as wheat and rice [14], and is a host for major cereal pathogens [32]. Additionally, RNAi is operational in Bd, with proven alteration of micro RNA (miRNA) expression patterns in response to abiotic stresses [33, 34]. In silico analyses revealed that the Bd genome, similar to other cereals, contains an expansion of DCL and AGO families [35].

Currently, the significance of cross-kingdom RNAi in mutualistic interactions is largely unknown. A recent in silico study predicted that the arbuscular mycorrhizal fungus Rhizophagus irregularis generates sRNAs, which have predicted targets in the host plant Medicago truncatula [36]. Moreover, tRNA-derived sRNA fragments from rhizobial bacteria were shown to regulate host nodulation-associated genes by utilizing the host's RNAi machinery [37]. To investigate the role of sRNAs in another agronomically relevant mutualistic interaction, we established a protocol for Si colonization of the model Pooideae Bd. Additionally, integrative high-throughput sequencing and transcriptome analysis assessed symbiosis-associated changes in the mRNA and sRNA expression patterns of both organisms. We discuss here possible sRNA-based regulation that might be critical for the establishment of the Sebacinalean symbiosis.

Results

Brachypodium distachyon Bd21-3 forms a mutualistic interaction with Serendipita indica

To investigate whether Bd can develop a beneficial interaction with Si, we established an inoculation protocol using one-week-old seedlings of Bd line Bd21-3, with dip-inoculation in 5 × 105 chlamydospores ml−1 for 3 h. Comparison of grain production in fully mature, colonized vs. non-colonized plants grown in soil showed that Si increased the number of filled grains/plant by 49.9% (Fig. 1a, Additional file 1), and total grain weight/plant increased by 38.1% (Fig. 1b, Additional file 1). Consistent with the observation that Si-colonized Bd21-3 plants flower several days earlier than control plants, they exhibited a 32.2% increase in the number of spikelets at 2 months after inoculation (Fig. 1c, Additional file 1). Concordantly, growth and biomass analyses of Bd21-3 seedlings revealed a significant 8.6% increase in shoot length (Fig. 1d) upon Si colonization.

Fig. 1
figure1

Root colonization by Serendipita indica (Si) increases growth and yield of Brachypodium distachyon (Bd) Bd21-3. a Number of full grains produced by non-colonized (control) vs. Si-colonized plants. Sample size n = 5. b Total grain weight of control vs. colonized plants. Sample size n = 5. c Number of spikelets of control vs. colonized plants. Sample size n = 5. d Shoot length of control vs. colonized plants. Sample size n = 20. The results are from three independent biological replicates that are represented by green, black, and red dots on the scatter graph, individual data values for 1a–c in Additional File 1. For a and b, 1-week-old seedlings were inoculated with 5 × 105 chlamydospores per ml and grown for approximately 3 months in F-E type LD 80 soil; for c, spikelets were counted on 2-month-old plants; for d, seedlings were grown for 3 weeks on a vermiculite:oil dri mixture (semi-sterile conditions). The significance threshold for p values, corrected for multiple testing (Benjamini-Hochberg) was set at 0.05 (*≤ 0.05) and the Si effect was calculated as ((Mean_Si-Mean_Control)/Mean_Control) × 102)

Further analysis of Bd21-3 seedlings indicated that Si colonization increased lateral root growth, as early as 4 days post inoculation (4 DPI, Additional file 2: Figure S1a). By 25 DPI, roots showed a more extensively branched structure (Additional file 2: Figure S1b). Microscopy of Si-inoculated Bd21-3 roots confirmed root surface colonization and proliferation of fungal spores after staining with chitin-specific WGA-AF 488 at 4 DPI (Fig. 2a–d) and further on at 7 DPI and 14 DPI (Additional file 2: Figure S2). Inter- and intracellular colonization of Bd21-3 cells in the root differentiation zone also was visible after WGA-AF 488 and propidium iodide staining (Fig. 2e–h). These results suggest that establishment of a mutualistic symbiosis correlates with observable phenotypic changes by 4 DPI; thus, this time point was used to further investigate the Bd21-3-Si system.

Fig. 2
figure2

Colonization pattern of Serendipita indica (Si) on Brachypodium distachyon Bd21-3 roots. ad Colonization at 4 DPI. a Fluorescence microscopy showing WGA-AF488 staining of Si cell walls (λexc494 nm, λem515). b Fluorescence control (λexc631 nm, λem642). c Bright-field microscopy to visualize Si chlamydospores. d Overlay showing Si chlamydospores (red arrows), which have germinated and formed a hyphal network on the root surface (blue arrows). eh Rhizodermal root colonization by Si at 4 DPI. e Fluorescence microscopy showing WGA-AF488 staining of Si cell walls (λexc494 nm, λem515). f Fluorescence microscopy to visualize propidium iodide staining of root cell walls (λexc535 nm, λem617). g Bright-field microscopy to visualize root cells. h Overlay showing extensive inter- and intracellular fungal growth on Bd21-3 roots. Imaging was done with a LEICA S8 confocal microscope (e-h: maximum projection; z-stack). For ad, 1-week-old seedlings were inoculated with 5 × 105 chlamydospores per ml and subsequently grown on a plastic mesh over 0.5X MS; for eh, Si-inoculated seedlings were grown on vermiculite:oil dri mix before harvesting at 4 DPI

Resequencing of the Si genome

To improve Si assembly, the genome was resequenced using MinION (25,167 reads, 324 Mb) and MiSeq (18,225,814 reads, 5.46 Gb); together they yielded approx. 6.0 Gb of sequence information. De novo assembly of the Nanopore sequence reads generated 57 contigs, accounting for a total length of 24.7 Mb and a N50 of 1.3 Mb. The draft genome sequence features GC content of 50.8%, similar to the first genome version with 2,359 contigs and GC content of 50.7% [22]. Analyses using the eukaryotic gene prediction tool Genemark-ES 4.33 [38] revealed 9441 predicted genes (75% of the genome; 59,045 exons), 9498 intergenic regions (25% of the genome), and a gene density of 380.68 genes/Mbp (Additional file 2: Table S1). Annotation of the Si genes using a GenDB version designed to process eukaryotic genomes possessing multi-exon genes [39, 40] revealed that 4756 have a predicted function. Comparison of predicted genes from the resequenced Si genome (Si-2020) vs. the 2011 assembly [22] indicated that the vast majority are shared (90.3%), with 915 genes unique to the Si-2020 genome (Additional file 2: Figure S3). There is a reduction in gene model numbers relative to the 2011 assembly, which can be attributed to improved gene prediction tools for eukaryotic organisms and a considerably reduced number of contigs. Additionally, Si shares 2585 genes with another member of the Sebacinales, Serendipita vermifera, while 156 genes are shared only with the ectomycorrhizal fungus Laccaria bicolor and 2729 genes are common in all three species (Additional file 2: Figure S4).

Establishment of the Si-Bd symbiosis is associated with extensive transcriptional reprogramming

To assess how the symbiotic interaction affects gene expression in both organisms, mRNA was sequenced and analyzed (Additional file 2: Table S2) from the roots of Si-colonized Bd21-3 seedlings (sample Bd-Si) and mock-treated plants (Bd-C) at 4 DPI, and from 4-week-old axenic Si cultures (Si-ax). Comparison of reads between Bd-Si and Si-ax identified 2963 differentially expressed fungal genes (DEGs, DESeq2: Wald test, Benjamini-Hochberg (BH) adjustment, padj < 0.05), which accounts for 31.4% of the 9441 predicted Si genes. Comparison of reads from Bd-Si and Bd-C revealed 317 plant DEGs (0.66% out of approximately 47,917 protein-coding transcripts disclosed in the JGI v1.1 annotation, padj < 0.05). The interaction-responsive DEGs in Si and Bd21-3, split into up- and downregulated groups are shown in Fig. 3.

Fig. 3
figure3

Volcano plots of colonization-associated, differentially expressed genes (DEGs). a Serendipita indica (Si) DEGs identified by comparing reads from colonized roots (Bd-Si) vs. axenic mycelium (Si-ax). b Brachypodium distachyon Bd21-3 DEGs identified by comparing colonized (Bd-Si) vs. mock-treated (Bd-C) roots. X-axis displays the log2 FoldChange and Y-axis displays the negative log10 of adjusted p values from DE analysis. The magnitude of up- or downregulation for the DEGs (represented by individual dots) is indicated by different colors, as designated in the legend for each plot

All significant DEGs were submitted to gene ontology term analysis against the reference background for Bd and a customized Si-specific background. The resulting enriched terms relate to metabolic, mainly redox processes and catalytic activity functions (Additional file 2: Table S3). Of the 25 highly downregulated Si DEGs, many encode proteins associated with metabolic reprogramming networks involved in nutrient exchange and adaptation to nutrient availability (Table 1). By contrast, highly upregulated Si DEGs encode for proteins involved in fungal catalytic and hydrolytic processes. This suggests that by 4 DPI, Si has entered a saprophytic-like growth phase similar to that detected in barley roots [22]. To investigate whether any Si DEG encodes effector proteins, a computational pipeline [41] was used to mine 982 genes identified in the Si-2020 genome that encode signal peptide-containing proteins, resulting in 480 putative protein effector genes. In total, 174 (36%) of these were significantly upregulated during colonization of Bd21-3 (Additional file 2: Table S4), including six DELD family proteins [22]. In Bd21-3, many of the highly downregulated interaction-responsive DEGs encode transcription factors or proteins associated with stress responses or circadian clock regulation. Those showing high levels of upregulation include genes linked to immune responses and hormone signaling networks (Table 2 and Additional file 2: Table S3). In order to further validate our sequencing data, we confirmed the expression of five Si and five Bd21-3 DEGs from Tables 1 and 2 by RT-qPCR. Generally, the qPCR results show a similar fold change for the selected genes between the colonized root and the respective controls, compared to the mRNA-seq results (Additional file 2: Figure S5). Together, these results show that both organisms utilize a complex enzymatic arsenal to establish and control the symbiosis.

Table 1 Top 25 Serendipita indica (Si) differentially expressed genes (DEGs) during colonization (4 DPI)
Table 2 Top 25 Brachypodium distachyon (Bd) differentially expressed genes (DEGs) during root colonization (4 DPI)

Prediction of Si RNAi genes

Since RNAi-mediated gene silencing has been documented in most but not all fungi [42], we assessed whether the Si-2020 genome encodes RNAi-related proteins with conserved domain architecture and homology to RNAi components in the model filamentous fungus Neurospora crassa [43]. Genes encoding predicted DCLs (G4U2H0, G4TBW9) with typical domains (dsRNA-binding, RNase III and helicase, [44]), QDE2-like proteins with PIWI domains typical of AGOs (G4TEK0, G4TLO4, [45]), an AGO-like protein (G4T5G9), and RdRPs (G4TNU7, G4TQP0) were identified and were expressed in axenic culture and Bd21-3-associated Si samples (Additional file 2: Table S5). Thus, the Si genome is predicted to contain genes encoding critical components of the RNAi machinery. Based on these new data, and the earlier discovery that AGO and DCL families are expanded in the Bd genome [35], we decided to sequence the sRNAs of both organisms, in order to assess the role of RNAi-based regulation and communication in symbiosis.

sRNA profiles undergo a substantial change at the onset of the Si-Bd symbiosis

To evaluate how the mutualistic interaction affects the sRNA profiles in the colonized root and respective Si and Bd controls, reads from Bd-C, Bd-Si, and Si-ax sRNA data sets were subjected to consecutive filtering steps (Additional file 2: Figure S6). This greatly reduced the number of raw reads to be analyzed and allowed us to distinguish between sRNAs with potential targets in the interacting organism (putative ck-sRNAs) and sRNAs with potential functions in the same organism (putative endogenous sRNAs) (Additional file 2: Table S6).

Analysis of sRNAs from the Bd-Si dataset revealed that the total number of putative ck-sRNAs exceeds that of endogenous sRNAs in both Si (786,732 vs. 261,478) and Bd21-3 (17 million vs. 1.6 million), but the converse is true for unique sRNAs (36,163 endogenous vs. 35,895 putative ck-sRNAs in Si and 483,352 endogenous vs. 286,198 putative ck-sRNAs in Bd21-3).

Size distribution profiles of Si and Bd21-3 sRNAs

Size distribution of sRNA reads from Bd21-3 and Si during colonization, and the respective controls was then assessed. For putative endogenous Bd21-3 sRNAs, peaks at 21 and 24 nt were identified, with the 24 nt sRNAs exhibiting greater diversity than those of 21 nt (Fig. 4a, b). These sizes are consistent with the expected peaks of RNAi-associated sRNAs in plants [46]. For putative endogenous Si sRNAs, a bimodal size distribution pattern was observed in the total fractions, with the first peak at 26 nt and second at 29–30 nt (Fig. 4c, d). A smaller peak of 21 nt long molecules was observed in the Bd-Si but not Si-ax samples, indicating that colonization affects the relative size distribution of Si sRNAs. Since previously identified ck-sRNAs range from 20 to 24 nt [12, 13], the size distribution of putative ck-sRNAs and corresponding reads in the control samples was assessed over a narrower window. Contrary to endogenous sRNAs, ck-sRNAs displayed no prominent peaks in the 20–24 nt window (Additional file 2: Figure S7).

Fig. 4
figure4

Size distribution of total and unique putative endogenous sRNAs. a Bd-C (mock-treated), b Bd-Si (colonized root), c Si-ax (axenic mycelium), and d Bd-Si (colonized root). All datasets represent three biological replicates and corresponding two technical replicates, merged together. sRNA length is displayed on the X-axis (nt) and number of total/unique sRNA counts on the Y-axis (× 103)

Before sRNAs can guide RNAi-mediated gene silencing, they must be loaded onto AGO proteins and assembled into a RISC. Previously, Arabidopsis AGO proteins were shown to preferentially recruit sRNAs with specific 5′ termini [47]. Hence, we analyzed the 5′ terminal nt composition of Bd-C, Si-ax, and Bd-Si sRNAs. For unique putative endogenous and ck-sRNAs, the 5′ nt composition was relatively consistent except for the 24 nt Bd21-3 sRNAs, which exhibited a strong bias towards 5′ A (Additional file 2: Figure S8, Figure S9). The total sRNA fractions exhibited somewhat greater variability in 5′ nt composition. Of the total Bd21-3 endogenous sRNAs, 24 nt molecules from colonized and non-colonized tissue showed a strong bias towards 5′ A, while 21 nt molecules were biased towards a terminal U (Additional file 2: Figure S10), and 20 nt ck-sRNAs had a higher percentage of 5′ Cs (Additional file 2: Figure S11). Of the total endogenous Si sRNAs, those from colonized samples generally had a stronger bias towards 5′ A than sRNA reads from Si-ax, especially at 26 nt and 21 nt (Additional file 2: Figure S10). A slightly higher percentage of 5′ As also was detected in total putative ck-sRNAs of 21 nt (Additional file 2: Figure S11).

Differentially expressed Si and Bd21-3 sRNAs

Analysis of unique plant sRNAs in Bd-C vs. Bd-Si revealed that 63% of the putative endogenous sRNAs were exclusively present in Bd-C, 30% were exclusively in Bd-Si and 7% were present in both (Fig. 5). For the reads from the putative ck-sRNA pipeline, 76% of the reads were exclusively present in Bd-C, 13% were exclusive to Bd-Si, and 11% were found in both. Comparison between the unique fungal sRNAs in Si-ax and Bd-Si indicated that 98.1% of the putative endogenous sRNAs were exclusively present in axenic culture, 0.98% were exclusively found in Bd-Si, and 0.92% were in both. Similarly, from the putative ck-sRNA pipeline, 98.2% of the sRNA reads were exclusive to the Si-ax sample, 0.3 % were exclusive to Bd-Si, and 1.5% were present in both. Considering only the Si sRNAs in the Bd-Si sample, 51.1% of the putative endogenous ones and 15% of the putative ck-sRNAs are exclusively present in the colonized sample. Among Bd21-3 sRNAs in the Bd-Si sample, there are 80.5% putative endogenous sRNAs and 54.3% ck-sRNAs exclusive for the colonized root. These data show that colonization induces many novel putative endogenous and ck-sRNAs in Bd21-3, and a smaller amount in Si, due to fungal quantity in the colonized roots.

Fig. 5
figure5

Venn diagrams showing the sample-exclusive or communal presence of unique putative endogenous or ck-sRNAs. a putative endogenous sRNAs in Bd-C (mock-treated) vs. Bd-Si (colonized root), b sRNA reads in Bd-C vs. putative ck-sRNAs in Bd-Si, c putative endogenous sRNAs in Si-ax (axenic mycelium) vs. Bd-Si (colonized root), and d sRNA reads in Si-ax vs. putative ck-RNAs in Bd-Si

Identification of Bd21-3 miRNAs during the Bd-Si interaction

Using the ShortStack analysis tool, we identified Bd21-3 loci that correspond to putative endogenous sRNAs expressed during Si colonization. The DicerCall function indicated loci whose predominant sRNAs are 20–24 nt. Comparison of these sRNAs with miRBase identified 16 sRNA-generating loci that correlate to known miRNAs (Table 3). These miRNAs belong to highly conserved plant miRNA families that regulate growth and development [33, 34]. We conducted the same analysis with Si sRNAs, but no predicted miRNA-like RNAs were identified in the colonized sample, possibly due to a lack of data about the fungal sRNA-generating loci.

Table 3 Predicted miRNA-generating loci identified in Bd21-3 roots colonized by Serendipita indica

In silico target prediction of putative Si and Bd21-3 ck-sRNAs

Since most examples of sRNA-based communication in plant-microbe interactions have the commonality of 21 nt long sRNAs that silence transcripts in the target organism [12, 13, 48], we predicted the targets of 21 nt colonization-induced Si and Bd putative ck-sRNAs and assessed their expression after colonization. Of 16,003 unique Bd21-3 targets predicted for 412 induced 21 nt Si sRNAs, 49 were confirmed as downregulated at 4DPI. This represents 15.4% of all significantly changed genes in Bd21-3 during Si colonization. Some 89% of these transcripts are predicted as targets of Si sRNAs that are expressed exclusively in colonized tissue or with log2FC > 1. A representative set of sRNA-mRNA duplexes, chosen based on target identity and expression of target and sRNA, indicates that putative ck-sRNA targets in Bd21-3 are associated with transcription factor families, signaling pathways, and basal plant defense (Table 4, Additional file 2: Table S7).

Table 4 Examples of deduced duplexes of Serendipita indica sRNAs and their downregulated targets in Brachypodium distachyon

To assess whether Bd21-3 generates ck-sRNAs that potentially target Si genes, we searched for predicted targets of 329 Bd21-3 sRNAs (21 nt long) induced in Si-colonized roots. Of 3,019 predicted unique Si targets, 358 were confirmed as downregulated after colonization. This represents 12% of all significantly changed Si genes. 35% of the 358 Si transcripts are predicted to be targeted by Bd21-3 sRNAs exclusive to colonized tissue and an additional 27.6 % are targeted by sRNAs that are highly upregulated in colonized tissue (log2FC > 1). A set of sRNA-mRNA duplexes, selected with the same criteria as for the Si sRNA – Bd21-3 targets (Table 4), shows that predicted Bd21-3 ck-sRNAs have putative targets in Si which include proteins associated with nutrient acquisition, development of cell walls, hyphal networks, pathogenic fungal activities, fungal starvation, and signaling (Table 5, Additional file 2: Table S7). In order to confirm the expression of some of these sRNAs, we conducted stem-loop PCR and sRNA-specific sequencing on 10 Si and 10 Bd21-3-originating sRNAs from our Bd-Si sample and Table S7. All Si and all Bd21-3 sRNAs were amplified in the stem-loop PCR. To verify the nature of the amplification products, a subset of four SisRNAs and three BdsRNAs were cloned and sent for Sanger sequencing, confirming the expected sRNA sequences (Additional file 2: Figure S12, original gel pictures in Additional file 3 and Additional File 4, sequencing results in Additional file 2: Table S8). Thus, predicted targets of putative ck-sRNAs within this system imply another layer of expression control within the mutualistic interaction.

Table 5 Examples of deduced duplexes of Brachypodium distachyon sRNAs and their downregulated targets in Serendipita indica

Discussion

We established and studied the interaction between Brachypodium distachyon—a model Pooideae plant with shared synteny to major cereal crops—and Serendipita indica—a beneficial endophyte with an exceptionally large host range. This particular combination of traits adorning the Bd-Si interaction has great translational value towards filling in the gaps in knowledge about plant symbioses, especially their transcriptomic and sRNA expression profiles and the significance of RNAi. We show that Si colonizes Bd, resulting in shoot growth promotion, earlier flowering, and improved grain development. In comparison, earlier studies characterizing the interaction between Si and barley have demonstrated that fungal hyphae establish an interface with the root cell plasma membrane at an early colonization stage, followed by expansion of an extracellular hyphal network, intercellular growth, and intracellular penetration of cortical and rhizodermal cells [26]. Around 3 to 5 DPI, Si starts the switch from a biotrophic to a saprophytic lifestyle [26, 27]. Although this change involves intracellular hyphae extensively colonizing dead host cells and gradual digestion of cortical cell walls, the plant still benefits from the fungal presence. Consistent with these findings, our microscopic analyses confirmed proliferation of Si chlamydospores and both inter- and intracellular colonization of Bd21-3 cells in the root differentiation zone from 4 to 14 DPI. Detection of proliferating hyphae that were not wrapped in plasma membrane further suggests that Si is colonizing dead surface plant root cells at 4 DPI [25, 26].

Transcriptional changes detected during the Bd-Si interaction

To investigate colonization of Bd21-3 by Si, we analyzed Si DEGs in colonized vs. axenic mycelium samples. Gene ontology analysis indicated enrichment in genes involved in various metabolic and catalytic processes. DEGs with the greatest changes in expression play roles in plant cell wall degradation, carbohydrate metabolism and nutrient acquisition. These changes in nutritional reprogramming are to be expected, considering the different nutritional content that Si has accessible in planta vs. axenic CM plates and a more detailed look into the roles of the changed genes unveils a typical switch of fungal lifestyle. Examples of upregulated Si genes involved in cell wall degradation include a probable Pectate lyase, Endo-1,4-beta-xylanases, Cellulose 1,4-beta-cellobiosidase, and Rhamnogalacturonan acetylesterase. The genes encoding these hydrolytic enzymes, which have undergone expansion in the Si genome [49], are similar to those upregulated in Si during saprophytic growth on autoclaved barley roots at 3 DPI and 5 DPI [22]. PiAMT1, encoding a high affinity ammonium transporter also was upregulated (logFC = 3.35; padj < 0.0001). Other upregulated genes encode enzymes involved in carbohydrate metabolism, including probable glucosidases, glucanase, and galactosidase. These proteins may modulate glucose concentration, which then regulates expression of some cell-wall degrading enzymes [22, 50]. Some of the 174 putative effector protein-encoding genes also are differentially expressed during Si colonization of barley and Arabidopsis [22, 27]. Six of these proteins (Additional file 2: Table S4) contain the Si-specific DELD domain, which suggests that Si utilizes a common protein effector arsenal to colonize various hosts. Considering highly downregulated Si DEGs, several encoded proteins are associated with adaption to nutrient availability (Accumulation of dyads protein 2, ADY2; Succinate-fumarate transporter) and nutrient acquisition (Acyl-CoA dehydrogenase; Carnitine acetyltransferase, CRAT; Phenylalanine ammonia-lyase, PAL [51, 52];). Their reduced expression suggests ample nutrient availability at 4 DPI. Given the similarities in the Si transcriptome during colonization of Bd and barley, we propose that these fungal-plant interactions follow a pattern, and that by 4 DPI, a network of plant-endophyte communication cues has initiated a tightly controlled transcriptional program, leading to a shift from biotrophic to saprophytic growth.

Roots of Bd21-3 plants also displayed substantial transcriptional reprogramming following Si colonization. Gene ontology term analysis indicated enrichment in genes involved in catalytic and oxidoreduction-associated processes. Bd21-3 DEGs exhibiting the greatest changes in expression between colonized and non-colonized plants are related to stress-response, defense, and plant development. Of the downregulated Bd21-3 genes, several encode proteins commonly associated with stress responses, including a peroxidase, a wound-induced protein, and a putative protein kinase. Additionally, members of the Heat-shock protein gene family [53] are commonly induced in Bd during abiotic stress and members of the Abscisic acid/water deficit stress (ABA/WDS)-induced protein and the Rare cold inducible (RCI2) gene families enhance abiotic stress tolerance in various plant species [54, 55]. Circadian clock and flowering regulation genes such as Pseudo-response regulator 7 (PRR7), Cold regulated protein 27 and Constans-like protein (CO8), also are downregulated during Si colonization. While members of the PRR and CO protein families work together to control flowering time [56, 57], any influence on early flowering in Si-colonized Bd21-3 is unclear. Circadian clock-associated genes also regulate lateral root development in Arabidopsis [58]; whether Si-induced changes in their expression influence root growth is unknown. Other downregulated development-associated DEGs include Fantastic four meristem regulator (FAF), which regulates shoot and root development [59], and putative Sulfoquinovosyltransferase (SQD2), which modulates seed setting and tiller development in rice [60]. Finally, several downregulated DEGs encode transcription factors, including MYB-related, GRAS, and bZIP.

In comparison, many of the upregulated Bd21-3 DEGs are associated with immune responses. Examples include genes encoding leucine-rich repeat (LRR) protein, a WRKY transcription factor, and thaumatin family protein. Increased expression of the defense gene Pathogenesis-related protein 1 (PR1) was similarly and transiently reported in Si-colonized Arabidopsis roots [61]. Upregulation of Glutathione S-transferase (GST) is consistent with the increased antioxidant capacity of Si-colonized plants, which provides protection against attack by necrotrophic pathogens [21, 62]. The upregulation of genes in other hormonal networks (PGP-like Phosphoglycoprotein auxin transporter) and redox processes (Multicopper oxidase) further suggests that Si colonization affects a range of signaling pathways.

Bd miRNAs detected in the Bd-Si sample

The role of miRNAs as regulators of gene expression in the Sebacinalean symbiosis is largely unexplored. One report showed that Si induces growth promotion-associated miRNAs in Oncidium orchid roots [63]. Analysis of putative endogenous Bd21-3 sRNAs expressed during Si colonization identified 16 miRNAs. Some of them have known targets in transcription factors associated with plant growth and development. For example, the bdi-MIR166 family targets mRNAs encoding Homeobox domain-leucine zipper transcription factors [64]. In Arabidopsis, repression of these transcription factors by the miR165/166 family modulates root growth, maintenance of the shoot apical meristem, and development of leaf polarity [65]. Plant-specific transcription factors encoded by Squamosa promoter-binding protein-like (SPL) genes are the presumed targets of bdi-MIR156 and bdi-MIR529 [66]. In Arabidopsis, miR156-mediated downregulation of SPLs modulates developmental timing, lateral root development, branching, and leaf morphology [65]. Members of the MYB superfamily of transcription factors, which regulate many aspects of development, are the predicted targets of bdi-MIR159 [34]. Interestingly, miRNAs belonging to the miR159 and miR166 families in cotton are known ck-sRNAs that target virulence genes in Verticillium dahliae [13].

Other miRNAs identified in Bd21-3 include bdi-MIR168, predicted to target AGO1 [64], and two miRNAs that regulate nutrition: bdi-MIR399, which is upregulated in Bd by phosphate starvation [64, 67], whereas bdi-MIR408 influences copper levels [34, 68]. Additionally, bdi-MIR408 (BdsRNA 10) has predicted ck targets in three Si transcripts: CCA72944, CCA72668, and CCA74115. Since various targets were predicted for bdi-MIR167 [34, 68] and no target was predicted for bdi-MIR9481, their endogenous functions in Bd are unclear. Interestingly, the miRNA families identified in our analysis, except bdi-MIR9481, also were detected in Si-colonized Oncidium [63]. Thus, this group of miRNAs may play an important role in reprogramming plant cells during Sebacinalean symbiosis establishment.

Putative Si and Bd21-3 ck-sRNAs and their predicted targets

To date, cross-kingdom RNAi has been demonstrated in pathogenic plant-fungal interactions [12, 13, 69], and while there are promising indications for its presence during plant-mycorrhiza interactions [36], whether it occurs in Si-plant associations is unknown. To investigate this possibility, we predicted targets for 21 nt putative ck-sRNAs from Si and Bd21-3 and confirmed their downregulation during colonization. This analysis uncovered 358 downregulated Si transcripts that are the predicted targets of 228 unique Bd21-3 sRNAs. Cross-kingdom RNAi-mediated downregulation of these targets might allow Bd21-3 to modulate Si growth during colonization. For example, PAL, Acetyl-CoA synthetase, Carnitine acetyl transferase, Isocitrate lyase, and Acyl-CoA dehydrogenase, which are targeted by BdsRNA 1, BdsRNA 2, BdsRNA 3, BdsRNA 6, and BdsRNA 16 (Table 5), are involved in fungal nutrient acquisition [22, 55, 70, 71]. Genes with important homologs in pathogenic fungi also are predicted targets, including Subtilisin-like serine protease (BdsRNA 18) [72], Alcohol dehydrogenase 1 (BdsRNA 7) [73], and Phosphoprotein phosphatase 2C (BdsRNA 9) [74]. Targeting of Hyphal anastomosis-2 (HAM-2) by BdsRNA 8 and BdsRNA 4 may provide another mechanism for controlling fungal growth, as HAM-2 is required for hyphal fusion in N. crassa [75]. Similarly, targeting of Chitinase (BdsRNA 14) may help control Si growth.

Concurrently, we identified 49 downregulated Bd21-3 mRNAs that are the predicted targets of 63 unique Si-generated ck-sRNAs. Downregulation of these target genes via cross-kingdom RNAi might facilitate Si growth during colonization. For example, Mannose-binding lectin (targeted by SisRNA 18) belongs to a family of defense-related genes whose products trigger immune responses following pathogen recognition [76]. SisRNA 8 and SisRNA 15 target a protein kinase domain/LRR gene (BdiBd21-3.4G0303000.1) that may belong to the LRR receptor kinase family, which regulates defense and developmental-related processes [77]. Transcripts encoding serine-carboxypeptidase-like (SCPL) proteins BdiBd21-3.2G0440200.1 and BdiBd21-3.1G0411900.1 (targeted by SisRNA 1 and SisRNA 15) are associated with defense against (a)biotic stresses in monocots [78]. Members of various transcription factors families also were identified as predicted targets (MYB by SisRNA 16, bZIP by SisRNA 2, and GRAS by SisRNA 10). These families are associated with (a)biotic stress responses, as well as plant growth and development [79,80,81]. Lastly, transcripts for proteins involved in circadian clock and flowering regulation (BdiBd21-3.1G0887100.1 and BdiBd21-3.3G0264400.1 [56];) are the presumed targets of SisRNA 1 and SisRNA 19. Together, these findings suggest that Si-derived ck-sRNAs may promote fungal colonization by targeting signaling processes associated with plant development and responses to (a)biotic stresses.

In combination with earlier studies on Bd RNAi proteins [35] and Bd interaction with the pathogen Magnaporthe oryzae [82], the in silico analyses presented here suggest that Si and Bd contain functional RNAi components and that both organisms generate ck-sRNAs, which potentially modulate this mutualistic interaction. However, further studies are necessary to validate cross-kingdom RNAi in a Sebacinalean symbiosis. Namely, degradome analysis is needed to confirm target degradation and evidence that Bd21-3 and Si AGOs associate with sRNAs expressed by the interacting organism is necessary for confirmation of cross-kingdom RNAi.

Conclusions

We report that Bd21-3 and Si form a mutualistic symbiosis with a promoting effect on plant yield and development, accompanied by changes in gene expression in both organisms, including putative protein Si effectors and RNAi-related genes. sRNA profiles of both organisms also changed, indicating that this model system will provide important insights into the multiple layers of regulation and interaction between beneficial fungi and cereal hosts. Within the broader scope of plant-mutualist interactions, we show that detection of putative RNAi-involved sRNAs in an interaction highly benefits from simultaneous transcriptome analysis and indicate an involvement of sRNA-based regulation in defense responses, nutritional reprogramming, and colonization maintenance. Alongside other experimental approaches in plant-microbe interactions (eg. sRNA uptake studies [83]), developing a deeper understanding of the communication mechanisms that modulate mutualistic interactions is highly relevant for establishing robust growth promotion and protection strategies in crops.

Methods

Bd and Si cultivation and inoculation

The seeds of Brachypodium distachyon (Bd) line Bd21-3 (gift from R. Sibout, INRA Versailles, France) were surface sterilized (3% active chlorine, sodium hypochlorite solution) for 15 min, washed three times, and placed on half-strength MS [84] medium in dark at 4 °C for 2 days and then 7 days at 24 °C and 16 h light/8 h dark cycle (47 μmol m− 2 s− 1 photon flux density). Serendipita indica (Si) (IPAZ-11827, Institute of Phytopathology, Giessen, Germany) was grown on complete media plates (CM [85]) at 23 °C in dark for 4–5 weeks.

For inoculation, Si mycelium was collected in 0.002% aqueous Tween 20 solution, filtered (Miracloth, Calbiochem), and pelleted by centrifugation (10 min/4000 rpm/20 °C) twice. Chlamydospore concentration of 5 × 105 conidia ml−1 in 0.002% Tween 20 solution was used to inoculate 7-day-old plant seedlings for 2–3 h. Control plants were mock treated with the 0.002% Tween 20 solution for the same time. Grain yield analyses were done on mature plants grown on soil (F-E type LD 80, Fruhstorfer Erde, Germany) under 16 h light (160 μmol m−2 s−1, 22 °C) and 8 h dark (18 °C) conditions at 60% relative humidity for 1 month, and then greenhouse conditions until seed maturity. Number of spikelets was assessed after 2 months. Shoot biomass was assessed 3 weeks after inoculation of seedlings grown on a mixture (2:1, v/v) of vermiculite (Deutsche Vermiculite GmbH) and oil dri (oil binder Typ III R Coarse grain, Damolin, Mettmann, Germany) under comparable conditions as for grain yield, and fertilized every 3 days with an aqueous solution of Wuxal Super NPK-8/8/6 (1:103 v/v; Haug, Düsseldorf, Germany). Samples for RNA-seq, RT-qPCR, stem-loop PCR, and microscopy were also grown under these conditions. To assess growth promotion in Si inoculated Bd relative to the control, we used the pairwise t test or the Mann-Whitney-Wilcoxon test on each of the three repetitions of experiments, after checking for normality and homogenous variances. Benjamini-Hochberg correction for multiple testing was used to correct the p values and the significance asterisks were assigned to the average p-value as follows: * for p ≤ 0.05, ** for p ≤ 0.001, and *** for p ≤ 0.0001.

Microscopy

Following Si inoculation, one-week-old seedlings were grown on plastic mesh (~ 90 μm) placed over half-strength MS medium or on vermiculite/oil dri prior to assessing root colonization. Si was visualized with the chitin-specific dye WGA-AF 488 (wheat germ agglutinin; Molecular Probes, Karlsruhe, Germany), as described in Deshmukh et al. (2006) [26], with boiling in KOH (10%) for 30 s, prior to incubation in phosphate-buffered saline (PBS, pH 7.4). Root cells were visualized by incubating with propidium iodide (10 μg ml−1) for 10 min and washing with sterile water. Confocal laser scanning microscopy was done (TCS SP8 microscope, Leica, Bensheim, Germany) and the Leica LAS X software was utilized for visualization and maximum (z-stack) projections.

Resequencing, assembly, and annotation of the Si genome

The MasterPure Yeast DNA Purification Kit (Epicentre, Illumina) was used to extract genomic DNA from 4-week-old axenic Si cultures. The Si genome was resequenced, assembled [86], and annotated as described [87], whereby a MinION sequencing library was prepared using the Nanopore Rapid DNA Sequencing kit. Sequencing was performed on an Oxford Nanopore MinION Mk1b sequencer using a R9.5 flow cell. Additionally, sequencing of an Illumina Nextra XT library was performed on the MiSeq platform (Illumina; 2 × 300 bp paired-end sequencing, v3 chemistry). Adapters and low-quality reads were removed by an in-house software pipeline prior to polishing [88]. MinKNOW (v1.13.1, Oxford Nanopore Technologies) was used to control the run with the 48 h sequencing run protocol, and base calling was performed offline using albacore (v2.3.1, Oxford Nanopore Technologies). The assembly was performed using Canu v1.6 ([89], default settings). The resulting contigs were polished with Illumina short read data using Pilon [90] for eight iterative cycles. BWA-MEM [91] was used for read mapping in the first four iterations and Bowtie2 v2.3.2 [92] in the second set. Gene prediction was performed with GeneMark-ES 4.3.3. ([38], default settings). Predicted genes were functionally annotated using a modified version of the genome annotation platform GenDB 2.0 [39] for eukaryotic genomes [40]. RNAi-associated proteins were predicted by searching the proteome [22] for typical domain structure and highest homology to Neurospora crassa RNAi proteins (NC12 genome assembly [93]). A modified version of the pipeline from Rafiqi et al. (2013) [41] was used to predict protein effectors. After identifying proteins with signal peptides (signalp-4.1 [94]), those predicted as transmembrane helix proteins (tmhmm [95]), mitochondrial proteins (target-1.1 [96]), and cell wall hydrolysis-associated proteins were removed. For comparative analysis of the Si (Si-2020 and DSM11827 ASM31354 v.1 [22]), Serendipita vermifera [97] and Laccaria bicolor [98] genomes, software platform EDGAR 2.3 [99] was used.

RNA extraction, library preparation, and mRNA/sRNA sequencing

Roots inoculated with Si (Bd-Si) or mock-inoculated (Bd-C), as described above, were grown for 4 days and pooled (three roots per sample). Si mycelium and spores were collected from 4-week-old axenic cultures grown on CM medium. All samples in triplicates were shock frozen, stored at − 80 °C, and ground in liquid N2. Total RNA was isolated using the ZymoBIOMICS RNA Mini Kit (Zymo Research, USA), quantified with DropSense16/Xpose (BIOKÉ, Netherlands), and analyzed with an Agilent 2100 Bioanalyzer Nano Chip (Agilent, Germany). RNA Clean and Concentrator 25 and 5 kits (Zymo Research) were utilized to separate total RNA into fractions: 17–200 nt and > 200 nt. 1.5 μg of the larger fractions were processed for mRNA library preparation (TruSeq Stranded mRNA protocol, Illumina, USA). Fragment Analyzer Automated CE System (Advanced Analytical Technologies, Austria) determined the quality of the generated polyA mRNA libraries. Quantity and quality of the smaller RNA fractions were assessed with the Qubit fluorometer (Invitrogen, Germany) and Agilent 2100 Bioanalyzer Pico Chip. sRNA library preparation was done with 50 ng of RNA (TruSeq Small RNA Library Prep, Illumina) and size selection with the BluePippin (Sage Science, USA) for fragments between 140 and 160 nt (15–35 nt without adapters) applied. Sequencing was accomplished on the Illumina HiSeq 1500.

Transcriptome analysis

Raw reads from mRNA sequencing [100] were submitted to quality check using FastQC [101] and aligned to the Bd21-3 v1.1 (DOE-JGI, http://phytozome.jgi.doe.gov/ [102]) or resequenced Si (Si-2020) genomes with HISAT2 [103]. An intron length of 20–2000 nt was allowed for Si [104] and 20–10,000 nt for Bd21-3 [105]. The reads were counted using HTSeq-count [106], differential gene expression was performed with DESeq2 [107], and gene enrichment analysis with AgriGO v.2 [108], with reference Bd21 setting for Bd21-3 (Bd 21 synonyms) and a customized background for Si. Volcano plots were generated using plotly [109] and ggplot2 [110] R [111] libraries. Gene descriptions were obtained from the organism annotations or Blast2GO [112].

sRNA analysis and prediction of putative endogenous and ck-sRNAs

Raw reads from sRNA sequencing [113] were submitted to FastQC [101] and adapter trimming [114]. Bowtie [115] was used for alignment as detailed in Additional file 2: Fig. S6. The resequenced Si genome was used for alignments of fungal origin. tRNA/rRNA sequences were downloaded from RNAcentral ([116], EMBL-EBI). Putative endogenous sRNA reads were submitted to ShortStack [117]. For filtering putative ck-sRNAs, a previously established pipeline [118] was utilized. Reads were normalized to the total number of mapped reads for a single genome and reads per million (RPM) and log2 (colonized/mock-treated) values calculated. Thus, a sRNA read was selected as a putative ck-sRNA if it was present exclusively or at a higher quantity (i.e., induced) in the colonized vs. control sample. Putative ck-sRNAs were submitted to psRNAtarget [119]. Since the separation of sRNA and mRNA fractions from each biological sample was facilitated, we checked for downregulation of mRNAs corresponding to predicted sRNA target genes within the DEGs. Transcriptomes used for these predictions were Bd 21-3 v1.1 (DOE-JGI, http://phytozome.jgi.doe.gov/ [102]) and Si DSM11827 ASM31354 v.1 [22]. Venn diagrams were generated using the VennDiagram R package [120].

Quantitative real-time PCR and stem-loop PCR for validation of sequencing results

To validate gene expression detected in the sequencing, we used quantitative real-time PCR (qRT-PCR). RNA extraction from mock treated and Si inoculated Bd21-3 roots, as well as Si axenic cultures, under the same conditions as explained above for the sequencing, was done with TRIzol (Thermo Fisher Scientific, Waltham, MA, USA), cDNA synthesized using qScriptTM cDNA kit (Quantabio, Beverly, MA, USA) and 10 ng of cDNA used as template in the QuantStudio 5 Real-Time PCR system (Applied Biosystems), with SYBR® green JumpStart Taq ReadyMix (Sigma-Aldrich, St. Louis, MO, USA). Each sample had three technical replicates. Primers used for these amplifications are listed in Table S9 (Additional file 2: Table S9). Transcript levels were calculated using the 2–ΔΔCt method [121], relatively to BdUbi4-3 for Bd21-3 and Si ITS sequence for Si.

For the identification of sRNAs in the interaction of Si with Bd21-3 stem-loop RT-PCR was employed [122]. cDNA was synthesized from DNase I-treated total RNA extracted from Si axenic culture or inoculated Bd21-3 roots. The folding of the hairpin primer was performed according to Kramer (2011) [123]. For each stem-loop reaction, six hairpin primers were multiplexed in a 20-μL reaction using the Revertaid RT enzyme according to the manufacturer’s instructions (Thermo Scientific). For primer annealing, the reaction was incubated for 30 min at 16 °C followed by an extension step at 42 °C for 30 min. The reaction was stopped at 85 °C for 5 min. cDNA was stored at − 80 °C until further use. Endpoint PCR was performed using an universal stem-loop primer and specific sRNA primer (Additional file 2: Table S10) under the following conditions: initial denaturation at 95 °C for 5 min followed by 35 cycles: 95 °C for 30 s, primer annealing at 60 °C for 30 s, and extension at 72 °C for 30 s. PCR products were separated by gel electrophoresis on a 2% (w/v) agarose gel. To obtain sequence information of the amplified sRNAs of the stem-loop reaction, PCR products were purified and cloned into the pGEM®-T Easy Vector Systems (Promega, Madison, WI, USA) following the manufacturer’s instructions. From each cloned sRNA, six colonies were further analyzed by Sanger sequencing using a M13 reverse primer.

Availability of data and materials

All data generated or analyzed during this study are included in this published article, its supplementary information files (Additional files 1, 2, 3 and 4) and publicly available repositories. mRNA and sRNA-seq of Brachypodium distachyon roots inoculated with Serendipita indica datasets are available in the ArrayExpress database at EMBL-EBI (www.ebi.ac.uk/arrayexpress) under accession numbers E-MTAB-10649 and E-MTAB-10650, respectively. The genome assembly data for Serendipita indica resequencing study is available in the European Nucleotide Archive (ENA) at EMBL-EBI under accession numbers: assembly GCA_910890315 and study PRJEB45884 (https://www.ebi.ac.uk/ena/browser/view/PRJEB45884).

Additional publicly available datasets have been accessed from the corresponding European Nucleotide Archive accessions and EnsemblFungi: Serendipita indica 2011 genome (GCA_000313545.1) [22], NC12 of Neurospora crassa (GCA_000182925.2) [93]. Serendipita vermifera subsp. bescii [97] has been accessed from JGI MycoCosm (https://mycocosm.jgi.doe.gov/Sebvebe1/Sebvebe1.home.html), as was the Laccaria bicolor genome ([98], https://mycocosm.jgi.doe.gov/Lacbi2/Lacbi2.home.html).

The Bd21-3 v1.1 genome [102] was accessed under early access conditions and these sequence data were produced by the US Department of Energy Joint Genome Institute - DOE-JGI, http://phytozome.jgi.doe.gov/.

Abbreviations

AGO:

Argonaute

Bd :

Brachypodium distachyon

Bd-C:

Mock-inoculated Bd root sample

Bd-Si:

Si-inoculated Bd root sample

BdsRNA:

Brachypodium distachyon-derived sRNA

ck-sRNA:

Cross-kingdom sRNA

DCL:

Dicer-like

DEG:

Differentially expressed gene

DPI:

Day(s) post inoculation

FC:

Fold change

miRNA:

MicroRNA

nt:

Nucleotide

RdRP:

RNA-dependent RNA polymerase

RNAi:

RNA interference

Si :

Serendipita indica

Si-ax:

Si axenic culture sample

SisRNA:

Serendipita indica-derived sRNA

sRNA:

Small RNA

References

  1. 1.

    Marschner P, Yang C-H, Lieberei R, Crowley DE. Soil and plant specific effects on bacterial community composition in the rhizosphere. Soil Biol Biochem. 2001;33(11):1437–45. https://doi.org/10.1016/S0038-0717(01)00052-9.

    CAS  Article  Google Scholar 

  2. 2.

    Lagunas B, Schäfer P, Gifford ML. Housing helpful invaders: the evolutionary and molecular architecture underlying plant root-mutualist microbe interactions. J Exp Bot. 2015;66(8):2177–86. https://doi.org/10.1093/jxb/erv038.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  3. 3.

    Kogel KH, Franken P, Hückelhoven R. Endophyte or parasite – what decides? Curr Opin Plant Biol. 2006;9(4):358–63. https://doi.org/10.1016/j.pbi.2006.05.001.

    Article  PubMed  Google Scholar 

  4. 4.

    Fesel PH, Zuccaro A. Dissecting endophytic lifestyle along the parasitism/mutualism continuum in Arabidopsis. Curr Opin Microbiol. 2016;32:103–12. https://doi.org/10.1016/j.mib.2016.05.008.

    Article  PubMed  Google Scholar 

  5. 5.

    Kloppholz S, Kuhn H, Requena N. A secreted fungal effector of Glomus intraradices promotes symbiotic biotrophy. Curr Biol. 2011;21(14):1204–9. https://doi.org/10.1016/j.cub.2011.06.044.

    CAS  Article  PubMed  Google Scholar 

  6. 6.

    Doehlemann G, Requena N, Schaefer P, Brunner F, O’Connell R, Parker JE. Reprogramming of plant cells by filamentous plant-colonizing microbes. New Phytol. 2014;204(4):803–14. https://doi.org/10.1111/nph.12938.

    Article  PubMed  Google Scholar 

  7. 7.

    Vishwakarma K, Kumar N, Shandilya C, Mohapatra S, Bhayana S, Varma A. Revisiting plant–microbe interactions and microbial consortia application for enhancing sustainable agriculture: a review. Front Microbiol. 2020;11:3195.

    Google Scholar 

  8. 8.

    Baulcombe D. RNA silencing in plants. Nature. 2004;431(7006):356–63. https://doi.org/10.1038/nature02874.

    CAS  Article  PubMed  Google Scholar 

  9. 9.

    Fang X, Qi Y. RNAi in plants: an Argonaute-centered view. Plant Cell. 2016;28(2):272–85. https://doi.org/10.1105/tpc.15.00920.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  10. 10.

    Uhse S, Djamei A. Effectors of plant-colonizing fungi and beyond. PLOS Pathog. 2018;14(6):e1006992. https://doi.org/10.1371/journal.ppat.1006992.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  11. 11.

    Akum FN, Steinbrenner J, Biedenkopf D, Imani J, Kogel KH. The Piriformospora indica effector PIIN_08944 promotes the mutualistic Sebacinalean symbiosis. Front Plant Sci. 2015;6:1–12.

    Article  Google Scholar 

  12. 12.

    Weiberg A, Wang M, Lin F-M, Zhao H, Zhang Z, Kaloshian I, et al. Fungal small RNAs suppress plant immunity by hijacking host RNA interference pathways. Science (80-). 2013;342(6154):118–23.

    CAS  Article  Google Scholar 

  13. 13.

    Zhang T, Zhao YL, Zhao JH, Wang S, Jin Y, Chen ZQ, et al. Cotton plants export microRNAs to inhibit virulence gene expression in a fungal pathogen. Nat Plants. 2016;2(10):1–6.

    Article  Google Scholar 

  14. 14.

    Kellogg EA. Brachypodium distachyon as a genetic model system. Annu Rev Genet. 2015;49(1):1–20. https://doi.org/10.1146/annurev-genet-112414-055135.

    CAS  Article  PubMed  Google Scholar 

  15. 15.

    Qiang X, Weiss M, Kogel KH, Schäfer P. Piriformospora indica-a mutualistic basidiomycete with an exceptionally large plant host range. Mol Plant Pathol. 2012;13(5):508–18. https://doi.org/10.1111/j.1364-3703.2011.00764.x.

    CAS  Article  PubMed  Google Scholar 

  16. 16.

    Stein E, Molitor A, Kogel KH, Waller F. Systemic resistance in Arabidopsis conferred by the mycorrhizal fungus Piriformospora indica requires jasmonic acid signaling and the cytoplasmic function of NPR1. Plant Cell Physiol. 2008;49(11):1747–51. https://doi.org/10.1093/pcp/pcn147.

    Article  PubMed  Google Scholar 

  17. 17.

    Deshmukh SD, Kogel KH. Piriformospora indica protects barley from root rot caused by Fusarium graminearum. J Plant Dis Prot. 2007;114(6):263–8. https://doi.org/10.1007/BF03356227.

    Article  Google Scholar 

  18. 18.

    Trzewik A, Maciorowski R, Klocke E, Orlikowska T. The influence of Piriformospora indica on the resistance of two rhododendron cultivars to Phytophthora cinnamomi and P. plurivora. Biol Control. 2020;140:104121.

    CAS  Article  Google Scholar 

  19. 19.

    Fakhro A, Andrade-Linares DR, von Bargen S, Bandte M, Büttner C, Grosch R, et al. Impact of Piriformospora indica on tomato growth and on interaction with fungal and viral pathogens. Mycorrhiza. 2010;20(3):191–200. https://doi.org/10.1007/s00572-009-0279-5.

    Article  PubMed  Google Scholar 

  20. 20.

    Baltruschat H, Fodor J, Harrach BD, Niemczyk E, Barna B, Gullner G, et al. Salt tolerance of barley induced by the root endophyte Piriformospora indica is associated with a strong increase in antioxidants. New Phytol. 2008;180(2):501–10. https://doi.org/10.1111/j.1469-8137.2008.02583.x.

    CAS  Article  PubMed  Google Scholar 

  21. 21.

    Waller F, Achatz B, Baltruschat H, Fodor J, Becker K, Fischer M, et al. The endophytic fungus Piriformospora indica reprograms barley to salt-stress tolerance, disease resistance, and higher yield. Proc Natl Acad Sci U S A. 2005;102(38):13386–91. https://doi.org/10.1073/pnas.0504423102.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  22. 22.

    Zuccaro A, Lahrmann U, Güldener U, Langen G, Pfiffi S, Biedenkopf D, et al. Endophytic life strategies decoded by genome and transcriptome analyses of the mutualistic root symbiont Piriformospora indica. PLoS Pathog. 2011;7(10) Serendipita indica genome assembly. 2011. GCA accession number: GCA_000313545.1 (https://www.ebi.ac.uk/ena/browser/view/GCA_000313545.1).

  23. 23.

    Zuccaro A, Basiewicz M, Zurawska M, Biedenkopf D, Kogel KH. Karyotype analysis, genome organization, and stable genetic transformation of the root colonizing fungus Piriformospora indica. Fungal Genet Biol. 2009;46(8):543–50. https://doi.org/10.1016/j.fgb.2009.03.009.

    CAS  Article  PubMed  Google Scholar 

  24. 24.

    Schäfer P, Pfiffi S, Voll LM, Zajic D, Chandler PM, Waller F, et al. Manipulation of plant innate immunity and gibberellin as factor of compatibility in the mutualistic association of barley roots with Piriformospora indica. Plant J. 2009;59(3):461–74. https://doi.org/10.1111/j.1365-313X.2009.03887.x.

    CAS  Article  PubMed  Google Scholar 

  25. 25.

    Jacobs S, Zechmann B, Molitor A, Trujillo M, Petutschnig E, Likpa V, et al. Broad-spectrum suppression of innate immunity is required for colonization of Arabidopsis roots by the fungus Piriformospora indica. Plant Physiol. 2011;156(2):726–40. https://doi.org/10.1104/pp.111.176446.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  26. 26.

    Deshmukh S, Hückelhoven R, Schäfer P, Imani J, Sharma M, Weiss M, et al. The root endophytic fungus Piriformospora indica requires host cell death for proliferation during mutualistic symbiosis with barley. Proc Natl Acad Sci U S A. 2006;103(49):18450–7. https://doi.org/10.1073/pnas.0605697103.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  27. 27.

    Lahrmann U, Ding Y, Banhara A, Rath M, Hajirezaei MR, Döhlemann S, et al. Host-related metabolic cues affect colonization strategies of a root endophyte. Proc Natl Acad Sci U S A. 2013;110(34):13965–70. https://doi.org/10.1073/pnas.1301653110.

    Article  PubMed  PubMed Central  Google Scholar 

  28. 28.

    Scholthof K-BG, Irigoyen S, Catalan P, Mandadi KK. Brachypodium: a monocot grass model genus for plant biology. Plant Cell. 2018;30(8):1673–94. https://doi.org/10.1105/tpc.18.00083.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  29. 29.

    Vogel J, Hill T. High-efficiency Agrobacterium-mediated transformation of Brachypodium distachyon inbred line Bd21-3. Plant Cell Rep. 2008;27(3):471–8. https://doi.org/10.1007/s00299-007-0472-y.

    CAS  Article  PubMed  Google Scholar 

  30. 30.

    Vogel JP, Garvin DF, Mockler TC, Schmutz J, Rokhsar D, Bevan MW, et al. Genome sequencing and analysis of the model grass Brachypodium distachyon. Nature. 2010;463(7282):763–8.

    CAS  Article  Google Scholar 

  31. 31.

    Hsia MM, O’Malley R, Cartwright A, Nieu R, Gordon SP, Kelly S, et al. Sequencing and functional validation of the JGI Brachypodium distachyon T-DNA collection. Plant J. 2017;91(3):361–70. https://doi.org/10.1111/tpj.13582.

    CAS  Article  PubMed  Google Scholar 

  32. 32.

    Fitzgerald TL, Powell JJ, Schneebeli K, Hsia MM, Gardiner DM, Bragg JN, et al. Brachypodium as an emerging model for cereal-pathogen interactions. Ann Bot. 2015;115(5):717–31. https://doi.org/10.1093/aob/mcv010.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  33. 33.

    Zhang J, Xu Y, Huan Q, Chong K. Deep sequencing of Brachypodium small RNAs at the global genome level identifies microRNAs involved in cold stress response. BMC Genomics. 2009;10(1):449. https://doi.org/10.1186/1471-2164-10-449.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  34. 34.

    Bertolini E, Verelst W, Horner DS, Gianfranceschi L, Piccolo V, Inzé D, et al. Addressing the role of microRNAs in reprogramming leaf growth during drought stress in Brachypodium distachyon. Mol Plant. 2013;6(2):423–43. https://doi.org/10.1093/mp/sss160.

    CAS  Article  PubMed  Google Scholar 

  35. 35.

    Šečić E, Zanini S, Kogel KH. Further elucidation of the Argonaute and Dicer protein families in the model grass species Brachypodium distachyon. Front Plant Sci. 2019;10:1332. https://doi.org/10.3389/fpls.2019.01332.

    Article  PubMed  PubMed Central  Google Scholar 

  36. 36.

    Silvestri A, Fiorilli V, Miozzi L, Accotto GP, Turina M, Lanfranco L. In silico analysis of fungal small RNA accumulation reveals putative plant mRNA targets in the symbiosis between an arbuscular mycorrhizal fungus and its host plant. BMC Genomics. 2019;20(1):1–18.

    Article  Google Scholar 

  37. 37.

    Ren B, Wang X, Duan J, Ma J. Rhizobial tRNA-derived small RNAs are signal molecules regulating plant nodulation. Science (80- ). 2019;365:919–22.

    CAS  Article  Google Scholar 

  38. 38.

    Ter-Hovhannisyan V, Lomsadze A, Chernoff YO, Borodovsky M. Gene prediction in novel fungal genomes using an ab initio algorithm with unsupervised training. Genome Res. 2008;18(12):1979–90. https://doi.org/10.1101/gr.081612.108.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  39. 39.

    Meyer F, Goesmann A, McHardy AC, Bartels D, Bekel T, Clausen J, et al. GenDB--an open source genome annotation system for prokaryote genomes. Nucleic Acids Res. 2003;31(8):2187–95. https://doi.org/10.1093/nar/gkg312.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  40. 40.

    Rupp O, Becker J, Brinkrolf K, Timmermann C, Borth N, Pühler A, et al. Construction of a public CHO cell line transcript database using versatile bioinformatics analysis pipelines. PLoS One. 2014;9(1):e85568. https://doi.org/10.1371/journal.pone.0085568.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  41. 41.

    Rafiqi M, Jelonek L, Akum NF, Zhang F, Kogel KH. Effector candidates in the secretome of Piriformospora indica, a ubiquitous plant-associated fungus. Front Plant Sci. 2013;4:228.

    Article  Google Scholar 

  42. 42.

    Billmyre RB, Calo S, Feretzaki M, Wang X, Heitman J. RNAi function, diversity, and loss in the fungal kingdom. Chromosome Res. 2013;21(6–7):561–72. https://doi.org/10.1007/s10577-013-9388-2.

    CAS  Article  PubMed  Google Scholar 

  43. 43.

    Dang Y, Yang Q, Xue Z, Liu Y. RNA interference in fungi: pathways, functions, and applications. Eukaryot Cell. 2011;10(9):1148 LP–1155.

    Article  Google Scholar 

  44. 44.

    Lau P-W, Guiley KZ, De N, Potter CS, Carragher B, MacRae IJ. The molecular architecture of human Dicer. Nat Struct Mol Biol. 2012;19(4):436–40. https://doi.org/10.1038/nsmb.2268.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  45. 45.

    Cenik ES, Zamore PD. Argonaute proteins. Curr Biol. 2011;21(12):R446–9. https://doi.org/10.1016/j.cub.2011.05.020.

    CAS  Article  PubMed  Google Scholar 

  46. 46.

    Meyers BC, Axtell MJ. MicroRNAs in plants: key findings from the early years. Plant Cell. 2019;31(6):1206–7. https://doi.org/10.1105/tpc.19.00310.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  47. 47.

    Mi S, Cai T, Hu Y, Chen Y, Hodges E, Ni F, et al. Sorting of small RNAs into Arabidopsis Argonaute complexes is directed by the 5’ terminal nucleotide. Cell. 2008;133(1):116–27. https://doi.org/10.1016/j.cell.2008.02.034.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  48. 48.

    Cai Q, Qiao L, Wang M, He B, Lin F-M, Palmquist J, et al. Plants send small RNAs in extracellular vesicles to fungal pathogen to silence virulence genes. Science (80- ). 2018;360:eaar4142.

    Article  Google Scholar 

  49. 49.

    Lahrmann U, Strehmel N, Langen G, Frerigmann H, Leson L, Ding Y, et al. Mutualistic root endophytism is not associated with the reduction of saprotrophic traits and requires a noncompromised plant innate immunity. New Phytol. 2015;207(3):841–57. https://doi.org/10.1111/nph.13411.

    CAS  Article  PubMed  Google Scholar 

  50. 50.

    Tonukari NJ, Scott-Craig JS, Waltonb JD. The Cochliobolus carbonum SNF1 gene is required for cell wall–degrading enzyme expression and virulence on maize. Plant Cell. 2000;12(2):237–47.

    CAS  PubMed  PubMed Central  Google Scholar 

  51. 51.

    Zhou H, Lorenz MC. Carnitine acetyltransferases are required for growth on non-fermentable carbon sources but not for pathogenesis in Candida albicans. Microbiology. 2008;154(2):500–9. https://doi.org/10.1099/mic.0.2007/014555-0.

    CAS  Article  PubMed  Google Scholar 

  52. 52.

    Hyun MW, Yun YH, Kim JY, Kim SH. Fungal and plant phenylalanine ammonia-lyase. Mycobiology. 2011;39(4):257–65. https://doi.org/10.5941/MYCO.2011.39.4.257.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  53. 53.

    Zhang M, Shen Z, Meng G, Lu Y, Wang Y. Genome-wide analysis of the Brachypodium distachyon (L.) P. Beauv. Hsp90 gene family reveals molecular evolution and expression profiling under drought and salt stresses. PLoS One. 2017;12(12):e0189187. https://doi.org/10.1371/journal.pone.0189187.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  54. 54.

    Brunetti SC, Arseneault MKM, Gulick PJ. Characterization of the Esi3/RCI2/PMP3 gene family in the Triticeae. BMC Genomics. 2018;19(1):898. https://doi.org/10.1186/s12864-018-5311-8.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  55. 55.

    Liang Y, Jiang Y, Du M, Li B, Chen L, Chen M, et al. ZmASR3 from the maize ASR gene family positively regulates drought tolerance in transgenic Arabidopsis. Int J Mol Sci. 2019;20(9):2278. https://doi.org/10.3390/ijms20092278.

    CAS  Article  PubMed Central  Google Scholar 

  56. 56.

    Hayama R, Sarid-Krebs L, Richter R, Fernández V, Jang S, Coupland G. PSEUDO RESPONSE REGULATORs stabilize CONSTANS protein to promote flowering in response to day length. EMBO J. 2017;36(7):904–18. https://doi.org/10.15252/embj.201693907.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  57. 57.

    Li Y, Xu M. CCT family genes in cereal crops: A current overview. Crop J. 2017;5(6):449–58. https://doi.org/10.1016/j.cj.2017.07.001.

    Article  Google Scholar 

  58. 58.

    Voß U, Wilson MH, Kenobi K, Gould PD, Robertson FC, Peer WA, et al. The circadian clock rephases during lateral root organ initiation in Arabidopsis thaliana. Nat Commun. 2015;6(1):7641. https://doi.org/10.1038/ncomms8641.

    Article  PubMed  Google Scholar 

  59. 59.

    Wahl V, Brand LH, Guo Y-L, Schmid M. The FANTASTIC FOUR proteins influence shoot meristem size in Arabidopsis thaliana. BMC Plant Biol. 2010;10(1):285. https://doi.org/10.1186/1471-2229-10-285.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  60. 60.

    Zhan X, Shen Q, Wang X, Hong Y. The Sulfoquinovosyltransferase-like enzyme SQD2.2 is involved in flavonoid glycosylation, regulating sugar metabolism and seed setting in rice. Sci Rep. 2017;7(1):4685.

    Article  Google Scholar 

  61. 61.

    Pedrotti L, Mueller MJ, Waller F. Piriformospora indica root colonization triggers local and systemic root responses and inhibits secondary colonization of distal roots. PLoS One. 2013;8(7):e69352. https://doi.org/10.1371/journal.pone.0069352.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  62. 62.

    Harrach B, Baltruschat H, Barna B, Fodor J, Kogel KH. The mutualistic fungus Piriformospora indica protects barley roots from a loss of antioxidant capacity caused by the necrotrophic pathogen Fusarium culmorum. Mol Plant Microbe Interact. 2013;26(5):599–605. https://doi.org/10.1094/MPMI-09-12-0216-R.

    CAS  Article  PubMed  Google Scholar 

  63. 63.

    Ye W, Shen C-H, Lin Y, Chen P-J, Xu X, Oelmüller R, et al. Growth promotion-related miRNAs in Oncidium orchid roots colonized by the endophytic fungus Piriformospora indica. PLoS One. 2014;9(1):e84920. https://doi.org/10.1371/journal.pone.0084920.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  64. 64.

    Jeong D-H, Schmidt SA, Rymarquis LA, Park S, Ganssmann M, German MA, et al. Parallel analysis of RNA ends enhances global investigation of microRNAs and target RNAs of Brachypodium distachyon. Genome Biol. 2013;14(12):R145. https://doi.org/10.1186/gb-2013-14-12-r145.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  65. 65.

    Šečić E, Kogel KH, Ladera-Carmona MJ. Biotic stress-associated microRNA families in plants. J Plant Physiol. 2021;263:153451. https://doi.org/10.1016/j.jplph.2021.153451.

    CAS  Article  PubMed  Google Scholar 

  66. 66.

    Morea E, Silva E, Silva G, Valente G, Barrera C, Vincentz M, et al. Functional and evolutionary analyses of the miR156 and miR529 families in land plants. BMC Plant Biol. 2016;16:40.

    Article  Google Scholar 

  67. 67.

    Hackenberg M, Shi B-J, Gustafson P, Langridge P. Characterization of phosphorus-regulated miR399 and miR827 and their isomirs in barley under phosphorus-sufficient and phosphorus-deficient conditions. BMC Plant Biol. 2013;13(1):214. https://doi.org/10.1186/1471-2229-13-214.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  68. 68.

    Baev V, Milev I, Naydenov M, Apostolova E, Minkov G, Minkov I, et al. Implementation of a de novo genome-wide computational approach for updating Brachypodium miRNAs. Genomics. 2011;97(5):282–93. https://doi.org/10.1016/j.ygeno.2011.02.008.

    CAS  Article  PubMed  Google Scholar 

  69. 69.

    Wang B, Sun Y, Song N, Zhao M, Liu R, Feng H, et al. Puccinia striiformis f. sp. tritici microRNA-like RNA 1 (Pst-milR1), an important pathogenicity factor of Pst, impairs wheat resistance to Pst by suppressing the wheat pathogenesis-related 2 gene. New Phytol. 2017;215(1):338–50. https://doi.org/10.1111/nph.14577.

    CAS  Article  PubMed  Google Scholar 

  70. 70.

    Lammers PJ, Jun J, Abubaker J, Arreola R, Gopalan A, Bago B, et al. The glyoxylate cycle in an arbuscular mycorrhizal fungus. Carbon flux and gene expression. Plant Physiol. 2001;127(3):1287–98. https://doi.org/10.1104/pp.010375.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  71. 71.

    Hynes MJ, Murray SL, Andrianopoulos A, Davis MA. Role of carnitine acetyltransferases in acetyl coenzyme A metabolism in Aspergillus nidulans. Eukaryot Cell. 2011;10(4):547–55. https://doi.org/10.1128/EC.00295-10.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  72. 72.

    Li J, Yu L, Yang J, Dong L, Tian B, Yu Z, et al. New insights into the evolution of subtilisin-like serine protease genes in Pezizomycotina. BMC Evol Biol. 2010;10(1):68. https://doi.org/10.1186/1471-2148-10-68.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  73. 73.

    Corrales Escobosa AR, Rangel Porras RA, Meza Carmen V, Gonzalez Hernandez GA, Torres Guzman JC, Wrobel K, et al. Fusarium oxysporum Adh1 has dual fermentative and oxidative functions and is involved in fungal virulence in tomato plants. Fungal Genet Biol. 2011;48(9):886–95. https://doi.org/10.1016/j.fgb.2011.06.004.

    CAS  Article  PubMed  Google Scholar 

  74. 74.

    Jiang J, Yun Y, Yang Q, Shim W-B, Wang Z, Ma Z. A type 2C protein phosphatase FgPtc3 is involved in cell wall integrity, lipid metabolism, and virulence in Fusarium graminearum. PLoS One. 2011;6(9):e25311. https://doi.org/10.1371/journal.pone.0025311.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  75. 75.

    Xiang Q, Rasmussen C, Glass NL. The ham-2 locus, encoding a putative transmembrane protein, is required for hyphal fusion in Neurospora crassa. Genetics. 2002;160(1):169–80. https://doi.org/10.1093/genetics/160.1.169.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  76. 76.

    Lannoo N, Van Damme EJM. Lectin domains at the frontiers of plant defense. Front Plant Sci. 2014;5:397.

    PubMed  PubMed Central  Google Scholar 

  77. 77.

    Torii KU. Leucine-rich repeat receptor kinases in plants: structure, function, and signal transduction pathways. Int Rev Cytol. 2004;234:1–46. https://doi.org/10.1016/S0074-7696(04)34001-5.

    CAS  Article  PubMed  Google Scholar 

  78. 78.

    Mugford ST, Osbourn A. Evolution of serine carboxypeptidase-like acyltransferases in the monocots. Plant Signal Behav. 2010;5(2):193–5. https://doi.org/10.4161/psb.5.2.11093.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  79. 79.

    Corrêa LGG, Riaño-Pachón DM, Schrago CG, dos Santos RV, Mueller-Roeber B, Vincentz M. The role of bZIP transcription factors in green plant evolution: adaptive features emerging from four founder genes. PLoS One. 2008;3(8):e2944. https://doi.org/10.1371/journal.pone.0002944.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  80. 80.

    Hirsch S, Oldroyd GED. GRAS-domain transcription factors that regulate plant development. Plant Signal Behav. 2009;4(8):698–700. https://doi.org/10.4161/psb.4.8.9176.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  81. 81.

    Ambawat S, Sharma P, Yadav NR, Yadav RC. MYB transcription factor genes as regulators for plant responses: an overview. Physiol Mol Biol Plants. 2013;19(3):307–21. https://doi.org/10.1007/s12298-013-0179-1.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  82. 82.

    Zanini S, Šečić E, Busche T, Galli M, Zheng Y, Kalinowski J, et al. Comparative analysis of transcriptome and sRNAs expression patterns in the Brachypodium distachyon- Magnaporthe oryzae pathosystems. Int J Mol Sci. 2021;22(2). https://doi.org/10.3390/ijms22020650.

  83. 83.

    Šečić E, Kogel KH. Requirements for fungal uptake of dsRNA and gene silencing in RNAi-based crop protection strategies. Curr Opin Biotech. 2021;70:136–42. https://doi.org/10.1016/j.copbio.2021.04.001.

    CAS  Article  PubMed  Google Scholar 

  84. 84.

    Murashige T, Skoog F. A revised medium for rapid growth and bio assays with tobacco tissue cultures. Physiol Plant. 1962;15(3):473–97. https://doi.org/10.1111/j.1399-3054.1962.tb08052.x.

    CAS  Article  Google Scholar 

  85. 85.

    Pontecorvo G, Roper JA, Hemmons LM, Macdonald KD, Bufton AW. The genetics of Aspergillus nidulans. Ad Genet. 1953;5:141–238. https://doi.org/10.1016/S0065-2660(08)60408-3.

    CAS  Article  Google Scholar 

  86. 86.

    Resequencing and genome assembly of Serendipita indica (syn. Piriformospora indica). European Nucleotide Archive. 2021. https://www.ebi.ac.uk/ena/browser/view/PRJEB45884

  87. 87.

    Wibberg D, Stadler M, Lambert C, Bunk B, Spröer C, Rückert C, et al. High quality genome sequences of thirteen Hypoxylaceae (Ascomycota) strengthen the phylogenetic family backbone and enable the discovery of new taxa. Fungal Divers. 2021;106:7–28.

  88. 88.

    Wibberg D, Andersson L, Tzelepis G, Rupp O, Blom J, Jelonek L, et al. Genome analysis of the sugar beet pathogen Rhizoctonia solani AG2-2IIIB revealed high numbers in secreted proteins and cell wall degrading enzymes. BMC Genomics. 2016;17(1):245. https://doi.org/10.1186/s12864-016-2561-1.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  89. 89.

    Koren S, Walenz BP, Berlin K, Miller JR, Bergman NH, Phillippy AM. Canu: scalable and accurate long-read assembly via adaptive k-mer weighting and repeat separation. Genome Res. 2017;27(5):722–36. https://doi.org/10.1101/gr.215087.116.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  90. 90.

    Walker BJ, Abeel T, Shea T, Priest M, Abouelliel A, Sakthikumar S, et al. Pilon: an integrated tool for comprehensive microbial variant detection and genome assembly improvement. PLoS One. 2014;9(11):e112963. https://doi.org/10.1371/journal.pone.0112963.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  91. 91.

    Li H. Aligning sequence reads, clone sequences and assembly contigs with BWA-MEM. ArXiv. 2013;1303.3997:1–3.

  92. 92.

    Langmead B, Salzberg SL. Fast gapped-read alignment with Bowtie 2. Nat Methods. 2012;9(4):357–9. https://doi.org/10.1038/nmeth.1923.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  93. 93.

    Galagan JE, Calvo SE, Borkovich KA, Selker EU, Read ND, Jaffe D, et al. The genome sequence of the filamentous fungus Neurospora crassa. Nature. 2003;422(6934):859–68 NC12 genome assembly. 2014. GCA accession number: GCA_000182925.2 (https://www.ebi.ac.uk/ena/browser/view/GCA_000182925.2).

    CAS  Article  Google Scholar 

  94. 94.

    Petersen TN, Brunak S, von Heijne G, Nielsen H. SignalP 4.0: discriminating signal peptides from transmembrane regions. Nat Methods. 2011;8(10):785–6. https://doi.org/10.1038/nmeth.1701.

    CAS  Article  PubMed  Google Scholar 

  95. 95.

    Krogh A, Larsson B, von Heijne G, Sonnhammer ELL. Predicting transmembrane protein topology with a hidden markov model: application to complete genomes. J Mol Biol. 2001;305(3):567–80. https://doi.org/10.1006/jmbi.2000.4315.

    CAS  Article  PubMed  Google Scholar 

  96. 96.

    Emanuelsson O, Nielsen H, Brunak S, von Heijne G. Predicting subcellular localization of proteins based on their N-terminal amino acid sequence. J Mol Biol. 2000;300(4):1005–16. https://doi.org/10.1006/jmbi.2000.3903.

    CAS  Article  PubMed  Google Scholar 

  97. 97.

    Ray P, Chi M-H, Guo Y, Chen C, Adam C, Kuo A, et al. Genome sequence of the plant growth promoting fungus Serendipita vermifera subsp. bescii: The First Native Strain from North America. Phytobiomes J. 2018;2(2):62–3 Serendipita vermifera ssp. bescii NFPB0129 v1.0. (https://mycocosm.jgi.doe.gov/Sebvebe1/Sebvebe1.home.html).

    Article  Google Scholar 

  98. 98.

    Martin F, Aerts A, Ahrén D, Brun A, Danchin EGJ, Duchaussoy F, et al. The genome of Laccaria bicolor provides insights into mycorrhizal symbiosis. Nature. 2008;452(7183):88–92 Laccaria bicolor v2.0. (https://mycocosm.jgi.doe.gov/Lacbi2/Lacbi2.home.html).

    CAS  Article  Google Scholar 

  99. 99.

    Blom J, Kreis J, Spänig S, Juhre T, Bertelli C, Ernst C, et al. EDGAR 2.0: an enhanced software platform for comparative gene content analyses. Nucleic Acids Res. 2016;44(W1):W22–8. https://doi.org/10.1093/nar/gkw255.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  100. 100.

    mRNA-seq of Brachypodium distachyon roots inoculated with Serendipita indica (syn. Piriformospora indica). ArrayExpress database at EMBL-EBI. 2021. https://www.ebi.ac.uk/arrayexpress/experiments/E-MTAB-10649

  101. 101.

    Andrews S, Krueger F, Segonds-Pichon A, Biggins L, Krueger C, Wingett S. FastQC: a quality control tool for high throughput sequence data. Babraham; 2016. http://www.bioinformatics.babraham.ac.uk/projects/fastqc [Accessed 1 Feb 2019]

  102. 102.

    Brachypodium distachyon Bd21-3 v1.1 DOE-JGI, http://phytozome.jgi.doe.gov/

  103. 103.

    Kim D, Paggi JM, Park C, Bennett C, Salzberg SL. Graph-based genome alignment and genotyping with HISAT2 and HISAT-genotype. Nat Biotechnol. 2019;37(8):907–15. https://doi.org/10.1038/s41587-019-0201-4.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  104. 104.

    Kupfer DM, Drabenstot SD, Buchanan KL, Lai H, Zhu H, Dyer DW, et al. Introns and splicing elements of five diverse fungi. Eukaryot Cell. 2004;3(5):1088–100. https://doi.org/10.1128/EC.3.5.1088-1100.2004.

    Article  PubMed  PubMed Central  Google Scholar 

  105. 105.

    Walters B, Lum G, Sablok G, Min XJ. Genome-wide landscape of alternative splicing events in Brachypodium distachyon. DNA Res. 2013;20(2):163–71. https://doi.org/10.1093/dnares/dss041.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  106. 106.

    Anders S, Pyl PT, Huber W. HTSeq--a Python framework to work with high-throughput sequencing data. Bioinformatics. 2015;31(2):166–9. https://doi.org/10.1093/bioinformatics/btu638.

    CAS  Article  PubMed  Google Scholar 

  107. 107.

    Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15(12):550. https://doi.org/10.1186/s13059-014-0550-8.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  108. 108.

    Du Z, Zhou X, Ling Y, Zhang Z, Su Z. agriGO: a GO analysis toolkit for the agricultural community. Nucleic Acids Res. 2010 Jul;38(suppl_2):W64–70. https://doi.org/10.1093/nar/gkq310.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  109. 109.

    Sievert C. Interactive web-based data visualization with R, plotly, and shiny: Chapman and Hall/CRC; 2020. Available from: https://plotly-r.com. https://doi.org/10.1201/9780429447273.

    Book  Google Scholar 

  110. 110.

    Wickham H. ggplot2: Elegant graphics for data analysis. New York: Springer-Verlag; 2016. Available from: http://ggplot2.org

    Book  Google Scholar 

  111. 111.

    R Core Team. R: A language and environment for statistical computing. Vienna; 2017. Available from: https://www.r-project.org/

  112. 112.

    Götz S, García-Gómez JM, Terol J, Williams TD, Nagaraj SH, Nueda MJ, et al. High-throughput functional annotation and data mining with the Blast2GO suite. Nucleic Acids Res. 2008;36(10):3420–35. https://doi.org/10.1093/nar/gkn176.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  113. 113.

    sRNA-seq of Brachypodium distachyon roots inoculated with Serendipita indica (syn. Piriformospora indica). ArrayExpress database at EMBL-EBI. 2021. https://www.ebi.ac.uk/arrayexpress/experiments/E-MTAB-10650

  114. 114.

    Martin M. Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet.journal. 2011;17(1):10–12. Next Gener Seq Data Anal - 1014806/ej171200.

  115. 115.

    Langmead B, Trapnell C, Pop M, Salzberg SL. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol. 2009;10(3):R25. https://doi.org/10.1186/gb-2009-10-3-r25.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  116. 116.

    The RNAcentral Consortium, Petrov AI, Kay SJE, Kalvari I, Howe KL, Gray KA, et al. RNAcentral: a comprehensive database of non-coding RNA sequences. Nucleic Acids Res. 2017;45(D1):D128–34. https://doi.org/10.1093/nar/gkw1008.

    CAS  Article  Google Scholar 

  117. 117.

    Johnson NR, Yeoh JM, Coruh C, Axtell MJ. Improved placement of multi-mapping small RNAs. G3 (Bethesda). 2016;6(7):2103–11.

    CAS  Article  Google Scholar 

  118. 118.

    Zanini S, Šečić E, Jelonek L, Kogel KH. A bioinformatics pipeline for the analysis and target prediction of RNA effectors in bidirectional communication during plant–microbe interactions. Front Plant Sci. 2018;9. https://doi.org/10.3389/fpls.2018.01212.

  119. 119.

    Dai X, Zhao PX. psRNATarget: a plant small RNA target analysis server. Nucleic Acids Res. 2011;39(suppl_2):W155–9. https://doi.org/10.1093/nar/gkr319.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  120. 120.

    Chen H, Boutros PC. VennDiagram: a package for the generation of highly-customizable Venn and Euler diagrams in R. BMC Bioinformatics. 2011;12(1):35. https://doi.org/10.1186/1471-2105-12-35.

    Article  PubMed  PubMed Central  Google Scholar 

  121. 121.

    Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2−ΔΔCT Method. Methods. 2001;25(4):402–8. https://doi.org/10.1006/meth.2001.1262.

    CAS  Article  Google Scholar 

  122. 122.

    Chen C, Ridzon DA, Broomer AJ, Zhou Z, Lee DH, Nguyen JT, et al. Real-time quantification of microRNAs by stem–loop RT–PCR. Nucleic Acids Res. 2005;33(20):e179. https://doi.org/10.1093/nar/gni178.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  123. 123.

    Kramer MF. Stem-Loop RT-qPCR for miRNAs. Curr Protoc Mol Biol. 2011;95:15.10.1.

    Article  Google Scholar 

Download references

Acknowledgements

We thank Elke Stein, Dagmar Biedenkopf, Ute Micknass and Christina Birkenstock for technical assistance. We thank Dr. John Vogel and the DOE-JGI for permission to use the Bd21-3 genome under early access conditions. Brachypodium distachyon Bd21-3 is a gift of R. Sibout, INRA Versailles. We are very thankful to D’Maris Dempsey for her assistance in preparation of the manuscript.

Funding

This work was supported by the Deutsche Forschungsgemeinschaft (DFG), Research Training Group (RTG) 2355, and Research Unit FOR5116 to KHK. We gratefully acknowledge the bioinformatics support of the BMBF-funded project Bielefeld-Gießen Center for Microbial Bioinformatics-BiGi (grant number 031A533) within the German Network for Bioinformatics Infrastructure (de.NBI). Open Access funding enabled and organized by Projekt DEAL.

Author information

Affiliations

Authors

Contributions

EŠ, SZ, DW, and KHK wrote the text and drafted the figures. EŠ and SZ designed and conducted the establishment experiments, EŠ, SZ, TB, DW, LJ, and JK conducted the RNAseq and genome resequencing experiments and analyzed the data. SN and JS conducted the stem-loop PCR experiments and sequencing analysis. JI and JT performed microscopy in sterile conditions. All authors read and approved the final manuscript.

Corresponding author

Correspondence to Karl-Heinz Kogel.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1.

Supporting individual data values displayed in Figure 1a-c and Additional File 2: Figure S5.

Additional file 2: Fig. S1.

Serendipita indica (Si) colonization has an effect on Brachypodium distachyon (Bd) root structure. Fig. S2. Progress of Serendipita indica (Si) spore proliferation during colonization of Brachypodium distachyon Bd21-3. Fig. S3. Comparison of the resequenced Serendipita indica (Si) genome to the 2011 assembly. Fig.S4. Comparison of the resequenced Serendipita indica (Si) genome to Serendipita vermifera and Laccaria bicolor. Fig. S5. qRT-PCR confirmation of DEGs identified during mRNA sequencing. Fig. S6. Filtering pipelines applied in the analysis. Fig. S7. Size distribution of total and unique putative ck (cross-kingdom) -sRNAs. Fig. S8. Percentage distribution of the 5′ terminal nucleotide in unique putative endogenous sRNAs. Fig. S9. Percentage distribution of the 5′ terminal nucleotide in unique putative ck-sRNAs. Fig. S10. Percentage distribution of the 5′ terminal nucleotide in total putative endogenous sRNAs. Fig. S11. Percentage distribution of the 5′ terminal nucleotide in total putative ck-sRNAs. Fig. S12. Stem-loop PCR (gel electrophoresis) of some Si and Bd21-3 sRNAs expressed in the Bd-Si sample. Table S1. Quantification of identified features in the resequenced genome of Serendipita indica (Si). Table S2. Total reads from the Bd-C (mock-treated), Bd-Si (colonized root) and Si-ax (axenic culture) samples and their alignment rate (HISAT2) to the corresponding genomes. Table S3. Significant gene ontology terms of molecular function in the differentially expressed genes (DEGs) datasets. Table S4. Predicted protein effectors identified in the resequenced Serendipita indica (Si) genome. Table S5. Candidate RNAi machinery proteins predicted from the resequenced Serendipita indica (Si) genome. Table S6. Total and unique reads for filtered sRNAs from Bd-C (mock-treated), Bd-Si (colonized root) and Si-ax (axenic culture). Table S7. Sequences of putative sRNAs in Tables 4 and 5. Table S8. sRNA sequences after stem-loop PCR amplification. Table S9. Primers used for qRT-PCR (Figure S5). Table S10. Primers used for stem-loop PCR and sequencing (Figure S12, Table S8).

Additional file 3.

Uncropped gel picture annotated in Additional file 2: Figure S12a.

Additional file 4.

Uncropped gel picture annotated in Additional file 2: Figure S12b.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Šečić, E., Zanini, S., Wibberg, D. et al. A novel plant-fungal association reveals fundamental sRNA and gene expression reprogramming at the onset of symbiosis. BMC Biol 19, 171 (2021). https://doi.org/10.1186/s12915-021-01104-2

Download citation

Keywords

  • Brachypodium distachyon
  • Genome sequencing
  • Sebacinalean symbiosis
  • Serendipita indica
  • Small RNAs