A High-Quality Genome Resource for the Oak Wilt Pathogen Bretziella fagacearum

Genome Announcement

Oak wilt, caused by the fungal pathogen Bretziella fagacearum (Bretz) Z.W. de Beer, Marinc., T.A. Duong & M.J. Wingf., (formerly Ceratocystis fagacearum [Bretz]), is a significant disease affecting oaks (Quercus spp.) in the United States (de Beer et al. 2017; Henry 1944). B. fagacearum spreads through belowground networks of grafted roots among neighboring trees, thereby forming centers or pockets of diseased host species (Blaedow and Juzwik 2010; Chahal et al. 2024). The pathogen is also transmitted aboveground by sap or nitidulid beetles (Coleoptera: Nitidulid) and oak bark beetles (Pityophthorus spp.) (Gibbs and French 1980; Rexrode and Jones 1970). B. fagacearum infections can disrupt both urban and natural forest ecosystems if not detected and managed in a timely manner. Currently, oak wilt has been confirmed in 24 states in the United States, and widespread epidemics are reported in the Great Lakes region and Texas (Appel 1986; Chahal et al. 2019b, 2021; Juzwik et al. 2011; Rexrode and Lincoln 1965). Oak wilt continues to spread, expanding into new geographical locations where it has not been previously reported (Gauthier et al. 2023; McLaughlin et al. 2022). Detection of oak wilt in Ontario, Canada, in June 2023 marked the first identification of B. fagacearum outside the United States (Canadian Food Inspection Agency 2023).

The epidemiology of oak wilt varies between the Great Lakes region and Texas disease epidemics. In the Great Lakes region, red oaks (section Lobatae) are prevalent and highly susceptible compared with white oaks (section Quercus) that exhibit moderate to high levels of resistance (Appel 1995; Juzwik et al. 2011). In the Midwestern United States, B. fagacearum infections have also been reported on orchard-grown chestnuts (Bretz and Long 1950; Chahal et al. 2024). In Texas, the widespread host Texas live oak (Quercus fusiformis, section Virentes) shows intermediate susceptibility compared with red and white oaks (Appel 1986, 1995). However, clonal propagation of Texas live oak through sprouts originating from a shared root system facilitates the spread of B. fagacearum, leading to expanding pockets of mortality (Appel 1995). The outcomes of oak wilt occurrences observed in Texas live oak populations have not been observed in any other host or state within its distribution range.

The sporulation structures of B. fagacearum, known as mycelial mats, predominantly occur on red oaks, rarely on white oaks and have never been reported on Texas live oak (Appel 1995; Juzwik et al. 2011). In the Great Lakes region, oak wilt incidence is high due to the homogenous stands of highly susceptible red oaks, prevalent fungal sporulation, sandy soils that favor root grafting, and the presence of active insect vectors (Juzwik et al. 2011). In the Great Lakes region, insect-vector activity and mycelial mat formation peak during the growing season, particularly from April to June (Chahal et al. 2021; Curl 1955; Juzwik et al. 2011; Morris et al. 2024). No mycelial mat formation has been reported in Michigan from mid-November through mid-April (Chahal et al. 2021), whereas in Texas, insect vectors remain active year-round, coinciding with the continuous formation of fungal mycelial mats (Appel 1995). In Texas, B. fagacearum infection centers typically span several hectares and expand at rates up to 50 m/year (Appel 2008). This contrasts with the slower growth and smaller sizes of foci in deciduous red oaks in the Great Lakes region, expanding at rates of 12 m/year in Michigan and 1.9 to 7.6 m/year in Minnesota (Bruhn et al. 1991; Gearman and Blinnikov 2019). Human-mediated activities, such as poorly timed pruning of oaks and the long-distance transmission of infected and mycelial mat-bearing firewood, contribute to vector-mediated infections across the distribution range of oak wilt (Gearman and Blinnikov 2019). Factors contributing to differences in oak wilt epidemiology include the land cover and distribution of different host species, climate, soil type, human population density (Gearman and Blinnikov 2019), and potentially the genotypes of B. fagacearum. The diversity of B. fagacearum genotypes within geographically different populations may contribute to observed variations in disease incidence and epidemiology (Chahal 2024).

