Abstract
Ethanol production from sugarcane is a key renewable fuel industry in Brazil. Major drivers of this alcoholic fermentation are Saccharomyces cerevisiae strains that originally were contaminants to the system and yet prevail in the industrial process. Here we present newly sequenced genomes (using Illumina short-read and PacBio long-read data) of two monosporic isolates (H3 and H4) of the S. cerevisiae PE-2, a predominant bioethanol strain in Brazil. The assembled genomes of H3 and H4, together with 42 draft genomes of sugarcane-fermenting (fuel ethanol plus cachaça) strains, were compared against those of the reference S288C and diverse S. cerevisiae. All genomes of bioethanol yeasts have amplified SNO2(3)/SNZ2(3) gene clusters for vitamin B1/B6 biosynthesis, and display ubiquitous presence of a particular family of SAM-dependent methyl transferases, rare in S. cerevisiae. Widespread amplifications of quinone oxidoreductases YCR102C/YLR460C/YNL134C, and the structural or punctual variations among aquaporins and components of the iron homeostasis system, likely represent adaptations to industrial fermentation. Interesting is the pervasive presence among the bioethanol/cachaça strains of a five-gene cluster (Region B) that is a known phylogenetic signature of European wine yeasts. Combining genomes of H3, H4, and 195 yeast strains, we comprehensively assessed whole-genome phylogeny of these taxa using an alignment-free approach. The 197-genome phylogeny substantiates that bioethanol yeasts are monophyletic and closely related to the cachaça and wine strains. Our results support the hypothesis that biofuel-producing yeasts in Brazil may have been co-opted from a pool of yeasts that were pre-adapted to alcoholic fermentation of sugarcane for the distillation of cachaça spirit, which historically is a much older industry than the large-scale fuel ethanol production.
Introduction
Fermentation of sugarcane extracts to distillate the “cachaça” spirit has been intertwined with the Brazilian culture since the sixteenth century. In contrast, the production in Brazil of the biofuel ethanol from sugarcane juices and molasses only reached economic prominence in the 1970s, when governmental stimulus propelled the establishment of a nationwide biofuel industry (). As a result, Brazil is currently the world’s second largest producer of ethanol, with the production of ∼35.5 billion liters in the 2019/20 period, only behind the United States of America (Brazilian Sugarcane Industry Association)1. Ethyl alcohol is derived from sugarcane sucrose via microbiological fermentation, primarily by S. cerevisiae strains. In Brazil, industrial alcoholic fermentation of sugarcane operates in a continuous manner (; ; Lopes et al., 2015). When a fermentation batch is over, yeast cells are centrifuged, washed, treated with sulfuric acid (to hamper bacterial contamination), and then reinoculated for a new round of fermentation (i.e., two per day across a season of up to 9 months). During the early period of 1970s–1990s, a common practice in the distilleries was to kick-off a fermentation season with starter cultures consisting of commercial S. cerevisiae strains (mainly baker’s yeasts) (; Lopes et al., 2015). However, it was soon observed that the inoculated yeasts were quickly displaced by S. cerevisiae strains that intruded the fermentation vats. Those contaminants were often referred to as “indigenous” or “wild-type” yeasts, although no phylogenetic evidence actually corroborated that notion (; ; ; ; ; Lopes et al., 2015). Regardless of their ambiguous provenance, “intruding” strains such as S. cerevisiae SA-1, BG-1, CAT-1, VR-1, and PE-2, were systematically isolated, characterized, and reintroduced into distilleries, becoming the standard bioethanol selected strains still in use today (; Lopes et al., 2015). Physiological experiments conducted on these strains revealed excellent fermentation profiles and exceptional capacity to withstand the harsh conditions of industrial fermentation that include high salt and sugar content [18–22% (w/w) total reducing sugars], high temperature (usually reaching 37°C), presence of bacterial contaminations (e.g., by Lactobacilli), high ethanol content (up to 12% v/v), and low pH (1.8–2.5) due to sulfuric acid treatment (; ; Walker and Basso, 2020). Overall, sugarcane fermentation represents an active system where contaminants (yeasts and bacteria) from the exogenous environment dynamically compete with inoculated strains. Even today, when traditional bioethanol strains (e.g., PE-2, BG-1, and CAT-1) are favored as starter cultures, they are often displaced over the year by exogenous S. cerevisiae (; Lopes et al., 2015). This begs the question about the phylogenetic provenance of strains that enter and prevail through the fermentation season.
Based on comparative genomic hybridization, Roche 454 and emerging Illumina technologies, the genome sequencing and structural analyses of the S. cerevisiae PE-2 and CAT-1 strains were reported about 10 years ago (; Stambuk et al., 2009; ). Those studies represent a milestone in our understanding of Brazilian bioethanol strains and substantiated the notion that these yeasts are often heterothallic diploids containing a highly polymorphic set of chromosomes. In particular, subtelomeric regions represent a fluidic part of the genome where frequent ectopic recombination events dynamically alter copy number of genes that encode industrial relevant traits, such as stress response and nutrient acquisition. This is the case for the SNO2(3)/SNZ2(3) gene cluster involved in pyridoxal 5’-phosphate (PLP, vitamin B6) and thiamin diphosphate (ThDP, vitamin B1) biosynthesis, which is represented by four copies in the S. cerevisiae PE-2 haploid derivative JAY291 (). Amplifications of the SNO2(3)/SNZ2(3) cluster, further observed in the CAT-1, VR-1, BG-1, and SA-1 strains, are thought to confer selective advantages to the yeasts when grown in conditions of thiamin repression under high sugar concentrations, possibly found in the industrial sugarcane fermentations (Stambuk et al., 2009).
However, the SNO2(3)/SNZ2(3) cluster and other gene families located at the repeat-rich subtelomeric regions cannot be adequately placed in genome assemblies using only short-read sequence data. To overcome this problem, long-read sequence data (e.g., generated using PacBio/Oxford Nanopore technologies) can better resolve these repetitive regions and have been used to improve the chromosomal positioning of telomeres, solo-long terminal repeats (solo LTRs), complete Ty transposable elements, and more importantly, of multi-paralog gene families (McIlwain et al., 2016; Istace et al., 2017; Salazar et al., 2017; Yue et al., 2017). High-throughput sequencing technologies are also facilitating multi-taxa population genomics studies that are clarifying the evolutionary relationships among key industrial and environmental groups of the S. cerevisiae sensu stricto complex (; ; Legras et al., 2018; Peter et al., 2018; Pontes et al., 2020). For example, some studies have hinted at the close evolutionary relationship of Brazilian bioethanol and cachaça strains to European wine yeasts, supporting that two sequential domestication events occurred in the evolutionary lineage of cachaça yeasts; i.e., the primary domestication of wild Mediterranean oak S. cerevisiae originating the European wine yeasts (), and a recent secondary adaptation in Brazil of wine strains for fermenting sugarcane extracts to produce the cachaça spirit (). To gain more insights into the evolutionary history of Brazilian sugarcane-fermenting yeasts we sequenced and assembled two haploid genomes derived from the diploid isolate S. cerevisiae PE-2. Incorporating 42 other Brazilian biofuel/cachaça strains in a comparative analysis, we systematically assessed gene content, polymorphisms, and copy-number variation (CNV) that represent features conserved among bioethanol yeasts. For the first time using a scalable alignment-free phylogenetic approach, we evaluated the evolutionary relationship among diverse industrial S. cerevisiae based on whole-genome sequences from 197 taxa. Our results provide unambiguous support to the notion that Brazilian bioethanol yeasts comprise a uniform group of domesticated and phylogenetically related yeasts.
Materials and Methods
Genome Sequencing of H3 and H4
A frozen stock derived from the original isolate Saccharomyces cerevisiae PE-2 was obtained from the collection of Prof. L. C. Basso at the University of São Paulo, Piracicaba, Brazil (). From this source, propagated cultures were sporulated and a single ascus was dissected to isolate four sibling spores. Two of them, denominated haploid H3 (MATα) and haploid H4 (MATa), were propagated and targeted for DNA extraction (DNeasy Plant Mini Kit, QIAGEN), DNA fragmentation (NEBNext dsDNA Fragmentase, New England BioLabs), and genomic library construction (Illumina TruSeq DNA PCR-Free Low Throughput Library Prep Kit, Illumina). The generated paired-end libraries were quantified by qPCR (KAPA Library Quantification Kit, Roche) and then sequenced at the Federal University of Rio Grande do Sul, Brazil, using the Illumina MiSeq platform (2 × 300 paired-end; MiSeq Reagent Kit v3 supporting 600-cycles, Illumina). The resulting sequence reads were filtered, retaining only those with Phred quality scores ≥ 30 and length ≥ 75 bases. This resulted in a final set of 11,535,554 paired reads from H3, and 12,366,234 paired reads from H4 for subsequent genome assembly (below). In parallel, genomic DNA samples (isolated from H3 and H4 as described above) were sent for PacBio long-read sequencing using the RS II platform at DNA Link (Seoul, Korea), for one SMRT cell each for H3 and H4. The resulting long-reads were processed using Canu v1.7 (Koren et al., 2017), resulting in 83,983 and 59,078 corrected reads from H3 and H4, respectively.
Assembly and Manual Curation of H3 and H4 Genomes
We assembled the H3 and H4 genomes independently using CLC Genomics Workbench version 8.0.1 (Qiagen)2, incorporating the Microbial Genome Finishing Module3. For each haploid sample, a four-step progressive assembly approach was adopted which included in-depth visual inspection and manual curation of each phase of the process. First, Illumina reads were assembled using CLC de novo assembler (a de Bruijn graph algorithm), resulting in 384 and 416 contigs for H3 and H4, respectively. Second, the assembled contigs were aligned and super-scaffolded against the chromosomal-scale genome sequences of the model reference S. cerevisiae S288C using BLASTN (), as implemented in the “contig aligner tool” of the CLC Genome Finishing Module. Overlapping and redundant contigs were merged and putative chromosomal structures were identified by concatenating contigs following the reference genome. A string of “Ns” was provisionally assigned to represent gaps of unknown length separating adjacent contigs. Third, a manual inspection of sequential contigs that overlaid on S288C genome sequences was conducted, chromosome-by-chromosome. When broken connections between contigs relative to S288C are represented by small (<500 bp) repetitive sequences (e.g., solo LTRs or regions of genes with highly similar paralogs), we mapped the H3 or H4 reads to the S288C genome to obtain a consensus sequence that bridge those specific gaps [i.e., a reference-guided assembly (; Lischer and Shimizu, 2017)]. After joining neighboring contigs, reads were realigned to the scaffold to verify that the number of reads mapped over the proposed connection was uniform (i.e., uniform read coverage, reads fully aligned, and unbroken read pairs). The CLC burrows-wheeler aligner was run at a stringent setting (minimal 0.9 length fraction and 0.9 similarity fraction) throughout the assembly process to ensure high accuracy of the alignments. Fourth, PacBio long-read data were used to resolve larger and/or more complex gapped regions that correspond to complete Ty elements (∼6,000 bp), genes with internal repetitions (e.g., the FLO family), or repeated subtelomeric regions (Koren and Phillippy, 2015). Corrected long reads or unitigs generated using Canu (Koren et al., 2017) were unambiguously assigned using BLASTN to unique (non-repeated) regions flanking both sides of any gapped regions. By realigning corrected reads and/or unitigs (CLC multiple sequence aligner, slow and accurate mode; gap open cost 10.0, gap extension cost 1.0), consensus sequence that bridge each gapped region was identified. In a similar approach, PacBio long reads and Canu unitigs assigned to chromosome extremities were used to extend core chromosome structures into subtelomeric regions and resolve segmental duplications involving paralog-rich gene families. Only one major Canu unitig (521) related to the H4 genome could not be properly assigned to one of three alternative positions related to the tips of chromosomes XIII, XV, and XVI (see Figure 1, bottom). Finally, where necessary, Illumina short reads were used to correct for incorrect bases inherent of PacBio assemblies (Koren and Phillippy, 2015). We tested the resulting H3 and H4 assemblies by thoroughly searching for inconsistencies in Illumina read-mapping profiles over these genomes. Diagnostic tests were implemented by the “analyze contigs tool” of the CLC Genome Finishing Module and included detection of sudden changes in read depth, regions with higher/lower coverage, unaligned read ends, and broken paired-end reads. Those features are usually associated with misassembled regions or the presence of structural polymorphisms (Suzuki et al., 2011; Nagarajan and Pop, 2013; Yang et al., 2019). Only a few problematic regions were retrieved that are related to highly repeated telomeric regions and Ty elements that were assembled relying almost exclusively on PacBio sequence reads of suboptimal quality. In addition, many short (<100 bases) unresolved polymeric structures, which cause errors for both Illumina and PacBio sequencing (Laehnemann et al., 2016), were assigned in our assembly with a string of 20 “N”s. Altogether, quality check reassured that our H3 and H4 genome assemblies generated reliable draft sequences for the 16 yeast chromosomes, mitochondrial genome, 2-micron plasmid, and could assign positions for highly repeated subtelomeric gene families that are key to this study. Genome sequence data for H3 and H4 are available at NCBI under the accession numbers HG994140-HG994157 and LR999872-LR999890, respectively. The BioProject accession number is PRJEB31792.
FIGURE 1
Genome Annotation and Gene-Content Comparison Against S288C Reference
We transferred the annotations en bloc from the reference S288C (GenBank accession number GCA_000146045.2) to the assembled H3 and H4 chromosomes using the “annotate from reference” tool from the CLC Genome Finishing Module. In parallel, we generated a genome-wide prediction of all possible coding regions (ORFs) larger than 150 bp. The CLC Graphical interface was used to visually inspect and manually reconcile the two types of information. In addition, ORF predictions allowed us to annotate (and catalog) genes or regions of H3 and H4 genomes that are not represented in S288C (below). The precise structure of each newly predicted feature was corroborated by comparisons against the S288C syntenic counterparts, specifically to correct any discrepancies in gene length, and identify possible cases of premature stop codons in H3 and H4 genomes. Finally, using custom in-house scripts, protein translations and index tags were assigned to each gene locus, and genomic features were further inspected, and corrected where necessary for any inconsistencies between gene and predicted protein lengths. In our analyses of gene content, we excluded uncharacterized and dubious ORFs (Lin et al., 2013) of short length (≤150 bp), and genes related to Ty elements and telomeric regions, given their highly repetitive nature (). Annotation of non-coding RNAs were transferred from the S288C genome and tRNA structures were predicted by using tRNAscan-SE 2.0 (Lowe and Chan, 2016).
Most cases of non-S288C genes were already discovered during our annotation process within regions in which no S288C annotations were transferred, or in cases where predicted ORFs found no counterpart in S288C. In addition, the H3 and H4 genomes were inspected by TBLASTN for the presence of previously described non-S288C genes (; McIlwain et al., 2016). Regions potentially introgressed from Saccharomyces paradoxus were identified by mapping H3 and H4 reads against the genome of either S. paradoxus strain [CBS432 (ASM207905v1) or UFRJ50816 (ASM207914v1)], and simultaneously to S288C. H3/H4 reads that mapped against S. paradoxus genome revealed chromosomal segments that are putatively introgressed into the genome S. cerevisiae PE-2. Those regions were further queried using BLASTN searches against the nucleotide collection database of GenBank to fulfill the criteria of returning top hits to S. paradoxus instead of any S. cerevisiae taxa. OrthoFinder v2.3.7 () was used to infer homologous protein families among H3, H4 and S288C, from which orthologs of H3/H4 (i.e., PE-2) not recovered among the S288C predicted proteins, or vice versa, were identified. Genes (and genomic regions) that are absent in PE-2 relative to S288C were alternatively identified by mapping of H3 and H4 Illumina reads against the S288C genome (“coverage analysis,” CLC Genome Finishing Module). Most cases of gene-copy number variation were identified during the earlier assembly of difficult-to-resolve regions, usually related to subtelomeres. For this reason, we used all subtelomeric regions as a query to search, using BLASTN, for further cases of segmental duplications and catalog CNVs between H3, H4 and S288C. In a complementary approach, OrthoFinder () was used to identify cases where protein families were expanded or contracted. CNVs were further probed based on analysis of read coverage (below).
Comparative Genomic Analyses Across S. cerevisiae Strains
For SNP comparison, we focus on H3, H4, JAY291 (another PE-2 haploid derivative), and S288C. Genomes of S. cerevisiae S288C (GCA_000146045.2) and JAY291 (ASM18231v2) were downloaded4. The first step in this analysis was to mask in the four genomes all repetitive regions such as telomers, transposable elements, rRNA and highly duplicated subtelomeric gene families. Second, unmasked regions corresponding to core chromosomes were aligned using nucmer from the MUMmer 3.0 package (Kurtz et al., 2004). Finally, tools from the same package were used for SNP count and characterization of aligned regions. GFF files were used to distinguish genic and intergenic regions in the analyzed genomes.
For more-comprehensive comparative genomic analysis, we focused on 44 bioethanol and cachaça strains (hereinafter referred to as the Brazilian bioethanol group; Supplementary Table 1). They correspond to isolates from ethanol plants, mostly located in the state of São Paulo, Brazil, and include the haploid genomes of H3, H4, JAY291 (), the diploid genome of BG-1 (), five strains sequenced by , plus 35 strains (mostly diploid genomes) from the 1,002 Yeast Genomes Project (Peter et al., 2018). For enrichment analysis of distinct genomic features (including genes and alleles) within those strains, a query specific for each feature was used to search by BLASTN against the genomes of the Brazilian bioethanol group (Supplementary Table 1) plus a dataset of 976 strains derived from the 1,002 Yeast Genomes Project (Peter et al., 2018). Query sequences, specific selection criteria, and BLAST settings for this analysis are described in Supplementary Table 2. A genetic feature was considered as present even in heterozygosis (i.e., both alternative alleles were present). The statistical significance of the results was assessed using Fisher’s exact test (implemented by GraphPad Prism 7.00, GraphPad Software, La Jolla, California, United States) comparing bioethanol to non-bioethanol strains and selecting for p < 0.05 as cut-off.
For analysis of CNVs, we extracted genes/regions of S288C for which copy-numbers were found to be amplified or contracted in H3 and H4, or in genomes of other bioethanol strains (Supplementary Table 3). Sequences from H3 and H4 were selected only in cases where no S288C counterpart was available. A nucleotide sequence (Supplementary Sequence File 1) was assembled by concatenating these genes/regions. Flanking regions of at least 250 bp (usually 500–1,000 bp), corresponding to promoter (upstream) or terminator (downstream) parts of genes, were added to separate the regions and to ensure full alignments of reads at the 5’- or 3’-end of the genes. Illumina reads from 40 Brazilian bioethanol strains (except H3 and H4) were downloaded from The European Nucleotide Archive (ENA)5, filtered by quality, and mapped against the S288C chromosomes to obtain the average genome coverage for each strain (Supplementary Table 4). For comparison, we included 10 non-bioethanol strains (Supplementary Table 5) randomly picked from the dataset of the 1002 Yeast Genomes Project (Peter et al., 2018). Read mappings were conducted using a burrows-wheeler aligner (BWA) (Li and Durbin, 2009). After obtaining strain-specific whole genome coverage, reads from each strain were mapped against the concatenated dataset (Supplementary Sequence File 1) to interrogate for CNVs over the selected genes/regions. CNVs were calculated by normalizing read coverage of the coding sequence to its length, and to the background read coverage of the genome. As an example, Supplementary Table 3 describes the numbers and calculations of CNVs for H3. For H3 and H4, we retrieved very similar CNV values estimated with our read depth approach to the numbers derived from our genome assembly [Supplementary Figure 1; For H3: Pearson’s r = 0.9976, Confidence Interval 95% (0.9956, 0.9987), p ≤ 0.0001; For H4: Pearson’s r = 0.9978, Confidence Interval 95% (0.9959–0.9988), p ≤ 0.0001]. Absolute copy numbers for probed genes/regions were subtracted from the numbers derived from the S288C genome annotation (GCA_000146045.2) to compose a matrix of differential copy numbers that was transformed in a heatmap representation using the GraphPad Prism 7.00 package (GraphPad Software, La Jolla, California, United States). A Mann-Whitney U test (implemented by GraphPad Prism 7.00) was applied to infer statistical significance of CNVs found in strains of the bioethanol group when compared to the non-bioethanol control. As a complementary approach, we included an analysis where key genes used in the previous read depth CNV analysis (Supplementary Table 3) were taken as queries in BLASTN searches to count gene copy numbers in 27 genomes sequenced by Illumina plus PacBio or Oxford Nanopore technologies (Istace et al., 2017; Yue et al., 2017). Those nearly complete genome assemblies are indicated on Supplementary Table 5 and mostly correspond to non-bioethanol strains, providing a further control group to evaluate significance of CNVs in the bioethanol group. BLASTN identity higher than 70% was selected as a cut-off of to imply the presence of the gene/region analyzed.
Alignment-Free Phylogenetic Analysis
We adopted an alignment-free (AF) approach (; ) for inferring phylogenetic relationships from whole-genome sequences, based on similarity of shorter, sub-sequences at length k (i.e., k-mers). For this analysis, we used a total of 197 genomes combining the dataset of bioethanol strains (Supplementary Table 1) and 157 genomes from a study of industrial yeasts (; Supplementary Table 6). Genomes of the Saccharomyces paradoxus strains CBS432 (ASM207905v1) and UFRJ50816 (ASM207914v1) were included as outgroup. The majority of these genomes are in chromosome-level resolution based on scaffolding against S288C, except for genomes of 35 bioethanol strains obtained from the 1002 Yeast Genomes Project, and BG-1. These 36-genomes data, although representing draft, fragmented assemblies that could bias enumeration of k-mers, were included to allow for comprehensive inclusion of available bioethanol strains. To ensure comparability of these genomes against the others in our dataset, we refined the 36 genome assemblies using Chromosomer v0.1.4a (Tamazian et al., 2016) to reconstruct chromosome-level assemblies, with S288C as reference. Where applicable, short (<1 kb) contigs were removed, and repeats were masked using RepeatMasker 4.0.96. Enumeration of k-mers for each set of genome sequences was conducted using Jellyfish (Marcais and Kingsford, 2011). The optimal k for phylogenetic analysis is based on a comprehensive assessment across all genomes over an increment of k (Greenfield and Roehm, 2012), the threshold at which the percentage of unique k-mers for each genome reach a plateau (Supplementary Figure 2). We chose k = 21 based on this assessment, as this value provides a saturated level of unique k-mers in each genome while providing sufficient information of shared k-mers among the genomes for calculation of pairwise distance. Following , for each possible genome pair, we calculated a pairwise distance based on the statistic (Wan et al., 2010). The pairwise distance matrix was used to infer a neighbor-joining tree, using neighbor implemented in PHYLIP v 3.697 (). To assess the support for each internal node of the resulting AF tree we adopted a jackknife approach (). Briefly, for the 197-genome set, we generated 100 pseudo-replicates, in which 40% of the bases for each set of genome sequences were randomly deleted, a 100 bp fragment at a time. A phylogenetic tree was then calculated from each pseudo-replicate (of a reduced genome-set) using the AF method (above). The support for each internal node within the reference tree was calculated based on the recovery of the node among the pseudo-replicate trees.
Results
H3, H4, and JAY291 Represent Alternative Polymorphic Genomes of the Bioethanol Strain S. cerevisiae PE-2
We isolated two alternative spores derived from a S. cerevisiae PE-2 single tetrad. For each haploid, we generated genome assemblies combining Illumina MiSeq short-read and PactBio RS II long-read data, guided by the reference genome of S. cerevisiae S288C at chromosomal resolution (see section “Materials and Methods”). We recovered genome sequences for the corresponding chromosomes of H3 and H4 (Table 1), including most of the subtelomeric repeats. These new genome assemblies of the Brazilian bioethanol strain S. cerevisiae PE-2 represent alternatives to the one from a different haploid isolate (JAY291) that was previously sequenced with Illumina (GAIIx) and Roche 454 technologies ().
TABLE 1
| Assembly information | H3 (MATα) | H4 (MATa) |
| Illumina coverage* | 190.8x | 213.6x |
| PacBio coverage* | 23.7x | 14.7x |
| Final assembly | 16 chromosomes; mitochondrion; 2μ plasmid | 16 chromosomes and 1 PacBio unitig; mitochondrion; 2μ plasmid |
| Nuclear genome | 11,930,681 bp | 11,932,493 bp |
| Mitochondrion | 83,189 bp | 83,189 bp |
| 2μ plasmid | 6,119 bp | 6,119 bp |
| CDSs | 5,732 | 5,748 |
| Non-coding RNAs | 97 | 97 |
| tRNAs | 298 | 298 |
| Ty elements | 36 | 35 |
| SNPs to H3 | N.A. | 0.93 SNPs/kb |
| SNPs to H4 | 0.93 SNPs/kb | N.A. |
| SNPs to JAY291 | 1.67 SNPs/kb | 1.61 SNPs/kb |
| SNPs to S288C | 4.41 SNPs/kb | 4.43 SNPs/kb |
H3 and H4 genome assembly and SNPs numbers.
*Calculated with the exclusion of chromosome XII (repetitive rRNA cluster), the mitochondrial genome, and the 2μ plasmid.
Direct comparison of non-repetitive genome regions between H3 and H4 estimated on average 0.93 SNPs/kb (Table 1 and Supplementary Table 7). Loss-of-heterozygosity (LOH) is observed across many chromosomal regions (with only 0.03 SNPs/kb identified in Chr. VI), consistent with the fact that H3 and H4 originated from sibling spores that were derived via meiosis from the same parental diploid genome. Patterns of LOH across the genome were also previously observed in segregants of JAY270, the parental diploid of JAY291 (Rodrigues-Prause et al., 2018; Sampaio et al., 2019). The higher polymorphism density observed between H3 and H4 along Chrs. IX and XIII (1.95 and 1.88 SNPs/kb, respectively; Supplementary Table 7) is close to an early estimation of about 2.0 SNPs/kb for the S. cerevisiae PE-2 diploid genome (). H3 and H4 are also more polymorphic when compared to S288C (about 4.42 SNPs/kb) and to JAY291 (about 1.64 SNPs/kb) (Supplementary Tables 8–11). Substantial heterozygosity among H3, H4, and JAY291 indicates that these monosporic isolates comprise three distinct genomic sets representing the bioethanol strain S. cerevisiae PE-2.
Besides heterozygosity at nucleotide level, H3 and H4 genomes further display substantial structural polymorphisms, especially at chromosome tips. For example, at the Chr. X left-end of H4, there is a 34.7 kb region (spanning 14 genes) that is absent in S288C and H3 (Figure 1 and Supplementary Figure 3). On the other hand, at the corresponding location on Chr. X, H3 carries within a 32.7 kb segment (Figures 1, 2A) the AQY3 cluster (THI5/AAD6/AGP3/AQY3/DAK2/ZNF1/IMA5/HXT8) and the MPR1/C2U11_6360 genes. The same 10-gene region is located at the Chr. VI left-end of both H3 and H4, and represents the long form of Chr. VI previously noticed in S. cerevisiae PE-2 (). Further structural variations at chromosome termini account for differences between H3 and H4 in copy numbers of some key genes (Figures 1, 2). Some examples are the SNO2(3)/SNZ2(3) cluster (five copies in H3 and seven in H4, Figure 2B), a family of SAM-dependent methyl transferases (SAM-mt, one in H3 and three in H4, Figure 2B), Region B (three in H3 and two in H4, Figure 2C), and the MAL cluster (three in H3 and two in H4, Figure 2D). Overall, genome assemblies of H3 and H4 corroborate the notion that in bioethanol strains (sub)telomeres are fluidic genomic regions, accounting for dynamic CNVs of potentially adaptive genes ().
FIGURE 2
Non-S288C Genetic Features Provide Hallmarks of Bioethanol Strains
In search for genetic features that may be linked to bioethanol production, we first compared H3 and H4 genomes to the reference S288C to identify genes exclusive to the Brazilian strains. We detected 38 and 46 genes in H3 and H4, respectively, that are not present in S288C (Supplementary Table 13). Many of them are embedded within gene clusters (Figures 1, 2 and Supplementary Table 14). A typical case of a non-S288C gene is provided by MRP1 (encoding L-azetidine-2-carboxylic acid acetyltransferases) which is present in H3 and H4 (Figures 1, 2A) and widespread among many yeasts (Figure 3). Another excellent example is the cluster of 14 genes (∼34.7 kb) located in the H4 genome at the Chr. X left-end (Figure 1). This region is rare among S. cerevisiae representatives and has been previously observed in JAY291 (
FIGURE 3

Genetic features enriched among yeasts of the bioethanol group. (A) The presence of key genetic features (represented on the x-axis by non-S288C genes/regions, loss-of-function alleles, and absence of specific genes) was determined in 44 strains (y-axis) by BLASTN searches. The color scheme depicts the presence (dark red) or absence (white) of the feature for each strain. S. cerevisiae isolates are described in the Supplementary Table 1. (B) Enrichment of each genetic feature in 44 strains from the bioethanol group (percentage of representation shown in red bars). For comparison, BLASTN searches were also conducted against a dataset of 976 strains and the fraction of presence for each genetic feature is shown in blue bars. AQY2 and OPT2 mutations combine the polymorphisms shown in panel (A). Asterisks (*) indicate significant enrichment of the feature among yeasts of the bioethanol group (Fisher’s exact test, p < 0.005).
Standing out within non-S288C gene families is the one encoding for SAM-dependent methyltransferases (SAM-mt, Pfam family Methyltransf_11, PF08241) that is found in H3 and H4 in one and three copies, respectively (Figures 1, 2B and Supplementary Tables 13, 14). This gene family is rare among S. cerevisiae and uncharacterized so far. Using BLASTN we assessed if SAM-mts and six other non-S288C genes/regions found in H3 and H4 genomes (green regions on chromosomes of Figure 1) are common among yeasts of the Brazilian bioethanol group (including 38 fuel-producing yeasts plus six cachaça isolates, listed on the Supplementary Table 1). For comparison with the Bioethanol group, we extended our BLASTN searches to interrogate a cluster of 976 genomes comprising 34 non-bioethanol clades from the 1,002 Yeast Genomes Project (see Supplementary Figure 5). As a result, we found that the SAM-mt gene is present in 42 out of 44 yeasts (i.e., 95.45% of taxa) from the Brazilian bioethanol group, including all biofuel strains (it is only absent from two cachaça strains) (Figures 3A,B). Conversely, SAM-mt genes are observed in only 95 out of 976 remaining non-bioethanol yeasts (i.e., 9.73% of the strains in the dataset) (Figure 3B and Supplementary Figure 5). This significant enrichment (p < 0.005, Fisher’s exact test) among the Brazilian strains makes the presence of SAM-mt a hallmark of bioethanol yeasts.
The Region B (Figures 1, 2C) encompasses five non-S288C genes that represent a signature of horizontal gene transfer from Zygosaccharomyces bailii to S. cerevisiae, postulated to have occurred at the basal lineage of the wine yeast clade (
There are other chromosomal regions in H3 and H4 that slightly differ from S288C (below 93% nucleotide identity), and are instead more similar to regions from S. paradoxus. These potentially introgressed regions are considered here as non-S288C features and were also selected for analysis by BLASTN. In H3 and H4, a four-gene region spanning ∼15.1 kb, distally placed on Chr. VI right-end (Figures 1, 2A) carries an ortholog of the S288C IRC7 that encodes a beta-lyase involved in the production of thiols (Roncoroni et al., 2011). The Chr. VI-encoded IRC7 from H3 and H4 has 82.3% nucleotide identity to the S288C counterpart, while carrying 99.8% identity to the S. paradoxus homolog. We found that the S. paradoxus IRC7-type is particularly associated (p < 0.005) to the bioethanol group [present in 34 (77.3%) taxa] when compared to non-bioethanol yeasts [only in 26 (2.66%) strains] (Figures 3A,B). Some bioethanol strains have both the S288C and S. paradoxus IRC7 paralogs. As recently reported (
Possible Fingerprints of Domestication in Bioethanol Strains: Aquaporins and Iron Homeostasis Related Genes
We further compared H3 and H4 genomes against the reference S288C to investigate cases of gene loss or defective alleles associated with bioethanol strains. Based on mapping of short reads, H3 and H4 lack 49 and 35 genes, respectively, that are encoded in the reference genome of S288C (Supplementary Table 15). Most of these genes correspond to subtelomeric gene clusters that are duplicated in S288C, but are completely absent or occur in single copies in H3 and/or H4. Some examples are the S288C regions MPH2/SOR2 (four genes on Chr. IV), VTH2/MPH3 (five genes on Chr. X), BDS1 (five genes on Chr. XV), and ENB1 (six genes on Chr. XV). Read depth analysis of the CNVs across multiple S. cerevisiae strains indicates that those regions are also under-represented among non-bioethanol yeasts (Supplementary Figure 6), suggesting that S288C genes not found at the syntenic position on H3 and H4 chromosomes may have arisen from gene gains or duplications specific to the S288C lineage.
We searched for gene models of H3 and H4 with ORFs shorter than their corresponding syntenic homologs in S288C, i.e., potentially having premature stop codons, as a proxy for loss of function. Using this approach, we identified 44 and 40 potential pseudogenes in the H3 and H4 genomes, respectively (Supplementary Table 16). Interestingly, we observed many instances of premature stop codons in genes within amplified regions [e.g., DDI2(3), SNZ2(3), AGP3, SEO1, and other examples on Figure 2], suggesting that the mutated gene is not the one under selection for amplification. Among the remaining putative defective alleles, some may be implicated in adaptation to industrial conditions. For example, most bioethanol strains have been described as heterothallic, as haploid cells cannot switch their mating-type (
Mutations that inactivate aquaporin genes are adaptive for fermentation of sugar-rich substrates in industrial processes because the loss-of-function alleles confer a significant selective advantage to withstand high osmolarity (
We observed many cases of defective or deleted alleles related to iron homeostasis. Ferric reductases (Fre) act at the cell periphery by reducing soluble or siderophore-bound ferric iron allowing the uptake of the ferrous form through the Fet3/Ftr1 complex (Philpott et al., 2002). We found in H3 and H4 an adenine insertion (frameshift) at position 1502 of the FRE3 gene that is also enriched among 17 (38.6%) strains of the Brazilian bioethanol group when compared to 178 (18.2%) of non-bioethanol yeasts (p < 0.005; Figure 3). The FRE3 gene is located within a cluster of iron regulated genes on the Chr. XV right-end that also harbors the genes FIT2 and FIT3 (encoding for mannoproteins that promote retention of siderophore-iron in the cell wall) and the ferric reductase homolog FRE5. We observed two types of deletion encompassing FIT2, FIT3, and FRE5 genes in the genomes of bioethanol strains (Supplementary Figure 7). A deletion of ∼5.36 kb is present in 10 strains and has breaking points surrounding the three genes. Interestingly, a second larger deletion that is found in 22 bioethanol strains is reminiscent of the ∼14 kb deletion described in flor yeasts (
FIGURE 4

Differential copy numbers of genes/regions in strains of the bioethanol group compared to S288C. Copy numbers of key genes/regions (labeled at the bottom) were measured in bioethanol yeasts by read depth based on Illumina sequence data. The estimated CNVs of probed genes/regions were discounted from the expected copy numbers according to the S288C genome annotation. A heatmap is shown depicting a color scheme in which extra copies in strains of the bioethanol group are quantitatively expressed in a red gradient (up to four copies). Equal number of copies found in both bioethanol and S288C strains is represented in white, while surplus of copy numbers in S288C is displayed as a blue gradient. Note that in some cases, hemizygous patterns account for non-integer gene copy numbers, since multiple diploid genomes are represented within the dataset (Peter et al., 2018). A group of ten randomly chosen non-bioethanol yeasts were included as external controls. Dark red colored boxes outside the red gradient range represent cases in which more than four extra copies of the probed gene/region were counted. Significant enrichment of CNV among the bioethanol clade is indicated with an asterisk (∗) above the respective column (p < 0.05, Mann-Whitney U test).
OPT2 (encoding the oligopeptide transporter 2) has been demonstrated to have genetic interactions with components of the iron homeostasis, such as FET3, SIT1, GEF1, and FRE3 (
Distinct Patterns of Gene CNVs Unite Bioethanol Yeasts
Many cases of gene CNVs in H3 and/or H4 genomes are evident relative to the S288C reference (Figures 1, 2 and Supplementary Tables 14, 15, 17). We further assessed if these CNVs are also observed in other bioethanol yeasts (Supplementary Table 4). Genome data (Illumina short reads) from 42 bioethanol strains (including H3 and H4) were independently mapped against a concatenated DNA sequence in which genes/regions probed for CNVs are represented (Supplementary Sequence File 1, see section “Materials and Methods”). For each strain, read depth within the analyzed sequence regions was normalized against the read depth over the genome as background (Supplementary Table 4), allowing estimation of differential CNVs relative to the S288C genome (Figure 4 and Supplementary Figure 6). CNVs for 10 non-bioethanol strains were also analyzed as external controls. Our analysis revealed that several cases of gene/region CNVs are enriched in the Brazilian bioethanol group (p < 0.05, Mann-Whitney U test; Figure 4). Remarkably, all bioethanol yeasts have increased copy numbers within the gene cluster for vitamin B1/B6 biosynthesis; the core genes of this cluster comprise SNO2(3) and SNZ2(3), and are usually flanked by DDI2(3), THI5(12), and AAD4. Based on molecular karyotyping and southern blot hybridization,
Besides the SNO2(3)/SNZ2(3) amplification, the occurrence of other gene CNVs in bioethanol yeasts has been described in previous comparative genomic hybridization analyses (
Among newly found cases of CNVs in the H3 and H4 genomes, the CUP1 gene (encoding for a metallothionein) stands out by occurring in tandem repeats with ∼18 estimated copies on Chr. VIII of both H3 and H4 (Supplementary Table 17 and Supplementary Figures 9, 10). CUP1 is known as a hotspot for adaptive and non-adaptive amplifications on S. cerevisiae genomes (Zhao et al., 2014), thus our observation of variable copy numbers of CUP1 among bioethanol and non-bioethanol yeasts (Supplementary Figure 10) is not surprising. Our gene CNVs analysis also confirmed that some of the non-S288C genes, such as the Region B, SAM-mt, and CSS3, often display more than one copy in strains from the Brazilian bioethanol group (Figure 4). For example, the AEG strain has four estimated copies of SAM-mt genes, while H4 has three paralogs. Likewise, Region B exists in three copies in H3, two located in the middle of Chrs. VI and XI, respectively, and one at the right tip of Chr. XI. H4 has a copy of Region B syntenic to H3 on the core of Chr. VI, and another at the subtelomeric part of Chr. VII right-end (Figures 1, 2C). These subtelomeric copies of Region B at the Chr. VII of H4 and Chr. XI of H3 are in fact part of a large (∼48.9 kb) cluster that also includes the regions MAL, ADH7(8), and SNO2(3)/SNZ2(3) (Figures 1, 2D). In H4, the ADH7(8) region (PAU3/ADH7/RDS1/AAD3) is also duplicated at the Chr. III right-end, while H3 displays at this corresponding position an extra copy of ECM34 (representing a fusion of the S288C ECM34 and YHL042W homologs) (Figures 1, 2G). We found that amplifications of ECM34 and the ADH7 region appear to be recurrent among bioethanol strains (Figure 4). Similarly, the FDH1 region (comprising FDH1/YOR389W/FEX1, Figures 1, 2F) occurs in three copies in the genomes of H3, H4, and other bioethanol strains, and only two in the S288C genome.
Previous comparative genomic hybridization analysis of the strain CAT-1 (
Of special interest is the amplification in H3 and H4 genomes of gene homologs annotated as quinone oxidoreductases (Figures 1, 2G). They are represented in S288C as YCR102C (Chr. III), YLR460C (Chr. XII) and YNL134C (Chr. XIV), where they are observed as monoallelic. In contrast to S288C, H3 and H4 have three homologs related to YCR102C on subtelomeric regions of Chrs. III, V, and XII, and a further paralog more similar to YLR460C at the Chr. IX left-end. In addition, H3 and H4 display a tandem duplication of YNL134C at the core of Chr. X (Figures 1, 2G), also previously observed in the strain CAT-1 (
In a complementary analysis, key gene/regions with differential CNVs in bioethanol strains were analyzed by BLASTN in 27 strains (2 bioethanol and 25 non-bioethanol) for which chromosomal-level (i.e., near finished) assemblies are available (Istace et al., 2017; Yue et al., 2017), allowing counting of gene copy numbers. The inclusion of H3 and H4 genomes enabled direct comparison of CNVs in four bioethanol strains to those found in the group of 25 non-bioethanol yeasts (Supplementary Figure 13). The analysis confirmed the significance for the presence and higher copy numbers in bioethanol strains of SAM-mt, Region B, the DDI2(3)/SNO2(3)/SNZ2(3) cluster, quinone oxidoreductases, and ECM34. In addition, the comparison demonstrated the absence in the 25 control strains of amplifications in GAP1, FCY1, NRT1, and HSP30 region, found in bioethanol yeasts by read depth analysis. Finally, amplifications of the regions ARR and MCH2/FRE2, widespread in bioethanol strains (Figure 4), were scarcely observed among the 25 non-bioethanol yeasts.
Whole-Genome Phylogeny Reveals a Single Origin of Bioethanol Yeasts
To assess the evolutionary origins of bioethanol yeasts in Brazil, we reconstructed a phylogenetic tree using an alignment-free (AF) approach, based on whole-genome sequences from 197 broadly sampled, representative yeast genomes (see Materials and Methods for details). This AF approach, based on analysis of short, sub-sequences of length k (k-mers), represents a faster, more-scalable alternative to standard phylogenetic inference based on multiple sequence alignment. Bypassing the unrealistic assumption of full-length contiguity among homologous sequences, and the computationally demanding step of multiple sequence alignment, AF approach allows for inferring phylogenetic relationships using whole-genome sequences quickly. Such an approach is known to be robust to distinct evolutionary scenarios, e.g., genetic rearrangement (
FIGURE 5

AF-based phylogenetic tree of 197 yeasts genomes. A distance matrix was calculated based on an AF approach and used to build a NJ phylogenetic tree. The Brazilian bioethanol group (red) was clustered in a single clade, with the exemption of the cachaça strain SP005. In addition, five major clades were resolved in accordance to
Discussion
Genetic Signatures of Bioethanol Yeasts and Adaptation to Industrial Fermentation
The present comparative genomics study is one of the first attempts to analyze genome structure, gene content, polymorphisms, gene CNVs, and phylogenetic relationships across diverse S. cerevisiae used for bioethanol production in Brazil. The Brazilian fuel ethanol production from sugarcane substrates represents a peculiar bioprocess in which yeast cells are recycled and reused in successive fermentation batches. Therefore, throughout the year, bioethanol yeasts are steadily challenged by multiple biotic and abiotic stresses (
FIGURE 6

Twenty genetic hallmarks of Brazilian bioethanol yeasts. Genetic features enriched (p < 0.05) in bioethanol yeasts are shown in red, including: (i) non-S288C genes, (ii) putative introgressions from S. paradoxus, (iii) large deletions, (iv) SNPs and small InDels, (v) amplified regions.
In Brazil, ethanol generation from sugarcane is traditionally coupled with sugar production (
Nutrient availability seems to be a key factor for industrial fermentation of sugarcane substrates. In most S. cerevisiae strains there are three copies for both SNO (encoding glutaminase) and SNZ (PLP synthase). These two enzymes, respectively, control successive catalytic steps for the synthesis of vitamin B6 (PLP) (Rodriguez-Navarro et al., 2002). SNO2(3)/SNZ2(3) are located in subtelomeric clusters together with the thiamin biosynthetic gene THI5(12), and their products physically interact in a complex that likely channels the PLP for thiamin production (Rodriguez-Navarro et al., 2002; Paxhia and Downs, 2019). In accordance with previous observations (
Mineral nutrition and/or toxicity may also be important for yeasts during fuel ethanol production, as suggested by amplification of ARR cluster for arsenic detoxification we observed in 20 bioethanol strains. The high iron content in molasses possibly represents an additional critical issue (
Other genetic markers of bioethanol strains could be implicated in varied stress responses. In the Brazilian sugarcane fermentation system temperatures frequently raise above 35°C (
Are Brazilian Bioethanol Yeasts Historically Derived From Cachaça Spirit Strains?
To understand the evolution of the Brazilian bioethanol strains we performed a phylogenetic analysis including taxa from the major groups of industrial S. cerevisiae, encompassing in total 197 yeast genomes. This represents one of the first attempts to use an AF method for inferring phylogenetic relationships among eukaryotic genomes. Our phylogenetic tree, largely topologically congruent with that generated by
It should be noted that our analyses were largely restricted to S. cerevisiae isolates from ethanol plants in the state of São Paulo. Therefore, a broader representation of industrial yeasts from different regions in Brazil is necessary to achieve a more-comprehensive overview of the phylogenetic origin and the common genetic features shared by Brazilian bioethanol strains. This is particularly relevant giving that some cachaça isolates display different phylogenetic affiliations from most others that are more-closely related to wine and bioethanol yeasts (
Considering our new data, we summarize the putative origin of Brazilian bioethanol S. cerevisiae in Figure 7. The present-day bioethanol yeasts were derived from “contaminants” that invaded the fermentation process and outcompeted the starter strains (
FIGURE 7

Model for the evolution of Brazilian bioethanol yeasts. Cachaça yeasts probably derive from an ancestral European yeast lineage domesticated for the wine-making process. Since cachaça and bioethanol yeasts have a mosaic genome structure, it is likely that at some point a hybridization event occurred between an European ancestor and a yet-unknown yeast (perhaps from a South American wild population?). Cachaça yeasts were latter co-opted for biofuel production becoming the modern-day bioethanol yeasts. The degree in which bioethanol yeasts may evolve within the process of fermentation is still unknown. This process-related evolution may account for the genetic divergence between cachaça and fuel ethanol strains.
A more likely scenario for the origin of bioethanol strains is through their European ancestry (Figure 7).
As bioethanol strains are phylogenetically closely related to C1 and C2 cachaça strains (
Finally, the degree to which bioethanol strains are evolving within the fermentation process should also be investigated by using genomic approaches (Figure 7). The industrial fermentation technology in Brazil is based on yeast biomass recycling and reuse in a continuous manner throughout the year, setting up an evolutionary experiment in which new mutations may arise and be selected according to their fitness contribution to withstand harsh industrial conditions. In fact, the frequently observed shifts in yeast karyotype patterns during industrial fermentation of sugarcane substrates indicate dynamic occurrences of chromosomal structural variations (
Statements
Data availability statement
All sequence data generated from this study are available at NCBI Short Read Archive (BioProject PRJEB31792). The two annotated genomes are available at GenBank, with nucleotide accessions HG994140-HG994157 for H3, and LR999872-LR999890 for H4. Accession numbers from other genomes analyzed in this work are available in Supplementary Tables 1, 4–6.
Author contributions
AJ, JG, and CXC: project conceptualization and manuscript writing. LB: strain/spore isolation characterization. AJ, LB, and JF: DNA preparation, genome library preparation, and genome sequencing. AJ and JG: genome assembly. AJ, JG, TS, PY, RG-P, YC, MC-C, KD, and CXC: genome annotation, diverse bioinformatic analyses, and statistical analyses. All authors contributed to the article and approved the submitted version.
Funding
This work was supported by research grants from the São Paulo Research Foundation, FAPESP: 2013/15743-9 and 2017/13972-1 to JG, and 2017/24453-5 to AJ. Travel funding and exchange activities were provided by the FAPESP Mobility Grant 2018/15159-9 to JG and CXC, supported by the São Paulo State University International Office (Brazil) and the University of Queensland Global Engagements (Australia).
Acknowledgments
This work was supported by the computational resources of the Australian National Computational Infrastructure (NCI) National Facility systems through the NCI Merit Allocation Scheme (Project d85) awarded to CXC.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmicb.2021.644089/full#supplementary-material
Footnotes
1.^http://www.unicadata.com.br/
2.^https://digitalinsights.qiagen.com
3.^https://digitalinsights.qiagen.com/plugins/clc-genome-finishing-module/
References
1
AlmeidaP.BarbosaR.ZalarP.ImanishiY.ShimizuK.TurchettiB.et al (2015). A population genomics insight into the Mediterranean origins of wine yeast domestication.Mol. Ecol.245412–5427. 10.1111/mec.13341
2
AltschulS. F.MaddenT. L.SchafferA. A.ZhangJ.ZhangZ.MillerW.et al (1997). Gapped BLAST and PSI-BLAST: a new generation of protein database search programs.Nucleic Acids Res.253389–3402. 10.1093/nar/25.17.3389
3
AntonangeloA. T.AlonsoD. P.RibollaP. E.ColombiD. (2013). Microsatellite marker-based assessment of the biodiversity of native bioethanol yeast strains.Yeast30307–317. 10.1002/yea.2964
4
ArguesoJ. L.CarazzolleM. F.MieczkowskiP. A.DuarteF. M.NettoO. V.MissawaS. K.et al (2009). Genome structure of a Saccharomyces cerevisiae strain widely used in bioethanol production.Genome Res.192258–2270. 10.1101/gr.091777.109
5
BabrzadehF.JaliliR.WangC.ShokrallaS.PierceS.Robinson-MosherA.et al (2012). Whole-genome sequencing of the efficient industrial fuel-ethanol fermentative Saccharomyces cerevisiae strain CAT-1.Mol. Genet. Genomics287485–494. 10.1007/s00438-012-0695-7
6
BarbosaR.AlmeidaP.SafarS. V.SantosR. O.MoraisP. B.Nielly-ThibaultL.et al (2016). Evidence of Natural Hybridization in Brazilian wild lineages of Saccharomyces cerevisiae.Genome Biol. Evol.8317–329. 10.1093/gbe/evv263
7
BarbosaR.PontesA.SantosR. O.MontandonG. G.De Ponzzes-GomesC. M.MoraisP. B.et al (2018). Multiple rounds o artificial selection promote microbe secondary domestics ion-the case of cachaca yeasts.Genome Biol. Evol.101939–1955. 10.1093/gbe/evy132
8
BassoL. C.BassoT. O.RochaS. N. (2011). “Ethanol production in Brazil: the industrial process and its impact on yeast fermentation,” in Biofuel Production – Recent Developments and Prospects, ed.BernardesM. A. S. (London: IntechOpen).
9
BassoL. C.De AmorimH. V.De OliveiraA. J.LopesM. L. (2008). Yeast selection for fuel ethanol production in Brazil.FEMS Yeast Res.81155–1163. 10.1111/j.1567-1364.2008.00428.x
10
BassoT. O.LinoF. S. O. (2019). “Clash of kingdoms: how do bacterial contaminants thrive in and interact with yeasts during ethanol production,” in Fuel Ethanol Production from Sugarcane, edsBassoT. P.BassoL. C. (London: IntechOpen).
11
BassoT. O.WalkerR. S. K.BassoL. C.WalkerG. M. (2019). “The future of bioethanol,” in Ethanol as a green alternative fuel: insight and perspectives. (Renewable Energy: Research, Development and Policies),edsTreichelH.Alves JúniorS. L.FongaroG.MüllerC. (Hauppauge, NY: Nova Science Publishers, Inc) 259–283.
12
BernardG.ChanC. X.ChanY. B.ChuaX. Y.CongY. N.HoganJ. M.et al (2019). Alignment-free inference of hierarchical and reticulate phylogenomic relationships.Brief. Bioinform.20426–435. 10.1093/bib/bbx067
13
BernardG.ChanC. X.RaganM. A. (2016). Alignment-free microbial phylogenomics under scenarios of sequence divergence, genome rearrangement and lateral genetic transfer.Sci. Rep.6:28970.
14
BernardG.GreenfieldP.RaganM. A.ChanC. X. (2018). k-mer similarity, networks of microbial genomes, and taxonomic rank.mSystems3:e00257-18. 10.1128/mSystems.00257-18
15
BonhiversM.CarbreyJ. M.GouldS. J.AgreP. (1998). Aquaporins in Saccharomyces. Genetic and functional distinctions between laboratory and wild-type strains.J. Biol. Chem.27327565–27572.
16
BornemanA. R.DesanyB. A.RichesD.AffourtitJ. P.ForganA. H.PretoriusI. S.et al (2011). Whole-genome comparison reveals novel genetic elements that characterize the genome of industrial strains of Saccharomyces cerevisiae.PLoS Genet.7:e1001287. 10.1371/journal.pgen.1001287
17
BrionC.AmbrosetC.DelobelP.SanchezI.BlondinB. (2014). Deciphering regulatory variation of THI genes in alcoholic fermentation indicate an impact of Thi3p on PDC1 expression.BMC Genomics15:1085. 10.1186/1471-2164-15-1085
18
CaspetaL.ChenY.GhiaciP.FeiziA.BuskovS.HallstromB. M.et al (2014). Biofuels. Altered sterol composition renders yeast thermotolerant.Science34675–78. 10.1126/science.1258137
19
ChanC. X.BernardG.PoirionO.HoganJ. M.RaganM. A. (2014). Inferring phylogenies of evolving sequences without multiple sequence alignment.Sci. Rep.4:6504.
20
ChenH.LiJ.WanC.FangQ.BaiF.ZhaoX. (2019). Improvement of inhibitor tolerance in Saccharomyces cerevisiae by overexpression of the quinone oxidoreductase family gene YCR102C.FEMS Yeast Res.19:foz055.
21
CostanzoM.VandersluisB.KochE. N.BaryshnikovaA.PonsC.TanG.et al (2016). A global genetic interaction network maps a wiring diagram of cellular function.Science353:aaf1420.
22
CoutouneN.MulatoA. T.Riano-PachonD. M.OliveiraJ. V. (2017). Draft genome sequence of Saccharomyces cerevisiae barra grande (BG-1), a Brazilian industrial bioethanol-producing strain.Genome Announc.5:e00111-17. 10.1128/genomeA.00111-17
23
da Silva-FilhoE. A.Brito Dos SantosS. K.Resende AdoM.De MoraisJ. O.De MoraisM. A.Jr.Ardaillon SimoesD. (2005). Yeast population dynamics of industrial fuel-ethanol fermentation process assessed by PCR-fingerprinting.Antonie Van Leeuwenhoek8813–23. 10.1007/s10482-005-7283-3
24
D’AngioloM.De ChiaraM.YueJ. X.IrizarA.StenbergS.PerssonK.et al (2020). A yeast living ancestor reveals the origin of genomic introgressions.Nature587420–425. 10.1038/s41586-020-2889-1
25
Della-BiancaB. E.BassoT. O.StambukB. U.BassoL. C.GombertA. K. (2013). What do we know about the yeast strains from the Brazilian fuel ethanol industry?Appl. Microbiol. Biotechnol.97979–991. 10.1007/s00253-012-4631-x
26
Elbaz-AlonY.MorganB.ClancyA.AmoakoT. N.ZalckvarE.DickT. P.et al (2014). The yeast oligopeptide transporter Opt2 is localized to peroxisomes and affects glutathione redox homeostasis.FEMS Yeast Res.141055–1067.
27
EldarovM. A.AvdaninaD. A.ShalamitskiiM. Y.IvanovaE. V.TanashchukT. N.KishkovskayaS. A.et al (2019). Polymorphism of the iron homeostasis genes and iron sensitivity in Saccharomyces cerevisiae Flor and Wine Strains.Microbiology88200–205. 10.1134/s0026261719020036
28
EldarovM. A.BeletskyA. V.TanashchukT. N.KishkovskayaS. A.RavinN. V.MardanovA. V. (2018). Whole-genome analysis of three yeast strains used for production of sherry-like wines revealed genetic traits specific to Flor Yeasts.Front. Microbiol.9:965. 10.3389/fmicb.2018.00965
29
EmmsD. M.KellyS. (2019). OrthoFinder: phylogenetic orthology inference for comparative genomics.Genome Biol.20:238.
30
EngelS. R.DietrichF. S.FiskD. G.BinkleyG.BalakrishnanR.CostanzoM. C.et al (2014). The reference genome sequence of Saccharomyces cerevisiae: then and now.G3 (Bethesda)4389–398.
31
FelsensteinJ. (1989). PHYLIP – phylogeny inference package (Version 3.2).Cladistics5164–166.
32
GaleoteV.BigeyF.BeyneE.NovoM.LegrasJ. L.CasaregolaS.et al (2011). Amplification of a Zygosaccharomyces bailii DNA segment in wine yeast genomes by extrachromosomal circular DNA formation.PLoS One6:e17872. 10.1371/journal.pone.0017872
33
GalloneB.SteenselsJ.PrahlT.SoriagaL.SaelsV.Herrera-MalaverB.et al (2016). Domestication and Divergence of Saccharomyces cerevisiae Beer Yeasts.Cell1661397–1410.e16.
34
GnerreS.LanderE. S.Lindblad-TohK.JaffeD. B. (2009). Assisted assembly: how to improve a de novo genome assembly by using related species.Genome Biol.10:R88.
35
GoncalvesM.PontesA.AlmeidaP.BarbosaR.SerraM.LibkindD.et al (2016). Distinct domestication trajectories in top-fermenting beer yeasts and wine yeasts.Curr. Biol.262750–2761. 10.1016/j.cub.2016.08.040
36
GreenfieldP.RoehmU. (2012). Answering biological questions by querying k-mer databases.Concurr. Comput.25497–509. 10.1002/cpe.2938
37
HeerD.HeineD.SauerU. (2009). Resistance of Saccharomyces cerevisiae to high concentrations of furfural is based on NADPH-dependent reduction by at least two oxireductases.Appl. Environ. Microbiol.757631–7638. 10.1128/aem.01649-09
38
IstaceB.FriedrichA.D’agataL.FayeS.PayenE.BelucheO.et al (2017). de novo assembly and population genomic survey of natural yeast isolates with the Oxford Nanopore MinION sequencer.Gigascience61–13.
39
KimS. K.JinY. S.ChoiI. G.ParkY. C.SeoJ. H. (2015). Enhanced tolerance of Saccharomyces cerevisiae to multiple lignocellulose-derived inhibitors through modulation of spermidine contents.Metab. Eng.2946–55. 10.1016/j.ymben.2015.02.004
40
KorenS.PhillippyA. M. (2015). One chromosome, one contig: complete microbial genomes from long-read sequencing and assembly.Curr. Opin. Microbiol.23110–120. 10.1016/j.mib.2014.11.014
41
KorenS.WalenzB. P.BerlinK.MillerJ. R.BergmanN. H.PhillippyA. M. (2017). Canu: scalable and accurate long-read assembly via adaptive k-mer weighting and repeat separation.Genome Res.27722–736. 10.1101/gr.215087.116
42
KurtzS.PhillippyA.DelcherA. L.SmootM.ShumwayM.AntonescuC.et al (2004). Versatile and open software for comparing large genomes.Genome Biol.5:R12.
43
LaehnemannD.BorkhardtA.MchardyA. C. (2016). Denoising DNA deep sequencing data-high-throughput sequencing errors and their correction.Brief. Bioinform.17154–179. 10.1093/bib/bbv029
44
LegrasJ. L.GaleoteV.BigeyF.CamarasaC.MarsitS.NideletT.et al (2018). Adaptation of S. cerevisiae to fermented food environments reveals remarkable genome plasticity and the footprints of domestication.Mol. Biol. Evol.351712–1727. 10.1093/molbev/msy066
45
LiH.DurbinR. (2009). Fast and accurate short read alignment with Burrows-Wheeler transform.Bioinformatics251754–1760. 10.1093/bioinformatics/btp324
46
LiJ.BissM.FuY.XuX.MooreS. A.XiaoW. (2015). Two duplicated genes DDI2 and DDI3 in budding yeast encode a cyanamide hydratase and are induced by cyanamide.J. Biol. Chem.29012664–12675. 10.1074/jbc.m115.645408
47
LinD.YinX.WangX.ZhouP.GuoF. B. (2013). Re-annotation of protein-coding genes in the genome of saccharomyces cerevisiae based on support vector machines.PLoS One8:e64477. 10.1371/journal.pone.0064477
48
LischerH. E. L.ShimizuK. K. (2017). Reference-guided de novo assembly approach improves genome reconstruction for related species.BMC Bioinform.18:474. 10.1186/s12859-017-1911-6
49
LopesM. L.PaulilloS. C. L.CherubinR. A.GodoyA.Amorim NetoH. B.AmorimH. V. (2015). Tailored Yeast Strains for Ethanol Production: Process-Driven Selection. Available online at: https://www.fermentec.com.br/downloads/eBook_-_Tailored_yeast_strains_selected_for_ethanol_production_-_Fermentec.pdf(accessed October 12, 2020).
50
LoweT. M.ChanP. P. (2016). tRNAscan-SE On-line: integrating search and context for analysis of transfer RNA genes.Nucleic Acids Res.44W54–W57.
51
MarcaisG.KingsfordC. (2011). A fast, lock-free approach for efficient parallel counting of occurrences of k-mers.Bioinformatics27764–770. 10.1093/bioinformatics/btr011
52
McIlwainS. J.PerisD.SardiM.MoskvinO. V.ZhanF. J.MyersK. S.et al (2016). Genome sequence and analysis of a stress-tolerant, wild-derived strain of Saccharomyces cerevisiae used in biofuels research.G361757–1766. 10.1534/g3.116.029389
53
MeenaR. C.ThakurS.ChakrabartiA. (2011). Regulation of Saccharomyces cerevisiae Plasma membrane H(+)-ATPase (Pma1) by Dextrose and Hsp30 during exposure to thermal stress.Indian J. Microbiol.51153–158. 10.1007/s12088-011-0137-y
54
NagamatsuS. T.CoutouneN.JoseJ.FiamenghiM. B.PereiraG. A. G.OliveiraJ. V. C.et al (2021). Ethanol production process driving changes on industrial strains.FEMS Yeast Res.21:foaa071.
55
NagarajanN.PopM. (2013). Sequence assembly demystified.Nat. Rev. Genet.14157–167. 10.1038/nrg3367
56
NavarreC.GhislainM.LetermeS.FerroudC.DufourJ. P.GoffeauA. (1992). Purification and complete sequence of a small proteolipid associated with the plasma membrane H(+)-ATPase of Saccharomyces cerevisiae.J. Biol. Chem.2676425–6428. 10.1016/s0021-9258(18)42713-5
57
PatelD. (2013). Characterization of Vacuole Aquaporin Function and Its Implication in Membrane Fission. Available online at: https://spectrum.library.concordia.ca/977797/1/Patel_MSc_F2013.pdf(accessed October 20, 2020).
58
PaxhiaM. D.DownsD. M. (2019). SNZ3 Encodes a PLP synthase involved in thiamine synthesis in Saccharomyces cerevisiae.G3 (Bethesda)9335–344.
59
PeterJ.De ChiaraM.FriedrichA.YueJ. X.PfliegerD.BergstromA.et al (2018). Genome evolution across 1,011 Saccharomyces cerevisiae isolates.Nature556339–344. 10.1038/s41586-018-0030-5
60
PhilpottC. C.ProtchenkoO.KimY. W.BoretskyY.Shakoury-ElizehM. (2002). The response to iron deprivation in Saccharomyces cerevisiae: expression of siderophore-based systems of iron uptake.Biochem. Soc. Trans.30698–702. 10.1042/bst0300698
61
PontesA.HutzlerM.BritoP. H.SampaioJ. P. (2020). Revisiting the taxonomic synonyms and populations of Saccharomyces cerevisiae-phylogeny, phenotypes, ecology and domestication.Microorganisms8:903. 10.3390/microorganisms8060903
62
Rodrigues-PrauseA.SampaioN. M. V.GurolT. M.AguirreG. M.SedamH. N. C.ChapmanM. J.et al (2018). A case study of genomic instability in an industrial strain of Saccharomyces cerevisiae.G3 (Bethesda)83703–3713.
63
Rodriguez-NavarroS.LlorenteB.Rodriguez-ManzanequeM. T.RamneA.UberG.MarchesanD.et al (2002). Functional analysis of yeast gene families involved in metabolism of vitamins B1 and B6.Yeast191261–1276. 10.1002/yea.916
64
RoncoroniM.SantiagoM.HooksD. O.MoroneyS.HarschM. J.LeeS. A.et al (2011). The yeast IRC7 gene encodes a beta-lyase responsible for production of the varietal thiol 4-mercapto-4-methylpentan-2-one in wine.Food Microbiol.28926–935. 10.1016/j.fm.2011.01.002
65
SalazarA. N.De VriesA. R. G.Van Den BroekM.WijsmanM.CortesP. D.BrickweddeA.et al (2017). Nanopore sequencing enables near-complete de novo assembly of Saccharomyces cerevisiae reference strain CEN.PK113-7D.FEMS Yeast Res.17:fox074.
66
SampaioN. M. V.WatsonR. A.ArguesoJ. L. (2019). Controlled reduction of genomic heterozygosity in an industrial yeast strain reveals wide cryptic phenotypic variation.Front. Genet.10:782. 10.3389/fgene.2019.00782
67
StambukB. U.DunnB.AlvesS. L.Jr.DuvalE. H.SherlockG. (2009). Industrial fuel ethanol yeasts contain adaptive copy number changes in genes involved in vitamin B1 and B6 biosynthesis.Genome Res.192271–2278. 10.1101/gr.094276.109
68
SteenselsJ.SnoekT.MeersmanE.Picca NicolinoM.VoordeckersK.VerstrepenK. J. (2014). Improving industrial yeast strains: exploiting natural and artificial diversity.FEMS Microbiol. Rev.38947–995. 10.1111/1574-6976.12073
69
SuzukiS.YasudaT.ShiraishiY.MiyanoS.NagasakiM. (2011). ClipCrop: a tool for detecting structural variations with single-base resolution using soft-clipping information.BMC Bioinform.12(Suppl. 14):S7. 10.1186/1471-2105-12-S14-S7
70
TamazianG.DobryninP.KrasheninnikovaK.KomissarovA.KoepfliK. P.O’brienS. J. (2016). Chromosomer: a reference-based genome arrangement tool for producing draft chromosome sequences.Gigascience5:38.
71
UemuraT.KashiwagiK.IgarashiK. (2005). Uptake of putrescine and spermidine by Gap1p on the plasma membrane in Saccharomyces cerevisiae.Biochem. Biophys. Res. Commun.3281028–1033. 10.1016/j.bbrc.2005.01.064
72
WalkerG. M.BassoT. O. (2020). Mitigating stress in industrial yeasts.Fungal Biol.124387–397. 10.1016/j.funbio.2019.10.010
73
WanL.ReinertG.SunF.WatermanM. S. (2010). Alignment-free sequence comparison (II): theoretical power of comparison statistics.J. Comput. Biol.171467–1490. 10.1089/cmb.2010.0056
74
WohlbachD. J.RovinskiyN.LewisJ. A.SardiM.SchackwitzW. S.MartinJ. A.et al (2014). Comparative genomics of Saccharomyces cerevisiae natural isolates for bioenergy production.Genome Biol. Evol.62557–2566. 10.1093/gbe/evu199
75
YangL. A.ChangY. J.ChenS. H.LinC. Y.HoJ. M. (2019). SQUAT: a sequencing quality assessment tool for data quality assessments of genome assemblies.BMC Genomics19:238. 10.1186/s12864-019-5445-3
76
YueJ. X.LiJ.AigrainL.HallinJ.PerssonK.OliverK.et al (2017). Contrasting evolutionary genome dynamics between domesticated and wild yeasts.Nat. Genet.49913–924. 10.1038/ng.3847
77
ZhaoX.TangJ.WangX.YangR.ZhangX.GuY.et al (2015). YNL134C from Saccharomyces cerevisiae encodes a novel protein with aldehyde reductase activity for detoxification of furfural derived from lignocellulosic biomass.Yeast32409–422. 10.1002/yea.3068
78
ZhaoY.StropeP. K.KozminS. G.MccuskerJ. H.DietrichF. S.KokoskaR. J.et al (2014). Structures of naturally evolved CUP1 tandem arrays in yeast indicate that these arrays are generated by unequal nonhomologous recombination.G3 (Bethesda)42259–2269. 10.1534/g3.114.012922
79
ZielezinskiA.GirgisH. Z.BernardG.LeimeisterC. A.TangK.DenckerT.et al (2019). Benchmarking of alignment-free sequence comparison methods.Genome Biol.20:144. 10.1186/s13059-019-1755-7
Summary
Keywords
bioethanol, yeasts, Saccharomyces cerevisiae, sugarcane, cachaça, genomics, alignment-free phylogenetics, evolution
Citation
Jacobus AP, Stephens TG, Youssef P, González-Pech R, Ciccotosto-Camp MM, Dougan KE, Chen Y, Basso LC, Frazzon J, Chan CX and Gross J (2021) Comparative Genomics Supports That Brazilian Bioethanol Saccharomyces cerevisiae Comprise a Unified Group of Domesticated Strains Related to Cachaça Spirit Yeasts. Front. Microbiol. 12:644089. doi: 10.3389/fmicb.2021.644089
Received
19 December 2020
Accepted
08 March 2021
Published
15 April 2021
Volume
12 - 2021
Edited by
Feng Gao, Tianjin University, China
Reviewed by
Xinqing Zhao, Shanghai Jiao Tong University, China; Jia-Xing Yue, Sun Yat-sen University Cancer Center (SYSUCC), China; Jose Paulo Sampaio, New University of Lisbon, Portugal
Updates

Check for updates
Copyright
© 2021 Jacobus, Stephens, Youssef, González-Pech, Ciccotosto-Camp, Dougan, Chen, Basso, Frazzon, Chan and Gross.
This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Jeferson Gross, jeferson.gross@unesp.br; jeferson.gross@googlemail.com
†Present Address: Timothy G. Stephens, Department of Biochemistry and Microbiology, Rutgers University, New Brunswick, NJ, United States Raul González-Pech, Department of Integrative Biology, University of South Florida, Tampa, FL, United States
This article was submitted to Evolutionary and Genomic Microbiology, a section of the journal Frontiers in Microbiology
Disclaimer
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.