Unravelling genomic diversity of Zygosaccharomyces rouxii complex with a link to its life cycle

Abstract

Zygosaccharomyces rouxii and the related species Zygosaccharomyces sapae (hereafter referred to as Z. rouxii complex) are protoploid hemiascomycete yeasts relevant in the elaboration and spoilage of foodstuff. Divergence of Z. rouxii complex before whole genome duplication, leading to the genus Saccharomyces, makes these yeasts very attractive for genome evolution study. Relatively little is known, however, about the diversity in this branch at the genetic and physiological levels. In this work, we investigated Z. rouxii complex, encompassing strains that in other works have been studied separately and comparing them in a comprehensive way. We showed that the majority of strains are unusually heterogeneous in their ribosomal DNA, a signal of relaxation of concerted evolution. Further analysis showed that they have hypervariable karyotypes, different levels of ploidy, and that housekeeping markers vary both in copy number and sequence. Overall, the results provide compelling evidence that the strains considered in this study are a complex of haploid, aneuploid and diploid mosaic lineages. The reproductive mode and life cycle of Zygosaccharomyces could lead to this unsuspected diversity.

Introduction

Genome structure is strictly dependent upon the life cycle and the breeding system of hemiascomycetous yeasts (Knop, 2006; Billiard et al., 2012). Different mechanisms of breeding (mating) are known that are implicated in individual heterozygosity and genetic diversity among yeast species (Knop, 2006). In the genus Saccharomyces, which includes the well-studied yeast Saccharomyces cerevisiae, the gene content and the rough organization of chromosomes are mostly preserved within isolates of the same species and among species (Kellis et al., 2003). This low variability is probably promoted by regular sexual cycles which help to preserve the genome organization (Fischer et al., 2000; Delneri et al., 2003). Natural S. cerevisiae strains are often diplontic, homothallic and heterozygous yeasts, which have the possibility of regenerating homozygous diploid cell from a haploid cell through self-fertilization, a mechanism interpreted as a way of genome renewal (Mortimer et al., 1994). The extent of intra- and interspecific genome variability is not well known for other hemiascomycetes. Recently, genome sequencing efforts have revealed an increasing number of yeast diploid mosaic species. In the GTG clade (so called because its constituent species translate CTG as serine instead of leucine), Candida albicans (Jones et al., 2004), Debaryomyces hansenii (Jacques et al., 2010), Pichia sorbitophila (Louis et al., 2012) and Millerozyma (Pichia) farinosa (Mallet et al., 2012) showed frequently complex diploid genomes, resulting from putative genetic exchanges between divergent strains and species at high frequency and thus referred to as ‘mosaic genomes’ (Zhaxybayeva et al., 2004).

Zygosaccharomyces clade encompasses industrially important halotolerant and osmotolerant yeasts that participate both in the elaboration and spoilage of foodstuff. Genome sequencing demonstrated that the type species of the genus, namely Zygosaccharomyces rouxii, is a protoploid lineage that did not undergo the whole genome duplication (WGD) leading to the genus Saccharomyces (Dujon et al., 2004). The genome of Z. rouxii strain CBS 732^T has a low number of duplicated genes (Souciet et al., 2009) and a haploid DNA content (Solieri et al., 2008). Given the tree topology of Kurtzman (2003) and the fact that the Vanderwaltozyma clade emerged after WGD (Scannell et al., 2007), Z. rouxii may be one of the closest relatives to the putative ancestral genome of S. cerevisiae (Souciet et al., 2009).

