Analysis of the genome sequence of Phomopsis longicolla: a fungal pathogen causing Phomopsis seed decay in soybean

Background Phomopsis longicolla T. W. Hobbs (syn. Diaporthe longicolla) is a seed-borne fungus causing Phomopsis seed decay in soybean. This disease is one of the most devastating diseases reducing soybean seed quality worldwide. To facilitate investigation of the genomic basis of pathogenicity and to understand the mechanism of the disease development, the genome of an isolate, MSPL10–6, from Mississippi, USA was sequenced, de novo assembled, and analyzed. Results The genome of MSPL 10–6 was estimated to be approximately 62 Mb in size with an overall G + C content of 48.6%. Of 16,597 predicted genes, 9866 genes (59.45%) had significant matches to genes in the NCBI nr database, while 18.01% of them did not link to any gene ontology classification, and 9.64% of genes did not significantly match any known genes. Analysis of the 1221 putative genes that encoded carbohydrate-activated enzymes (CAZys) indicated that 715 genes belong to three classes of CAZy that have a direct role in degrading plant cell walls. A novel fungal ulvan lyase (PL24; EC 4.2.2.-) was identified. Approximately 12.7% of the P. longicolla genome consists of repetitive elements. A total of 510 potentially horizontally transferred genes were identified. They appeared to originate from 22 other fungi, 26 eubacteria and 5 archaebacteria. Conclusions The genome of the P. longicolla isolate MSPL10–6 represented the first reported genome sequence in the fungal Diaporthe-Phomopsis complex causing soybean diseases. The genome contained a number of Pfams not described previously. Information obtained from this study enhances our knowledge about this seed-borne pathogen and will facilitate further research on the genomic basis and pathogenicity mechanism of P. longicolla and aids in development of improved strategies for efficient management of Phomopsis seed decay in soybean. Electronic supplementary material The online version of this article (doi:10.1186/s12864-017-4075-x) contains supplementary material, which is available to authorized users.


