Leveraging single-cell genomics to expand the fungal tree of life

Abstract

Environmental DNA surveys reveal that most fungal diversity represents uncultured species. We sequenced the genomes of eight uncultured species across the fungal tree of life using a new single-cell genomics pipeline. We show that, despite a large variation in genome and gene space recovery from each single amplified genome (SAG), ≥90% can be recovered by combining multiple SAGs. SAGs provide robust placement for early-diverging lineages and infer a diploid ancestor of fungi. Early-diverging fungi share metabolic deficiencies and show unique gene expansions correlated with parasitism and unculturability. Single-cell genomics holds great promise in exploring fungal diversity, life cycles and metabolic potential.

Fig. 1. Phylogenetic tree with assembly and annotation statistics.

a, RAxML tree constructed from MCL clustering from across the fungi and deep-branching eukaryote outgroups. Support values are based on 1,000 bootstrap replicates. Bootstrap values (<100%) are indicated on branches. Fungal species are annotated with simplified lifestyle icons. Species sequenced using single-cell methods are coloured in magenta. Enrichment cultures for R. allomycis and C. protostelioides and a related species of M. bicuspidata are denoted with an asterisk. b, Estimated completeness (CEGMA (%)), assembly size and gene model statistics for fungal genomes. Metrics for the most complete individual libraries are indicated on the respective graphs, either as ‘single’ (referring to 1-cell libraries), ‘multiple’ (referring to 100-cell libraries in all cases except for T. sphaerospora, in which it refers to a 50-cell library) or ‘CoASM’ (referring to co-assemblies). ‘Isolate’ refers to assemblies derived from genomic material sequenced using non-single-cell methods.

Fig. 2. CEGMA correlation with assembly size and read mapping.

a, Scatter plot showing the correlation between estimated completeness and assembly size for individual libraries and co-assemblies for a given species. b, Scatter plot showing the correlation between estimated completeness and the percentage of reads from a given library mapped to the co-assembly for that organism. Individual library sizes are reflected in the size of the points.

Fig. 3. Ploidy analysis.

a–g, Analysis of diversity in population sequencing (a–d) and single-cell sequencing (e–g). R. allomycis (a,b) and C. protostelioides (c,d) population variability is shown in the k-mer distribution plot (a,c) and in the allele frequency spectra (b,d). In panel e, recovery of isolate SNPs from a single cell (1) and multiple cells (1–100) are shown. In panel f, the amount of heterozygous SNPs (allele frequencies above 25%) in single-cell and multiple-cell libraries is indicative of higher ploidy. In panel g, high numbers of shared heterozygous SNPs in non-haploid organisms are shown. For e–g, n = 4 biologically independent samples, in which each dot corresponds to one library. The horizontal black lines indicate the mean values.

Fig. 4. Core metabolism and individual pathways.

a, Metabolic pathway map with enzymes found in ~70% of ‘free-living’ fungi and absent in at least one of the single-cell fungi shown as a gradient: darker red shades indicate more single-cell species missing the same enzyme in a given pathway. b–f, Certain pathways with a high degree of such common losses: assimilatory sulfate reduction (sulfate to sulfide) (b), thiamine biosynthesis (HMP and TZE to thiamine phosphate) (c), spermidine synthesis (S-adenosylmethionine to spermidine) (d), biotin metabolism (dethiobiotin to biotin (and biotinyl-5′-AMP)) (e), and citrate synthesis (citrate to/from oxaloacetate and acetyl-CoA) (f). APS, adenylyl sulfate; HMP, 4-amino-5-hydroxymethyl-2-methylpyrimidine; HMP-P, HMP phosphate; HMP-PP, HMP diphosphate; PAPS, phosphoadenylyl sulfate; TZE, 5-(2-hydroxyethyl)-4-methylthiazole; TZP, TZE phosphate. The coloured boxes for b–f indicate the presence (solid) or the absence (shaded) of a given enzyme in a genome: Ra, R. allomycis (red); Bh, B. helicus (blue); Cp, C. protostelioides (green); Dc, D. cristalligena (purple); Pc, P. cylindrospora (orange); Ts, T. sphaerospora (yellow); Sp, S. pseudoplumigaleata (brown); Mb, M. bicuspidata (pink).

Fig. 5. CAZY-to-protease ratios.

CAZymes and proteases (defined as in the CAZY database and the MEROPS database) were predicted for each species in the taxon list. Total counts of each were plotted to describe the overall ratio. Correlation was determined using the linear model function (lm) in R, which implements QR decomposition. A trend line for saprobe species (n = 19) was omitted owing to a low r² value. Similarly, a trend line for mycorrhizal species was omitted owing to a low species count (n = 2).

Fig. 6. Subtilases.

a, Network figure highlighting the relationships between fungal subtilases. The distance reflects sequence similarity. Groups marked as Li-1 through to Li-4 refer to clusters previously reported by Li et al., from which the data set originated. Bh, B. helicus; Cp, C. protostelioides; Dc, D. cristalligena; Mb, M. bicuspidata; Pc, P. cylindrospora; Ra, R. allomycis; Sp, S. pseudoplumigaleata; Ts, T. sphaerospora. b, Subtilase domain architecture in Zoopagomycota genomes, for proteins in a given homologous cluster. Mycoparasites sequenced using single-cell methods from this study (labelled with *) have generally more subtilase proteins. For S. pseudoplumigaleata and T. sphaerospora, these are predominately single-domain (PF00082) proteins. Domain architecture colours: blue, PF00082; light grey, PF02225; dark grey, PF06280. Lifestyles are indicated with simplified icons. Species abbreviations are consistent with panel a and also include: Cc, Conidiobolus coronatus; Cr, C. reversa; Lp, L. pennispora.