Zygosaccharomyces rouxii and its close relatives display an important phylogenetic position, yet the extent of genome variability and the life cycle of species in this clade have not been studied extensively. These yeasts were generally regarded as spending their vegetative life in the haploid phase and meiosis follows immediately after mating of two cells (Wickerham & Burton, 1960). Some genetic surveys have reported that Z. rouxii strains had unusual hypervariable karyotypes (De Jonge et al., 1986; Török et al., 1993; Oda & Tonomura, 1995). Other works have also suggested that Z. rouxii and phylogenetically closely related species exhibit mosaic genome structure with respect to some nuclear and mitochondrial phylogenetic markers (James et al., 2005; Solieri et al., 2007; Gordon & Wolfe, 2008; Suezawa et al., 2008; James & Stratford, 2011). This allowed us to refer the group of these species as Z. rouxii complex. For instance, a new diploid species distinct from haploid Z. rouxii, namely Zygosaccharomyces sapae, has been recently delineated on the basis of multi-gene phylogeny and physiological evidence (Solieri et al., 2008, 2013). Intragenomic heterogeneity in the internal transcribed spacer region and the intervening 5.8S rRNA gene and the distribution of housekeeping gene markers partially supported the hypothesis of reticulate speciation (Solieri et al., 2013). Additionally, strains extensively used in Z. rouxii genetics and molecular biology, such as ATCC 42981 and CBS 4837, harbour mosaic genomes with two copies of many genes. Sequence analysis of individual genes indicated that the parental strains contributing to the mosaic genome closely resemble Z. rouxii type strain CBS 732^T (t-subgenome) and Zygosaccharomyces pseudorouxii nom. inval. (p-subgenome) (James et al., 2005; Gordon & Wolfe, 2008). However, strain CBS 4837 was previously described as being a Z. rouxii haploid and heterothallic strain with the mating type MTa (Wickerham & Burton, 1960), a condition incongruent with the hypothesis of allodiploidization (James et al., 2005). Furthermore, strain CBS 4837 is able to mate with the sibling strain CBS 4838 (mating type MTα) (Wickerham & Burton, 1960), but the genomic features of the latter strain have not yet been investigated.

To cover these gaps and to unravel the genetic structure of species of Z. rouxii complex, we compared Z. rouxii CBS 732^T with six relatives on the basis of a polyphasic analysis. These strains included two that were retrieved from sweet food-grade and recently ascribed to Z. sapae; the remaining four, isolated from salt-rich foodstuff, are currently classified as Z. rouxii or of doubt assignment. The results provided evidence that the strains considered in this study are a complex of haploid, aneuploid and diploid mosaic lineages. The life cycle of Zygosaccharomyces strains was discussed to explain the high diversity within Z. rouxii complex.

Materials and methods

Strains and culture conditions

The strains used in this study are listed in Table 1. Strains were obtained as lyophilized stocks from CBS (Centralbureau voor Schimmelcultures, Delft, the Netherlands) and ATCC (American Type Culture Collection). Strain OUT7136 was kindly provided by Prof. Y. Kaneko (Osaka University, Japan). Cells were routinely cultured at 28 °C in YPD medium (1% yeast extract, 2% bactopeptone, 2% glucose) and 1-mL stocks containing 25% glycerol were stored at −80 °C. For ploidy estimation and flow cytometry (FCM) analysis, yeast nitrogen base (Difco) medium supplemented with 5% glucose (YNB5%G) was used.

DNA manipulations and PCR reactions

Genomic DNA extraction was carried out as described by Sambrook et al. (1989). Randomly amplified polymorphic DNA (RAPD) with M13 primer and microsatellite primed-PCR with synthetic primer (GTG)₅ [(GTG)₅-PCR] were performed and analyzed, as previously reported (Solieri et al., 2013). D1/D2 domains of the large subunit (LSU) of rRNA gene (LSU D1/D2) were PCR-amplified according to O'Donnell (1993). ITS regions were amplified as described by Kurtzman & Robnett (1998). The HIS3 and ZrSOD2 genes, encoding imidazole-glycerol-phosphate dehydrate gene and plasma-membrane Na⁺/H⁺ antiporter gene, respectively, were amplified, as previously reported (Solieri et al., 2007; Harrison et al., 2011). The mitochondrial cytochrome C oxidase II gene (COX2) was amplified with primers COII-5 and COII-3 as reported by Belloch et al. (2000).

Sequencing, subcloning of PCR-multiplied loci, and sequence analysis

