Aspergillus parasiticus SU-1 Genome Sequence, Predicted Chromosome Structure, and Comparative Gene Expression under Aflatoxin-Inducing Conditions: Evidence that Differential Expression Contributes to Species Phenotype

ABSTRACT The filamentous fungi Aspergillus parasiticus and Aspergillus flavus produce the carcinogenic secondary metabolite aflatoxin on susceptible crops. These species differ in the quantity of aflatoxins B1, B2, G1, and G2 produced in culture, in the ability to produce the mycotoxin cyclopiazonic acid, and in morphology of mycelia and conidiospores. To understand the genetic basis for differences in biochemistry and morphology, we conducted next-generation sequence (NGS) analysis of the A. parasiticus strain SU-1 genome and comparative gene expression (RNA sequence analysis [RNA Seq]) analysis of A. parasiticus SU-1 and A. flavus strain NRRL 3357 (3357) grown under aflatoxin-inducing and -noninducing culture conditions. Although A. parasiticus SU-1 and A. flavus 3357 are highly similar in genome structure and gene organization, we observed differences in the presence of specific mycotoxin gene clusters and differential expression of specific mycotoxin genes and gene clusters that help explain differences in the type and quantity of mycotoxins synthesized. Using computer-aided analysis of secondary metabolite clusters (antiSMASH), we demonstrated that A. parasiticus SU-1 and A. flavus 3357 may carry up to 93 secondary metabolite gene clusters, and surprisingly, up to 10% of the genome appears to be dedicated to secondary metabolite synthesis. The data also suggest that fungus-specific zinc binuclear cluster (C6) transcription factors play an important role in regulation of secondary metabolite cluster expression. Finally, we identified uniquely expressed genes in A. parasiticus SU-1 that encode C6 transcription factors and genes involved in secondary metabolism and stress response/cellular defense. Future work will focus on these differentially expressed A. parasiticus SU-1 loci to reveal their role in determining distinct species characteristics.

