Aggressiveness, growth rate and OA production
The typical symptoms caused by S. rolfsii strains ZY and GP3 on the peanut stems included unrestricted lesions on infection sites followed by tissue maceration, finally partial plant even whole plant wilting. Disease severity was scored at 14 days after inoculation (dpi) and disease index showed a significant difference between these two strains. The disease index of ZY was 82.34, which was classified as highly aggressive. The disease index of GP3 was 32.2, which was regarded as weakly aggressive (Figure 1a, 1b). The growth rate of these two strains was similar on PDA plate and showed no significant difference (Figure 1c, 1d). There was no significant difference in the amount of oxalic acid (OA) produced by these two strains either by haloes revealing on the PDA plate containing bromophenol blue, or by OA amount in the culture filtrate as analyzed by KMnO4 titration after grown in PDB medium without shaking for 7 days (Figure 1e, 1f).
Genome sequence and assembly
A total of 9, 97 Gb with 8.80 kb average subreads was generated for ZY and 6.34 Gb with 10.68 kb mean subreads for GP3 by SMRT sequencing. After polishing with Illumina data, the assembled genomes of GP3 and ZY were 70.51 Mb and 70.61 Mb, respectively, containing 27 contigs with an N50 length of 3.67 Mb for GP3, and 23 contigs with an N50 length of 3.71 Mb for ZY (Table 1). The two strains had genome assemblies of a similar size, both slightly smaller than that of S. rolfsii Indian strain M10 (73.18 Mb) [20]. The completeness of the genome assemblies was assessed using BUSCO [21]. About 97.5% (1301/1335) and 97.2% (1298/1305) of gene groups required for the correct assembly of Basidiomycota were present in GP3 and ZY, respectively (Figure S1). The average GC contents of the resulting S. rolfsii genomes of GP3 (46.27%) and ZY (46.29%) were comparable to S. rolfsii M10 (46.16%) (Table 1). Gene candidates in the S. rolfsii GP3 and ZY genomes were predicted by a combined approach, and 17,097 and 16,743 genes with an average gene length of 2,013.91 bp and 2,039.76 bp were identified (Table 1). Approximately 93.27% (15,947) of GP3 genes and 93.93% (15,727) of ZY genes could be annotated by non-redundant nucleotide and protein sequences in the Cluster of Orthologous Groups (KOG), Gene Ontology (GO), Kyoto Encyclopedia of Genes and Genomes (KEGG), Non-redundant Protein (NR), and Swiss-Prot databases (Figure S2–5). The number of genes predicted in S. rolfsii strains GP3 and ZY was similar to that in S. rolfsii strain M10 (16,830 genes) (Table 1). In this study, we identified 356 tRNAs, 48 rRNAs and 32 snRNAs in the S. rolfsii GP3 genome, and 415 tRNAs, 55 rRNAs and 32 snRNAs in the S. rolfsii ZY genome (Table S1). Comparison of gene orthologous with nine Agricomycetes fungi by OrthoMCL [22], GP3 and ZY shared a similarly low number of unique genes with 75 for GP3 and 37 for ZY distributed in 62 and 19 gene families (Table S2), respectively. Sequence comparison between contigs of two whole-genome assemblies indicated good macrosynteny between GP3 and ZY. Especially, contigs3, 7, 10, 15, 16, and 17 of GP3 corresponded well with contigs 1, 6, 14, 18, 15, and 16 of ZY (Figure 2).
Repetitive element analysis
De novo and homology approaches were combined to identify repetitive sequences in the S. rolfsii GP3 and ZY genomes. A total of 14.75% and 14.66% repetitive sequences were generated for GP3 and ZY, respectively (Table 1, Table S3). The abundance of repetitive sequences was similar between two strains and much more than that of S. rolfsii strain MR10, which had a repetitive sequence content of 3.73% (Table 1). GP3 and ZY contained repetitive elements including DNA transposons, retroelements, and satellites. Retroelements were abundant in the studied genomes, accounting for 10.28% and 10.79% in GP3 and ZY. LTR was abundant in the retroelements, accounting for 9.85% and 10.38% in GP3 and ZY (Figure 3, Table S4). Both abundance of LTR elements and retroelements in repetitive sequences were also found in S. rolfsii MR10 genome (Table S4).
Orthology analysis and phylogenetic analysis
The entire sets of predicted proteins of S. rolfsii GP3 and ZY were clustered with the OrthoMCL program [22] to identify gene families. Comparative analysis of the genomes of related species of Agaricomycetes, Basidiomycota showed that S. rolfsii strains had larger genomes but fewer total genes in comparison to most of the other species (Figure S6). Of gene families, the unclustered genes number of GP3 and ZY were the least among fungi in Agaricomycetes. A Venn diagram of the OrthoMCL revealed that S. rolfsii strains shared 4813 genes with other four Agaricomycetes species (Figure 4a).
To understand the genetic relationship of GP3 and ZY to the related Agaricomycetes species, we generated a phylogenetic tree of single-copy genes based the orthologous gene family analysis between two S. rolfsii strains and other Agriocomycetes fungi, including Armillaria gallica, Auricularia subglabra [23], Exidia glandulosa [24], Galerina marginata, Gymnopus luxurians [25], Hydnomerulius pinast [25], Psilocybe Cyanescens [25], Scleroderma citrinum [25], and Piloderma croceum [25]. The phylogenetic tree indicated that S. rolfsii strains were more closely related to E. glandulsa and A. subglabra, which belonged to Auriculariale, than to P. croceum, which belonged to Atheliaceae, the same as S. rolfsii (Figure 4b).
Genes involved in pathogenicity
Homologs in PHI base
In total, we identified 4,600 and 4,603 potential pathogen-host interaction (PHI) genes by searching the PHI base (Figure 5). Among them, 24 genes were annotated as effector and 172 genes were annotated as increased virulence in GP3, while ZY had 25 effectors and 138 genes related to increased virulence. Compared with S. rolfsii GP3, a total of 45 genes were unique in ZY, two of which annotated as effector and one annotated as increased virulence. We also found 58 genes in GP3 were not present in the ZY genome, 12 and 18 of which was annotated as loss of pathogenicity and reduced virulence, respectively (Table S5).
CAZymes
The genomes of S. rolfsii GP3 and ZY contained 957 and 925 genes encoding putative CAZymes, distributed in 118 and 119 CAZyme families. Glycoside hydrolases (GH) were dominant in the GP3 and ZY genomes (51.62% and 52.54%), followed by carbohydrate-binding modules (CBMs) and glycosyltransferases (GTs) (Figure 6a). The CAZyme content of GP3 was slightly larger than that of ZY, and CAZyme content of both GP3 and ZY was more than that of S. rolfsii MR10 (902) (Figure 6b).
Comparison of CAZyme content of S. rolfsii strains with other plant pathogens including six necrotrophic fungi (Aspergillus niger, Bortytis cinerea, Pecillium digitatum, Sclerotinia sclerotium, Rhizcotonia solani, and Verticillum dahliae), three hemibiotrophic fungi (Colletotrichum higginsianum, Fusarium granimearum, and Magnaporthe oryzae), and three biotrophic fungi (Puccinia graminis, Peronospora effusa, and Ustilago maydis) showed that the CAZyme content of S. rolfsii genome was the highest among above analyzed pathogens (Figure 6b). Necrotrophic fungi had more CAZymes than biotrophic and hemibiotrophic fungi. In comparison with other necrotrophic plant pathogens with a broad host range, such as S. sclerotinia, B. cinema, and V. dahliae, the CAZyme content of S. rolfsii was much more than these fungi. Compared to Basidiomycota plant pathogens, CAZyme content of S. rolfsii was three times as much as R. solani and P. graminis, and four times as much as U. maydis (Figure 6b). Besides differences in CAZyme content, the number of CAZymes involved in cellulose, hemicellulose, and pectin degradation of S. rolfsii strains GP and ZY was also different from that of other analyzed pathogens (Tables S6–S8). The number of those CAZymes in S. rolfsii strains was noticeably larger than that in other analyzed pathogens, especially in the pectin degrading capacity.
Glycoside hydrolases are known to catalyze the hydrolysis of glycosidic bonds in carbohydrate molecules. S. rolfsii was enriched in one glycosyl hydrolase family, GH28, a class of polygalacturonases involved in pectin degradation. The amount of GH28 was the same in GP3 and ZY (62 vs 62) and was significantly larger than the amount found in the other analyzed pathogen species (Table S8). The expansion of GH28 was not found in the biotrophic and some hemibiotrophic pathogens, such as U. maydis, P. graminis, and M. oryzae. In comparison with the other analyzed necrotrophic pathogenic fungi, S. rolfsii strains had three times more GH28 than those pathogens. Besides GH28, some other glycoside hydrolases involved in pectin degradation in S. rolfsii, such as GH35, GH51, and GH78, also had higher number in comparison to other pathogen species (Table S8).
Secretome and effector
The putative secreted proteins of S. rolfsii GP and ZY were identified based on a comprehensive pipeline (Figure S7). The genomes of GP3 and ZY were predicted to encode 536 (3.14%) and 551(3.29%) secreted proteins, respectively. Among the secreted proteins candidates, there are 151 and 30 secreted CAZymes genes for ZY and GP3, including 113 GH, 20 CE, 15 CBM, and 3 PL genes for ZY, while 22 GH, 6 CE, and 2 CBM genes for GP3 (Figure 6c, Table S9). In compare to secreted CAZymes involved in cellulose, hemicelluloses and pectin degradation, ZY had more of these genes than GP3 (Figure 6d).
A total of 50 and 46 putative effector candidates for GP3 and ZY, respectively, were predicted by Effector P.1. After manual inspection with the criteria of 50≤ molecular weight ≤ 300 kDa, 0–1 predicted trans-membrane domain, and ≥ 4 cysteine residues, a total of 30 and 27 effector candidates for GP3 and ZY were identified (Table 2). Most of the putative effector candidates were small (average length of 146 and 152 amino acids, ranging from 52 to 278, and 58 to 291 amino acids for GP3 and ZY). These candidates were rich in cysteines (the average cysteine composition was 8.5% (GP3) and 8.6% (ZY)). The functions of most effector candidates (73.33% and 44.44% of GP3 and ZY) were unknown. Comparison of putative effectors with PHI and CAZymes candidate genes showed that the number of genes, for functional effector, loss of pathogenicity, reduced virulence, GH, and CBM, differed between these two strains. ZY had two effectors and five GH genes, while GP3 had one GH gene and no effector overlapping with PHI and CAZymes candidate genes (Table S10). The function of these predicted effectors needs to be further verified in future research.
Secondary metabolites
The antiSMASH 4.0 software was used to identify the secondary metabolite gene clusters in the genome sequence of S. rolfsii ZY and GP3. A total of 46 and 31 gene clusters were predicted to be related to secondary metabolism in ZY and GP3, respectively (Figure 7). In ZY, two clusters containing genes encoding non-ribosomal peptide synthase (NPRS) were identified. Three, one, and 12 clusters were predicted as Type I polyketide synthase (T1 PKS), NPRS/ T1 PKS, and terpene, respectively. Besides, 28 clusters were predicted as others. Compared to ZY, GP3 contained no NPRS cluster, the same number of NPRS/ T1 PKS clusters, two fewer T1 PKS clusters, three fewer terpene clusters, and 8 fewer other clusters (Figure 7).