The PCR products were purified using the DNA Clean and Concentrator kit (Zymo Research, Orange, CA) according to the supplier's instructions. Single-extension sequencing from both strands of the PCR products from genomic DNA (isolated from single colonies) was performed by a commercial DNA sequencing service provider (MWG-Biotech, Ebersberg, Germany) under BigDye terminator cycling conditions. For 26S D1/D2 amplicons, additional internal primers NL2a and NL3a were employed according to Kurtzman & Robnett (1998). Full-length sequences were assembled with seqman (Lasergene software package, DNAStar Inc., Madison, WI) and compared with reference sequences available in GenBank database using the NCBI blast program (Altschul et al., 1997). The sequences were aligned with clustal x program (Thompson et al., 1997) and used for phylogenetic analysis by the neighbour-joining (NJ) method (Saitou & Nei, 1987) with 1000 bootstrap replicates (Felsenstein, 1985).

When the sequencing chromatograms showed evidence of two different sequences for one locus, we subcloned the PCR products into the pGEM-T Easy Vector, following the manufacturer's instructions (Promega, Madison, WI). Eight to 12 plasmids per fragment were isolated using a Qiaprep Spin Miniprep kit (Qiagen, Westburg BV, Leusden, the Netherlands), and each individual plasmid was separately amplified, sequenced and analyzed as described above. Cloned sequences were screened for potential PCR chimeras analysis with the software package bellerophon (Huber et al., 2004).

FCM

For ploidy estimation by FCM, the strains were grown to the mid-log phase (OD_{600 nm} of 0.5 to 0.6) in YNB5%G at 28 °C with shaking at 180 r.p.m., and the cell cultures were then stained with propidium iodide (PI; Sigma) at a final concentration of 4 μg mL⁻¹, as previously reported (Solieri et al., 2008). To remove clumps and avoid clogging problems, two consecutive ultrasound pulses were applied, before FCM analysis, at 70% of total output for 6 s with an interval of 1–2 s between the two pulses using a Microson Ultrasonic Cell Disruptor XL (Misonix Inc., New Highway Farmingdale, NY). Samples were stored on ice and protected from light before FCM analysis.

All FCM experiments were performed in duplicate on an EPICs XL flow cytometer (Beckman Coulter Inc., Miami, FL) equipped with a 15-mW 488-nm air-cooled argon laser. Forward scatter (FSC) and side scatter were analyzed on logarithmic scale, and red fluorescence intensity (FL3) was analyzed on linear scale. For DNA content analysis, doublets and aggregates were discarded by electronic analysis of integral and height signals from the particles analyzed. The fluorescence histograms of 30 000 cells were generated using the gated data. Data were analyzed using EXPO 32 ADC Analysis (Advanced Cytometry Systems). Saccharomyces cerevisiae haploid strain BY4742, Euroscarf Accession No. Y10000 (MATα; his3-Δ1; leu2Δ0; lys2Δ0; ura3Δ0), and diploid strain BY4743, Euroscarf Accession No. Y23146 (MATa/α; his3-Δ1/his3-Δ1; leu2Δ0/ leu2Δ0; lys2Δ0/LYS2; met15-Δ0/MET15; ura3-Δ0/ura3-Δ0) were used for calibration, as previously reported (Solieri et al., 2008).

Chromosomal DNA preparation and pulse-field gel electrophoresis

Chromosomal DNA was isolated as described by Pribylová et al. (2007). The 1% agarose gels were prepared with Seakem GTG agarose (FMC Bioproducts) and run in 0.5× TBE buffer at 12 °C on a Bio-Rad contour-clamped homogeneous electric field (CHEF) apparatus for 98 h, according to running conditions reported by Solieri et al. (2008). The chromosomal DNA size marker of S. cerevisiae S288C (Bio-Rad Laboratories) was used to estimate the chromosome size. Pictures were digitally captured and analyzed using the BioDocAnalyze gel imaging and analysis system (Biometra, Göttingen, Germany).

Nucleotide sequence accession numbers

The DNA sequences described here were deposited in the GenBank data library under the accession numbers HE687308 and HE687309, and the accession numbers ranging from HE664087 to HE664107.

Results

Structuring Z. rouxii complex in two lineages

