Differential regulation of mycelial growth and aflatoxin biosynthesis by Aspergillus flavus under different temperatures as revealed by strand-specific RNA-Seq.

Although several regulatory pathways have been reported for Aspergillus flavus, the regulation of aflatoxin production and mycelial growth under different temperatures remains unclear. In this study, A. flavus differentially expressed genes (DEGs) and regulatory pathways were analyzed under three temperatures, by strand-specific RNA-Seq. Results show that a total of 2,428 and 1,474 DEGs were identified in fungal mycelia cultured at 20°C and 37°C, respectively, as compared with the control (28°C). Approximately ~ 79% of DEGs in the 37°C samples were up-regulated genes, while ~ 63% of DEGs in the 20°C samples were down-regulated genes. Most of the DEG pathways enriched by lower temperatures differed from those enriched by higher temperatures, while only a small portion of the pathways were shared by A. flavus grown under different temperatures. Aflatoxin biosynthesis, Butanoate metabolism, oxidation-reduction process, and benzene-containing compound metabolic process were the shared down-regulated pathways, while steroid biosynthesis, oxidoreductase activity, cellular protein modification process, DNA binding, protein complex were the shared up-regulated pathways between lower and higher temperatures. The shared genes and pathways are the key regulatory candidates for aflatoxin biosynthesis with changes of temperature. In addition, the identification of both up-regulated and down-regulated genes provides a useful gene set for further investigation of the aflatoxin biosynthesis among Aspergillus.

INTRODUCTION

Aflatoxins are a highly toxic secondary metabolite produced by certain molds (e.g., Aspergillus flavus, A. parasiticus) (Chalivendra, DeRobertis, Chang, & Damann, 2017; Klich, 2007). Many important agricultural crops such as maize, rice, wheat, peanuts, and cotton, are prone to contamination by aflatoxins both pre‐ and postharvest (Bai et al., 2015; Yu et al., 2011). Humans and animals are commonly exposed to aflatoxins via the ingestion of contaminated food or products. Several studies have shown that the consumption of food contaminated with aflatoxins can be harmful to the liver, kidney, epididymis, testis, heart, brain, nervous system and can lead to immune suppression and carcinogenic effects, or even death (Kumar, Mahato, Kamle, Mohanta, & Kang, 2016; Nierman et al., 2015). Although a large number of techniques including physical, chemical, and biological methods have been developed to prevent or reduce the occurrence of aflatoxins in foods, these strategies are not always effective in eliminating grain contamination (Gressel & Polturak, 2018; Klich, 2007; Kumar et al., 2016). To date, the regulatory mechanisms for mycelial growth and aflatoxin biosynthesis under different environmental conditions remain unclear. Studies have also shown that the occurrence of aflatoxins in crops is influenced by various environmental factors, such as nitrogen levels, light, temperature, water, redox status, and pH level (Georgianna & Payne, 2009; Khlangwiset, Shephard, & Wu, 2011). Of these factors, high temperature and drought stress are commonly reported precursors to aflatoxin outbreaks in corn and other crops (Klich, 2007).

Previous investigations have shown that the genes for aflatoxin biosynthesis are clustered in the genome of Aspergillus (Sarma, Bhetaria, Devi, & Varma, 2017; Yu et al., 2004). Twenty‐five genes were identified from the 70 kb gene cluster in A. flavus and 82 kb gene cluster in A. parasiticus (Alkhayyat & Yu, 2014; Yu et al., 2004). Of these, the functions of 19 genes have been assigned, while the structures of at least 15 intermediate gene products have been defined (Alkhayyat & Yu, 2014; Yu et al., 2004). Enzymes convert the initial substrate acetate into the four major aflatoxins (B1, B2, G1, and G2) (Yu et al., 2004). It has been also confirmed that the expression of genes in the cluster is regulated by aflR and aflS. The product of aflR is a Gal 4‐type 47‐kDa polypeptide which binds to the palindromic sequence 5'‐TCGN5CGA‐3' in the cluster gene promoters, resulting in transcriptional activation of aflatoxin biosynthesis genes (Alkhayyat & Yu, 2014; Yu et al., 2004). AflS can bind AflR inhibitors and act as a transcriptional enhancer to optimize AflR activity (Alkhayyat & Yu, 2014).

The influence of temperature on aflatoxin biosynthesis by A. flavus has been investigated in many studies (Bai et al., 2015; Medina et al., 2017; Yu et al., 2011). The optimum temperature for aflatoxin formation is ~30°C, while the optimum temperature for mycelial growth is ~37°C (Yu et al., 2011). To understand the effects of varying temperatures on aflatoxin biosynthesis, several high‐throughput technologies, (e.g., transcriptomic and proteomic analyses) have been introduced. Yu et al., (2011) observed a large number of differentially expressed genes between 30°C and 37°C, in mycelia that were harvested 24 hr after inoculation. The average transcription level for the 30 aflatoxin biosynthesis genes increased by ~3,300‐fold at 30°C, as compared with 37°C (Yu et al., 2011). Bai et al., (2015) applied transcriptomic and proteomic analyses to identify changes in A. flavus at 37°C for 1.5 days and 28°C for 3 days, showing that post‐transcriptional processes play a critical role in regulating the protein level between the two temperatures. Lind, Smith, Saterlee, Calvo, and Rokas (2016) found that 11 temperature‐regulated gene clusters, associated with secondary metabolites, were required VeA at 37°C, and LaeA at both 30°C and 37°C. In addition, a large number of gene ontology and KEGG pathways have been identified in maize kernels colonized by A. flavus under different water activities (aw; 0.99 and 0.91) and temperatures (30°C, 37°C) after 10 days (Medina et al., 2017). With the help of high‐throughput technologies, the underlying mechanisms are increasingly being established, although little information is available on the changes in A. flavus mycelial growth and aflatoxin production under low temperature conditions and conflicting results exist on the feasibility of aflatoxin production by A. flavus at 37°C.

