The evolution of drug resistance in clinical isolates of Candida albicans

Abstract

Candida albicans is both a member of the healthy human microbiome and a major pathogen in immunocompromised individuals. Infections are typically treated with azole inhibitors of ergosterol biosynthesis often leading to drug resistance. Studies in clinical isolates have implicated multiple mechanisms in resistance, but have focused on large-scale aberrations or candidate genes, and do not comprehensively chart the genetic basis of adaptation. Here, we leveraged next-generation sequencing to analyze 43 isolates from 11 oral candidiasis patients. We detected newly selected mutations, including single-nucleotide polymorphisms (SNPs), copy-number variations and loss-of-heterozygosity (LOH) events. LOH events were commonly associated with acquired resistance, and SNPs in 240 genes may be related to host adaptation. Conversely, most aneuploidies were transient and did not correlate with drug resistance. Our analysis also shows that isolates also varied in adherence, filamentation, and virulence. Our work reveals new molecular mechanisms underlying the evolution of drug resistance and host adaptation.

Keywords: Candida albicans; drug resistance; evolution; evolutionary biology; genomics; infectious disease; microbiology; virulence.

Figure 1.. Overview of study design.

(A) Background, persistent, transient, recurrent, and driver mutations in patient time courses. Shown is a schematic illustration of the genomes of isolates (gray bars) from two patient time courses (Patient A and B, left and right panels, respectively), ordered from the first isolate (progenitor, top) to the last (evolved, bottom). Background mutations (purple) exist in the all isolates; persistent mutations (yellow) are not in the progenitor, but found in all subsequent isolates after their first occurrence; transient mutations (pink) are not in the progenitor and only in some later isolates; recurrently polymorphic genes contain persistent mutations that occur in the same gene in more than one patient (black box). LOH events were also evaluated for persistence (light teal bar). Driver mutations, where a new persistent homozygous allele appears (e.g., G/T > A/A), are annotated in association with persistent LOH events (dark teal) and independent of these events (not shown). Each of these can be associated with a change in phenotype, such as drug resistance (boxes, right). (B) Sampling in the context of de novo mutation and selection bottlenecks. Each strain is a single clone (circle) isolated from an evolving population (represented by a phylogenetic tree). The population evolves and undergoes selective sweeps (dashed lines), with phenotypic changes occurring during the course of infection and treatment (i.e., drug resistance, black: high, white: low; gray scale at bottom). Persistent mutations (yellow lightning bolt) have likely swept through the population, whereas transient mutations (pink lightning bolt) have not. (C) Sampling in the context of selection on existing variation. Selection acts to vary the frequency of different pre-existing genotypes in the population. Persistent mutations (yellow lightning bolt) have risen in the population to a frequency that they are repeatedly sampled (large circles) whereas transient mutations (pink lightning bolt) have not (small circle). DOI: http://dx.doi.org/10.7554/eLife.00662.004

Figure 1—figure supplement 1.. Analysis of discordant sites.

(A) Degree of concordance (Y axis) with Sequenom iPlex genotyping for 1973 SNP X strain combinations overall (leftmost red bar; 93.9%) and in each tested strain (X axis). (B) Shown are the classes of discordant sites by genotype as defined by Illumina (orange) or Sequenome (teal) (X axis) and the prevalence (Y axis) of that genotype call in Sequenom (blue) and Illumina (orange) based discordant calls. The most common discrepancies arose when Sequenom typing classified a site as homozygous, but Illumina sequencing identified it as heterozygous. (C–G) Comparison on distributions of quality features between concordant (blue bars) and discordant (green bars) sites: (C) depth of coverage, (D) RMS Mapping Quality (MQ) score, (E) PHRED scaled quality score for each base call, shown as log-normalized ‘QUAL’ scores, (F) quality by depth (QD) score for each variant site, and (G) the allele balance ratio (AB Score) for each variant site. DOI: http://dx.doi.org/10.7554/eLife.00662.005

Figure 2.. Most isolates from the same patient are clonal.

(A) Two possible models of infection may underlie serial isolates. In the ‘clonal model’ (top) each subsequent sample (circle) is related to the other isolates. In the non-clonal model (bottom) isolates in a series are un-related. (B) The phylogenetic relationship of the isolates (black) from 11 patients (blue) was inferred based on 201,793 informative SNP positions using maximum parsimony in PAUP*. Isolates from the same patient separated by a branch distance greater than 20,000 were considered non-clonal (3281, 2823, 3184, 1691, red). Most nodes were supported by 100% of 1000 bootstrap replicates (indicated by *), expect as indicated (in gray). Clade identifiers were included as appropriate. DOI: http://dx.doi.org/10.7554/eLife.00662.006

