Basidioascus undulatus: genome, origins, and sexuality

Abstract

Basidioascus undulatus is a soil basidiomycete belonging to the order Geminibasidiales. The taxonomic status of the order was unclear as originally it was only tentatively classified in the class Wallemiomycetes. The fungi in Geminibasidiales have an ambiguously defined sexual cycle. In this study, we sequenced the genome of B. undulatus to gain insights into its sexuality and evolutionary origins. The assembled genome draft was approximately 32 Mb in size, had a median nucleotide coverage of 24X, and contained 6123 predicted genes. Previous morphological descriptions of B. undulatus relied on interpretation of putative sexual structures. In this study, nuclear staining and confocal microscopy showed meiosis occurring in basidia and genome analysis confirmed the existence of genes involved in meiosis and mating. Using 35 protein-coding genes extracted from genomic information, phylogenomic and molecular dating analyses confirmed that B. undulatus indeed belongs to a lineage distantly related to Wallemia while retaining a basal position in Agaricomycotina. These results, combined with differences in septal pore morphology, led us to move the order Geminibasidiales out of the Wallemiomycetes and into the new class Geminibasidiomycetes cl. nov. Finally, the concept of Agaricomycotina is emended to include both Wallemiomycetes and Geminibasidiomycetes.

Introduction

Using morphological characters, Matsushima (2003) described Basidioascus undulatus as the only species in the new genus Basidioascus based on a single strain isolated from the tropical rainforest soil in Cape Tribulation National Park, Queensland, Australia. Matsushima (2003) photographed fertile structures of B. undulatus, which he interpreted as asci giving rise to ascospores, bearing hooks that could be either croziers or clamps. The genus was named Basidioascus after what he interpreted as “basidialike asci” and the species epithet undulatus was given to recall “ascospores with a wavy wall” (Matsushima 2003). Also, B. undulatus produced a geotrichum-like asexual morph in culture, characterized by chains of aseptate arthroconidia. MycoBank (Robert et al. 2013) and Ainsworth & Bisby’s Dictionary of the Fungi (Kirk et al. 2008) classified B. undulatus in Saccharomycetes because the asexual morph was assumed to be a Geotrichum (a genus typified by an asexual morph and usually associated with sexual morphs in Galactomyces or Dipodascus).

During a survey of heat resistant fungi in Canadian soils, Nguyen et al. (2013) isolated nine additional strains of B. undulatus, and a second species of Basidioascus, named B. magus. A third yeast-like species, B. persicus, was recently described from soil in Iran (Nasr et al. 2014). Soil appears to be the main habitat for Basidioascus species and their distribution is probably broad. However, their ecological role is currently unknown, but they are presumably saprobic as are many soil inhabiting fungi (Domsch et al. 1980).

Surprisingly, phylogenetic analyses with rDNA sequences showed that Basidioascus was related to Wallemiomycetes (Basidiomycota) rather than Saccharomycetes (Ascomycota) (Nguyen et al. 2013). This finding initiated a revision of its taxonomy and a re-interpretation of its morphology as a basidiomycete. The structures identified as asci and ascospores by Matsushima (2003) are reinterpretted as thick-walled basidiospores, and the subtending cell as a basidium that usually produces a single basidiospore. Most unusual was that the basidia appeared to be forcibly discharged, leaving them collapsed with the basidiospore still attached by a long, cylindrical sterigma (Nguyen et al. 2013). The species of Basidioascus, and of its sister genus Geminibasidium (G. donsium and G. hirsutum), were classified in the new order Geminibasidiales (Nguyen et al. 2013). The Geminibasidiales are a phylogenetic sister group to Wallemiales and were placed tentatively under the class Wallemiomycetes (Nguyen et al. 2013). Wallemiales currently includes a single genus Wallemia with three species: W. sebi, W. muriae, and W. ichthyophaga (Zalar et al. 2005). The phylogenetic placement of Wallemiomycetes in the fungal kingdom was at first ambiguous (Matheny et al. 2006) because only a few protein coding genes were used in phylogenetic analyses and because ribosomal genes did not provide robustly supported conclusions. However, a few recent studies, through phylogenomic analyses with a large number of protein coding genes, demonstrate that this lineage is an early diverging one within Agaricomycotina (Padamsee et al. 2012, Zajc et al. 2013).

In this study, our first objective was to gain further insight into the sexuality of B. undulatus because the structures referred to as basidia and basidiospores were only putatively identified as such (Nguyen et al. 2013). For this purpose, we performed nuclear staining on these presumed sexual structures and observed them with laser confocal microscopy. Further, we sequenced the genome of B. undulatus and looked for genes involved in meiosis and mating to support our findings. Our second objective was to resolve the tentative placement of Geminibasidiales in Wallemiomycetes. We conducted phylogenomic analysis using 35 single copy protein-coding genes from the B. undulatus genome and we performed a molecular clock analysis to date the divergence of B. undulatus from Wallemia species and other fungi. The third objective was to investigate the septal pore morphology, which has proved significant in basidiomycete systematics, especially at class rank and particularly in lineages of Agaricomycotina (van Driel et al. 2009). We imaged the septal pore of B. undulatus, G. donsium, and W. sebi using transmission electron microscopy to support our interpretation of the higher classification of the Geminibasidiales.

