Alternative splicing regulates targeting of malate dehydrogenase in Yarrowia lipolytica.

Alternative pre-mRNA splicing is a major mechanism contributing to the proteome complexity of most eukaryotes, especially mammals. In less complex organisms, such as yeasts, the numbers of genes that contain introns are low and cases of alternative splicing (AS) with functional implications are rare. We report the first case of AS with functional consequences in the yeast Yarrowia lipolytica. The splicing pattern was found to govern the cellular localization of malate dehydrogenase, an enzyme of the central carbon metabolism. This ubiquitous enzyme is involved in the tricarboxylic acid cycle in mitochondria and in the glyoxylate cycle, which takes place in peroxisomes and the cytosol. In Saccharomyces cerevisiae, three genes encode three compartment-specific enzymes. In contrast, only two genes exist in Y. lipolytica. One gene (YlMDH1, YALI0D16753g) encodes a predicted mitochondrial protein, whereas the second gene (YlMDH2, YALI0E14190g) generates the cytosolic and peroxisomal forms through the alternative use of two 3'-splice sites in the second intron. Both splicing variants were detected in cDNA libraries obtained from cells grown under different conditions. Mutants expressing the individual YlMdh2p isoforms tagged with fluorescent proteins confirmed that they localized to either the cytosolic or the peroxisomal compartment.

Introduction

The hemiascomycetous yeast Yarrowia lipolytica has relatively few introns compared with other Opisthokonts such as vertebrates or even compared with basidiomycota and filamentous ascomycota. However, with introns in ∼15% of its genes,[1] Y. lipolytica has many more introns than any other hemiascomycete yeast, e.g. four times more than Saccharomyces cerevisiae.[2] Despite this paucity of introns, the S. cerevisiae transcriptome is proving to be more complex than previously appreciated.[3,4] Alternative transcription start sites have been reported and, for some intron-containing genes, transcript variants result from either alternative splicing (AS) or intron retention. Overall, alternative transcripts leading to different proteins remain uncommon in S. cerevisiae and only one case of AS (SRC1[5]) and one of intron retention (PTC7[6]) have been described, along with several cases of alternative start sites (e.g. YCAT[7] and SUC2[8]). The large majority of alternative transcripts in yeasts are presumed to be untranslated, as most of them contain premature termination codons and therefore activate the nonsense-mediated mRNA decay pathway that degrades defective mRNAs prior to translation (for review, see Behm-Ansmant et al.[9] and Stalder and Muhlemann[10]).

AS has been investigated at the genome scale in Y. lipolytica and all known modes of AS have been observed, i.e. exon skipping and alternative 5′- or 3′-splice sites, as well as intron retention and alternative start sites.[1] Several of these genes are involved in central carbon metabolism. The malate dehydrogenase (MDH) gene (YlMDH2, YALI0E14190g) was particularly attractive, as AS does not generate a premature termination codon, but leads to putative functional variants of the protein.

The MDH isoenzymes catalyse the conversion of malate into oxaloacetate with a concomitant reduction of NAD+ (see Minarik et al.[11] for a review). This reaction, which is reversible, represents an important step in the tricarboxylic acid (TCA) cycle, a central metabolic pathway occurring in mitochondria and critical for cellular respiration and ATP production.[12,13] The reaction also takes place in the glyoxylate cycle, which is a variant of the TCA cycle that shares three of its five enzymes with the TCA cycle, including MDH. In plants, nematodes and yeasts, the glyoxylate cycle is partly localized in peroxisomes, and thus contributes to the degradation of free fatty acids by providing coenzyme A and/or reoxidizing NADH for the β-oxidation cycle. MDH is also involved in gluconeogenesis, which takes place in the cytosol and plays a significant role in the malate/aspartate shuttle across the mitochondrial membrane.

In S. cerevisiae, three MDH genes (MDH1, MDH2 and MDH3) have been identified and each encodes an enzyme targeted for a different subcellular compartment. Mdh1p is localized in mitochondria,[14,15] Mdh2p in the cytosol[16,17] and Mdh3p in peroxisomes.[18,19] Large sequence variations have been observed in the MDH regions that encode the putative domains involved in compartmental targeting. For example, MDH1 encodes a 17-amino acid N-terminal extension that is absent from the other isozymes and is removed upon mitochondrial import.[15] Mdh3p has a unique C-terminal tripeptide sequence, Ser-Lys-Leu, characteristic of peroxisomal targeting sequences (PTSs).[20] This targeting sequence, called PTS1, is conserved in S. cerevisiae, but other conservative variants, such as Ala-Lys-Ile, exist among yeasts.[21-23]

Homologues of all three S. cerevisiae MDH genes are present in many other hemiascomycetous species.[24] In contrast, the oleagineous yeast Y. lipolytica contains only two MDH genes: YALI0D16753g (YlMDH1) and YALI0E14190g (YlMDH2). Whereas the first gene encodes a protein predicted to have a mitochondrial location, the second gene encodes a putatively cytosolic form. The presence of only two MDH genes in Y. lipolytica compared with the three genes possessed by other hemiascomycetes, and the fact that YlMDH2 is potentially subject to AS,[1] suggest a possible role for AS in the regulation of MDH compartmentalization in this yeast. We thus focused on these Y. lipolytica MDH genes to decipher the roles and regulation of the different isoenzymes. In this study, we show that YlMDH2 encodes the cytosolic and peroxisomal forms of MDH. This dual localization is due to the presence of an alternative 3′-splice site in the second intron of the gene, located at the 3′-end of the coding sequence. The use of this alternative splice site creates an mRNA that encodes a carboxyl-terminal peroxisomal targeting sequence (PTS1) and allows specific peroxisomal localization, which was revealed by colocalization studies.

