101 Dothideomycetes genomes: A test case for predicting lifestyles and emergence of pathogens.

Dothideomycetes is the largest class of kingdom Fungi and comprises an incredible diversity of lifestyles, many of which have evolved multiple times. Plant pathogens represent a major ecological niche of the class Dothideomycetes and they are known to infect most major food crops and feedstocks for biomass and biofuel production. Studying the ecology and evolution of Dothideomycetes has significant implications for our fundamental understanding of fungal evolution, their adaptation to stress and host specificity, and practical implications with regard to the effects of climate change and on the food, feed, and livestock elements of the agro-economy. In this study, we present the first large-scale, whole-genome comparison of 101 Dothideomycetes introducing 55 newly sequenced species. The availability of whole-genome data produced a high-confidence phylogeny leading to reclassification of 25 organisms, provided a clearer picture of the relationships among the various families, and indicated that pathogenicity evolved multiple times within this class. We also identified gene family expansions and contractions across the Dothideomycetes phylogeny linked to ecological niches providing insights into genome evolution and adaptation across this group. Using machine-learning methods we classified fungi into lifestyle classes with >95 % accuracy and identified a small number of gene families that positively correlated with these distinctions. This can become a valuable tool for genome-based prediction of species lifestyle, especially for rarely seen and poorly studied species.

Introduction

Dothideomycetes represent the largest and one of the most important classes of ascomycete fungi, currently encompassing more than 23 orders, 110 families, 1 261 genera and 19 000 species (Wijayawardene ). Numerous species remain undescribed, however, and many well-known taxa actually represent species complexes waiting to be resolved, suggesting that there are likely still thousands of species to be discovered (Lücking & Hawksworth 2018). Dothideomycetes are cosmopolitan fungi of diverse lifestyles, which can associate with a wide range of hosts/substrates. The class Dothideomycetes is well supported as sister to Arthroniomycetes, a class of 1 500-plus lichenized species (Wijayawardene ) within the subphylum Pezizomycotina, which includes most filamentous ascomycetes. Dothideomycetes includes two subclasses, Pleosporomycetidae and Dothideomycetidae, based on phylogeny and the presence or absence of pseudoparaphyses (Schoch et al., 2006, Schoch et al., 2009). Members of the order Capnodiales (Dothideomycetidae) seem to have adapted to (hemi-)biotrophy, while Pleosporales (Pleosporomycetidae) species tend to have adapted to necrotrophy. A previous study of the Dothideomycetes by Ohm showed that Capnodiales have fewer genes involved in carbohydrate degradation, proteolysis and secondary metabolism. Due to limited sampling, however, it is not known if this holds true for all orders in each subclass, or if further evolutionary specialization has taken place to coincide with their differing ecologies.