In this study, an RNA‐Seq approach was used to identify differentially expressed genes (DEGs) and regulatory pathways in A. flavus, under three temperature conditions (20°C, 28°C, and 37°C). Furthermore, the results of bioinformatic analysis were verified by experiments. The result of this study assists our understanding of the regulatory mechanisms of aflatoxin formation under different temperature conditions and can help develop strategies to control the production of aflatoxins in the food chain.

MATERIALS AND METHODS

Fungal strain and cultivation conditions

An aflatoxin‐producing strain, A. flavus NRRL 3357, was provided by Dr. Zhumei He (Sun Yat‐sen University, China). This strain can produce high amounts of aflatoxin B1 on YES agar (20 g/L yeast extract, 150 g/L sucrose, 15 g/L agar) under permissive conditions. In order to analyze the mycelial growth and toxin production abilities of A. flavus NRRL 3357 under different temperatures, agar plates were overlaid with sterile 8.5 cellophane sheets and single point inoculated (centrally) with 10 μl of spore suspension (106 spores in sterile water) and incubated at 20°C, 28°C, or 37°C. Three biological replicates were used for all subsequent analyses.

Fungal growth and aflatoxin analyses

Fungal mycelia were collected from the cellophane surface using a scraper, for weight and aflatoxin analyses. All fungal mycelia from each plate were transferred into a 50 ml tube containing 5 ml of methanol at room temperature, then incubated with continual agitation at 150 rpm, for 30 min. The supernatant was collected by centrifugation at 3,000 g and filtered through a syringe filter (RC 0.22 μm, Alltech). The presence of aflatoxin B1 was determined by HPLC with fluorescence detection, using a Waters 600 series HPLC equipped with a 600 pump, a 2,707 autosampler, and a 600 column thermostat set at 30°C. Detection was performed using a 2,475 Multi λ fluorescence detector set at 365 nm (λex) and 465 nm (λem), with a Waters Empower Windows xp operating system (Waters). The analytical column was a Luna 3u C18 (2) (150 × 4.6 mm, 3 μm) (Phenomenex) preceded by a SecurityGuard TM precolumn (C18, 4 × 3.0 mm, Phenomenex). The mobile phase consisted of methanol: water (55:45), eluted at a flow rate of 0.6 ml/min, with 20 μl of filtered extract injected into the HPLC per run. Aflatoxin B1 production was measured in μg/g of mycelia.

RNA isolation

Total RNA of A. flavus NRRL 3357 was isolated using a RNAiso plus (Takara) according to manufacturer's instructions. The quality and quantity of total RNAs were characterized using an Agilent 2100 and a Nano‐Drop 2000c instrument (Thermo Scientific).

RNA sequencing

Crude RNA was digested using 10 U DNase I (TaKaRa) at 37°C for 30 min. Ribosomal RNAs was removed using TruSeq Stranded mRNA Sample Preparation Kit (Illumina) according to the manufacturer's instructions. Strand‐specific RNA‐sequencing libraries were constructed using the TruSeq RNA Sample Prep Kit v2 (Illumina) according to the manufacturer's instructions. The fragments of an expected size were purified and amplified by PCR, with the purified PCR products sequenced using an Illumina Hiseq 4000 platform (BGI‐shenzhen).

Bioinformatics analyses

