Comparative Transcriptomic Analysis of Virulence Factors in Leptosphaeria maculans during Compatible and Incompatible Interactions with Canola

Leptosphaeria maculans is a hemibiotrophic fungus that causes blackleg of canola (Brassica napus), one of the most devastating diseases of this crop. In the present study, transcriptome profiling of L. maculans was performed in an effort to understand and define the pathogenicity genes that govern both the biotrophic and the necrotrophic phase of the fungus, as well as those that separate a compatible from an incompatible interaction. For this purpose, comparative RNA-seq analyses were performed on L. maculans isolate D5 at four different time points following inoculation on susceptible cultivar Topas-DH16516 or resistant introgression line Topas-Rlm2. Analysis of 1.6 billion Illumina reads readily identified differentially expressed genes that were over represented by candidate secretory effector proteins, CAZymes, and other pathogenicity genes. Comparisons between the compatible and incompatible interactions led to the identification of 28 effector proteins whose chronology and level of expression suggested a role in the establishment and maintenance of biotrophy with the plant. These included all known Avr genes of isolate D5 along with eight newly characterized effectors. In addition, another 15 effector proteins were found to be exclusively expressed during the necrotrophic phase of the fungus, which supports the concept that L. maculans has a separate and distinct arsenal contributing to each phase. As for CAZymes, they were often highly expressed at 3 dpi but with no difference in expression between the compatible and incompatible interactions, indicating that other factors were necessary to determine the outcome of the interaction. However, their significantly higher expression at 11 dpi in the compatible interaction confirmed that they contributed to the necrotrophic phase of the fungus. A notable exception was LysM genes whose high expression was singularly observed on the susceptible host at 7 dpi. In the case of TFs, their higher expression at 7 and 11 dpi on susceptible Topas support an important role in regulating the genes involved in the different pathogenic phases of L. maculans. In conclusion, comparison of the transcriptome of L. maculans during compatible and incompatible interactions has led to the identification of key pathogenicity genes that regulate not only the fate of the interaction but also lifestyle transitions of the fungus.

Introduction

Blackleg disease (stem canker) caused by Leptosphaeria maculans (Desm.) Ces. & De Not. is one of the major constraints to canola (Brassica napus L.) production worldwide (Fitt et al., 2006). Infection by the fungus is known to cause more than 50% yield losses in canola (Kutcher et al., 2013). The major difficulty for combating the pathogen lies in the understanding of its complex lifestyle, which includes alternative biotrophic, and necrotrophic, phases, along with a symptomless endophytic phase (Howlett et al., 2001; Van de Wouw et al., 2016). Management of blackleg disease includes crop rotations, seed treatment and fungicide applications, and preferably, disease-resistant cultivars, arguably the most effective approach (Delourme et al., 2006).

Canola shows two types of resistance against L. maculans: qualitative and quantitative. Single or few major genes that are known to be involved in a gene for gene interaction govern qualitative resistance. Major genes provide resistance particularly at the seedling stage whereas quantitative resistance involves many small effect genes that are mostly expressed during the adult plant stage (Raman et al., 2012; Huang et al., 2014). To date, several resistance genes with major effect have been identified in Brassica species, but only two, LepR3 and Rlm2, have been cloned and well characterized (Delourme et al., 2006; Long et al., 2011; Larkan et al., 2013, 2015; Van de Wouw et al., 2014). On the other hand, better progress has been achieved with L. maculans where 14 avirulence genes have been identified, and seven of them, namely AvrLm1, AvrLm4-7, AvrLm6, AvrLm11, AvrLmJ1, AvrLm2, and AvrLm3 have been cloned (Gout et al., 2006; Fudal et al., 2007; Parlange et al., 2009; Balesdent et al., 2013; Van de Wouw et al., 2014, 2016; Ghanbarnia et al., 2015; Plissonneau et al., 2016). Interestingly, some of these avirulence genes have been found to be clustered, with clusters AvrLm1-2-6 and AvrLm3-4-7-9-AvrLepR1 being the notable examples (Balesdent et al., 2002; Ghanbarnia et al., 2012). For the most part, avirulence genes, including L. maculans Avrs, are small-secreted proteins (SSPs) with several cysteine residues, and often referred to as effectors (Stukenbrock and McDonald, 2009; Rouxel et al., 2011).

Effectors are key elements in fungal virulence against plants and particularly important during the biotrophic phase of infection (Kloppholz et al., 2011). L. maculans, being a hemibiotroph, will initially rely on effectors to suppress plant defenses, and then will subsequently use effectors to kill plant cells. In L. maculans, most putative or candidate effectors are localized in transposon-rich repetitive DNA and are affected by a repeat-induced point mutation (Rouxel et al., 2011). The putative effector genes are mostly over-expressed during primary leaf infection (Soyer et al., 2014). Such information about the genomic localization, gene organization, and expression dynamics is helpful to understand the host-pathogen interaction and more particularly for the identification of bona fide effectors. Similarly, transcription factors (TFs), and carbohydrate active enzymes (CAZymes) are known to play a pivotal role in host-pathogen interactions, and are, along with effectors, prime targets for studying virulence factors in fungi (Guo et al., 2011; Lombard et al., 2014; Lowe et al., 2014; Malinovsky et al., 2014).