Figure 2—figure supplement 1.. SNP heterozygosity profiles for each strain.

The heterozygosity profiles shows, in chromosomal order (top), each variant locus that exists in at least one strain in the series (white is a heterozygous SNP, blue is homozygous for the SC5314 allele, red is a homozygous SNP relative to SC5314). (A) Patient 1; (B) Patient 7; (C) Patient 9; (D) Patient 14; (E) Patient 15; (F) Patient 16; (G) Patient 30; (H) Patient 42; (I) Patient 43; (J) Patient 59. Only Patient 9 (C), Patient 16 (F), Patient 42 (H), and Patient 64 (not shown**) contain un-related isolates. ** Patient 64 contained an isolate (4380) whose genome aligned poorly to the C. albicans reference, but aligned well to C. dubliniensis. DOI: http://dx.doi.org/10.7554/eLife.00662.007

Figure 3.. LOH events were often persistent while aneuploidies were often transient.

For each time series shown are the genomes of all isolates (rows) from a patient, ordered from the first isolate (progenitor, top) to the last (evolved, bottom). Boxes on right indicate the MIC of the respective strain (black: high, white: low; gray scale at bottom). Persistent LOHs: blue, transient LOHs: pink; trisomies (all transient): green. The sequence coverage along each chromosome is indicated by black tickmarks. (A) Patient 1 has four LOH events, each coinciding with an increase in MIC (gray scale boxes, right). One LOH is transient (isolate 2, chromosome R, pink) and three are persistent (isolate 3, chromosome 3; isolate 13, chromosome 5; and isolate 16, chromosome 5, blue). The ploidy changes (isolates 6, 8, 13) are all transient. (B) Patient 7 has one LOH event (isolate 2307, chromosome 3, blue) which coincides with an increase in MIC. DOI: http://dx.doi.org/10.7554/eLife.00662.010

Figure 4.. Persistent and transient LOH and aneuploidies.

For each time series shown are the genomes of all isolates (rows), ordered from the first isolate (progenitor, top) to the last (evolved, bottom). Boxes on right indicate the MIC of the respective strain (black: high, white: low, gray scale at bottom). Persistent LOHs: light blue, transient LOHs: pink; trisomies (all transient): green. The coverage along each chromosome is indicated by black tickmarks. (A) Patient 9; (B) Patient 14; (C) Patient 15; (D) Patient 16; (E) Patient 30; (F) Patient 42; (G) Patient 43; (H) Patient 59. Several LOHs are recurrent (right arm of chromosome 3, left arm of chromosome 5, and chromosome 1). Please note: data in Figure 3—source data 1 also applies to this figure. DOI: http://dx.doi.org/10.7554/eLife.00662.012

Figure 5.. Co-occurrence of nonsynonymous substitutions across isolates reveals functional clusters.

(A) For each of the recurrently mutated 240 genes (genes in which nonsynonymous persistent SNPs appear in more than three patients and are not within an LOH region), we constructed a patient-by-gene binary vector. We clustered the resulting patient-by-gene matrix using NMF clustering to reveal five coherent clusters (correlation matrix of the clusters left; red: positive correlation; blue: negative correlation; white: no correlation). (B) Co-occurrence clusters. For the genes in each cluster (rows), shown are their mutated occurrences in each patient (columns); green: gene is persistently mutated in patient, white: no persistent mutation, yellow circle: driver mutation. Functional enrichment of clusters was revealed using gene ontology, and genes matching the enriched cluster function are bolded. We have overlaid recurrent driver mutations (e.g., G/T > A/A) (n = 17) occurring outside of LOH regions (yellow circle, green box) and inside LOH regions (yellow circle, white box). DOI: http://dx.doi.org/10.7554/eLife.00662.014

Figure 5—figure supplement 1.. Co-occurrence of nonsynonymous SNPs occurring in conjunction with a shift in MIC.