The quality of 150‐bp reads was assessed using FASTQC software (http://www.bioinformatics.babraham.ac.uk/projects/fastqc/). The paired‐end raw reads from RNA‐Seq were trimmed with low quality base‐calls and adaptor sequences using the pipeline Trimmomatic (v0.33) tool (Bolger, Lohse, & Usadel, 2014). Cleaned reads were mapped to the genome of A. flavus NRRL 3357 via HISAT2 (v2.1.0) (Kim et al., 2013). Uniquely mapped reads were used to quantify the raw counts using HTSeq (v0.9.1) (Anders, Pyl, & Huber, 2015). DEGs were calculated via DESeq2 using the parameters: p < 0.05 and a fold change >2) (Anders & Huber, 2010). Gene functions were annotated via the BLAST pipeline against references of the protein‐encoding sequence from the Nr of GenBank, Gene Ontology (GO) and the Kyoto Encyclopedia of Genes and Genomes (KEGG) (Ashburner et al., 2000; Kanehisa & Goto, 2000). Fisher's exact test was used to obtain enriched functional terms.

Real‐time quantitative PCR

Real‐time quantitative PCR was used to verify the gene expression level calculated from transcriptomic data. DEGs that may regulate fungal growth and aflatoxin production were verified and selected for further investigation. Crude RNA was used to synthesize cDNA using a transScript® first‐strand cDNA synthesis superMix kit (Transgen), where the 20 μL reaction system consisted of 10 μl SYBR® Fast qPCR Mix (2x), 0.5 μl of each primer (10 μmol/l) and 1 μl cDNA. The real‐time quantitative PCR program was set to the following sequence: 95°C for 30 s, followed by 40 cycles of 95°C for 5 s and finally, 60°C for 10 s. The β‐tubulin gene was used as an endogenous control, with three biological replicates assessed for each sample. Relative expression levels were calculated using the 2−ΔΔCT method.

RESULTS

Effects of temperature on mycelial growth and aflatoxin production

To assess the effects of temperature on fungal growth and aflatoxin production, three different culture temperatures were assessed. Figure 1a shows that optimal growth of fungal mycelia was observed at 37°C, while the slowest growth occurred at 20°C, as compared with the 28°C control. Cultivation at 28°C resulted in the highest level of aflatoxin B1 production on YES agar, while no aflatoxin B1 was detected by HPLC in all samples grown at 37°C, from 2 days to 7 days (Figure 1b). Cultivation at 20°C formed significantly less aflatoxin B1 than with that of the 28°C control (p < 0.05) (Figure 1b). Large differences were observed in both the fungal biomass and aflatoxin B1 content of A. flavus cultured at varying temperatures at 4 days and as changes in gene expression occurs in advance of changes in fungal growth and aflatoxin production, samples at 3 days were selected for RNA‐Seq.

Figure 1

Effects of temperature on A. flavus fungal growth (a) and toxin production (b)

Impact of different temperatures on gene expression

The quality and quantity of total RNAs were examined using Agilent 2100 and Nano‐Drop 2000c instruments, confirming that the isolated total RNAs were of a good enough quality for cDNA library construction (Table A1 in Appendix 1). Sequencing of all samples yielded a total of 155,177,520 raw paired‐end 150‐bp reads (Table A2 in Appendix 1). Assessment of the quality of raw reads by FASTQC showed that an overwhelming majority of reads had quality scores above Q30 (Appendix Figure A1), indicating the quality of the raw reads of all samples could be used for further analyses. 95.37% of reads (147,998,105 paired reads) remained as clean reads, following removal of adaptor, unknown, low quality and rRNA sequences (Table A2 in Appendix 1). The clean reads were used for mapping onto the genome of A. flavus NRRL 3357 for measurement of gene expression levels.

Table A1

The quality and quantity of isolated total RNAs

Sample	Concentration (ng/μL)	Volume (μL)	Total amount (μg)	OD260/280	OD260/230	RIN	28S/18S	Class
20T_1	882	40	35.28	1.99	2.34	7.1	2	A
20T_2	786	40	31.44	1.98	2.35	7.1	1.9	A
20T_3	1,416	40	56.64	2.03	2.39	7.1	2	A
28T_1	1,014	100	101.4	1.96	2.12	9.9	1.7	A
28T_2	750	20	15	2.01	1.95	9.2	2	A
28T_3	1,464	100	146.4	1.85	2.09	7.2	1.6	A
37T_1	1521	150	228.15	2.05	2.27	8.5	1.7	A
37T_2	1,370	100	137	1.95	2.14	8.9	1.7	A
37T_3	1834	100	183.4	1.81	2.02	8.3	1.7	A

Table A2

Statistical analysis of RNA‐Seq reads mapping results

Sample	Number of raw Seq	Number of clean Seq	Overall read mapping rate (%)	Uniquely mapped reads (%)
20T_1	17,197,487	16,355,506	76.22	68.14
20T_2	17,154,192	16,259,858	75.41	66.97
20T_3	17,383,718	16,616,781	77.96	69.77
28T_1	17,251,932	16,515,704	71.4	66.14
28T_2	17,221,767	16,446,529	68.74	62.68
28T_3	17,255,279	16,484,140	67.13	61.59
37T_1	17,221,828	16,449,309	67.45	62.45
37T_2	17,188,718	16,316,984	67.95	62.58
37T_3	17,302,599	16,553,294	68.17	63.16

The number of reads was expressed in pairs.

Figure A1

Quality of raw reads of two arbitrarily selected samples

Identification of DEGs

The genome sequence and related annotation files of A. flavus NRRL 3357 were obtained from the J. Craig Venter Institute (https://www.jcvi.org/). More than 61% of paired clean reads could be uniquely mapped to the A. flavus NRRL 3357 genome via pipeline HISAT2 (Table A2 in Appendix 1). 81.22% (10,953 out of 13,486) of putative protein‐coding genes could be detected throughout all samples and at the cut‐off count of ≧ 10 in at least one sample. Further analysis between samples showed a high correlation of the gene expression levels in the replicates of each treatment (Figure 2).

Figure 2

Correlation of gene expression levels between samples cultured under different temperatures

Due to A. flavus NRRL 3357 producing more aflatoxin at 28°C than at 20°C or 37°C, gene expression following mycelial growth at 28°C was selected as the control. Figure 3 and Table A3 in Appendix 1 show that a total of 2,428 and 1,474 more DEGs were identified in fungal mycelia grown at 20°C and 37°C, respectively, as compared with that of the 28°C control. It is of note that ~79% of the DEGs in 37°C samples belonged to up‐regulated genes, while only ~21% of the DEGs were down‐regulated. Conversely, ~ 63% of the DEGs in 20°C samples belonged to down‐regulated genes (Table A3 in Appendix 1). The Venn diagram shows that 7.1% of DEGs (137 genes) were shared between up‐regulated genes, and 4.4% of DEGs (77 genes) were shared between down‐regulated genes (Figure 4).

Figure 3

Volcano plots displaying DEGs for two different samples. The y‐axis shows the mean expression value of log10 (adjusted p value) and the x‐axis displays the log2 fold change value. The significantly up‐ and down‐regulated genes are shown as red dots and blue dots, respectively (p < 0.05, fold change >2). 20‐vs.‐28, mycelia at 20°C compared with control; 37‐vs.‐28, mycelia at 37°C compared with control

Table A3

Number of up‐regulated and down‐regulated genes

Sample comparison	Up	Down	Total
20°C vs. 28°C	889	1,539	2,428
37°C vs. 28°C	1,171	303	1,474

Figure 4

Venn diagram depicting the overlap in DEGs shared between up‐regulated genes or between down‐regulated genes

GO term and KEGG pathway enrichment analyses of the down‐regulated DEGs

The enriched GO terms of the down‐regulated DEGs impacted by high temperatures were apparently less than that of lower temperatures. The GO terms of the down‐regulated DEGs impacted by lower temperatures were enriched in serine‐type carboxypeptidase activity, oxidoreductase activity, nitrate metabolic process, fatty acid catabolic process, catalytic complex, cytoskeletal part, benzene‐containing compound metabolic process, etc, while oxidoreductase activity, aflatoxin biosynthetic process, austinol biosynthetic process, benzene‐containing compound metabolic process, protein complex, among others, were the enriched terms under higher temperature conditions (Figure 5).

Figure 5

GO term enrichments of down‐regulated DEGs. 20‐down, mycelia at 20°C compared with control; 37‐down, mycelia at 37°C compared with control. BP, Biological process; CC, Cellular component; MF, Molecular function

KEGG analysis showed that the down‐regulated DEGs enriched by lower temperatures were involved in carbon metabolism, nitrogen metabolism, amino acid metabolism, fatty acid degradation, peroxisome, among others (Figure 6). The enriched pathways of down‐regulated the DEGs by higher temperatures were involved in aflatoxin biosynthesis, butanoate metabolism, C5‐branched dibasic acid metabolism, biosynthesis of amino acids, etc (Figure 6). Aflatoxin biosynthesis, butanoate metabolism, oxidation‐reduction process, and benzene‐containing compound metabolic process were the shared down‐regulated pathways between lower and higher temperatures.

Figure 6

KEGG pathway enrichment of down‐regulated DEGs. 20‐down, mycelia at 20°C compared with control; 37‐down, mycelia at 37°C compared with control

GO term and KEGG pathway enrichment analyses of up‐regulated DEGs

The enriched GO terms of the up‐regulated DEGs impacted by high temperatures were apparently more than that at lower temperatures. The GO terms of the up‐regulated DEGs impacted by lower temperatures were enriched in iron ion homeostasis, amino acid transport, oxidoreductase activity, monooxygenase activity, carboxylic acid transmembrane transporter activity, steroid metabolic process, cellular protein modification process, oxidation−reduction process, cellular protein modification process, protein complex, etc, while nitrate metabolic process, melanin biosynthetic process, glycosaminoglycan catabolic process, asexual spore wall assembly, O‐methyltransferase activity, galactosidase activity, among others were the enriched terms under higher temperature (Figure 7).

Figure 7

GO term enrichments of up‐regulated DEGs. 20‐up, mycelia at 20°C compared with control; 37‐up, mycelia at 37°C compared with control. BP, Biological process; CC, Cellular component; MF, Molecular function

KEGG analysis showed that the up‐regulated DEGs enriched by lower temperature were involved only in steroid biosynthesis (Figure 7). The enriched pathways of up‐regulated the DEGs by higher temperatures were involved in steroid biosynthesis, glycosphingolipid biosynthesis, Nitrogen metabolism, Amino sugar and nucleotide sugar metabolism, etc (Figure 8). Steroid biosynthesis, oxidoreductase activity, cellular protein modification process, DNA binding, protein complex were the shared up‐regulated pathways between lower and higher temperatures.

Figure 8

KEGG pathway enrichment of up‐regulated DEGs. 20‐up, mycelia at 20°C compared with control; 37‐up, mycelia at 37°C compared with control

Aflatoxin biosynthesis processes

A. flavus aflatoxin biosynthesis genes were first analyzed using the SMURF informatics tool (Khaldi et al., 2010). Thirty genes were annotated in the aflatoxin biosynthetic cluster, with twenty‐two aflatoxin biosynthetic genes down‐regulated by lower temperatures, while all 30 genes were down‐regulated by higher temperatures, as compared to the control conditions (Table 1).

Table 1

The expression levels of genes in the Aflatoxin biosynthetic cluster

CDS_ID	Protein_id	Annotation	20°C (FPKM)	28°C (FPKM)	37°C (FPKM)
AFLA_139200	EED51155	aflQ/ordA/ord‐1/oxidoreductase/cytochrome P450 monooxigenase	20.68	57.68	10.21
AFLA_139210	EED51156	aflP/omtA/omt‐1/O‐methyltransferase A	55.92	162.99	25.22
AFLA_139220	EED51157	aflO/omtB/dmtA/O‐methyltransferase B	336.20	592.08	85.29
AFLA_139230	EED51158	aflI/avfA/cytochrome P450 monooxygenase	15.42	21.13	4.11
AFLA_139240	EED51159	aflLa/hypB/hypothetical protein	41.56	102.93	21.14
AFLA_139250	EED51160	aflL/verB/desaturase/P450 monooxygenase	81.25	192.50	36.76
AFLA_139260	EED51161	aflG/avnA/ord‐1/cytochrome P450 monooxygenase	19.68	46.12	10.72
AFLA_139270	EED51162	aflNa/hypD/hypothetical protein	1517.30	1,451.73	303.44
AFLA_139280	EED51163	aflN/verA/monooxygenase	57.34	77.57	19.94
AFLA_139290	EED51164	aflMa/hypE/hypothetical protein	121.35	320.65	59.98
AFLA_139300	EED51165	aflM/ver‐1/dehydrogenase/ketoreductase	551.98	1,105.51	163.80
AFLA_139310	EED51166	aflE/norA/aad/adh‐2/NOR reductase/dehydrogenase	320.36	450.69	64.53
AFLA_139320	EED51167	aflJ/estA/esterase	368.31	631.40	94.41
AFLA_139330	EED51168	aflH/adhA/short‐chain alcohol dehydrogenase	163.41	274.31	38.44
AFLA_139340	EED51169	aflS/pathway regulator	215.12	275.53	207.27
AFLA_139360	EED51170	aflR/apa‐2/afl‐2/transcription activator	24.44	22.74	20.89
AFLA_139370	EED51171	aflB/fas‐1/fatty acid synthase beta subunit	32.57	57.55	11.88
AFLA_139380	EED51172	aflA/fas‐2/hexA/fatty acid synthase alpha subunit	11.35	40.30	10.40
AFLA_139390	EED51173	aflD/nor‐1/reductase	76.12	348.68	47.04
AFLA_139400	EED51174	aflCa/hypC/hypothetical protein	29.85	195.80	35.97
AFLA_139410	EED51175	aflC/pksA/pksL1/polyketide synthase	30.04	278.49	61.46
AFLA_139420	EED51176	aflT/aflT/transmembrane protein	54.71	78.62	48.04
AFLA_139430	EED51177	aflU/cypA/P450 monooxygenase	0.03	0.06	0.04
AFLA_139440	EED51178	aflF/norB/dehydrogenase	8.88	7.63	6.94
AFLA_139450	EED51179	Conserved hypothetical protein	0.00	0.15	0.10
AFLA_139460	EED51180	MFS multidrug transporter putative	495.81	344.31	289.30
AFLA_139470	EED51181	FAD‐dependent oxidoreductase putative	1108.93	570.24	233.71
AFLA_139480	EED51182	Dimethylallyl tryptophan synthase putative	1118.70	682.57	248.59
AFLA_139490	EED51183	Hybrid PKS/NRPS enzyme putative	168.29	103.42	45.36
AFLA_139500	EED51184	Conserved hypothetical protein	0.50	0.29	0.10

Oxidoreductase activity

Oxidoreductase activity is driven by laccase or multicopper oxidase. Although the expression levels of the two genes encoding oxidoreductase activity were significantly higher at 37°C than at 28°C, the total expression levels were ~9% lower at 37°C than at 28°C (Table 2). It is of note, that the total expression level was highest at 20°C, suggesting that oxidation–reduction (redox) reactions are influenced by temperature.

Table 2

The expression levels of members of oxidoreductase activity

CDS_ID	Protein_id	Seq. Description	20°C (FPKM)	28°C (FPKM)	37°C (FPKM)
AFLA_084170	EED57720	Laccase	2.78	2.34	4.34
AFLA_089660	EED48886	Iron transport multicopper oxidase fet3	0.00	0.09	0.12
AFLA_000890	EED47448	Laccase	0.00	0.00	0.02
AFLA_006620	EED48018	Iron transport multicopper oxidase fet3	241.78	96.82	80.56
AFLA_120890	EED45860	Extracellular dihydrogeodin oxidase laccase	3.77	4.06	8.24

The DEGs shared by lower and higher temperature

It can be seen that the 77 down‐regulated genes shared by lower and higher temperatures were involved in aflatoxin biosynthetic process, and sterigmatocystin biosynthetic process (Table A4 in Appendix 1, Figure 9). The 137 up‐regulated genes shared by lower and higher temperatures were involved in cellular macromolecule biosynthetic process, gene expression, nucleic acid metabolic process, and gliotoxin biosynthetic process (Table A5 in Appendix 1, Figure 9). The functions of many DEGs are still unclear.

Table A4

Annotation of down‐regulated genes shared by lower and higher temperature

CDS Id	Protein Id	Seq. description
AFLA_072350	EED56541	Conserved hypothetical protein
AFLA_073450	EED56651	Amine oxidase
AFLA_074340	EED56740	Uncharacterized oxidoreductase Ygbj
AFLA_074510	EED56757	Short‐chain dehydrogenase family
AFLA_082440	EED57547	Methyltransferase domain‐containing protein
AFLA_084250	EED57728	Ketose‐bisphosphate aldolase class‐ii family protein
AFLA_087620	EED58063	nmra family transcriptional regulator
AFLA_024510	EED55180	MFS quinate transporter
AFLA_033970	EED56125	Conserved hypothetical protein
AFLA_034390	EED56167	Conserved hypothetical protein
AFLA_013220	EED54070	Beta‐lactamase family protein
AFLA_016350	EED54383	Nadph‐dependent fmn reductase
AFLA_016620	EED54410	Tryptophanyl‐trna synthetase
AFLA_107880	EED53413	Facilitator of iron transport 3
AFLA_109300	EED53553	fad synthetase
AFLA_110220	EED53645	mgs207 protein
AFLA_038420	EED52141	Chitosanase
AFLA_039270	EED52225	Carboxylesterase family protein
AFLA_039280	EED52226	Carboxylesterase family protein
AFLA_039650	EED52263	Nadh‐cytochrome b5 reductase
AFLA_041280	EED52426	Dimethylaniline monooxygenase
AFLA_061120	EED51849	Polyamine transporter
AFLA_061760	EED51913	Conserved hypothetical protein
AFLA_062700	EED52007	Mitochondrial carrier
AFLA_062890	EED52026	Hypothetical protein AFLA_062890
AFLA_136790	EED50916	Acetolactate synthase
AFLA_138150	EED51051	Hypothetical protein AFLA_138150
AFLA_138870	EED51122	Cyanovirin‐n family protein
AFLA_138920	EED51127	Dimethylaniline monooxygenase
AFLA_139170	EED51152	Sterigmatocystin biosynthesis monooxygenase stcw
AFLA_139200	EED51155	Cytochrome p450 monooxygenase
AFLA_139210	EED51156	o‐ partial
AFLA_139240	EED51159	Partial
AFLA_139250	EED51160	afll verb desaturase p450 monooxygenase
AFLA_139260	EED51161	aflg avna ord‐1 cytochrome p450 monooxygenase
AFLA_139290	EED51164	Hypothetical e partial
AFLA_139380	EED51172	Fatty acid synthase alpha subunit
AFLA_139390	EED51173	Norsolorinic acid partial
AFLA_139400	EED51174	duf1772‐domain‐containing protein
AFLA_139410	EED51175	Polyketide synthase
AFLA_064900	EED49669	Nadh‐ubiquinone oxidoreductase kda mitochondrial
AFLA_066480	EED49826	Cyclase
AFLA_066970	EED49875	Endonuclease exonuclease phosphatase family protein
AFLA_066980	EED49876	Polyketide synthase
AFLA_067740	EED49952	Conserved hypothetical protein
AFLA_067880	EED49966	MFS transporter
AFLA_069030	EED50081	Conserved hypothetical protein
AFLA_093890	EED49308	Conserved hypothetical protein
AFLA_096130	EED49531	clavata3 esr‐like protein
AFLA_124290	EED48206	c6 transcription
AFLA_124300	EED48207	MFS general substrate transporter
AFLA_124970	EED48274	Ankyrin domain protein
AFLA_125000	EED48277	MFS multidrug transporter
AFLA_125330	EED48310	Ankyrin repeat‐containing protein
AFLA_128540	EED48631	Proline oxidase
AFLA_005040	EED47863	Phenazine biosynthesis‐like
AFLA_005570	EED47915	Short‐chain dehydrogenase reductase family
AFLA_052610	EED47207	Succinyl‐ :3‐ketoacid‐coenzyme a transferase
AFLA_053540	EED47299	fad‐dependent oxidoreductase
AFLA_054060	EED47351	atp gtp‐binding protein
AFLA_054520	EED47397	1‐aminocyclopropane‐1‐carboxylate synthase
AFLA_054530	EED47398	Synaptic vesicle transporter svop
AFLA_054550	EED47400	myo‐inositol 2‐dehydrogenase
AFLA_097300	EED46067	Metal‐activated pyridoxal enzyme
AFLA_097340	EED46071	Transmembrane amino acid transporter protein
AFLA_097530	EED46090	duf1857 domain‐containing protein
AFLA_098140	EED46151	Monocarboxylate
AFLA_101930	EED46529	Succinate‐semialdehyde dehydrogenase
AFLA_101990	EED46535	Hexose carrier protein
AFLA_116550	EED45426	Glycoside hydrolase family 24 protein
AFLA_117000	EED45471	RNA
AFLA_118420	EED45613	Hypothetical protein AFLA_118420
AFLA_118450	EED45616	Six‐hairpin glycosidase
AFLA_118740	EED45645	Xylose isomerase tim barrel
AFLA_120990	EED45870	o‐methyltransferase
AFLA_122140	EED45985	Acetyltransferase
AFLA_009040	EED45020	3‐Isopropylmalate dehydrogenase

Figure 9

GO term enrichments of shared up‐regulated DEGs and shared down‐regulated DEGs. Same_down, The enriched GO terms by same down‐regulated genes; Same_up, The enriched GO terms by same up‐regulated genes

Table A5

Annotation of up‐regulated genes shared by lower and higher temperature

CDS Id	Protein Id	Seq. description
AFLA_078540	EED57159	Immediate early response protein ier
AFLA_082090	EED57512	Fungal‐specific transcription factor domain‐containing protein
AFLA_082720	EED57575	Fatty acid elongase
AFLA_083800	EED57683	Endonuclease exonuclease phosphatase family protein
AFLA_087280	EED58029	Alpha‐ketoglutarate‐dependent taurine dioxygenase
AFLA_022820	EED55021	Conserved hypothetical protein
AFLA_022830	EED55022	Conserved hypothetical protein
AFLA_023350	EED55064	Beta‐glucosidase m
AFLA_023870	EED55116	Transmembrane domain of the epidermal growth factor receptor family of protein tyrosine kinase
AFLA_023880	EED55117	Hsp70 family chaperone
AFLA_023890	EED55118	Conserved hypothetical protein
AFLA_025190	EED55248	Lipase
AFLA_027990	EED55527	Conserved hypothetical protein
AFLA_028640	EED55592	Cytochrome p450 61
AFLA_028890	EED55617	Tartrate dehydrogenase
AFLA_030430	EED55771	Fatty acid oxygenase
AFLA_032440	EED55972	Conserved hypothetical protein
AFLA_034810	EED56209	Hypothetical protein AFLA_034810
AFLA_035890	EED56317	Acyl‐ n‐acyltransferase
AFLA_036370	EED56365	Phosphoenolpyruvate carboxykinase
AFLA_037250	EED56451	Cyanide hydratase
AFLA_014030	EED54151	Hypothetical protein AFLA_014030
AFLA_014040	EED54152	Hypothetical protein AFLA_014040
AFLA_014300	EED54178	Sodium bile acid symporter family protein
AFLA_016530	EED54401	Beta‐galactosidase
AFLA_017480	EED54496	Sun domain protein
AFLA_017760	EED54524	Phenol 2‐
AFLA_018350	EED54583	4‐coumarate‐‐ ligase‐like 7
AFLA_018740	EED54622	Cora family metal ion transporter
AFLA_020960	EED54844	Copper resistance‐associated p‐type atpase
AFLA_104700	EED53096	Monooxygenase
AFLA_105270	EED53153	Conserved hypothetical protein
AFLA_106900	EED53315	Major facilitator superfamily general substrate transporter
AFLA_109230	EED53546	2‐hydroxyacid dehydrogenase
AFLA_109460	EED53569	Family taurine
AFLA_113430	EED53966	Transcription factor subunit 5
AFLA_038530	EED52152	Elastinolytic metalloproteinase mep
AFLA_039540	EED52252	Conserved hypothetical protein
AFLA_040140	EED52312	Aquaporin
AFLA_040580	EED52356	Serine threonine protein kinase
AFLA_041410	EED52439	Aldehyde dehydrogenase family protein
AFLA_046740	EED52972	nadp‐dependent malic enzyme
AFLA_057600	EED51497	Heat shock protein
AFLA_057680	EED51505	Beta‐n‐
AFLA_057710	EED51508	Oxaloacetate acetylhydrolase
AFLA_057810	EED51518	Glutamine‐serine‐proline rich
AFLA_059990	EED51736	flavin‐dependent halogenase o‐methyltransferase bifunctional protein
AFLA_060020	EED51739	Polyketide synthase
AFLA_060770	EED51814	Protein alcs
AFLA_062630	EED52000	hyp effector
AFLA_062820	EED52019	Polyketide synthase
AFLA_063040	EED52041	Glycosyl hydrolase family 3 n terminal domain‐containing protein
AFLA_063240	EED52061	Hypothetical protein AFLA_063240
AFLA_063250	EED52062	Glutaminyl cyclase
AFLA_063260	EED52063	Lwamide neuropeptide partial
AFLA_063270	EED52064	Hypothetical protein AFLA_063270
AFLA_063280	EED52065	Conserved hypothetical protein
AFLA_063290	EED52066	Conserved hypothetical protein
AFLA_063300	EED52067	Conserved hypothetical protein
AFLA_063310	EED52068	tat pathway signal sequence protein
AFLA_132470	EED50485	Toxin biosynthesis protein
AFLA_133640	EED50602	Cell wall cysteine‐rich protein
AFLA_133810	EED50619	Conserved hypothetical protein
AFLA_137770	EED51013	umta methyltransferase family protein
AFLA_138060	EED51042	c−24 sterol reductase
AFLA_138400	EED51076	nad‐dependent epimerase dehydratase
AFLA_138760	EED51111	trx2p
AFLA_064380	EED49617	radh flavin‐dependent halogenase
AFLA_064390	EED49618	Cytochrome p450
AFLA_064400	EED49619	Cytochrome p450
AFLA_064440	EED49623	Heavy metal tolerance protein
AFLA_064450	EED49624	1‐aminocyclopropane‐1‐carboxylate synthase
AFLA_064460	EED49625	Toxin biosynthesis protein
AFLA_064470	EED49626	Cytochrome p450
AFLA_064480	EED49627	Thioredoxin reductase
AFLA_064490	EED49628	Methyltransferase
AFLA_064510	EED49630	Thioredoxin reductase
AFLA_064530	EED49632	Glutathione s‐transferase
AFLA_064540	EED49633	Cytochrome p450 monooxygenase
AFLA_064550	EED49634	Membrane dipeptidase
AFLA_064560	EED49635	Nonribosomal peptide synthase ‐like
AFLA_064570	EED49636	ncs1 nucleoside transporter
AFLA_064580	EED49637	Dimeric dihydrodiol
AFLA_064590	EED49638	o‐methyltransferase
AFLA_064600	EED49639	Major facilitator superfamily domain
AFLA_066050	EED49783	DNA repair family protein
AFLA_066710	EED49849	Oxoglutarate iron‐dependent dioxygenase
AFLA_070280	EED50206	Siderophore esterase ‐like protein
AFLA_070400	EED50218	aaa family atpase
AFLA_070420	EED50220	Siderochrome‐iron transporter
AFLA_090590	EED48978	Alpha‐ ‐ subfamily
AFLA_093580	EED49277	Integral membrane protein
AFLA_095310	EED49450	Conserved hypothetical protein
AFLA_096180	EED49536	duf636 domain protein
AFLA_096650	EED49583	Conserved hypothetical protein
AFLA_096660	EED49584	Conserved hypothetical protein
AFLA_122840	EED48082	Conserved hypothetical protein
AFLA_123700	EED48147	Extracellular proline‐rich protein
AFLA_125620	EED48339	dj‐1 ‐type
AFLA_125760	EED48353	Squalene cyclase
AFLA_126510	EED48428	Copper‐transporting atpase
AFLA_127490	EED48526	Hypothetical protein AOR_1_770164
AFLA_128040	EED48581	Major facilitator superfamily transporter
AFLA_128050	EED48582	Serine hydrolase fsh
AFLA_128110	EED48588	Aquaglyceroporin
AFLA_129750	EED48752	mt‐A70 family
AFLA_000850	EED47444	Isoamyl alcohol
AFLA_003960	EED47755	Hypothetical protein AFLA_003960
AFLA_004270	EED47786	Protein kinase‐like domain
AFLA_004870	EED47846	Cytochrome p450
AFLA_005760	EED47934	Conserved hypothetical protein
AFLA_007170	EED48073	Pumilio‐family rna binding repeat protein
AFLA_048390	EED46786	S‐Adenosyl‐l‐methionine‐dependent methyltransferase
AFLA_049160	EED46863	Cyclopentanone ‐monooxygenase
AFLA_049210	EED46868	Integral membrane protein
AFLA_049520	EED46899	Integral membrane protein pth11
AFLA_052520	EED47198	Hypothetical protein AFLA_052520
AFLA_054360	EED47381	Methyltransferase type 11
AFLA_099110	EED46248	Fibronectin type iii domain‐containing protein
AFLA_099750	EED46311	Epoxide hydrolase
AFLA_100260	EED46362	t5orf172 domain protein
AFLA_101540	EED46490	Protein
AFLA_116390	EED45410	Amino acid transporter
AFLA_120630	EED45834	Formate dehydrogenase
AFLA_121190	EED45890	Zinc‐binding oxidoreductase
AFLA_121730	EED45944	Alpha‐galactosidase c
AFLA_121740	EED45945	Hypothetical protein AFLA_121740
AFLA_122040	EED45975	Oleate delta‐12 desaturase
AFLA_007600	EED44877	Oligopeptide transporter opt superfamily
AFLA_009910	EED45107	Membrane fusion mating protein Figure 1
AFLA_010590	EED45175	Siderophore biosynthesis lipase
AFLA_010610	EED45177	enoyl‐ hydratase isomerase family protein
AFLA_010620	EED45178	amp‐dependent synthetase ligase
AFLA_010630	EED45179	abc multidrug transporter
AFLA_010640	EED45180	Siderophore iron transporter
AFLA_010740	EED45190	Carboxypeptidase s1
AFLA_011540	EED45270	Multiple drug resistance protein

Real‐time PCR verification and analysis of several DEGs

To verify the reliability of DEGs identified by RNA‐Seq, the relative expression levels of several DEGs were further investigated via real‐time PCR. Results show that a similar expression pattern was observed between RNA‐Seq and real‐time PCR data (Appendix Figure A2), indicating that the relative expression level identified via RNA‐Seq was reliable.

Figure A2

Relative expression of several DEGs via Real‐time PCR

DISCUSSION

The results of the present study demonstrate the complexity of aflatoxin biosynthesis regulation and fungal mycelial growth, under different temperatures. As compared to the control, 1,539 genes were significantly down‐regulated by the reduced temperature of 20°C, while only 303 genes were significantly down‐regulated by the higher temperature of 37°C, in mycelia at 3 days (Table A3 in Appendix 1). It is also very interesting that a majority of the down‐regulated genes at higher temperatures related to secondary metabolic processes (Figures 5 and 6), while a majority of genes up‐regulated by higher temperatures were related to primary metabolic processes (Figures 7 and 8) related to fungal growth (Wisecaver, Slot, & Rokas, 2014). It also has been established that sterigmatocystin compounds are the precursor substances for aflatoxin synthesis (Georgianna & Payne, 2009). Majority of the genes related to sterigmatocystin biosynthesis and aflatoxin synthesis were down‐regulated by the higher temperature of 37°C, which was consistent with the observation that no aflatoxin B1 was detected by HPLC in all samples grown at 37°C.

Majority of the enriched pathways at lower temperatures and higher temperature were different. Only a small portion of the pathways were shared by mycelium grown under different temperatures. Aflatoxin biosynthesis, butanoate metabolism, oxidation–reduction process, and benzene‐containing compound metabolic process were the shared down‐regulated pathways, while steroid biosynthesis, oxidoreductase activity, cellular protein modification process, DNA binding, protein complex were the shared up‐regulated pathways between lower and higher temperatures. Secondary metabolic process (toxin metabolic process) and oxidoreductase activity were also shared by the absence of two key members (VeA and LaeA) of the Velvet protein complex (Lind et al., 2016). The results of the present study are in accordance with previous investigation, which implied that the enriched pathways in this study should be reliable and can be used for further investigation. The results of the present study also demonstrate that the total expression level of genes related to oxidoreductase activity was lowest at 37°C and highest at 20°C. These results indicate that fast growth mycelia may be less affected by oxidative stress at 37°C, while slow growth mycelia might be affected by oxidative stress at 20°C. Previous studies have shown that oxidative stress can cause changes in cytosolic and mitochondrial calcium concentrations in A. nidulans (Greene, Cao, Schanne, & Bartelt, 2002). Several investigations have also implied that antioxidants can significantly inhibit aflatoxin production, while oxidants enhanced aflatoxin production (Kim et al., 2008; Narasaiah, Sashidhar, & Subramanyam, 2006; Reverberi et al., 2005). According to our results and previously reported studies, oxidative stress resulting from reactive oxygen species might be involved in instigating aflatoxin biosynthesis. The formation of aflatoxin should be regulated via different pathways while being independent from fungal growth.

Previous studies have observed that spore and pigment production is accompanied by the formation of toxins (Chang, 2008; Georgianna & Payne, 2009). In the present study, genes for asexual spore wall assembly were significantly up‐regulated under higher temperatures (Figure 7), suggesting that spores are associated or involved in the process of syntheses. Some previous studies have not detected aflatoxins in growth medium at 37°C, while other studies have detected aflatoxins under the same temperatures. According to the enriched pathways, it may be speculated that higher temperatures may block the formation of aflatoxin by delaying the synthesis of its precursor substances during the fast growth stage of fungal mycelia. According to this theory, we cultivated the fungus again as described in the materials and methods section, incubating samples at 37°C for a longer period of time and at 10 days, aflatoxin concentrations reached as high as 10 μg/g. These results support the proposed theory and explain the contrasting results reported by previous studies.

Several genes have been reported to be involved in the production of aflatoxin. Members of the Velvet protein complex coding genes (VeA, LaeA) and homeobox gene (hbx1) were proved to be required for the formation of aflatoxin (Cary et al., 2019; Lind et al., 2016). All of them were not the down‐regulated genes shared by lower and higher temperatures (Table A4 in Appendix 1). Only the expression of VeA was down‐regulated by lower temperature, which suggested that the regulation of aflatoxin biosynthesis by different temperatures is very complex. It is worth noting that a few down‐regulated genes shared by lower and higher temperatures were related to toxin biosynthesis, several shared DEGs belong to transcription factor and methyltransferase. In addition, the functions of many DEGs are unclear and still need further investigation. Among the 77 down‐regulated genes shared by lower and higher temperature (Table A4 in Appendix 1), and the 137 up‐regulated genes shared by lower and higher temperature (Table A5 in Appendix 1), majority of the genes are likely involved in the direct aflatoxin biosynthesis or indirect regulation of aflatoxin biosynthesis. The shared DEGs list provides a useful gene set for further investigation of the aflatoxin biosynthesis among Aspergillus.

CONFLICT OF INTERESTS

The authors declare no conflict of interest.

AUTHOR CONTRIBUTIONS

Fang Tao conceived and designed experiments. Guomin Han carried out the bioinformatics studies, and Kai Zhao conducted experiments. Guomin Han and Fang Tao wrote the manuscript. Xiaodan Yan, Fangzhi Xiang, and Xuede Li helped draft the manuscript.

ETHICS STATEMENT

None required.