We applied two typing methods, such as M13-RAPD and (GTG)₅-PCR, to detect clonality and assess genetic relatedness among Z. rouxii CBS732^T and related strains listed in Table 1. Both methods generated complex fingerprinting patterns, with PCR products ranging from 650 to 1600 bp in M13-RAPD and from 300 to 3000 bp in (GTG)₅-PCR profiles. The strains fell into two groups, as identified by the cluster analyses of M13-RAPD patterns (Supporting Information, Fig. S1). (GTG)₅-PCR analysis also placed the strains into two clusters, with the exception of strain CBS 732^T, which did not fall into any of the clusters (Fig. S1). In both cases, one group was delineated that includes Z. sapae strains isolated from high sugary foodstuff, whereas the other one grouped strains from salt-rich food, including strains previously described as being allodiploid (James et al., 2005; Gordon & Wolfe, 2008). For this reason we referred to the latter cluster as the mosaic lineage.

Analysis of mtDNA

To address the question of the mtDNA inheritance within the mosaic lineage, sequence polymorphism of mitochondrial COX2 marker was investigated in strains ATCC 42981, CBS 4838 and OUT7136. Phylogenesis based on COX2 sequences showed that strains of the mosaic lineage shared the COX2 allele identical to that of Z. rouxii CBS 732^T, consistent with the hypothesis that these allodiploids derived from crosses between clades, one of which could be Z. rouxii. Moreover, the mosaic lineage differed from Z. sapae strains, which carried a unique allele at the mitochondrial gene COX2, divergent both from Z. rouxii (2.1% divergence) and Zygosaccharomyces mellis (2.5% divergence), a close relative of Z. rouxii (Fig. S2). These results support the differentiation of strains into two groups obtained by PCR fingerprinting.

Sequence variation within intragenomic rRNA gene arrays

The phylogenetic markers LSU D1/D2 domains and ITS region are members of the structural RNA components, which are highly repetitive and conserved owing to homogenization among unit copies through unequal crossing-over and gene conversion (Dover, 1982). Their analysis is described as a useful tool for studying genome origin and stability (James et al., 2009), so these markers were amplified and sequenced in strains CBS 4837, CBS 4838 and OUT7136 and compared with those available for Z. sapae ABT301^T and ABT601 and allodiploid strain ATCC 42981.

Sequencing of LSU D1/D2 domains revealed intra-individual heterogeneity in strains CBS 4838 and CBS 4837. Mixed populations of PCR products were generated, which need to be cloned before sequencing. The results showed that both strains harbour at least two types of LSU D1/D2 domains, suggesting that they could be either diploid heterozygotes (variation between rRNA gene clusters located on two homeologous/homologous chromosomes) or heterogeneous for this marker (variation harboured by repeating units arranged in the same rRNA gene cluster). In both strains, we arbitrarily referred to the two D1/D2 types as copy 1 (100% Z. rouxii similarity in both strains) and copy 2 (100% and 99.7% Z. sapae similarity for CBS 4837 and CBS 4838, respectively). The latter copy was slightly different in strains CBS 4837 and CBS 4838 at five variable sites (divergence 0.87% in 572 nt analyzed). The LSU D1/D2 copies 1 and 2 of strains CBS 4837 differed from each other in one single indel and 15 variable sites, including 10 transitions and five transversions, resulting in an overall divergence of 2.8%. The occurrence of two divergent copies of D1/D2 domain in CBS 4837 partially disagreed with that previously reported by another study, which described only one LSU D1/D2 sequence occurring in this strain (James et al., 2005). Peculiar PCR conditions have been demonstrated to favour preferential amplification of a few repeat types (Fenton et al., 1998) and this could explain the homogeneity in CBS 4837 LSU D1/D2 domains previously observed. The difference between the LSU D1/D2 copies 1 and 2 of strain CBS 4838 was of 1.9% and consisted of one single indel and 10 variable sites, including five transitions and five transversions. Similarly, sequence variation within intra-individual LSU D1/D2 genes was previously described in strain ATCC 42981 (Gordon & Wolfe, 2008). Although a member of the mosaic lineage, strain OUT7136 showed a single LSU D1/D2 domain sequence, 99% similar to that of Z. sapae.

Phylogenetic relationships based on alignment of variable LSU D1/D2 sequences both among and within individual yeast LSU D1/D2 domains are shown in Fig. 1. This analysis separated the strains into two groups, namely Z. rouxii and Z. sapae. The Z. sapae clade consists of six strains, of which Z. sapae ABT301^T and ABT601 and strain OUT7136 share a unique copy of LSU D1/D2 domain. The remaining three, ATCC 42981, CBS 4838 and CBS 4837, harbour two copies, one resembling that of Z. sapae and the other that of Z. rouxii.