The aim of this study was to generate a high-quality genome sequence of B. fagacearum and annotate it using RNA sequencing data from the same isolate, providing a reference genome for future studies. B. fagacearum isolate C519 was isolated from infected Q. rubra in Sherburne County, MN, in 1992. The isolate was confirmed as B. fagacearum using real-time PCR (Bourgault et al. 2022) and by PCR amplifying and sequencing the 60S, LSU, and MCM7 gene regions (de Beer et al. 2017). Mycelium was grown on potato dextrose broth and harvested following the established protocol (Chahal et al. 2019a, 2022). Genomic DNA was extracted from harvested mycelium using a phenol-chloroform protocol described in Parada-Rojas and Quesada-Ocampo (2018); 5 µg of genomic DNA was sheared to >10 kb using Covaris g-Tubes. The sheared DNA was treated with exonuclease to remove single-stranded ends and a DNA damage repair mix, followed by end repair and ligation of blunt adapters using the SMRTbell Template Prep Kit 1.0 (Pacific Biosciences). The library was purified with AMPure PB beads. PacBio sequencing primer was annealed to the SMRTbell template library, and sequencing polymerase was bound using the Sequel II Binding Kit 1.0. The prepared SMRTbell template libraries were sequenced on a Pacific Biosciences Sequel II sequencer using 8M v1 SMRT cells and version 1.0 sequencing chemistry with 1 × 900 min sequencing movie run times. Sequence data were processed with the JGI QC pipeline to remove artifacts. The genome was assembled with Falcon version pb-assembly = 0.0.2|falcon-kit = 1.2.3|pypeflow = 2.1.0 (Chin et al. 2016), improved with FinisherSC version 2.1 (Lam et al. 2015), and polished with Arrow version SMRT Link v7.0.1.66975. Contigs shorter than 1,000 bp were excluded, and those identified as contaminants were removed. Contaminants were detected by performing basic local alignment search tool (BLAST) searches of the contigs separately against bacterial/viral and fungal UniProt accessions. Contigs with a higher number of bacterial/viral hits in length compared with fungal hits were discarded. To better assemble B. fagacearum, which had substantially higher coverage than the contaminant, preassembled reads (preads) were subsampled using BBTools v38.79 (reformat.sh samplerate = 0.3). The subsampled reads were assembled with Flye v2.6 (flye –pacbio-corr -g 40m –asm-coverage 50) and polished with GCpp (–algorithm arrow, SMRT Link v8.0.0.80529). The genome features are presented in Circos plots (Krzywinski et al. 2009) that were created using SAMtools (Danecek et al. 2021), bedtools (Quinlan and Hall 2010) was used for data extraction, and bedtops (Neph et al. 2012) was used for converting between formats (Fig. 1). The mitochondrial genome was assembled separately from the Falcon preads using the joint genome institute tool, assemblemito.py (https://code.jgi.doe.gov/gaa/jgi-fungal-asm/-/tree/dev/bin). The preads were filtered and polished with Arrow (SMRT Link v7.0.1.66975) (Fig. 2).

Download as PowerPoint

Download as PowerPoint

For transcriptome analysis, total RNA from harvested mycelium of B. fagacearum C519 was extracted using the EZNA Fungal RNA Mini Kit (Omega Bio-tek Inc., Norcross, GA). Plate-based RNA sample preparation was performed on the PerkinElmer Sciclone NGS robotic liquid handling system using Illumina's TruSeq Stranded mRNA HT sample prep kit, following the protocol outlined by Illumina for poly-A selection of mRNA. Total RNA starting material was 1 µg/sample, and eight cycles of PCR were used for library amplification. The prepared library was quantified using KAPA Biosystem's next-generation sequencing library qPCR kit (Roche) and run on a Roche Light Cycler 480 real-time PCR instrument. The quantified libraries were multiplexed and prepared for sequencing on the Illumina NovaSeq 6000 platform using NovaSeq XP v1.5 reagent kits and an S4 flow cell, following a 2 × 150-bp indexed run recipe. After sequencing, read artifacts (k-mer = 25 bp, 1 mismatch) were detected with BBDuk (https://sourceforge.net/projects/bbmap/). Detected artifacts were trimmed from the 3′ end of reads. General quality trimming was performed with the Phred trimming method set at Q6. Reads shorter than 25 bp or one-third of the original read length were removed, as well as RNA spike-in reads, PhiX reads, and reads containing any ambiguous characters (Ns). Filtered transcriptome sequencing (RNA-Seq) reads were assembled using Trinity v2.8.5 (Grabherr et al. 2011). RNA-Seq capture was performed by aligning a 1% subsample of RNA reads to the final DNA assembly using BBTools v38.67. The genome was annotated using the JGI Annotation pipeline (Grigoriev et al. 2014). Gene clusters of interest were searched in the clustering table available at the JGI MycoCosm portal (Grigoriev et al. 2014). Clusters were computed following the TRIBE-Markov cluster (MCL) clustering method (Enright et al. 2002) from all-versus-all BLAST analysis of the proteins in the set of organisms included in a clustering run.

Following joint genome institute genome annotation pipeline, we predicted 7,445 full length transcripts, with an average transcript length of 1, 722 bp and an average gene length of 1,928 bp. The assembled genome consisted of nine scaffolds or pseudo-chromosomes, with N50 of 3 Mbp and L50 value of 4.39 Mbp. Except for scaffold 3, all scaffolds have telomere-to-telomere assembly, starting with repeats of “ACCCTA” and ending in complementary repeats of “TAGGGT” (Fig. 1) (Diotti et al. 2021). The benchmarking universal single-copy orthologs (BUSCO) (Manni et al. 2021) analysis of the whole genome with the Ascomycota fungal group revealed a completeness of 96.7%, with 1,650 complete BUSCOs found out of 1,706 total BUSCO groups searched. Of this, 1,649 were found to be complete and single copy, one was complete and had a duplicated BUSCO, while 52 were missing, and four fragmented BUSCO were reported. The BUSCO analysis of the transcriptome revealed a completeness of 97.9%, consisting of 1,669 complete BUSCO out of a total of 1,706. Among these, there are 1,668 complete single-copy BUSCO, one complete duplicated BUSCO, and four fragmented BUSCO, with 33 BUSCO reported as missing. The mitochondrial genome was also sequenced, yielding an assembly length of 1.7 Mbp with a single scaffold. A total of 54 open reading frames were annotated, including 14 core genes, 28 tRNAs, 4 rRNAs, and 8 hypothetical proteins.

Although a draft genome of B. fagacearum is available from Wingfield et al. (2016), this genome assembly is significantly improved (Table 1). Our genome assembly is gapless, with contigs and scaffolds being synonymous and extending telomere-to-telomere in all cases except for scaffold 3, in contrast to the 2,949 contigs and 1,257 scaffolds in Wingfield et al. (2016). Moreover, this assembly includes 270 Kbp more sequence than the previous genome.

TABLE 1 Parameters of available genomes of Bretziella fagacearum

Assembly parameters	Wingfield et al. (2016)	This study
Genome assembly size (Mbp)	26.8	27.07
Sequencing read coverage depth	123.2×	1,464.36×
Number of contigs	2,949	9
Number of scaffolds	1,257	9
Scaffold L50 (Mbp)	0.0422	4.39
Scaffold N50	190	3
Number of gaps (Mbp)	0.2	0
Number of proteins	8,713	7,554
Assembly level	Scaffold	Pseudo-chromosome
Sequencing technology	Illumina	PacBio

View as image HTML

This publicly available genome resource can serve as a reference for future studies, including population and comparative genomics, to explore molecular mechanisms underlying the epidemiology and biology of B. fagacearum. Insights gained from this genome will enhance the understanding of genetic diversity, evolutionary history, and pathogenicity, aiding in the development of effective disease management strategies.