Materials and methods

Strains and culture conditions

The yeast and bacterial strains used in this study are listed in Table 1. The bacterial strains, Mach1 T1® (Invitrogen, Cergy Pontoise, France) and DH5α (Gibco BRL, Rockville, MD, USA) used for the amplification of recombinant plasmids, were grown at 37°C in Luria-Bertani medium supplemented with 100 μg/l ampicillin or 40 μg/l kanamycin, if required.

Table 1.

Yarrowia lipolytica and bacterial strains used in this study

Strain	Genotype	Reference	Comments
Yeasts
PO1d	Mat A, leu2–270, ura3–302, xpr2–322	25
JMY1699	PO1d, URA3-YALI0E14190g Δintron 1335-1422	This study	Cytosolic form of YlMDH2
JMY1711	PO1d, URA3-YALI0E14190g Δintron 1335-1426	This study	Peroxisomal form of YlMDH2
JMY1685	PO1d, URA3-YALI0E14190g	This study	Wild-type form of YlMDH2
JMY2416	JMY1699, LEU2	This study	Cytosolic form of YlMDH2
JMY2426	JMY1711, LEU2	This study	Peroxisomal form of YlMDH2
JMY2428	JMY1685, LEU2	This study	Wild-type form of YlMDH2
JMY2451	PO1d, URA3-pPOX-eYFP-YALI0E14190g Δintron 1335-1422	This study	Cytosolic form of YlMDH2
JMY2457	PO1d, URA3-pPOX-eYFP-YALI0E14190g Δintron 1335-1426	This study	Peroxisomal form of YlMDH2
JMY2459	PO1d, URA3-pPOX-eYFP-YALI0E14190g	This study	Wild-type form of YlMDH2
JMY2499	JMY2451, LEU2-pTEF-REDSTAR2	This study	Cytosolic form of YlMDH2
JMY2501	JMY2457, LEU2-pTEF-REDSTAR2	This study	Peroxisomal form of YlMDH2
JMY2500	JMY2459, LEU2-pTEF-REDSTAR2	This study	Wild-type form of YlMDH2
JMY2493	JMY2451, LEU2-pTEF-REDSTAR2-SKL	This study	Cytosolic form of YlMDH2
JMY2496	JMY2457, LEU2-pTEF-REDSTAR2-SKL	This study	Peroxisomal form of YlMDH2
JMY2495	JMY2459, LEU2-pTEF-REDSTAR2-SKL	This study	Wild-type form of YlMDH2
Bacteria	Plasmid/genotype
JME802	JMP62 LEU2ex, pPOX2 expression vector with the excisable LEU2ex marker	66
JME803	JMP62 URA3ex, pPOX2 expression vector with the excisable URA3ex marker	56
JME1018	JMP62 URA3 pPOX eYFP N-Ter	This study
JME1392	JMP62 LEU2ex pTEF Redstar2-SKL	This study
JME1394	JMP62 LEU2ex pTEF Redstar2	This study

Growth media and conditions used for Y. lipolytica have been previously described.[25] The yeasts were grown on rich medium YPD and minimum medium YNB[26] or YNBE (YNB supplemented with 0.1% yeast extract), supplemented with 0.2% (w/v) casa-amino acids or 0.05 g/ml uracil, if required. Carbon sources were dextrose (D), oleic acid (OA), tributyrin (TB), triolein (TO), alkane (hexadecane), acetate (ammonium or sodium), glycerol or ethanol at a final concentration of 2%. The hydrophobic substrates were emulsified by sonication in a mixture containing 20% OA, TB or TO and 0.625% (v/v) Tween 40. For solid media, 2% agar was added. Yeast strains were grown at 28°C. Growth on 96-well plates in 100 µl was performed under constant agitation with glucose, malate, succinate or citrate at 0.5% as the carbon source. Growth was monitored by measuring the optical density at 600 nm every 10 min using a microtiter plate reader (Biotek, Colmar, France). For each strain and condition, four biological replicates were performed. Calculations of average OD, blank reduction, the lag phase (µ) and the maximum OD were performed using the grofit package in R.[27] Differences between growth curves were evaluated by analysis of variance with R statistical software.[28]

Molecular biology techniques

Standard molecular biology techniques were used throughout this study.[29] Restriction enzymes were obtained from Ozyme (Saint-Quentin-en-Yvelines, France). Genomic DNA from yeast was prepared as previously reported.[30] Polymerase chain reaction (PCR) amplifications were performed using an Eppendorf 2720 thermocycler, Pyrobest DNA polymerase (Lonza, Levallois-Perret, France) and the primers are listed in Table 2. PCR fragments were purified using the QIAquick PCR Purification Kit (Qiagen, Courtaboeuf, France). DNA fragments were recovered from agarose gels using the QIAquick Gel Extraction Kit (Qiagen).

Table 2.

Oligonucleotides used in this study