In addition to the analysis of LSU D1/D2 domain variation, we applied a similar strategy to examine the patterns of variations in the ITS regions among and within strains CBS 4837, CBS 4838 and OUT7136. We detected a high degree of heterogeneity in strains CBS 4838 and CBS 4837 but not in OUT7136. ITS PCR products from strains CBS 4837 and CBS 4838 failed in sequencing and were subcloned. Approximately 12 resulting plasmids from each subcloning experiment were verified for the insert sequence. Comparisons of sequence variation patterns from conserved (5.8S) and less constrained (ITS 1 and ITS 2) regions showed that CBS 4837 and CBS 4838 harbour three divergent ITS sequences (arbitrarily referred to as copy 1–3). The set of ITS variants was only partially identical between the two strains. Strain CBS 4838 showed an ITS copy 1 identical to that of CBS 732^T and very similar to ITS copy 1 of Z. sapae ABT301^T and t-subgenome copy of ATCC 42981 (Fig. 2). The remaining copies 2 and 3 from strain CBS 4838 showed an identity of 88.4% and 87.2% compared with Z. rouxii CBS 732^T, respectively. Strain CBS 4837 displayed three variants in ITS regions, but all of them were significantly divergent from Z. rouxii (identity ranging from 91.7% to 87.2%). Conversely, sequencing of ITS domain from strain OUT7136 resulted in a single sequence identical to that of Z. rouxii CBS732^T (Fig. 2). This result disagreed with that reported for LSU D1/D2 domain, which was 100% identical to Z. sapae, and suggested that OUT7136 genome could harbour a patchy combination of Z. rouxii and Z. sapae sequences in different loci of rRNA gene arrays. The tree topology based on ITS sequences showed that clones from single individuals were scattered across the tree into three lineages, called S-group, R-group and recombinant group (Fig. 3). In particular, ITS copy 1 of CBS 4838 and copies 2 of CBS 4837 and CBS 4838 belonged to the R-group, whereas copy 1 of CBS 4837 and copy 3 from CBS 4838 and CBS 4837 belonged to the recombinant-group and displayed an nt identity higher than 98% with the strain IFO 0494, a uncharacterized Zygosaccharomyces yeast isolated from miso (Suezawa et al., 2008).

Overall, the data suggested that the rRNA gene composition of CBS 4838 and CBS 4837 is similar to that of ATCC 42981, but there were some differences. Strain ATCC 42981 has a standard allodiploid rearrangement with two types of ITS regions and two types of LSU D1/D2 domains (Gordon & Wolfe, 2008), whereas strains CBS 4837 and CBS 4838 have two copies of LSU D1/D2 domain and three copies of ITS regions (Table 2). Strains CBS 4837 and CBS 4838 differed also from Z. sapae strains, which showed heterogeneity in ITS regions but not in LSU D1/D2 domains (Solieri et al., 2007). Three possible scenarios can explain the observed data: (1) these strains may be either aneuploid or diploid, with each chromosome of the pair of homeologous/homologous chromosomes bearing one rRNA gene variant; (2) there are different tandem repeat variants arranged in their rRNA gene arrays located on the same chromosome; (3) there are divergent rRNA gene arrays dispersed among different chromosomes.

Sequencing of housekeeping markers

Previous studies have demonstrated that both Z. sapae and allodiploid strains ATCC 42981 and CBS 4837 possess divergent copies of the housekeeping markers ZrSOD and HIS3 (Iwaki et al., 1998; Solieri et al., 2007). Conversely, a single copy of ZrSOD was found in CBS 732^T (Kinclová et al., 2001). We performed allele-specific PCR assays and sequencing of ZrSOD and HIS3 genes for strains OUT7136 and CBS 4838 and results are shown in Table 2.