A spergillus parasiticus and Aspergillus flavus are the predominant aflatoxin producers on susceptible crops (1). The A. flavus NRRL 3357 (3357) complete genome sequence is available (http://www.aspergillusflavus.org/genomics/), and the sequence and location of the aflatoxin biosynthetic gene cluster were determined as part of this analysis. The nucleotide sequence of the A. parasiticus SU-1 aflatoxin biosynthetic gene cluster is also available (2), and A. flavus 3357 and A. parasiticus SU-1 exhibit a high degree of sequence identity between the aflatoxin genes as well as between other genes that have been sequenced in both organisms.
Although it is clear that A. flavus and A. parasiticus are closely related species, several phenotypic characteristics have been used to differentiate these species; these phenotypic characteristics include the quantity and type of mycotoxins synthesized in culture and on plants. For example, A. flavus strains synthesize predominantly B aflatoxins, while A. parasiticus strains produce both Band G-type aflatoxins (1,3). Less than 50% of A. flavus isolates are reported to produce aflatoxins, while almost all strains of A. parasiticus are toxigenic (4). A. parasiticus strains generally synthesize greater quantities of aflatoxin than A. flavus strains do (1,3). A. flavus strains are also reported to synthesize cyclopiazonic acid (CPA), while A. parasiticus strains are not known to synthesize this mycotoxin (1).
A. parasiticus and A. flavus also exhibit morphological differences, and gross mycelial structure and conidiospore morphology have been used to help assign species designations to isolates in the genus Aspergillus (5). A. flavus and A. parasiticus can grow and produce aflatoxin on a wide variety of susceptible crops. A. flavus is the most widely reported fungus associated with food, and it appears to have a broader host range because it is strongly associated with contamination of corn, cotton, peanuts, and tree nuts, while A. parasiticus is primarily associated with contamination of peanuts (4). Both species share an ecological niche in the soil, but it is not clear how well these organisms compete there.
To identify the genetic basis for reported differences in species morphology and biochemistry, we conducted Illumina next-generation sequence (NGS) analysis on the A. parasiticus SU-1 genome. We compared the A. parasiticus SU-1 genome sequence and predicted chromosome structure with the available A. flavus 3357 genome and compared expression (RNA sequence analysis [RNA Seq]) of all detectable A. parasiticus SU-1 and A. flavus 3357 genes under identical aflatoxin-inducing (YES growth medium, 40 h of growth; standard conditions, see Materials and Methods) and aflatoxin-noninducing (YEP medium, 40 h) growth conditions. We observed that the genome sequences of A. parasiticus SU-1 and A. flavus 3357 exhibit a high degree of identity, and the chromosome structure and gene distribution are predicted to be similar. Despite these similarities, we observed significant differences in the presence of specific secondary metabolism genes and entire gene clusters as well as differential expression of secondary metabolism and other genes between the two species that may contribute to species differences in phenotype. In addition, we identified three sets of uniquely expressed genes in A. parasiticus that encode fungus-specific zinc binuclear cluster (C6) transcription factors as well as enzymes involved in secondary metabolism and stress response/cellular defense. Identification of these differentially expressed genes prompts further study of the genetic basis for unique phenotypes observed in these closely related fungal species.

MATERIALS AND METHODS
Microorganisms and culture conditions. Aspergillus parasiticus strain SU-1 (ATCC 56775) and A. flavus strain NRRL 3357 (3357) were used as the wild-type aflatoxin-producing strains. YES liquid medium (2% yeast extract, 6% sucrose [pH 5.8]) (6) was used as an aflatoxin-inducing rich growth medium, and YEP liquid medium (2% yeast extract, 6% peptone [pH 5.8]) (7) was used as an aflatoxin-noninducing medium. An aliquot of frozen conidiospores from each strain was used to inoculate 100 ml of these growth media (250-ml Erlenmeyer flask carrying 5-or 6-mm glass beads) at 10 4 spores per ml, and these cultures were incubated for 40 h in the dark at 30°C with shaking (150 rpm) as described previously (standard growth conditions) (7).
DNA preparation and Illumina NGS analysis. Genomic DNA was prepared from A. parasiticus SU-1 grown for 44 h in YES medium under standard conditions using the method of Cihlar and Sypherd (8) as modified by Horng et al. (9). Genomic DNA was subjected to Illumina nextgeneration sequence (NGS) analysis by Functional Biosciences (Madison, WI), and bioinformatics analysis was conducted by ContigExpress, LLC (New York, NY). Quality control checks on the raw sequencing data (80 million reads in the forward and reverse directions) were performed using FastQC version 0.10.1. K-mer error correction was performed on the raw reads using Quake version 0.3.4. Paired-end reads were extracted from the corrected read pool, and the remaining reads were deposited as single-end reads. Both paired-end and single-end reads were used in de novo genome assembly. SOAPDenovo version 1.05 was used to perform assembly of the error-corrected reads. Standard parameters for paired-end reads were used. The K-mer that generated scaffold sequences with a maximal N 50 was 67. The resulting scaffolds were gap filled using corrected paired-end reads. The final scaffold sequences consisted of 375 gap-filled scaffold sequences larger than 1 kb (see A. parasiticus FASTA scaffolds in the supplemental file in the supplemental material).
Comparison of chromosome structure and gene distribution in A. parasiticus SU-1 and A. flavus 3357. The 31 largest scaffolds from the A. parasiticus SU-1 de novo genome sequence assembly were compared to the 16 largest A. flavus 3357 scaffolds that composed the arms of the eight chromosomes in this filamentous fungus (http://www.aspergillusflavus .org/genomics/). The relative position of each of the 31 A. parasiticus SU-1 scaffolds was assigned by direct sequence comparison using the 16 largest A. flavus 3357 scaffolds as a reference. These position assignments were supported by examination of the relative locations of 43 secondary metabolite gene clusters that were shared by A. flavus 3357 and A. parasiticus SU-1. The locations of A. parasiticus SU-1 open reading frames on A. parasiticus SU-1 scaffolds and chromosomes also were completed based on direct sequence comparison of the largest A. flavus 3357 and A. parasiticus SU-1 scaffolds.
antiSMASH (antibiotics and secondary metabolite analysis shell) software (10) was used to identify secondary metabolite loci within the A. flavus 3357 and A. parasiticus SU-1 genomes. FASTA formatted files of the A. flavus 3357 genome sequence (NCBI accession number AAIH00000000) and the A. parasiticus SU-1 de novo genome assembly were uploaded onto the web version of antiSMASH (http://antismash .secondarymetabolites.org/) with default settings. antiSMASH output displays a list of identified clusters based on the chemical backbones of their predicted products.
The predicted function of the uniquely expressed genes was determined by ContigExpress, LLC, New York, NY. Raw FASTA reads of 213 uniquely expressed transcripts in A. parasiticus SU-1 were translated in six frames and submitted as BLAST queries using NCBI nonredundant protein sequences. Function and Gene Ontology (GO) assignment was conducted using BLAST2GO. A subsequent BLAST of the translated transcript sequences was compared against a Kyoto Encyclopedia of Genes and Genomes (KEGG) database and assigned KEGG gene identification (ID) and pathway analysis using KOBAS. Function, GO, and functional pathway analysis were compiled in an Excel file.
RNA preparation and RNA Seq analysis. Total RNA was isolated from two biological replicates of A. parasiticus SU-1 and two biological replicates of A. flavus NRRL 3357 grown under standard conditions in YES medium for 40 h using the TRIzol method (TRIzol reagent; Invitrogen, Carlsbad, CA) as described by the manufacturer. Total RNA was also isolated from a single sample of A. parasiticus SU-1 grown in YEP medium under standard conditions for 40 h. The integrity and purity of the RNA samples were assessed using an Agilent 2200 TapeStation bioanalyzer system (Agilent, Santa Clara, CA). All samples had an RNA integrity number (RIN) of 9.0 or higher with an optical density at 260 or 280 nm (OD 260/280 ) between 1.78 and 1.80. Illumina RNA sequence (RNA Seq) analysis was  conducted on purified total RNA samples by Otogenetics (Norcross, GA), and bioinformatics analysis of the resulting data was conducted by ContigExpress, LLC (New York, NY). cDNA library construction was conducted as follows. One hundred nanograms of total RNA was reverse transcribed into cDNA using a Clontech SmartPCR cDNA kit (Clontech, Mountain View, CA). A Fragment Analyzer (Advanced Analytical Technologies, Inc. [AATI], Ames, IA) and Qubit (Life Technologies, Grand Island, NY) were used to assess cDNA quality. Restriction enzyme digestion end repair was conducted by removing adaptor sequences, and fragmentation of the resulting cDNA was conducted using a Covaris M220 focused ultrasonicator (Woburn, MA). Fragmented cDNA was subjected to Illumina library preparation with NEBNext-based reagents (New England BioLabs [NEB], Ipswich, MA). The quality, quantity, and size distribution of Illumina libraries were determined using an Agilent Bioanalyzer 2100. cDNA clusters were generated on a cBot automated station, and the Illumina libraries were submitted for Illumina HiSeq 2000 100-bp paired-end sequencing (Illumina, San Diego, CA).
Data analysis workflow was developed by ContigExpress, LLC (New York, NY). FASTQC (Babraham Institute, Cambridge, United Kingdom) was used to assess quality of paired-end 100-bp Illumina raw sequence reads. The numbers of raw reads per RNA sample were as follows: A. parasiticus SU-1 in YES medium, replicate 1 (35,870,008); SU-1 in YES medium, replicate 2 (23,352,926); A. flavus 3357 in YES medium, replicate 1 (21,887, 60); 3357 in YES medium, replicate 2 (28,694,452); SU-1 in YEP medium (22,730,224). Processed reads were then mapped to the de novo scaffold sequences of the A. parasiticus SU-1 genome with paired-end configuration using Bowtie2, version 2.0.0-beta7 (11). The A. flavus 3357 genome (http://www.aspergillusflavus.org/genomics/) was used to obtain A. flavus cDNA gene designations which were compared to the scaffold sequences using BLAST, version 2.2.24ϩ (12). The matched scaffold region of the top hit of each A. flavus gene, when available, was used to annotate that region with an A. flavus gene ID. The annotation information was assembled in the GTF format for downstream analysis. The mapping files of three A. parasiticus samples and two A. flavus samples were used to assemble the master gene annotations. The BLAST-based gene annotation was submitted to Cufflinks (version 2.0.2) as the reference annotation (13). Differential gene expression analysis was performed with the master gene annotations using Cufflinks.
Identification and function of "uniquely expressed loci" in A. parasiticus. We compared the level of gene expression in the two A. flavus 3357 biological replicates grown for 40 h on YES medium under standard conditions versus the two A. parasiticus SU-1 biological replicates grown for 40 h under the same conditions. We compiled a list of 1,409 genes that exhibited significant levels of differential expression between the two fungal species grown for 40 h on YES medium (q value Ͻ 0.05). Of these 1,409 genes, 318 were not identified in the A. flavus 3357 genome (http://www .aspergillusflavus.org/genomics/). We hypothesized that either these genes were truly absent in A. flavus 3357, or they were omitted due to gaps in coverage in the A. flavus 3357 genome sequence. We also observed that of the 318 genes not present in the A. flavus 3357 genome, 213 were expressed at significant levels in A. parasiticus SU-1, but the transcript was not detected in A. flavus 3357. Because these 213 genes and their transcript were absent in A. flavus 3357 but present in A. parasiticus SU-1, we propose that they represent genes that are uniquely expressed in A. parasiticus SU-1 (compared to A. flavus 3357). For the remaining 105 genes, transcripts were detected in A. flavus 3357, suggesting that they were not identified previously due to gaps in construction of the genome sequence.
Nucleotide sequence accession number. The A. parasiticus SU-1 scaffold sequences have been submitted to NCBI and assigned accession number JMUG00000000.

RESULTS
A. parasiticus SU-1 genome sequence and predicted chromosome structure. As part of the bioinformatics analysis, contigs from Illumina sequence analysis of the A. parasiticus SU-1 genome were assembled into 375 scaffolds greater than 1 kb in size (see Table S1 in the supplemental material). The resulting scaffolds were then listed by decreasing length. The top 31 A. parasiticus SU-1 scaffolds by length represented approximately 95% of the predicted SU-1 genome ( Table 1). The locations of all 375 A. parasiticus SU-1 scaffolds were assigned based on direct sequence comparison with 16 A. flavus 3357 scaffolds presented in the A. flavus 3357 genome browser (http://www.aspergillusflavus.org /genomics/) (see Table S2 in the supplemental material). The locations of the top 31 A. parasiticus SU-1 scaffolds by length were assigned to A. flavus 3357 chromosomes by a similar analysis, and these data were used to construct a schematic of the A. parasiticus SU-1 genome (Fig. 1A). We observed that the predicted chromosome structures of the A. parasiticus SU-1 and A. flavus 3357 genomes are quite similar ( Fig. 1A and B). The A. parasiticus SU-1 genome is predicted to consist of 8 chromosomes that have sizes similar to those of the 8 homologous A. flavus 3357 chromosomes. A total of 13,290 A. parasiticus SU-1 genes were identified based on the A. flavus 3357 genome annotation (carries 13,487 genes), and the relative locations of these loci were assigned to A. parasiticus SU-1 scaffolds and chromosomes by direct sequence comparison to the homolog in A. flavus 3357 scaffolds and chromosomes (see Table S3 in the supplemental material).
A. parasiticus SU-1 genes were provided A. flavus 3357 gene designations based on identity, and these genes were located on the A. parasiticus SU-1 scaffolds and chromosomes using the A. flavus 3357 genome as a template. We observed that 96% of A. flavus 3357 open reading frames (ORFs) have at least one homolog in the A. flavus 3357 genome, and 4% of the A. flavus 3357 ORFs have no detectable homolog in the A. parasiticus SU-1 genome ( Fig. 2A). A. parasiticus SU-1 and A. flavus 3357 share greater than 90% sequence identity over greater than 90% of genome sequence (Fig. 2B).
Fifty-five secondary metabolite gene clusters were identified and located on A. flavus 3357 chromosomes as part of the original   1 and Table 2). The exception is cluster 34, which appears to be located on the left arm of chromosome 2 (C2L) rather than the left arm of chromosome 6 (C6L). On the basis of this labor-intensive manual analysis, we noted that A. parasiticus SU-1 lacks 11 of the secondary metabolite clusters observed in A. flavus 3357, including clusters 4, 22, 28, 33, 36, 38, 40, 45, 48, 49, and cluster 55, which mediates synthesis of cyclopiazonic acid (CPA) ( Table 2). In contrast to manual analysis, antiSMASH software (default parameters) detected 93 secondary metabolite clusters in A. parasiticus SU-1, including 30 polyketide synthetase (PKS) clusters, 22 nonribosomal peptide (NRPS) clusters, 14 terpene clusters, 2 siderophore clusters, and 25 clusters designated "putative" or "other." antiSMASH detected 87 secondary metabolite clusters in A. flavus 3357, including 28 PKS clusters, 22 NRPS clusters, 22 terpene clusters, 1 siderophore cluster, and 14 clusters designated "other." We recently expanded this analysis to perform direct sequence comparison at the individual cluster level (J. Wee and J. E. Linz, unpublished data). Bioinformatics analysis support our initial observation that the proportions of the genome devoted to secondary metabolism are highly similar in A. parasiticus SU-1 and A. flavus 3357.
Comparison of A. parasiticus SU-1 and A. flavus 3357 gene expression under aflatoxin-inducing conditions. The level of differential expression during growth on YES medium (40 h) was calculated for all genes in the A. parasiticus SU-1 genome (see Table S4 in the supplemental material). Of 13,290 total transcripts expressed by A. parasiticus SU-1, 1,408 (11%) exhibited significantly different levels of gene expression between A. parasiticus SU-1 and A. flavus 3357 during growth in YES medium at 40 h. This medium induces aflatoxin gene expression in both species, and aflatoxin gene expression increases at maximum rates between 30 and 40 h in YES medium (14). The number of genes exhibiting statistically significant differences in expression are much lower between two biological replicates of A. flavus 3357 (63 genes) or two biological replicates of A. parasiticus SU-1 (34 genes) analyzed in this study (see Table S5 in the supplemental material), suggesting that the observed differences in expression between A. parasiticus SU-1 and A. flavus 3357 are not due to variability in gene expression between biological replicates.
Six aflatoxin genes (cypA, aflT, nor-1, estA, hypC, and ordA) and one additional secondary metabolism gene (aroM) are expressed at significantly higher levels (YES medium, at 40 h) in A. parasiticus SU-1 than in A. flavus 3357 (Table 3). In support of this observation, 13 additional genes in the aflatoxin gene cluster (cluster 54) are expressed at a minimum of 2-fold-higher levels in A. parasiticus SU-1 than in A. flavus 3357 including fas-2, a conserved hypothetical protein, omtB, avfA, hypB, verB, avnA, hypE, ordB, omtA, vbs, cypX, and moxY (Table 3 and Fig. 3); these differences in expression were not statistically significant. In contrast, the gene encoding the hybrid polyketide synthase/nonribosomal peptide synthase in the CPA gene cluster 55 and three stress response genes (Hsp70, MnSOD, and xylulose reductase) are expressed at significantly higher levels in A. flavus 3357 than in A. parasiticus SU-1. These data may help explain differences in the quantity of aflatoxin synthesized by the two species under our laboratory conditions (Fig. 4) and differences in CPA produced by A. flavus 3357. In support of these data, we did not detect the CPA NRPS/PKS (gene designation AFLA_13949) in the A. parasiticus de novo genome sequence. We also did not detect transcripts for adhA and, aflS (in cluster 54) as well as major facilitator superfamily (MFS) multidrug transporter, flavin adenine dinucleotide (FAD)-dependent oxidoreductase, and dimethylallyl tryptophan synthase (in cluster 55, the CPA cluster) in A. parasiticus SU-1.
Of the 44 secondary metabolite gene clusters detected manually in A. parasiticus, 26 were not expressed in either A. flavus 3357 or A. parasiticus SU-1 under the aflatoxin-inducing or -noninduc- The first value represents the expression level of the first member of the comparison pair, and the second number represents the expression level of the second member of the comparison pair. If the difference in expression between the first and second member of the pair is statistically significant, "yes" follows the values (based on the q value). If the difference in expression is not statistically significant, "no" follows the values. In some cases, the differences in expression were at least 2-fold, but this was not statistically significant; in this case, "no*" follows the values. ϩ, up-regulation; Ϫ, down-regulation.
ing conditions utilized in this study ( Table 2), suggesting that these gene clusters are either not expressed at significant levels or are expressed under growth conditions that these fungi encounter in the wild. The NRPS in clusters 37 and 45 as well as Arp1 in cluster 10 were expressed at significantly higher levels after 40 h of growth of A. parasiticus SU-1 in YES medium, while the PKS in cluster 40 was expressed at significantly higher levels in A. flavus 3357 under these same conditions ( Fig. 3 and Table 3). We also All values reported in these two columns are significant (q Ͻ 0.05) except for the values with asterisks. The values with asterisks display trends in differential expression, but they are not statistically significant (q Ͼ 0.05). Symbols; ϩ, upregulated in YES or YEP medium; Ϫ, downregulated in YES or YEP medium. b The level of expression of A. parasiticus SU-1 grown for 40 h in YES medium and the level of expression of A. flavus 3357 grown for 40 h in YES medium (SU-1/3357 YES) are shown. The first value represents the raw expression data for strain SU-1 grown in YES medium at 40 h, and the second value represents the raw expression data for strain 3357 grown in YES medium for 40 h. The ratio of these values determines the fold upregulation (ϩ) or downregulation (Ϫ), which is shown after the first two values. All values reported in this column are statistically significant.  (Table 3), providing evidence that these species differ greatly in the type and quantity of secondary metabolites synthesized under the specified growth conditions. Gene expression in A. parasiticus SU-1 and A. flavus 3357 was also compared during growth in YES medium (aflatoxin-inducing growth conditions) and YEP medium (aflatoxin-noninducing growth conditions) (see Table S6 in the supplemental material) with a focus on aflatoxin, secondary metabolism, and cellular defense/stress response genes (Table 3). A total of 1,284 A. parasiticus SU-1 genes and 1,802 A. flavus 3357 genes exhibited significant differences in gene expression when grown in YES medium or in YEP medium at 40 h (see Table S6 in the supplemental material). We analyzed only a single biological replicate of strain SU-1 on YEP medium for these comparisons, so these data are preliminary in nature. However, most genes within the aflatoxin gene cluster are upregulated in YES medium as expected for both A. flavus 3357 and A. parasiticus SU-1 ( Fig. 3 and Table 3), so it is likely that these trends in expression will be confirmed upon further analysis.
Identification of "uniquely expressed genes" in A. parasiticus. We observed that 860 A. parasiticus SU-1 loci have no homolog in the A. flavus 3357 genome sequence based on direct DNA sequence comparison between the genomes and RNA Seq analysis of gene expression on YES medium and/or YEP medium (see Tables S4 and S6 in the supplemental material). The absence of these loci in the A. flavus 3357 genome may imply that these loci are unique A. parasiticus SU-1 gene sequences. Alternatively, these loci may be located within gaps in the A. flavus 3357 genome sequence. Of these 860 loci, the difference in the level of expression between A. flavus 3357 and A. parasiticus SU-1 for 318 of these loci is statistically significant. For 213 of these loci, RNA transcripts (expression) could not be detected in A. flavus 3357, while transcripts were detected at significant levels in A. parasiticus SU-1 ( Fig. 5 and Table 4). We identified gene function for 68 of these genes by BLAST analysis of the translated protein products and then assigned predicted function by KOBAS GO analysis. Of particular importance to the current study, we identified 9 fungusspecific transcription factors, 4 enzymes involved in secondary metabolism, and 14 enzymes involved in stress response/cellular defense that are unique to A. parasiticus SU-1 and likely contribute to phenotypic and biochemical differences observed in the two species. Surprisingly, no known function was assigned to 145 of the 213 uniquely expressed genes (68%) based on BLAST analysis of the NCBI database (Fig. 5), suggesting that they may contribute unique functions in the biology of A. parasiticus SU-1.

DISCUSSION
The focus of this ongoing work was to identify the genetic basis for the observed species differences between A. parasiticus SU-1 and A. flavus 3357. Overall, the number of genes identified in our analysis of the A. parasiticus SU-1 genome (13,290) was similar to that identified in A. flavus 3357 (13,487). Even though their genome sequences and chromosome structure are similar, we observed statistically significant levels of differential expression between A. flavus 3357 and A. parasiticus SU-1 during growth in YES medium for 40 h and during growth on YES medium or YEP medium for 40 h. Furthermore, expression of 2,930 A. parasiticus SU-1 genes could not be detected in YES medium, and expression of 3,836 genes could not be detected on YEP growth medium. Since we analyzed only a single strain of each species, the results we present provide a starting point for a broader analysis of many different strains of each species to determine whether differences we observed in the current work extend beyond the two strains we analyzed.
One short-term goal of the current work was to explain specific species differences in mycotoxin biosynthesis, including the quantity and type of aflatoxin produced in culture as well as the synthesis of cyclopiazonic acid (CPA). We observed significantly higher levels of expression of 6 key genes (cypA, aflT, nor-1, estA, hypC, and ordA) in A. parasiticus SU-1 as well as a minimum 2-fold increase in expression of 13 other aflatoxin genes that could indicated by the distance above (upregulation) and below (downregulation) the zero difference horizontal line. The number of transcripts is shown on the x axis. (B) Differential expression of specific genes in secondary metabolism clusters 54 and 55. Genes that are upregulated (blue diamonds) and genes that are downregulated (red diamonds) are indicated. The common gene name(s) is presented next to the colored diamond and is based on annotations in the A. flavus genome browser. In specific cases, RNA Seq analysis presented aggregate expression data for closely spaced genes within cluster 54 (e.g., omtB, avfA, hypB, verB, avnA, and hypE [omtB/avfA/hypB/verB/avnA/hypE]). Genes included within these aggregates are closely spaced and encoded on the same strand of DNA. Presumably RNA SEQ was unable to resolve independent data for these groups of genes. (C) Bar graph presentation of differential expression data in panel B. in part explain species differences in aflatoxin biosynthesis, and future work can focus on expanding the focus of this analysis to additional strains of A. flavus and A. parasiticus. Upregulation of all 19 of these aflatoxin genes in A. parasiticus SU-1 could account for higher levels of aflatoxin synthesized by this species, while higher levels of ordA and cypA specifically (which mediate conversion of B aflatoxins to the corresponding G aflatoxins) could explain why most strains of A. parasiticus synthesize both B and G aflatoxins, while most strains of A. flavus synthesize predominantly B aflatoxins. We also detected significantly higher expression of the hybrid PKS/NRPS in A. flavus 3357 that catalyzes the initial steps in the synthesis of CPA in many A. flavus strains. The CPA PKS/NRPS was expressed at up to 100-fold-higher levels in A. flavus 3357 than in A. parasiticus SU-1. In support of this observation, this hybrid gene and three other genes in cluster 55 were not detected in the assembled A. parasiticus SU-1 genome.
We also compared the presence and location of 55 secondary metabolism clusters present in the A. flavus 3357 genome. While 42 of these clusters were detected in both organisms at the same approximate location, 1 cluster was present on a different chromosome, perhaps providing evidence for intragenomic recombination leading to translocation of this cluster. During analysis of the 318 loci that were present in the A. parasiticus SU-1 genome but not A. flavus 3357, we detected additional enzymes like PKS and NRPS, so it appears that both species are able to synthesize unique secondary metabolites. Further analysis of the A. flavus 3357 genome must be conducted to determine whether these genes were present but detected as part of the original sequence analysis or if they are indeed unique to A. parasiticus SU-1.
antiSMASH is a simple and rapid bioinformatics analysis for identification and functional annotation for secondary metabolite gene clusters in bacterial and fungal genome sequences (10). antiSMASH algorithms are based on specific hallmarks (functional motifs) and sequence identity in clusters of secondary metabolism Clusters of Orthologous Groups (smCOGs). This software program was reported to have 97% reliability when it was used to compare 484 cloned fungal and bacterial secondary metabolism gene clusters in the GenBank database. While we manually detected 44 clusters based on a single key enzyme present in a specific pathway, antiSMASH expanded upon this analysis and detected 93 clusters based on a collection of enzymes present in the same cluster. antiSMASH detected similar numbers of PKS and NRPS clusters in A. parasiticus SU-1 (30 PKS and 22 NRPS) and A. flavus 3357 (28 PKS and 22 NRPS), and approximately 60% of the secondary metabolite clusters are committed to synthesis of PKS and NRPS in both species. Of particular interest, there was a clear difference in the number of terpene clusters in A. parasiticus SU-1 (14 clusters) and A. flavus 3357 (22 clusters). Analyzing how the observed differences in the number and predicted function of terpene and siderophore clusters impact species phenotype will be one focus of future research.
Finally, we detected 213 genes that appear to be "uniquely expressed" in A. parasiticus SU-1. These genes and their transcripts were not identified in A. flavus 3357 as part of the current study or in the A. flavus 3357 genome available online. KOBAS GO Pathway analysis enabled us to place 27 of these 213 genes into functional groups, including secondary metabolism (4 genes), cellular defense/stress response (14 genes), and C6 (zinc binuclear cluster) transcription factors (9 genes) (15,16).
Members of the C6 transcription factor family are particularly interesting, because they are specific to fungi and are strongly associated with many secondary metabolism gene clusters in both Of 213 uniquely expressed A. parasiticus SU-1 genes, we were unable to assign predicted functions for 145 genes (68%) based on BLAST analysis of the NCBI genome database. The remaining 68 genes were assigned predicted functions and assigned to functional groups. Functional group annotations appear on the y axis as follows: Unassigned, transcripts with no known or assigned function; Other, transcripts with predicted functions other than SR, TF, and SM; SM, secondary metabolism; TF, transcription factors; SR, oxidative stress response. Each bar represents the number of genes in each functional group.