Phylogenomics and Comparative Genomics Highlight Specific Genetic Features in Ganoderma Species

The Ganoderma species in Polyporales are ecologically and economically relevant wood decayers used in traditional medicine, but their genomic traits are still poorly documented. In the present study, we carried out a phylogenomic and comparative genomic analyses to better understand the genetic blueprint of this fungal lineage. We investigated seven Ganoderma genomes, including three new genomes, G. australe, G. leucocontextum, and G. lingzhi. The size of the newly sequenced genomes ranged from 60.34 to 84.27 Mb and they encoded 15,007 to 20,460 genes. A total of 58 species, including 40 white-rot fungi, 11 brown-rot fungi, four ectomycorrhizal fungi, one endophyte fungus, and two pathogens in Basidiomycota, were used for phylogenomic analyses based on 143 single-copy genes. It confirmed that Ganoderma species belong to the core polyporoid clade. Comparing to the other selected species, the genomes of the Ganoderma species encoded a larger set of genes involved in terpene metabolism and coding for secreted proteins (CAZymes, lipases, proteases and SSPs). Of note, G. australe has the largest genome size with no obvious genome wide duplication, but showed transposable elements (TEs) expansion and the largest set of terpene gene clusters, suggesting a high ability to produce terpenoids for medicinal treatment. G. australe also encoded the largest set of proteins containing domains for cytochrome P450s, heterokaryon incompatibility and major facilitator families. Besides, the size of G. australe secretome is the largest, including CAZymes (AA9, GH18, A01A), proteases G01, and lipases GGGX, which may enhance the catabolism of cell wall carbohydrates, proteins, and fats during hosts colonization. The current genomic resource will be used to develop further biotechnology and medicinal applications, together with ecological studies of the Ganoderma species.