Analysis of ZrSOD locus showed that strain OUT7136 harbours two ZrSOD variants, one 99.6% identical to ATCC 42981 ZrSOD2 and the other one 100% identical to ATCC 42981 ZrSOD22. Similarly, strain CBS 4838 showed two copies, both 100% identical to ATCC 42981 ZrSOD2 and ZrSOD22 genes, respectively. The same combination of variants has been reported for strain CBS 4837 (James et al., 2005), suggesting that all members of the mosaic lineage share the identical pattern of ZrSOD variants. Differently from the mosaic lineage, Z. sapae possesses partially divergent ZrSOD variants: one copy was 98.6% identical to ZrSOD2-22, a variant previously found only in CBS 732^T, whereas the other one was similar to ATCC 42981 ZrSOD22 (99.8 and 98.7% identity in ABT301^T and ABT601, respectively) (Solieri et al., 2007). Finally, when we constructed a phylogenetic tree, two divergent lineages were delineated (Fig. 4, panel I). With the exception of CBS 732^T, each strain was represented in the two different lineages (designated as A and B, Fig. 4, panel I), suggesting that copies from different strains, orthologues (sequences within a single lineage), are more similar to each other than are the two copies within the same strain.

Similarly, allele-specific amplification and sequencing of HIS3 gene were performed using two different primer pairs targeting CBS 732^T (copy 1) and ATCC 42981 (copy 2) HIS3 genes, respectively. A successful amplification resulted from both primer pairs in strains OUT7136 and CBS 4838, confirming the presence of two different HIS3 variants in their genomes. Similar results were found in all the other strains considered in this study, but not in Z. rouxii CBS 732^T (Table 2). The topology of HIS3-derived NJ tree was very similar to that of the ZrSOD marker and confirmed that two HIS3 variants in each strain genome belong to two divergent phylogenetic lineages (Fig. 4, panel II). The presence of orthologues more similar to each other than the two paralogues within the same strain suggested that the ancestor of all strains belonging to both Z. sapae and the mosaic lineage had at least two variants of each tested housekeeping marker.

Genome size and ploidy level

Analysis of housekeeping markers indicates that these yeasts could be not simple haploids (the classical state of Z. rouxii strains) and that at least some parts of the genome are duplicated Accordingly, strain ATCC 42981 and Z. sapae have been determined previously to be diploid (Solieri et al., 2008). The occurrence of different variants in the genomes of the remaining strains CBS 4837, CBS 4838 and OUT7136 could suggest that they can be either diploid or aneuploid. We addressed the question by measuring total DNA content using an FCM approach. For CBS 732^T we obtained two peaks corresponding to two subpopulations, one with a 1n content (G₀/G₁ phases) and the other with a 2n content (G₂/M phases) of the haploid genome (see Fig. S3). A clear shift of the two peaks towards a double amount of DNA was observed for strains OUT7136, CBS 4837 and CBS 4838 (Fig. S3, panels I, II and III, respectively). In all cases, the first peak (G₀/G₁) was close to the second peak (G₂/M) of the haploid CBS 732^T, consistent with a diploid or at least an aneuploid status of these strains. The genome sizes were estimated from the G₀/G₁ peak median fluorescence intensity (MFI), using the genome size of S. cerevisiae as the calibration (Table 3). Zygosaccharomyces sapae strains ABT301^T and ABT601, and allodiploid ATCC 42981 have genome sizes ranging from 21.9 to 28.1 Mb (Solieri et al., 2008). Strains OUT7136, CBS 4837 and CBS 4838 have estimated genome sizes of 19.57–21.94 Mb, significantly higher than the estimated 9.8-Mb genome of CBS 732^T (Souciet et al., 2009).

We used PFGE to separate chromosomes and compare haploid and diploid strains in order to get an insight into the diversity of chromosome structure in Z. rouxii complex. Karyotyping showed that strains contained from six to 11 chromosomes. In some strains, a few bands showed a higher intensity (see Fig. S4), suggesting that some chromosomes overlap in size. The chromosome separation of strain CBS 732^T displays, as expected, a total of six bands, the 1.5-Mb band being a doublet made of chromosomes C and D. With the exception of a few chromosome length polymorphisms, the strains belonging to the mosaic lineage tend to harbour a similar number of chromosomal bands, amounting to eight chromosome bands per strain (Table 3). The chromosome sizes ranged from 1.2 to 2.2 Mb for strain OUT7136 and from 1.2 to 2.3 Mb for strains CBS 4837 and CBS 4838. Combination of FCM results and karyotypic data resulted in a ploidy level of 1.74 for OUT7136, 1.96 for CBS 4837, and 1.90 for CBS 4838, confirming the aneuploid/diploid status of these strains (Table 3). On the other hand, the chromosome separation of Z. sapae strains gave a higher number of chromosomes ranging from 10 to 11 bands, consistent with their higher genome sizes (Solieri et al., 2008). The heterogeneity in karyotype within Z. rouxii complex suggests that the genome has rearranged very fast upon the separation of single lineages.

