- Research article
- Open Access
From yeast to hypha: defining transcriptomic signatures of the morphological switch in the dimorphic fungal pathogen Ophiostoma novo-ulmi
BMC Genomicsvolume 17, Article number: 920 (2016)
Yeast-to-hypha transition is a major morphological change in fungi. Molecular regulators and pathways that are involved in this process have been extensively studied in model species, including Saccharomyces cerevisiae. The Mitogen-Actived Protein Kinase (MAPK) cascade, for example, is known to be involved in the yeast-to-pseudohypha switch. Yet the conservation of mechanisms regulating such morphological changes in non-model fungi is still poorly understood. Here, we investigate cell remodeling and transcriptomic modifications that occur during this morphological switch in the highly aggressive ascomycete fungus Ophiostoma novo-ulmi, the causal agent of Dutch elm disease.
Using a combination of light microscopy, scanning electron microscopy and flow cytometry, we demonstrate that the morphological switch occurs in less than 27 h, with phenotypic cell modifications being detected within the first 4 h. Using RNAseq, we found that over 22% of the genome of O. novo-ulmi is differentially expressed during the transition. By performing clustering analyses of time series gene expression data, we identified several sets of genes that are differentially expressed according to distinct and representative temporal profiles. Further, we found that several genes that are homologous to S. cerevisiae MAPK genes are regulated during the yeast-to-hypha transition in O. novo-ulmi and mostly over-expressed, suggesting convergence in gene expression regulation.
Our results are the first report of a time-course experiment monitoring the morphological transition in a non-model Sordariomycota species and reveal many genes of interest for further functional investigations of fungal dimorphism.
Ophiostoma novo-ulmi (Ascomycota, Sordariomycetes) is the highly aggressive dimorphic pathogen that is responsible for the ongoing pandemic of Dutch elm disease (DED) . This fungus is capable of assuming two distinct forms: a unicellular yeast stage and a multicellular hyphal form. In the context of DED, both forms are thought to be involved in pathogenesis. The yeast cells are considered to be responsible for the passive colonization of the xylem vessels as they are carried by the sap flow and found rapidly after invasion both in the shoots and in the roots . The hyphal form (also known as the invasive phase) is suggested as being responsible for fungal dispersion between adjacent vessels by passing through the vessel pits. This pluricellular form is also very important during the saprophytic stage of O. novo-ulmi for the colonization of the principal egg galleries that are dug by female elm bark beetles on dead trees and further into secondary galleries that are excavated by the beetle larvae .
Many species within the genus Ophiostoma are known to be dimorphic, including the other two DED fungi, O. ulmi and O. himal-ulmi as well as the sapwood-staining fungi O. piceae [4, 5] and O. quercus . Few studies have been devoted to the identification of key factors that regulate the yeast-to-hypha (Y-to-H) transition within this genus. Nutritional factors, such as nitrogen sources [7, 8], pyridoxine  or linoleic acid , have been shown to be involved, together with other molecular factors, such as Ca2+-calmodulin interaction [11, 12] or cyclic Adenosine MonoPhosphate (cAMP) . Inoculum size is also implicated in the regulation of the Y-to-H transition in DED fungi and effects have been confirmed for O. ulmi, both subspecies of O. novo-ulmi (americana and novo-ulmi), and O. himal-ulmi [8, 14–16]. Moreover, quorum sensing has been shown in O. ulmi, O. piceae and O. floccosum [15, 17, 18]. Until now, only two genes have been identified as being involved in the dimorphism in O. novo-ulmi, namely COL1  and an as yet uncharacterized gene that also affects asexual sporulation (synnematospores) and pathogenicity [20, 21]. The knockout col1 mutant exhibits reduced mycelial growth, whereas yeast growth is barely affected by the mutation . However, pathways and genes regulating the Y-to-H transition in O. novo-ulmi are still largely unknown.
Dimorphism in fungi is a morphological characteristic that has been studied for many years in model species, such as Saccharomyces cerevisiae, the human pathogens Candida albicans and Histoplasma capsulatum, and the plant pathogen Ustilago maydis. The hyphal state is very variable among species, and is represented by pseudo-hyphae in S. cerevisiae (not common in nature) or by true septate hyphae in filamentous fungi. In all of these species, several conserved pathways have been reported to be linked to the Y-to-H switch, including the Mitogen-Activated Protein Kinase (MAPK) cascade [22–24], the Protein Kinase A (PKA) pathway [25–27], and the pH-dependent RIM pathway [28, 29].
In S. cerevisiae, the MAPK cascade that is involved in morphological changes is well described (for a simplified schematic, see Fig. 1). The initial signal is perceived by the osmoreceptors Sho1 and Msb2, which are located on the plasma membrane. The signal is then transduced to a first kinase Ste20, and a cascade of phosphorylations is activated through Ste11 (MAPKKK), Ste7 (MAPKK), and Kss1 (MAPK). The unphosphorylated Kss1 acts as a repressor of the morphological switch since it sequestrates the transcription factors Tec1 and Ste12 in the nucleus. The phosphorylation of Kss1 by Ste7 induces the release of Tec1 and Ste12 [23, 30, 31]. An alternative pathway involves the activation by Ste11 of another MAPKK, Pbs2, and then by Hog1, a MAPK that responds to hyperosmolar conditions [32, 33]. Hog1 activates Msn4, which is a stress-response transcriptional factor that is linked to Ste12 . In C. albicans, both cascades (via Ste7 or Hog1) are conserved and involved in the control of morphogenesis [34, 35].
The PKA protein of S. cerevisiae is a tetramer composed of a homodimer of regulator sub-units (Bcy1) and a homo- or heterodimer of catalytic sub-units (two of the three possible Tpk1, Tpk2 or Tpk3 sub-units). The activation of the PKA pathway depends upon the presence of cAMP  (for a simplified schematic, see Fig. 1). Of the three Tpk proteins, only Tpk2 exerts a positive effect on pseudohyphal growth in S. cerevisiae, since deletion of Tpk2 represses filamentous growth. In contrast, Tpk1 deletion has no effect, while knockout mutants of Tpk3 exhibit hyperfilamentous growth [36, 37]. Both Tpk1 and Tpk2 regulate genes that are involved in the formation of pseudohyphae. Tpk1 regulates the activity of the dual-specificity tyrosine-regulated kinase Yak1, which controls the expression of Flo11. Flo11 is a gene encoding a cell surface flocculin through the transcription factors Sok1 and Phd1 . Through phosphorylation, Tpk2 activates the transcriptional factor Flo8 which, in turn, activates the expression of filamentation target genes, such as Flo11 . Also through Tpk2, the PKA pathway induces the activation of the transcriptional factor Phd1. PKA is highly conserved in C. albicans and the Phd1 S. cerevisiae homolog, Efg1, is also activated through Tpk2 . Both Phd1 and Efg1 are important regulators of dimorphism [22, 39].
Finally, the RIM pathway, which is also called the PAL cascade in filamentous fungi, is dependent on variation in pH within the environment. In C. albicans and S. cerevisiae a change from acidic to neutral pH induces the Y-to-H switch (Y-to-pseudohyphae in S. cerevisiae) [40, 41], whereas hyphal form in U. maydis is favoured by acidic pH . The PAL/RIM pathway involves at least six proteins which regulate the activation of a zinc-finger-like transcriptional factor called PacC or RIM101 . RIM101 regulates, in turn, the expression of downstream genes, such as PHR1 and PHR2 in C. albicans (GAS1 in S. cerevisiae) .
A few large-scale transcriptomic analyses have been conducted to determine the molecular regulation of the Y-to-H transition in model species. RNAseq was recently used in Penicillium marneffei where 2718 genes (28% of the gene content) were differentially expressed during the morphological change . DNA microarrays in Candida albicans were used under different conditions to identify sets of genes that are regulated during the transition [43–45]. Yet genome-wide analyses of this major morphological change, which would describe the global cell response to Y-to-H transition-promoting stimuli, are still lacking. Indeed, most of the studies focus on a particular protein or pathway.
The O. novo-ulmi genome was recently sequenced  and almost 75% of the 8640 predicted genes have since been annotated , thereby facilitating further genomic studies. Comparative analyses of the transcriptomes of yeast and mycelial forms of O. novo-ulmi strain H327  show a clear difference in gene expression profiles between the two life stages and reveal processes that are specific to each form. Based on this study, we investigate the regulation of gene expression during the Y-to-H transition using a time-series approach, which allows the description of dynamic biological processes. We used RNAseq technology to characterize major molecular changes that are associated with the morphological switch. We also focused attention on the regulation of gene expression in O. novo-ulmi MAPK, PKA and RIM pathway homologs.
Microscopy and flow cytometry
We first defined the Y-to-H transition in O. novo-ulmi by recording the modifications that were observed at the population and cell levels using light microscopy. During the Y-to-H transition that was induced by incubation in still-liquid culture in OCM (complemented with proline), the percentage of yeast cells decreases within the first 10 h following transfer to the induction medium, while the percentage of mycelium increases (Fig. 2a). Further, cell size increases with time (Fig. 2b).
The Y-to-H transition was then followed using flow cytometry. We compared conditions that do not favor transition (shake-liquid OMM complemented with proline [7, 8], Fig. 2c) with the condition that was selected in this study (Fig. 2c). We observe that both cell size (FSC) and cell internal granularity (SSC) vary through time under conditions that induce the Y-to-H transition. Within the first 10 h, peaks of density shift towards higher FSC and SSC. Finally, scanning electron micrographs (SEM) that were taken at each time point confirm the Y-to-H transition and highlight cell shape modification within the first 4 h, as cells may switch from an ovoid (OY, Fig. 2d) to a spherical shape (SY).
Global view of transcriptomes of yeast-to-hypha transition
For each of the six time points that were previously defined by microscopy and flow cytometry (i.e., 0, 2, 4, 6, 10 and 27 h after transfer to the induction medium), samples were collected and total RNA was extracted to sequence the whole transcriptome. After filtering and trimming, RNAseq samples contained between 7 and 12 million reads (Table 1). By alignment on the complete set of exons of the O. novo-ulmi genome, we are able to achieve from 25.4 up to 63.6 X exon coverage depth.
We retained a total of 7605 genes expressed under our experimental conditions out of the 8640 genes that are predicted in O. novo-ulmi for downstream analyses. Each sample contains at least 7256 genes that are represented by more than one read (Table 1, Additional file 1: Figure S2, Additional file 2: Table S3). Multidimensional scaling (MDS) analysis shows that the samples cluster according to the time after transfer to OCM and, within-time points, the three biological replicates cluster together (Fig. 3). The length of incubation of cells in non-agitated OCM (0–27 h) explains most of the variability in the data (first axis). Molecular variability is small between 10 and 27 h. In the second dimension, there is a large difference between samples that were collected at the beginning of the experiment and those that were collected as early as 2 h after transfer to non-agitated OCM.
Identification and analysis of differentially Expressed Genes (DEGs)
The Next-maSigPro package in R (v3.0.1) applies generalized linear models (GLM) to regression models in order to determine differentially expressed genes (DEGs) in RNAseq time-course analyses [49, 50]. Here, when the package was used on normalized read counts, a total of 1897 genes out of the 7605 expressed genes were identified as DEGs during the Y-to-H transition (0–27 h), when a false rate discovery (FDR) threshold  of 0.05 was employed. (Additional file 2: Table S3). This set of DEGs was then analyzed by clustering methods with STEM software [52, 53]. Here, we present results for the assignment of DEGs to 50 predefined gene expression model profiles (0–49) (Fig. 4). A total of five profiles are identified as containing more genes than expected by chance (significant P-values that were corrected for FDR, ranging from 72 to 815 genes, colored profiles). Some genes are present in more than one profile (Additional file 2: Table S4). The five profiles were then consolidated into four clusters that gather together similar model profiles, which were represented by different colors, with one cluster including two profiles (red, Fig. 4). We performed a more in depth investigation of profiles 39 (815 genes, 9.4% of the total gene content) and 8 (256 genes, 3% of the total gene content), since they are the most significant (P-values of 0 and 1e−144, respectively) and exact opposites in terms of gene expression regulation (Figs. 4 and 5a, b).
Genes with increasing expression during Y-to-H transition (profile 39)
Profile 39 comprises genes that are increasingly expressed during the Y-to-H transition. The most notable gene in this profile is Onu3773, which is annotated as coding for an extracellular serine-threonine-rich protein, and also related to the Adhesin protein Mad1 (Additional file 2: Table S4). This gene is over-expressed by 7819 times at 27 h compared to 0 h. Another notable gene is Onu6790.1, which is annotated as coding for the Woronin body major protein, and is homologous to Hex1 in O. floccosum. In O. novo-ulmi, this gene is barely expressed at the beginning of the switch, but its expression increases rapidly within the 27 first hours (959 times more strongly expressed at 27 h than at 0 h; Additional file 2: Table S4). Genes coding for actin, tubulin and septin are also present in profile 39. Of the 67 genes encoding carbohydrate-active enzymes (CAZymes) that were found within the total set of DEGs, 31 genes are grouped in this profile. This set includes glycosyltransferases (GT, n = 21), glycoside hydroxylases (GH, n = 8) and carbohydrate binding modules (CBM, n = 2) (Additional file 2: Table S4). Finally, 90 genes found in profile 39 have been already described as being over-expressed in the mycelium growth phase of O. novo-ulmi  (Additional file 2: Table S3). By gene ontology (GO) term-enrichment analysis, we determined that 112 biological processes (BPs) containing more than four genes are significantly enriched (P-value ≤ 0.01) in profile 39 (Additional file 2: Table S5). After using REVIGO software  to reduce redundancy, 28 terms were retained (Additional file 3: Figure S3). Among the enriched BPs, many are related to cellular organization, including localization (GO:0051179, 127 genes), cellular component organization (GO:0016043, 67 genes) and microtubule-based movement (GO:0007018, 10 genes) (Additional file 2: Table S5). Cell division, vesicle transport and protein modifications are also dominant processes.
Genes with decreasing expression during Y-to-H transition (profile 8)
Profile 8 comprises genes for which expression decreases during the Y-to-H transition (Additional file 2: Table S6). Only four genes encoding CAZymes (three GHs and one GT) are found. Gene Onu5171 is one of the top genes in the profile (40 times more repressed at 27 h compared to 0 h), and encoding a siderophore iron transporter 1 (Ferrioxamine B permease). A total of 15 genes that were associated with yeast growth conditions by Nigg et al.  are down-regulated during the Y-to-H transition and are present in profile 8 (Additional file 2: Table S3). In this profile, we further found gene Onu8210, which was annotated as encoding a DNA-binding protein. This protein has a zinc-finger DNA binding domain, which is homologous to the domain of S. cerevisiae Mig1, Nrg1 and Msn4 proteins (data not shown). Gene Onu8210 is down-regulated during the Y-to-H transition (5.75 times less expressed at 27 h after transfer to OCM).
Finally, GO term enrichment analysis that was followed by the application of REVIGO highlighted 22 BPs and revealed over-representation of terms that were associated with non-coding RNA (ncRNA) metabolism and ribosome biogenesis (Additional file 4: Figure S4 and Additional file 2: Table S7).
Other genes of interest
Gene Onu6217-CAT1 homolog, which has already been described as being over-expressed in the yeast form compared to the mycelium phase , is highly down-regulated along the Y-to-H transition (profile 0, 192 times less expressed after 2 h) (Additional file 2: Table S3). Another gene, Onu7070, which codes for a type 1 endochitinase, was previously reported as being over-expressed in the mycelium phase . This gene is also regulated during the Y-to-H switch and ends up being 14 times more strongly expressed at 27 h compared to 0 h (profile 27) (Additional file 2: Table S3). Gene Onu4296 codes for Cerato-ulmin (CU), a hydrophobin that has been well studied in the context of DED. This gene is down-regulated through time and belongs to profiles 0 and 1 (Additional file 2: Table S3).
Identification and expression analysis of S. cerevisiae MAPK, PKA, and RIM pathway homologs
S. cerevisiae MAPK protein homologs
We identified genes coding for potential proteins that were homologous to S. cerevisiae MAPK proteins. Of the 11 S. cerevisiae proteins with homologs in O. novo-ulmi, four have more than one potential homolog (Ste20, Kss1, Sho1 and Msn4; Additional file 2: Table S2). Eight genes encoding these homologous proteins are actually orthologs: Onu0581-Ste20, Onu2999-Ste11, Onu4018-Pbs2, Onu1002-Hog1, Onu2279-Cdc42, Onu7491-Ste12, Onu1392-Ste7, and Onu5199-Tec1 (Additional file 2: Table S2). We found the following eight genes to be differentially expressed during the Y-to-H transition: Onu1002-Hog1, Onu2992-Ste20, Onu1392-Ste7, Onu7491-Ste12, Onu0070-Msn4, Onu8210-Msn4, Onu1780-Kss1 and Onu8395-Kss1 (Fig. 6b). Hog1, Ste20, Ste7, Kss1 and Ste12 homologs are all over-expressed during the Y-to-H transition and are present in STEM profile 39 (Fig. 6b).
S. cerevisiae cAMP-PKA protein homologs
We found homologous genes for the most important genes that were involved in the cAMP-PKA pathway of S. cerevisiae (Bcy1, Tk1-3, Yak1, Phd1, Sok2, Flo8 or Pde2; Additional file 2: Table S2). Only Onu6515-Bcy1 and Onu7499-Pde2 are orthologous genes. Phd1 and Sok2 are putative homologs of the same protein in O. novo-ulmi, Onu0977. Surprisingly, we did not find a homologous protein for Flo11. The gene encoding the best homologous match for Yak1, i.e., Onu0836, is the only one that is differentially expressed under our set of experimental conditions. This gene is part of profile 39 and is already over-expressed almost 4-fold at 4 h following induction of the transition. For genes Onu0977-Phd1/Sok2 and Onu7708-Tpk, large differences among biological replicates might explain why these genes do not pass the threshold set for DEG identification.
S. cerevisiae RIM protein homologs
All six well-described proteins that are involved in the RIM pathway have a homolog in O. novo-ulmi (Additional file 2: Table S2) and five of them are actually encoded by orthologous genes. However, Onu0751-RIM101 is the only gene that is differentially expressed during the Y-to-H transition. Here, this gene is down-regulated (12-fold lower at 27 h compared to 0 h) and is present in profile 0 (Additional file 2: Table S3).
Quantitative Reverse Transcription-PCR (qRT-PCR) analysis of selected genes
In order to confirm RNAseq data, we chose genes that were detected as being DE and tested their expression using qRT-PCR. We selected genes from the MAPK pathway (Onu1392-Ste7, Onu7491-Ste12, Onu8395-Kss1, Onu1002-Hog1 and Onu8210-Msn4) and genes with contrasting expression profiles: over-expressed genes (Onu0980 and Onu6970.1-Hex1) and down-regulated genes (Onu0303 and Onu5171). For all nine genes, we confirmed a differential expression tendency seen in RNAseq by using qRT-PCR with three biological replicates (Additional file 5: Figure S5A-D). Indeed, all genes that were shown to be down-regulated in RNAseq are also found to be down-regulated in qRT-PCR (around a 2-fold decrease at 10 h compared to 0 h). In contrast, all genes that were shown to be upregulated in RNAseq have their expression increased through time in qRT-PCR with comparable fold changes between the two methods.
Large-scale gene expression regulation during yeast-to-hypha (Y-to-H) transition in dimorphic fungi has been studied in model fungal species, such as Saccharomyces cerevisiae and Candida albicans, using RNAseq, DNA microarrays or phenotypic screening of mutants [33, 45, 55, 56]. It also has been recently investigated in Penicillium marneffei . In Ophiostoma novo-ulmi, large-scale transcriptomic analyses have been performed to characterize molecular signatures of yeast and mycelial forms via Expressed Sequence Tags analyses  or, more recently, by RNAseq . Here, we present the first whole-transcriptome study of a time-course experiment in which we monitored the Y-to-H transition in O. novo-ulmi using Ion Torrent RNAseq technology.
Microscopy and flow cytometry
We first characterized the morphological switch process using various quantitative and qualitative techniques (light microscopy, scanning-electron microscopy, and flow cytometry) to assess the relevance of growth conditions that had been selected for the time-course experiment. Results show a greater degree of evidence for Y-to-H transition under the conditions that we have tested. We confirm the emergence of hyphae in culture by showing increasing cell size and the modification of single cell shapes in accordance with previous reports . These changes are seen already within the first 4 h after transfer in a medium that induces Y-to-H. These results could be explained by the phenomenon of cell swelling at 2 h after transfer to static liquid OCM, prior to germination, which has been described earlier .
Global view of transcriptomes of yeast-to-hypha transition
Second, in order to characterize transcriptome modifications during this major morphological change, we used mRNA sequencing. Interestingly, we observed large changes in gene expression among samples that were taken at different times during the Y-to-H transition, with almost 22% of O. novo-ulmi total gene content being differentially expressed. In particular, we found large transcriptomic variation between samples that were collected at the beginning of the experiment and those that were collected 2 h after transfer to static liquid OCM. While this difference agrees with the cell swelling and reshaping observed through microscopy and flow cytometry, it may also be related to the perception of the signal triggering the Y-to-H switch and the regulation of many processes that are required to achieve the transition. In contrast, between 10 and 27 h after transfer to static liquid OCM, molecular variability is relatively small, which could be explained by the fact that the percentage of hyphae that are produced is not increasing. Only the hyphae that were already present in the medium are growing between 10 and 27 h.
By analyzing processes that were enriched in the set of up-regulated genes during the morphological switch, we highlight biological processes that were related to major cell modification and remodeling. For instance, the activation of the cell cycle has already been shown in S. cerevisiae and C. albicans as being critical during morphological changes [58–60]. Furthermore, CAZyme genes are more likely to be actively transcribed than repressed, suggesting an active process of modifications in carbohydrate compounds, which is consistent with cell wall remodeling. By more in depth exploration of genes that are over-expressed during the morphological switch, we found genes that were related to specific hyphal structures, supporting the formation of hyphae. One of them is a gene that is predicted to encode a major protein (Hex1) that is associated with the formation of Woronin bodies. These organelles are associated with septal pores that are present in filamentous fungi [61, 62]. The Hex1 gene appears to be conserved across species and has been shown to be over-expressed in hyphae of Coccidioides spp.. In Neurospora crassa, Woronin bodies have been reported as being responsible for protecting hyphae against cellular damage [64, 65]. Even though their role has been clearly associated with stress response, Woronin bodies are found associated with septal pores but also in the cytoplasm of non-damaged hyphae in A. fumigatus . Finally, in Magnaporthe grisea, the Hex1 protein is apparently required for efficient pathogenesis and in responses to nitrogen starvation, since partial deletion induced a reduction of the virulence and made the fungus unable to survive nitrogen starvation . The Hex1 gene is included in the Pathogen-Host interaction database (PHI-Base, www.phi-base.org) . The high expression of Hex1 in O. novo-ulmi at 27 h after induction of the transition suggests that the expression of this gene is correlated with hyphal formation. Putative actin-, tubulin- and septin-coding genes are also over-expressed through time, consistent with major cellular reorganization due to the formation of pluricellular structures and septa. Septins, in particular, have been reported to play an important role in the filamentous growth of other dimorphic fungi such as Aspergillus nidulans, C. albicans and U. maydis [69–72].
Interestingly, the most highly over-expressed gene, Onu3773, is related to the gene coding for the adhesin protein Mad1. In Metarhizium robertsii (previously M. anisopliae ), which is a well-known insect pathogen that can also colonize plant roots, Mad1 is involved in the adhesion of blastospores to the insect cuticle . Also, the Mad1 gene is highly expressed in the insect hemolymph . It is also responsible for cytoskeleton organization, cell division and actin polymerization. All of these processes are over-represented in putative functions of genes that are over-expressed during the Y-to-H transition in O. novo-ulmi. Deletion of Mad1 (ΔMad1) in M. robertsii induces a delay in conidial germination, a reduction in blastospore formation, and the production of large multicellular hyphal bodies. Moreover, whereas Mad1 is not expressed in freshly collected conidia, its expression increases as soon as conidia start swelling prior to germination . Accordingly, in O. novo-ulmi, Onu3773 is not expressed 2 h after incubation under conditions promoting Y-to-H transition. However, after 4 h of incubation in static OCM, its expression level increases rapidly. This response could be connected with the cell swelling that was observed earlier in this study. Finally, M. robertsii ΔMad1 mutants are less virulent towards insects than the wild-type strain . If the functions of the Mad1 homolog in O. novo-ulmi are conserved, this gene could be a candidate for the link between morphology and virulence, although O. novo-ulmi itself is not entomopathogenic.
In contrast, gene Onu4296, which encodes Cerato-ulmin (CU), is down-regulated under our experimental conditions. Cerato-ulmin is a hydrophobin  that has been studied extensively in the context of DED [77–80]. Its role as a pathogenicity factor remains controversial, but it has been suggested that like Mad1 in M. robertsii, CU promotes the adherence of infectious spores to the cuticle of elm bark beetles . The CU expression level was previously reported to be 20–60% higher in mycelium than in yeast cells that were grown in vitro . Yet the CU gene was equally expressed in RNAseq data of yeast and mycelial forms . Taken together, these contrasting results may reflect medium-dependent regulation and show the lack of conservation of expression for the CU gene among different studies, which in turn explains the controversy regarding the specific role of CU in pathogenicity.
Even though growth conditions and sequencing technologies differ from those that we used to compare yeast and mycelium growth phases , we found convincing overlaps between the two studies. For instance, genes such as Onu4877, Onu5171, Onu6217 and Onu7070, which are already described as being over-expressed in yeast or mycelium growth phases, show the same expression pattern and are over- or down-regulated during the Y-to-H transition. In total, we found 105 genes for which the regulation of gene expression appears to be conserved under variable growth conditions.
Identification and expression analysis of S. cerevisiae MAPK, PKA, and RIM pathway homologs
We had previously shown that transcriptomes of yeast and mycelium phases of O. novo-ulmi are distinct from those of H. capsulatum and C. albicans . However, we also found O. novo-ulmi homologs of genes encoding proteins involved in MAPK, PKA and RIM pathways. The role of the MAPK cascade in the Y-to-H transition has been studied extensively in model dimorphic fungal species such as S. cerevisiae [30, 82] (Fig. 1). We show that 11 key proteins have at least one homolog in O. novo-ulmi, of which four have two putative homologs (Ste20, Kss1, Sho1 and Msn4). More interestingly, the expression of eight genes encoding these proteins is differentially regulated during the Y-to-H transition in O. novo-ulmi. Homologs of Ste12, Kss1, Ste7 and Ste20 are all over-expressed through time, consistent with data that were collected for S. cerevisiae . Ste12 is a transcriptional factor that regulates the expression of multiple genes that are related to the formation of pseudohyphae in S. cerevisiae . The Ste12 homolog in Candida albicans is also involved in hyphal formation and its knockout mutation suppresses filamentation and virulence [39, 84]. During activation of pseudohyphal formation in S. cerevisiae, Ste12 acts as a heterodimer when combined together with Tec1 protein . Here, the Tec1 homologous gene (Onu5199) is not differentially expressed during the Y-to-H transition. This result might suggest that Tec1 is not required in the filamentation process in O. novo-ulmi. In contrast, upregulation of the Ste12 homolog in O. novo-ulmi is consistent with a conserved role in the activation of hypha formation.
Moreover, the Hog1 homolog is also over-expressed during Y-to-H transition in O. novo-ulmi. The MAPK encoded by this gene is described as a key regulator of hyper-osmolarity stress perception and is also related to morphological changes and virulence in C. albicans . In S. cerevisiae, the Hog1 protein regulates invasive growth and contributes to the repression of pseudohyphal growth [33, 86]. Protein Msn4 of S. cerevisiae is a stress-responsive transcriptional activator activated by Hog1. This protein has two putative homologs in O. novo-ulmi (Onu8210 and Onu0070). Interestingly, a paralog of Msn4, known as Msn2, is found in the genome of S. cerevisiae . These two paralogs are largely functionally redundant and regulate stress response . In S. cerevisiae, Msn4 expression is induced by stress, whereas Msn2 expression is constitutive . Under our experimental conditions, the two homologs of Msn4 have contrasting gene expression regulation. Gene Onu8210 is downregulated, whereas Onu0070 is upregulated. The latter is consistent with the upregulation of Hog1 and of the entire MAPK cascade.
Gene Onu8210 encodes a putative transcriptional factor containing a zinc-finger DNA binding domain. Interestingly, this domain is not only homologous with the Msn4 domain but also with those from S. cerevisiae Mig1 and Nrg1 proteins. Both Mig1 and Nrg1 are negative regulators of genes that are involved in multiple processes, including filamentous growth of S. cerevisiae under glucose-limiting conditions [90, 91]. The downregulation of Onu8210 is consistent with the repression of a negative regulator of filamentatous growth, suggesting that this gene might encode a protein that acts more like Mig1 or Nig1 than Msn4.
Taken together, these results constitute the first report revealing potential involvement of the MAPK cascade in the process of morphological change in O. novo-ulmi. Further, they suggest a conserved mechanism across fungal families. Roles and activities for each of the proteins that are encoded by the identified genes must be further confirmed by functional analyses, such as gene knockout experiments, for which efficient protocols are still lacking in O. novo-ulmi.
Under the growth conditions that were selected for Y-to-H induction in O. novo-ulmi, gene Onu0836-Yak1 is the only gene-encoding homolog of proteins involved in the PKA pathway in S. cerevisiae which is differentially expressed. Both genes encoding homologs of Phd1/Sok2 (Onu0977) and Tkp1-3 (Onu7708) showed variation in the number of reads between 0 and 4 h. Phd1 and Sok2 are actually paralogous proteins in S. cerevisiae, which explains why the homologous protein in O. novo-ulmi is the same for both proteins. The fact that all three putative catalytic sub-units of the PKA (Tpk1-3) have a unique homolog in O. novo-ulmi suggests that there is only one catalytic sub-unit in this species. In previous studies , cAMP was implicated in the dimorphism of O. ulmi NRRL6404. Over-expression of the gene encoding the downstream transcription factor of the PKA pathway, Phd1 (S. cerevisiae) or Efg1 (C. albicans), was sufficient to enhance filamentous growth in both species [22, 39, 92]. The slight modification in the expression of Phd1 and Tpk homologs in O. novo-ulmi could reveal either no obvious role for these genes in the Y-to-H transition or that a very small increase in gene expression is sufficient for activation of the positive functions of both genes in the process of Y-to-H transition. In the latter case, the modification of Tpk expression might induce Yak1 expression, which then positively regulates Phd1 expression. The absence of a Flo11 homolog in O. novo-ulmi suggests that there likely is an as yet unknown gene that is targeted by Phd1. Further studies monitoring cAMP concentrations during the Y-to-H transition would aid in demonstrating whether there truly is a relationship between gene expression and cAMP levels.
Finally, analyses of the six putative candidates for the RIM pathway in O. novo-ulmi revealed overexpression of the gene encoding the homolog of RIM101, which is the downstream transcriptional factor regulated by the RIM pathway in S. cerevisiae and C. albicans [28, 40]. We did not monitor changes in pH of the growth medium over the entire experiment. However, we noted that the pH of the 3-day-old pre-culture of yeast cells in OMM was 3.0, whereas fresh OCM that was used to induce transition to mycelium had a pH of 5.5 (data not shown). Therefore, modification in the expression of Onu0751-RIM101 could be the result of the change in pH that was induced by the transfer to OCM, even though this modification is not as drastic as a switch from an acidic to a basic pH. Nevertheless, this also suggests that this gene is involved in the switch and, thus, has a conserved role among dimorphic fungal species.
Our survey of morphological and transcriptomic changes during the yeast to hypha transition in Ophiostoma novo-ulmi highlights modifications in both the structure of cells and the regulation of gene expression within 10 h after induction of the switch. Analyses of gene expression regulation in a time-course experiment provided the first insights into the processes that were activated and those that were repressed during the Y-to-H transition in O. novo-ulmi. These analyses further confirmed the occurrence of major rearrangement of the cell wall and cellular structures. We pointed out both genes and pathways that are either specific to O. novo-ulmi or not described in model dimorphic fungi, (Hex1, Mad1 or CU), or which converge in gene expression regulation for genes encoding proteins that are homologous to Saccharomyces cerevisiae proteins involved in the MAPK cascade related to the yeast to pseudohypha switch. Our experimental conditions also suggest the involvement of PKA and RIM pathways, which should be confirmed in future studies. Taken together, our results indicate many genes of interest for further investigations of the Y-to-H transition in O. novo-ulmi and highlight the possible convergence of molecular regulation of dimorphism between very divergent fungal species.
Materials and methods
Fungal growth conditions and sample collection
Preparation of RNAseq samples
Yeast cells (1 × 106 spores mL−1) of Ophiostoma novo-ulmi ssp. novo-ulmi strain H327 (Centre d’Etude de la Forêt, fungal collection, http://www.cef-cfr.ca/index.php?n=CEF.Collections ChampignonsPathogenes) were first grown in liquid minimum medium (OMM)  that was supplemented with 1.15 g L−1 proline as the sole nitrogen source on an orbital shaker (Infors HT Ecotron) at 150 rpm (22 °C). After 3 days, cultures were filtered through 16 layers of sterile cheesecloth and diluted to 1 × 106 spores mL−1 in flasks containing 100 mL of complete medium (OCM) with 1.15 g L−1 proline. Flasks were incubated under still conditions in darkness at room temperature (22 °C). At each point during the time-course experiment (0, 2, 4, 6, 10 and 27 h after transfer of cells to static liquid OCM), five photographs per flask were taken with a microscope (Olympus BX41, 400 × magnification) to measure the proportion of yeasts and hyphae. Cells were considered hyphal structures when they were at least twice as long as a typical yeast cell. About 60 structures were counted and measured at each time. One hundred milliliters of culture from each flask were transferred each time into two sterile Falcon 50 mL conical centrifuge tubes and centrifuged for 5 min at 4 °C and 3000 g. Pellets were resuspended in sterile water and transferred to 1.5 mL sterile microtubes. Samples were again centrifuged for 5 min at 4 °C and 3000 g, after which the pellets were flash-frozen in liquid N2. Samples were stored at −80 °C until RNA extraction. Three time series replicates were collected (18 samples in total).
Time points for flow cytometry measurements were 0, 2, 4, 6 and 10 h after transfer of cells to static liquid OCM. Samples that were incubated for 27 h in OCM were not analyzed since cell size and shape were no longer compatible with the flow cytometer. Cells were grown in OCM under the same conditions as for the RNAseq samples. At each time tested, cultures were diluted with sterile water in order to have 50–500 cells per μL in 300 μL of solution that was transferred to 96-well plates for analysis in a Guava® easyCyte HT Sampling flow cytometer (EMD Millipore Corp.). OMM that was supplemented with 1.15 g L−1 proline was used as negative control for the transition [7, 8].
Scanning electron microscopy (SEM)
Cells were grown and collected following the same procedure as for RNA extractions. Pellets of cells were resuspended in fixation buffer (2.5% glutaraldehyde in 0.1 M sodium cacodylate buffer, pH 7.3). Samples were fixed for 20 min at room temperature and stored at 4 °C for at least 12 h. Further processing of samples was performed at the Plateforme de Microscopie (IBIS/Université Laval) following the same procedure as described previously . Samples were examined with a JEOL JSM6360LV microscope (JEOL Inc., Tokyo, Japan).
RNA extraction, cDNA library production and RNA sequencing
Total RNA from each of the 18 samples was extracted using the Ambion® Ribopure™ Yeast kit (Life Technologies) (for schematic see Additional file 6: Figure S1). The quality was verified using a spectrophotometer (NanoDrop ND-1000, Thermo Scientific) and a bioanalyser (BioAnalyzer RNA 6000 n Kit, Agilent Technologies).
Complementary DNA (cDNA) libraries were constructed at the Plateforme d’Analyses Génomiques (IBIS/Université Laval) using the kapa stranded RNA-seq kit (Kapa Biosystems) with Y-adaptors that were custom made to be compatible with Ion Proton™ (Thermo Fisher Scientific) sequencing. Single-end RNA sequencing was performed following the manufacturer’s recommendations.
RNAseq data preprocessing
RNAseq read quality was verified using FastQC software (http://www.bioinformatics.babraham.ac.uk/projects/fastqc/). Given the poor quality of 5’ and 3’ tails, 15 nt were removed in 5’ and sequences were trimmed to a length of 150 nt using Prinseq v0.20.4 . Finally, sequences shorter than 25 nt were removed.
RNAseq data analysis
Filtered reads were mapped onto the exon sequences of the O. novo-ulmi H327 genome  with TopHat2 (v.2.0.10) using default parameters . For greater readability in the manuscript, genes of O. novo-ulmi are named “Onu” followed with the number of the gene. All further analyses were performed in the R statistical environment (v3.0.1; ). Raw data were filtered and transformed with the Bioconductor package EdgeR (correction by the library size and filtering of genes for which the row sum was less than 3 counts per million [CPM]) . A multidimensional scaling (MDS) plot was used for the visualization of the variability among samples (function implemented in EdgeR), following the procedure described previously . The MDS plot was produced using the R package ggplot2 . The updated version of maSigPro R package , Next-maSigPro , was used for the detection of differentially expressed genes (DEGs) in R. P-values were corrected using the False Discovery Rate (FDR)  with a threshold of 0.05.
Clustering analysis was performed using Short Time-series Expression Miner (STEM) software  on the set of DEGs with our own annotation file. The software assigns genes to a set of model profiles predefined in order to capture all the potential distinct patterns that can be expected from the experiment . We chose a maximum number of 50 for model profiles with a significance threshold of 0.01 and False Rate Discovery (FDR) correction. Determination of differentially expressed genes was based on ‘Maximum − Minimum’ that is the largest difference in the gene’s value between any two time points, which are not necessarily consecutive . A minimum fold change of two was required over time between normalized count values for a gene to be retained. Differential expression for each gene was calculated on the normalized count values, as:
where v(i) is the average number of normalized counts at a time point of interest, and v(0) is the average number of normalized counts at the beginning of the experiment. Values less than 0 were inverted to ease of interpretation; in this case, the equation was written as:
Gene ontology (GO) term enrichment analyses were performed as described in Nigg and collaborators (2015) using the Goseq package in Bioconductor 3.0  and GO Trimming . REVIGO was used to reduce redundancy and produce treemaps (space-filling visualization of hierarchical structures) using the treemap package in R . Default parameters were used and the redundancy threshold (allowed similarity) was set at 0.4 (very small).
Protein homology and gene orthology assessment
Sequences of MAPK, PKA and RIM proteins from S. cerevisiae were downloaded from NCBI and the yeast genome database (www.yeastgenome.com). A standalone version of BLASTp  was used (parameters: word-size = 3; e-value = 0.01) to assess homology with O. novo-ulmi proteins.
Orthologous genes among the identified homologs were highlighted by reciprocal best blast hits (RBH) using tblastx . Exon sequences of the two species (O. novo-ulmi and S. cerevisiae) were compared using local databases and the RBH was conducted using the following parameters: e-value = 0.001 and word size = 5.
Reverse transcription, PCR and quantitative Reverse Transcription-PCR (qRT-PCR)
Samples for qRT-PCR were collected following the procedure that was described above. Total RNA was extracted and samples were conserved at −20 °C for a maximum of 1 month prior to complementary DNA (cDNA) synthesis.
The cDNAs for qRT-PCR were synthesized by reverse transcription using Superscript II enzyme (Invitrogen) and according to the following protocol: 500 ng of mRNA (final volume of 8 μL, completed with nuclease-free water) were mixed with 1 μL of dNTPs 10 mM and 1 μL of oligo-dT 0.5 μg/μL. The mixture was incubated at 65 °C for 5 min. Samples were then cooled on ice for 1 min and quickly centrifuged. A mixture of 10X Buffer RT (2 μL), MgCl2 (25 mM, 4 μL), DTT (0.1 M, 2 μL), RNase Out (1 μL, Invitrogen) and Superscript II (200 U, 0.25 μL) was added to each sample. Preparations were gently mixed, quickly centrifuged and incubated for 50 min at 42 °C. Reactions were completed by incubating samples for 15 min at 72 °C. Samples were finally cooled down on ice for 5 min and 80 μL of nuclease-free water were added. cDNA samples were stored at −20 °C.
For each gene that was tested in qRT-PCR, primers for fragments that were 120–200 base pairs in length were designed using Primer3plus (http://primer3plus.com/cgi-bin/dev/primer3plus.cgi) and properties were verified with Oligocalc (http://www.basic.northwestern.edu/biotools/oligocalc.html) (Additional file 6: Table S1). Primers were then synthesized at the Plateforme de séquençage et de génotypage des génomes (Research Centre, Centre Hospitalier de l’Université Laval, Quebec QC, CA).
Prior to qRT-PCR, primers were tested on cDNA (9 ng μL−1) in PCR. Reaction volume of each sample was 25 μL (12.2 μL Premix Ex Taq TaKaRa, 1.25 μL of each primer, 5 μL of sterile water and 5 μL of cDNA). PCR conditions for amplification were: 10 s at 98 °C followed by 30 cycles of 10 s at 98 °C, annealing for 30 s at Tm (usually 55 °C), amplification for 12 s at 72 °C and a final step of 5 min at 72 °C. PCR products were transferred to 2% (m/v) agarose gel for electrophoresis (migration at 130 V) followed by staining in ethidium bromide bath for 15 min. DNA fragments were visualized under UV light.
Quantitative RT-PCR (qRT-PCR) reactions were performed with the PerfeCta™ SYBR® Green FastMix™, Low ROX kit (Quanta BioSciences) following the manufacturer’s instructions, with 10 μL of PerfeCta™ SYBR® Green FastMix™, Low ROX mix, 2.5 μL of 5’ and 3’ 1 μM primers for each gene that was tested and 5 μL of 1 ng μL−1 cDNA (or water for negative control) for a final volume of 20 μL in MicroAmp® Fast Optical 96-well plates (Applied Biosystems™). Gene expression was quantified using the 7500 Fast Real-Time PCR system (Applied Biosystems™). Genes Onu3626 (Ubiquitin conjugating enzyme), Onu1683 (Ubiquitin conjugating enzyme E2 6) and Onu0623 (RING finger ubiquitin ligase) were defined as reference genes using the R script of NormFinder software . For each of the nine genes that were tested (Onu1392, Onu7491, Onu8395, Onu1002, Onu8210, Onu0980, Onu6970.1, Onu0303 and Onu5171), three biological replicates with three technical replicates were used to quantify expression with a comparative method (2-ΔΔCt). The three reference genes were used to calculate the level of transcripts at each time point (2 h, 4 h, 6 h and 10 h) relative to 0 h (beginning of the Y-to-H transition experiment).
Brasier CM. Ophiostoma novo-ulmi sp.nov., causative agent of current Dutch elm disease pandemics. Mycopathologia. 1991;115:151–61.
Campana RJ. Inoculation and fungal invasion of the tree. In: Sinclair WA, Campana RJ, editors. Dutch elm Dis. Perspect. after 60 years. Ithaca: Cornell Un; 1978. p. 17–20.
Webber JF, Brasier CM. The transmission of Dutch elm disease: a study of the process involved. In: Anderson JM, Rayner ADM, Walton DWH, editors. Invertebrate-microbial interactions. Cambridge: Cambridge University Press; 1984. p. 271–306.
Dogra N, Breuil C. Suppressive subtractive hybridization and differential screening identified genes differentially expressed in yeast and mycelial forms of Ophiostoma piceae. FEMS Microbiol Lett. 2004;238:175–81.
Krokene P, Solheim H. Pathogenicity of four blue-stain fungi associated with aggressive and nonaggressive bark beetles. Phytopathology. 1998;88:39–44.
Halmschlager E, Messner R, Kowalski T, Prillinger H. Differentiation of Ophiostoma piceae and Ophiostoma quercus by morphology and RAPD analysis. Syst Appl Microbiol. 1994;17:554–62.
Naruzawa ES, Bernier L. Control of yeast-mycelium dimorphism in vitro in Dutch elm disease fungi by manipulation of specific external stimuli. Fungal Biol. 2014;118:872–84.
Kulkarni RK, Nickerson KW. Nutritional control of dimorphism in Ceratocystis ulmi. Exp Mycol. 1981;5:148–54.
Dalpé Y. L’influence de la carence en pyridoxine sur la morphologie et l'ultrastructure cellulaire de Ceratocystis ulmi. Can J Bot. 1983;61:2079–84.
Naruzawa ES, Malagnac F, Bernier L. Effect of linoleic acid on reproduction and yeast-mycelium dimorphism in the Dutch elm disease pathogens. Botany. 2016;96:31–9.
Muthukumar G, Kulkarni RK, Nickerson KW. Calmodulin levels in the yeast and mycelial phases of Ceratocystis ulmi. J Bacteriol. 1985;162:47–9.
Muthukumar G, Nickerson KW. Ca(II)-calmodulin regulation of fungal dimorphism in Ceratocystis ulmi. J Bacteriol. 1984;159:390–2.
Brunton AH, Gadd GM. The effect of exogenously-supplied nucleosides and nucleotides and the involvement of adenosine 3’:5'-cyclic monophosphate (cyclic AMP) in the yeast mycelium transition of Ceratocystis (= Ophiostoma) ulmi. FEMS Microbiol Lett. 1989;60:49–53.
Hornby JM, Jacobitz-Kizzier SM, McNeel DJ, Jensen EC, Treves DS, Nickerson KW. Inoculum size effect in dimorphic fungi: extracellular control of yeast-mycelium dimorphism in Ceratocystis ulmi. Appl Environ Microbiol. 2004;70:1356–9.
Berrocal A, Navarrete J, Oviedo C, Nickerson KW. Quorum sensing activity in Ophiostoma ulmi: effects of fusel oils and branched chain amino acids on yeast-mycelial dimorphism. J Appl Microbiol. 2012;113:126–34.
Wedge MÈ, Naruzawa ES, Nigg M, Bernier L. Diversity in yeast–mycelium dimorphism response of the Dutch elm disease pathogens: the inoculum size effect. Can J Microbiol. 2016;62:1–5.
Berrocal A, Oviedo C, Nickerson KW, Navarrete J. Quorum sensing activity and control of yeast-mycelium dimorphism in Ophiostoma floccosum. Biotechnol Lett. 2014;36:1503–13.
de Salas F, Martínez MJ, Barriuso J. Quorum sensing mechanisms mediated by farnesol in Ophiostoma piceae: its effect on the secretion of sterol esterase. Appl Environ Microbiol. 2015;81:4351–7.
Pereira V, Royer JC, Hintz WE, Field D, Bowden C, Kokurewicz K, et al. A gene associated with filamentous growth in Ophiostoma novo-ulmi has RNA-binding motifs and is similar to a yeast gene involved in mRNA splicing. Curr Genet. 2000;37:94–103.
Richards WC, Takai S, Lin D, Hiratsuka Y, Asina S. An abnormal strain of Ceratocystis ulmi incapable of producing external symptoms of Dutch elm disease. Eur J For Pathol. 1982;12:193–202.
Richards WC. Nonsporulation in the Dutch elm disease fungus Ophiostoma ulmi: evidence for control by a single nuclear gene. Rev Can Bot. 1994;72:461–7.
Gancedo JM. Control of pseudohyphae formation in Saccharomyces cerevisiae. FEMS Microbiol Rev. 2001;25:107–23.
Sánchez-Martínez C, Pérez-Martín J. Dimorphism in fungal pathogens: Candida albicans and Ustilago maydis – similar inputs, different outputs. Curr Opin Microbiol. 2001;4:214–21.
Nadal M, García-Pedrajas MD, Gold SE. Dimorphism in fungal plant pathogens. FEMS Microbiol Lett. 2008;284:127–34.
Martínez-Espinoza AD, Ruiz-Herrera J, León-Ramírez CG, Gold SE. MAP kinase and cAMP signaling pathways modulate the pH-induced yeast-to-mycelium dimorphic transition in the corn smut fungus Ustilago maydis. Curr Microbiol. 2004;49:274–81.
Agarwal C, Aulakh KB, Edelen K, Cooper M, Wallen RM, Adams S, et al. Ustilago maydis phosphodiesterases play a role in the dimorphic switch and in pathogenicity. Microbiology. 2013;159:857–68.
Sonneborn A, Bockmühl DP, Gerads M, Kurpanek K, Sanglard D, Ernst JF. Protein kinase A encoded by TPK2 regulates dimorphism of Candida albicans. Mol Microbiol. 2000;35:386–96.
Selvig K, Alspaugh JA. pH response pathways in fungi: adapting to host-derived and environmental signals. Mycobiology. 2011;39:249–56.
Peñalva MA, Tilburn J, Bignell E, Arst HN. Ambient pH gene regulation in fungi: making connections. Trends Microbiol. 2008;16:291–300.
Madhani HD, Fink GR. The control of filamentous differentiation and virulence in fungi. Trends Cell Biol. 1998;8:348–53.
Hamel L-P, Nicole M-C, Duplessis S, Ellis BE. Mitogen-activated protein kinase signaling in plant-interacting fungi: distinct messages from conserved messengers. Plant Cell. 2012;24:1327–51.
Ferrigno P, Posas F, Koepp D, Saito H, Silver PA. Regulated nucleo/cytoplasmic exchange of Hog1 MAPK requires the importin beta homologs Nmd5 and Xpo1. EMBO J. 1998;17:5606–14.
Shively CA, Eckwahl MJ, Dobry CJ, Mellacheruvu D, Nesvizhskii A, Kumar A. Genetic networks inducing invasive growth in Saccharomyces cerevisiae identified through systematic genome-wide overexpression. Genetics. 2013;193:1297–310.
Alonso-Monge R, Navarro-García F, Molero G, Diez-Orejas R, Gustin M, Pla J, et al. Role of the mitogen-activated protein kinase Hog1p in morphogenesis and virulence of Candida albicans. J Bacteriol. 1999;181:3058–68.
Lengeler KB, Davidson RC, D’souza C, Harashima T, Shen WC, Wang P, et al. Signal transduction cascades regulating fungal development and virulence. Microbiol Mol Biol Rev. 2000;64:746–85.
Pan X, Heitman J. Cyclic AMP-dependent protein kinase regulates pseudohyphal differentiation in Saccharomyces cerevisiae. Mol Cell Biol. 1999;19:4874–87.
Robertson LS, Fink GR. The three yeast A kinases have specific signaling functions in pseudohyphal growth. Proc Natl Acad Sci U S A. 1998;95:13783–7.
Malcher M, Schladebeck S, Mösch HU. The Yak1 protein kinase lies at the center of a regulatory cascade affecting adhesive growth and stress resistance in Saccharomyces cerevisiae. Genetics. 2011;187:717–30.
Lo HJ, Köhler JR, Didomenico B, Loebenberg D, Cacciapuoti A, Fink GR. Nonfilamentous C. albicans mutants are avirulent. Cell. 1997;90:939–49.
Davis D. Adaptation to environmental pH in Candida albicans and its relation to pathogenesis. Curr Genet. 2003;44:1–7.
Penãlva MA, Arst HN. Regulation of gene expression by ambient pH in filamentous fungi and yeasts. Microbiol Mol Biol Rev. 2002;66:426–46.
Yang E, Chow W, Wang G, Woo PCY, Lau SKP, Yuen K, et al. Signature gene expression reveals novel clues to the molecular mechanisms of dimorphic transition in Penicillium marneffei. PLoS Genet. 2014;10, e1004662.
Nantel A, Dignard D, Bachewich C, Harcus D, Marcil A, Bouin A-P, et al. Transcription profiling of Candida albicans cells undergoing the yeast-to-hyphal transition. Mol Biol Cell. 2002;13:3452–65.
Kadosh D, Johnson A. Induction of the Candida albicans filamentous growth program by relief of transcriptional repression: a genome-wide analysis. Mol Biol Cell. 2005;16:2903–12.
Carlisle PL, Kadosh D. A genome-wide transcriptional analysis of morphology determination in Candida albicans. Mol Biol Cell. 2012;24:246–60.
Forgetta V, Leveque G, Dias J, Grove D, Lyons R, Genik S, et al. Sequencing of the Dutch elm disease fungus genome using the Roche/454 GS-FLX Titanium System in a comparison of multiple genomics core facilities. J Biomol Tech. 2013;24:39–49.
Comeau AM, Dufour J, Bouvet GF, Jacobi V, Nigg M, Henrissat B, et al. Functional annotation of the Ophiostoma novo-ulmi genome: insights into the phytopathogenicity of the fungal agent of Dutch elm disease. Genome Biol Evol. 2015;7:410–30.
Nigg M, Laroche J, Landry CR, Bernier L. RNAseq analysis highlights specific transcriptome signatures of yeast and mycelial growth phases in the Dutch elm disease fungus Ophiostoma novo-ulmi. G3. 2015;5:2487–95.
Nueda MJ, Tarazona S, Conesa A. Next maSigPro: updating maSigPro bioconductor package for RNA-seq time series. Bioinformatics. 2014;30:2598–602.
Conesa A, Nueda MJ, Ferrer A, Talon M. maSigPro: a method to identify significantly differential expression profiles in time-course microarray experiments. Bioinformatics. 2006;22:1096–102.
Benjamini Y, Hochberg Y. Controlling the false discovery rate: a practical and powerful approach to multiple testing. J R Statistical Soc. 1995;57:289–300.
Ernst J, Bar-Joseph Z. STEM: a tool for the analysis of short time series gene expression data. BMC Bioinformatics. 2006;7:191.
Ernst J, Nau GJ, Bar-Joseph Z. Clustering short time series gene expression data. Bioinformatics. 2005;21:159–68.
Supek F, Bošnjak M, Škunca N, Šmuc T. REVIGO summarizes and visualizes long lists of gene ontology terms. PLoS One. 2011;6, e21800.
Desai PR, van Wijlick L, Kurtz D, Juchimiuk M, Ernst JF. Hypoxia and temperature regulated morphogenesis in Candida albicans. PLoS Genet. 2015;11, e1005447.
Jin R, Dobry CJ, McCown PJ, Kumar K. Large-scale analysis of yeast filamentous growth by systematic gene disruption and overexpression. Mol Biol Cell. 2008;19:284–96.
Hintz W, Pinchback M, de la Bastide P, Burgess S, Jacobi V, Hamelin R, et al. Functional categorization of unique expressed sequence tags obtained from the yeast-like growth phase of the elm pathogen Ophiostoma novo-ulmi. BMC Genomics. 2011;12:431.
Lew DJ, Reed SI. Morphogenesis in the yeast cell cycle: regulation by Cdc28 and cyclins. J Cell Biol. 1993;120:1305–20.
Wang Y. CDKs and the yeast-hyphal decision. Curr Opin Microbiol. 2009;12:644–9.
Moseley JB, Nurse P. Cdk1 and cell morphology: connections and directions. Curr Opin Cell Biol. 2009;21:82–8.
Woronin M. Zur Entwicklungsgeschichte der Ascobolus pulcherrimus Cr. und einer Pezizen. Abh. Senkend. Naturforsch. 1964;5:333–44.
Reichle RE, Alexander JV. Multiperforate septations, woronin bodies, and septal plugs in Fusarium. J. Cell Biol. 1965;24:498–96.
Whiston E, Zhang Wise H, Sharpton TJ, Jui G, Cole GT, Taylor JW. Comparative transcriptomics of the saprobic and parasitic growth phases in Coccidioides spp. PLoS One. 2012;7, e41034.
Jedd G, Chua NH. A new self-assembled peroxisomal vesicle required for efficient resealing of the plasma membrane. Nat Cell Biol. 2000;2:226–31.
Tenney K, Hunt I, Sweigard J, Pounder JI, McClain C, Bowman EJ, et al. Hex-1, a gene unique to filamentous fungi, encodes the major protein of the Woronin body and functions as a plug for septal pores. Fungal Genet Biol. 2000;31:205–17.
Beck J, Ebel F. Characterization of the major Woronin body protein HexA of the human pathogenic mold Aspergillus fumigatus. Int J Med Microbiol. 2013;303:90–7.
Soundararajan S, Jedd G, Li X, Ramos-pamplon M, Chua NH, Naqvi NI. Woronin body function in Magnaporthe grisea is essential for efficient pathogenesis and for survival during nitrogen starvation stress. Plant Cell. 2004;16:1564–74.
Winnenburg R, Urban M, Beacham A, Baldwin TK, Holland S, Lindeberg M, et al. PHI-base update: additions to the pathogen-host interaction database. Nucleic Acids Res. 2008;36:572–6.
Gladfelter AS. Control of filamentous fungal cell shape by septins and formins. Nat Rev Microbiol. 2006;4:223–9.
Lindsey R, Ha Y, Momany M. A septin from the filamentous fungus A. nidulans induces atypical pseudohyphae in the budding yeast S. cerevisiae. PLoS One. 2010;5:1–9.
Boyce KJ, Chang H, Souza CAD, Kronstad JW. An Ustilago maydis septin is required for filamentous growth in culture and for full symptom development on maize. Eukaryot Cell. 2005;4:2044–56.
Warenda AJ, Konopka JB. Septin function in Candida albicans morphogenesis. Mol Biol Cell. 2003;13:2732–46.
Bischoff JF, Rehner SA, Humber RA. A multilocus phylogeny of the Metarhizium anisopliae lineage. Mycologia. 2009;101:512–30.
Wang C, St Leger RJ. The MAD1 adhesin of Metarhizium anisopliae links adhesion with blastospore production and virulence to insects, and the MAD2 adhesin enables attachment to plants. Eukaryot. Cell. 2007;6:808–16.
Wang C, Hu G, St Leger RJ. Differential gene expression by Metarhizium anisopliae growing in root exudate and host (Manduca sexta) cuticle or hemolymph reveals mechanisms of physiological adaptation. Fungal Genet. Biol. 2005;42:704–18.
Stringer MA, Timberlake WE. Cerato-ulmin, a toxin involved in Dutch elm disease, is a fungal hydrophobin. Plant Cell. 1993;5:145–6.
Bowden CG, Hintz WE, Jeng R, Hubbes M, Horgen PA. Isolation and characterization of the cerato-ulmin toxin of the Dutch elm disease pathogen. Ophiostoma ulmi Curr Genet. 1994;75:323–9.
Takai S, Richards WC, Stevenson KJ. Evidence for the involvement of cerato-ulmin, the Ceratocystis ulmi toxin, in the development of Dutch elm disease. Physiol Plant Pathol. 1983;23:275–80.
Temple B, Horgen PA, Bernier L, Hintz WE. Cerato-ulmin, a hydrophobin secreted by the causal agents of Dutch elm disease, is a parasitic fitness factor. Fungal Genet Biol. 1997;22:39–53.
Sherif S, Jones AMP, Shukla MR, Saxena PK. Establishment of invasive and non-invasive reporter systems to investigate American elm-Ophiostoma novo-ulmi interactions. Fungal Genet Biol. 2014;71:32–41.
Tadesse Y, Bernier L, Hintz WE, Horgen PA. Real time RT-PCR quantification and Northern analysis of cerato-ulmin (CU) gene transcription in different strains of the phytopathogens Ophiostoma ulmi and O. novo-ulmi. Mol. Genet. Genomics. 2003;269:789–96.
Cook JG, Bardwell L, Thorner J. Inhibitory and activating functions for MAPK Kss1 in the S. cerevisiae filamentous-growth signalling pathway. Nature. 1997;390:85–8.
Liu H, Styles CA, Fink GR. Elements of the yeast pheromone response pathway required for filamentous growth of diploids. Science. 1993;262:1741–4.
Liu H, Köhler J, Fink GR. Suppression of hyphal formation in Candida albicans by mutation of a STE12 homolog. Science. 1994;266:1723–6.
Madhani HD, Fink GR. Combinatorial control required for the specificity of yeast MAPK signaling. Science. 1997;275:1314–7.
O’ Rourke SM, Herskowitz I. The Hog1 MAPK prevents cross talk between the HOG and pheromone response MAPK pathways in Saccharomyces cerevisiae. Genes Dev. 1998;12:2874–86.
Estruch F, Carlson M. Two homologous zinc finger genes identified by multicopy suppression in a SNF1 protein kinase mutant of Saccharomyces cerevisiae. Mol Cell Biol. 1993;13:3872–81.
Martínez-Pastor MT, Marchler G, Schüller C, Marchler-Bauer A, Ruis H, Estruch F. The Saccharomyces cerevisiae zinc finger proteins Msn2p and Msn4p are required for transcriptional induction through the stress response element (STRE). EMBO J. 1996;15:2227–35.
Gasch AP, Spellman PT, Kao CM, Carmel-Harel O, Eisen MB, Storz G, et al. Genomic expression programs in the response of yeast cells to environmental changes. Mol Biol Cell. 2000;11:4241–57.
Kuchin S, Vyas VK, Carlson M. Snf1 protein kinase and the repressors Nrg1 and Nrg2 regulate FLO11, haploid invasive growth, and diploid pseudohyphal differentiation. Mol Cell Biol. 2002;22:3994–4000.
Karunanithi S, Cullen PJ. The filamentous growth MAPK pathway responds to glucose starvation through the Mig1/2 transcriptional repressors in Saccharomyces cerevisiae. Genetics. 2012;192:869–87.
Gimeno CJ, Fink GR. Induction of pseudohyphal growth by overexpression of PHD1, a Saccharomyces cerevisiae gene related to transcriptional regulators of fungal development. Mol Cell Biol. 1994;14:2100–12.
Bernier L, Hubbes M. Mutations in Ophiostoma ulmi induced by N-methyl-N’-nitro-N-nitrosoguanidine. Can J Bot. 1990;68:225–31.
Aoun M, Rioux D, Simard M, Bernier L. Fungal colonization and host defense reactions in Ulmus americana callus cultures inoculated with Ophiostoma novo-ulmi. Phytopathology. 2009;99:642–50.
Schmieder R, Edwards R. Quality control and preprocessing of metagenomic datasets. Bioinformatics. 2011;27:863–4.
Kim D, Pertea G, Trapnell C, Pimentel H, Kelley R, Salzberg SL. TopHat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions. Genome Biol. 2013;14:R36.
R Development Core Team. R: A language and environment for statistical computing. Vienna: R Foundation for Statistical Computing; 2015. http://www.r-project.org/.
Robinson MD, McCarthy DJ, Smyth GK. EdgeR: a bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2009;26:139–40.
Wickham H. ggplot2: elegant graphics for data analysis. New York: Springer; 2009.
Young MD, Wakefield MJ, Smyth GK, Oshlack A. Gene ontology analysis for RNA-seq: accounting for selection bias. Genome Biol. 2010;11:R14.
Jantzen SG, Sutherland BJ, Minkley DR, Koop BF. GO Trimming: Systematically reducing redundancy in large Gene Ontology datasets. BMC Res Notes. 2011;4:267.
Altschul SF, Gish W, Miller W, Myers EW, Lipman DJ. Basic local alignment search tool. J Mol Biol. 1990;215:403–10.
Andersen CL, Jensen JL, Ørntoft TF. Normalization of real-time quantitative reverse transcription-PCR data: a model-based variance estimation approach to identify genes suited for normalization, applied to bladder and colon cancer data sets. Cancer Res. 2004;64:5245–50.
We thank Brian Boyle and staff from the Plateforme d’Analyses Génomiques (IBIS, Université Laval, Quebec, QC, CA) for the production of RNAseq samples and Richard Janvier from the Plateforme de Microscopie (IBIS, Université Laval, Quebec, QC, CA) for the SEM samples and pictures. We are indebted to Pr. Christian R. Landry and members of his group (IBIS, Université Laval, Quebec, QC, CA) for help with the design of the experiment and technical support during flow cytometry experiments. We also thank Jean-Guy Catford (CEF/IBIS, Université Laval, Quebec, QC, CA) for assistance with RT-qPCR experiments and Isabelle Giguère (CEF/IBIS, Université Laval, Quebec, QC, CA) for technical support in molecular biology. We thank François-Olivier Hébert Gagnon (IBIS, Université Laval, Quebec, QC, CA) and Jean-Baptiste Leducq (Université de Montréal, Montreal, QC, CA) for their useful comments on the manuscript. Finally, we acknowledge Dr. William F. J. Parsons (CEF, Université de Sherbrooke, Sherbrooke, QC, CA) who edited the English.
Our work was funded by the Natural Sciences and Engineering Research Council of Canada (NSERC, Discovery Grant 105519 to L.B.).
Availability of data and materials
The data set supporting the results of this article is included within the article and its additional files. The raw data from the 18 RNAseq samples have been submitted to the National Center for Biotechnology Information (NCBI) in the BioProject PRJNA325932, study accession SRP078831.
MN and LB planned the experiments. MN did the laboratory work and the analyses. MN wrote the paper. LB provided the funding and helped with the writing. Both authors have read and approved the final version of the manuscript.
The authors declare that they have no competing interests.
Consent for publication
Ethics approval and consent to participate
Distribution of the number of reads (log10 scale) per gene, per condition (mean of 3 replicates). (PNG 76 kb)
Table S1: qRT-PCR primers for reference genes and tested genes. Table S2: Results of BLASTp analysis of Saccharomyces cerevisiae dimorphism related proteins with Ophiostoma novo-ulmi genome. The last column indicates if the homologous proteins are orthologs (results of reciprocal best blast hits (RBH) using tblastx on exon sequences for each gene of the two species, e-value = 1x10^3 and word size = 5). Table S3: Normalized read counts for the 7605 genes analyzed using EdgeR and maSigPro in R. 1-3: 0 h; 4-6: 2 h; 7-9: 4 h; 10-12: 6 h; 13-15: 10 h; 16-18: 27 h. DEGs : Differentially Expressed Genes detected using maSigPro, FDR ≤ 0.05. For each DEG, the number of the STEM profile associated is indicated. Nigg et al. 2015: genes overexpressed in Yeast or Mycelium growth phases in the study of Nigg and collaborators in 2015 . Table S4: Description of the 815 genes found in the STEM profile 39. Values are means of differential gene expressions for each time compared to 0 h.Table S5: Enriched biological processes in the STEM profile 39 containing at least 4 DEGs genes after filtration using GO Trimming and REVIGO. GO term in bold are parent term retained as name for a regroupment. Table S6: Description of the 256 genes found in the STEM profile 8. Values are means of differential gene expressions for each time compared to 0 h. Table S7: Enriched biological processes in the STEM profile 8 containing at least 4 DEGs genes after filtration using GO Trimming and REVIGO. GO term in bold are parent term retained as name for a regroupment. (XLSX 1029 kb)
REVIGO Gene Ontology treemap for the 112 biological processes that were enriched in profile 39. The size of boxes represents the absolute P-value for enrichment of each GO term in the gene set of the profile (log10). (PNG 39 kb)
REVIGO Gene Ontology treemap for the 22 biological processes enriched in profile 8. The size of boxes represents the absolute P-value for enrichment of each GO term in the gene set of the profile (log10). (PNG 27 kb)
Quantitative Reverse Transcriptase-Polymerase Chain Reaction (qRT-PCR) results at 0, 2, 4, 6, and 10 h following transfer to fresh complete growth medium (OCM) for A, downregulated genes; B, C and D, over-expressed genes. Differential expression: 2-ΔΔCt, fold increase in expression level compared to the start of the experiment (0 h). Standard deviations are calculated on three biological replicates for each gene. Three technical replicates per gene were run. Transcript levels were normalized with three control genes (Onu3626, Onu1683, Onu0623) that had the most stable expression and were used to calculate a reliable normalization factor according to Normfinder software . (PNG 161 kb)
Workflow for RNAseq data production and analysis. (PNG 115 kb)