Evaluation of rice cultivars against bakanae disease
Out of 12 genotypes screened against bakanae disease, Selenio resulted the most resistant (disease index 17.0 %), while Dorella turned out to be the most susceptible (disease index 82.5 %). A sporadic presence of yellowish plants, without typical symptoms of bakanae disease was observed in cv. Selenio. On the contrary, typical bakanae-diseased symptoms with elongated and thin internodes which progressively led to plant death were observed in plants of cv. Dorella (69 %). The other tested genotypes showed intermediate levels of resistance (Fig. 1).
The F. fujikuroi colonization was quantified in the Selenio and Dorella genotypes at 1 and 3 weeks post germination (wpg) (Fig. 2). Dorella showed abundant infection with F. fujikuroi at 3 wpg (15 times more than Selenio), confirming its sensitivity, while at 1 wpg a low F. fujikuroi infection level with no significant differences between the two cultivars was shown. Phenotypic and quantitative molecular differences in Selenio and Dorella inoculated with F. fujikuroi allowed to select the two cultivars for further transcriptomic studies.
Differentially expressed genes in bakanae-infected rice
RNA was isolated from leaves of both cvs. Selenio and Dorella, inoculated with F. fujikuroi or with the mock, at 1 and 3 wpg, and subjected to whole transcriptome sequencing via RNA-seq. Selenio was constituted in 1987, while Dorella was constituted in 1998, both in Italy. Both cultivars showed similar growth rate and developmental stage during the sampling period. The growing cycle from sowing to spike lasts 95 days for both cultivars and the growing cycle from sowing to maturity lasts, respectively, 145 days for Selenio (plant height: 76.0 cm; ear length: 17.0 cm) and 140 days for Dorella (plant height: 96.0 cm; ear length: 21.0 cm). Both cultivars are not pigmented and the grains are not aromatic and not glutinous. Three biological replicates were sequenced for each genotype (Selenio vs Dorella), disease condition (mock vs bakanae-inoculated), and infection time point (1 vs 3 wpg). Raw reads generated from the Illumina Genome Analyzer GAIIx were filtered by Illumina passed-filter call. Subsequently, adapters identified by fastQC and low-quality regions were filtered out by cutadapt application [17]. Filtered, 51-base reads for each biological replicate (on average, 17 million reads, Additional file 1: Table S1) were mapped with Bowtie/2TopHat2 to the rice genome sequence (Oryza sativa Japonica - Nipponbare IRGSP-1.0.20).
Read counts were generated from Bam alignment files with HTSeq software [18]. Data normalization and call of differentially expressed genes (DEGs) was implemented with DESeq2 R package [19, 20].
Pearson correlation coefficients for normalized expression values of samples are shown in Additional file 2: Table S2. All biological replicates sharing cultivar type (Selenio vs Dorella), infection time point (1 vs 3 wpg) and treatment (mock vs bakanae-inoculated) showed correlation coefficients above 0.9 indicating a good reproducibility between biological replicates. Expression values of all detected transcribed genes and non-coding RNAs are reported as DESeq-normalized read counts and log2 fold changes in Additional file 3: Table S3.
The number of DEGs at 1 wpg was 80 in Selenio (Selenio mock 1 wpg vs Selenio inoculated 1 wpg) of which 24 were Selenio-specific expressed loci, 5 were common DEGs, 7 were genes called as DEGs in Selenio but expressed loci also in Dorella, and 44 were DEGs common in both cultivars at 3 wpg; and 1,285 in Dorella (Dorella mock 1 wpg vs Dorella inoculated 1 wpg) of which 797 were Dorella-specific expressed loci, 7 were expressed loci in Dorella but called as DEGs also in Selenio, 5 were common DEGs, and 476 were DEGs common in Selenio and Dorella at 3 wpg.
The number of DEGs was higher at 3 wpg in both cultivars: 3,119 in Selenio (Selenio mock 3 wpg vs Selenio inoculated 3 wpg), of which 1,560 were Selenio-specific expressed loci, 19 were common DEGs, 1,346 were genes called as DEGs in Selenio but expressed loci also in Dorella, and 194 were DEGs common in Selenio and Dorella at 1 wpg; and 5,095 in Dorella (Dorella mock 3 wpg vs Dorella inoculated 3 wpg), of which 3,347 were Dorella-specific expressed loci, 1,346 were expressed loci in Dorella but called as DEGs also in Selenio, 19 were common DEGs, and 383 were DEGs common in both cultivars at 1 wpg (Fig. 3).
Modulated genes for both genotypes were therefore more abundant at 3 wpg than at 1 wpg, indicating a higher infection level in mature rice plants compared to younger plants [21]. Rice plants germinated one week after sowing, and the first bakanae symptoms occurred at 1 wpg. Susceptible cultivars of rice, inoculated with pathogenic F. fujikuroi, started to die at 3 wpg, indicating that this is another crucial time point [22]. The higher number of DEGs identified in the susceptible cv. Dorella compared to the resistant cv. Selenio is consistent with a highly visible infection of F. fujikuroi in the tissues of Dorella, resulting in a numerically higher number of transcriptionally modulated genes. This finding is similar to what has been described for the interaction between rice and other rice pathogens, where the susceptible cultivar showed much more DEGs than the resistant one [14, 16]. However, the defence mechanisms involved in the bakanae disease responses were largely divergent between Selenio and Dorella since the two genotypes shared only 5 and 19 DEGs at 1 and 3 wpg, respectively.
The proportion of up-regulated DEGs was higher in the susceptible genotype (62.0 %) with respect to the resistant one (38.8 %) at 1 wpg, whereas it was higher in the resistant genotype (75.9 %) compared to the susceptible one (41.2 %) at 3 wpg, as indicated in the MA-plots (Fig. 4). Indeed, a higher proportion of down-regulated genes was observed in Dorella at 3 wpg (58.8 %) when systemic infection, characterized by specific symptoms of bakanae in all the green tissues of rice, was well established.
Functional annotations of the rice transcripts available from RAP-DB (http://rapdb.dna.affrc.go.jp/) were used to investigate the functions and the role of the selected genes differentially expressed in response to F. fujikuroi infection. The criteria used for the description of genes in the following sections were based either on the higher fold changes (FC) in cv. Selenio in the case of cytochrome P450 monooxygenases, or on the defence role of genes already reported in rice infected with other pathogens [14, 15, 23–26].
Cytochrome P450 monooxygenases
Cytochrome P450 monooxygenases (P450s) are universal enzymes that catalyse the oxidation of many substrates by the activation of molecular oxygen. Plant P450s catalyse several enzymatic steps for various plant metabolites such as defence-related metabolites, phytohormones, pigments, fatty acids, and lignin [27, 28]. Here, we found a putative cytochrome P450 (Os12g0582700; family 94) as the most up-regulated gene in Selenio (3 wpg) with a log2FC 6.60 (Additional file 4: Table S4). Other highly up-regulated P450s in Selenio (3 wpg), also located on chromosome 12, were Os12g0443000 and Os12g0150200 (both genes belonging to family 94). The three P450s were not DEGs either in Selenio or in Dorella at 1 wpg. Interestingly, P450 Os12g0443000 was the highest down-regulated gene in Dorella with log2FC -3.66 at 3 wpg, while the other two P450s were not present among the DEGs in this cultivar.
PR genes
Pathogenesis-related (PR) proteins have antifungal activity against many phytopathogenic fungi, such as Phytophthora infestans, Phytophthora parasitica, Uromyces fabae, and Erysiphe graminis [29, 30]. PR1 plant proteins have homology to the superfamily of cysteine-rich proteins, and when act as protease inhibitors they can suppress the enzymatic activities in response to proteinase attack produced by plant pathogens [31, 32]. PR1 genes were induced in wheat infected by Fusarium pseudograminearum and in rice infected by Magnaporthe oryzae [15, 23, 33].
We found, respectively five and six up-regulated PR1 genes in Selenio and in Dorella at 3 wpg (Additional file 5: Table S5). Furthermore, there were no down-regulated proteinase inhibitors in Selenio, while in Dorella five of them were down-regulated. Interestingly, out of five up-regulated genes in Selenio (all proteinase inhibitors), three of them were down-regulated in Dorella (Os01g0615050, Os03g0429000 and Os01g0127600). At 1 wpg, there were no induced PR1 genes in Selenio, whereas only one PR1 gene was transcriptionally activated in Dorella.
Few other PR genes, including thaumatin-like genes, were identified as modulated by F. fujikuroi infection mainly at 3 wpg in Selenio and Dorella (4 and 8 genes, respectively). Out of four up-regulated PR genes in Selenio, three were down-regulated in Dorella (Os01g0693400, Os12g0555300, and Os12g0569500; Additional file 5: Table S5). Furthermore, the four up-regulated genes in Selenio were also induced in rice during M. oryzae infection (Os01g0693400 and Os12g0569500 in compatible interactions, and Os04g0398000 and Os12g0555300 in incompatible interactions) [15].
Glycoside hydrolases
Glycoside hydrolases (GH) catalyse the hydrolysis of glycosidic bonds in cell wall polymers. The rice GH are classified into 34 families encoding 437 GH genes [26]. The largest GH17 family contains PR2 genes producing β-1,3-glucanases that hydrolyse β-1,3-D-glucosidic linkages in β-1,3-glucans [34]. Fungal cell walls usually contain both chitins and β-1,3-glucans, and the contemporary expression of chitinases and β-1,3-glucanases indicate their mutual function in fungal cell wall degradation [35]. Sharma et al. [26] found 128 up-regulated and 82 down-regulated GH genes during response of rice to four rice pathogens, whose 19 up-regulated in response to M. oryzae, the causal agent of blast disease. Kawahara et al. [15] reported up-regulation of 23 GH17 genes in early-stage of M. oryzae infection.
In the present work, an higher up-regulation of GH genes, comprising β-1,3-glucanases, was observed at 3 wpg; 16 genes in Selenio and 13 genes in Dorella (Additional file 6: Table S6). Five up-regulated GH genes (Os01g0713200, Os01g0860800, Os01g0940800, Os03g0792800, and Os11g0704600) in Selenio were in common with the rice-M. oryzae interaction [15].
At 1 wpg, Selenio did not induce any GH gene, while Dorella up-regulated 9 glucanase genes. Dorella shared 4 up-regulated GH genes during both infection stages. A higher activation of glucanase genes in the susceptible genotype compared to the resistant one at 1 wpg could be attributed to the F. fujikuroi colonization of Dorella plants and to the consequent activation of the mechanisms of fungal cell-wall degradation.
Chitinases
Chitinases belong to the PR protein families PR3, 4, 8, and 11, and are involved in the plant defence response against pathogens and pests by hydrolysing the chitin in the cell wall of fungi and in the skeleton of insects [36, 37]. Three up-regulated chitinases were common in Dorella at 1 and 3 wpg, while Selenio did not show activation of chitinases at 1 wpg (Additional file 7: Table S7). Interestingly, among the up-regulated chitinases in Selenio at 3 wpg, Os02g0605900 was down-regulated in Dorella during both infection time points. Os02g0605900 was also up-regulated during incompatible interaction between rice and M. oryzae [15]. At 3 wpg, we found 10 up-regulated chitinase genes in the susceptible genotype, whereas Selenio showed only 2 up-regulated chitinases. Higher up-regulation of chitinase genes in Dorella compared to Selenio at 3 wpg might be a consequence of fungal colonization of plant cells and subsequent activation of fungal cell-wall degradation, as in the case of glucanase genes at 1 wpg.
Peroxidases
Peroxidases are PR9 proteins induced in plant tissues upon pathogen infection and they are known as defence-related proteins. Peroxidases are produced to prevent cellular diffusion of pathogens by massive production of reactive oxygen and nitrogen species, and subsequent creation of an highly toxic environment, or by development/reinforcement of structural barriers [38–40]. In our study, we found 12 up-regulated peroxidase genes in Dorella and 7 in Selenio at 3 wpg (Additional file 8: Table S8). Activation of a higher number of peroxidases was also found at 1 wpg in susceptible genotype (9) as compared to the resistant genotype (2). Selenio did not share activated peroxidases between 1 and 3 wpg, while one induced peroxidase was common in Dorella during both infection time points (Os02g0236600). Lower induction of peroxidases in Selenio could be attributable to a reduced spreading of F. fujikuroi in this genotype. Three peroxidases (Os01g0787000, Os03g0285700, and Os07g0677600) up-regulated in Selenio were also induced in rice infected with M. oryzae [15].
Germin-like proteins
Oxalate oxidase-like genes, recently known as germin-like protein (GLP) genes, are included in plant defence mechanisms, where some of them increase their expression upon pathogen infection or insect attack [41–44]. GLP genes are known to contribute to rice resistance. Manosalva et al. [24] demonstrated that when GLP genes of chromosome 8 were suppressed, rice plants were more susceptible to rice blast and sheath blight (caused by Rhizoctonia solani).
At 3 wpg, 4 up-regulated GLP genes were identified in Selenio, while 5 up-regulated and 2 down-regulated GLP genes were found in Dorella (Additional file 9: Table S9). Interestingly, 2 up-regulated GLPs (Os08g0189600 and Os11g0537350) in Selenio were down-regulated in Dorella. Out of four up-regulated GLPs in Selenio, two genes (Os08g0189600 and Os08g0189400) are located on chromosome 8, where other GLP genes involved in rice fungal resistance are present [24]. On the other hand, Selenio did not activate GLP genes at 1 wpg, whereas Dorella had all GLPs genes up-regulated at this stage.
MAP kinases
Mitogen-activated protein kinases (MAPKs) play an important role in many resistance (R)-mediated defence responses to plant pathogens [45]. MAPK cascades are involved in a range of signalling pathways downstream of receptor kinases, such as the biosynthesis of phytoalexins in plants [46]. MAPK cascades are composed of three functionally linked protein kinases: MAPK, phosphorylated and activated by a MAPK kinase (MAPKK), which itself is activated by another protein kinase, a MAPKK kinase (MAPKKK) [47].
A MAPK cascade, including the MAPKK OsMKK4 and two MAPKs, OsMPK3 and OsMPK6, is involved in chitin elicitor-induced biosynthesis of diterpenoid phytoalexins in rice [25]. Two of these genes were up-regulated in Selenio at 3 wpg (OsMKK4 or Os02g0787300 with a log2FC 3.27, and OsMPK3 or Os03g0285800 with a log2FC 1.64), whereas one of them was the most down-regulated MAP-kinase in Dorella (OsMKK4 with a log2FC -2.83) (Additional file 10: Table S10). Furthermore, Selenio had five further MAP kinases activated while the majority of MAP-kinases were down-regulated in Dorella (8 genes in total). No MAP kinases were DEGs either in Selenio or in Dorella at 1 wpg.
WRKY transcriptional factors
WRKY transcription factors are a superfamily of zinc-finger transcription factors which contain a highly conserved DNA-binding WRKY domain. Most WRKY proteins bind to the conserved W-box, which is present in the promoters of many defence-related genes [48, 49].
At 3 wpg, Selenio showed an induction of 22 WRKY genes with FC ranging from 2 to 15 and a down-regulation for only two WRKY genes (Additional file 11: Table S11). Conversely, at the same time point Dorella up-regulated only 4 WRKY genes (FC from 2 to 4), and down-regulated 15 WRKY genes (Additional file 11: Table S11). Selenio showed no activation of WRKY genes at 1 wpg, while Dorella showed 4 WRKY up-regulated genes with only one of them shared between 1 and 3 wpg. Interestingly, 8 up-regulated WRKY genes in Selenio at 3 wpg were down-regulated in Dorella at the same infection time point. On the contrary, one of the down-regulated genes in Selenio (WRKY29) was highly up-regulated in Dorella.
Comparison of WRKY gene induction identified in the present study with analyses carried out on M. oryzae-rice interactions [14, 15], highlighted a common activation of WRKY-62, −76, −19, and −50 genes in resistant rice cultivars.
Expression of defence-related genes in Selenio and Dorella in response to F. fujikuroi
Comparison of the transcript profiles for Selenio and Dorella at 1 and 3 wpg reveals a different modulation of the plant response to F. fujikuroi. In detail, Selenio activated, albeit to a moderate extent, only peroxidase genes at 1 wpg (Os07g0677200 and Os01g0326300), while many genes related to glucanases, germin-like proteins, peroxidases, MAP kinases, WRKY transcriptional factors, and PR genes were activated at 3 wpg (Fig. 5). Overall, Selenio up-regulated the highest number WRKY transcription factors, MAPK, and glucanases, and down-regulated the majority of peroxidases at 3 wpg. Furthermore, Dorella activated all groups of defence-related genes at both infection time points (with exception of MAP-kinases at 1 wpg), and it activated the highest number of chitinases, glucanases, and peroxidases and down-regulated the majority of WRKY transcription factors at 3 wpg (Fig. 5). Further analyses focused on the results at 3 wpg, because this infection time point with a higher transcriptomic response in both cultivars, showed to be crucial in the F. fujikuroi-rice interactions.
GO enrichment analyses
As GO terms covering genes with longer transcripts are more likely to be determined as enriched due to higher read counts and increased statistical power for DEG call [50], the goseq R package was implemented in order to limit such length-related biases [51]. Test sets for goseq inputs were DEGs assessed by DESeq2 Goseq output, resulting in 79 and 97 enriched GO terms for Selenio and Dorella, respectively. Specific and common GO terms at 3 wpg, i.e. 46 common, 33 Selenio-specific, and 51 Dorella-specific GO terms, are detailed in Fig. 6.
Among 46 common enriched terms in both rice genotypes, some GO terms were associated to general plant response to fungus attack such as ‘defence response’ (GO:0006952), ‘MAP kinase kinase kinase activity’ (GO:0004709), ‘oxidation-reduction process’ (GO:0055114), ‘response to stress’ (GO:0006950), and ‘response to wounding’ (GO:0009611). On the other hand, some GO terms associated with fungal infection were specific.
Among the GO terms found solely in Selenio, but not in Dorella, there are: ‘response to chitin’ (GO:0010200), ‘brassinosteroid mediated signaling pathway (GO:0009742), ‘response to abscisic acid stimulus’ (GO:0009737), ‘jasmonic acid biosynthetic process’ (GO:0009695), ‘plant-type hypersensitive response (GO:0009626), and ‘transmembrane receptor protein tyrosine kinase signalling pathway’ (GO:0007169).
Although Dorella displayed a higher number of specific GO groups, many of them were not associated with a specific plant response to biotic stress. Furthermore, Dorella activated a different battery of response mechanisms to F. fujikuroi as compared to Selenio, e.g. ‘response to salicylic acid stimulus’ (GO:0009751), and ‘gibberellin metabolic process’ (GO:0009685).
In the following sections, some Selenio-specific GO terms most likely involved in bakanae resistance (hypersensitive response, jasmonic acid biosynthetic process and response to chitin) and two Dorella-specific GO terms associated with bakanae disease (response to salicylic acid and gibberellin metabolism) are described.
Hypersensitive response
The interaction of avirulence (Avr) genes in pathogens and disease resistance (R) genes in plants can cause a localized programmed cell death named hypersensitive response (HR) [52]. The protein Avr9 of Cladosporium fulvum provokes various defence responses in tomato plants bearing the Cf9 gene, and after Avr9 recognition the cells close to the infection site die [53, 54]. The GO term ‘plant-type hypersensitive response’ (GO:0009626), enriched only in Selenio (at 3 wpg), includes defence genes encoding Avr9/Cf-9 elicited protein, C2 domain containing protein, Pto kinase interactor 1 and chitinase 6, many of them related with HR (Additional file 12: Table S12). On the contrary, all genes comprised in this GO term were down-regulated in Dorella with the exception of putative chitin-binding allergen Bra r 2. Thus, activation of a strong hypersensitive response could represent a major differential mechanism conferring resistance to bakanae disease in Selenio.
Jasmonic acid biosynthetic process
Jasmonic acid (JA), salicylic acid (SA), and ethylene (ET) are phytohormones involved in defence processes. Defence against necrotrophic pathogens and leaf-chewing insects is regulated by ET pathway and JA dependent signalling, while biotrophic pathogens mainly activate SA-dependent signalling pathway [55, 56]. In fact, based on recent studies, gibberellins (GAs) stimulate colonization by necrotrophic fungi, such as F. fujikuroi, through suppression of the JA signalling pathway [57].
Genes associated to the GO term ‘jasmonic acid biosynthetic process’ (GO:0009695) were found up-regulated in Selenio (3 wpg) and down-regulated in Dorella (Additional file 13: Table S13). This GO included 12-oxophytodienoic acid reductase, 3-ketoacyl-CoA thiolase, and allene oxide cyclase which participate in the biosynthesis of JA [58–60]. Furthermore, phospholipase D (Os10g0524400 and Os01g0172400) and lipoxygenase-like (Os08g0508800, Os04g0447100, Os08g0509100, and Os03g0700700) genes, necessary for the initial steps of JA synthesis in Arabidopsis thaliana [61], were all up-regulated among DEGs in Selenio at 3 wpg. Thus, activation of JA signalling pathway genes may also play a crucial role in the resistance of Selenio to F. fujikuroi, while fungal infection and down-regulation of JA-biosynthesis related genes probably induced the suppression of this defence mechanism in Dorella and could be favoured by high production of gibberellins during infection and colonization processes.
Response to chitin
The GO term ‘response to chitin’ (GO:0010200) was specifically enriched in Selenio; genes associated to this group include various transcription factors, such as heat-stress, GRAS, and ET-responsive genes which were found up-regulated in Selenio, but not in Dorella.
Heat-stress transcription factors are reported to play a role in fungal resistance in other plant species [62, 63]. Selenio up-regulated the heat-stress transcription factors B-2b, B-2c, and A-2a, possibly as a result of a F. fujikuroi chitin stimulus (Additional file 14: Table S14). Three transcription factors similar to ET-responsive factors (Os05g0572000, Os05g0420300, and Os03g0860100) were up-regulated in Selenio, but not in Dorella. These factors could indicate that ET-mediated signalling pathway may also be involved in Selenio resistance to F. fujikuroi. GRAS transcriptional regulators participate in plant disease resistance and mechanical stress response [64], and the GRAS transcriptional factor Os03g0723000 included in this GO term was found up-regulated in Selenio.
The zinc-finger proteins are known to be involved in rice resistance to pathogens, such as M. oryzae, and Xanthomonas oryzae pv oryzae [65–67]. In this GO group, two zinc-finger RING/FYVE/PHD-type genes (Os03g0240600 and Os01g0755700) were both up-regulated in Selenio and down-regulated in Dorella, and they may contribute to resistance response of Selenio to F. fujiuroi (Additional file 14: Table S14).
Response to salicylic acid stimulus
SA plays an important role in resistance and plant defences against pathogen attacks. R2R3-MYB genes are involved in signalling pathways of salicylic acid [68] and in antagonistic activation of SA defence mechanisms and repression of JA defence mechanisms in A. thaliana [69]. In this study, the GO term ‘response to salicylic acid stimulus’ (GO:0009751) was specifically enriched in Dorella. Dorella activated at 3 wpg two P-type R2R3 MYB genes (Os08g0437300 and Os04g0594100), and other putative MYB genes which were either down-regulated in Selenio or not found among DEGs (Additional file 15: Table S15). This indicated that salicylic signalling pathway could be employed in response of Dorella to F. fujikuroi attack, but not in Selenio.
Gibberellin metabolic process
The GA role in stimulating stem growth in rice plants was discovered for the first time during the studies of bakanae disease [70]. The GAs were reported to be produced not only by the causal agent of bakanae disease inducing the stem elongation, reduced root growth, and inhibition of chlorophyll synthesis, but the plant itself is able to synthesize GAs [71]. In our studies, the enzymes involved in the GA synthesis belonging to the GO term ‘gibberellin metabolic process’ (GO:0009685) were up-regulated in Dorella and down-regulated in Selenio at 3 wpg: gibberellin 2-beta-dioxygenase (Os05g0560900), gibberellin 20 oxidase 1 (Os03g0856700), and syn-copalyl diphosphate synthase (Os04g0178300) (Additional file 16: Table S16). Moreover, putative cytochrome P450 dwarf3 gene (Os06g0110000) belonging to this GO term was up-regulated in Dorella, and not included as DEG in Selenio. It was reported that dwarf3 gene encoding a cytochrome P450-type enzyme (or OsKAO gene) is involved in early steps of gibberellin biosynthesis in rice [72]. Interestingly, the dwarf3 gene is involved also in leaf longevity in rice during the stress conditions, and its expression in Dorella is possibly related with pronounced leaf senescence and cell death in response to F. fujikuroi beside the gibberellin metabolism [73].
The activation of gibberellin metabolic genes in the susceptible rice genotype, where typical bakanae symptoms and active growth of pathogen were observed, was in accordance with GA3 biosynthesis measured through high performance liquid chromatography (HPLC) analysis. As shown in Fig. 7, Dorella inoculated by F. fujikuroi had twelve times more GA3 compared to inoculated Selenio at 3 wpg. On the basis of these results, it is evident that the gibberellin group of growth hormones plays important role in response of rice to F. fujikuroi. By modulating the production of gibberellins, rice plants are coordinating their level of tolerance to bakanae disease possibly through interaction of the GA-signalling molecules with components of the JA signalling pathway, as already reported under other stress conditions in plants [74], and confirmed in our work by enriched JA metabolism in Selenio.
KEGG and MapMan maps
In order to compare and summarize the response of the two genotypes to infection, the genes were mapped to KEGG rice plant-pathogen interaction diagram (osa04626; http://www.genome.ad.jp/kegg/) [75]. Fold change values for the two contrasts Dorella mock vs Dorella inoculated (3 wpg) and Selenio mock vs Selenio inoculated (3 wpg) were integrated in the diagram by colour-coding, implementing the pathview R package [76]. This allowed a prompt comparative evaluation of expression responses (Fig. 8). Boxes are placeholders for one or more genes assigned to the same KEGG orthology (KO) group. When more than one gene is mapped to the same group, expression fold-changes were summed up according to default pathview settings.
Based on the plant-pathogen interaction map, when compared to Dorella, Selenio exhibited higher up-regulation of genes in pathways involved in hypersensitive response, response to chitin and JA-dependent signalling confirming the results of the GO enrichment analyses. As an example, following fungal Avr9 and plant Cf9 interaction, Ca2+-dependent protein kinase gene (CDPK, e.g. Os01g0808400) involved in hypersensitive response and signal transduction was found up-regulated in Selenio, and down-regulated in Dorella (Fig. 8).
Chitin is a typical microbe-associated molecular pattern (MAMP) molecule from fungal cell walls which elicits plant immune responses. Two plasma membrane proteins, OsCEBiP and OsCERK1 were found necessary to regulate chitin as elicitor of signalling in rice [77, 78]. As summarized in the KEGG map, CEBiP and CERK1 genes were induced in Selenio, but not in Dorella inducing subsequently up-regulation of RIN4 (e.g. Os03g0848600) and HSP90 (e.g.Os06g0716700)-like genes (similar genes are included in these KO groups and fold changes are summed up in the boxes), downstream components of the hypersensitive response (Fig. 8).
Furthermore, JASMONATE ZIM-DOMAIN-like (JAZ; e.g. Os03g0180800) and MYC2 (Os10g0575000), two genes involved in JA signalling, were found by KEGG annotation induced exclusively in Selenio (Fig. 8) which emphasizes the possible role of JA-dependent signalling in Selenio resistance to F. fujikuroi.
Additionally, employing MapMan software [79], we presented a summary of the main expression changes putatively involved in biotic stress (F. fujikuroi) in Selenio and Dorella, as shown by Fig. 9. In contrast to Dorella, Selenio up-regulated the majority of DEGs involved in the candidate resistance pathways already described above, including jasmonic- and ethylene-mediated signalling.
Validation of RNA-Seq technology
To validate the RNA-Seq technique, seven DEGs were selected based on their expression patterns at 3 wpg for quantitative RT-PCR (qRT-PCR) by using the same RNA extracts as for RNA-seq experiments (Fig. 10). The results of the selected DEGs showed that the qRT-PCR was consistent with the RNA-Seq results showing similar expression pattern of up- and down-regulated genes by using both, RNA-Seq and qRT-PCR, analyses.