Members of Dothideomycetes are known to infect most major crops of mono- and dicotyledons, including species that serve as staple food crops, biomass, and biofuel production (Ohm et al., 2012, Condon et al., 2013). Well-known examples include southern corn leaf blight caused by Bipolaris maydis (Turgeon and Baker, 2007, Marin-Felix et al., 2017), the Sigatoka complex of banana by Pseudocercospora spp. (Chang ), and multiple wheat pathogens including Parastagonospora nodorum (Quaedvlieg ), Zymoseptoria tritici (Stukenbrock ), and Pyrenophora tritici-repentis (Ellwood ). Diseases of woody hosts caused by Dothideomycetes include Sphaerulina canker of poplar by Sphaerulina spp. (Foster ), stem cankers by species of Botryosphaeria, Diplodia and Neofusicoccum (Phillips ), apple scab by Venturia inaequalis (Marin-Felix , Citrus Black Spot by Phyllosticta citricarpa (Guarnaccia ), and Teratosphaeria leaf blight of Eucalyptus by Teratosphaeria spp. (Crous ), to name a few.

In addition to numerous plant-pathogenic species, Dothideomycetes also include an incredible diversity of non-plant-pathogenic life styles, many of which have evolved multiple times. Members of the class are distributed globally, and have important roles as foliar endophytes (Alternaria, Botryosphaeria, Plenodomus, Phyllosticta), ectophytes (Phaeothecoidiella, Schizothyrium, Uwebraunia), and epiphytes (sooty molds) (Capnodium, Scorias). They can grow on other fungi as mycoparasites (Dissoconium, Pseudoveronaea), on scale insects (Myriangium), tolerate extreme climatic conditions, and on or in a diversity of substrates, from soil to rocks (Cryomyces, Coniosporium, Friedmanniomyces). The rock-inhabiting fungi are found in several unrelated lineages in Dothideomycetes, and could represent an ancestral trait (Ruibal ). Several genera are also found in fresh water (Jahnula, Massiosphaeria) or salt water (Halomassarina, Morosphaeriaceae) environments, where they occur in the intertidal zones. A small number of species are lichenised and Dothideomycetes appears to have the highest number of transitions to the lichenised lifestyle within a class (Nelsen ), as it occurs in three orders. Dothideomycetes also harbour several examples of mycoparasitism of lichens, especially in the Phoma complex (Lawrey et al., 2012, Valenzuela-Lopez et al., 2018).

Some taxa have been shown to be mycorrhizal, such as Cenococcum geophilum (Peter ), which is one of the more commonly encountered ectomycorrhizae, especially in more xeric conditions. Dothideomycetes saprobes are commonly found in plant litter, or in more extreme ecological niches. For example, the black yeasts are highly tolerant to extreme temperatures, solar radiation and desiccation. Some species are also tolerant to ethanol vapour, such as Baudoinia spp. that grow in the vicinity of whiskey distilleries, or Zasmidium cellare, growing in wine cellars (Goodwin et al., 2016, Videira et al., 2017). Several members are also prominent pathogens of humans and animals, such as Hortaea werneckii causing tinea nigra of skin, Piedraia hortae on hair (De Hoog ), or are allergenic, such as species of Alternaria and Cladosporium (Bensch ).

Understanding the ecology and evolutionary biology of Dothideomycetes is therefore fundamental to interpreting fungal evolution and adaptation to abiotic and biotic stress and host specificity. With a phylogenetically interspersed mix of several lifestyles, including a large number of plant pathogens responsible for economic losses, a major outstanding question is whether there are identifiable genomic signatures that correlate with diverse dothideomycete ecological niches. Here, we test this hypothesis by first generating a robust phylogeny (which additionally allowed us to reclassify 25 Dothideomycetes) and predicting the ancestral state of Dothideomycetes, then explore relationships between genomic content and lifestyle. Using this approach, we were able to identify several gene family expansions and contractions consistent with changes in lifestyle across the Dothideomycetes. Furthermore, our large data set enabled the use of machine learning techniques and demonstrate that high quality classification of species can be achieved using a small number of easily identifiable genes. With increasing genomic data becoming available, the ideas presented here suggest that machine-learning techniques can be used to gain insights into ecological speciation and provides clues to several genes that may play a role in lifestyle switches across the fungal tree of life.

Materials and methods

Genome sequencing, assembly and annotation

We used 101 Dothideomycetes genomes and eight other Pezizomycotina genomes as outgroups representing the Pezizomycetes, Eurotiomycetes, Leotiomycetes and Sordariomycetes (Supplementary Table 1). The 55 newly sequenced Dothideomycetes in this study were assembled and annotated between January 2012 and November 2016. The rapid progress in sequencing and assembly techniques during this period is reflected in the various techniques used to assemble these genomes. Two primary sequencing technologies were used, namely Illumina for earlier genomes and Pacific Biosciences more recently. For the oldest, Cercospora zeae-maydis, the consensus assembly was from a mixture of Roche (454) pyrosequencing, Sanger fosmids, and Illumina data using Velvet (Zerbino & Birney 2008). The Illumina-based assemblies fall into two main classes: 1) minimal draft assemblies from single 270-bp insert size libraries using 2 × 150-bp paired-end reads assembled first with Velvet to create an in-silico long mate-pair library with insert sizes of 3 000 +/- 300 bp, which was then assembled together with the original Illumina library using AllPaths-LG (Gnerre ); and 2) standard draft Illumina assembles where 2 × 100-bp paired-end reads from a 4-kb-long mate-pair library was used in addition to the reads from the 270-bp insert size libraries and assembled using AllPaths-LG. For PacBio sequencing, filtered subread data were assembled together with Falcon (https://github.com/PacificBiosciences/FALCON), improved with finisherSC (Lam ), and polished with Quiver (https://github.com/PacificBiosciences/GenomicConsensus). A summary of technologies used and assembly methods is shown in Supplementary Table 5. Genomes were annotated using the JGI Annotation pipeline and made available via the JGI fungal genome portal MycoCosm (Grigoriev ) https://mycocosm.jgi.doe.gov. An overview of the JGI Fungal Annotation Pipeline is available at https://mycocosm.jgi.doe.gov/programs/fungi/FungalGenomeAnnotationSOP.pdf. Functional annotation was performed as in the overview and has also been previously described (Ohm ). Previously published genomes were functionally re-annotated using the JGI annotation pipeline, but no changes to gene models were made.

Phylogenetic analyses

The phylogenetic tree is based on 738 single-copy orthologs identified in all 109 genomes (101 Dothideomycetes + eight outgroup) using OrthoMCL v. 2.0.9 (Li ). The genes were aligned using Mafft v. 7.123b (Katoh, 2002, Li et al., 2003) using the --auto option, concatenated and filtered using Gblocks v. 0.91b (Castresana 2000) using block options -b4 = 5 -b5 = h. To determine the optimal substitution model for tree building, we used PartitionFinder v2.1.1 (Lanfear ) and analysed the alignment in blocks of 300, 500 and 1 000 bp with options --no-ml-tree and --rcluster-max 100. In all cases, we obtained aamodelpr = fixed(wag) and rates = invgamma as optimal parameters. The alignment was used as input to RAxML v. 8.2.2 (Stamatakis 2014) under the GAMMA model of rate heterogeneity, ML estimate of alpha-parameter, and WAG substitution matrix.

To infer ancestral ecological character states, ecologies of each species was coded based on known ecologies states and using FUNGuild (Nguyen ). Known ecological states were coded either as multistate: acidophile, aquatic, human-associated, insect pathogen, lichen, mycoparasite, mycorrhiza, plant pathogen, and saprobe, or as binary: plant pathogen and non-plant pathogen. Ancestral ecological character states were reconstructed in Mesquite (http://www.mesquiteproject.org). Multistate ecological character states were analysed using an unordered parsimony model, while binary ecological character states were analysed using the Mk1 (Markov k-state 1 parameter model) likelihood reconstruction model. Similar analysis was also performed using BayesTraits analysis using maximum likelihood (http://www.evolution.rdg.ac.uk/SoftwareMain.html).

Enrichment analysis

To identify functional annotation differences among the various lifestyles in the Dothideomycetes, the two-tailed Fisher exact test for functional annotation was performed using the scipy.stats package in Python with a p-value threshold of 0.05. In order to further limit false positives, annotations with fewer than 100 genes in background distribution were ignored.

Machine learning

We used the Support Vector Machine method as implemented in Python's scikit-learn library (http://scikit-learn.org) to identify the most informative features (e.g., OrthoMCL gene clusters, PFAM annotations) for differentiating between pathogen and saprobe genomes. The dataset was cleaned using the nearZeroVar function in R to reduce noise in the analysis. We used a subset of known pathogens and saprobes (Supplementary Table 4) to train and identify features that were differentially distributed between the pathogens and saprobes. Features were ranked on ability to predict using C-Support Vector Classification within scikit-learn (sklearn.svm.SVC). Strength of prediction was validated using sklearn.cross_validation and jackknife, and best features were selected using sklearn.feature_selection.SelectKBest. Features that were highly correlated (>0.8) to lifestyle were used in all combinations to improve prediction accuracy on the training set. The best-performing combinations of features were used on the entire study set to assess the accuracy of predictions.

Results and discussion

Whole-genome data enables reclassification of 25 species

Our study set included 55 newly sequenced Dothideomycetes genomes and 45 previously published genomes (Supplementary Table 1) in this class. Genomes were sequenced and assembled as described in Methods, then annotated using the JGI Annotation pipeline (Grigoriev ). Using these 101 Dothideomycetes and 8 outgroups, we were able to generate a high confidence whole-genome-based phylogenetic tree (Fig. 1). The Gblocks filtered alignment of 738 single-copy orthologs consisted of 238 889 unique sites. The species tree provided strong support for the two subclasses, Pleosporomycetidae and Dothideomycetidae, proposed by Schoch to accommodate the pseudoparaphysate Pleosporales (order Pleosporales) and the aparaphysate Dothideales (orders Dothideales, Capnodiales and Myriangiales) (Fig. 1).

Fig. 1

A phylogenetic tree of the Dothideomycetes and eight outgroups using whole-genome data. The two sub-classes in the Dothideomycetes (Pleosporomycetidae and the Dothideomycetidae) are well resolved and expanded in Fig. 2, Fig. 3.

Fig. 2

A phylogenetic tree of 76 species from the Pleosporomycetidae used in this study. All bootstrap values are 100 % except for those shown. Well resolved orders with multiple species are indicated by brackets on the right. The three icons left of species names represent lifestyle classification based on organism data, FunGuild classification and Machine Learning predictions, respectively. *Represents newly sequenced species for this study. Arrows point to newly reclassified species.

Fig. 3

A phylogenetic tree of 25 species from the Dothideomycetidae used in this study. All bootstrap support values are 100 percent. Well resolved orders with multiple species are indicated by brackets on the right. The three icons left of species names represent lifestyle classification based on organism data, FunGuild classification and Machine Learning predictions, respectively. *Represents newly sequenced species for this study. Arrows point to newly reclassified species.

In the detailed phylogenetic analyses (Fig. 2, Fig. 3), the Pleosporomycetidae is further expanded to include Asterinales, Coniosporiales, Eremomycetales, Hysteriales, Lineolatales, Lizoniales, Microthyriales, Mytilinidiales, Patellariales, Phaeotrichales, and Venturiales, while the Trypetheliales is added to the Dothideomycetidae, with all branches having 99 % or higher bootstrap support. The original division proposed by Schoch still applies, with the newly placed orders listed above exhibiting these characters, except for the Pleosporomycetidae genera Aulographum (Asterinales), Eremomyces (Eremomycetales), and Microthyrium (Microthyriales) that lack pseudoparaphyses. The latter examples are, however, embedded in a subclade exhibiting pseudoparaphyses and represent secondary losses. All branches within the Dothideomycetidae tree had 100 % bootstrap support (Fig. 3) and indicated that Acidomyces richmondensis and Hortaea acidophila both belong in the family Teratosphaeriaceae of the Capnodiales. Based on this phylogenetic analysis, we have proposed new classifications for 25 species (Supplementary Table 2, Supplementary text) including four new orders (Aulographales, Coniosporiales, Eremomycetales and Lineolatales) and three new families (Coniosporiaceae, Lineolataceae and Rhizodiscinaceae).

Saprobes show larger genomes and proteomes than plant pathogens

The Dothideomycetes genomes varied more than tenfold in size (from <17 Mbp to >177 Mbp) and encoded from 7 572 to 21 730 protein-coding genes, with neither genome size nor gene count correlated with the lifestyle (Supplementary Figs S1–S3) or phylogenetic position of the species (Fig. 4). However, pathogens generally showed smaller genomes as seen previously, including the recently sequenced species of Geosmithia where pathogenic species had smaller genomes than non-pathogenic Geosmithia species (Schuelke ). Piedraia hortae, a non-pathogenic opportunist that infects human hair causing black hardened nodules, had the smallest genome assembly at 16.95 Mbp while the two root-associated fungi Zopfia rhizophila, an opportunistic pathogen, and Cenococcum geophilum, an ectomycorrhizal species, had the largest genomes at 152.8 Mbp and 177.6 Mbp, respectively. Repeat content also varied greatly between the genomes, but this is influenced by the different sequencing and assembly technologies. As expected, smaller genomes had fewer repeats compared to larger assembles (R2 = 0.89). The two largest assemblies were C. geophilum and Z. rhizophila, which consisted of 65.4 % and 49.6 % repeat content, respectively. This made the non-repeat genome space of Z. rhizophila (77.1 Mbp) slightly larger than that of C. geophilum (61.5 Mbp). The median assembly was just under 36.5 Mbp or about 33.1 Mbp without repeats. The sizes of the non-repeat genome spaces also did not show any correlation related to lifestyle (Supplementary Fig. S2) or phylogenetic position.

Fig. 4

Genome sizes and gene content of 101 Dothideomycetes. A. Phylogenetic tree of the Dothideomycetes showing plant pathogens in blue and presumed saprobes in green. B. Genome sizes showing proportion of repeat content in these genomes. C. Gene content showing predicted proteome sizes and content. Core genes are found in all Dothideomycetes in this study while common genes are found in more than one genome. “Unique to the Pleosporomycetidae” and “Unique to the Dothideomycetidae” are found in all genomes in the Pleosporomycetidae and the Dothideomycetidae, respectively.

The number of predicted genes was related to non-repeat genome size (R2 = 0.69), but not to functional annotation, lifestyle or phylogenetic position, although saprobes generally had a slightly higher gene count than pathogens (Supplementary Fig. S3). Piedraia hortae, which had the smallest genome, also had the fewest genes at 7 572, and was the only one in our study group which had fewer than 10 000 genes. The average Dothideomycetes genome had about 12 750 genes of which about 2 100 had homologs in all Dothideomycetes (core genes) in our analysis and an additional 8 700 had homologs in at least one other Dothideomycetes. Even at our study depth, we were able to identify about 1 950 gene clusters unique to each newly sequenced species (Fig. 5). As expected, the vast majority of core genes had some functional annotation, with over 91 % encoding a PFAM domain. On the other hand, only 8.6 % of “unique genes” encoded a PFAM domain showing the vast gap in our knowledge about key genes that could account for species-specific adaptations. While the average unique gene was smaller than the average core gene (394 vs 462 codons), it is well above the 100–150 aa threshold below which gene prediction is often fraught with errors (Haridas ).

Fig. 5

Identifying core, group-specific (common to two or more genomes), and unique genes in OrthoMCL data. Error bars represent 100 random subsamples from the study group.

We annotated the gene models to identify PFAMs, transporters, proteases, KEGG Orthology (KO), Carbohydrate-Active Enzymes (CAZy), the secretome and small secreted proteins. We were able to assign some form of functional annotation to about 79 % (average; 67 % median) of the genes in these Dothideomycetes genomes (Fig. 4C). The most common was PFAM assignment in over 73 % (average; 62 % median) of the genes. As expected, the number of genes with PFAM annotations was highly correlated with the number of genes (R2 = 0.85). Other functional annotations are also correlated with the number of genes (proteases R2 = 0.66, transporters R2 = 0.67 and secretome R2 = 0.68). In general, saprobes had a higher number of genes, and consequently proportionally higher numbers of PFAMs, proteases and transporters. However, small, secreted cysteine-rich proteins (R2 = 0.46) showed several outliers with the plant pathogen Venturia inaequalis, the causal agent of apple scab disease, and Zymoseptoria ardabiliae, a hemibiotrophic pathogen of certain grasses, showing the highest relative numbers. We also found that saprobes had a proportionally higher number of CAZy enzymes compared to pathogens. Ohm reported a reduction in genes involved in carbohydrate degeneration in the Capnodiales, but this group is primarily composed of plant pathogens even at our increased sampling of 16 Capnodiales species compared to seven in Ohm . Ancestral state reconstruction (see below) suggests that the switch to parasitism occurred early in the Capnodiales (Supplementary Fig. S4, S5) and the reduction in CAZy enzymes may be concurrent with this switch.

The ancestral dothideomycete was a saprophyte

Parsimony reconstruction of multistate character coding using Mesquite inferred the ancestral ecological character state of Dothideomycetes as saprobe (Supplementary Fig. S4), and likelihood reconstruction of binary character coding significantly favoured the common ancestor of Dothideomycetes as non-plant pathogen (NPP = 0.9987) as compared to plant pathogen (PP = 0.0013) (node 102: Supplementary Fig. S5, Supplementary Table T1). A similar result was obtained using BayesTraits (saprobe: 0.999698, Pathogen: 0.000277, other: 0.000025). Plant pathogenicity is supported as evolving six to eight times within Dothideomycetes with at least four major transitions from NPP to PP including in the most recent common ancestor (MRCA) of Mycosphaerellaceae of Dothideomycetidae (node 158: NPP = 0.1713; PP = 0.8287), and the MRCAs of Venturia (node 96: NPP = 0.0033, PP = 0.9967), Botryosphaeriales (node 172: NPP = 0.0770, PP = 0.9230), and the Setomelanomma-Bipolaris clade (node 67: NPP = 0.0177, PP = 0.9823) of Pleosporomycetidae. The increased sampling depth in the Pleosporales from nine in Ohm to 49 in this study showed that the transition from saprotrophy occurred late, and is restricted to the Pleosporineae (Supplementary Fig. S4, S5). The transition from NPP to PP at the Setomelanomma-Bipolaris clade could possibly have occurred one node earlier (MRCA of Dothidotthia-Bipolaris clade at node 23), but the parsimony reconstruction was equivocal and the likelihood reconstructions showed increased proportional support for NPP (NPP = 0.1815, PP = 0.8185). Within Dothideomycetidae, acidophiles (Acidomyces richmondensis, Hortaea acidophila) and human-associated (Piedraia hortae) species were restricted to a subclade of Capnodiales, and lichen (Viridothelium virens, Trypetheliales) and insect pathogen (Myriangium duriaei, Myriangiales) each displayed unique origins. Within Pleosporomycetidae, mycoparasitic species Ampelomyces quisqualis and Pseudovirgaria hyperparasitica arose independently; the seven aquatic species sampled represented six independent transitions to aquatic environments, but ectomycorrhizae had a single origin in Cenococcocum geophilum. These results are consistent with findings of Schoch , who inferred numerous transitions from saprobic life histories to plant pathogens and lichens, and multiple transitions from terrestrial to aquatic habitats.

Some gene family expansions are correlated with lifestyles

Most of the species in our study set of 101 Dothideomycetes were plant pathogens or presumed saprobes, but ecologies can be difficult to code unequivocally for some taxa. The saprobes sampled here also showed a wide variety of ecological niches such as growth in aquatic habitats (e.g., Amniculicola lignicola, Lineolata rhizophorae), on dung (e.g., Delitschia confertaspora, Eremomyces bilateralis, Sporormia fimentaria and Trichodelitschia bisporula), and at relatively high concentrations of alcohol in the air (Baudoinia panamericana and Zasmidium cellare). Other lifestyles in our study group included mycorrhizal (Cenococcum geophilum), human hair pathogen (Piedraia hortae), mycoparasites (Ampelomyces quisqualis and Pseudovirgaria hyperparasitica), and fungal partners of lichen symbiosis (Viridothelium virens). As summarized above, the various lifestyles are not grouped phylogenetically, but appear interspersed among different branches of the tree and most likely arose from saprobic ancestors. We used Fisher's exact test with Bonferroni correction to adjust p-values for multiple comparisons and identified gene family expansions and contractions in convergent lifestyle evolution on various branches of the Dothideomycetes tree (Supplementary data - Enrichment.xls).

Among the Dothideomycetes, Baudoinia panamericana and Zasmidium cellare are found in wine caves or cellars or proximal to distilleries, environments associated with significant levels of alcohol vapours. These genomes showed a >30 % increase (p value = 0.008) in transport genes related to the Major Facilitator SuperFamily (MFS). The ability to assimilate alcohol from vapours presumably requires transporters, but analysis of their genomes did not identify any unique genes or gene families. To identify specific MFS transporters that may be involved in alcohol transport, we looked at the distribution of OrthoMCL gene clusters containing MFS genes. Baudoinia panamericana and Z. cellare did not show an expansion in specific MFS clusters compared to other Dothideomycetes in this study. This could suggest that genes in particular MFS families could have been co-opted for alcohol transport, or this ability was independently gained in these two genomes, which would remain lineage-specific, and were therefore omitted in our analysis.

The two animal-associated fungi in our study group, Piedraia hortae, which forms dark, hardened nodules (hence piedra, stone) firmly attached to human hair in humid tropical environments, and Myriangium duriaei, which is associated with scale insects and is presumably parasitic, showed a 25 % reduction (p-value = 0.001) in MFS and a 16-fold increase (p-value = 0.0) in PF05922 (protease inhibitor). These differences are presumably due to the narrow ecological niche and relatively stable environment inhabited by these two fungi.

An analysis of the genomes of Pseudovirgaria hyperparasitica and Ampelomyces quisqualis showed a two-fold increase (p-value = 0.0) in PF01636 (Phosphotransferase enzyme family), a family of bacterial antibiotic resistance proteins. These mycoparasites may be using these genes to overcome antibiotics produced by their hosts during infection. However, this class of proteins also showed high numbers in insect-associated fungi like Tolypocladium inflatum (Bushley ) and Metarhizium robertsii (Gao ), and also some plant pathogens such as Neonectria ditissima (Gómez-Cortecero ). In addition, several well-known mycoparasites like Trichoderma gamsii (Baroncelli ), T. asperellum and T. harzianum (Druzhinina ) do not show an increase in PF01636, suggesting that they either use alternative strategies or the lack of sufficient sampling in our study set introduced a statistical false positive in this case. Our analysis also showed PF01636 enrichment in aquatic and saprophytic lineages, suggesting that this family has diverse functions and testing individual genes using experimental approaches for roles in mycoparasitism might yield valuable insights.

Several gene families are reduced in plant pathogens compared to saprobes

Plant pathogens and saprobes were the largest group of lifestyles in our study set. They split evenly into 42 pathogens and 42 presumed saprobes, with other lifestyles making up the rest of the study set. Among the functionally annotated genes (pfams, CAZy, etc), none were unique to either the pathogens or saprobes. Because about a third of the genes in these genomes have no functional annotation, we looked for OrthoMCL gene clusters that were unique to either the pathogens or saprobes. A cluster consisting of 35 genes from 26 pathogen genomes with no representative from other lifestyle classes was identified (Supplementary Table 6). This group consisted of small (<200 aa) proteins, most with signal peptides, but no other functional annotation. We deleted one of these genes in Bipolaris maydis strain C4, where only one copy was identified, but the virulence phenotype on the maize host was unaltered when compared to inoculation with the wild type strain (Supplementary Fig. S7). Although we would not expect a single gene to have a large effect on pathogenicity of diverse species to their hosts, a gene conserved only in pathogens likely plays a quantitative role in pathogenicity and additional higher resolution virulence assays in a range of species may be warranted.

Fisher's exact test identified several gene families that were expanded or contracted in pathogens and saprobes. To reduce errors due to misclassification of species into various lifestyles, we used a curated list of 39 well-known saprobes and 39 plant pathogens (removing additional ambiguous lifestyle species) to identify enrichment of gene families in these two groups. We used a cutoff of p < 0.005 and at least 100 genes in the test subset to identify 14 GO annotations and 28 PFAM annotations that showed different enrichment between the saprobes and pathogens. Only one PFAM (PF01328, Peroxidase) and one GO (0032440, 2-alkenal reductase activity) were enriched in pathogens compared to saprobes. While we observe a strong association between presence of these functional annotations and lifestyle, further experimental work is necessary to elucidate their potential roles in pathogenicity.

Among the PFAM domains that were reduced in pathogens were several tetratricopeptide repeat (TPR) domains (Supplementary Fig. S6). The TPR structural motif was originally identified in the fission yeast Schizosaccharomyces pombe and appears in all branches of life in various unrelated protein families; it is linked to virulence in some human pathogens like Candida albicans (Kaneko ) and Porphyromonas gingivalis (Kondo ). The TPR motif has a role in protein-protein interactions and is found in proteins involved in the anaphase-promoting complex, NADPH oxidases, hsp90-binding proteins, transcription factors, protein kinase inhibitors, and peroxisomal and mitochondrial import proteins, among others (Das ). TPR proteins are also part of the plant hormone signal-transduction pathway (Schapire ). Our dataset showed that fungal saprobes had more than double the number of TPR proteins compared to plant pathogens. Much of this expansion was driven by root-associated saprobes (TPE 3, 4, 7, 8, 10 and 16), but TPR 12 was expanded (68 %) in the pathogens compared to saprobes (p-value = 0.0). It is interesting to speculate whether this reduction is due to a decrease in signalling-related TPR proteins in pathogenic fungi, where they may be affected by plant hormone signalling induced by the fungal attack. Other interesting groups of PFAMs depleted in plant pathogens are PF17111 (Fungal N-terminal domain of STAND proteins), PF05729 (NACHT domain) and PF14479 (involved in prion-inhibition and propagation). These domains have a role in heterokaryon incompatibility (Paoletti and Saupe, 2009, Greenwald et al., 2010, Daskalov et al., 2012). In Colletotrichum lindemuthianum, fusions between incompatible strains during colony initiation survived and formed heterokaryotic colonies from conidial anastomosis tubes (Ishikawa ). The authors hypothesized that such fusions allow asexual fungi to increase their genetic diversity and acquire new pathogenic traits. The reduction in PFAMs related to heterokaryon incompatibility suggests that other strategies to reduce heterokaryon incompatibility may exist in pathogenic fungi that may allow for an increase in fitness related to pathogenic traits.

Similarly, among the GO annotations, protein kinases (GO:0004672), and transcription factors (GO:0003700) were under-represented in pathogens compared to saprobes. In addition, GO:0006808 (regulation of nitrogen utilization), and GO:0008270 (Zn-binding) also showed a reduction in pathogens compared to saprobes. The availability of zinc and nitrogen is a critical factor in infection by fungi (Staats et al., 2013, Mur et al., 2017), which can be explained in part by the reduction of these genes in pathogens compared to saprobes. Moreover, several hosts have developed nutritional immunity strategies by withholding essential micronutrients from invading microbes. In such a scenario, an increase in uptake efficiency may be preferable to an increase in genes related to assimilation. An evolved reliance on uptake of plant amino acids would be energy sparing for the pathogens given the complexity and energy demands of amino acid biosynthesis pathways.

Machine learning distinguishes saprobes from plant pathogens

Our knowledge of lifestyles and ecologies of many species sequenced in this study are limited and we attempted to predict the major lifestyles based on genome sequencing alone. Therefore, we used FUNGuild (Nguyen ) to predict lifestyles of poorly studied taxa and compared this with our current understanding (Fig. 2, Fig. 3). In several cases, we found a discrepancy between what we expected and what was predicted by FUNGuild. For example, we sequenced Plenodomus tracheiphilus and Ophiobolus disseminans, which are known plant pathogens, but FUNGuild classifies them as probable saprobes. On the other hand, we assumed Massarina eburnea and Periconia macrospinosa are saprobic and endophytic respectively based on their known lifestyles, but are classified as pathogenic by FUNGuild (Supplementary Table 3). Since FUNGuild results are based on database queries, the prediction accuracy is affected by the quality of data present in the database. To reconcile these discrepancies, we used machine-learning methods to identify hidden correlations between lifestyle and gene content. Large data sets are critical to machine learning, and only plant pathogens and saprobes had sufficient sampling depth in our dataset to be adequately separated using these methods. We selected a training set of 34 well known pathogens and 24 saprobes across the tree to identify genes that could discriminate between pathogens and saprobes in this study. Ten PFAM annotations provided accuracy >0.8 in their ability to separate pathogens from saprobes in the machine learning analysis. The best PFAM discriminator was PF13350 (Protein tyrosine phosphatase), which showed an average of three versus four gene copies in pathogens versus saprobes. Other PFAMs showing reduction between pathogens and saprobes included PF06544 (DUF1115), PF00017 (SH2 domain), PF11700 (Vacuole effluxer Atg22 like), PF07522 (DNA repair metallo-beta-lactamase), PF05057 (DUF676), PF03537 (Glycoside hydrolase), PF01183 (Glycoside hydrolase family 25), PF00413 (Matrixin), and PF00076 (RNA recognition motif) (Supplementary data - Enrichment.xls). Of these, enrichment analysis using Fisher’s Exact Test had identified PF05057 (p = 0.03) and PF03537 (p = 0.0) as reduced in pathogens.

As seen in Fig. 4C, over a third of the genes in the Dothideomycetes genomes do not have functional annotation. On the other hand, OrthoMCL gene clusters (genes with orthologs) encompass 85 % of the gene space in these genomes (Fig. 5). Using the OrthoMCL gene clusters dataset (which ignores lineage-specific genes), 27 OrthoMCL gene clusters (presumed orthologs) produced accuracy scores >0.8 in their ability to distinguish between saprobes and pathogens. The best discriminator among the OrthoMCL gene clusters was cluster #9837 (C009837 in Supplementary Table 4) which was found in 80 % of the pathogens and less than 5 % of saprobes and produced a prediction accuracy of 0.89. This cluster consisted of small (<200 aa) proteins with no other functional annotation. To improve prediction accuracy, we used all 27 clusters in all possible combinations. Our tests indicated that high prediction accuracy (>95 %) could be achieved by using a small handful of gene clusters (Supplementary Fig. S8).

Among the 50 highest-ranked combinations, six clusters, none of which had any functional annotation, appeared more than 30 times. We used these six clusters to classify the Dothideomycetes into pathogens and saprobes as shown in Fig. 2, Fig. 3, Fig. 6. These six orthologous gene clusters were able to differentiate between pathogens and saprobes with >95 % accuracy in our study set. While the vectors underlying machine-learning remain obscure within the algorithm, the data shows that two of these gene clusters are primarily found in pathogens (PRG, pathogen related genes) while the other four are mostly in the presumed saprobes (SRG, saprobe related genes) (Supplementary Table 4). The data suggest that saprobes have SRG at least two times the number of PRG except for Macroventuria anomochaeta and Stagonospora sp. SRC1lsM3a. Macroventuria anomochaeta has both the PRG and none of the SRG and is strongly predicted to be a pathogen (Macan1 in Fig. 6). The sexual morph of Macroventuria anomochaeta was isolated from decayed canvas in a desert in South Africa and its ecology is presumed to be saprotrophic, but the underlying enzymatic mechanisms are unknown. Stagonospora sp. SRC1lsM3a (Stasp1 in Fig. 6) has two of the SRG and one PRG and is predicted to be a saprobe. While several Stagonospora species are plant pathogens, Stagonospora sp. SRC1lsM3a was isolated and studied for its role in the remediation of metal-polluted environments due to the oxidation of manganese compounds. Species-level determination of Stagonospora sp. SRC1lsM3a in phylogenetic analysis have been unsuccessful, and the absence of spores or reproductive structures in culture prevents morphological identification (Zeiner ) (https://mycocosm.jgi.doe.gov/Stasp1). Additional genomic data from closely related species are required to resolve this.

Fig. 6

Support Vector Machine (SVM)-based prediction of lifestyle based on 6 gene clusters showed a >95 % accuracy in correctly predicting plant pathogens (blue) vs saprobes (green). (Other lifestyles are indicated in black). Some poorly studied saprobes are predicted to be pathogens suggesting that these may be weak pathogens or have recently diverged from pathogenic ancestors. Short names used are shown in Supplementary Table 4.

Another organism isolated and studied for its ability to oxidize Mn compounds, Pyrenochaeta sp. DS3sAY3a (Pyrsp1), had two PRG and two SRG and therefore has a near equal probability of being classified as a pathogen or saprobe. This ratio (2PRG : 2SRG) is also seen in Ophiobolus disseminans (Ophdi1), which was isolated from the highly poisonous plant dog’s mercury (Mercurialis perennis) in Sweden, but whose ecology remains unknown. Five species in our Dothideomycetes study set did not have any of these six gene clusters in their genomes. This included Dissoconium aciculare, a hyperparasite of powdery mildew, the alcohol vapour-loving sooty mould Baudoinia panamericana, the extremophile Hortaea acidophila, the aquatic saprophyte Lineolata rhizophorae, and the human hair pathogen Piedraia hortae. We generated MEME profiles from the MAFFT alignments of these orthoMCL clusters using the MEME suite (Bailey ) and scanned 374 previously published fungal genomes (including 47 Dothideomycetes) using MAST (part of the MEME suite) for genes matching these profiles in their proteomes. We were able to identify homologs for one or more of these MEME profiles in 134 proteomes, all in the Pezizomycotina. Like in the Dothideomycetes, the presence or absence of the six gene clusters in these 134 proteomes was able to differentiate the pathogens and saprobes with >95 % accuracy. This can become a valuable tool to predict lifestyles for rarely seen and poorly studied species beyond Dothideomycetes, across Pezizomycotina as our efforts to understand the breadth of the fungal tree of life continues.

Conclusions

This study represents the largest genome comparison of fungal plant pathogens and introduces the genomes of 55 newly sequenced Dothideomycetes. The availability of whole-genome data also produced high-confidence phylogenetic trees with which we have been able to reclassify 25 organisms including the clarification of phylogenetic placement of several species that were previously considered to be “incertae sedis” and the addition of new taxonomic units. Ancestral character state analyses support a terrestrial saprobic lifestyle as being ancestral with subsequent derivations of plant pathogen, lichen, mycorrhiza, mycoparasite, aquatic and extremophile. Genomic traits were identified as enriched or reduced for selected ecologies, but often these were unannotated orthologous clusters or PFAM domains with broad high-level annotations, and subsequent experimentation is needed to test potential function related to ecology. The identification of gene family expansion/contraction is just one piece of the puzzle, with additional factors like changes in gene expression regulation, SNPs, etc. also playing important roles in determining the proteome and consequently the ecological adaptations of a particular species. Growth specific gene regulation is likely a key factor in understanding the saprophyte-pathogen continuum of several species and closely related lineages. Nevertheless, this large-scale comparison combined with machine learning was able to identify six gene clusters, most of which contain single-copy genes that can classify plant pathogens and saprobes in Dothideomycetes with greater than 95 % accuracy and we are encouraged by the ability of machine learning to differentiate pathogen vs saprobe based on genomic features alone. With increasing genomic data becoming available, a fundamental requirement for machine learning techniques, we can envision the use of these methods to identify genes and gene families driving other phenotypic changes, for example identifying genes related to fungal morphology (filamentous vs yeast vs dimorphic growth). This study encompassed 48 of the 110 known families in Dothideomycetes, but represents an important step in our efforts to explore fungal diversity and identify sets of genes determining fungal lifestyle.

Data availability

Genome data is available on the JGI MycoCosm portal using the DBID shown in Supplementary Table 4 in the form https://mycocosm.jgi.doe.gov/DBID; eg: https://mycocosm.jgi.doe.gov/Aaoar1 for Aaosphaeria arxii. All Dothideomycetes genomes available on MycoCosm can be seen at https://mycocosm.jgi.doe.gov/Dothideomycetes. OrthoMCL Clustering run data for 101 Dothideomycetes is available at https://mycocosm.jgi.doe.gov/clm/run/101-Dothideomycetes.3460. Data for newly sequenced genomes has been submitted to GenBank and accession numbers are available in Supplementary Table 5.