(A) For each of the 166 recurrently mutated genes associated with a change in MIC, we constructed a patient-by-gene binary vector. We clustered the resulting patient by gene matrix using NMF clustering to reveal 5 coherent clusters (correlation matrix of the clusters left; red: positive correlation; blue: negative correlation; white: no correlation). (B) Co-occurrence clusters. For the genes in each cluster (rows), shown are their mutated occurrences in each patient (columns); green: gene is persistently mutated in patient, white: no persistent mutation. Functional enrichment of clusters was revealed using gene ontology, and genes matching the enriched cluster function are bolded. DOI: http://dx.doi.org/10.7554/eLife.00662.015

Figure 6.. Filamentation, adhesion and virulence increase concurrently with fitness.

For each pair of consecutive isolates (green preceding blue), shown are the fitness, adhesion, filamentation, and virulence in a worm model of infection (each described in ‘Materials and methods’). A subset of fitness values are duplicated from Figure 7A, with selection coefficient (s) shown on the Y-axis. A subset of adhesion values are plotted from Figure 6—source data 1, with Abs590 nm on the Y-axis. A subset of images showing filamentation on spider media are shown, with the full set found in Figure 6—figure supplement 1. For virulence, shown are Kaplan–Meier plots of survival rates from C.elegans infection with the specified C. albicans isolates (‘Materials and methods’). For each isolate pair, significant changes in virulence were observed between the two isolates (in all cases, p < 0.001, log-rank test), with three of the four evolved isolates being more virulent than their corresponding progenitor. (A) Patient 30 isolates 5106 and 5108; (B) Patient 43 isolates 1649 and 3034; (C) Patient 1 isolates 12 and 13; (D) Patient 59 isolates 3917 and 4617. DOI: http://dx.doi.org/10.7554/eLife.00662.017

Figure 6—figure supplement 1.. Filamentation increases in many patient series.

For several patient series, shown are the filamentation assay results after 7 days of growth on Spider Media (‘Materials and Methods’). These data, a subset of which is shown in Figure 6, demonstrate the heterogeneity seen between strains, as well as the general trend for filamentation to increase over time. DOI: http://dx.doi.org/10.7554/eLife.00662.018

Figure 7.. Emergence of increased drug resistance often coincides with reduction in fitness in the absence of drug, but an increase in the presence of drug.

(A) For each patient (panel) shown is the fitness (‘Materials and methods’) of each strain (Y axis, mean ± STDV), ordered from the progenitor to evolved isolates (left to right, X axis). Fitness is calculated relative to an ENO1::YFP SC5314 reference isolate. The MIC of each strain is shown in the gray boxes on top (white: low; black: high, color bar at bottom). Green: isolates with aneuploidies; Blue: euploid isolates. (B) Shown is the mean difference between fitness in the absence and presence of drug (Y axis, error bars are ± STDV; n > 3) for isolates (X axis) that showed a decrease in fitness (Figure 7A) in the absence of drug concomitant with an increase in MIC (asterisks), and flanking isolates in Patient 1 and 59 (ordered from the progenitor to evolved isolates, left to right, X axis). The difference in fitness is calculated as the difference in selection coefficient (s, Y axis) between matching competition experiments in RPMI and those in RPMI with one half the MIC for fluconazole (Table 1) for each isolate tested (X axis). Negative values indicate that the strain had higher fitness in the presence of fluconazole vs assays without fluconazole. For each assay, the fluconazole-resistant isolate 4639 ENO1::YFP was used as the reference strain. DOI: http://dx.doi.org/10.7554/eLife.00662.020

Figure 8.. Deletion mutants of recurrently mutated genes reveal changes in relative fitness.

Shown is the fitness (‘Materials and methods’) for each deletion mutant strain and the corresponding wild-type strain (Y axis, mean ± STDV). The wild-type parental strain (SN250) is on the far left (red bar and dashed line). Fitness is calculated relative to an ENO1::YFP SC5314 reference isolate. Locus names are given for the mutant isolates (X axis). Asterisks denote statistical significance (* < 0.05, ** < 0.01, *** < 0.001, **** < 0.0001) by one-way ANOVA with Holm–Sidak correction for multiple comparisons. DOI: http://dx.doi.org/10.7554/eLife.00662.021

Author response image 1.

Null distributions of (a) the expected number of persistent, recurrent non synonymous SNPs outside of LOH regions, and (b) persistent non-synonymous SNPs associated with MIC changes outside of LOH regions.

Author response image 2.

The results of three separate MCMC analyses based on differing priors for mutation rate are shown. This analysis is based only on SNPs from isolates of Patient 1.