Discussion

By evaluating genomic properties and nuclear and mitochondrial markers, we recognized two lineages well-separated from the Z. rouxii strain CBS 732^T, namely Z. sapae and the mosaic lineage (Table 4). With the exception of strain OUT7136, the strains considered in this study have an unusual level of heterogeneity in rRNA gene regions, but the pattern of this heterogeneity varies significantly between Z. sapae and the mosaic lineage (Table 4). The co-occurrence of rRNA gene variants in the genome of a single individual suggests relaxation of concerted evolution, a recombination-driven process that is responsible for homogenizing rRNA gene repeats (Birky, 1996). Relaxation was demonstrated to occur in transition stages of concerted evolution, when location of rRNA gene loci on nonhomologous chromosomes potentially disrupts concerted evolution, when organisms are polyploid, as a result of interspecific hybridization, or when the mutation rate exceeds the rate of concerted evolution, as in length variants in the intergenic spacer (as quoted by Paun et al., 2007; Albertin & Marullo, 2012). In particular, the concerted evolution model assumed that the co-evolving genes have to be organized in tandem arrays. A previous study showed that the repeats of the 5S ribosomal genes of filamentous fungi species are dispersed among the genomes and are escaping the concerted evolution model (Rooney & Ward, 2005). The Birth-and-Death evolution model was more appropriate to describe the dynamics of rRNA genes not organized in tandem arrays (Nei & Rooney, 2005). In that model, the genetic evolution is regulated by a balance between gene duplication, turnover and maintenance, which allows the occurrence of several haplotypes within the same genome. As of today, nothing is known about the organizational pattern of the ribosomal genes in Z. sapae and the mosaic lineage and therefore we cannot exclude that the Birth-and-Death evolution model could explain the heterogeneity detected in the present study.

Intra-individual variations in rRNA gene regions combined with the presence of divergent pairs of ZrSOD and HIS3 paralogues may be the expression either of a duplicated gene/genome structure, which preceded the divergence of the lineages, or of a hybridization event (Birky, 1996). Analysis of these nuclear housekeeping markers demonstrated that two homologues for each gene are present in diploid and aneuploid strains and that the orthologous sequences originating from one of the paralogous types are more similar in different strains to each other than to the sequences of the other paralogous type. Moreover, multiple lines of evidence demonstrated that Z. sapae and the mosaic lineage include diploid and aneuploid strains, in contrast to the haploid state reported for strains CBS 4837 and CBS 4838 in the early literature (Wickerham & Burton, 1960). These findings led us to suppose an evolutionary model in which the presumptive diploid-like status originated prior to sorting of Z. sapae and the mosaic lineage into different groups. Genome polyploidization has been detected in many asexual taxa (Otto & Whitton, 2000) and it is an effective way to buffer the genome against the effects of accumulating mutations (Kondrashov, 1994). Interestingly, putative asexual fungal species seem to have a propensity for within-individual rRNA gene variation (Simon & Weiß, 2008). In predominantly asexual yeasts, such as C. albicans (Rustchenko, 2007), Candida glabrata (Muller et al., 2009), Dekkera bruxellensis (Hellborg & Piškur, 2009), and D. hansenii (Jacques et al., 2010), pronounced karyotype variability is common and a relaxed control over the chromosome structure has been supposed to increase the genome variability and competitiveness (Poláková et al., 2009). Significantly, chromosome number and genome size vary significantly also within Z. rouxii complex. This chromosome polymorphism makes it difficult to believe that these strains regularly undergo meiotic recombination. The variability in genome organization may ‘isolate’ the strains and prevent them from successful recombination, if meiosis occurs. Furthermore, in Zygosaccharomyces, mating between heterothallic haploid vegetative cells should precede meiosis and sporulation (Wickerham & Burton, 1960). The reproductive mode, outlined in Fig. 5a and b, precludes the process of genome renewal by self-diploidization occurring in S. cerevisiae, and allows for reticulation. Under this hypothetical evolutionary scenario, mating between two divergent haploid cells would have resulted in a diploid ancestor, which lost meiosis ability and gave rise to stable diploid lineages (namely Z. sapae and the mosaic lineage), which clonally reproduce and independently evolve (Fig. 5c, I). The sequence divergence between variants of rRNA gene gene and housekeeping marker within an individual, as well as chromosome variability, supported such suggestion.