Background
Phomopsis longicolla T. W. Hobbs (syn. Diaporthe longicolla) is a seed-borne fungus primarily causing Phomopsis seed decay (PSD) in soybean, Glycine max (L.) Merrill [1][2][3][4]. This disease decreases seed quality and has been found in most soybean production areas, worldwide [2,4,5]. The common symptoms of PSD include discolored seed that are both shriveled and elongated or cracked seed coats and are chalky-white in appearance. Soybean seed infected by P. longicolla often lack of visible symptoms or signs at harvest [6]. It has been reported that soybean seed infected by P. longicolla, whether symptomless or having symptoms, could have very low seed germination, reduced seedling vigor, and poor stands [4,7]. Poor seed quality of soybean could be due to the alteration of seed composition, reduction of oil quality, or moldy and/or split seed caused by P. longicolla [8]. PSD is one of the most economically important diseases of soybean. This disease has caused significant soybean yield loss. [9,10]. If the environment is warm and humid during the late growing season from pod fill through harvest, it favors pathogen growth and PSD development [11].
Management of PSD is very challenging. Inconsistent reductions of PSD have been reported when common agronomic practices were used. Practices included crop rotation with non-host or non-legume crops, conventional tillage to reduce pathogen inoculum, and prompt harvest when soybeans matured to avoid late season wet weather [2]. In addition, fungicide treatments could be used as an option to reduce PSD and other soybean diseases. However, they were not always effective in controlling PSD [12][13][14]. Planting cultivars with resistance to PSD is a long-term strategy to manage PSD. In past decades, most research conducted on the host resistance, such as identifying resistance sources by screening soybean germplasms, commercial cultivars, and breeding lines [5,[15][16][17], breeding for resistant lines and cultivars [18], and investigating inheritance of resistance to PSD [19][20][21]. In addition, genetic mapping of resistance to PSD was reported [22]. However, information about the genomic features and mechanisms underlying the pathogenicity of P. longicolla on soybean were lacking. It is well-known that plant cell walls are the primary barrier against pathogen invasions. In order to infect plants, a plant pathogen should have the ability to pass through the plant cell wall. Plant cell wall degrading enzymes (PCWDEs) are a subset of carbohydrate-activated enzymes (CAZy) that are produced by plant pathogens to degrade plant cell walls. There was no information about PCWDEs in P. longicolla. Further, horizontal gene transfer or lateral gene transfer has been inferred to be the movement of genetic material between different organisms [23,24]. If true it is a major force driving the evolution of both bacteria and eukaryotes [25,26]. To date there was no report inferring the possibility of HGT in P. longicolla. Understanding the nature of the pathogen and mechanisms of PSD development in host plants will help us develop better disease management strategies. In recent years, genomic studies have made important contributions to research and disease management in plant pathology [27]. The next-generation sequencing technology has facilitated the genomics-based approached to both improve disease resistance in crops and enhance our understanding the mechanism of pathogenicity. The genomic approaches could provide an alternative way to identify host resistances.
To facilitate investigation of the genomic basis of pathogenicity in P. longicolla and to understand the mechanism of the disease development, the genome of isolate MSPL10-6 was sequenced and de novo assembled [28]. This research was conducted to analyze the genome sequences of the P. longicolla isolate MSPL10-6. The aims here were to understand the genome features of P. longicolla, identify genes encoding plant cell wall degrading enzymes, discover and classify the repeat elements in the genome, and investigate the potentially horizontally transferred genes in the P. longicolla genome.

General genome features
The genome of the P. longicolla isolate MSPL 10-6 was assembled from both the paired-end and mate-pair libraries with the short oligonucleotide assembler package (SOAP), a denovo assembler. The oligonucleotides formed 108 scaffolds of 500 bases or larger. As reported in previous studies, the N50 length was 1,039,102 bp, and the largest scaffold contained 6,247,470 bp. The genome size was estimated to be approximately 62 Mb with an overall G + C content of 48.6% [28]. Statistics of genome sequencing and assembly are summarized in Table 1.

Gene prediction and annotation
Gene prediction analysis yielded a total of 16,597 genes (Average Length was 1704 bp, Total Length was 28,287,360 bp, Total Coding Length was 24,840,981 bp), of which 4334 genes where found to consist of a single exon (Average Length = 1219 bp). The total number of exons in all predicted genes was 47,213 (Average Length was 3622 bp, Total Length = 4,435,952 bp).
Of 16,597 genes predicted, 9866 genes (59.45%) had significant matches to genes in the NCBI nr database, while 18.01% of them did not link to any gene ontology (GO) classification. Further, 9.64% of the genes did not significantly match any known genes. Enzyme codes were assigned to 15.45% of the genes. The gene prediction statistics are summarized in Table 2. Functional categorization and distribution of potential genes in the P. longicolla genome are shown in Fig. 1.

Plant cell-wall degrading enzymes
The enzyme classification code (EC number) class distributions (level 3) are summarized in Fig. 2. Of 2674 EC enzyme-like orthologs identified, 1184 (44.3%) were related to hydrolases or hydrolytic enzyme, while 63 (2.4%) sequences were related to isomerases. An abundance of genes encoding plant cell-wall degrading enzymes (PCWDEs) were found in the P. longicolla genome. Of 1221 putative genes that encode carbohydrate-activated enzymes (CAZys), (Additional file 1: Table S1), 199 genes encoded carbohydrate esterases (CE), 471 encoded glycoside hydrolases (GH), and 45 encoded polysaccharide lyases (PL) ( Table 2). In the CAZy family, enzymes that have the same substrate and description could have different "domain" structures and coding sequences, such as CE1 -CE5, GH 10 and GH11 (Table 3). Three other classes of CAZys with indirect roles on degrading carbohydrates were auxiliary activities (AAs), carbohydrate-binding modules (CBMs), and glycosyl-transferases (GTs). The number of putative genes identified in the AA, CBM, and GT classes were 259, 113, and 134, respectively ( Table 4). Comparisons of the numbers of CAZys in P. longicolla with other Ascomycete fungi are shown in Fig. 3.

Horizontal gene transfers
A total of 510 potential horizontal gene transfers (HGTs) were identified in the genome of the MSPL 10-6 isolate (Additional file 2: Table S2). They were originally from 53 species including 22 fungi, 26 eubacteria and 5 archaebacteria ( Table 6). The majority of HGTs were from fungal origins (85.3%), while 13.3% and 1.4% of the HGTs were from eu-and archae-bacterial origins, respectively. Results of annotation of the HGTs based on gene ontology analysis are shown in Fig. 4. Over 70% of the HGTs were related to molecular functions.

Discussion
In this study, the general genome feature of P. longicolla isolate MSPL10-6 from Mississippi, USA was described. The assembly of the genome was the result of combining the output from analyzing the sequences from both paired end and mate pair libraries. A de novo nuclear genome assembly was generated and characterized. The P. longicolla genome was estimated to be approximately 62 Mb, using both kmer and read coverage analyses. The size of the genome of P. longicolla MSPL 10-6 isolate appeared larger than other reported ascomycete soybean pathogens, such as Diaporthe aspalathi (55 Mb) [29], Fusarium virguliforme (50.5 Mb) [30] and Macrophomina phaseolina (49.3 Mb) [31]. The overall number of predicted gene in P. longicolla was 16,597, while it was 14, 962, 14,845, and 14,249 in D. aspalathi [29], F. virguliforme [30], and M. phaseolina [31], respectively. It is unknown whether the larger size of the genome and the bigger number of the predicted gene contributed in part to the fungal specialization and the pathogenicity of P. longicolla on soybean. Both P. longicolla and D.  aspalathi are the members of Diaporthe-Phomopsis complex causing soybean diseases. They have similar culture morphology and very close relationships in taxonomy. However, P. longicolla is the primary cause of Phomopsis seed decay, while D. aspalathi primarily causes canker on soybean stem. It has been reported that soybean seeds, instead of stems and other soybean tissues, are more susceptible to P. longicolla than to other Diaphorthe species [14]. A comparative genome analysis is underway to address the questions about the pathogenicity mechanisms of those two species (Li et al., unpublished).
To investigate the genetic basis of pathogenicity, plant cell wall degrading enzymes (PCWDEs) in P. longicolla were identified and annotated. Results here indicated that P. longicolla contained abundance of genes encoding PCWDEs, which include all six classes of CAZys. For enzymes in the class of carbohydrate esterases (CE), the genome of P. longicolla contained 105 GH10 like genes that are absent in other pathogens, like F. virguliforme [32]. Plant cell walls contain abundant cellulose. Most of cellulose-degrading enzymes are classified into the glycoside hydrolase (GH) class. As indicated in Table 3, the substrates for the GH1, GH3, and GH5 degrading enzymes included not only cellulose, but also hemicellulose and pectin. It has been reported that most of biotrophic fungi do not have GH1 [33], while the genome of P. longicolla encodes six GH1 genes. Notably, there are 23 GH3 and 27 GH5 genes in the P. longicolla genome. Both GH3 and GH5 were common in pathogenic oomycetes, hemibiotrophic and necrotrophic fungi, which have more genes than biotrophic fungi and Pythium species [32]. The polysaccharide lyases (PL) are one of the important classes of PCWDEs. The PL class specializes in pectin degradation. It is wellknown that pectin is the most divergent component of plant cell walls with different modifications on the side chains. The microbial pectinolytic enzymes have been studied and reviewed [34]. PL1 and PL3 were the most common pectin lyases found in plant pathogens. In this study, there were 21 and 10 copies of PL1 and PL3 in P. longicolla, respectively. Three other PL members, PL9, PL 20, and PL22, have been thought previously to be unique to N. haematococca and F. virguliforme, and were not found in other plant fungal pathogens or oomycetes [32], but both PL9 and PL22 domains were present in the P. longicolla genome. F. graminearum [35], F. fujikuroi [36], and F. virguliforme [32], P. longicolla had PL1, PL3, PL4, PL9, and PL22. They all had PL11 and PL24, except F. virguliforme [32]. However, P. longicolla did not have PL20 as F. virguliforme had.
Significantly, a novel fungal ulvan lyase (PL24; EC 4.2.2.-) was found in the P. longicolla genome. Ulvan lyases degrade ulvan, an anionic polysaccharide. This enzyme has only been reported recently, as the first members of a new polysaccharide lyase family in bacteria [37]. Ulvan is the most abundant component of the green algal cell wall. The role of ulvan lyases in the pathogenicity of P. longicolla is unknown. It will require further investigation. The PCWDE is one of the most important factors associated with pathogenicity of fungal pathogens. It could play a crucial role in infecting plants and influencing host resistance. The list of PCWDEs identified in P. longicolla laid the foundation for dissecting the mechanisms of the fungal pathogenicity through further functional analyses of genes encoding PCWDEs. Those analyses should aid in   developing new strategies for breeding for resistance to Phomopsis seed decay in soybean.
To examine the genome architecture of P. longicolla, repetitive elements were analyzed. As noted, members of both class I and class II repetitive elements were found in the genomes of filamentous fungi [38]. Approximately 13% of the P. longicolla genome consists of repetitive elements. This is greater than the 10% and 1% in other ascomycetes Magnaporthe oryzae and M. poae genomes, respectively [39][40][41]. The majority (57.1%) of repetitive elements in the P. longicolla genome are categorized as Class I elements (retrotransposons). They are transcribed from DNA to RNA, and the RNA produced is then reverse transcribed into DNA. Moreover, retrotransposons mobilize via a "copy-and-paste", which allow for many copies to be inserted throughout the genome. Thus, retrotransposons are the most common transposon in eukaryotes [42,43] including P. longicolla. In addition, the long terminal repeat (LTR) is one of the main groups of retrotransposons [44]. The two main superfamilies of LTR retrotransposons found in fungi are Gypsy and Copia. In P. longicolla there were 30.7% of LTR/Copia and 24.4% of LTR/Gypsy among the repetitive elements. Since there are abundances of LTR/Copia present in the genomes of plants, fungi, animals, algae and several protists, it has been proposed that the ancestors of the LTR/Copia family probably co-existed with the ancestors of LTR/Gypsy before the separation between plants and other kingdoms [45].
The DNA transposons (class II elements) mobilize via a cut-and-paste mechanism that use a DNA intermediate, in which, the DNA itself is excised from the genome and integrated elsewhere. Mariner-like elements are one of prominent classes of the DNA transposons found in multiple species, including humans. This Class II transposable element is known for its uncanny ability to be transmitted horizontally between many species [46,47]. It estimated that there were 14,000 copies of mariner in the human genome encompassing 2.6 million base pairs [48]. The first mariner-element transposons outside of animals were found in Trichomonas vaginalis, the most common pathogenic protozoan infection of humans [49]. Interestedly, the most major transposon in the P. longicolla genome was DNA/TcMar-Fot1 (41.8%).  The function of this transposon in P. longicolla is unknown.
In the past, many documented HGT events inferred to fungi involved bacterial donors. For example, in a search for eubacterial-derived HGTs in 60 fully sequenced fungal species, Marcet-Houben and Gabaldon detected 713 transfer genes from bacteria [50]. Gene transfer between fungi has already been reported [51]. In a comparative genomic study of Fusarium species, four of F. oxysporum's 15 chromosomes inferred to have been acquired through HGT from a fungal source [52]. Notably, chromosome 14, which is essential for pathogenicity of tomato, could be transferred between pathogenic and non-pathogenic strains of F. oxysporum resulting in conversion of non-pathogenic strains into pathogenic strains. Here, the majority of HGTs in the P. longicolla genome were of fungal origin (85.3%), while only 13.3% of the HGTs were from bacteria. Almost half of the HGT genes in the P. longicolla genome were related to molecular function. Further research will be necessary to address many open questions such as the impact of HGTs on the genome structure, gene function, and pathogenicity of P. longicolla.
The genome of the isolate MSPL10-6 was the first reported genome sequence in the fungal Diaphothe-Phomopsis complex causing soybean diseases. Our study represents the first genomic effort to discover the genome structure of P. longicolla. The genome data provide new insights into the gene repertoire and physiological potential of seed-borne pathogens. Additionally, the genomic resources presented here, including the genome sequences and annotations, detail lists of cell wall degrading enzymes, repetitive elements and horizontal transferred genes, enhance our knowledge about the biology and genetics of P. longicolla. These discoveries will facilitate further research on the genomic basis and pathogenicity mechanism of P. longicolla, and aid in development of improved strategies for efficient management of this pathogen.

Conclusions
Phomopsis seed decay of soybean is one of the most devastating diseases affecting soybean seed quality worldwide. However, genomic basis and mechanism of the pathogenicity of P. longicolla on soybean was lacking. The draft genome of the P. longicolla isolate MSPL10-6 represents the first reported genome sequence in the fungal  Diaporthe-Phomopsis complex causing soybean diseases. The MSPL 10-6 genome contains a number of unique genomic features, including an abundance of genes encoding cell-wall degrading enzymes, numerous repetitive elements, as well as horizontal transferred genes from eubacteria, fungi and archaebacteria. Information obtained from this study enhances our knowledge about the biology and genetics of the seed-borne pathogen will facilitate further research on the genomic basis and pathogenicity mechanism of P. longicolla. The study will aid in development of improved strategies for efficient management of this pathogen.

Methods
Isolation, identification, and cultivation of P. longicolla isolate A P. longicolla isolate MSPL10-6 was isolated from field-grown soybean seed in Mississippi, USA in 2010 using the seed plating method. Briefly, over 100 randomly chosen soybean seeds that were harvested from the field were surface-disinfected in 0.5% sodium hypochlorite for 3 min, rinsed in sterile distilled water 3 times (3 min each time), and then placed on potato dextrose agar (Difico Laboratories, Detroit, MI) that was acidified (pH 4.8) with 25% lactic acid after autoclaving (APDA). Five seeds were plated on each 100 mm-diameter Petri dish. After 4 days of incubation at 24°C in the dark, putative/potential P. longicolla was isolated, streaked to the new APDA plates and incubated under 12-h lightand dark cycles. After 4-6 days, monoconidial cultures were obtained. Identification of P. longicolla was first based on morphological characteristics according to Hobbs et al. [1] and then was confirmed by analysis of the ITS region of rDNA amplified by PCR with primers ITS1, 5′-TCCG TAGGTGAACCTGCGG-3'and ITS4, 5′-TCCTCCGCT  Total TATTGATATGC-3′ [53] and the translation elongation factor 1-α gene primer set EF1-728F, 5′-CAT CGA GAA GTT CGA GAA GG -3′, and EF1-986R, 5′-TAC TTG AAG GAA CCC TTA CC -3′ [54,55]. Pathogenicity tests were performed using a cut-seedling inoculation method as described by Li et al. [56]. Isolate MSPL 10-6 was one of the most aggressive isolates causing severe soybean stem lesion in the greenhouse tests (data not shown). This isolate has also been used to screen soybean germplasm and successfully identified 23 new sources of resistance to PSD [17,57].

Genomic DNA extraction and sequencing
For DNA extraction, mycelial plugs (3-mm in diameter) from the margin of a 10-day old culture of MSPL 10-6 on APDA were cut and placed in potato dextrose broth (Difico Laboratories, Detroit, MI). After 4 days of incubation at 24°C under 12-h light-and dark cycles, mycelia were collected on sterile cheesecloth, washed with sterile water, immediately frozen with liquid nitrogen, and lyophilized with a freeze-drier (IMC Instruments, Inc., Wisconsin, USA). Fungal mycelia were ground with a mortar and pestle and pulverized in liquid nitrogen. The Genomic DNA of the P. longicolla MSPL10-6 isolate was used to generate sequencing libraries as previous described [28]. Briefly, paired-end libraries were made with the TruSeq DNA PCR-Free Sample Preparation kit (Illumina San Diego, CA), while the no-gel mate-pair libraries were generated with the Nextera Mate-Pair Sample Preparation kit (Illumina San Diego, CA) according to the manufacturer's protocols. All libraries were sequenced in separate lanes on an Illumina HiSeq 2500 sequencer using a TruSeq SBS sequencing kit (version 3, Illumina) at the Genomics Core Facility, Purdue University, West Lafayette, IN.

De novo genome assembly
Adapter sequences and poor quality bases (Phred score < 20) for each sequence read were trimmed using the FASTX-Toolkit (http://hannonlab.cshl.edu/ fastx_toolkit/index.html). A total of 72,216,734 mate-pair reads with a total of 8.2 billion bp representing 128-fold coverage, and 63,763,666 paired-end reads with a total of 6.2 billion bp, representing 97-fold coverage, were generated. The P. longicolla genome was assembled from both libraries using the software SOAPdenovo assembler version 2.04 [58]. Raw sequence data was deposited into NCBI's SRA database, under accession number: AYRD00000000 (1) (http://www.ncbi.nlm.nih. gov/nuccore/AYRD00000000).

Gene prediction and annotation
Gene prediction analysis was performed using a combination of homology searching and de novo prediction using Augustus web server [59,60] with complete gene option enabled and default for the rest of the parameters. F. graminearum [31] was used as the reference species due to the relatively close phylogenetic relationship to P. longicolla MSPL10-6 among the genome sequences available in GenBank.
Predicted genes were functionally annotated using Blast2GO [61]. The gene models were BLAST-ed (BLASTx) [62] against the NCBI non-redundant protein database. Then domain finding searches were done using InterProScan [63]. Enzyme codes and GO ontologies were then assigned to the gene models as described by the Gene Ontology Consortium (http://geneontology.org).
Identification of carbohydrate-activated enzymes (CAZys) in the P. longicolla genome To identify CAZys in the genome of P. longicolla isolate MSPL 10-6, Augustus program (http://bioinf.unigreifswald.de/webaugustus/prediction/create) [59] trained with the parameters of the species Fusarium graminearum was used to predict putative proteins of the CAZy family. Using the web resources of dbCAN, CAZy domains in the genome of P. longicolla were identified with a cutoff E value of 10 −3 [64]. Classification of CAZy was conducted as described by Chang et al., 2016 [32].
Identification and classification of repeat elements in the P. longicolla genome Repeat elements were identified using two well-cited software packages, RepeatMasker (v. 4.0.5) and Censor (v. 4.2.29) [65,66]. The analysis was partially carried out at the Bioinformatics Core at Purdue University. For Repeat-Masker, the optimized default parameters were utilized with the '-lib' option to find repeats associated with the RepBase file (v. 20.03) above. Additionally, the '-species Fungi' option was also used in a separate analysis to find fungal repeats based on RepBasederived libraries provided by RepeatMasker. Similarly, for Censor, with the 'censor.ncbi' module, the optimized default parameters were utilized and the '-lib' option to find repeats.
The predicted protein sequences were compared against the NCBI non redundant database using BLASTP to find the top hits. The genes from contigs larger than 500 bp were used since the same contigs were used for repeat prediction. Similarly, the predicted protein sequences were compared against protein sequences of the closest relative (Diaporthe ampelina) by creating a BLAST-able database. The combined annotation file containing the annotation of each predicted genes and D. ampelina database was generated using in house scripts.

Additional files
Additional file 1: Table S1. A list of putative genes encoding cell wall degrading enzymes. (TXT 73 kb) Additional file 2: Table S2. A list of horizontal transferred genes identified in the genome of Phomopsis longicolla isolate MSPL 10-6. (XLSX 75 kb)