Introduction
Ganoderma species (Ganodermataceae, Polyporales, Basidiomycota) are both ecologically and economically relevant fungi found in forest ecosystems. As wood-decay fungi, Ganoderma species efficiently decompose the components of plant cell walls, including lignin, cellulose and hemicellulose [1]. G. boninense Pat. causes a basal stem rot (BSR) on oil palm trees [2,3], while G. lucidum and G. australe have been used for degradation of environmental pollutants [4,5]. As medicines, Ganoderma fruiting bodies have been used over the last 2000 years to produce drugs used for improving immunity, and in anti-aging and anti-cancer treatments in humans [6][7][8]. For example, spore powder or basidiocarp slices of G. lingzhi Sheng H. Wu, Y. Cao, and Y.C. Dai and G. tsugae Murrill are used in Asian traditional medicine to improve health. On the other hand, polysaccharides and glycans extracted from G. sinense J.D. Zhao, L.W. Hsu and X.Q. Zhang have been used as clinical drugs [9].
Ganoderma is the largest genus in Ganodermataceae including 461 taxa recorded in Index Fungorum (http://www.indexfungorum.org/, accessed on 16 January 2022) as of 15 September 2021. However, only four complete genomes of Ganoderma species have been reported at this date. The genome of the well-known medicinal polypore, G. lucidum, was published by Chen et al. [10]. This study focused on the analyses of genes encoding cytochrome P450s (CYPs), transporters and regulatory proteins which are related to secondary metabolism (SM) and wood degradation. Then, Binder et al. [11] used the genome of an unknown North American Ganoderma species to unravel phylogenetic relationships within the Polyporales. The sequencing and analyses of G. sinense, another well-known medicinal fungus, provided new highlights on the genome methylation patterns, small RNA transcriptome, SM and defense processes [12]. Utomo et al. [13] sequenced and analyzed a pathogenic strain of G. boninense isolated from an oil palm tree with severe symptoms of BSR disease. Additional genomes are required to explore the full diversity of Ganoderma gene sets involved secondary metabolism and wood white-rot decay and identify the genetic traits specific to the different species. Therefore, we sequenced and analyzed three additional genomes of Ganoderma species and compared them to other previously sequenced taxa, i.e., G. australe, G. leucocontextum and G. lingzhi.
Ganoderma australe is widely distributed on deciduous trees in central, eastern and southern China, and G. lingzhi is widely cultivated thanks to its high medicinal value, whereas G. leucocontextum was recently described from southwestern China and is now cultivated in Yunnan and Tibet for medicinal production. In addition to these three new genomes, we included in our analysis seven publicly available genomes of Ganoderma species. Finally, we compared the Ganoderma gene sets to those obtained on 51 Basidiomycota species, including white-rot and brown-rot wood decayers, and biotrophic fungi. Here, we provided new insights on the evolutionary relationships between Ganoderma species and other Polyporales, and we identified idiosyncrasies in the gene repertoire involved in SM and secreted enzymes.

Strain and Culture Conditions
Fruiting bodies of Ganoderma leucocontextum Dai 12418 was collected on Quercus, and the strain was cultured on Potato Dextrose Agar (Potato Extract 200 g, Agar 20 g, Dextrose 20 g, water 1 L); G. australe Cui 17254 and G. lingzhi Cui 9166 were collected on other angiosperm logs., and the strains were cultured on Malt Extract Agar (Malt Extract 20 g, Agar 18 g, KH 2 PO 4 3 g, Glucose 10 g, water 1 L). The strains were cultured for 7-14 days in the dark at 25 • C. The taxonomic affiliation of the three strains was confirmed by phylogenetic analyses based on ITS sequences. Isolates were deposited in the culture collection of the Institute of Microbiology, Beijing Forestry University and are available upon request.

DNA Extractions
For genomic DNA extraction, fresh 7-day-old vegetative agar mycelial cultures grown on cellophane sheets were harvested, snap frozen in liquid nitrogen and stored at −80 • C. High molecular weight genomic DNA of Ganoderma australe was extracted by using the sodium dodecylsulfate (SDS) method [14], whereas a modified cetyltrimethyl ammonium bromide (CTAB) method [15] was used for genomic DNA extraction of G. leucocontextum and G. lingzhi.

Genome Sequencing and Assembly
The genome of the strain Dai 12418 from Ganoderma leucocontextum was sequenced using the Pacific Biosciences platform (PacBio, Menlo Park, CA, USA); 1.5 µg of genomic DNA was sheared to 10 kb using Covaris g-Tube or Diagenode megaruptor tube. The sheared DNA was treated with DNA prep to remove single-stranded ends and DNA damage repair mix followed by end repair/A Tail and ligation of barcoded overhang adapters using SMRTbell Express Template Prep 2.0 Kit (PacBio, Menlo Park, CA, USA). The library was purified with AMPure PB beads and libraries with different barcodes were pooled at equimolar for up to four maximum fungal genomes with a total sum of genome sizes of 200 Mb. A BluePippin size selection (Sage Science, Beverly, MA, USA) was then performed on the pool to remove the shorter peaks physically. PacBio Sequencing primer was then annealed to the SMRTbell template library and sequencing polymerase was bound to them using Sequel II Binding kit 2.0. The prepared SMRTbell template libraries were then sequenced on a Pacific Biosystems' Sequel II sequencer using 8 M v1 SMRT cells and Version 2.0 sequencing chemistry with 1 × 1800 sequencing movie run times. Filtered subread data was processed to remove artifacts. Mitochondria was assembled separately with the CCS reads using an in-house tool (assemblemito.py), used to filter the CCS reads, and polished with gcpp-algorithm arrow version SMRTLink v8.0.0.80529 (https://www.pacb.com/support/software-downloads, accessed on 16 January 2022). The mitochondria-filtered CCS reads were then assembled with Flye version 2.7.1-b1590 (https://github.com/fenderglass/Flye, accessed on 16 January 2022) (-g 40 M-asmcoverage 50-pacbio-corr) and polished with gcpp-algorithm arrow version SMRTLINK v8.0.0.80529. Contigs less than 1000 bp were excluded.
The Ganoderma australe Cui 17254 and G. lingzhi Cui 9166 genomes were sequenced using PacBio and Illumina platform. Sequencing libraries for Illumina were generated using NEBNext ® Ultra™ DNA Library Prep Kit (NEB, Ipswich, MA, USA). One µg of genomic DNA was fragmented by sonication to 350 bp (400 bp for G. lingzhi DNA) using Covaris g-Tube. The DNA fragments were end-polished, A-tailed, and ligated with the full-length adaptor for Illumina sequencing with further PCR amplification. The libraries were purified by AMPure XP system (Beckman Coulter, Brea, IN, USA) and were analyzed for size distribution by Agilent2100 Bioanalyzer and quantified using real-time PCR. The libraries for single-molecule real-time (SMRT) on PacBio platform was constructed with an insert size of 20 kb using the SMRT bell TM Template kit (version 1.0, Pacific Biosciences, Menlo Park, CA, USA). The DNA fragments were repaired DNA damage and ends and prepared blunt ligation reaction. The library was purified with 0.45X AMPure PB beads and size-selection using the BluePippin System (Sage Science, Beverly, MA, USA). The libraries were analyzed for size distribution by Agilent 2100 Bioanalyzer (Agilent, Santa Clara, CA, USA). The whole genome of G. australe was sequenced using PacBio Sequel and Illumina NovaSeq PE150 at Beijing Novogene Bioinformatics Technology Co., Ltd. (Beijing, China). The whole genome of G. lingzhi was sequenced using PacBio RS II platform and Illumina MiSeq platform at Shanghai Personalbio Technology Co., Ltd. (Shanghai, China). The lowquality reads were filtered (less than 500 bp) to obtain clean data. Preliminary assembly was conducted with SMRTLink v5.0.1 (https://www.pacb.com/support/software-downloads, accessed on 16 January 2022), and long reads (more than 6000 bp) were selected. By the variant Caller module of the SMRT Link software, the arrow algorithm was used to correct and count the variant sites in the preliminary assembly results.

Synteny Analyses
The synteny analyses were conducted on the 10 largest scaffolds of each of the selected Ganoderma genomes. Pair-wise comparisons and identification of syntenic blocks were performed by using the R package DECIPHER [22] with default parameters. The synteny blocks between every two species were visualized with the R package Circlize [23]. Data management, integration, and visualization were as described in Hage et al. [24].

Phylogenetic Analyses
The phylogenetic analyses of 12 Ganoderma species were conducted by Maximum Likelihood (ML) with 23 ITS sequences. Sanguinoderma sp. Cui 17238 was used as the outgroup. The ML analyses were performed in RAxML-HPC v. 8.2.3 [25] involving 1000 ML searches under the GTRGAMMA model and 1000 rapid bootstrap replicates with the GTRCAT model to obtain the best tree and ML bootstrap. All trees were viewed in FigTree 1.4.2 (http://tree.bio.ed.ac.uk/software/figtree/, accessed on 16 January 2022). The ML bootstrap ≥ 50% were presented on topologies from ML analyses.

Main Genome Features
Three Ganoderma species were newly sequenced in this study: Ganoderma australe strain Cui 17254, G. leucocontextum strain Dai 12418 and G. lingzhi strain Cui 9166. The size of the G. australe assembly was 84.27 Mb and 20,460 protein-coding genes were predicted ( Table 2). This is the largest Ganoderma genome sequenced to date. The size of the genome assemblies for G. leucocontextum and G. lingzhi was lower at 60.34 and 60.56 Mb, and 15,007 and 16,592 protein-coding genes were predicted on the assemblies, respectively. Between 73.1% to 99.8% of a benchmark set of conserved fungal genes (BUSCO) were found in genome assemblies, indicating that assembled genomes captured most of the coding gene space, although the gene annotation for G. lingzhi appeared to be more fragmented ( Table 2). We compared the three newly sequenced Ganoderma genomes to four published Ganoderma genomes and 51 other Basidiomycota genomes. The genomes size of the 58 investigated species ranged from 19.66 to 109.88 Mb (Figure 1), with 6785 to 26,226 predicted genes (Table S2). Except for Melampsora larici-populina, Ganoderma species displayed the largest genomes, i.e., G. australe and G. boninense were two-to three-fold larger than other wood decayers. No significant differences (p > 0.05) in the average genome size were found between white-rot fungi, brown-rot fungi and ectomycorrhizal fungi. The number of duplicated BUSCO genes was higher (23.4%) in G. australe. The synteny analyses showed no evidence for whole-genome or segmental duplications ( Figure S1), it may indicate a polymorphic dikaryon (Table S2). found between white-rot fungi, brown-rot fungi and ectomycorrhizal fungi. The number of duplicated BUSCO genes was higher (23.4%) in G. australe. The synteny analyses showed no evidence for whole-genome or segmental duplications ( Figure S1), it may indicate a polymorphic dikaryon (Table S2).

Macrosynteny between Ganoderma Genomes
The top 10 scaffolds of each Ganoderma genomes, covering 13% to 75% of the whole assemblies, were selected to perform a macrosynteny analysis ( Figure S2, Table S3A). We observed the highest percentage of syntenic segments (72%) between G. lingzhi and G. lucidum, while G. boninense and G. leucocontextum showed a lower rate when compared to other species (Table S3B), reflecting a higher genome divergence.

Phylogenetic Analyses of Ganoderma Strains
The phylogenetic relationship of 12 Ganoderma species was conducted based on 23 ITS sequences and one Sanguinoderma sp. as outgroup (Figure 2). Two G. lucidum strains were not clustered together indicating that the identification of G. lucidum G.260125-1 should be considered as G. lingzhi actually. Ganoderma sp. 10597 SS1 clustered with G. sessile suggesting that this strain may pertain to this sessile species in this phylogeny.

Macrosynteny between Ganoderma Genomes
The top 10 scaffolds of each Ganoderma genomes, covering 13% to 75% of the whole assemblies, were selected to perform a macrosynteny analysis ( Figure S2, Table S3A). We observed the highest percentage of syntenic segments (72%) between G. lingzhi and G. lucidum, while G. boninense and G. leucocontextum showed a lower rate when compared to other species (Table S3B), reflecting a higher genome divergence.

Phylogenetic Analyses of Ganoderma Strains
The phylogenetic relationship of 12 Ganoderma species was conducted based on 23 ITS sequences and one Sanguinoderma sp. as outgroup (Figure 2). Two G. lucidum strains were not clustered together indicating that the identification of G. lucidum G.260125-1 should be considered as G. lingzhi actually. Ganoderma sp. 10597 SS1 clustered with G. sessile suggesting that this strain may pertain to this sessile species in this phylogeny.

Phylogenomic Analyses of the Ganoderma Species and Related Polyporales
Our phylogenetic analyses, based on 143 single-copy conserved protein sequences, was in agreement with previous phylogenetic analyses using single or multi-locus approaches [24,63]. We also identified five major clades in Polyporales: the core polyporoid clade, the antrodia clade, the gelatoporia clade, the phlebioid clade, and the residual polyporoid clade. Ganoderma species clustered in the core polyporoid clade with a high bootstrap value (Figure 3, 100% ML bootstrap). The phylogenomic status of Ganoderma species was consistent with the macro-synteny conservation results.
The taxonomic status of Rigidoporus microporus has been changed to Hymenochaetales [64], and here, it was confirmed again by phylogenomic analyses. Among the sampled species, the brown-rot fungi and the ectomycorrhizal fungi formed monophyletic groups. The present analyses also confirmed that the white-rot lifestyle is evolutionary polyphyletic.

Phylogenomic Analyses of the Ganoderma Species and Related Polyporales
Our phylogenetic analyses, based on 143 single-copy conserved protein sequences, was in agreement with previous phylogenetic analyses using single or multi-locus approaches [24,63]. We also identified five major clades in Polyporales: the core polyporoid clade, the antrodia clade, the gelatoporia clade, the phlebioid clade, and the residual polyporoid clade. Ganoderma species clustered in the core polyporoid clade with a high bootstrap value (Figure 3, 100% ML bootstrap). The phylogenomic status of Ganoderma species was consistent with the macro-synteny conservation results.
The taxonomic status of Rigidoporus microporus has been changed to Hymenochaetales [64], and here, it was confirmed again by phylogenomic analyses. Among the sampled species, the brown-rot fungi and the ectomycorrhizal fungi formed monophyletic groups. The present analyses also confirmed that the white-rot lifestyle is evolutionary polyphyletic.

Transposable Element Identification
The dominant TEs in Ganoderma species belonged to the Gypsy and Copia families of long terminal repeats (LTR) retrotransposons. Besides, the proportion of unknown TEs in Ganoderma species was large, especially in G. australe (10.5% of the total assembly) and G. lingzhi (11.1% of the total assembly), indicating they likely played a key role in genome rearrangements. G. leucocontextum displayed the largest TE coverage (18.32% of the total assembly), and the most diverse TE distribution, including simple repeats (Repetitive), IS3EU and Dada of DNA transposon, RTEX and L1 of non-LTR retrotransposon which are unique repeat elements in Ganoderma species ( Figure 4A). The number of IS3EU and Dada sequences in G. leucocontextum was significantly larger comparing with other species (Figure 4B). Although occurring at moderate copy numbers, the G. boninense genome contained unique TE families, such as the DIRS LTR retrotransposons and several DNA repeated elements. G. australe also contained more copy of Tad1 non-LTR retrotransposon and Helitron DNA transposons by comparing to other Ganoderma species.

Transposable Element Identification
The dominant TEs in Ganoderma species belonged to the Gypsy and Copia families of long terminal repeats (LTR) retrotransposons. Besides, the proportion of unknown TEs in Ganoderma species was large, especially in G. australe (10.5% of the total assembly) and G. lingzhi (11.1% of the total assembly), indicating they likely played a key role in genome rearrangements. G. leucocontextum displayed the largest TE coverage (18.32% of the total assembly), and the most diverse TE distribution, including simple repeats (Repetitive), IS3EU and Dada of DNA transposon, RTEX and L1 of non-LTR retrotransposon which are unique repeat elements in Ganoderma species ( Figure 4A). The number of IS3EU and Dada sequences in G. leucocontextum was significantly larger comparing with other species (Figure 4B). Although occurring at moderate copy numbers, the G. boninense genome contained unique TE families, such as the DIRS LTR retrotransposons and several DNA repeated elements. G. australe also contained more copy of Tad1 non-LTR retrotransposon and Helitron DNA transposons by comparing to other Ganoderma species.
TE coverage in the 58 analyzed genomes ranged from 0.37% (Phlebiopsis gigantea) to 41.69% (Fomitiporia mediterranea) (Table S4). Comparing to other lifestyle fungi, white-rot fungi have low TE coverage except F. mediterranea (p < 0.01), among which Ganoderma species have relatively higher TE coverage (p < 0.05). Ectomycorrhizal fungi showed the larger repeat element coverage, in which unknown TE elements were the most abundant (>12%). TE coverage in the 58 analyzed genomes ranged from 0.37% (Phlebiopsis gigantea) to 41.69% (Fomitiporia mediterranea) (Table S4). Comparing to other lifestyle fungi, white-rot fungi have low TE coverage except F. mediterranea (p < 0.01), among which Ganoderma species have relatively higher TE coverage (p < 0.05). Ectomycorrhizal fungi showed the larger repeat element coverage, in which unknown TE elements were the most abundant (>12%).

Biosynthetic Gene Clusters
We identified a total of 16 types of SM biosynthetic clusters ( Figure 5, Table S5). Except for the Rhizopogon vinicolor, the pathogens, ectomycorrhizal fungi and endophytes generally had a lower number of genes involved in SM than saprotrophic species. This difference was mainly driven by a lower content in t1pks, terpenes related genes and miscellaneous genes tagged as "putative" and "other" by the antiSMASH software. The antrodia clade composed of brown-rot species was segregated from other Polyporales clades by an enrichment in t1pks, indole, and fatty acid associated genes, and a depletion in NRPS. Among white-rot fungi in Polyporales, the phlebioid clade displayed a higher content in t3pks and a lower content in terpenes associated genes. The core polyporoid clade contained the species with diverse sets of biosynthetic genes clusters. Gene clusters involved in terpene synthesis and unknown metabolites (i.e., putative SM gene clusters) were enriched in Ganoderma species, especially in G. australe which contained 31 terpene associated gene clusters and 80 putative biosynthetic gene clusters. Noteworthy, G. lingzhi and G. lucidum as the main medically relevant species showed the lowest number of terpenes associated gene clusters, whereas G. australe and the pathogenic G. boninense had the highest content. G. lucidum, G. leucocontextum, Trametes versicolor, and T. pubescens genomes encoded a cluster related to lantipeptide production, and the first three were identified as associated with a t1pks.

Pfam Protein Domains Found in the Genomes
More than 5000 Pfam protein domains were identified in the 58 selected genomes. A total of 32 Pfam categories were sorted according to their gene copy number (>100) (Figure 6, Table S6). Ganoderma australe showed the higher number of Pfam protein domains, i.e., five-fold larger than Ustilago maydis (with 864 protein domains). In G. australe, 13 protein domains were prominent, including cytochrome P450s (PF00067: p450) involved in SM, heterokaryon incompatibility protein (PF06985: HET) often related to vegetative incompatibility (VI) and membrane transporters of the major facilitator superfamily (PF07690: MFS_1).

The Predicted Secretome
Ganoderma species have the largest secretome among the 58 selected fungal species (Figure 7, Table S7). Small secreted proteins (SSPs) represented 41% to 81% of the total secreted proteins and 4% to 28% SSPs were annotated as CAZymes, lipases, or proteases. Most known SSPs were annotated as CAZymes, especially enriched in G. australe and G. lingzhi (Figure 8, Table S8).

The Predicted Secretome
Ganoderma species have the largest secretome among the 58 selected fungal species (Figure 7, Table S7). Small secreted proteins (SSPs) represented 41% to 81% of the total secreted proteins and 4% to 28% SSPs were annotated as CAZymes, lipases, or proteases. Most known SSPs were annotated as CAZymes, especially enriched in G. australe and G. lingzhi (Figure 8, Table S8).
As secreted CAZymes plays a key role in wood degradation, the evolution of their content through Polyporales was investigated. The different clades in Polyporales have a distinctive CAZyme content and the Ganoderma species also contained a distinctive CAZyme repertoire within the core polyporoid clade. The antrodia clade composed of brown-rot fungi have a two-fold lower content of secreted CAZymes compared to the white-rot species in Polyporales (141 ± 24 vs. 237 ± 71 respectively, Table S9). Among the 148 CAZymes sub-families (grouped in 101 families) with secreted genes representatives, 48 (grouped in 38 families) exhibited a lower number of genes in brown-rots compared to white-rots in Polyporales (BM test, FDR padj < 0.01). It included cellulolytic enzymes (GH6, GH7), LPMOs (AA9), ligninolytic PODs (AA1_1), heme-associated PODs (AA2), the carbohydrate-binding module CBM1, CAZymes involved in bacterial cell wall degradation (GH25, GH79), fungal cell wall degradation (GH20, GH76, GH92, GH128, GH135, GH152), and other CAZymes involved in plant cell wall degradation. CAZymes was globally enriched in the core polyporoid clade compared to the phlebioid clade. This enrichment was associated to the expansion of AA1, AA14, CE1, EXPN, GH16, GH17, GH18, GH25, and GH30 families. Of note, the AA7 family was absent from the phlebioid clade and AA12 was absent from the core polyporoid clade. Among the core polyporoid clade, the Ganoderma species have enriched CAZymes, and 11 families were expanded, including GH18, GH16, AA1, GH43, CE16, GH3, GH128, GH47, GH115, GH25, and GH1 families. Whereas only PL4 was enriched (and only found) in the other species of the core polyporoid clade. Compared to other species, G. australe had the largest set of secreted CAZymes, in which GH18 and EXPN families were highly enriched in this species. It also contained the largest repertoire of PCWDEs, MCWDEs, and enzymes acting on pectin, peptidoglycans, and chitin ( Figure 9, Table S10).   As secreted CAZymes plays a key role in wood degradation, the evolution of their content through Polyporales was investigated. The different clades in Polyporales have a distinctive CAZyme content and the Ganoderma species also contained a distinctive CAZyme repertoire within the core polyporoid clade. The antrodia clade composed of brown-rot fungi have a two-fold lower content of secreted CAZymes compared to the white-rot species in Polyporales (141 ± 24 vs. 237 ± 71 respectively, Table S9). Among the 148 CAZymes sub-families (grouped in 101 families) with secreted genes representatives, 48 (grouped in 38 families) exhibited a lower number of genes in brown-rots compared to white-rots in Polyporales (BM test, FDR padj < 0.01). It included cellulolytic enzymes (GH6, GH7), LPMOs (AA9), ligninolytic PODs (AA1_1), heme-associated PODs (AA2), the carbohydrate-binding module CBM1, CAZymes involved in bacterial cell wall degradation (GH25, GH79), fungal cell wall degradation (GH20, GH76, GH92, GH128, GH135, GH152), and other CAZymes involved in plant cell wall degradation. CAZymes was globally enriched in the core polyporoid clade compared to the phlebioid clade. This enrichment was associated to the expansion of AA1, AA14, CE1, EXPN, GH16, GH17, GH18, GH25, and GH30 families. Of note, the AA7 family was absent from the phlebioid clade and AA12 was absent from the core polyporoid clade. Among the core polyporoid clade, the Ganoderma species have enriched CAZymes, and 11 families were expanded, including GH18, GH16, AA1, GH43, CE16, GH3, GH128, GH47, GH115, GH25, and GH1 families. Whereas only PL4 was enriched (and only found) in the other species of the core polyporoid clade. Compared to other species, G. australe had the largest set of secreted CAZymes, in which GH18 and EXPN families were highly enriched in this species. It also contained the largest repertoire of PCWDEs, MCWDEs, and enzymes acting on pectin, peptidoglycans, and chitin ( Figure 9, Table S10).  A total of 3880 genes encoding for secreted proteases (in 73 MEROPS families) were identified in 58 fungal genomes (Table S11). Eight protease families (A01A, G01, M28E, M36, S08A, S10, S28, and S53) have more than 100 members, and the largest represented protease subfamily was A01A with 1332 proteins. No significant differences were found between saprotrophic and symbiotic species (p > 0.05), but five protease families (C40, M24X, M57, M77, S08B) were only found in white-rot fungi, although with low gene copy numbers. Among these families, the M57 protease family including four proteins, was found in Ganoderma australe, G. leucocontextum, G. lingzhi, and G. lucidum. In addition, Ganoderma species have more G01, M28E, and M35 proteases.
Comparing to other lifestyles, white-rot fungi have the higher set of secreted lipases (p < 0.001) and among the 58 analyzed fungal secretomes, Ganoderma australe contained the highest number of secreted lipases, the most abundant being GGGX lipases (Table S12). Most genomes displayed less than seven GX lipases, except for G. boninense and G. leucocontextum with ten and nine GX lipases respectively. They are prominent enzymes catalyzing a wide range of reactions on various cellular substrates [65,66].

Discussion
In this study, we provided three newly sequenced genomes of ecologically and economically relevant Ganoderma species. A phylogeny, based on 23 rDNA ITS sequences from 12 Ganoderma species, allowed us to determine the phylogenetic status of the newly sequenced species in Ganoderma. This phylogeny concurred with previous studies which were carried out by using a smaller set of Ganoderma species [11,63,67], except for G. lucidum. Indeed, this species displayed inconsistence between the two available strains, suggesting a misidentification of G. lucidum G.260125-1 strain which was purchased from a company. According to our study, this strain may pertain to the G. lingzhi which is widely cultivated for its medicinal usage in China [7]. This phylogeny also indicated that Ganoderma sp. 10597-SS1 probably pertained to G. sessile. The phylogenomic analyses based on 58 genomes, including seven Ganoderma species, supported the ITS phylogenetic result and confirmed that Ganoderma belong to the core polyporoid clade of Polyporales as defined by Justo et al. [63].
Ganoderma species have been used for centuries in traditional medicine thanks to their well-known arsenal of antimicrobial, anti-aging, antioxidant, anti-inflammatory, and immunomodulating compounds (e.g., polysaccharides, triterpenoids, and peptides) [68,69]. The drastic reduction of sequencing cost in the last decades allowed genome-wide mining of medicinal compound. As a result, the number of new genomes rapidly increased over the year, such as genomes of G. tsugae CCMJ4178 [67] and G. leucocontextum DH-8 [70]. The strain DH-8 of G. leucocontextum was newly sequenced and its genome was 50.05 Mb with 58 scaffolds [70], the genome of strain Dai 12418 of G. leucocontextum sequenced in our study was 60.34 Mb with 843 scaffolds. In our study, 21 terpene gene clusters were predicted in strain Dai 12418 against the 10 terpene gene clusters predicted in strain DH-8. The two-fold enrichment in strain Dai 12418 compared to strain DH-8 was unexpected, despite taking into differences in annotation tools and the potential bias of assembly quality. It may suggest a substantial intraspecific genome polymorphism. Direct comparisons of other functional gene categories were hampered by major differences in annotation methodologies. For example, Liu et al. [70] predicted 614 CAZymes genes using HMM, while we only found 291 CAZymes genes, using the in-house pipeline from the CAZymes database followed by expert manual curation.
One of the main targets for medically-relevant products are the secondary metabolites associated genes, especially the terpenes in Ganoderma species. Compared to other species, the terpene genes of Ganoderma species were indeed expanded, however, substantial differences were observed among this genus. G. lingzhi and G. lucidum as the main medically species have the lowest number of terpenes associated clusters whereas the widespread G. australe and the pathogenic G. boninense have the highest content. Terpene related gene expansion in G. australe and G. boninense might explain their ecological ability to develop on a broader set of substrates. On the other hand, a shift of the terpenome composition in G. lingzhi, G. tsugae, and G. leucocontextum might have led to the production of terpenes with beneficial properties for human health. This study also allowed us to identify a polyketide synthase complex involved in the synthesis of the antibacterial lantipeptide [71] in G. lucidum, G. leucocontextum, Trametes versicolor, and T. pubescens, likely explaining the use of formulations based on these fungi as antibiotics [10,[72][73][74]. Besides, the large content of gene clusters encoding unknown biochemical function(s) were found in Polyporales species, especially G. australe, suggested an outstanding ability to synthesize a large set of secondary metabolites of yet unknown function. This confirmed that the economically and medically relevant secondary metabolites of these fungi represent an untapped resource.
As Polyporales is an important group of wood decayers, we investigated the decay ability of these species by analyzing their secretome. The current evolution was posited that early Agaricomycetes were saprotrophic and different lifestyles were derived from it [75,76]. It also supported that the white-rot lifestyle is linked to the ability to degrade lignin. The acquisition of lignin degradation in the Agaricomycetes was estimated in the late Carboniferous and further evolved with the evolution of lignin complexity in Polyporales and Agaricales taxa [77]. The monophyletic antrodia clade, only composed of brow-rot species, was derived from white-rot lineages in current phylogenomic analyses. According to Baldrian and Valášková [75], brown-rot fungi, despite having several independent origins, were associated with the loss of ligninolytic PODs, heme dye-decolorizing PODs, heme-thiolate POD/peroxygenases (HTPs), cellulolytic enzymes (GH6, GH7, LPMO), and the carbohydrate-binding module CBM1 reduced their ability to degrade lignin into its partial modification by releasing Fenton-generated hydroxyl radicals in the colonized material. This statement was confirmed in the antrodia clade, in which those losses were also accompanied by losses in other PCW-degrading enzymes, and FCW-and BCWdegrading enzymes. The enhanced ability of white-rot fungi to decompose substrates was also marked by a secreted lipase enrichment. White-rots specific proteases (C40, M24X, M57, M77, S08B) were identified, and they were related to genetic regulation, bacterial cell-wall modification, nutrient transformation, and other cellular physiological functions.
The various white-rot Polyporales clades also displayed divergent wood degrading abilities. The potential ligninolytic ability and more generally PCW degradation ability of the core polyporoid clade was the highest among Polyporales due to the expansion of ligninolytic PODs (AA1) and hemicellulose degrading enzyme (AA14, CE1, GH30). The difference among white-rot clades also showed in the enzyme arsenal to compete with other microbes for substrates and copy with biotic threats. The specific expansion in BCW-degrading enzyme (GH25) and FCW-degrading enzyme (GH16, GH17, GH18) was also observed in the core polyporoid clade. Ganoderma species deepened their lignocellulose degradation ability with further expansion of the ligninolytic PODs (AA1), hemicellulose (GH43, CE16, GH115), BCW-(GH25), and FCW-(GH16, GH18, GH128, GH47) degrading enzymes. Aminopeptidase Ap1 (M28E), deuterolysin (M35) and scytalidoglutamic peptidase (G01) were also enriched in Ganoderma species, but further functional analyses are needed to clarify their role(s). Noteworthy, G. australe showed an even wider degrading gene repertoire including higher PCW-degrading CAZyme, lipase and protease content. Similarly, numerous cytochrome P450 monooxygenases (PF00067: p450) were found in G. australe. These enzymes are known for their role in lignin and xenobiotic degradation [78]. This pattern supported the known ability of G. australe to more efficiently decompose lignocellulose than other Ganoderma species. An increased CAZyme repertoire in BCW-and FCW-degrading enzymes was also found in G. australe, which could perform the higher competitiveness than other microorganism for the substrate during colonization on the hosts. Associated with its enriched number of biosynthetic genes and clusters, such as lytic transglycolase (PF03330) and cerato-platanins (PF07249), this repertoire could be used to combat with miscellaneous biotic threats.
The abundance of HET genes related to heterokaryon incompatibility, MFS transporters, and epigenetic regulation in Ganoderma australe indicated that these mechanisms play a key role in the stable reproduction and evolution of this wood decayer [79][80][81]. With 168 protein domains detected among seven Ganoderma species (PF00249, PF08914, PF11831, PF12776, PF13837, PF13921, PF15963), MYB transcription factors were also abundant in Ganoderma species. These transcription factors are the largest transcription factor families in eukaryotic organisms and play key role in variable development and physiological activities [82,83]. Wang et al. identified 75 MYB transcription factors in five Ganoderma species after manually curation, and of the gene copy number found in each species was lower than the detected results in our study. This is likely due to the differences between strains and technical methods used here. Further exploration of MYB genes can help to clarify its potential function during the growth and development of Ganoderma fungi.
The present comparative analysis of the publicly available Ganoderma genomes revealed a series of genetic features specific to this lineage of wood decayers. This study provided foundational information to characterize further ecological traits of this important group of decomposers. In addition, this information will be used to characterize the regulation of genes involved in SM biosynthesis pathways at the transcriptomic level, including antimicrobial compounds and medicinal drugs.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/jof8030311/s1, Figure S1: Synteny analyses of Ganoderma australe itself. The scatter plot shows correlation between all contigs, Figure S2: Macrosynteny between Ganoderma species. The percentage of syntenic hits occurring in the same order on the compared block of sequences between two genomes are shown inside each Circos plot, Table S1: Model list for each orthologous group of concatenate protein alignments, Table S2: Main features of the 58 genomes used in this study, Table S3: (A) Size and genome coverage of the 10 largest scaffolds in Ganoderma genomes, (B) synteny summary statistic of the 10 largest scaffolds in Ganoderma species, Table S4: Summary statistics of components and coverage for transposable elements (TEs) in the 58 genomes, Table S5: Number of the gene clusters involved in secondary metabolites biosynthesis identified in the 58 genomes, Table S6: Number of the filtered protein domains (abundance > 100) in Pfam database among the 58 genomes, Table S7: Count of each category (CAZymes, lipases, SSP, proteases) in 58 fungal secretomes, Table S8: Number of small secreted proteins in 58 secretomes, Table S9: Number of secreted CAZymes in 58 secretomes, Table S10: Composition and ratios of secreted CAZymes for plant and microbial cell wall degradation in 58 genomes, Table S11: Number of secreted proteases in 58 secretomes, Table S12