Primer	Sequence	Purpose
E14190F1	GCCTACATCTACCTTGAC	cDNA sequencing
MDH1-Start	CAGATCGATTAACTGACTAGGTTAAAGCTGTCGTTGCCGGAGC	MDH disruption cassette
MDH1-Sa	AGGCCTAGGTTTAGATCCTAGCTAGAATGGTTAGTGATCGTGTAGTTCAAATGG	MDH disruption cassette
MDH1-04	AGGCCTAGGTTTAGATCCTAGTTGGCAGGAGGAGGGTTAACAATG.	MDH disruption cassette
MDH1-09	AGGCCTAGGTTTAGATCTTGGCAGGAGGAGGGTTAACAATG.	MDH disruption cassette
MDHS1	AACAGGCAACAATGGCATGC	Probe for Southern blot
MDHS2	AACTGGATTCGCTTGACGAG	Probe for Southern blot
MDHgfp1	CGCTGATCAAgtgagttatcatggtgggag	Fusion cassette
MDHgfp2	CGCTGATCAATGGTTAAAGCTGTCGTTGCC	Fusion cassette
RedSKL1	CACGGATCCCACAATGAGTGCTTCTTCTGAAGATG	Localization cassette
RedSKL2	TGGTGCCTAGGCTGCTTAAAGCTTGGACAAGAACAAGTGGTGTCTACC	Localization cassette
RedNoSKL2	TGGTGCCTAGGCTGCTTACAAGAACAAGTGGTGTCTACC	Localization cassette

Restriction sites are underlined: ClaI ATCGAT, AvrII CCTAGG, BclI TGATCA, BamHI GGATCC. Stop codons and methionine codons (or the first nucleotide of the codon) are in bold. Nucleotides complementary to the S-K-L codons are in bold italic. Intron nucleotide sequence is in lower case.

Genetic modifications

Three constructions were made to determine the phenotypes of the strains containing deletions of either the short or long introns of YlMDH2. The primers used for PCR amplification are listed in Table 2. The forward primer contained a ClaI restriction site and stop codons in the three phases. The reverse primers contained an AvrII site and sequences corresponding to the wild-type gene and to the two variants. The corresponding ClaI-AvrII PCR fragments were cloned into the JMP62 URA3ex vector[31] previously digested with ClaI-AvrII. The resulting plasmids were purified, digested with BamHI (a BamHI site is present in the middle of the gene), and used to transform the PO1d strain by the lithium acetate method.[32] After single crossing-over, the mutated version had integrated at the YlMDH2 locus under its own promoter and the genomic version was invalidated by the introduction of stop codons at the 5′-end of the wild-type gene (Supplementary Fig. S1). The resulting strains JMY1699, JMY1707 and JMY1711 (Table 1) expressed the predicted cytosolic, peroxisomal and wild-type versions of the YlMDH2 gene products and were confirmed by southern blot (data not shown). Prototrophic strains were obtained by transforming the different mutants with the LEU2 marker, which gave rise to the strains JMY2416, JMY2426 and JMY2428, respectively (Table 1).

We monitored localization of the MDH proteins in different subcellular compartments by expressing recombinant genes that encoded the peroxisomal or cytosolic mutant or the wild-type protein fused to the eYFP fluorescent protein tag at their N-termini. The three variant constructs were PCR amplified as BclI-AvrII fragments using specific primers (Table 2) and then cloned into the JMP62 URA3ex-YFP-N vector previously digested by BamHI and AvrII. This vector was derived from the JMP62 URA3ex vector[31] and allows N-terminal, eYFP-tagged proteins to be produced under control of the POX2 promoter. The three different expression cassettes were excised from their respective vectors by NotI digestion, gel purified, and subsequently used to transform PO1d to create the JMY2451, JMY2457 and JMY2459 strains, respectively (Table 1). These strains were subsequently transformed with constructions containing the sequence of the fluorescent protein RedStar2. The RedStar2 gene was PCR amplified from the plasmid pYM43[33] using primers containing BamHI and AvrII restriction sites (Table 2) and cloned into the JMP62 LEU2ex pTEF vector previously digested by BamHI and AvrII. The JMP62 LEU2ex pTEF plasmid corresponds to JMP62 LEU2ex in which the pPOX2 promoter has been exchanged with the pTEF promoter.[34] Primers were designed to add or not add a PTS1 consensus sequence (Ser-Lys-Leu; SKL) onto the C-terminus of RedStar2 and allow expression of cytosolic or peroxisome-targeted RedStar2 protein. Each RedStar2 expression cassette was released by NotI digestion, gel purified, and subsequently used to transform the JMY2451, JMY2457 and JMY2459 strains. The resulting strains harboured one of the three different eYFP-tagged MDH protein expression constructs controlled by the POX2 promoter and an addition construct for expressing either the RedStar2 or the RedStar2-SKL protein under control of pTEF. All the strains generated are listed in Table 1.

Microscopy