Materials and Methods

Growth, DNA extraction and sequencing

The ex-type strain of Basidioascus undulatus (DAOM 241956) was inoculated in 2% malt extract broth in an Erlenmeyer flask on an orbital shaker at 25 °C for 2 wk. The broth culture was transferred to two 50 mL Falcon tubes and centrifuged at 10000 × g for 5 min. The liquid was decanted, leaving only the fungal tissue. The fungal tissue was frozen in liquid nitrogen and crushed with a sterile pestle. DNA was extracted with the OmniPrep kit (G-Biosciences, St Louis, MO) following the manufacturer’s instructions. DNA quality and quantity were verified with Qbit (Life Technologies, Burlington, Canada). Whole-genome sequencing (101 base pairs (bp) paired-end) was performed on an Illumina HiSeq 2500 with TrueSeq V3 chemistry at the National Research Council Canada facilty in Saskatoon (Saskatchewan).

Genome assembly and annotation

The quality of the reads was checked with the program FastQC v. 0.10.1 (http://www.bioinformatics.babraham.ac.uk/projects/fastqc/). Using fastx_trimmer (part of the FASTX-Toolkit v. 0.0.13; http://hannonlab.cshl.edu/fastx_toolkit/), eight bases from the 5′ end were trimmed to yield reads of 93 bp in length of higher quality. De novo assembly was performed using SPAdes v. 3.0 (Bankevich et al. 2012) with the BayesHammer error correction (Nikolenko et al. 2013) and mismatch correction enabled (parameters: —careful and k=21, 31, 41, 51, 61, 71, 81, 91). Final contigs were assembled into scaffolds with SSPACE v2.0 (Boetzer et al. 2011) (parameters: -x 1 -m 45 -o 10 -t 0 -r 0.7 -k 5 -a 0.7 -n 10 -z 2000 -T 16 -p 1) and contigs shorter than 2000 bp were discarded. Assembly statistics were generated with QUAST v. 2.3 (Gurevich et al. 2013). The assembly was checked by aligning the corrected reads onto the scaffolds using Bowtie2 v2.1.0 (Langmead & Salzberg 2012). Alignments produced by Bowtie2 in SAM format were converted to sorted BAM format by SAMtools v. 0.1.19 (Li et al. 2009) and statistics for coverage were generated with Qualimap v. 0.8.1 (Garcia-Alcalde et al. 2012). To benchmark the completeness of our genome assembly, CEGMA v. 2.5 (Parra et al. 2007) was run on the scaffolds to detect the percentage of conserved eukaryotic genes (CEGs).

Genome annotation was performed following established guidelines (Haas et al. 2011). Repeats in scaffold sequences were masked with RepeatMasker v. 4.0.5 (http://www.repeatmasker.org) (parameters: -no_is -species fungi) using the Repbase libraries (http://www.girinst.org/). The masked scaffolds were used as input for the MAKER2 v. 2.10 (Holt & Yandell 2011) genome annotation pipeline. In the MAKER2 pipeline, the GeneMark-ES v. 2.3e (Borodovsky & Lomsadze 2011) ab initio gene prediction tool was enabled and the NCBI RefSeq protein sequences were aligned to the genome using exonerate v. 2.2.0 (Slater & Birney 2005). Predicted gene models exhibiting strong evidence by exon alignment were exported as protein sequences and coding nucleotide sequences (CDS). Predicted gene models lacking evidence from exon alignment were discarded in downstream analyses. To determine function, the protein sequences were used as input for InterProScan 5RC6 (Jones et al. 2014) (parameters: -dp -f -t p -iprlookup -pa -goterms) and were also compared to the manually curated protein data set from UniProt/Swiss-Prot by blastp v. 2.2.28+. The results in XML format from blastp v. 2.2.28+ and InterProScan were loaded into Blast2GO v. 2.7.1 (Conesa et al. 2005) and merged to create an annotation table (available from the first author on request). The gene models with BLAST hits having e-value of less than 1.0E⁻¹⁰⁰ and mean similarity hit of ≥ 70% were assumed to be orthologs and they were given names following recommended conventions (http://www.uniprot.org/docs/proknameprot). Ribosomal RNA’s were predicted by RNAmmer v. 1.2 (Lagesen et al. 2007). Data files are publicly available at NCBI (Genome Accession No. JTLS00000000 version JTLS01000000; BioProject Accession No. PRJNA247992) and JGI MycoCosm portal (Grigoriev et al. 2014).

Identification of meiosis and mating genes

BLAST was used for finding the Basidioascus undulatus mating and meiosis genes. The protein sequences predicted from evidence-supported gene models, determined above, were formatted into a local BLAST database with makeblastdb v. 2.2.28+ (Camacho et al. 2009).

Genes previously determined to be involved in meiosis in Saccharomyces cerevisiae and Cryptococcus neoformans (Halary et al. 2011) were chosen as input queries (e-value cut off < 1.0E⁻⁰⁵) for the B. undulatus protein BLAST database using blastp v. 2.2.28+. Saccharomyces cerevisiae and C. neoformans were chosen because they are well studied genetically and all meiosis-specific proteins are present (Halary et al. 2011).

Mating genes were located with the protein domains identified by InterProScan, and by blastp v. 2.2.28+ using known mating genes in Saccharomyces cerevisiae as blastp input queries (e-value cut off < 1.0E⁻⁰⁵).

Phylogenomics and molecular dating

Protein sequences of selected fungi (Agaricus bisporus var. bisporus (H97) v. 2.0 (Morin et al. 2012), Alternaria brassicicola (Ohm et al. 2012), Arthrobotrys oligospora ATCC 24927 (Yang et al. 2011), Aureobasidium pullulans var. pullulans EXF-150 (Gostincar et al. 2014), Auricularia subglabra v. 2.0 (Floudas et al. 2012), Botrytis cinerea v. 1.0 (Staats and van Kan 2012), Coccidioides immitis RS (Sharpton et al. 2009), Coprinopsis cinerea (Stajich et al. 2010), Cryptococcus neoformans var. neoformans JEC21 (Loftus et al. 2005), Dacryopinax sp. DJM 731 SSP1 v. 1.0 (Floudas et al. 2012), Fomitiporia mediterranea v. 1.0 (Floudas et al. 2012), Fomitopsis pinicola FP-58527 SS1 v. 3.0 (Floudas et al. 2012), Malassezia globosa (Xu et al. 2007), Mixia osmundae IAM 14324 v. 1.0 (Toome et al. 2014b), Monacrosporium haptotylum CBS 200.50 (Meerupati et al. 2013), Neurospora crassa OR74A v. 2.0 (Galagan et al. 2003), Penicillium chrysogenum Wisconsin 54-1255 (van den Berg et al. 2008), Puccinia striiformis f. sp. tritici PST-130 (Cantu et al. 2011), Pyronema confluens CBS 100304 (Traeger et al. 2013), Rhizophagus irregularis DAOM 181602 v. 1.0 (Tisserant et al. 2013), Saccharomyces cerevisiae S288C (Goffeau et al. 1996), Sclerotinia sclerotiorum v. 1.0 (Amselem et al. 2011), Taphrina deformans (Cisee et al. 2013), Tilletiaria anomala (Toome et al. 2014a), Tremella mesenterica Fries v. 1.0 (Floudas et al. 2012), Trichoderma atroviride v. 2.0 (Kubicek et al. 2011), Tuber melanosporum (Martin et al. 2010), Ustilago maydis (Kamper et al. 2006), Wallemia ichthyophaga EXF-994 (Zajc et al. 2013), and W. sebi v. 1.0 (Padamsee et al. 2012)) were downloaded from the JGI MycoCosm portal (Grigoriev et al. 2014) and formatted into separate BLAST databases with makeblastdb v. 2.2.28+. The 246 reliable single copy ortholog protein data set from FUNYbase (Marthey et al. 2008) was downloaded. Only amino acid sequences coming from nuclear genes yielding a topological score of > 90% (Marthey et al. 2008) were considered for our phylogenomic analysis (information on exact genes chosen are available from the first author on request). These amino acid sequences were used as blastp v. 2.2.28+ search queries (e-value threshold < 1.0E⁻⁰⁵) against the protein databases built from data downloaded from JGI MycoCosm described above. Protein sequences were aligned with T-Coffee v10.00.r1613 (Notredame et al. 2000) (parameters: t_coffee sequence.fasta -output score_ ascii, aln) and poorly aligned regions and columns containing gaps were automatically discarded (parameters: t_coffee -other_pg seq_reformat -in sequence.aln -struc_in sequence. score_ascii -struc_in_f number_aln -action +use_cons +keep ‘[8-9]’ +rm_gap 1 > sequence.best.aln). The alignments were concatenated and converted to PHYLIP format with SeaView v. 4.5.3 (Gouy et al. 2010). Three independent phylogenomic analyses were performed with the parallelized version of PhyloBayes 3 (pb_mpi v1.4) (Lartillot et al. 2013, Lartillot et al. 2009) using the CAT-GTR model (Lartillot & Philippe 2004). Analyses were stopped when convergence was attained (effective size > 100 and maxdiff < 0.1 determined with the programs bpcomp (parameters: -x 1000 50 run1 run2 run3) and tracecomp (parameters: -x 1000 run1 run2 run3), which are part of pb_mpi software package). A lognormal ‘relaxed clock’ molecular dating analysis was performed with the non-parallelized version of PhyloBayes 3 (pb v3.2e) with a birth-death prior using the tree topology generated from the converged Bayesian analysis above (parameters: -d combined.phy -T bpcomp.con.tre -r outgroup.txt -cal calib. txt -ln -bd md1). Rhizophagus irregularis (Glomeromycota) was specified as the outgroup. Date constraints previously determined were used as calibrations (Hibbett et al. 1997, Smith et al. 2004, Taylor & Berbee 2007, Berbee & Taylor 2010, Prieto & Wedin 2013); the exact calibrations are available from the first author on request. Chronograms and statistics were obtained with readdiv (parameters: -x 1000 50 md1), which is part of PhyloBayes 3.

Confocal laser scanning microscopy

To study nuclear behaviour and to look for indicators of meiosis, Basidioascus undulatus DAOM 241956 was grown on corn meal agar (CMA, Acumedia Manufacturers, Lansing, MI) for 1 wkand mounted in DNA stains: DAPI-Fluoromount-G™ mounting medium (EMS, Hatfield, PA) or aqueous SYTO 9 (25 µM) (Life Technologies, Burlington, ON). Samples were visualized under confocal laser scanning microscopy using an LSM 510 DUO (Carl Zeiss MicroImaging, Göttingen, Germany) with a Plan-Apochromat 40×/1.4 Oil DIC objective and electronic zoom 4. An excitation diode laser (405 nm) and emission light (420–700 nm) were used for DAPI. An excitation Argon laser (488 nm) and emission light (505–550 nm) were used for SYTO 9. Images were captured using ZEN 2009 Imaging Software (Carl Zeiss MicroImaging).

Transmission electron microscopy

Actively growing hyphae of Basidioascus undulatus DAOM 241956, Geminibasidium donsium DAOM 241966, and Wallemia sebi CBS 633.66 were prepared for transmission electron microscopy using cryo-preparation methods. Hyphae were grown on thin, sterile, deionized dialysis membrane segments overlaying appropriate media at 23 °C. The leading edge of growing mycelia and supporting membranes were trimmed with a sharp razor blade to approximately 5 × 5 mm and after 30–40 min (time to recover from trimming) were removed from the agar surface and immediately cryofixed by rapid plunging into liquid propane cooled to −186 °C with liquid nitrogen (Hoch 1986, Roberson & Fuller 1988, McDaniel & Roberson 2000). After rapid freezing, the samples were freeze-substituted in 1% glutaraldehyde (w/v) and 1% tannic acid (w/v) in anhydrous acetone at −85°C for 72 h. After washing in cold acetone (−85 °C), the samples were warmed slowly to room temperature in 1% OsO₄ (w/v) in acetone, washed in acetone, and infiltrated and flat embedded on glass slides in Spurr’s resin (Spurr 1969). Using phase contrast optics (100×), we examined the slides for well-preserved hyphae and mounted the selected cells on resin blocks (Howard & O’Donnell 1987) then hand-trimmed them. Selected hyphae were sectioned using a Leica Ultracut R ultramicrotome (Leica Microsystems, Bannockburn, IL), collected on copper grids, and post-stained for 10 min in 2% uranyl acetate in 50% ethanol and for 5 min in Sato’s lead citrate (Hanaichi et al. 1986). Sections were then examined using a JEOL 1200EX (JEOL, Tokyo,) transmission electron microscope equipped with a SIA L3C CCD camera (SIA, Duluth, GA). Measurements from captured images were made with ImageJ (Schneider et al. 2012).

Results

Genome sequencing, assembly and annotation

Short-read Illumina sequencing generated approximately 12 million paired end reads (six million reads in the forward (R1) direction, and six million reads in the reverse (R2) direction) of 101 bp in length each. After trimming to 93 bp, about 1.1 Gb of data was assembled de novo to yield a genome assembly size of about 32 Mb. The GC content in Basidioascus undulatus was 58%. The final assembly contained 2992 scaffolds and the longest scaffold was 97 Kb. The nucleotide coverage varied for each scaffold but it was 28X on average and the median nucleotide coverage was 24X. The N50 statistic was 15 Kb. According to CEGMA, 81% and 88% of complete and partial CEGs, were detected respectively.

A total of 13935 gene models were detected ab initio using GeneMark-ES, but only 6123 gene models were supported by evidence from protein alignment to the NCBI RefSeq fungal protein data. We kept the set of 6123 gene models with evidence for further manual annotation and used it in downstream analyses. Only 3681 of these (60%) were considered complete because they contained a start and stop codon while the remaining 2442 (40%) lacked either a start codon, a stop codon, or both.

Meiosis and meiosis specific genes

To determine whether the putative sexual structures identified by Nguyen et al. (2013) truly represented basidia, we followed the fungus’s ontogeny and performed nuclear staining by laser confocal microscopy (Fig. 1). Three nuclear divisions normally occur during basidiospore maturation: meiosis I, meiosis II, followed by four different patterns of post-meiotic mitosis (Duncan & Galbraith 1972). We observed dikaryotic nuclei in the basidia (Fig. 1A), karyogamy (Fig. 1B), anaphase I (Fig. 1C), and the telophase I stage where the basal lateral projection was collapsed (Fig. 1D). Four nuclei from telophase II (Fig. 1E) were seen in both non-discharged and discharged basidia (Fig. 1F). The migration of a nucleus through the sterigma (Fig. 1G) into the basidiospore (Fig. 1H) was observed as the basidiospore matured. There were three remaining nuclei in the collapsing basidium (Fig. 1I) and they eventually degenerated at a later stage when the basidium completely collapsed (Fig. 1J). Whole arthroconidia became stained instead of just their nuclei when DAPI was used. A single nucleus could be spotted in each arthroconidium when stained with SYTO 9 (Fig. 1K–L).

Fig. 1

Basidioascus undulatus (DAOM 241956) sexual and asexual structures stained with DAPI (A–J) and SYTO 9 (K–L) and imaged with confocal microscopy. A. Dikaryotic basidium (2 nuclei). B. Karyogamy. C. Anaphase I. D. Telophase I and collapsed basal lateral projection. E. Telophase II (4 nuclei). F. Ejected basidium (4 nuclei). G–H. Maturation of a basidiospore on a basidium and the migration ofa nucleus through the sterigma and into the basidiospore. I. Collapsing basidium with the three remaining nuclei. J. A totally collapsed basidium and probably one nucleus inside the mature basidiospore. K–L. Single nuclei inside arthroconidia. Bar = 5 µm.

Our microscopic observations correlated with our analysis of the genome to detect meiotic genes (Table 1). All meiosis-specific genes, as defined in Malik et al. (2008) and Halary et al. (2011), were found in the genome.

Table 1 List of known and putative meiosis genes in the genome of Saccharomyces cerevisiae, Cryptococcus neoformans and Basidioascus undulatus. Meiosis specific proteins are highlighted in red.

Mating genes

We looked for mating genes such as those encoding for homeodomain proteins, G-protein coupled pheromone receptor, high mobility group (HMG) DNA binding proteins, mitogen-activated protein kinases (MAPK, MAPKK, MAPKKK) and the subunits of the trimeric GTPase protein (Gα, Gβ, Gγ) (James et al. 2013). We followed the Saccharomyces cerevisiae standard for gene names. Most of the mating response genes and the pheromone processing genes (Table 2) were detected in the B. undulatus genome. All of the mating genes detected were located on different scaffolds, except for KSS1 and FUS3, which co-occurred on scaffold 2837. The presence of most mating genes suggests that a mating type locus exists, but its structure and gene order could not be determined because our assembled genome was too fragmented.

Table 2 List of mating genes.

Septal pore morphology

We imaged the septal pores of Wallemia sebi, Basidioascus undulatus, and Gemmibasidium donsium, all currently classified in Wallemiomycetes, using identical fixation methods to strengthen comparisons (Fig. 2). Basidioascus undulatus, G. donsium and W. sebi all have a dolipore septum characteristic of Agaricomycotina. The pore swelling of B. undulatus and G. donsium is electron dense, but is not in W. sebi. The adseptal tubular extensions that arise from sheets of endoplasmic reticulum in W. sebi were obvious (Fig. 2E–F) but they were either absent or unclear in B. undulatus (Fig. 2A–B) and G. donsium (Fig. 2C–D). Electron-dense regions were evident near the septal pore in all three fungi, but this area was denser in W. sebi than in B. undulatus and G. donsium. An electron-dense septal pore occlusion extending across the septal pore was observed in all three fungi. This occlusion was non-membranous in W. sebi, which agrees with previously published findings (Terracina 1974, Padamsee et al. 2012) and the W. sebi data on the AFTOL Structural and Biochemical Database (Celio et al. 2007). The septal pore occlusion in B. undulatus and G. donsium had a membrane, making it more clearly defined compared to that of W. sebi. Furthermore, this occlusion was almost twice as wide in B. undulatus and G. donsium (∼0.12 µm) compared with W. sebi (∼0.06 µm). Also, the septal pore of W. sebi had striations that appeared to be fine fibrils vertically arranged at the pore opening, but these were absent in B. undulatus and G. donsium. We did not see a septal pore cap in our sections of W. sebi, B. undulatus or G. donsium.

Fig. 2

Transmission electron micrographs showing septal pore morphology. A–B. Basidioascus undulatus (DAOM 241956 ex-type). C–D. Geminibasidium donsium (DAOM 241966 ex-type). E–F. Wallemia sebi (CBS 633.66). Bar = 20 µm.

Phylogenomic analysis and molecular dating

Because Basidioascus is currently classified in Wallemiomycetes, and the phylogenetic position of the fungi classified in this class has been unstable in different analyses as a result of sparse taxon sampling, low gene sampling and conflicts in certain genealogies (Zalar et al. 2005, Matheny et al. 2006), we performed a phylogenomic analysis using protein sequences from 35 single copy protein-coding genes to study its phylogenetic position in a more robust fashion (Fig. 3A). Our results were similar to those obtained from the rDNA phylogenetic analysis in Nguyen et al. (2013), showing Basidioascus as a distant lineage to Wallemia, with a posterior probability of 1.00. The results were also similar to those of Padamsee et al. (2012) and Zajc et al. (2013), where the Wallemiomycetes were shown to be an early diverging lineage of Agaricomycotina.

Fig. 3

Phylogenetic trees resulting from phylogenomic analysis and molecular dating. A. Consensus topology and branch lengths from analyses of concatenated amino acid sequences from 35 single copy genes with a total of 10129 data columns. Gapped and poorly aligned sites were removed. Analyses were performed with the CAT-GTR model. Posterior probabilities are shown at the nodes of the tree. Scale bar indicates expected changes per site. B. Chronogram resulting from a lognormal relaxed molecular clock analysis with the birth-death prior. The mean divergence time with standard error are shown at each node. Circled nodes were pre-calibrated before the analysis. The paleontological periods, in million years ago (Mya), are shown as a scale at the bottom. Rhizophagus irregularis (Glomeromycota) was used as the outgroup.

We dated the divergence of Basidioascus from other lineages using some reliable fossil data of basidiomycetes (Hibbett et al. 1997, Smith et al. 2004) and the latest calculated calibration points (Berbee & Taylor 2010, Prieto & Wedin 2013) (Fig. 3B). Under these assumptions, the split between B. undulatus and Wallemia was estimated at 250±29 Mya, and that between W. sebi and W. ichthyophaga was about 33±6 Mya.

Discussion

Genome sequencing and annotation

Basidioascus undulatus is the first fungus in Geminibasidiales with a sequenced genome. The 1.1 Gb of sequence data gave a median nucleotide coverage of 24X and a final assembly that contained 2992 scaffolds, which is similar to what has been produced in many recent whole genome sequencing studies (e.g. Sims et al. 2014). Wallemia sebi and W. ichthyophaga were the closest known relatives to B. undulatus. Despite being considered close relatives, the genome assembly size of B. undulatus (32 Mb) was larger than W. sebi (9.8 Mb) (Padamsee et al. 2012) and W. ichthyophaga (9.6 Mb) (Zajc et al. 2013). Also, the GC content of B. undulatus was higher (58%) compared to W. sebi (40%) and W. ichthyophaga (45%) (Table 3). Based on these results, the fungi of these two genera are more distantly related than initially suspected.

Table 3 Genome assembly and annotation statistics of Basidioascus undulatus compared to Wallemia sebi and W. Ichthyophaga.

We were able to detect 88% of partial core eukaryotic genes (CEGs). We assumed that downstream annotation with this genome assembly should, in theory, reveal close to the same proportion of the total number of genes contained in B. undulatus.

Many gene prediction tools are available, but GeneMark-ES was selected because it was shown previously to be accurate at detecting genes in fungal genomes (Ter-Hovhannisyan et al. 2008). GeneMark-ES considers branch point sequences in the intron model, which guides lariat formation during splicing, providing greater accuracy at locating intron boundaries. Because RNA data were not acquired in our study, we chose to use GeneMark-ES, which does not require RNAdata for algorithm training. Nevertheless, genome annotation based on genome sequences alone has limitations. First, without RNA data, CDS predictions cannot be validated and genes from splice sites that use donor and acceptor sequences other than the canonical GT–AG introns cannot be predicted. The frequencies of the non-canonical GC–AG introns are 1.0–1.2% in some ascomycetes (Rep et al. 2006) and could be as high as 3% in basidiomycetes (Misiek & Hoffmeister 2008). To compensate for the absence of RNA data, the fungal NCBI RefSeq protein sequence data set (13 March 2014 release), which contains curated protein sequences from completed fungal genome projects, was aligned to our B. undulatus genome using the program “exonerate”. In the MAKER2 pipeline, aligning protein sequences to predicted genes on the genome generates ‘evidence’ expressed as the Annotation Evidence Distance (AED) score for exons (eAED) (Eilbeck et al. 2009). The eAED score is a metric that measures how a predicted gene determined from GeneMark-ES agrees with the protein alignment evidence from exonerate while accounting for protein reading frames shifts. This score is helpful for assessing annotation quality of a predicted gene (Eilbeck et al. 2009). Overall, 13 935 genes were predicted by GeneMark-ES but only 6123 of those genes were supported by protein alignment evidence. For each of the 6123 predicted genes, we made further functional annotation with a BLAST search against the UniProt/Swiss-Prot manually curated data set and an InterProScan analysis. Using this information, 452 genes were manually annotated with confidence. Future releases of the B. undulatus genome with RNA sequence data would increase the number of annotated genes and allow discovery and validation of completely novel genes or splice variants.

Meiosis

We used genomic information in combination with confocal microscopy to gain further understanding of sexuality in Basidioascus undulatus. Fungal taxonomists infer function of meiosporangia in newly discovered fungi based on morphological similarities to proven meiosporangia in similar or related fungi. Proof of meiosis by nuclear staining is not usually required. In our case, B. undulatus is different morphologically, ontogenetically, and phylogenetically from known basidiomycetes and it was difficult to assess whether previous morphological interpretations by Nguyen et al. (2013) were correct without nuclear staining experiments. We suggested previously that the clavate structures found on somatic hyphae were basidia (Nguyen et al. 2013). These basidia were deciduous, forcibly discharged, and had a basal lateral projection that eventually collapsed. Sterigmata grew on the basidium, and a basidiospore developed at the tip of each sterigma. To verify these interpretations, we visualized these putative sexual structures with nuclear staining and confocal microscopy. Our observations suggested that meiosis occurred in the clavate structures (Fig. 1) providing the evidence for previous interpretation of them as basidia. Post-meiotic mitosis in the basidium was not seen but it was difficult to visualize nuclei in the mature basidiospores. Although we saw evidence of meiosis with microscopy, we also confirmed the presence of meiotic genes in the genome (Table 1). We did not detect HOP2 among the 10 meiosis specific genes. The HOP2 gene product is involved in preventing synapsis between non-homologous chromosomes during meiotic double-strand break repair (Malik et al. 2008). The function of HOP2 could be served by another unknown gene or gene product in B. undulatus, or perhaps our genome assembly was not complete enough to detect it. With 88% of CEGs detected, about 1 out of every 10 genes would not be detected during the annotation procedure.

Basidioascus undulatus produces cells attached in chains, which are presumed to be arthroconidia (Nguyen et al. 2013). These arthroconidia were indeed single-nucleated (Fig. 1K–L) and therefore represent an asexual morph of B. undulatus. Given these results, B. undulatus indeed exhibits both sexual and asexual morphs in culture.

Mating

The genes located at the MAT locus orchestrate the fungal mating process and determine the sex of individuals (Lee et al. 2010). In the fungal kingdom, mating is most comprehensively studied in Saccharomyces cerevisiae. Therefore, we used what is known about mating in S. cerevisiae as a model to guide our interpretations of mating genes in the B. undulatus genome.

Homeodomain (HD) transcription factors control the expression of pheromone and pheromone receptor genes. Normally, dimerization occurs between two paralogous HD transcription factors, HD1 and HD2, to form HD1-HD2 complex. In Ustilago maydis, this is termed the bE/bW heterodimer, and in S. cerevisiae the a1/α2 complex (Lee et al. 2010). We identified the HD1 (8946F71C4C) and HD2 (AF24C22F02) proteins in B. undulatus.

Pheromone production and processing are also important in mating. All orthologs for pheromone processing (KEX1, KEX2, STE13, RAM1, RAM2, STE14, AXL11, STE6) were detected in B. undulatus except for the RCE1/STE24 genes responsible for cleavage of the pheromone peptide (Table 2).

Upon release into the environment, pheromones must bind to pheromone receptors to initiate a signaling cascade involving the mitogen-activated protein kinase (MAPK) pathway to turn on mating genes (Jones & Bennett 2011). MAPK pathways regulate the activity of high mobility group (HMG) DNA binding proteins, which are transcription factors regulating pheromone responsive genes in the MAT loci (Hartmann et al. 1996). Most of the mating response genes were detected (STE3, KSS1, FUS3, STE7, STE11, GPA1, STE4, STE20) (Table 2) except the G-protein γ subunit (STE18), the transcription factor STE12, the coordinator of MAPK pathway STE5, and the cell cycle control FAR1 protein. We identified four HMG DNA binding proteins (B013112742, C8B77D01D8, C1AD40EA22, and 9CC40B0A26).

We were unable to isolate single arthroconidia or successfully germinate basidiospores so mating tests were not possible. Thus, it is unknown whether B. undulatus is homothallic or heterothallic. Wallemia sebi apparently has a bipolar mating system (Padamsee et al. 2012) but no discernable mating type locus could be identified in the genome of W. ichthyophaga (Zajc et al. 2013). Our genome assembly was too fragmented to determine the structure of the mating type locus conclusively and therefore the mating system in B. undulatus remains uncharacterized until more sequencing is done. The detection of a near complete cellular machinery for the mating response, pheromone processing, and meiosis suggests that B. undulatus is capable of outcrossing. However, future experiments with gene knockout mutants will be needed to validate these identified genes.

Higher classification and divergence of Geminibasidiales

Basidioascus undulatus and Geminibasidium donsium have been classified in Geminibasidiales, a sister order to Wallemiales. Geminibasidiales was placed tentatively in the class Wallemiomycetes (Nguyen et al. 2013). The main reason for this initial inconclusive placement was the absence of septal pore ultrastructure data and that only ribosomal genes could be used to make a phylogenetic analysis with enough taxon sampling to reach a conclusion. More genomes are available now and phylogenomic analysis with protein coding genes is more feasible and informative.

Our phylogenomic analysis shows the same relationships previously determined with rDNA, where Wallemia and Basidioascus are distantly related sister groups that occupy a basal phylogenetic position in Agaricomycotina (Fig. 3A). The distance of this relationship is also reflected in differences in GC content and genome size between B. undulatus and Wallemia discussed above. According to the sampling and assumptions made in our analyses, molecular dating suggests that the split between B. undulatus and W. sebi occurred 250 Mya (Fig. 3B), which is around the same time that Sordariomycetes split from Leotiomycetes (253 Mya) and Exobasidiomycetes from the Ustilaginomycetes and Malasseziomycetes (262 Mya), but we also note that not all nodes corresponding to class ranks in our analysis represent divergences at around 250 Mya. If our dates are accurate, it means Wallemia and B. undulatus shared a common ancestor until a divergence event occurred between the Permian and Triassic period. Interestingly, it was during this period that Earth experienced its most severe extinction event (Sahney & Benton 2008). The split between W. sebi and W. ichthyophaga was dated at 33±6 Mya, an estimate that differs from previous analyses; a 11.9 Mya split was proposed by Zajc et al. (2013), but with different methods, taxon sampling, and calibrations used for the molecular dating analysis. These estimates should all be considered approximate because the fungal fossil record is rather limited and many studies generate radically different age estimates for the same divergence events (Lücking et al. 2009).

Before sequence data began to dominate fungal systematics, morphological characters of the septa and septal pore caps were considered important for classifying fungi at the class or higher taxonomic levels (Khan & Kimbrough 1982, Oberwinkler & Bandoni 1982, Müller et al. 2000, van Driel et al. 2009) and septal pore morphology is still considered essential for delimiting classes in Basidiomycota. However, the overall septal pore morphologies of B. undulatus and G. donsium are similar and there are significant differences when compared to W. sebi (see p. 219 above; Fig. 2). The septal pore cap is probably not a reliable or practical feature for the classification of Wallemia, Basidioascus and Geminibasidium because it is sometimes present or absent in Wallemia (Padamsee et al. 2012, D.J. McLaughlin pers. comm.) and was not seen in B. undulatus or G. donsium.

Given the results obtained from phylogenomics analysis, molecular dating, and septal pore ultrastructure, we propose to classify Geminibasidiales in a new class Geminibasidiomycetes. Although Geminibasidium could not be included in the phylogenomic analysis and molecular dating because its genome is not yet sequenced, we are confidentfrom previous rDNA phylogenetic analysis (Nguyen et al. 2013), morphological characters, and similarities in the septal pore morphology (Fig. 2) that Geminibasidium will remain a close sister group to Basidioascus. Previous phylogenies position Wallemiomycetes at the base of Basidiomycota (Matheny et al. 2006, Zalar et al. 2005), as a sister group to or as the earliest diverging lineage of Agaricomycotina (Padamsee et al. 2012, Nguyen et al. 2013, Zajc et al. 2013). The classification for the Wallemiomycetes is shown as incertae sedis in the Basidiomycota in the NCBI taxonomy (Federhen 2012), MycoBank (Robert et al. 2013), and the Encyclopedia of Life (http://www.eol.org). In the taxonomic section below, we emend the description of the subphylum Agaricomycotina to include the classes Geminibasidiomycetes and Wallemiomycetes along with the Agaricomycetes, Tremellomycetes, and Dacrymycetes, which are well delimited classes of Agaricomycotina (Hibbett 2007). By formally connecting the basal Geminibasidiomycetes and Wallemiomycetes to the Agaricomycotina, we eliminate the uncertain status of these taxa. The subphylum name Wallemiomycotina (MycoBank no. MB 550364) was published recently to accommodate Wallemia (Doweld 2014). There was no discussion of the rationale behind this proposal and we consider it unnecessary to include Wallemiomycetes (with or without the Geminibasidiomycetes) in its own subphylum.

Taxonomy

Agaricomycotina Doweld, Prosyllabus: lxxxvii (2001).

MycoBank MB560553

Homonym: Agaricomycotina R. Bauer et al., Mycol. Prog. 5: 45 (2006).

Description: Members of the Basidiomycota that have a cell wall carbohydrate composition with a dominance of glucose and presence of xylose and having a type B secondary structure of the 5S RNA (Bauer et al. 2006). Fungi that belong to this subphylum are classified in the classes Agaricomycetes, Dacrymycetes, Tremellomycetes, Wallemiomycetes or Geminibasidiomycetes.

Notes: Although Agaricomycotina Dowell (2001) is the oldest name, it is the concept of proposed by Bauer et al. (2006) (i.e. Hymenomycetes sensu Swann & Taylor 1995) that is widely used today. We modified it above only to explicitly state all of the accepted classes within Agaricomycotina. The cell wall of B. persicus was composed mostly of glucose and xylose (Nasr et al. 2014), and assuming that this is the same for other species of Basidioascus and Geminibasidium, this placement for the class Geminibasidiomycetes is appropriate. The cell wall carbohydrate composition of Wallemia species is unknown.

Wallemiomycetes Zalar et al., Antonie van Leeuwenhoek 87: 322 (2005).

MycoBank MB501496

Description: Class of xerophilic basidiomycetes belonging to the subphylum Agaricomycotina. These fungi produce basauxic anamorphs and do not produce basidiomata in culture. Species have dolipore septa with adseptal tubular extensions that arise from sheets of endoplasmic reticulum that form the septal pore cap. The septal pore cap is sometimes absent. The septal pore has an electron-dense non-membranous septal pore occlusion and striations that are oriented vertically.

Type order Wallemiales.

Notes: This description is altered to exclude characters of Geminibasidiales previously added to the concept of Wallemiomycetes by Nguyen et al. (2013). Information about the septal pore morphology was added, as detailed in Padamsee et al. (2012).

Geminibasidiomycetes H.D.T. Nguyen & Seifert, cl. nov.

MycoBank MB811680

Description: Class of xerotolerant basidiomycetes belonging to the subphylum Agaricomycotina. Basidiomata not produced in culture. Basidia arising from somatic hyphae or from swollen basidium-bearing cells (primary cells) with a basal lateral projection occurring either on the basidium or the swollen primary cell. Basidiospores symmetrical on sterigma, not forcibly discharged, and brown at maturity. Arthroconidial and/or yeast-like asexual morphs sometimes produced. Species have a dolipore septum that is electron-dense at the pore swelling with an electron-dense membranous septal pore occlusion. Some species are heat resistant.

Type order: Geminibasidiales.