However, evidence of rare conjugated asci in Z. sapae (Solieri et al., 2013) and successful mating in laboratory conditions between sibling heterothallic strains CBS 4837 and CBS 4838 (Mori & Onishi, 1967) has been reported. These observations gave rise two hypotheses. The first is that mating between diploid cells results either in a transient tetraploid zygote that undergoes meiosis or in a massive chromosome loss in order to restore the diploid status (Fig. 5, panel C, II), a mechanism similar to what occurs during the parasexual cycle of C. albicans (Hull et al., 2000). The second is that the resulting transient tetraploid zygote is dikaryotic and the nuclei segregate into mitotic spores without karyogamy both in conjugated and unconjugated asci (Fig. 5, panel C, III and IV, respectively). In Zygosaccharomyces bailii, a close relative of Z. rouxii, absence of nuclear fusion and, as a result, the formation of mitotic spores have been reported (Rodrigues et al., 2003). The occurrence of binucleate cells in late stages of the Z. bailii life cycle was further confirmed by Dato et al. (2008). Similarly, M. farinosa allodiploid strains did not undergo meiosis, producing mainly diploid spores (Mallet et al., 2012). The abundance of Z. sapae cells with conjugation tubes and without zygote supported the second rather than the first hypothesis: mating between diploid cells might occur without nuclear fusion and meiosis.

Finally, the extent of rRNA gene variants within Z. rouxii complex, has evolutionary and diagnostic implications and points emphasizes the limit of classical phylogenetic analysis, which may lead eventually to erroneous trees. These phenomena create a network of paralogous sequence relationships potentially confounding accurate phylogenetic reconstruction. The question has been addressed by other authors (Mallet, 2007; Liti et al., 2006; Wu et al., 2008; Casaregola et al., 2011), who reported how ploidy and heterozygosity assessments, as well as phylogenetic networks and incorporation of population genetics into phylogenetic analysis, are pivotal tools to reveal mosaic species.

Overall, our results indicate that the group of strains included in this work is a complex of haploid and diploid heterogeneous species, including: (1) the haploid species Z. rouxii, e.g. CBS 732^T; (2) the diploid species Z. sapae, isolated from high sugar environments; (3) a diploid mosaic lineage that includes strains retrieved from salt environment. We hypothesize that the alternative mode of propagation to classical sexuality described in Saccharomyces species can account for this genome complexity and makes Z. rouxii complex prone to genome mosaicism and reticulate evolution. An increasing amount of evidence has shown that genome mosaicism is widespread in protoploid hemiascomycetes other than Z. rouxii complex, such as D. hansenii (Jacques et al., 2010), P. sorbitophila (Louis et al., 2012) and M. farinosa (Mallet et al., 2012). Notably, all these yeasts are able to survive to harsh conditions, such as a high extracellular concentration of osmolytes, suggesting that this genome complexity may occur in response to environmental cues (Gompert et al., 2006).

References

Supporting Information

Fig. S1.UPGMA dendrograms of Zygosaccharomyces rouxii CBS 732T and six relatives by M13-RAPD and (GTG)5-PCR fingerprinting.

Fig. S2. Phylogenetic relationships of Zygosaccharomyces rouxii complex, according to the COX2 mitochondrial gene neighbour-joining tree (Saitou & Nei, 1987).

Fig. S3. Fluorescence histograms of various Zygosaccharomyces rouxii complex strains after propidium iodide staining.

Fig. S4. Chromosome separation in Zygosaccharomyces rouxii complex.

Author notes