Evolutionary divergence in the fungal response to fluconazole revealed by soft clustering

Abstract

Background: Fungal infections are an emerging health risk, especially those involving yeast that are resistant to antifungal agents. To understand the range of mechanisms by which yeasts can respond to anti-fungals, we compared gene expression patterns across three evolutionarily distant species - Saccharomyces cerevisiae, Candida glabrata and Kluyveromyces lactis - over time following fluconazole exposure.

Results: Conserved and diverged expression patterns were identified using a novel soft clustering algorithm that concurrently clusters data from all species while incorporating sequence orthology. The analysis suggests complementary strategies for coping with ergosterol depletion by azoles - Saccharomyces imports exogenous ergosterol, Candida exports fluconazole, while Kluyveromyces does neither, leading to extreme sensitivity. In support of this hypothesis we find that only Saccharomyces becomes more azole resistant in ergosterol-supplemented media; that this depends on sterol importers Aus1 and Pdr11; and that transgenic expression of sterol importers in Kluyveromyces alleviates its drug sensitivity.

Conclusions: We have compared the dynamic transcriptional responses of three diverse yeast species to fluconazole treatment using a novel clustering algorithm. This approach revealed significant divergence among regulatory programs associated with fluconazole sensitivity. In future, such approaches might be used to survey a wider range of species, drug concentrations and stimuli to reveal conserved and divergent molecular response pathways.

Figure 1

Differentially expressed genes. (a) Number of differentially expressed (up- and down-regulated) genes by species versus the number of cell doublings. (b) Venn diagram showing the overlap in the sets of differentially expressed genes selected in each species at a false discovery rate of q ≤ 0.1. The number of differentially expressed genes in each region of the Venn diagram is not identical across species, since the number of genes that a species contributes to an orthologous group (that is, number of paralogs) can vary. Ratios in parentheses indicate the number of differentially expressed orthologs by the total number of differentially expressed genes (not all genes possess orthologs).

Figure 2

Soft clustering method. (a) Standard clustering based on expression only: two sets of orthologs are depicted (color represents orthology, shape represents species) where orthologs are split between clusters 1 and 2. For illustrative purposes, only two time points (t and t + 1) are shown. (b) Soft clustering based on expression and orthology: dashed circles denote regions where orthologs will be co-clustered. Since the purple square has no orthologs in cluster 1, it remains assigned to cluster 2. (c) Effect of number of clusters k and orthology weight W on GO term enrichment. (d) The number of enriched GO terms, variance, and fraction of co-clustered orthologs for k = 17 as a function of W in comparison to randomized paralogs/orthologs. Randomization was performed as described in Additional file 1: Randomizing the Orthology Mapping. (e) Since k-means is non-deterministic, to ensure robustness we performed 50 runs of the algorithm recording the fraction of times each gene pair was co-clustered (including all genes from all species). This matrix was hierarchically clustered.

Figure 3

Cluster structure and dynamics. (a) Each of the 17 clusters appears as a bubble containing up to three colored nodes whose sizes represent the number of genes contributed by each species. Edge thickness denotes the percent of gene orthology shared within or between clusters, measured using the size of the intersection divided by the size of the union of the sample sets. Only significant edges (P < 0.01) are shown. Several clusters show conserved orthology but not dynamics (for example, cluster 10 Sc, Cg with cluster 15 Kl). Note that clusters were ordered to minimize orthology edge crossings. (b) Expression dynamics of the 17 soft clusters over time following fluconazole exposure. Separate plots for each species can be found in Additional file 1. The width of each band corresponds to ± one standard deviation about the mean. A selection of enriched GO terms are shown for different clusters; see Figure S11 in Additional file 1 for full GO enrichment results. The number of genes for each species in each cluster is also shown.

Figure 4

Pathway expression conservation and divergence. (a) Top conserved and (b) diverged pathway responses as revealed by the soft clustering approach. Each pathway is represented by a pie with four slices - green, yellow, red, and black - denoting the percentage of orthologs in that pathway for which all three species co-clustered, two species co-clustered, no two species co-clustered, and no species' orthologs were differentially expressed, respectively. Pathways were defined using GO biological process annotations. (c) Schematic of ergosterol biosynthesis, the most conserved pathway response. Interestingly, this pathway includes isoprenoid biosynthesis, for which the response was one of the most divergent. (d) mRNA expression responses of ergosterol pathway genes are shown in order of occurrence in the pathway. Expression levels of genes 3 to 8 (boxed, and red) corresponding to isoprenoid biosynthesis are strikingly divergent. The fluconazole target Erg11 is boxed. (e) Hierarchically clustered mRNA expression responses of methionine biosynthesis genes show extensive divergence across species. Grey expression values denote a gene for which the species lacks an ortholog.

Figure 5

Divergence in transporter usage. Cross-species expression profiles of (a) ATP-binding cassette (ABC) and (b) major facilitator superfamily (MFS) transporters are shown. Grey expression values denote a gene for which the species lacks an ortholog. (c) Change in cell density with addition of exogenous ergosterol at the fluconazole 50% inhibitory concentration across different mutant backgrounds. Sc.bpt1Δ is a gene knockout unrelated to fluconazole response and is included as a control. Error bars indicate one standard deviation. (d) Model for differential usage of transporters among Sc, Cg, and Kl.