Images were acquired using a Zeiss Axio Imager.M2 microscope (Zeiss, Le Pecq, France) using a 100× objective and Zeiss filters 45 and 46 for fluorescent microscopy. Image acquisition was performed with the Axiovision 4.8 software (Zeiss). All images were post-processed similarly (background reduction) and merged using ImageJ software (http://rsbweb.nih.gov/ij/). Strains were grown 12h on YNBE 2% OA, YNBE 2% D or YNBE with both 1% D and 2% OA.

RNA-seq analysis of YlMDH2 transcripts

RNA-seq data produced in the lab were screened for the different splicing variants of YlMDH2. Total RNAs from cells cultivated in six different conditions (glucose, OA, TO, TB, alkane, and glycerol) were prepared using the Qiagen RNeasy kit (Qiagen). mRNAs were purified by the selection of poly(A)+ transcripts, which were then sequenced by the Illumina Solexa sequencing technology with either a Genome Analyzer IIX or a HiSeq sequencing system. Thirteen to 29 million single-end reads were generated per sample with a read length of 36 nt, 50 nt or 100 nt (accession numbers E-MTAB-939 and E-MTAB-940). The reads were aligned to both forms of YlMDH2 spliced transcripts (long or short intron) using SOAPaligner version 2.20.[35] Only reads aligned to the exon–exon junction (34 nt, 48 nt or 98 nt from each exon, depending on the read length) of intron 2 were counted and the ratio between the two isoforms was determined. We used Fisher's exact tests to determine whether the expression levels of YlMDH2 differed in the different growth conditions and to determine the significance of the differences observed among the AS ratios.

Bioinformatics

Homologues of S. cerevisiae MDH were identified in fully sequenced Y. lipolytica genomes from Génolevures (http://www.genolevures.org/) using BLASTp.[36] Homologues were also found in Pichia pastoris (Pichia_GS115.pep_0509; https://bioinformatics.psb.ugent.be/gdb/pichia/) and in Schizosaccharomyces pombe (pompep_17022010; http://www.sanger.ac.uk/Projects/S_pombe/). Amino acids sequences were then aligned using MUSCLE[37] or Multalin[38] and columns of gap-containing residues were removed manually. Phylogenetic trees were constructed with the Neighbour-Joining algorithm using ClustalX.[39]

Predictions for mitochondrial protein targeting were performed using MITOPROT II version 1.101 (http://ihg.gsf.de/ihg/mitoprot.html).[40] PTS1 targeting sequence predictions were performed using the PTS1 predictor using default cut-off (http://mendel.imp.ac.at/mendeljsp/sat/pts1/PTS1predictor.jsp).[41]

Results

Malate dehydrogenase protein families in hemiascomycetous yeasts

In S. cerevisiae, the three genes that encode MDH are: YKL085W (MDH1), which encodes the mitochondrial enzyme;[14,15] YOL126C (MDH2), which encodes the cytosolic enzyme;[16,17] and YDL378C (MDH3), which encodes the peroxisomal enzyme.[18,19] The three isozymes show amino acid identities ranging from 43 to 50%.[19] A Blast search for S. cerevisiae MDH homologues in the Y. lipolytica genome[42] identified only two genes, YALI0D16753g (YlMDH1) and YALI0E14190g (YlMDH2). Both genes exhibit the highest levels of amino acid conservation with the S. cerevisiae MDH1 gene, with 63% identity and 79% similarity for YALI0D16753g and 47% identity and 62% similarity for YlMDH2. The computed proteome families of the Génolevures database[24] revealed that both Y. lipolytica MDH genes belong to the GL3R0092 family, which includes all of the MDH proteins encoded by the nine fully sequenced yeast genomes in the Génolevures database. All of these genomes contain three MDH genes. The only exceptions are Y. lipolytica (two genes) and Lachancea (Kluyveromyces) thermotolerans, which possesses three additional genes. We found these additional genes (KLTH0D00440g, KLTH0G19536g and KLTH0G19558g) to be closely related to each other (64–66% identity and 80–83% similarity) but very divergent from all other hemiascomycetous MDH genes. Due to their mutual similarity and their localization (one is subtelomeric and two are repeated in tandem), they likely represent species-specific gene duplications. Given their divergence from the other hemiascomycete proteins, we did not consider them in our phylogenetic analysis. As Y. lipolytica is relatively isolated among the hemiascomycetes from a phylogenetic point of view, we added two other yeasts, not included in the previous yeast set, more closely related to Y. lipolytica: Arxula adeninivorans (recently sequenced and annotated by the Génolevures consortium, to be published) and P. pastoris.[43] Like Y. lipolytica, these two species contain only two MDH genes. We used these data to construct a phylogenetic tree of MDHs in the hemiascomycetous yeasts with the unique MDH gene from S. pombe as the outgroup (Fig. 1A). Simultaneously, we performed in silico sequence analysis to predict mitochondrial and peroxisomal localization for the different genes using MITOPROT and PTS1 prediction online software (Fig. 1B). The phylogenetic tree highlights clear distinctions between the three types of MDH. The most conserved branch corresponds to the mitochondrial MDH, with high levels of amino acid conservation among the different proteins: 61–85% identity and 76–92% similarity (Fig. 1A). The peroxisomal branch is less conserved, with minimums of 34% identity and 44% similarity for the most distantly related proteins. The cytosolic MDHs segregated into two branches, with similarity and identity values close to those obtained for the peroxisomal branch (35–62% identity and 45–77% similarity). According to the localization predictions, only one gene, DEHA2F09020g from Debaryomyces hansenii, which has a PTS1 peroxisome targeting sequence (C-terminal sequence SKL and PTS1 predictor score of −8.044), is mislocalized in the cytosolic branch (Fig. 1). Eremothecium (Ashbya) gossypii has three MDH genes, but the localization prediction failed to identify any PTS1 sequence in these three genes. The ERGO0D15114g protein belonging to the peroxisomal branch harbours an unexpectedly long C-terminal extension of ∼150 amino acids. In a recently corrected version of the Ashbya Genome Database (http://agd.vital-it.ch/index.html[44]; version 11 of 26 May 2011), the gene AGOS_ADR252W encodes a protein of 339 amino acids (vs 486 aa in ERGO0D15114p) and has a putative PTS1 with a significant score (C-terminal sequence ARL and PTS1 predictor score 3.865).

Figure 1.

(A) Phylogenetic tree of the MDHs from 11 fully sequenced hemiascomycetous yeast species using the unique MDH gene from S. pombe (SPCC306.08) as an outgroup. The yeasts are Candida glabrata (CAGL), D. hansenii (DEHA), Kluyveromyces lactis (KLLA), Lachancea thermotolerans (KLTH), Lachancea kluyveri (SAKL), Zygosaccharomyces rouxii (ZYRO), Eremothecium gossypii (ERGO), P. pastoris (Pipas), A. adeninivorans (ARAD), Y. lipolytica (YALI) and S. cerevisiae. The blue zone indicates genes with a typical mitochondrial targeting sequence and the pink zone indicates genes with a potential PTS. The percentages of amino acid identity and similarity among the proteins of each group (peroxisomal, mitochondrial and cytosolic) are indicated next to each group. (B) In silico predictions of protein targeting. The MITOPROT score, PTS1 predictor score and C-terminal sequence are indicated for each MDH.

The product of the Y. lipolytica YlMDH1 gene is predicted to be targeted to the mitochondria (MITOPROT probability 0.9891), whereas the product of the YlMDH2 gene is not (MITOPROT probability 0.3528). Neither of these genes harbours any PTSs. This is in line with their positions in the phylogenetic tree, i.e. YlMDH1 is positioned within the mitochondrial branch, and YlMDH2 is within the cytosolic branch of the tree. Similarly, A. adeninivorans and P. pastoris both have a mitochondrial form (ARAD1B14036g and Pipas_chr2-1_0238, respectively) and a cytosolic form (ARAD1D04884g and Pipas_chr4_0815, respectively) according to phylogeny and signal sequence prediction. Consequently, for these three closely related species, sequence analysis of their complete genomes did not identify any peroxisomal forms of MDH.

The YlMDH2 gene is alternatively spliced

The peroxisomal form of MDH has been reported to be mandatory for fatty acid degradation in S. cerevisiae.[45] Yarrowia lipolytica is an oleagineous yeast, with an affinity for fatty acids as a carbon source. Thus, the absence of a peroxisomal form of MDH was unlikely. A deeper analysis of the sequence of YlMDH2 revealed a putative alternative 3′-splice site for the second intron, defined by comparison with the splicing sequence pattern of the genosplicing database (http://genome.jouy.inra.fr/genosplicing2). The gene model of YlMDH2 is represented in Fig. 2A, with the positions of the putative intron boundaries (5′- and 3′-splice sites) and the branch point (BP) at the exon 2–intron 2–exon 3 junctions in Fig. 2B. The putative AS event leads to a very minor change in the length of the mRNA sequence (4 nt) but introduces a 3′-frameshift and thus generates different stop codons. If the intron is spliced at the upstream 3′-splice site (short intron), the amino acid sequence of the C-terminus encoded by the YlMDH2 mRNA is P-A-N. If splicing occurs at the downstream 3′-splice site (long intron), the C-terminal amino acid sequence ends with A-K-I (Fig. 2B), which is a characteristic PTS1 sequence identified in several yeasts.[21,23]

Figure 2.

Gene model for YlMDH2. (A) Schematic representation of alternative transcripts from the multi-intronic MDH gene YlMDH2. Exons are represented by grey rectangles and introns are symbolized by thin black articulated lines. Vertical bars on each of the three phases (0, +1 and +2) represent in-frame stop codons. (B) Sequence representation of the 3′ regions centred on the second intron. Coloured circles indicate the 5′-splice site, the BP and the two 3′-splice sites used to generate the two mRNA variants. Exon parts are represented by grey rectangles and the two putative C-terminal protein sequences are indicated.

We generated three cDNA libraries from Y. lipolytica cells cultivated in media containing either glucose (exponential and stationary phases) or oleate (exponential phase) as the carbon source[1] and obtained 30 cDNA clones corresponding to YlMDH2 transcripts. All of these clones derived from cells harvested from glucose medium: 10 clones from cells harvested in the exponential phase (expo library) and 20 clones from cells harvested in the stationary phase (stat library). We sequenced the 3′-ends of 11 of these YlMDH2 cDNA clones and found that seven corresponded to mRNAs spliced at the upstream 3′-splice site (short intron 2 splice) and four corresponded to mRNAs spliced at the downstream 3′-splice site (long intron 2). As expected from our model, the two groups of cDNAs differed by only four nucleotides (Fig. 2B) and thus validated AS for YlMDH2. The different numbers of cDNAs for the two spliced forms (five cDNAs with a short intron and two cDNAs with a long intron in the expo library, and only four cDNAs with a long intron, none with a short intron, in the expo library) was a first indication that the AS was regulated by the environment of the cells, i.e. the growth conditions. However, due to the low number of cDNAs sequenced, the statistical significance of this difference could not be confirmed. This prompted us to investigate larger data sets using RNA-seq experiments.

In a first experiment (Illumina Solexa with single reads of 36 bp; Table 3) involving six different growth conditions, we retrieved a total of 1829 RNA-seq reads containing the exon–exon junction of the second intron of YlMDH2: 1181 reads corresponded to the short intron variant, and 648 reads corresponded to the long intron variant. The short intron/long intron usage ratio varied from 0.9 to 5.05 according to the carbon source used for growth (Table 3). These results strongly suggest that environmental conditions influence the regulation of YlMDH2 AS. For example, for Y. lipolytica grown on TB medium, the number of transcripts containing the upstream 3′-splice site is five times the number of transcripts containing the downstream 3′-splice site. This ratio is significantly different from those found in all other conditions (P < 0.0001, for all comparisons to TB). In contrast, the levels of both variants were equivalent in cells grown on alkane. Additional RNA-seq experiments with two replicates on glucose (HiSeq with single reads of 50 bp; Table 3) and OA (HiSeq with paired-end reads of 2 × 100 bp) media confirmed that both ratios are conserved. However, these promising results on splicing regulation are presented as preliminary results since they came from a single time point, with a single concentration of carbon source and using RNA-seq as a single method of quantification. They obviously should be confirmed by additional experiments and alternative methods to understand the kinetic and the carbon source concentration-dependent regulations. This will be further investigated and will constitute the aim of a separate study.

Table 3.

Number of reads specific for the short or long intron 2 of YlMDH2 obtained by RNA-seq analysis under various growth conditions

Carbon source	Illumina Solexa sequencing technology and read length	Total reads	Reads mapping short intron 2	Reads mapping long intron 2	Ratio short/long
Triolein	GAIIX single reads 36 nt	19 589 043	94 (57%)	71 (43%)	1.33
Tributyrin	GAIIX single reads 36 nt	15 996 762	298 (83%)	59 (17%)	5.05a
Glycerol	GAIIX single reads 36 nt	22 935 285	308 (62%)	186 (38%)	1.66
Alkane	GAIIX single reads 36 nt	23 903 839	99 (47%)	110 (53%)	0.9b
Glucose	GAIIX single reads 36 nt	13 022 908	121 (68%)	56 (32%)	2.16
Glucose	Hiseq single reads 50 nt	12 965 094	202 (67%)	86 (33%)	2.21
Glucose	Hiseq single reads 50 nt	10 592 928	203 (63%)	116 (37%)	1.70
Oleic acid	GAIIX single reads 36 nt	19 270 870	261 (61%)	166 (39%)	1.57
Oleic acid	Hiseq paired-end 100 nt	25 974 957	118 (60%)	78 (40%)	1.51
Oleic acid	Hiseq paired-end 100 nt	13 316 230	260 (67%)	130 (33%)	2.00

aThe ratio on tributyrin is statistically different from that on all other media (P < 0.0001).bThe ratio on alkane is statistically different from that on all other media, except on triolein (P = 0.0765).

Absence of either the peroxisomal or the cytosolic form of the MDH does not affect growth rate, irrespective of the carbon source

In order to evaluate the biological roles of the two MDH isoforms encoded by YlMDH2, we constructed strains that expressed cDNAs encoding either the cytosolic (cyto), peroxisomal (pero) or wild-type (wt) form of YlMdh2p at the YlMDH2 locus in the Y. lipolytica PO1d strain (JMY2416, JMY2426 and JMY2428 strains, respectively). After complementation to restore prototrophy, strains were screened for growth on various substrates known to affect the growth of S. cerevisiae MDH mutants, i.e. glucose, acetate, ethanol and OA.[17-19,46,47] Surprisingly, we did not observe any growth rate differences among the three strains in either liquid or solid culture using these substrates as the sole carbon sources (see Fig. 3 for growth on solid substrates and Fig. 4A for glucose on liquid medium). We also screened these strains for growth on other lipid substrates, i.e. TB and TO, but failed to identify any specific growth phenotypes (Fig. 3). We also performed growth rate comparisons in a 96-well plate system using various carbon sources. Under these experimental conditions, we observed significant increases in the growth rate and maximum cell growth for the cyto strain (JMY2416) for all carbon sources tested, but observed no differences in the lag phase. For example, the P-values calculated for four replicates grown on glucose are 1.544e−2 and 1.184e−4 for the growth rate and maximum cell growth, respectively (Fig. 4B). We speculated that this result may be due to low oxygenation of the culture in a microtiter plate. In order to simulate low oxygen levels in a flask culture, we complemented the growth medium with antimycin A, a drug that impairs the respiratory chain, at a concentration that allows Y. lipolytica to grow but with a reduced rate (0.1 µg/ml). However, the growth rates of all strains were similarly reduced (data not shown), which indicated that the increased growth rate of the cyto strain is not directly linked to respiration efficiency. This phenomenon could be linked to oxygen access for other pathways or could reflect a different morphological state that affects optical density and may not rely on growth rate directly.

Figure 3.

Expression of alternatively spliced forms of YlMDH2 and growth on different substrates. Serial dilutions (serial dilution factor of five) of cultures of the wild-type (wt—JMY2428) strain, the cytosolic variant (cyto—JMY2416) and the peroxisomal variant (pero—JMY2426) were inoculated on YNB medium supplemented with different carbon sources. No growth differences between the mutants were detected; both of them were able to grow on all the media.

Figure 4.

Comparative growth of YlMDH2 mutants. (A) Growth curves in flasks with agitation in YNBE medium with glucose 2%. (B) Growth curves on 96-well plates in YNBE with glucose 0.5%. Coloured curves represent the different splicing mutants of YlMDH2: the cytosolic form is in black (MDHc), the peroxisomal form is in red (MDHp) and the wild-type is in green (MDHwt). OD, optical density measured at 600 nm.

Malate dehydrogenase localization

In order to locate the two MDH isoforms produced by AS, we constructed vectors to express the eYFP fluorescent protein fused to the N-termini of the cytosolic and peroxisomal MDH isoforms and to the N-terminus of MDH expressed from the wt YlMDH2 gene. Yarrowia lipolytica strains containing these constructs were subsequently transformed with a construction encoding the RedStar2 fluorescent protein without or with a PTS1 sequence (Ser-Lys-Leu; SKL) and used for colocalization studies. The native RedStar2 protein was cytosolic, whereas the RedStar2-SKL protein was targeted to the peroxisomes, which were visible as fluorescent dots in the cells (Fig. 5). In the strain expressing the predicted cytosolic form of eYFP-MDH, green fluorescence colocalized with red fluorescence produced by cytosolic RedStar2 (Fig. 5; JMY2499), but did not colocalize with the peroxisomal form of RedStar2 (Fig. 5; JMY2493). Conversely, the predicted peroxisomal form of eYFP-MDH colocalized with the peroxisomal form of RedStar2 (colocalization of fluorescent dots; Fig. 5; JMY2496) but not with cytosolic RedStar2 (Fig. 5; JMY2501). These colocalization experiments confirmed that the upstream splice variant mRNA encodes a cytosolic form of MDH, whereas the downstream splice variant encodes MDH with a C-terminal PTS1 and is peroxisomal. For all these constructions, we observed similar localizations and colocalizations with either glucose or OA as the carbon source. For the eYFP-tagged wild-type copy of YlMdh2p, which encoded both MDH isoforms, we observed clear cytosolic fluorescence (Fig. 5; JMY2500) but no clear labelling of peroxisomes (Fig. 5; JMY2495). However, overlapping fluorescence might have prevented distinguishing the peroxisomal and cytosolic localizations, especially since the expression level of the cytosolic isoform is higher than that of the peroxisomal form, as suggested by the cyto/pero RNA ratio (Table 3).

Figure 5.

Colocalization of the two YlMdh2p isoforms with cytosolic or peroxisomal forms of the RedStar2 protein. eYFP-tagged peroxisomal YlMdh2p (eYFP-MDHp), cytosolic YlMdh2p (eYFP-MDHc) or wild-type YlMdh2p (eYFP-MDHwt) were co-expressed with either the peroxisomal (redstar2p) or cytosolic (redstar2c) forms of the RedStar2 protein. For each strain, both proteins (MDH and RedStar2) were visualized simultaneously by fluorescence microscopy. Cells were imaged after 12h of growth in YNBE 2% OA using differential interference contrast (DIC) for eYFP fluorescence (eYFP, green) and for RedStar2 fluorescence (redstar2, red). eYFP and RedStar2 images were merged (right panel). Yellow colour indicates overlapping fluorescence and evidences colocalization.

AS of MDH is unique to Y. lipolytica

A comparative study of the MDH gene family revealed that P. pastoris and A. adeninivorans also contain only two genes coding for MDH, and in each case, a peroxisomal form is absent. Whereas the P. pastoris MDH genes do not contain any introns, model predictions for the gene that encodes the cytosolic form of MDH in A. adeninivorans revealed the presence of two introns. The second intron is located at exactly the same position as the second intron of YlMDH2 (Supplementary Fig. S2). Unexpectedly, the terminal exon has undergone a large expansion, probably by fusion with a downstream CDS, which is conserved in synteny in Y. lipolytica. Despite the presence of the same positioned intron, no alternative 3′-splice site was detected in A. adeninivorans.

Discussion

Yarrowia lipolytica genes contain more introns than those of S. cerevisiae (15 and 4%, respectively) or other hemiascomycetes, and Y. lipolytica exhibits all known forms of AS.[1,2] This suggests possible roles for splicing in the post-transcriptional regulation of gene expression and/or in generating additional proteome complexity in Y. lipolytica, and that splicing might have a much higher impact in Y. lipolytica than in S. cerevisiae. As an example of protein targeting regulation, we found clear evidence that the MDH gene YlMDH2 is alternatively spliced and leads to the production of two proteins destined for distinct cellular compartments: cytosol and peroxisomes. The HAC1 gene in Y. lipolytica[48] is also subject to functional AS; however, this gene possesses a non-spliceosomal intron and the splicing mechanism is mediated under stress conditions by Ire1p, a transmembrane kinase/endonuclease.[49] Thus, the YlMDH2 gene constitutes the first case of functional AS involving a spliceosomal intron in Y. lipolytica.

Overall, few cases of AS have been described in hemiascomycete yeasts, which agrees with their low intron content. The only examples that contribute to the diversity of both the transcriptome and proteome were found in S. cerevisiae. In two reported cases, the protein isoforms encoded by the alternatively spliced transcripts are targeted to different subcellular compartments: PTC7 proteins localize to the nuclear envelope and the mitochondria[6] and YCAT proteins to the mitochondria and peroxisomes.[7] In the case of YlMDH2, AS also leads to dual localization. We clearly demonstrated by cDNA sequencing that two transcript variants exist, which result from the use of two alternative 3′-splice sites separated by 4 nt. The predicted protein isoforms differ by only two C-terminal amino acids (AN or AKI), the latter of which constitutes a typical PTS1. Visualization of the eYFP-labelled cyto and pero isoforms showed that both proteins localized to their respective predicted compartments.

In S. cerevisiae, the distinct phenotypes of mdh2- and mdh3-deleted strains have been extensively studied. The mdh2 mutant, which is equivalent to the Y. lipolytica pero strain (deleted for the cytosolic form), exhibits no growth defect on medium with glucose as a carbon source but is unable to grow on acetate or ethanol in minimum medium.[17-19,46,47] We did not observe a similar phenotype for the Y. lipolytica pero strain. In S. cerevisiae, a functional glyoxylate cycle is mandatory for growth on C2 carbon sources and the cytosolic oxaloacetate-malate conversion step by Mdh2p is likely to be required for the glyoxylate cycle.[16] Our finding that a Y. lipolytica strain deleted for cytosolic MDH is able to grow on ethanol and acetate suggests that the peroxisomal form can perform the malate–oxaloacetate conversion of the glyoxylate cycle. If we assume that the Y. lipolytica glyoxylate cycle follows the S. cerevisiae model, this result has two implications. First, the cytosolic and peroxisomal MDH proteins in Y. lipolytica must share the same enzymatic activities, which is entirely likely as the two isoforms differ by only two C-terminal amino acids. Moreover, in S. cerevisiae, cytosolic versions of Mdh3p (lacking the PTS1) or Mdh1p (lacking the mitochondrial targeting sequence) are able to complement an mdh2 deletion for growth on ethanol, indicating that all isoforms share similar enzymatic activities despite their amino acid sequence differences.[18,50] Second, in cases of growth on C2 compounds, the peroxisomal form may be retained in the cytosol independent of the splicing event and the intermediates of the glyoxylate cycle, i.e. malate and oxaloacetate, may be able to cross the peroxisomal membrane in both directions. Examples of mislocalization of peroxisomal proteins have been reported in S. cerevisiae: Mls1p is sequestered in the cytosol when the yeast is cultivated on ethanol[51] and the Mdh3p–GFP fusion protein remained in the cytosol in an mdh2Δ mutant.[52] However, an alternative model to that of S. cerevisiae is that the entire glyoxylate cycle takes place in the peroxisome. Indeed, in Y. lipolytica, isocitrate lyase, one of the two enzymes specific to the glyoxylate cycle, is located in peroxisomes, whereas it is located into the cytosol in S. cerevisiae.[53] Thus, peroxisomal MDH may be involved in the malate–oxaloacetate conversion. However, the deletion of the peroxisomal form, as exemplified by our cyto strain, had no effect on growth on C2 compounds, suggesting that the S. cerevisiae model also reflects the situation in Y. lipolytica, at least in the absence of peroxisomal MDH.

van Roermund et al. reported that growth of the S. cerevisiae mdh3 mutant (peroxisomal form disrupted) on OA was impaired. They showed that β-oxidation is also impaired in this mutant, which suggested an indirect role for Mdh3p in this process through reoxidation of the NADH generated during the β-oxidation of fatty acids.[45] This relationship has also been demonstrated in Arabidopsis thaliana.[54] In Y. lipolytica, strains expressing only the cytosolic MDH (cyto strain) grow as well as the wild-type strain on various lipid substrates and we did not observe any other distinct phenotypes. Thus, the peroxisomal MDH is not essential for NADH reoxidation in the peroxisome. This suggests that reoxidation is carried out by an alternative pathway in peroxisomes and might take the form of shuttles based on either lactate dehydrogenase (LDH) or glycerol-3-phosphate dehydrogenase (G3PDH), which are present in mammalian peroxisomes that do not contain MDH.[55] This alternative activity remains to be discovered in Y. lipolytica.

Several less convincing alternative hypotheses can also be addressed. First, cytosolic MDH may be targeted to the peroxisomes by a PTS distinct from PTS1 or may be imported through a peroxisomal protein complex. Second, the Y. lipolytica peroxisomal membrane may be permeable to NADH metabolites, which would allow NADH recycling through a cytosolic reoxidation. However, this is clearly not the case in S. cerevisiae.[45] Third, another β-oxidation pathway allowing growth on fatty acids in the absence of peroxisomal MDH may exist in Y. lipolytica. In line with this concept, vestigial mitochondrial β-oxidation may exist in Y. lipolytica.[56,57]

From the sequencing of our cDNA libraries and the RNA-seq data, it appears that YlMDH2 splicing is probably regulated or affected by the carbon source and/or the growth conditions. Several examples of splicing regulation in yeast, principally in S. cerevisiae, are known to exist. The unspliced form of SUS1 mRNA (intron retention) accumulates after a temperature shift.[58,59] Splicing of the MER genes in response to sporulation conditions and/or meiosis has also been well documented.[3,60-62] AS autoregulation by the protein product of the gene itself, through negative or positive feedback, has been reported for RPL30,[63] YRA1[64] and SUS1.[58] In Y. lipolytica, regulation of YlMDH2 splicing could thus control the level of MDH in the peroxisome or the cytosol. Dual localization of a protein can be regulated at the level of transcription, usually involving alternative start sites, and at the post-transcriptional level. The differential targeting of YCAT is controlled at the transcriptional level by means of alternative transcription initiation sites, whose frequency of use depends on the carbon source available.[7] In the case of PTC7 (intron retention), splicing is not regulated in response to carbon source availability, but the ratio of the two protein isoforms produced by the gene depends on whether the available carbon sources are fermentable or not, indicating that PTC7 isoform production is regulated post-transcriptionally.[6]

Thus, splicing regulation, as well as dual protein localization, appears to be more prevalent in yeasts than previously believed.[65] Whether this is reminiscent of ancestral attributes or corresponds to evolutionary changes that led to selective advantages remains to be investigated. Our data illustrate distinct differences in the metabolic behaviour of Y. lipolytica and S. cerevisiae in terms of their glyoxylate and TCA pathways and their MDH enzymatic activities. This underscores a possible role in the adaptation for growth on lipid substrates and for maintaining efficient β-oxidation across evolution. The AS mechanism that regulates MDH localization in Y. lipolytica remains so far an exclusive innovation of this yeast. YlMDH2 AS regulation and its impact on protein localization and metabolism are under investigation and will provide new insights into the role of splicing on the proteome complexity in Y. lipolytica.