Transcription factors are essential players in the signal transduction pathways. In L. maculans, TF LmStuA is found to be required for normal growth, perithecium formation, pathogenicity on oilseed rape leaves, and expression of effectors (Soyer et al., 2015). The silencing of LmStuA triggers drastic effects on the morphogenesis and pathogenicity of L. maculans, indicating that it may affect a large number of genes and pathways (Soyer et al., 2015). Similarly, several CAZymes in L. maculans genome have been predicted to have a functional role in pathogenesis (Lowe et al., 2014). CAZymes are important to break down the polysaccharides of plant cell walls, to establish infection, and also, to facilitate access to nutrients during the necrotrophic and saprophytic growth phases. For instance, global transcriptomic analyses of the hemibiotroph Colletotrichum higginsianum revealed that genes encoding secreted proteins without a functional annotation are expressed predominantly during the initial biotrophic phase, whereas expression of secreted lytic enzymes (including CAZymes) was higher in the subsequent necrotrophic phase (O'Connell et al., 2012). A similar finding was observed in Leptosphaeria biglobosa, a necrotroph expressing more cell wall degrading genes than L. maculans (Lowe et al., 2014). However, L. maculans expressed many genes in the carbohydrate binding module (CBM) class of CAZymes, particularly CBM50 genes, during early infection, and cell wall degrading enzymes at later stages of growth (Lowe et al., 2014). This suggests that expression of secreted proteins without functional annotation is a general feature of biotrophy, whereas expression of cell wall degrading enzymes is generally associated with necrotrophy. Other important necrotrophy-related genes code for sirodesmin PL (Sir), a phytotoxin that belongs to the class of epipolythiodioxopiperazine (ETP). The production of sirodesmin by L. maculans is thought to be suppressed by brassinin, a phytoalexin of canola (Pedras et al., 1993). In L. maculans, a cluster of 23 genes including polyketide synthase (PKS), non-ribosomal peptide synthase (NRPs) genes, and 18 Sir genes have been identified (Gardiner et al., 2004).

The infection process is highly dependent on host recognition and molecular cross-talk between the host and the pathogen where pathogenicity-related genes play an important role. However, any given host-pathogen interaction is a very complex phenomenon, which makes it difficult to understand the factors dictating compatibility or incompatibility. In Arabidopsis, Huibers et al. (2009) were able to discriminate genes induced during compatible from incompatible interactions with the downy mildew pathogen Hyaloperonospora arabidopsidis (downy mildew). Similarly, other studies have been conducted to compare gene expression profiling under compatible and incompatible interactions (Wang et al., 2010; Sestili et al., 2011; Li et al., 2015).

Expression profiling during the development of disease is an effective approach to better understand the pathogenesis process. Current improvements in sequencing technologies have provided new opportunities to evaluate gene expression under different conditions by sequencing the entire transcriptome (Wang et al., 2009). Compared to other hybridization based transcriptome profiling platforms like microarrays, RNA-seq provides access to simultaneous transcript discovery and abundance estimation, identification of differentially-expressed genes (DEGs) and associated molecular cellular pathways, and alternative splicing variants (Wang et al., 2009; Trapnell et al., 2012). However, RNA-seq analysis requires scalable, fast, and statistically relevant software tools that can handle complex and large sequence data. Fortunately, considerable efforts have been devoted to the development of specialized software tools to perform effective RNA-seq analysis (Garber et al., 2011; Trapnell et al., 2012; Seyednasrollah et al., 2015; Sonah et al., 2016).

The L. maculans-canola compatible interaction has been previously addressed using the RNA-seq approach (Lowe et al., 2014; Haddadi et al., 2016). These studies made significant efforts toward the understanding of susceptible reaction in canola after L. maculans infection and provides the first in-depth look into the transcriptomic profile of this interaction. In this study, our objective was to compare the pathogen responses against susceptible and resistant host genotypes in order to obtain a precise definition of the virulence factors expressed by L. maculans during its biotrophic and necrotrophic phase. For this purpose, we performed the RNA-seq transcriptome profiling of L. maculans inoculated on susceptible and resistant canola lines at four developmental stages over five biological replications. A particular emphasis was placed on the identification of DEGs including effectors, CAZymes, and other pathogenesis-related genes during blackleg disease development in canola.

Materials and Methods

Plant Material and L. maculans Inoculation

Canola (B. napus) breeding lines, Topas DH16516 (Topas-wild), a double haploid line susceptible to L. maculans and Topas-Rlm2, an introgression line resistant to L. maculans isolates carrying AvrLm2 (Larkan et al., 2016a,b) were used as plant material. Seven-day old seedlings of both resistant and susceptible canola lines grown under controlled environment were point inoculated with pycnidiospores suspension of L. maculans isolate D5. L. maculans isolate D5 contains known avirulence effectors AvrLm1, AvrLm4-7, AvrLm2, and AvrLmJ1 but lacks AvrLm6 (Raman et al., 2012). Five biological replicates each with leaves from eight plants per sample were collected at 3, 5, 7, and 11 days post inoculation (dpi). The L. maculans isolate D5 was cultured in V8^® agar plates at 25°C. For RNA extraction, two fungal agar plugs were inoculated in V8^® liquid media and incubated at 22°C. After 7 days, the mycelium was harvested for RNA extraction.

Sample Preparation and Illumina Sequencing

Total RNA for all samples was extracted using a combined Trizol/Qiagen RNeasy mini kit (Qiagen, Mississauga, ON, Canada). The RNA quality and quantity was accessed using NanoDrop ND-1000 spectrophotometer (NanoDrop technologies) and further verified by agarose gel electrophoresis. Four micrograms of total RNA from all samples were used to make individual barcoded cDNA library using TruSeq RNA library preparation kit v2 with some modifications. Individual library was assessed using the Agilent 2100 Bioanalyzer (Agilent Technologies). The products resulted in a smear with an average fragment size of approximately 260 bp. A total of six individual libraries were pooled using the uniform amount of each library and the quality of the final library pools was also assessed. Single end, 100 bp sequencing was performed using an Illumina HiSeq 2000 sequencer at the Génome Québec Innovation Centre (McGill University, Quebec, Canada).

RNA-seq Data Analysis

RNA-seq reads were quality-checked with fastqc, which performs various quality checks for the raw reads. Read processing was performed by using Trimmomatic software (Bolger et al., 2014). Reads were trimmed from both ends until the average of all 5 bp sliding windows reached a Phred score of 25 or higher and all the sequences shorter than 35 bases were discarded. Processed reads were aligned to the L. maculans genome and transcriptome with Tophat2 (Trapnell et al., 2009). Most of the parameters in Tophat was set as default. A mismatch of two bases were allowed for the alignment. The minimum and maximum intron length was set to 50 and 500,000 respectively. Reads aligned to multiple sites were removed prior to further analysis.

Novel transcripts that did not overlap with any annotated transcripts were identified using Cufflink tools. The gene expression level for each annotated as well as non-annotated novel transcripts were estimated as the number of Fragments (reads) per kilobase of transcript per million mapped reads (FPKM), considering only uniquely mapped reads in exonic region by using the Cufflink software. The differentially expressed genes (DEGs) were identified using four different tools including Cuffdiff (Trapnell et al., 2012), EdgeR (Robinson et al., 2010), DESeq2 (Love et al., 2014), and CLC Genomics Workbench. We used FDR < 0.0001 and the absolute value of log₂ (Fold-change) >1.5 relative to axenic culture as the threshold for the identification of DEGs. We used quartile normalization that excludes the top 25% of expressed genes to improve detection of less abundant genes. With the Cufflink software, we used -M option to mask rRNA, -b, and -u option for bias correction and option to normalize for aligned tags instead of total tags. HTseq tool was used to count reads prior to DEGs identification. The DESeq2 package was used to estimate sample quality (PCA) and the expression level of transcripts. The regularized rlog transformation and variance stabilizing transformation were used for data visualization. Time course analysis for all the four time points were carried out to find the genes that reacted in a time-specific manner using DESeq2 package.

Functional Annotation and Gene Ontology

Standard gene ontology (GO) was used to describe DE gene functionality, a hypergeometric test and the p < 0.05 of Pearson Chi-Square test between the gene numbers of the two input dataset were used to map the DT genes to GO terms based on the BGI WEGO (Web Gene Ontology Annotation Plot, http://wego.genomics.org.cn/cgi-bin/wego/index.pl). Single enrichment analysis (SEA), a function of AgriGO was used to examine GO term enrichment using Magnaporthe grisea as a background reference using Fisher statistical test and 0.05 significance level and other default parameters. Transcription factors were identified using Fungal Transcription Factor Database (FTFD) (Park et al., 2008). The FTFD pipeline sorts fungal TFs initially based on the relevant InterPro terms like DNA-binding motifs, and then the false-positive TFs are filtered with different criteria. CAZymes were identified using dbCAN server (Yin et al., 2012). The dbCAN hosted an analytical pipeline compiled with CDD (conserved domain database) search, family specific hidden Markov model and literature curation (Yin et al., 2012). For the classification of putative secreted peptidases, the sequences for the secreted proteins predicted by WoLF PSORT (cutoff score = 15) were submitted to MEROPS Batch Blast analysis (http://merops.sanger.ac.uk, Rawlings et al., 2014).

Identification of Candidate Secreted Effector Proteins (CSEPs)

Predicted protein sequences from the L. maculans genome were retrieved from the JGI MycoCosm (http://genome.jgi.doe.gov/Lepmu1/Lepmu1.home.html) (Rouxel et al., 2011; Grigoriev et al., 2014). SignalP (cutoff probability = 0.8), TargetP (cutoff probability = 0.8), Psort (cutoff score = 15), BlastP (cutoff e = 1e⁻⁰⁵), and TMHMM software tools along with Secretool pipeline were used to predict small-secretory proteins (Petersen et al., 2011; Cortázar et al., 2014). To prioritize candidate effector genes, the entire secretome was analyzed by EffectorP a machine learning method optimized for the prediction of fungal effectors (Sperschneider et al., 2015). Crinkler type effectors were also searched in the L. maculans genome with the FEMO software tool implemented in MEME suit by using conserved domain LFLAK (http://meme-suite.org/tools/fimo). Proteins with LFLAK-domain within the initial 100 AA were considered for further analysis. To verify the search parameters, Crinkler effector search was also performed in Phytophthora sojae genome and compared with earlier report by Haas et al. (2009). After an initial search with FIMO, candidate Crinkler effectors were analyzed with secretool pipeline and subsequently with effectorP software tool.

Results

Disease Progress on Compatible and Incompatible Host

On plants of either the compatible host Topas-wild or the introgression line Topas-Rlm2 carrying a major resistance gene, disease symptoms were not observed until 7 dpi. After this asymptomatic early growth stage, symptoms became visible exclusively on Topas-wild plants expanding into clear lesions and zones of chlorotic tissues at 11 dpi. In the case of Topas-Rlm2 plants, there was no visible lesion beyond the site of inoculation (Figure 1). The disease progress on susceptible Topas-wild cotyledons caused by L. maculans was similar to that previously described on Westar (Lowe et al., 2014).

FIGURE 1

Figure 1. Disease symptoms at 11 days post inoculation on leaves of compatible (Topas-wild) and incompatible (Topas-Rlm2) canola host inoculated with Leptosphaeria maculans D5 isolate or water.

Transcriptome Sequencing with the Time Course of Disease Progress

A total of 1.6 billion single-end reads consisting of approximately an average of 33 million reads for each cDNA library were obtained. The entire data set was submitted to SRA NCBI database and can be accessed with accession number SRP078092. The raw reads obtained were uniform across the different libraries. After performing quality assessment and read processing, about 0.1% of the reads with poor quality and shorter length were discarded. Mapping of the processed reads to the L. maculans genome showed a very high percentage of mapping for the axenic samples (Supplementary Table 1). The percentage of reads from in planta samples mapped to the L. maculans reference genome increased over time from <0.5% at 3 dpi to about 12% at 11 dpi in compatible host Topas-wild (Figure 2, Supplementary Table 1). By contrast, mapped fungal reads remain stable in the incompatible interaction with Topas-Rlm2 and were well below 1% even at 11 dpi (Figure 2, Supplementary Table 1), a result well in line with the phenotype observed in Figure 1. Comparing the percent of mapped reads in the early asymptomatic growth stage in Topas-wild plants, only a marginal increase was observed between 3 and 5 dpi. This percentage more than doubled between 5 and 7 dpi and registered a 10-fold increase from 7 to 11 dpi (Supplementary Table 1).

FIGURE 2

Figure 2. Percentage of mapped reads to the Leptosphaeria maculans genome at 3, 5, 7, and 11 days post inoculation (dpi). RNA-seq was performed for L. maculans isolate D5 inoculated to Topas-wild (compatible) and Topas-Rlm2 (incompatible) canola genotypes.

The principal components (PCs) analysis highlighted a clear differential effect of the treatments along with the uniformity of the five biological replications within a treatment (Figure 3). The first principal component explained 58% of the expression variation supporting the large phenotypic differences between the conditions over the time period.

FIGURE 3

Figure 3. Principal components (PCs) of RNA-seq expression profile obtained for Leptosphaeria maculans in planta samples during compatible and incompatible interactions. A total of five biological replications for 3, 5, 7, and 11 days post inoculation (dpi) were used for RNA-seq expression profiling of L. maculans inoculated to compatible host Topas-wild (blue) and incompatible host Topas-Rlm2 (Red).

Comparison of Gene Expression Profiling during In vitro and In planta Growth of L. maculans

The number of DEGs in L. maculans between in vitro axenic samples and samples at different in planta growth stages showed an interesting pattern of gene expression turnover over time. The overall pattern of DEGs identified with different software tools was similar (Supplementary Figure 1). The highest number of DEGs for both the compatible (Topas-wild) and incompatible (Topas-Rlm2) interactions was recorded in the early events (3 dpi) and was similar in both cases (Figure 4). This number was reduced by around three-fold at 5 dpi in both interactions and remained fairly level over the next sampling times in Topas-Rlm2 plants. On the other hand, it increased steadily in Topas-wild plants to exceed by roughly three times the number of DEGs found in Topas-Rlm2 plants (Figure 4).

FIGURE 4

Figure 4. Differentially expressed Leptosphaeria maculans genes identified at 3, 5, 7, and 11 dpi during compatible and incompatible interactions. Differentially expressed genes were identified by comparing in planta samples with axenic cultures. Analyses were performed with five biological replicates and a threshold value of FDR < 0.0001 and log₂ fold-change > 1.5.

When looking at the top 20 upregulated DEGs at 3 dpi in the compatible interaction (Topas-wild), most of them were linked to uncharacterized proteins for which functional annotation was not available (Table 1). Since the uncharacterized proteins are unique to L. maculans, homology-based annotation failed to characterize them. Genes involved in chitin binding, fasciclin, and related adhesion glycoproteins known to play roles in early stages of infection were found to be highly upregulated (Table 1). Interestingly, of the top 20 upregulated genes of L. maculans identified during the compatible interaction with Topas-wild, 18 were also found associated with the incompatible interaction (Topas-Rlm2) (Table 1, Supplementary Table 2). However, the number of common genes gradually decreased over time to 11, 6, and 4 at 5, 7, and 11 dpi, respectively (Table 1, Supplementary Table 2). Time series analysis showed high expression of Avr genes during the biotrophic phase at 7 dpi while most of the highly expressed genes during the necrotrophic phase at 11 dpi were associated with molecular functions involved in catalase activity, hydrolases, CAZymes, peptidases, and transporters (Table 1, Figure 5).

TABLE 1

Table 1. List of the top 20 upregulated Leptosphaeria maculans genes during the compatible interaction with canola cultivar Topas-wild at each of the four sampling points.

FIGURE 5

Figure 5. Heat map of time series analysis showing expression pattern of the top 22 genes in Leptosphaeria maculans during compatible interaction with host Topas-wild. Gene_2397, gene_10780, and gene_11099 represents AvrLm4-7, AvrLm1, and Lm5LysM, respectively. Analyses were performed with five biological replicates.

Expression of CSEPs in L. maculans during Compatible and Incompatible Interactions

A total of 552 classically secreted proteins were identified using computational pipeline comprising SignalP and TMHMM software tools along with Secretool pipeline. Following further analyses, 134 genes were prioritized as high confidence CSEPs based on the results obtained with EffectorP software. In the case of the compatible interaction with Topas-wild, an important increase in upregulated CSEPs was observed at 7 and 11 dpi, while very limited differences over time were noted in the incompatible interaction (Figure 6).

FIGURE 6

Figure 6. Number of upregulated candidate small secretory effector proteins (CSEPs) identified by RNA-seq analyses performed in Leptosphaeria maculans at different stages of disease development during compatible and incompatible interactions. Upregulated genes were identified by comparing in planta samples at different growth stages with axenic cultures. dpi, days post inoculation. Analyses were performed with five biological replicates and a threshold value of FDR < 0.0001 and log₂ fold-change > 1.5.

Out of the 134 CSEPs, 35 genes were not expressed in either the compatible or incompatible host. However, a total of 28 genes showed the highest expression level at 7 dpi in the compatible interaction, which would link them to the biotrophic phase of the fungus (Table 2). This profile was similar to that of known Avr genes, which were mostly found to be highly expressed at 7 dpi in the compatible interaction, with the obvious exception of AvrLm6 that remained unchanged throughout the sampling periods (Figures 5, 7). The expression patterns of AvrLm1, AvrLm2, AvrLm4-7, AvrLm11, and AvrLmJ1 were similar in the compatible interaction as their expression increased over time to reach a peak at 7 dpi, and subsided at 11 dpi. As for the incompatible interaction, limited expression was observed for all known Avr genes with no distinctive pattern (Table 2, Figure 7). Of additional significance, an expression pattern similar to Avr genes was found for 23 other effectors in the compatible interaction, including eight specific to L. maculans and not previously reported as effectors (Table 2). Of these, gene_6114 and gene_2728 are particularly interesting for their differentially higher expression at 7 dpi in Topas-wild. Finally, gene_9004 showed the highest expression at 7 dpi during the compatible interaction at a level of almost 24,000 FPKM.

TABLE 2

Table 2. Expression pattern of 28 effector proteins, including known Avr genes, during Leptosphaeria maculans compatible and incompatible interactions with canola.

FIGURE 7

Figure 7. Expression pattern of known Avr genes at different disease developmental stages in Leptosphaeria maculans during compatible and incompatible interactions. FPKM, Fragments per kilo-base of transcript per million mapped reads; dpi, days post inoculation. Analyses were performed with five biological replicates.

In addition to the effectors listed in Table 2, 15 other genes identified as CSEPs showed higher levels of expression exclusively during the compatible interaction at 11 dpi when the fungus had entered its necrotrophic phase (Table 3). Furthermore, the crinkler effector search performed initially in P. sojae identified exactly the same set of proteins reported earlier by Haas et al. (2009) (Supplementary Table 3). This validated the method used here to identify crinkler effectors in L. maculans genome. The initial search identified 63 proteins with LFLAK-like domain present within the initial 100 AA (Supplementary Table 4). Further analysis with secretool revealed only five proteins and EffectorP confirmed only one candidate (Gene_2728, Lema_T080290.1) as a crinkler effector (Gene_2728, Lema_T080290.1). Gene_2728 showed the highest expression at 7 dpi during the compatible interaction (Supplementary Figure 2).

TABLE 3

Table 3. Expression pattern of 15 effector proteins highly expressed during the necrotrophic phase of a compatible interaction of Leptosphaeria maculans with canola.

Transcription Factors Associated with Effector Expression

TFs regulate the expression of a number of genes simultaneously, and their upregulation during in vivo conditions was indicative of their importance in the outcome of the interaction. Both compatible and incompatible interactions were regulated by a large number of upregulated TFs in L. maculans in the early events (3 dpi) followed by a sharp decline at 5 dpi (Figure 8). As with the phenotypes, a clear distinction between compatible and incompatible interactions appeared at 7 dpi where nearly 100 TFs were upregulated in the former case compared to 40 in the latter, and these differences were even more manifest at 11 dpi (Supplementary Table 5). Among TFs, the LmStuA TF (gene_1191), a member of APSES domain containing proteins was found to be strongly upregulated at 7 dpi in the compatible interaction with Topas-wild, thus suggesting its importance in the infection process. Mostly Zn2Cys6 TFs were upregulated during the compatible interaction (Supplementary Table 5).

FIGURE 8

Figure 8. Number of transcription factors (TFs) upregulated in Leptosphaeria maculans at 3, 5, 7, and 11 days post inoculation (dpi) on compatible host Topas-wild and incompatible host Topas-Rlm2. Upregulated genes were identified by comparing in planta samples at different growth stages with axenic cultures. Analyses were performed with five biological replicates and a threshold value of FDR < 0.0001 and log₂ fold-change > 1.5.

Expression of CAZymes in L. maculans during Compatible and Incompatible Interactions

Plant-fungal interactions involve a variety of CAZymes that are used by the fungal pathogen to infect its host. RNA-seq data for L. maculans infecting compatible host Topas-wild and incompatible host Topas-Rlm2 showed a high number of differentially expressed CAZymes at 3 dpi (Figure 9). At 5 dpi, this number dropped drastically. It remained fairly constant thereafter in the case of the incompatible interaction (Figure 9). By contrast, we observed a rapid increase in differentially expressed CAZymes from 5 to 7 dpi in the compatible interaction, and this number exceeded 200 at 11 dpi, supporting their role in the necrotrophic phase of the fungus.

FIGURE 9

Figure 9. Number of carbohydrate active enzymes (CAZymes) upregulated at 3, 5, 7, and 11 days post inoculation (dpi) of Leptosphaeria maculans inoculated to compatible host Topas-wild and in-compatible host Topas-Rlm2. Upregulated genes were identified by comparing in planta samples with axenic cultures. Analyses were performed with five biological replicates and a threshold value of FDR < 0.0001 and log₂ fold-change > 1.5.

Among the differentially expressed CAZymes, glycosyl hydrolase (GH), carbohydrate esterase (CE) domain containing genes were more prevalent followed by those belonging to auxiliary activities (AAs) class (Supplementary Table 6). During the early stages (3 dpi) of L. maculans infection, cellulose and pectin-degrading enzymes such as PL3, CE4, GH43, and GH3 were upregulated in both interactions compared to axenic growth. At the later stage (11 dpi), CE1, PL1, GT34, and GH28 were only upregulated in the compatible interaction. Carbohydrate binding molecules (CBM) were another prevalent group among DEGs. In L. maculans, enzymes with LysM motifs are a well-studied class of CAZymes. Among the LysM domain containing genes present in L. maculans genome, Lm2LysM (Gene_4592) and Lm5LysM (Gene_11099) genes showed a higher level of expression at 7 dpi during the compatible interaction (Supplementary Figure 3). However, Lm4LysM (Gene_7646) gene was not expressed during either the compatible or incompatible interaction. The expression pattern for some of the CAZymes at different growth stages of L. maculans during compatible and incompatible interaction is shown in Supplementary Figure 4.

Expression of Important Peptidases, Secondary Metabolites, and Necrosis-Inducing Proteins (NIPs)

At the early stages of infection, a similar number of differentially expressed peptidases were observed in both interactions. However, as the fungus transitioned from the biotrophic phase at 7 dpi to the necrotrophic phase at 11 dpi, the difference in the number of differentially expressed peptidases between the compatible and incompatible interaction steadily increases to a ratio of over 60 to 1 at 11 dpi (Figure 10, Supplementary Table 7). Among the different classes of peptidases, serine proteases and carboxypeptidases were highly upregulated during L. maculans compatible interaction (Figure 11, Supplementary Table 7).

FIGURE 10

Figure 10. Number of upregulated peptidases observed at 3, 5, 7, and 11 days post inoculation (dpi) of Leptosphaeria maculans inoculated to compatible host Topas-wild and in-compatible host Topas-Rlm2. Upregulated genes were identified by comparing in planta samples with axenic cultures. Analyses were performed with five biological replicates and a threshold value of FDR < 0.0001 and log₂ fold-change > 1.5.

FIGURE 11

Figure 11. Functional categorization of upregulated genes at 7 days post inoculation (dpi) during compatible and incompatible canola-Leptosphaeria maculans interactions. The right y-axis indicates the number of genes in a category. The left y-axis indicates the percentage of a specific category of genes in the main category. Upregulated genes were identified by comparing in planta samples with axenic cultures. Analyses were performed with five biological replicates and a threshold value of FDR < 0.0001 and log₂ fold-change > 1.5.

In L. maculans, a cluster of 23 genes including PKS, NRPs genes and 18 Sir genes have been identified (Gardiner et al., 2004). Most of these genes showed a similar pattern of expression at all the stages of infection with a higher expression at 11 dpi compared to early growth stages of L. maculans in planta (Supplementary Figure 5A). By contrast, Sir O (oxidoreductase) and Sir T (thioredoxin reductase) showed higher level of expression at both 3 dpi and 11 dpi (Supplementary Figure 5A). However, comparative analyses indicated that there was no difference in expression at 3 dpi between the interactions (Supplementary Figure 5B). Most of the NRPs showed high level of expression at 3 dpi and 11 dpi but NRP gene_904 showed high expression at 7 dpi (FPKM 876) that receded at 11 dpi (FPKM 112.2). For the two necrosis-and ethylene inducing proteins (Neps) in L. maculans, Nep1-like proteins (NPP1, gene_11090) showed higher expression during the pathogenic phases of the fungus (Supplementary Table 8).

Functional Categorization and Gene Ontology (GO) Enrichment of the DEGs

Differential molecular responses in L. maculans under compatible and incompatible interactions were observed based on the functional categorization of DEGs (Figure 11, Supplementary Figure 6). The DEGs were grouped into three main classes, molecular function, cellular component, and biological process, and several subclasses based on gene ontology terms assigned using WEGO tool. During early asymptomatic stage at 3 dpi, upregulated genes with similar functional categories were observed under both conditions (Supplementary Figure 6). However, the molecular responses of L. maculans between the compatible and incompatible hosts differed significantly over time, and were more prominent at 7 and 11 dpi (Figure 11; Supplementary Figure 6). The most striking differences were found in the number of genes expressed in each sub-category and in the prevalence of functional categories linked to a transition between the biotrophic and necrotrophic phases such as carboxypeptidases, ligase, peroxidases, cell-wall biogenesis, and stimuli (Figure 11).

Hierarchical clustering and GO enrichment of DEGs during the compatible interaction was also performed using AgriGO tool. Most of the DEGs were found to belong to binding and catalytic functions (Supplementary Figure 7). Localization (GO:0051179, 3.21e-62), cellular process (GO:0009987, 2.37e-258), signal transduction (GO:0007165, 3.48e-13), response to stress (GO:0006950, 8.44e-11), and transcription regulator activity (GO:0030528, 3.3e-15), were the most significantly enriched GO terms at 3 and 5 dpi. At 7 dpi, localization (GO:0051179, 2.35e-61), transport (GO:0006810, 1.68e-61), catabolic process (GO:0009056, 1.03e-11), metabolic process (GO:0008152, 0), gene expression (GO:0010467, 3.93e-38), binding (GO:0005488, 6.84e-319), and transferase activity (GO:0016741, 8.4e-09) were linked to the highest number of genes, and at 11 dpi, hydrolytic activity, mostly monooxygenase activity (GO:0051179, 5.86e-09), transcription regulator activity (GO:0030528, 2.49e-21), transcription factor activity (GO:0051179, 3.14e-14), and regulation of biological process (GO:0050789, 7.09e-26) were the most common enriched GO terms (Supplementary Figure 7).

Discussion

In the present study, we performed a comprehensive analysis of the L. maculans transcriptome profile during compatible and incompatible interactions with canola. Based on comparative analyses, key genes that dictate both the interaction between canola and L. maculans and the different pathogenic stages of the fungus were highlighted. Among the genes of particular significance, our results have identified candidate effectors, TFs, CAZymes, peptidases, and other pathogenesis-related genes that are specifically upregulated as L. maculans initiates a biotrophic interaction with the plant, and transitions to a necrotrophic phase. The differential expression of genes during compatible and incompatible interactions offers a precise insight into the mechanisms of pathogenesis in the L. maculans-canola interactions.

The effectiveness of transcriptome analyses in plant-pathogen interactions depends primarily on the approach of expression quantification, statistical methods avoiding possible errors and reliance on normalized comparisons, experimental design with sufficient replications and appropriate plant or pathogen material to address relevant biological questions (Williams et al., 2014). In this work, introgression lines Topas-wild and Topas-Rlm2 were used to study L. maculans molecular responses during disease development. Topas-wild is a common cultivar from Canada well known and exploited for its susceptibility to L. maculans isolate D5 (Larkan et al., 2013). On the other hand, Topas-Rlm2 is a recently developed cultivar that carries Rlm2, a major resistance gene that prevents infection from L. maculans isolate D5 (Larkan et al., 2015). As such, this cultivar provided a unique opportunity to investigate the subtle elements that distinguish the ability of L. maculans to infect or not its host. The importance of obtaining good reproducible phenotypes for transcriptomic analyses cannot be overstated as it remains the reference basis for all analyses. On the basis of visual observations, this condition was clearly met as Topas-wild plants exhibited clear symptoms of infection over the course of the experiment that culminated with the presence of extensive necrotic tissues at 11 dpi while the infection never extended beyond the point of inoculation in Topas-Rlm2 plants. Other studies have reported similar symptom progression at varying time points (Lowe et al., 2014; Haddadi et al., 2016), thus suggesting that experimental conditions can influence the rapidity with which L. maculans can infect its host. For this reason, it is critical to ensure that sampling procedures include enough replications that will capture an accurate biological variability within a condition. In this work we have used five biological replications to enhance both the biological and statistical power to compare gene expression across stages and conditions. We have also collected the whole infected cotyledons for analysis to achieve a thorough understanding of pathogen transcriptome activities at cotyledon level during infection. This approach was validated by our PC analysis (see Figure 3) where the clustering of samples confirmed a uniformity within a given treatment and a variability among treatments thereby supporting that observed phenotypes were indeed associated with differential gene expression. This further supported subsequent statistical analyses of DEGs.

Expression of Known Avrs and CSEPs during Biotrophic and Necrotrophic Phases

Several Avr and R genes have been identified and/or proposed to play a role in the L. maculans-canola interaction (Balesdent et al., 2002; Ghanbarnia et al., 2012). Following our analyses, all known Avr genes observed by Lowe et al. (2014) were also found to be highly expressed at 7 dpi, and reduced in their expression at 11 dpi in the compatible interaction. Based on a comparison with Topas-Rlm2, it is apparent that the period 5–7 dpi harbors the biotrophic phase of L. maculans and 11 dpi is more consistent with the necrotrophic phase. As such the differential expression of Avrs at 7 dpi indicates that some of them are clearly involved primarily in the establishment of biotrophy, and possibly the transition to necrotrophy but are no longer relevant when L. maculans has entered its necrotrophic phase.

With the availability of full genome sequences and more advanced computational tools and pipelines, we were able to identify and characterize effectors in L. maculans and better address their potential role/functionality. From strict computational predictions, there are 552 classically secreted proteins in L. maculans, which represents an unrealistic number in terms of functional effectors (Sonah et al., 2016). Recently developed tools based on machine learning were used here to further prioritize 134 CSEPs that showed evidence of upregulation in RNA-seq data. At 3 dpi, which corresponded to an asymptomatic growth phase, a common set of effectors highly expressed in both the compatible and incompatible conditions was observed, which suggests they play a minimal role in the fate of the interaction. By contrast, the expression pattern of effectors varied drastically between the compatible and incompatible conditions at 7 dpi. This approach allowed to narrow down the list of possible functional effectors to 28 that were uniquely upregulated at 7 dpi only under compatible conditions and thus presumed to be important for the establishment and maintenance of the biotrophic phase. As a matter of fact, many of the identified effectors here were either known Avr, or genes for which a role in pathogenicity was suggested (Haddadi et al., 2016). In addition, eight new effectors are proposed on the basis of their features and expression, and should be interesting candidates in future functionality assays. The studies by Haddadi et al. (2016) and Lowe et al. (2014) used different strains of L. maculans, which explains the differences in the number of differentially expressed effectors. For instance, AvrLm6 was not expressed in our study since the L. maculans strain D5 does not have AvrLm6. This suggests that it should be important to properly assess resistant germplasm with specific L. maculans strains present in a given region.

Compared to biotrophic effectors, very few effectors responsible for the necrotrophic phase are known (Lo Presti et al., 2015). In the present study, on the basis of comparative expression in compatible and incompatible interactions, we have identified 15 effectors that were distinct from Avr genes (and associated effectors) in their chronology of expression. This suggests that these effectors are specific to the necrotrophic phase of the fungus and that L. maculans, as a hemibiotroph, has indeed evolved different mechanisms to support its biotrophic and necrotrophic phase. At the same time, high expression of a single crinkler-type effector at 7 dpi in L. maculans suggests its involvement in the transition from a biotrophic to a necrotrophic phase. However, the occurrence of a single crinkler in L. maculans compared to 74 in P. sojae suggests that L. maculans does not have the crinkler-mediated mechanism leading to the necrotrophic phase found in P. sojae and other oomycetes.

A biotrophic-necrotrophic effector system that has been well studied is that of Phytophthora infestans in its interaction with potato (Whisson et al., 2007). During the biotrophic phase, P. infestans secretes AVR3a from its haustoria to suppress cell-death; as the oomycete moves to a necrotrophic stage, AVR3a is downregulated (Whisson et al., 2007). Similarly, our results showed that all 28 genes, including known Avr genes, were highly expressed during the biotrophic phase at 7 dpi and their expression reduced during the necrotrophic phase at 11 dpi. The same genes showed limited or no expression in the incompatible interaction. This suggests that these genes produce functional effectors and are key factors responsible for cross-talk between L. maculans and its host. Another gene found here, Lm5LysM, is highly expressed at 7 dpi and shows homology with the SLP1 gene of M. grisea that is expressed at the interface between the fungal cell wall and host cell plasma membrane during biotrophic invasion. The LysM-type effectors have been previously associated with biotrophy and shown to be involved in plant-fungus interactions (Gust et al., 2012; Kombrink and Thomma, 2013; Lowe et al., 2014).

Expression Pattern of NLPs Confirming the Necrotrophic Stage in L. maculans

In the present study, the highest level of expression for Lm-NLP (gene_11090) gene at 11 dpi, as L. maculans entered into its necrotrophic phase, is an observation consistent with the results of Haddadi et al. (2016). Indeed, these authors used the NLP gene expression profile as a means to distinguish genes related to the biotrophic or necrotrophic phase in L. maculans; the induction of necrosis by Lm-NLP was confirmed with a transient assay in tobacco (Haddadi et al., 2016). Our results do confirm that the highest expression of NLPs is synchronized with the necrotrophic phase, while that of known Avrs with the biotrophic phase. Avrs and NLPs can thus be considered as valid markers to characterize CSEPs as biotrophic or necrotrophic effectors. Similar findings were observed in P. infestans where INF1 and Nep1-like effectors are secreted at later stages of infection that correspond to P. infestans transitioning from a biotrophic to a necrotrophic stage (Kanneganti et al., 2006).

Expression Dynamics of CAZymes and TFs

The increasing number of upregulated CAZymes from 5 to 11 dpi in the compatible interaction can be related to the biphasic life style of L. maculans as previously suggested by Lowe et al. (2014). This is particularly relevant when compared to the incompatible interaction, where this number remained relatively unchanged over the course of infection. Our results showed that GH, CE, and AA were the most prevalent groups of CAZymes among upregulated genes. During the early stages (3 dpi) of L. maculans infection, mostly cellulose and pectin-degrading enzymes were the most prevalent groups of CAZymes among upregulated genes. On the other hand, GH, CE and AA families were predominantly upregulated during the necrotrophic stage. Pathogenic fungi face the plant cell wall as a first barrier to establish infection. Plant cell walls are mainly composed of carbohydrates and glycoproteins. To breakdown this barrier, plant pathogenic fungi need to secrete a diverse range of carbohydrate-active enzymes (CAZymes). However, our results clearly showed that the release of CAZymes in the early stages is not a key determinant of the interaction since their number and expression were similar in both interactions.

CAZymes are also involved in nutrient uptake and those prominently expressed at later stages of infection are thought to be involved in uptake of amino acid and sugars from the host. Our results bring stronger support to the concept suggested by Lowe et al. (2014) that CAZymes play an important role in L. maculans-canola interaction, notably during the establishment of necrotrophy.

Biphasic expression turnover similar to that of CAZymes was also observed with TFs. As the infection progressed in the compatible interaction, a much higher number of differentially expressed TFs was noted especially from 5 to 7 dpi. This period is critical in the fate of the interaction since it is clearly synchronized with the establishment of the infection in the compatible interaction. Of particular importance, the LmStuA TF, a member of APSES domain-containing proteins, showed its highest expression at 7 dpi in L. maculans interaction with Topas-wild. The expression profile of APSES domain-containing genes observed in this study is well aligned with an earlier report by Soyer et al. (2015). The StuA TF was found to be involved in morphogenesis, metabolites production and effector regulation (Baeza-Montañez et al., 2015; Soyer et al., 2015). Recently, StuA was suggested to play a key role in the L. maculans-canola interaction since its silencing led to a reduced expression of AvrLm1, AvrLm6, and AvrLm4-7, at 7dpi (Soyer et al., 2015). Moreover, the Zn2Cys6 TFs were upregulated in the compatible interaction supporting its role in regulating pathogenicity-related genes during the biotrophic phase. The Zn2Cys6 TFs have been reported to be involved in different regulatory functions. They are also unique to fungi and have been extensively studied in Saccharomyces cerevisiae and Aspergillus nidulans (Shimizu et al., 2003; Vienken et al., 2005).

Secondary Metabolites Involved in the Pathogenesis of L. maculans

Secondary metabolites like NRPs and phytotoxins such as sirodesmin were found to have a differential expression pattern between the compatible and incompatible interactions. Differences in gene expression pattern between the compatible and incompatible interaction were observed for NRPs, known to be involved in the production of phytotoxins, siderophores, and pigments, mostly at 11 dpi. Another important phytotoxin involved in the establishment of necrosis is sirodesmin, which is regulated by Sir genes (Gardiner et al., 2004; Haddadi et al., 2016). Most Sir genes were found to be highly expressed at 11 dpi compared to in planta L. maculans early growth stages in this study. A notable exception was SirO and SirT, which were highly expressed at 3 dpi and are thought to be involved in the production of sirodesmin early in the infection process. Indeed, Gardiner et al. (2004) have also reported significant amount of sirodesmin production at 4 dpi with L. maculans. On the other hand, the fact that SirT was equally expressed in the compatible and incompatible interactions at 3 dpi in our work would suggest that additional factors must complement their activity for infection to occur.

Based on our WEGO functional annotation, the respective vigorous and restricted growth of L. maculans on the compatible and incompatible hosts correlated very well with enrichment of genes in functional categories. Toward the later stages of the experiment, very few classes of functional categories were represented in the incompatible interaction, an observation that clearly confirms the inability of L. maculans to infect Topas-Rlm2. However, it is not possible to conduct robust statistical analyses because of the disparity in the number of genes between the compatible (572) and incompatible interactions (10). Nevertheless, these analyses give a good qualitative visualization of the genes and functions, and their relative importance in the case of compatibility. Interestingly, these analyses also highlighted in the compatible interaction how specific functional classes were indicative of the biotrophic or necrotrophic stage, and the transition from one to the other.

Validation of DEGs identified with transcriptome profiling is always a concern to build a confidence about the significance of results. The validation of DEGs by qPCR is expected for the microarray studies or transcriptomic studies with no or limited biological replications (Deshmukh et al., 2010; Xie et al., 2016). For instance, microarray studies mostly relied on qPCR to validate the results since there is a limit for the number of probes hybridized to the target in microarrays. Therefore, the highest level of gene expression is never truly represented for highly expressed genes. This also leads to a bias of data normalization and reduced linear range. By contrast, RNA-seq provides digital expression for the entire transcripts without any minimum or maximum count limit. In our study, we have used five replications, in line with the robust standards accepted for RNA-seq experiments (Fang and Cui, 2011). This is further supported by the fact that all previously identified Avr genes (AvrLm1, AvrLm4-7, AvrLm11, AvrLmJ1, AvrLm2) had an expression pattern consistent with expectations and previous reports (Lowe et al., 2014; Haddadi et al., 2016).

Conclusions

Hemibiotrophic fungi such as L. maculans have a distinctive life style that involves complex molecular processes and the expression turnover of thousands of genes. The atlas of gene expression provided here will be helpful to understand the molecular crosstalk between canola and L. maculans as it relates to compatibility or incompatibility, and could be exploited toward the deployment of novel strategies to overcome blackleg disease. The comparison made between the compatible and incompatible interaction highlighted the role of specific CSEPs, CAZymes, TFs, and secondary metabolites involved in the infection process. The differential expression pattern observed for these classes of pathogenicity-related genes can serve as a valuable resource to differentiate compatible and incompatible interactions in the context of developing resistant canola germplasm. In addition, our time series analysis based on a comprehensive and comparative differential gene expression as L. maculans infects a susceptible or resistant host has brought to light elements that define the biotrophic and necrotrophic phases of the fungus, as well as some of the mechanisms involved in the transition between the two phases.

Author Contributions

HS performed RNA extractions, library preparations, data analysis, and compilation; XZ and WF performed all bio-assays, maintained the fungal cultures, and performed RNA extractions, HB developed introgression line and provided seed material, RD and RB contributed to data analysis and overall conclusions, HS and RB wrote initial draft of the MS with final contributions from XZ, RD, HB, and WF. RB designed and directed the project.

Funding

This study is funded by the Agri-Innovation program Growing Forward 2, SaskCanola and Agriculture and Agri-Food Canada, and the Canada Research Chair in plant protection to RB.

Conflict of Interest Statement

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

Acknowledgments

The authors would like to thank the Agri-Innovation Program Growing Forward 2, SaskCanola, Agriculture and Agri-Food Canada and the Canada Research Chair in Plant Protection to RRB for financial support. The authors also thank Brian Boyle from Institut de Biologie integrative et des Systemes at Laval University for his valuable suggestions regarding RNA-seq library preparations, Caroline Labbé for technical assistance in many aspects of this work, and Paula Parks for help with the greenhouse assays.

Supplementary Material

The Supplementary Material for this article can be found online at: http://journal.frontiersin.org/article/10.3389/fpls.2016.01784/full#supplementary-material

References

Baeza-Montañez, L., Gold, S. E., Espeso, E. A., and García-Pedrajas, M. D. (2015). Conserved and distinct functions of the “stunted” (StuA)-homolog Ust1 during cell differentiation in the corn smut fungus Ustilago maydis. Mol. Plant Microbe Interact. 28, 86–102. doi: 10.1094/MPMI-07-14-0215-R

PubMed Abstract | CrossRef Full Text | Google Scholar

Balesdent, M. H., Attard, A., Kühn, M. L., and Rouxel, T. (2002). New avirulence genes in the phytopathogenic fungus Leptosphaeria maculans. Phytopathology 92, 1122–1133. doi: 10.1094/PHYTO.2002.92.10.1122

PubMed Abstract | CrossRef Full Text | Google Scholar

Balesdent, M. H., Fudal, I., Ollivier, B., Bally, P., Grandaubert, J., Eber, F., et al. (2013). The dispensable chromosome of Leptosphaeria maculans shelters an effector gene conferring avirulence towards Brassica rapa. New Phytol. 198, 887–898. doi: 10.1111/nph.12178

PubMed Abstract | CrossRef Full Text | Google Scholar

Bolger, A. M., Lohse, M., and Usadel, B. (2014). Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics 30, 2114–2120. doi: 10.1093/bioinformatics/btu170

PubMed Abstract | CrossRef Full Text

Cortázar, A. R., Aransay, A. M., Alfaro, M., Oguiza, J. A., and Lavín, J. L. (2014). SECRETOOL: integrated secretome analysis tool for fungi. Amino Acids 46, 471–473. doi: 10.1007/s00726-013-1649-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Delourme, R., Chevre, A., Brun, H., Rouxel, T., Balesdent, M., Dias, J., et al. (2006). Major gene and polygenic resistance to Leptosphaeria maculans in oilseed rape (Brassica napus). Eur. J. Plant Pathol. 114, 41–52. doi: 10.1007/s10658-005-2108-9

CrossRef Full Text | Google Scholar

Deshmukh, R., Singh, A., Jain, N., Anand, S., Gacche, R., Singh, A., et al. (2010). Identification of candidate genes for grain number in rice (Oryza sativa L.). Funct. Integr. Genomics 10, 339–347. doi: 10.1007/s10142-010-0167-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Fang, Z., and Cui, X. (2011). Design and validation issues in RNA-seq experiments. Brief. Bioinform. 12, 280–287. doi: 10.1093/bib/bbr004

PubMed Abstract | CrossRef Full Text | Google Scholar

Fitt, B. D. L., Brun, H., Barbetti, M. J., and Rimmer, S. R. (2006). “World-wide importance of phoma stem canker (Leptosphaeria maculans and L. biglobosa) on oilseed rape (Brassica napus),” in Sustainable Strategies for Managing Brassica napus (Oilseed Rape) Resistance to Leptosphaeria maculans (Phoma Stem Canker), eds B. D. L. Fitt, N. Evans, B. J. Howlett, and B. M. Cooke (Springer), 3–15.

Google Scholar

Fudal, I., Ross, S., Gout, L., Blaise, F., Kuhn, M. L., Eckert, M. R., et al. (2007). Heterochromatin-like regions as ecological niches for avirulence genes in the Leptosphaeria maculans genome: map-based cloning of AvrLm6. Mol. Plant Microbe Interact. 20, 459–470. doi: 10.1094/MPMI-20-4-0459

PubMed Abstract | CrossRef Full Text | Google Scholar

Garber, M., Grabherr, M. G., Guttman, M., and Trapnell, C. (2011). Computational methods for transcriptome annotation and quantification using RNA-seq. Nat. Methods 8, 469–477. doi: 10.1038/nmeth.1613

PubMed Abstract | CrossRef Full Text | Google Scholar

Gardiner, D. M., Cozijnsen, A. J., Wilson, L. M., Pedras, M. S. C., and Howlett, B. J. (2004). The sirodesmin biosynthetic gene cluster of the plant pathogenic fungus Leptosphaeria maculans. Mol. Microbiol. 53, 1307–1318. doi: 10.1111/j.1365-2958.2004.04215.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Ghanbarnia, K., Fudal, I., Larkan, N. J., Links, M. G., Balesdent, M. H., Profotova, B., et al. (2015). Rapid identification of the Leptosphaeria maculans avirulence gene AvrLm2 using an intraspecific comparative genomics approach. Mol. Plant Pathol. 16, 699–709. doi: 10.1111/mpp.12228

PubMed Abstract | CrossRef Full Text | Google Scholar

Ghanbarnia, K., Lydiate, D. J., Rimmer, S. R., Li, G., Kutcher, H. R., Larkan, N. J., et al. (2012). Genetic mapping of the Leptosphaeria maculans avirulence gene corresponding to the LepR1 resistance gene of Brassica napus. Theor. Appl. Genet. 124, 505–513. doi: 10.1007/s00122-011-1724-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Gout, L., Fudal, I., Kuhn, M. L., Blaise, F., Eckert, M., Cattolico, L., et al. (2006). Lost in the middle of nowhere: the AvrLm1 avirulence gene of the Dothideomycete Leptosphaeria maculans. Mol. Microbiol. 60, 67–80. doi: 10.1111/j.1365-2958.2006.05076.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Grigoriev, I. V., Nikitin, R., Haridas, S., Kuo, A., Ohm, R., Otillar, R., et al. (2014). MycoCosm portal: gearing up for 1000 fungal genomes. Nucleic Acids Res. 42, D699–D704. doi: 10.1093/nar/gkt1183

PubMed Abstract | CrossRef Full Text | Google Scholar

Guo, M., Chen, Y., Du, Y., Dong, Y., Guo, W., Zhai, S., et al. (2011). The bZIP transcription factor MoAP1 mediates the oxidative stress response and is critical for pathogenicity of the rice blast fungus Magnaporthe oryzae. PLoS Pathog. 7:e1001302. doi: 10.1371/journal.ppat.1001302

PubMed Abstract | CrossRef Full Text | Google Scholar

Gust, A. A., Willmann, R., Desaki, Y., Grabherr, H. M., and Nürnberger, T. (2012). Plant LysM proteins: modules mediating symbiosis and immunity. Trends Plant Sci. 17, 495–502. doi: 10.1016/j.tplants.2012.04.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Haas, B. J., Kamoun, S., Zody, M. C., Jiang, R. H., Handsaker, R. E., Cano, L. M., et al. (2009). Genome sequence and analysis of the Irish potato famine pathogen Phytophthora infestans. Nature 461, 393–398. doi: 10.1038/nature08358

PubMed Abstract | CrossRef Full Text | Google Scholar

Haddadi, P., Ma, L., Wang, H., and Borhan, M. H. (2016). Genome-wide transcriptomic analyses provide insights into the lifestyle transition and effector repertoire of Leptosphaeria maculans during the colonization of Brassica napus seedlings. Mol. Plant Pathol. 17, 1196–1210. doi: 10.1111/mpp.12356

PubMed Abstract | CrossRef Full Text | Google Scholar

Howlett, B. J., Idnurm, A., and Pedras, M. S. C. (2001). Leptosphaeria maculans, the causal agent of blackleg disease of Brassicas. Fungal Genet. Biol. 33, 1–14. doi: 10.1006/fgbi.2001.1274

PubMed Abstract | CrossRef Full Text | Google Scholar

Huang, Y.-J., Qi, A., King, G. J., and Fitt, B. D. (2014). Assessing quantitative resistance against Leptosphaeria maculans (phoma stem canker) in Brassica napus (oilseed rape) in young plants. PLoS ONE 9:e84924. doi: 10.1371/journal.pone.0084924

PubMed Abstract | CrossRef Full Text | Google Scholar

Huibers, R. P., De Jong, M., Dekter, R. W., and Van Den Ackerveken, G. (2009). Disease-specific expression of host genes during downy mildew infection of Arabidopsis. Mol. Plant Microbe Interact. 22, 1104–1115. doi: 10.1094/MPMI-22-9-1104

PubMed Abstract | CrossRef Full Text | Google Scholar

Kanneganti, T. D., Huitema, E., Cakir, C., and Kamoun, S. (2006). Synergistic interactions of the plant cell death pathways induced by Phytophthora infestans Nep1-like protein PiNPP1. 1 and INF1 elicitin. Mol. Plant Microbe Interact. 19, 854–863. doi: 10.1094/MPMI-19-0854

PubMed Abstract | CrossRef Full Text | Google Scholar

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

PubMed Abstract | CrossRef Full Text | Google Scholar

Kombrink, A., and Thomma, B. P. (2013). LysM effectors: secreted proteins supporting fungal life. PLoS Pathog. 9:e1003769. doi: 10.1371/journal.ppat.1003769

PubMed Abstract | CrossRef Full Text | Google Scholar

Kutcher, H., Brandt, S., Smith, E., Ulrich, D., Malhi, S., and Johnston, A. (2013). Blackleg disease of canola mitigated by resistant cultivars and four-year crop rotations in western Canada. Can. J. Plant Pathol. 35, 209–221. doi: 10.1080/07060661.2013.775600

CrossRef Full Text | Google Scholar

Larkan, N., Yu, F., Lydiate, D., Rimmer, S. R., and Borhan, M. H. (2016a). Single R gene introgression lines for accurate dissection of the Brassica - Leptosphaeria pathosystem. Front. Plant Sci. 7:1771. doi: 10.3389/fpls.2016.01771

CrossRef Full Text

Larkan, N. J., Lydiate, D. J., Parkin, I. A. P., Nelson, M. N., Epp, D. J., Cowling, W. A., et al. (2013). The Brassica napus blackleg resistance gene LepR3 encodes a receptor-like protein triggered by the Leptosphaeria maculans effector AVRLM1. New Phytol. 197, 595–605. doi: 10.1111/nph.12043

PubMed Abstract | CrossRef Full Text | Google Scholar

Larkan, N. J., Ma, L., and Borhan, M. H. (2015). The Brassica napus receptor-like protein Rlm2 is encoded by a second allele of the LepR3/Rlm2 blackleg resistance locus. Plant Biotechnol. J. 13, 983–992. doi: 10.1111/pbi.12341

PubMed Abstract | CrossRef Full Text | Google Scholar

Larkan, N. J., Raman, H., Lydiate, D. J., Robinson, S. J., Yu, F., Barbulescu, D. M., et al. (2016b). Multi-environment QTL studies suggest a role for cysteine-rich protein kinase genes in quantitative resistance to blackleg disease in Brassica napus. BMC Plant Biol. 16:183. doi: 10.1186/s12870-016-0877-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, X., Wu, J., Yin, L., Zhang, Y., Qu, J., and Lu, J. (2015). Comparative transcriptome analysis reveals defense-related genes and pathways against downy mildew in Vitis amurensis grapevine. Plant Physiol. Biochem. 95, 1–14. doi: 10.1016/j.plaphy.2015.06.016

PubMed Abstract | CrossRef Full Text | Google Scholar

Lo Presti, L., Lanver, D., Schweizer, G., Tanaka, S., Liang, L., Tollot, M., et al. (2015). Fungal effectors and plant susceptibility. Annu. Rev. Plant Biol. 66, 513–545. doi: 10.1146/annurev-arplant-043014-114623

PubMed Abstract | CrossRef Full Text | Google Scholar

Lombard, V., Golaconda Ramulu, H., Drula, E., Coutinho, P. M., and Henrissat, B. (2014). The carbohydrate-active enzymes database (CAZy) in 2013. Nucleic Acids Res. 42, D490–D495. doi: 10.1093/nar/gkt1178

PubMed Abstract | CrossRef Full Text | Google Scholar

Long, Y., Wang, Z., Sun, Z., Fernando, D. W., McVetty, P. B., and Li, G. (2011). Identification of two blackleg resistance genes and fine mapping of one of these two genes in a Brassica napus canola cultivar 'surpass 400'. Theor. Appl. Genet. 122, 1223–1231. doi: 10.1007/s00122-010-1526-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Love, M., Anders, S., and Huber, W. (2014). Differential analysis of count data–the DESeq2 package. Genome Biol. 15, p550. doi: 10.1186/s13059-014-0550-8

CrossRef Full Text

Lowe, R. G., Cassin, A., Grandaubert, J., Clark, B. L., Van De Wouw, A. P., Rouxel, T., et al. (2014). Genomes and transcriptomes of partners in plant-fungal-interactions between canola (Brassica napus) and two Leptosphaeria species. PLoS ONE 9:e103098. doi: 10.1371/journal.pone.0103098

PubMed Abstract | CrossRef Full Text | Google Scholar

Malinovsky, F. G., Batoux, M., Schwessinger, B., Youn, J. H., Stransfeld, L., Win, J., et al. (2014). Antagonistic regulation of growth and immunity by the Arabidopsis basic helix-loop-helix transcription factor homolog of brassinosteroid enhanced expression2 interacting with increased leaf inclination1 binding bHLH1. Plant Physiol. 164, 1443–1455. doi: 10.1104/pp.113.234625

PubMed Abstract | CrossRef Full Text | Google Scholar

O'Connell, R. J., Thon, M. R., Hacquard, S., Amyotte, S. G., Kleemann, J., Torres, M. F., et al. (2012). Lifestyle transitions in plant pathogenic Colletotrichum fungi deciphered by genome and transcriptome analyses. Nat. Genet. 44, 1060–1065. doi: 10.1038/ng.2372

PubMed Abstract | CrossRef Full Text | Google Scholar

Park, J., Park, J., Jang, S., Kim, S., Kong, S., Choi, J., et al. (2008). FTFD: an informatics pipeline supporting phylogenomic analysis of fungal transcription factors. Bioinformatics 24, 1024–1025. doi: 10.1093/bioinformatics/btn058

PubMed Abstract | CrossRef Full Text | Google Scholar

Parlange, F., Daverdin, G., Fudal, I., Kuhn, M. L., Balesdent, M. H., Blaise, F., et al. (2009). Leptosphaeria maculans avirulence gene AvrLm4-7 confers a dual recognition specificity by the Rlm4 and Rlm7 resistance genes of oilseed rape, and circumvents Rlm4-mediated recognition through a single amino acid change. Mol. Microbiol. 71, 851–863. doi: 10.1111/j.1365-2958.2008.06547.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Pedras, M. S. C., Taylor, J. L., and Nakashima, T. T. (1993). A novel chemical signal from the “blackleg” fungus: beyond phytotoxins and phytoalexins. J. Org. Chem. 58, 4778–4780. doi: 10.1021/jo00070a002

CrossRef Full Text | Google Scholar

Petersen, T. N., Brunak, S., Von Heijne, G., and Nielsen, H. (2011). SignalP 4.0: discriminating signal peptides from transmembrane regions. Nat. Methods 8, 785–786. doi: 10.1038/nmeth.1701

PubMed Abstract | CrossRef Full Text | Google Scholar

Plissonneau, C., Daverdin, G., Ollivier, B., Blaise, F., Degrave, A., Fudal, I., et al. (2016). A game of hide and seek between avirulence genes AvrLm4-7 and AvrLm3 in Leptosphaeria maculans. New Phytol. 209, 1613–1624. doi: 10.1111/nph.13736

PubMed Abstract | CrossRef Full Text | Google Scholar

Raman, R., Taylor, B., Marcroft, S., Stiller, J., Eckermann, P., Coombes, N., et al. (2012). Molecular mapping of qualitative and quantitative loci for resistance to Leptosphaeria maculans causing blackleg disease in canola (Brassica napus L.). Theor. Appl. Genet. 125, 405–418. doi: 10.1007/s00122-012-1842-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Rawlings, N. D., Waller, M., Barrett, A. J., and Bateman, A. (2014). MEROPS: the database of proteolytic enzymes, their substrates and inhibitors. Nucleic Acids Res. 42, D503–D509. doi: 10.1093/nar/gkt953

PubMed Abstract | CrossRef Full Text | Google Scholar

Robinson, M. D., McCarthy, D. J., and Smyth, G. K. (2010). edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics 26, 139–140. doi: 10.1093/bioinformatics/btp616

PubMed Abstract | CrossRef Full Text | Google Scholar

Rouxel, T., Grandaubert, J., Hane, J. K., Hoede, C., Van De Wouw, A. P., Couloux, A., et al. (2011). Effector diversification within compartments of the Leptosphaeria maculans genome affected by repeat-induced point mutations. Nat. Commun. 2, 202. doi: 10.1038/ncomms1189

PubMed Abstract | CrossRef Full Text | Google Scholar

Sestili, S., Polverari, A., Luongo, L., Ferrarini, A., Scotton, M., Hussain, J., et al. (2011). Distinct colonization patterns and cDNA-AFLP transcriptome profiles in compatible and incompatible interactions between melon and different races of Fusarium oxysporum f. sp. melonis. BMC Genomics 12:122. doi: 10.1186/1471-2164-12-122

PubMed Abstract | CrossRef Full Text | Google Scholar

Seyednasrollah, F., Laiho, A., and Elo, L. L. (2015). Comparison of software packages for detecting differential expression in RNA-seq studies. Brief. Bioinform. 16, 59–70. doi: 10.1093/bib/bbt086

PubMed Abstract | CrossRef Full Text | Google Scholar

Shimizu, K., Hicks, J. K., Huang, T. P., and Keller, N. P. (2003). Pka, Ras and RGS protein interactions regulate activity of AflR, a Zn (II) 2Cys6 transcription factor in Aspergillus nidulans. Genetics 165, 1095–1104.

PubMed Abstract | Google Scholar

Sonah, H., Deshmukh, R. K., and Bélanger, R. R. (2016). Computational prediction of effector proteins in fungi: opportunities and challenges. Front. Plant Sci. 7:126. doi: 10.3389/fpls.2016.00126

PubMed Abstract | CrossRef Full Text | Google Scholar

Soyer, J. L., El Ghalid, M., Glaser, N., Ollivier, B., Linglin, J., Grandaubert, J., et al. (2014). Epigenetic control of effector gene expression in the plant pathogenic fungus Leptosphaeria maculans. PLoS Genet. 10:e1004227. doi: 10.1371/journal.pgen.1004227

PubMed Abstract | CrossRef Full Text | Google Scholar

Soyer, J. L., Hamiot, A., Ollivier, B., Balesdent, M. H., Rouxel, T., and Fudal, I. (2015). The APSES transcription factor LmStuA is required for sporulation, pathogenic development and effector gene expression in Leptosphaeria maculans. Mol. Plant Pathol. 16, 1000–1005. doi: 10.1111/mpp.12249

PubMed Abstract | CrossRef Full Text | Google Scholar

Sperschneider, J., Gardiner, D. M., Dodds, P. N., Tini, F., Covarelli, L., Singh, K. B., et al. (2015). EffectorP: predicting fungal effector proteins from secretomes using machine learning. New Phytol. 210, 743–761. doi: 10.1111/nph.13794

PubMed Abstract | CrossRef Full Text | Google Scholar

Stukenbrock, E. H., and McDonald, B. A. (2009). Population genetics of fungal and oomycete effectors involved in gene-for-gene interactions. Mol. Plant Microbe Interact. 22, 371–380. doi: 10.1094/MPMI-22-4-0371

PubMed Abstract | CrossRef Full Text | Google Scholar

Trapnell, C., Pachter, L., and Salzberg, S. L. (2009). TopHat: discovering splice junctions with RNA-Seq. Bioinformatics 25, 1105–1111. doi: 10.1093/bioinformatics/btp120

PubMed Abstract | CrossRef Full Text | Google Scholar

Trapnell, C., Roberts, A., Goff, L., Pertea, G., Kim, D., Kelley, D. R., et al. (2012). Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks. Nat. Protoc. 7, 562–578. doi: 10.1038/nprot.2012.016

PubMed Abstract | CrossRef Full Text | Google Scholar

Van de Wouw, A. P., Lowe, R. G. T., Elliott, C. E., Dubois, D. J., and Howlett, B. J. (2014). An avirulence gene, AvrLmJ1, from the blackleg fungus, Leptosphaeria maculans, confers avirulence to Brassica juncea varieties. Mol. Plant Pathol. 15, 523–530. doi: 10.1111/mpp.12105

PubMed Abstract | CrossRef Full Text | Google Scholar

Van de Wouw, A. P., Marcroft, S. J., and Howlett, B. J. (2016). Blackleg disease of canola in Australia. Crop Pasture Sci. 67, 273–283. doi: 10.1071/CP15221

CrossRef Full Text | Google Scholar

Vienken, K., Scherer, M., and Fischer, R. (2005). The Zn(II)2Cys6 putative Aspergillus nidulans transcription factor repressor of sexual development inhibits sexual development under low-carbon conditions and in submersed culture. Genetics 169, 619–630 doi: 10.1534/genetics.104.030767

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, X., Liu, W., Chen, X., Tang, C., Dong, Y., Ma, J., et al. (2010). Differential gene expression in incompatible interaction between wheat and stripe rust fungus revealed by cDNA-AFLP and comparison to compatible interaction. BMC Plant Biol. 10:9. doi: 10.1186/1471-2229-10-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, Z., Gerstein, M., and Snyder, M. (2009). RNA-Seq: a revolutionary tool for transcriptomics. Nat. Rev. Genet. 10, 57–63. doi: 10.1038/nrg2484

PubMed Abstract | CrossRef Full Text | Google Scholar

Whisson, S. C., Boevink, P. C., Moleleki, L., Avrova, A. O., Morales, J. G., Gilroy, E. M., et al. (2007). A translocation signal for delivery of oomycete effector proteins into host plant cells. Nature 450, 115–118. doi: 10.1038/nature06203

PubMed Abstract | CrossRef Full Text | Google Scholar

Williams, A. G., Thomas, S., Wyman, S. K., and Holloway, A. K. (2014). RNA-seq Data: challenges in and recommendations for experimental design and analysis. Curr. Protoc. Hum. Genet. 83, 11.13.1–11.13.20. doi: 10.1002/0471142905.hg1113s83

PubMed Abstract | CrossRef Full Text | Google Scholar

Xie, J., Li, S., Mo, C., Xiao, X., Peng, D., Wang, G., et al. (2016). Genome and transcriptome sequences reveal the specific parasitism of the nematophagous Purpureocillium lilacinum 36-1. Front. Microbiol. 7:1084. doi: 10.3389/fmicb.2016.01084

PubMed Abstract | CrossRef Full Text | Google Scholar

Yin, Y., Mao, X., Yang, J., Chen, X., Mao, F., and Xu, Y. (2012). dbCAN: a web resource for automated carbohydrate-active enzyme annotation. Nucleic Acids Res. 40, W445–W451. doi: 10.1093/nar/gks479

PubMed Abstract | CrossRef Full Text | Google Scholar