1 Transcriptome analysis of non-ochratoxigenic Aspergillus carbonarius strains and 2 interactions with other ochratoxigenic black aspergilli 3 4 Gemma Castellá1, M. Rosa Bragulat1, Riccardo Aiese Cigliano2, F. Javier Cabañes1* 5 1 6 Veterinary Mycology Group, Department of Animal Health and Anatomy, Universitat 7 Autònoma de Barcelona, Bellaterra, Catalonia, Spain. 2 8 Sequentia Biotech SL, Barcelona, Catalonia, Spain. 9 10 * Corresponding author: 11 F. Javier Cabañes 12 Grup de Micologia Veterinària 13 Departament de Sanitat i d'Anatomia Animals 14 Facultat de Veterinària, Universitat Autònoma de Barcelona 15 08193 Bellaterra (Barcelona), Spain. 16 Tel: +34 93 5811749 17 Fax: +34 93 5812006 18 e-mail: javier.cabanes@uab.es 1 19 Abstract 20 Aspergillus carbonarius consistently produces large amounts of ochratoxin A (OTA), a 21 mycotoxin with nephrotoxic effects on animals and humans. In the present study, we 22 analyzed the transcriptional changes associated to OTA production in three atypical 23 non-ochratoxigenic strains of A. carbonarius. In addition, in vitro interactions between 24 ochratoxigenic strains of A. carbonarius and A. niger and non-ochratoxigenic strains of 25 A. carbonarius and A. tubingensis were studied in order to evaluate their potential for 26 controlling OTA production. RNA-seq analysis revealed that there are 696 differentially 27 expressed genes identified in the three non-OTA producing strains, including 280 up- 28 regulated and 333 down-regulated genes. A functional and gene ontology enrichment 29 analysis revealed that the processes related to metabolic and oxidation processes, 30 associated with functions such oxidoreductase and hydrolase activity were down 31 regulated. All the genes related with OTA biosynthesis in A. carbonarius were the most 32 down-regulated genes in non-ochratoxigenic strains. We also showed that these strains 33 possess a deleterious mutation in the AcOTApks gene required for OTA biosynthesis. 34 Moreover, one of these strains gave the best control of OTA production resulting in an 35 OTA reduction of 98-100% in co-inoculation with an ochratoxigenic strain of A. niger 36 and an OTA reduction of 79-89% with an ochratoxigenic strain of A. carbonarius. 37 Results of this study provided novel insights into the knowledge of the OTA 38 biosynthetic pathway in these non-ochratoxigenic wild strains, and showed the 39 biocontrol potential of these strains. 40 41 Keywords: Aspergillus carbonarius, ochratoxin A, transcriptome, fungal interactions, 42 biological control agents. 43 2 44 1. Introduction 45 Ochratoxin A (OTA) is a nephrotoxic mycotoxin that has been related to human 46 and animal diseases. It is classified as a 2B carcinogen by the International Agency for 47 Research on Cancer (Pfohl-Leszkowicz and Manderville, 2012) and can be found in a 48 variety of common foods and beverages (Pit and Hocking, 2009). This mycotoxin is 49 produced by several species of Penicillium and Aspergillus, being Aspergillus 50 carbonarius the main OTA source on grapes and derived products. This species is also 51 considered a potential source of OTA in coffee and cocoa, and it has been also reported 52 from other foods such as figs, maize, paprika, and peanuts among others (Cabañes and 53 Bragulat, 2018). 54 The genes involved in the OTA biosynthetic pathway have not been completely 55 characterized, even if many studies have been carried out in several OTA-producing 56 fungi (Wang et al., 2016). Recently, a consensus OTA biosynthetic pathway has been 57 proposed in Aspergillus ochraceus (Wang et al., 2018). Key genes of the OTA pathway 58 have been identified by comparative genomic analysis with other OTA-producing fungi 59 as A. carbonarius, A. niger, A. steynii, A. westerdijkiae and P. nordicum. 60 In A. carbonarius a hypothetical OTA gene cluster (cluster 38) has been 61 described (Gerin et al., 2016). This cluster contained 21 genes and only some have been 62 correlated with OTA biosynthesis. Three genes have been directly related to OTA 63 biosynthesis: the nonribosomal peptide synthetase AcOTAnrps (Gallo et al., 2012), the 64 polyketide synthase AcOTApks (Gallo et al., 2014) and the halogenase AcOTAhal 65 (Ferrara et al., 2016). Two other genes located in the same genomic region have been 66 described: the cytochrome p450 monooxygenase AcOTAp450 and the transcription 67 factor AcOTAbZIP. It has been shown that under OTA induction conditions all these 68 genes were up-regulated and co-expressed (Gerin et al., 2016). 3 69 OTA production is a very consistent property of A. carbonarius. While A. 70 carbonarius consistently produces large amounts of ochratoxin A, the reported 71 percentages of OTA-producing strains in the A. niger aggregate is much lower (Cabañes 72 and Bragulat, 2018). In non-ochratoxigenic strains of A. niger, a 21-kb deletion in a 73 remnant of the polyketide sinthase gene of the putative OTA cluster has been identified 74 (Andersen et al., 2011). In the same way, the genome analysis of A. tubingensis 75 revealed that it does not contain an orthologue of this gene (Choque et al., 2018). In 76 previous studies, we resequenced the genome of three atypical and unique non-OTA 77 producing strains of A. carbonarius (Cabañes et al., 2015; Castellá et al., 2018). 78 Although no large deletions in functional genes related with OTA production were 79 found, some private missense variants of non-ochratoxigneic strains in AcOTApks gene 80 were detected (Castellá et al., 2018). The rest of the OTA cluster genes (AcOTAnrps, 81 AcOTAp450, AcOTAhal, and AcOTAbZIP) showed no privative mutations in the non- 82 OTA producing strains. 83 On the other hand, one of the more promising strategies for control of aflatoxins 84 in crops involves the use of atoxic strains of A. flavus as biological control agents 85 (BCAs) (Ehrlich, 2014). These carcinogenic mycotoxins contaminate worldwide food 86 and feed goods (e.g. maize, cotton seed, and peanuts), which pose a serious health risk 87 to humans and domestic animals, and are the cause of crop losses. It is not difficult to 88 isolate this kind of atoxigenic strains in this species because a high percentage of 89 nonaflatoxigenic strains may be present in some natural environments (Alaniz, 2018; 90 Ehrlich, 2008). On the contrary, OTA production is a very consistent property of A. 91 carbonarius and for this reason atoxigenic isolates of this species are very rarely found 92 (Cabañes et al., 2013). 4 93 In this study, we applied the High-throughput Next Generation Sequencing- 94 based RNA Seq method to carry out a global transcriptional analysis on four A. 95 carbonarius strains, one OTA producer and the three atypical non-OTA producing 96 strains. The aim of this study was to analyze the differentially expressed genes between 97 the OTA producer strain and the three non-ochratoxigenic strains to improve knowledge 98 on OTA-related gene expression and to identify genes directly or indirectly related to 99 OTA biosynthetic pathway. Besides, in vitro interactions between ochratoxigenic strains 100 of A. carbonarius and A. niger and non-ochratoxigenic strains of A. carbonarius and A. 101 tubingensis were assessed in order to evaluate their potential for control of OTA 102 production. 103 104 2. Material and methods 105 106 2.1. Strains and growing conditions 107 The OTA-producing strain A-1137 and three non-OTA-producing strains of A. 108 carbonarius A-2160, A-2579 and A-2594 from our fungal collection were analyzed in 109 this study. Strains were grown in triplicate in Czapek Yeast extract broth medium in the 110 dark at 25ºC for 48 h without shaking. After incubation, six plugs were removed from 111 different points of the colony for OTA determination. The rest of the mycelium was 112 removed and stored frozen at -80 ºC prior to total RNA extraction. 113 114 2.2. OTA quantification 115 OTA production was confirmed using a previously described HPLC screening 116 method designed in our laboratory (Bragulat et al., 2001). Plugs removed were extracted 5 117 with 0.5 ml of methanol. The extracts were filtered and maintained at 4ºC until their 118 analysis. 119 OTA quantification was made by a Waters 2695 chromatograph with a 120 fluorescence detector Waters 2475 (excitation wavelength: 330 nm/emission 121 wavelength: 460 nm), and with a Sunfire C18 column, 150×4.6 mm, i.d., 3.5 µm. 122 Twenty μl of each extract were applied. The mobile phase was acetonitril/water/acetic 123 acid (57/41/2, v/v/v) eluted at a flow rate of 1 ml/min. The extracts with the same 124 retention time as OTA (around 4.8 min), were considered positive. The limit of 125 quantification of the HPLC technique with the extraction procedure was 0.045 µg/ml for 126 OTA. 127 128 2.3. RNA extraction, RNA-Seq library preparation and sequencing 129 Frozen mycelium (100 mg) was powdered in liquid nitrogen and total RNA was 130 isolated using QIAzol lysis reagent (Qiagen, Madrid, Spain) and purified by using 131 RNeasy Plant Mini Kit (Qiagen, Madrid, Spain), following the manufacturer's 132 protocol.Contaminating genomic DNA was removed with DNAse (DNAse I, 133 amplification grade, Invitrogen, Carlsbad, CA, USA). The amount and quality of total 134 RNA was estimated by a Nanodrop 2000 spectrophotometer (Thermo Fisher Scientific 135 Inc., Wilmington, Delaware, USA) and a Bioanalyzer 2100 (Agilent Technologies, 136 Santa Clara, CA, USA). Next-generation sequencing including quality control was 137 performed by Sequentia Biotech SL (www.sequentiabiotech.com). Indexed libraries 138 were prepared from 1g of purified RNA with TruSeq Stranded mRNA Sample Prep 139 Kits (Illumina, San Diego, CA, USA) according to the manufacturer’s instructions. 140 Library size and integrity were determined using the TapeStation 4200 / Agilent 141 Bioanalyzer 2100 (Santa Clara, CA, USA). Libraries were pooled so that each index- 6 142 tagged sample was present in equimolar amounts, with a final concentration of 2nM. 143 The pooled samples were subject to cluster generation and sequencing using an Illumina 144 NextSeq 500 System in paired-end mode (2x150 bp). 145 146 2.4. RNA-Seq data analysis 147 Sequence analysis was performed using the A.I.R. software 148 (https://transcriptomics.sequentiabiotech.com/) developed by Sequentia Biotech. 149 Briefly, raw sequence files were first subjected to quality control analysis by using 150 FastQC v0.10.1 before trimming and removal of adapters with BBDuk 151 (https://jgi.doe.gov/data-and-tools/bbtools/). Reads were then mapped against the 152 Aspergillus carbonarius genome (ITEM 5010 v3 - JGI Genome Portal) with STAR v2.6 153 (Dobin et al., 2013). FeatureCounts v1.6.1 (Liao et al., 2014) was then used to obtain 154 raw expression counts for each annotated gene. The differential expression analysis was 155 conducted with the R package edgeR (Robinson et al., 2010) using the Trimmed mean 156 of M-values (TMM) normalization method. Fragments Per Kilobase Million (FPKM) 157 were obtained with edgeR. Gene Ontology Enrichment Analysis was performed using 158 in-house scripts based on the AgriGO publication (Tian et al., 2017). The analysis of 159 transcription factor binding motifs was performed with the software DMINDA2 (Li et 160 al., 2010) and the database JASPAR (Khan et al., 2017). 161 162 2.5. Validation of RNA-Seq results by RT-qPCR 163 The transcription profiles of five genes of the ocratoxin A biosynthesis cluster of 164 A. carbonarius (AcOTApks, AcOTAnrps, AcOTAhal, AcOTAp450 and AcOTAbZIP) 165 were analyzed in all strains by using real-time quantitative reverse transcription-PCR 166 (qRT-PCR) as described previously(Castellá et al., 2018). The RNA samples for each 7 167 replication and cDNA synthesis were run in triplicate. The comparative 2-Ct method 168 was used to calculate relative gene expression and all data was normalized to -tubulin 169 and ubiquitin. Expression difference among strains was assessed for statistical 170 significance using the REST 2009 v2.0.13 software (Pfaffl, Horgan and Dempfle, 171 2002). The calibrator sample corresponded to the value of expression of the OTA- 172 producing strain A-1137. 173 174 2.6. Confirmation of AcOTApks gene mutations 175 The OTA-producing strain (A-1137) and three non-OTA-producing strains of A. 176 carbonarius A-2160, A-2579 and A-2594 were grown on malt extract broth medium in 177 the dark at 25ºC for 48 h. Mycelium was recovered and DNA was extracted using the 178 FastDNA Spin kit (MP Biomedicals, Biolink, Barcelona, Spain). The DNA was kept at 179 -20 °C until used as template for PCR amplification. Primers for partial amplification of 180 AcOTApks gene (173482) were designed with Primer3Plus tool (Untergasser et al., 181 2007). Primers PKS1F (5'-GGATCATCCATGGAGGCCAG-3') and PKS1R (5'- 182 TATTGTGCGGGCCAGATCTC-3') and PKS2F (5'-GCCACAGGGGGAATAGTCTT- 183 3') and PKS2R (5'- CCAGGGAACAAGTCAGACCA-3') were designed to confirm the 184 five common mutations found in non-OTA producing strains and located at nucleotide 185 positions 948986, 949202, 950405, 950501 and 950904 of the AcOTApks gene. The 186 Big-Dye Terminator v3.1 Cycle Sequencing kit (Applied Biosystems, Foster City, CA, 187 USA) and primers PKS1F/R and PKS2F/R were used for sequencing as specified by the 188 manufacturer. Sequencing was performed using an Applied Biosystems 3730 analyzer. 189 The Clustal X v2.0.12 program was used to perform sequence alignments (Larkin et al., 190 2007). 8 191 To predict if the amino acid substitutions have an impact on the biological 192 function of the AcOTApks gene, we used the PROVEAN software tool (Choi and Chan, 193 2015). (http://provean.jcvi.org/index.php). 194 195 2.7. Assay of the in vitro interactions between ochratoxigenic strains of A. carbonarius 196 and A. niger and non-ochratoxigenic strains of A. carbonarius and A. tubingensis and 197 their effects on OTA production 198 In total, four strains of black aspergilli isolated from grapes were selected from 199 our culture collection. All the strains had been confirmed previously for identity by 200 DNA sequencing. An OTA-producing strain of A.carbonarius (A-1136) (Ac+) and an 201 OTA-producing strain of A. niger (A-1241) (An+) were used in this study. Two non- 202 OTA-producing strains, A. carbonarius (A-2579) (Ac-) and A. tubingensis (A-1747) 203 (At), used as potential biocontrol agents, were screened for their ability to inhibit or 204 reduce OTA production in co-inoculation with the OTA-producing strains (Ac+ and 205 An+). We designed an adaptation of a technique previously described for 206 ecophysiological studies of toxigenic fungal species using microtiter plates (Abarca et 207 al., 2014; 2019). Briefly, the inoculum suspensions of these strains were prepared in 208 sterile saline (0.85%) containing 0.05% Tween 80 from 7-day-old cultures on malt 209 extract agar (MEA) at 25ºC and spore concentration was adjusted to around 106 210 conidia/ml by using a haemocytometer chamber. The inoculum size was confirmed by 211 quantitative colony counts (CFU/ml) on MEA plates. 212 Sterile 96-well flat-bottom microtiter plates were used. The liquid culture 213 medium Yeast Extract Sucrose broth (YESB) (per liter, yeast extract, 20 g; sucrose, 150 214 g; FeSO4·7 H2O, 0.5 g; pH adjusted to 6.5) (Samson et al., 2000) was used. The 215 adjusted fungal suspensions were diluted 1:100 in YESB. 9 216 The effect of the potential biological agent on OTA production was evaluated 217 using different mixed spore suspensions of OTA-producing strain:non-OTA-producing 218 strain ratios of 100:0, 75:25, 50:50, 25:75 and 0:100 respectively in YESB, following a 219 modification of the technique described by Samsudin and Magan (2016) which used a 220 solid culture medium. In each microplate column, seven wells were inoculated with 200 221 µl of the diluted suspension of each inoculum ratio. One well, used as a blank, was 222 filled with 200 µl of the un-inoculated YESB. OTA production was determined after 5 223 and 11 days of incubation at 25ºC. The entire experiment was repeated twice on 224 different days. 225 OTA production was detected using a previously described HPLC screening 226 method developed in our laboratory for fungi growing in microtiter wells (Abarca et al., 227 2014). On each sampling occasion, the seven replicate wells inoculated for inoculum 228 ratio were removed and extracted with 0.5 ml of methanol. The extracts were filtered 229 and injected into the HPLC. The limit of quantification was 0.045 µg/ml for this 230 mycotoxin. 231 Data obtained from the different conditions tested were statistically analyzed by 232 means of one-way analysis of variance test. All statistical analyses were performed 233 using Minitab 17 statistical software (Minitab Inc., State College, Pennsylvania, USA). 234 235 3. Results 236 237 3.1. OTA production of the A. carbonarius strains studied 238 All strains presented good growth with proper sporulation forming typical black 239 colonies. A. carbonarius A-1137 produced OTA at detectable levels in the three 10 240 replicates, with a mean value of 5.73 ± 0.41 g/ml. The non-OTA producing strains A- 241 2160, A-2579 and A-2594 were not able to produce OTA in none of the three replicates. 242 243 3.2. RNA-Seq analysis 244 A summary of RNA-Seq data is reported in Table 1. In total, 18,510,000, 245 17,383,036, 17,438,795, and 18,276,380 high-quality (QS_30) paired-end reads (150 246 bp) were obtained from A-1137, A-2160, A-2579 and A-2594 transcriptome, 247 respectively, and deposited in the NCBI Sequence Read Archive (SRA) database under 248 study accession number PRJNA550023. Approximately 81% of short reads were 249 successfully mapped on the A. carbonarius genome. No significant differences were 250 observed in terms of total number of aligned reads among the replicates and the strains. 251 Most of the mapped reads were uniquely mapped, while the remaining reads were 252 unmapped (12–25%) or showed multi-position matches (6–10%). 253 254 3.3. Differential gene expression and functional analysis 255 Based on the TMM normalized values, differentially expressed genes (DEG) 256 were identified (with FDR ≤ 0.05, log2 fold change ≥1 or ≤-1) between the OTA 257 producing strain A- 1137 and each of the non-OTA producing strains. Analysis of the 258 transcriptional profiles revealed a total of 1,526 DEGs (710 up-regulated and 816 down- 259 regulated), 2,188 DEGs (1159 up-regulated and 1029 down-regulated), and 1,676 DEGs 260 (860 up-regulated and 816 down-regulated) in A-2160, A-2579, and A-2594, 261 respectively. A total of 696 differentially expressed genes were identified comparing the 262 OTA producing strain vs. the three non-OTA producing strains. Among these DEGs, 263 280 genes were up-regulated (Figure 1A) and 333 genes were down-regulated in the 11 264 three non-OTA producing strains (Figure 1B). The rest of differentially expressed genes 265 (83) were up-regulated or down-regulated depending on the strain. 266 These differentially expressed genes were subjected to GO functional 267 enrichment analysis (Supplementary Table S1). The GO terms over-represented in DEG 268 common in the three non-ochratoxigenic strains were related to metabolic process, 269 oxidation-reduction process and transport. Molecular function was related to catalytic 270 activity, oxidoreductase activity, heme binding, iron ion binding, and monooxygenase 271 activity. 272 Among up-regulated DEGs (Table 2, Figure 2), the GO enrichment analysis 273 revealed that several categories for biological processes were enriched among the DEGs 274 of non-OTA producing strains, including ribosome biogenesis, rRNA and tRNA 275 processing, polysaccharide catabolic process, ion transport, and oxidation-reduction 276 process. The oxidation-reduction processes were the most represented. As to molecular 277 function is concerned, monooxygenase activity, helicase activity, phosphopantetheine 278 binding, oxidoreductase activity, iron ion binding, heme binding, catalytic activity, 279 transaminase activity, and FAD binding (flavin adenine dinucleotide) were over- 280 represented. Among down-regulated DEGs (Table 3, Figure 2), over-represented 281 biological processes categories were oxidation-reduction and metabolic processes. As to 282 molecular function is concerned, oxidoreductase and hydrolase activities were over- 283 represented. 284 Within the most down-regulated genes (Supplementary Table S2), we found a 285 polyketide synthase (ID 173482) and one nonribosomal peptide synthetase (ID 132610). 286 The polyketide synthase ID 173482 corresponded to AcOTApks gene and the 287 nonribosomal peptide synthetase ID 1326010 corresponded to AcOTAnrps gene, both 288 included in the putative OTA cluster (cluster 38) of A. carbonarius. Additional down 12 289 regulated genes were putatively involved in biosynthesis of mycotoxins and other 290 secondary metabolites such as cytochrome P450 monooxygenases, monoxigenases, 291 dehydrogenases, hydrolases and methyltransferases. Different cytochrome P450 292 monoxigenases were down-regulated, but the most down-regulated in the three non- 293 ochratoxigenic strains was the cytochrome P450 monooxygenase ID 517149 which is 294 AcOTAp450 gene included in the OTA cluster. Moreover, the bZip transcription factor 295 AcOTAbZip (ID 7821) and a pepsin B (ID 7823) located in the same cluster were also 296 down-regulated. The halogenase AcOTAhal (ID 209543) was also down-regulated 297 although no significant differences were observed (Table 4). 298 Within the up-regulated genes (Supplementary Table S3), we also find genes 299 putatively involved in biosynthesis of secondary metabolites as a polyketide synthase 300 (ID 166447), some nonribosomal peptide synthetases (ID 133906, ID 503574), 301 cytochrome P450 monooxygenases, monoxigenases, dehydrogenases, hydrolases and 302 methyltransferases. 303 We analyzed the expression of 146 highly up-regulated genes describe in A. 304 carbonarius strains under OTA induction conditions (Gerin et al., 2016). In the non- 305 ochratoxigenic strains an oxalacetate acetylhidrolase (ID 158065) was down-regulated. 306 Also, three genes related with transport and secretion activities were down-regulated in 307 the non-ochratoxigenic strains. These genes were the MFS transporter (FLU1) (ID 308 135554), the ABC multidrug transporter (ID 171384) and the Formate/nitrate family 309 transporter (ID 132418). One gene related to transferase activity (methyltranferase ID 310 133219) and two genes related to oxidation-reduction processes (cytochrome P450 311 monoxygenase ID 517149 and multicopper oxidase ID 519260) were also down- 312 regulated in these strains. 13 313 We also analyzed expression of gene clusters co-regulated with cluster 38. In the 314 non-ochratoxigenic strains, no DEGs were found in clusters 24 and 48. In cluster 40 315 only 1 gene was down-regulated, a binding factor (ID 133215). In cluster 42 which 316 contained 17 genes, only two genes were DEG, a down-regulated glutation transferase 317 (ID 175671) and an up-regulated molybdenum cofactor sulfurase (ID 175665). 318 319 3.4. Validation of RNA-Seq results by RT-qPCR 320 RT-qPCR was performed on AcOTApks, AcOTAnrps, AcOTAhal, AcOTAp450 321 and AcOTAbZIP. Two genes, β-tubulin and ubiquitin, were used as reference genes. 322 Gene expression patterns of the selected DEGs were consistent with those obtained by 323 RNA-Seq (Figure 3), confirming the reliability and accuracy of the NGS analysis. 324 325 3.5. Transcription factors 326 By using Gene Ontology annotations and information in the OrthoMCL database 327 (https://orthomcl.org/orthomcl/), a total of 341 transcriptions factors were identified in 328 the Aspergillus carbonarius genome. Only 14 of them were differentially expressed, of 329 which five were down-regulated and nine were up-regulated (Table 5). Only three 330 transcription factors showed conserved motifs known to act as binding sites specific for 331 the putative OTA cluster genes AcOTApks (ID 173482), AcOTAnrps (ID 132610), 332 AcOTAp450 (ID 517149), and pepsin B (ID 7823) and all they were all down-regulated. 333 These transcriptions factors were the bZIP transcription factor AcOTAbZIP (ID 7821), 334 located in the OTA cluster, and two Zn2Cys6 transcription factors, (ID 131099, ID 335 514084) outside the OTA cluster. 336 337 3.6. AcOTApks mutations confirmation by Sanger sequencing 14 338 In order to confirm the five common missense variants found in AcOTApks gene 339 of these atoxigenic strains, we designed the primer pairs PKS1F/R and PKS2F/R. We 340 confirm the mutations located at positions 948986, 949202, 950405, 950501 and 341 950904 in scaffold 12. These missense variants produce amino acid substitutions that 342 have a neutral effect on biological functions of the AcOTApks with the exception of 343 mutation at position 949202 (Figure 4, Supplementary Figure S1). This mutation 344 produced an amino acid substitution (Y728H) that can have a deleterious impact on the 345 biological function of the PKS_AT domain of the AcOTApks. 346 347 3.7. Assay of the in vitro interactions between ochratoxigenic strains of A. carbonarius 348 and A. niger and non-ochratoxigenic strains of A. carbonarius and A. tubingensis and 349 their effects on OTA production 350 Non-OTA-producing strains of A. carbonarius and A. tubingensis (Ac- and At) 351 were screened for their ability to inhibit OTA production by ochratoxigenic strains of A. 352 carbonarius and A. niger (Ac+ and An+). Inocula enumerated with a cell counting 353 haemocytometer provided suspensions from 1.3 x 106 to 1.4 x 106 conidia/ml and 354 colony counts were between 4.0 x 105 and 6.6 x 105 CFU/ml. No statistically significant 355 differences were observed between replicates (p > 0.05). 356 The effects of OTA-producing strain:non-OTA-producing strain inoculum ratios 357 on OTA production are detailed in Table 6. The results were expressed as a mean of the 358 14 replicates. No statistically significant differences were observed in OTA production 359 values neither between replicates nor between experiments (p > 0.05). In most cases, no 360 significant differences were observed in OTA production values between the two 361 incubation times (5 and 11 days). 15 362 Regarding OTA production, in all cases of co-inoculation the presence of Ac- 363 and At decreased the production at different percentages, depending on their inoculum 364 load. As more load of the non-ochratoxigenic strains assayed, the higher percentages of 365 OTA reduction were detected. 366 Of the two non-ochratoxigenic species assayed, Ac- gave the best control 367 resulting in practically complete inhibition of OTA production in co-inoculation with 368 An+, at all inoculum ratios at the two incubation times. Figure 5A shows some selected 369 chromatograms of extracts of the interactions between these species at different 370 inoculum ratios. Positive control of OTA production (100 An+:0 Ac-) presented a clear 371 peak of OTA (retention time of 5.2 minutes). A small peak with the same retention time 372 is observed at ratio 75:25. The remaining inoculum ratios showed no signals at the same 373 retention time of OTA. 374 For Ac-, OTA reduction was higher than 74% in co-inoculation with Ac+, at all 375 strain ratios at the two incubation times (see Figure 5B). Percentages of OTA reduction 376 ranged from 96% to 100% when At was co-inoculated with An+ at all inoculum ratios 377 at the two incubation times. On the contrary, co-inoculation of At with Ac+ showed 378 lower percentages of OTA reduction. 379 380 4. Discussion 381 Despite some studies have been carried out on molecular aspects of OTA 382 biosynthesis, the OTA cluster remains not completely defined and most of the 383 regulatory aspects underlying OTA production remain unclear. In A. ochraceus a cluster 384 that contains four biosynthetic genes and two regulators have been related to OTA 385 production (Wang et al., 2018). The biosynthetic genes were a polyketides sinthase 386 (OtaA), a nonribosomal peptide synthetase (Ota B), a cytochrome P450 (Ota C) and a 16 387 halogenase (OtaD). Two regulators have been described, a bZIP transcription factor 388 (otaR1), which probably functions as an OTA pathway-specific regulator, and a 389 Zn2Cys6 binuclear DNA-binding protein (ota R2), which is adjacent to biosynthetic 390 genes. All the biosynthetic genes of the OTA cluster of A. ochraceus were present in A. 391 carbonarius in cluster 38, and corresponded to AcOTApks, AcOTAnrps, AcOTAp450, 392 AcOTAhal genes. In A. carbonarius the transcription factor AcOTAbZIP which is 393 located between AcOTAp450 and AcOTAhal is homologous to otaR1 found in A. 394 ochraceus. The transcriptome analysis of A. carbonarius strains under OTA inducing 395 conditions showed that AcOTApks, AcOTAnrps, AcOTAp450, and AcOTAhal were 396 highly up-regulated in OTA inducing conditions and all were co-expressed in addition 397 to the AcOTAbZIP (Gerin et al., 2016) supporting the fact that the cluster 38 is likely the 398 OTA gen cluster in A. carbonarius. 399 In the present study, we applied the RNA-seq to study the whole transcriptional 400 changes associated to OTA production in three non-OTA producing strains of A. 401 carbonarius. Functional analysis and enrichment GO terms in DEG revealed that the 402 processes related to metabolic and oxidation processes, associated with functions such 403 oxidoreductase and hydrolase activity were down regulated, indicating that changes in 404 the fungal metabolism occurred in the three non-OTA producing strains. Within the 405 most down-regulated genes, we found the putative OTA genes described in A. 406 carbonarius in cluster 38. The expression level of all the genes related to OTA 407 biosynthesis was down-regulated. These genes corresponded to AcOTApks (ID 173482), 408 AcOTAnrps (ID 132610), AcOTAp450 (ID 517149), AcOTAhal (ID 209543), and 409 AcOTAbZIP (ID 7821). 410 The AcOTAbZip has a high similarity to otaR1 and could function as a cluster 411 specific regulator also in A. carbonarius. In fact, this transcription factor showed 17 412 conserved motifs known to act as binding sites of some of the putative OTA genes in A. 413 carbonarius. AcOTAbZip down-regulation may explain the lack of OTA production in 414 the three non-ochratoxigenic strains. The bZIP transcriptions factors typically form 415 homo or hetero-dimers and bind to the promoter regions to regulate the transcription of 416 their target genes. Some bZIP transcriptions factors have been described to regulate 417 mycotoxin production directly by binding at promoters of genes encoding key elements 418 of the biosynthetic pathway, as shown for RsmA regulating sterigmatocystin in A. 419 nidulans (Yin et al., 2012) and gliotoxin in A. fumigatus (Xiao et al., 2010) and AtfB 420 regulating aflatoxin in Aspergillus parasiticus (Roze et al., 2011). In A. ochraceus, the 421 bZIP otaR1 control the expression of the four OTA biosynthetic genes, particularly 422 otaA but did not regulate the expression of otaR2, a second regulator described in the 423 OTA cluster of this species (Wang et al., 2018). The regulatory mechanisms of otaR1 424 are not known. The second regulator found in OTA cluster, ota R2, can modulate the 425 expression of otaA, otaC and otaD but did not regulate the expression of otaB, otaE or 426 otaR1. In A. carbonarius, toward the edges of the OTA cluster, two genes encoding 427 transcription factors have been described (Gallo et al., 2017). These transcription factors 428 (ID 173470, ID 132614) contain Cys –rich motives which are present in zinc-cluster 429 proteins (Zn2Cys6), the most common regulators of fungal secondary metabolism 430 (Brakhage, 2013). In our study, both of them were not differentially expressed in the 431 non-ochratoxigenic strains. We also found two Zn2Cys6 transcriptions factors outside 432 the genomic region of the OTA cluster and with binding sites to AcOTApks, AcOTAnrps 433 and AcOTAp450. Both of them were down-regulated in non-ochratoxigenic strains. 434 These transcription factors (ID 131099 and ID 514084) were similar to other 435 transcription factors present in different Aspergillus species and their function is 436 unknown. 18 437 The genome analysis of our non-ochratoxigenic strains revealed five private 438 missense variants in AcOTApks gene (Castellá et al., 2018). In the present study these 439 mutations of AcOTApks have been confirmed by Sanger sequencing. One of the 440 predicted non-conservative amino acid change was deleterious to enzymatic activity of 441 the AcOTApks. OTA biosynthesis begins with the AcOTApks gene, utilizing acetyl 442 coenzime A and malonyl-CoA to synthesize 7-methylmellein. This first step could be 443 disrupted in the non-ochratoxigenic strains due to the deleterious mutation found in the 444 acyltransferase domain of the AcOTApks. The nucleotide polymorphism (A->G) at nt 445 949202 in AcOTApks of non-ochratoxigenic strains is sufficient to account for the lack 446 of OTA production by these strains. A similar situation has been reported in some 447 atoxigenic strains of A. flavus, where a nucleotide polymorphism near the beginning of 448 the code sequence of pksA gene prevent enzyme production and aflatoxin accumulation 449 (Ehrlich and Cotty, 2004). 450 In the non-ochratoxigenic strains the FLU1 transporter (ID 135554) was also 451 significantly down-regulated. This MFS transporter has been putatively involved in 452 OTA transport in A. carbonarius (Crespo-Sempere et al., 2010). It has been shown to be 453 significantly up-regulated under OTA induction conditions (Crespo-Sempere et al., 454 2010; Gerin et al., 2016) although its role in OTA transport has not been elucidated. 455 Moreover, the transcript ID158065, a homologous to the A. niger ANI_1_92174 gene 456 encoding oxalacetate acetylhydrolase, was also significantly down-regulatd in non- 457 ochratoxingenic strains. This enzyme converts oxalacetate to oxalate and acetate, which 458 is a potential source of polyketide precursors (Kubicek et al., 1988). 459 A. carbonarius is the main responsible source of OTA in wine or dried vine 460 fruits from main viticultural regions worldwide and it is also considered a potential 461 source of OTA in coffee and cocoa, among other foods. Besides, it is very consistent in 19 462 producing OTA, and non-OTA-producing strains in this species are very rare. On the 463 contrary, A. niger and A. tubingensis have a worldwide distribution and are very 464 frequently isolated from grapes. However, the reported percentage of OTA-producing 465 strains in A. niger is very low whereas A. tubingensis is considered a non- 466 ochratoxigenic species (Cabañes and Bragulat, 2018; Frisvad et al., 2011). 467 In our study, Ac- gave the best control resulting in practically complete 468 inhibition of OTA production in co-inoculation with An+, and high percentages of OTA 469 reduction with Ac+. On the other hand, almost complete inhibition of OTA production 470 was also detected by At in co-inoculation with An+. 471 These preliminary results obtained with these no-toxigenic wild strains are very 472 promising for the OTA control. An approach would be the use of these wild atoxigenic 473 black aspergilli to outcompete toxigenic A. carbonarius strains in the field and minimise 474 OTA contamination. Biological control using microbial antagonists either alone or as 475 part of an integrated control strategy to reduce pesticide inputs, has emerge as a 476 promising approach for control of mycotoxins in crops, both pre- and post-harvest. 477 However, other studies are needed in order to know if these strains are able to work as 478 effective BCAs in the field. Despite the fact that the use of BCAs has often been 479 successful in the lab, their commercialization largely depends on whether they can 480 consistently control mycotoxigenic fungi in different locations and cultivars and reduce 481 mycotoxin levels below the legislative limits (Ehrlich, 2014). These strains are available 482 and our research group is open for collaboration through a transfer agreement of the 483 material to industry. 484 485 5. Conclusions 20 486 The RNA-seq analysis of the three atypical non-ochratoxigenic strains of A. 487 carbonarius showed that all the genes described to be related with OTA biosynthesis 488 were down-regulated. Also, the AcOTApks gene showed a deleterious mutation that it is 489 sufficient to prevent OTA production. Moreover, one of these strains gave the best 490 control resulting in practically complete inhibition of OTA production in co-inoculation 491 with an ochratoxigenic strain of A. niger and high percentages of OTA reduction with 492 an ochratoxigenic strain of A. carbonarius. This study represents an important step 493 forward to understand the OTA biosynthetic pathway in these non-ochratoxigenic wild 494 strains, which can lead to the development of improved control strategies to reduce the 495 risk of OTA contamination in food products. 496 497 Acknowledgments 498 We thank C. Gómez for her valuable technical assistance. 499 500 Funding 501 This research was supported by the Ministerio de Economía y Competitividad of 502 the Spanish Government (AGL2014- 52516-R). 503 504 Supplementary material 505 Supplementary information accompanies this paper at 506 507 Conflict of Interest Statement 508 RAC was employed by company Sequentia Biotech SL. All other authors 509 declare no competing interests. 510 21 511 References 512 Abarca, M. L., Bragulat, M. R., Cabañes, F.J., 2014. A new in vitro method to detect 513 growth and ochratoxin A-producing ability of multiple fungal species commonly 514 found in food commodities. Food Microbiol. 44, 243-248. doi: 515 10.1016/j.fm.2014.06.014 516 Abarca, M. L., Bragulat, M. R., Castellá, G., Cabañes, F.J., 2019. Impact of some 517 environmental factors on growth and ochratoxin A production by Aspergillus niger 518 and Aspergillus welwitschiae. Int. J. Food Microbiol. 291, 10-16. doi: 519 10.1016/j.ijfoodmicro.2018.11.001 520 Alaniz Zanon, M.S., Clemente, M.P., Chulze, S.N., 2018. Characterization and 521 competitive ability of non-aflatoxigenic Aspergillus flavus isolated from the maize 522 agro-ecosystem in Argentina as potential aflatoxin biocontrol agents. Int. J. Food 523 Microbiol. 277, 58-63. doi: 10.1016/j.ijfoodmicro.2018.04.020. 524 Andersen, M. R., Salazar, M. P., Schaap, P. J., van de Vondervoort, P. J. I., Culley, D., 525 Thykaer, J., Frisvad, J.C., Nielsen, K.F., Albang, R., Albermann, K., Berka, R.M., 526 Braus, G.H., Braus-Stromeyer, S.A., Corrochano, L.M., Dai, Z., van Dijck, P.W., 527 Hofmann, G., Lasure, L.L., Magnuson, J.K., Menke, H., Meijer, M., Meijer, S.L., 528 Nielsen, J.B,, Nielsen, M.L., van Ooyen, A.J., Pel, H.J., Poulsen, L., Samson, R.A., 529 Stam, H., Tsang, A., van den Brink, J.M., Atkins, A., Aerts, A., Shapiro, H., 530 Pangilinan, J., Salamov, A., Lou, Y., Lindquist, E., Lucas, S., Grimwood, J., 531 Grigoriev, I.V., Kubicek, C.P., Martinez, D., van Peij, N.N., Roubos, J.A., Nielsen, 532 J., Baker, S.E., 2011. Comparative genomics of citric-acid-producing Aspergillus 533 niger ATCC 1015 versus enzyme-producing CBS 513.88. Genome Res. 21, 885- 534 897. doi: 10.1101/gr.112169.110 22 535 Bragulat, M. R., Abarca, M. L., Cabañes, F. J., 2001. An easy screening method for 536 fungi producing ochratoxin A in pure culture. Int. J. Food Microbiol. 71, 139-144. 537 Cabañes, F.J., Bragulat, M.R., 2018. Black aspergilli and ochratoxin A-producing 538 species in foods. Curr. Opin. Food Sci. 23, 1-10. doi: 10.1016/j.cofs.2018.01.006 539 Cabañes, F. J., Bragulat, M. R., Castellá, G., 2013. Characterization of 540 nonochratoxigenic strains of Aspergillus carbonarius from grapes. Food Microbiol. 541 36, 135–141. doi: 10.1016/j.fm.2013.05.004 542 Cabañes, F. J., Sanseverino, W., Castellá, G., Bragulat, M. R., Cigliano, R. A., Sanchez, 543 A., 2015. Rapid genome resequencing of an atoxigenic strain of Aspergillus 544 carbonarius. Sci. Rep. 5, 9086. doi: 10.1038/srep09086 545 Castellá, G., Bragulat, M.R., Puig, L., Sanseverino, W., Cabañes, F.J., 2018. Genomic 546 diversity in ochratoxigenic and non ochratoxigenic strains of Aspergillus 547 carbonarius. Sci. Rep. 8, 5439. doi: 10.1038/s41598-018-23802-8 548 Choi, Y., Chan, A.P., 2015. PROVEAN web server: a tool to predict the functional 549 effect of amino acid substitutions and indels. Bioinformatics 31, 2745-2747. doi: 550 10.1093/bioinformatics/btv195 551 Choque, E., Klopp, C., Valiere, S., Raynal, J., Mathieu, F., 2018. Whole-genome 552 sequencing of Aspergillus tubingensis G131 and overview of its secondary 553 metabolism potential. BMC Genomics 19: 200. doi: 10.1186/s12864-018-4574-4 554 Crespo-Sempere, A., González-Candelas, L., Martínez-Culebras, P.V., 2010. Genes 555 differentially expressed by Aspergillus carbonarius strains under ochratoxin A 556 producing conditions. Int. J. Food Microbiol. 142, 170-179. doi: 557 10.1016/j.ijfoodmicro.2010.06.019 558 Brakhage, A. A., 2013. Regulation of fungal secondary metabolism. Nature Rev. 11, 559 21-32. doi: 10.1038/nrmicro2916 23 560 Dobin, A., Davis, C.A., Schlesinger, F., Drenkow, J., Zaleski, C., Jha, S., Batut, P., 561 Chaisson, M., Gingeras, T.R., 2013. STAR: ultrafast universal RNA-seq aligner. 562 Bioinformatics 29, 15-21. doi: 10.1093/bioinformatics/bts635 563 Ehrlich, K. C., 2008. Genetic diversity in Aspergillus flavus and its implications for 564 agriculture. In: Varga, J. and Samson, R.A. (Eds.), Aspergillus in the Genomic Era, 565 Wageningen Academic Publishers, Wageningen, pp. 233–247. 566 Ehrlich, K.C., 2014. Non-aflatoxigenic Aspergillus flavus to prevent aflatoxin 567 contamination in crops: advantages and limitations. Front. Microbiol. 5, 50. doi: 568 10.3389/fmicb.2014.00050. doi: 10.3389/fmicb.2014.00050 569 Ehrlich K.C., Cotty, P.J., 2004. An isolate of Aspergillus flavus used to reduce aflatoxin 570 contamination in cottonseed has a defective polyketide synthase gene. Appl. 571 Microbiol. Biotechnol. 65, 473-478. 572 Ferrara, M., Perrone, G., Gambacorta, L., Epifani, F., Solfrizzo, M., Gallo, A., 2016. 573 Identification of a halogenase involved in the biosynthesis of ochratoxin A in 574 Aspergillus carbonarius. Appl. Environ. Microbiol. 82, 5631-5641. doi: 575 10.1128/AEM.01209-16 576 Frisvad, J.C., Larsen, T.O., Thrane,U., Meijer, M., Varga, J., Samson, R.A., Nielsen, 577 K.F., 2011. Fumonisin and ochratoxin production in industrial Aspergillus niger 578 strains. PLoS ONE 6, e23496. doi: 10.1371/journal.pone.0023496 579 Gallo, A., Bruno, K. S., Solfrizzo, M., Perrone, G., Mulè, G., Visconti, A., Baker, S. E., 580 2012. New insight in the ochratoxin A biosynthetic pathway by deletion of an nrps 581 gene in Aspergillus carbonarius. Appl. Environ. Microbiol. 78, 8208-8218. doi: 582 10.1128/AEM.02508-12 24 583 Gallo, A., Ferrara, M., Perrone, G., 2017. Recent advances on the molecular aspects of 584 ochratoxin A biosynthesis. Curr Opin Food Sci. 17, 49-56. doi: 585 10.1016/j.cofs.2017.09.011 586 Gallo, A., Knox, B.P., Bruno, K.S., Solfrizzo, M., Baker, S.E., Perrone, G., 2014. 587 Identification and characterization of the polyketide synthase involved in 588 ochratoxin A biosynthesis in Aspergillus carbonarius. Int. J. Food Microbiol. 179, 589 10–17. doi: 10.1016/j.ijfoodmicro.2014.03.013 590 Gerin, D., De Miccolis Angelini, R.M., Pollastro, E., Faretra, F., 2016. RNA-Seq 591 reveals OTA-related gene transcriptional changes in Aspergillus carbonarius. PLoS 592 One 11, e0147089. doi: 10.1371/journal.pone.0147089 593 Khan, A., Fornes, O., Stigliani, A., Gheorghe, M., Castro-Mondragon, J.A., van der 594 Lee, R., Bessy, A., Chèneby, J., Kulkarni, S.R., Tan, G., Baranasic, D., Arenillas, 595 D.J., Sandelin, A., Vandepoele, K., Lenhard, B., Ballester, B., Wasserman, W.W., 596 Parcy, F., Mathelier, A., 2018. JASPAR 2018: update of the open-access database 597 of transcription factor binding profiles and its web framework. Nucleic Acids Res. 598 46, D260–D266. doi: 10.1093/nar/gkx1126 599 Kubicek CP, Schreferl-Kunar G, Wöhrer W, Röhr M., 1988. Evidence for a cytoplasmic 600 pathway of oxalate biosynthesis in Aspergillus niger. Appl. Environ. Microbiol. 54, 601 633-637. 602 Larkin, M. A. Blackshields, G., Brown, N.P., Chenna, R., McGettigan, P.A., 603 McWilliam, H., Valentin, F., Wallace, I.M., Wilm, A., Lopez, R., Thompson, J.D., 604 Gibson, T.J., Higgins, D.G., 2007. Clustal W and Clustal X version 2.0. 605 Bioinformatics 23, 2947–2948. 25 606 Li, G., Liu, B., Xu, Y., 2010. Accurate recognition of cis -regulatory motifs with the 607 correct lengths in prokaryotic genomes. Nucleic Acids Res. 38 e12. doi: 608 10.1093/nar/gkp907 609 Liao, Y., Smyth, G.K., Shi W., 2014. featureCounts: an efficient general purpose 610 program for assigning sequence reads to genomic features. Bioinformatics 30, 923– 611 930. doi:10.1093/bioinformatics/btt656 612 Pfaffl, M. V., Horgan, G. W., Dempfle, L., 2002. Relative expression software tool 613 (REST©) for group-wise comparison and statistical analysis of relative expression 614 results in real-time PCR. Nucleic Acids Res. 30, e36. doi: 10.1093/nar/30.9.e36 615 Pfohl-Leszkowicz, A., Manderville, R. A., 2012. An update on direct genotoxicity as a 616 molecular mechanism of ochratoxin a carcinogenicity. Chem. Res. Toxicol. 25, 617 252-262. doi: 10.1021/tx200430f 618 Pitt, J.I., Hocking, A.D., 2009. Fungi and Food Spoilage, Springer, New York. 619 Robinson, M.D., McCarthy, D.J., Smyth, G.K., 2010. edgeR: a Bioconductor package 620 for differential expression analysis of digital gene expression data. Bioinformatics 621 26, 139-140. doi: 10.1093/bioinformatics/btp616 622 Roze, L.V., Chanda, A., Wee, J., Awad, D., Linz, J.E., 2011. Stress-related transcription 623 factor AtfB integrates secondary metabolism with oxidative stress response in 624 aspergilli. J. Biol. Chem. 286, 35137-35148. doi: 10.1074/jbc.M111.253468 625 Samson, R. A., Hoekstra, E.S., Frisvad, J. C., Filtenborg, O., 2000. Introduction to 626 Food- and Airborne Fungi, Centraalbureau voor Schimmelcultures, Utrech. 627 Samsudin, N.P.I., Magan, N., 2016. Efficacy of potential biocontrol agents for control 628 of Fusarium verticillioides and fumonisin B1 under different environmental 629 conditions. World Mycotoxin J. 9, 205-213. doi: 10.3920/WMJ2015.1886 26 630 Tian, T., Liu, Y., Yan, H., You, Q., Yi, X., Du, Z., Xu, W., Su, Z., 2017. agriGO v2.0: a 631 GO analysis toolkit for the agricultural community, 2017 update. Nucleic Acids 632 Res. 45, W122–W129. doi:10.1093/nar/gkx382 633 Untergasser, A., Cutcutache, I., Koressaar, T., Ye,J., Faircloth, B.C., Remm, M., and 634 Rozen, S.G., 2012. Primer3 - new capabilities and interfaces. Nucleic Acids Res. 635 40, e115. doi: 10.1093/nar/gks596 636 Wang, Y., Wang, L., Liu, F., Wang, Q., Selvaraj, J.N., Xing, F., Zhao, Y., Liu, Y., 637 2016. Ochratoxin A producing fungi, biosynthetic pathway and regulatory 638 mechanisms. Toxins 8, 83. doi: 10.3390/toxins8030083 639 Wang, Y., Wang, L., Wu, F., Liu, F., Wang, Q., Zhang, X., Selvaraj, J.N., Zhao, Y., 640 Xing, F., Yin, W.B., Liu, Y., 2018. A consensus Ochratoxin A biosynthetic 641 pathway: insights from the genome sequence of Aspergillus ochraceus and a 642 comparative genomic analysis. Appl. Environ. Microbiol. 84, e01009-18. doi: 643 10.1128/AEM.01009-18 644 Xiao, P., Shin, K.-S., Wang, T., Yu, J.-H., 2010. Aspergillus fumigatus flbB encodes 645 two basic leucine zipper domain (bZIP) proteins required for proper asexual 646 development and gliotoxin production. Eukaryot. Cell 9, 1711-1723. doi: 647 10.1128/EC.00198-10 648 Yin, W. B., Amaike, S., Wohlbach, D.J., Gasch, A.P., Chiang, Y. M., Wang, C. C. C., 649 Bok, J.W., Rohlfs, M., Keller, N.P., 2012. An Aspergillus nidulans bZIP response 650 pathway hardwired for defensive secondary metabolism operates through aflR. 651 Mol. Microbiol. 83, 1024-1034. doi: 10.1111/j.1365-2958.2012.07986.x 652 27 653 Figure captions 654 655 Fig. 1. Venn diagram of (A) up-regulated differentially expressed genes (DEGs) (log2 656 FC ≥1) and (B) down-regulated DEGs (log2 FC ≤-1) comparing the ochratoxin A 657 (OTA) producing strain A-1137 vs. the three non-OTA producing strains A-2160, A- 658 2579, A-2594. 659 660 Fig. 2. Most represented GO terms in functional enrichment analyses for differentially 661 expressed genes down- and up-regulated. 662 663 Fig. 3. Relative expression analysis by real time PCR of AcOTApks, AcOTAnrps, 664 AcOTAhal, AcOTAbZIP and AcOTAp450 genes in A. carbonarius strains A-2160, 665 A2579 and A-2594 grown on Czapek Yeast extract broth. (*P < 0.05). 666 667 Fig. 4. Alignment of a portion of the deduced amino acid sequences of AcOTApks gene 668 of Aspergillus carbonarius strains assayed. Identical residues are indicated by dots. a 669 Amino acid sequence of the reference strain A. carbonarius ITEM 5010. 670 671 Fig. 5. Selected HPLC-FLD chromatograms of the extracts of (A) the OTA-producing 672 strain A. niger A-1241 in co-inoculation with the non-OTA-producing strain A. 673 carbonarius A-2579, and (B) the OTA-producing strain A. carbonarius A-1136 in co- 674 inoculation with the non-OTA-producing strain A. carbonarius A-2579 (ND: Not 675 detected. No signals at the same retention time of OTA). 28 676 Table 1. 677 Statistics of the Illumina 150 bp paired-end reads and mapping on A. carbonarius genome sequence. 678 Genome Biological No. total read N. mapped reads Strain replicates pairs Total mapped Unique match Multi-position matches reads (%) A-1137 1 17,286,120 82,34% 12,949,113 1,283,919 2 17,540,106 82,88% 13,375,126 1,161525 3 20,703,776 83,23% 15,643,222 1,589,570 A-2160 1 18,876,503 87,42% 15,277,532 1,224,674 2 18,190,280 87,50% 14,731,510 1,184,859 3 15,082,325 80,67% 11,247,505 919,146 A-2579 1 18,880,195 82,10% 14,141,545 1,360,033 2 15,033,474 78,00% 10,732,980 992,571 3 18,402,717 74,30% 11,702,543 1,971,291 A-2594 1 21,905,338 81,87% 16,362,498 1,572,395 2 15,774,483 76,41% 10,979,552 1,073,423 3 17,149,319 78,14% 12,258,516 1,142,589 679 29 680 Table 2. 681 Enrichment analysis of GO terms up-regulated in differentially expressed genes common in the three non-ochratoxigenic strains. 682 GO ID Class Description Query Item p-value FDR Enrichment Score GO:0042254 biological_process ribosome biogenesis 5 1,01E-03 5,77E-03 4,28 GO:0006364 biological_process rRNA processing 5 1,42E-03 7,97E-03 4,01 GO:0000272 biological_process polysaccharide catabolic process 3 2,67E-03 1,17E-02 5,02 GO:0008033 biological_process tRNA processing 3 1,64E-02 3,63E-02 3,04 GO:0006811 biological_process ion transport 3 2,12E-02 4,34E-02 2,82 GO:0055114 biological_process oxidation-reduction process 42 2,13E-02 4,34E-02 1,33 GO:0005730 cellular_component nucleolus 4 1,09E-02 2,71E-02 2,96 GO:0005886 cellular_component plasma membrane 3 1,50E-02 3,40E-02 3,12 GO:0004497 molecular_function monooxygenase activity 13 3,24E-04 4,21E-03 2,70 GO:0004386 molecular_function helicase activity 7 6,96E-04 4,26E-03 3,59 GO:0031177 molecular_function phosphopantetheine binding 5 7,88E-04 4,71E-03 4,47 GO:0016705 molecular_function oxidoreductase activity, acting on paired donors, with 10 2,11E-03 9,59E-03 2,52 incorporation or reduction of molecular oxygen GO:0005506 molecular_function iron ion binding 11 2,21E-03 9,86E-03 2,39 GO:0020037 molecular_function heme binding 11 3,79E-03 1,44E-02 2,24 GO:0003824 molecular_function catalytic activity 30 5,69E-03 1,79E-02 1,55 GO:0008483 molecular_function transaminase activity 3 1,24E-02 2,96E-02 3,30 GO:0071949 molecular_function FAD binding 5 1,86E-02 3,94E-02 2,38 683 684 30 685 Table 3. 686 Enrichment analysis of GO terms down-regulated in differentially expressed genes common in the three non-ochratoxigenic strains. 687 Query Query Background Backgroun Enrichment GO ID Class Description p-value FDR Total Item Total d Item Score oxidation-reduction GO:0055114 biological_process process 283 57 9428 1218 1,94E-04 1,84E-03 1,56 GO:0070941 biological_process eisosome assembly 283 3 9428 17 1,39E-03 7,43E-03 5,88 GO:0008152 biological_process metabolic process 283 27 9428 590 1,06E-02 3,00E-02 1,52 GO:0032126 cellular_component eisosome 283 3 9428 18 1,74E-03 9,12E-03 5,55 GO:0005576 cellular_component extracellular region 283 4 9428 49 1,53E-02 3,99E-02 2,72 GO:0016407 molecular_function acetyltransferase activity 283 3 9428 6 1,14E-05 1,17E-04 16,66 GO:0016491 molecular_function oxidoreductase activity 283 44 9428 894 3,50E-04 3,08E-03 1,64 GO:0016810 molecular_function hydrolase activity, acting 283 4 9428 36 4,13E-03 1,52E-02 3,70 on carbon-nitrogen (but not peptide) bonds GO:0016787 molecular_function hydrolase activity 283 38 9428 858 5,55E-03 1,77E-02 1,48 GO:0010181 molecular_function FMN binding 283 4 9428 48 1,40E-02 3,71E-02 2,78 688 31 689 Table 4. 690 Down-regulated genes in ochratoxin A (OTA) cluster in the three non-OTA producing strains. 691 Log2 FC Transcript ID Description A2160 A2579 A2594 173482 Polyketide synthase -8,33 -7,55 -8,52 132610 Nonribosomal peptide synthetase -8,07 -6,30 -7,57 209543 Halogenase -8,97 -6,47 -7,15 517149 Cytochrome P450 -6,87 -6,12 -6,46 7821 BZIP transcription factor -5,46 -5,62 -4,48 7823 Pepsin B -2,60 -2,12 -1,68 692 32 693 Table 5. 694 Transcriptions factors differentially expressed in non-ochratoxigenic strains. 695 Transcript ID OrthoGroup Description A2160 A2579 A2594 7821 OG5_147188 BZIP transcription factor -5,46 -5,62 -4,48 131099 OG5_155600 Transcription factor -2,77 -2,05 -2,08 154243 OG5_188739 Transcription factor -1,78 -1,76 -1,58 13310 OG5_159084 Transcription factor -1,4 -1,21 -1,87 514084 OG5_175802 C6 zinc finger domain protein -1,15 -0,92 -0,78 143623 OG5_127437 C2H2 transcription factor 0,5 0,48 0,75 146632 OG5_139332 C6 zinc finger domain protein 0,61 0,57 0,6 144278 OG5_140254 C6 transcription factor (Gal4) 0,99 1,46 1,44 33782 OG5_188739 Transcription factor 1,22 1,59 1,11 513323 OG5_152733 Xylanolytic transcriptional activator xlnR 1,41 1,76 1,62 132388 OG5_176156 BZIP transcription factor 1,56 3,08 2,29 125349 OG5_149802 BZIP family transcription factor 1,66 2,09 2,02 130875 OG5_187896 Transcription factor 1,81 1,81 1,84 145697 OG5_159360 Transcription factor 3,09 1,5 1,91 696 33 697 Table 6. 698 Mean values of ochratoxin A (OTA) concentration in µg / ml produced by the strains assayed in the two experiments. Co-inoculated OTA days strains 100 / 0 † 75/25 (%R) 50 / 50 (%R) 25 / 75 (%R) 0 / 100 a b b,c b,c c Ac + / Ac - 5 21.25 ± 11.0 5.42 ± 2.78 (74%) 3.3 ± 2.1 (84%) 2.4 ± 2.6 (89%) ND a b c c d 11 30.5 ± 4.7 7.64 ± 3.63 (75%) 3.6 ± 3.2 (88%) 3.4 ± 3.4 (89%) ND a b b b b An + / Ac - 5 6.49 ± 2.13 0.10 ± 0.24 (98%) ND (100%) ND (100%) ND a b b b b 11 3.42 ± 2.04 0.014 ± 0.027 (99%) ND (100%) ND (100%) ND a b b b b An + / At 5 10.54 ± 3.68 0.45 ± 0.57 (96%) 0.023 ± 0.028 (99.8%) ND (100%) ND a b b b b 11 5.14 ± 2.61 0.13 ± 0.16 (97%) 0.009 ± 0.02 (99.8%) ND (100%) ND a a a,b b c Ac + / At 5 8.87 ± 4.46 7.86 ± 3.43 (11%) 6.72 ± 3.34 (24%) 4.96 ± 3.35 (44%) ND a b b c d 11 19.76 ± 4.91 10.89 ± 4.59 (45%) 9.51 ± 2.25 (52%) 6.55 ± 3.15 (67%) ND †, spore suspension ratios of each strain (in %). (%R), denotes de percentage of OTA reduction. a,b,c,d In files, values with the same superscript are not significantly different (p > 0.05). Ac, A. carbonarius; An, A. niger; At, A. tubingensis; +, OTA-producing strain; -, non-OTA-producing strain 34 Figure 1 A-2160 vs. A-1137 A-2579 vs. A-1137 A A-2594 vs. A-1137 A-2160 vs. A-1137 A-2579 vs. A-1137 B A-2594 vs. A-1137 35 Figure 2 36 Figure 3 0 -2 log2 fold change -4 A2160 * A2579 * * A2594 -6 * * * -8 * * * * * * -10 * * * AcOTAbZIP AcOTAhal AcOTAnrps AcOTAp450 AcOTApks 37 Figure 4 671 770 ....|....|....|....|....|....|....|....|....|....|....|....|....|....|....|....|....|....|....|....| PKS173482a QISRQAAWKVAFCRGQVCARRTDGQGRMLAAAMPVTQLEQLVARVNKGQSTAVKVGCYNSPKNLTLTGRAEDILRAKLELDDVGALNRLLPVKVAYHSDY A-1137 .................................................................................................... A-2160 .........................................................H.......................................... A-2579 .........................................................H.......................................... A-2594 .........................................................H.......................................... 771 836 ....|....|....|....|....|....|....|....|....|....|....|....|....|. PKS173482a MRDAAPEYLDLLGDLDFGDSIHADAGIKMVSSVTGRAVSAGEAQQPSYWVDNLVSPVRFSTALLAS A-1137 .................................................................. A-2160 .............................I.................................... A-2579 .............................I.................................... A-2594 .............................I.................................... 38 Figure 5 (a) 5500.00 700.00 800.00 100:0 75:25 700.00 50:50 25:75 0:100 OA - 5.238 5000.00 700.00 900.00 600.00 4500.00 6.39 g/ml 600.00 800.00 4000.00 600.00 500.00 700.00 500.00 3500.00 500.00 400.00 600.00 3000.00 400.00 400.00 EU 500.00 EU EU EU EU 2500.00 300.00 0.098 g/ml300.00 400.00 300.00 2000.00 ND ND300.00 ND 200.00 1500.00 200.00 200.00 200.00 1000.00 5.240 100.00 100.00 3.570 4.651 4.637 3.574 5.801 4.654 100.00 3.572 4.656 5.802 4.267 5.795 5.888 500.00 100.00 0.00 0.00 0.00 0.00 0.00 0.00 1.00 2.00 3.00 4.00 5.00 6.00 0.00 1.00 2.00 3.00 4.00 5.00 6.00 0.00 1.00 2.00 3.00 4.00 5.00 6.00 0.00 1.00 2.00 3.00 4.00 5.00 6.00 0.00 1.00 2.00 3.00 4.00 5.00 6.00 Minutes Minutes Minutes Minutes Minutes (b) 3200.00 13000.0 900.00 3000.00 100:0 75:25 OA - 4.908 800.00 50:50 25:75 750.00 0:100 4.908 12000.0 OA - 4. 908 2800.00 800.00 700.00 2600.00 11000.0 14.69 g/ml 3.55 g/ml 700.00 1.01 g/ml 650.00 2400.00 700.00 10000.0 600.00 2200.00 600.00 9000.0 550.00 2000.00 600.00 8000.0 500.00 1800.00 500.00 500.00 450.00 7000.0 1600.00 OA - 4.904 EU EU 0.40 g/ml 400.00 EU EU 400.00 1400.00 EU 6000.0 400.00 350.00 1200.00 5000.0 300.00 300.00 300.00 4000.0 1000.00 250.00 ND 800.00 3000.0 200.00 200.00 200.00 600.00 150.00 3.292 2000.0 3.286 400.00 4.372 3.293 100.00 4.375 100.00 100.00 4.373 4.375 3.427 3.297 5.250 5.367 4.376 1000.0 200.00 50.00 0.0 0.00 0.00 0.00 0.00 0.00 1.00 2.00 3.00 4.00 5.00 6.00 0.00 1.00 2.00 3.00 4.00 5.00 6.00 0.00 1 .00 2.00 3.00 4.00 5.00 6.00 0.00 1.00 2.00 3.00 4.00 5.00 6.00 0.00 1.00 2.00 3.00 4.00 5.00 6.00 Minutes Minutes Minutes Minutes Minutes 39