The transcriptome and proteome change dynamically as cells respond to environmental stress; however, prior proteomic studies reported poor correlation between mRNA and protein, rendering their relationships unclear. To address this, we combined high mass accuracy mass spectrometry with isobaric tagging to quantify dynamic changes in ∼2500 Saccharomyces cerevisiae proteins, in biological triplicate and with paired mRNA samples, as cells acclimated to high osmolarity. Surprisingly, while transcript induction correlated extremely well with protein increase, transcript reduction produced little to no change in the corresponding proteins. We constructed a mathematical model of dynamic protein changes and propose that the lack of protein reduction is explained by cell‐division arrest, while transcript reduction supports redistribution of translational machinery. Furthermore, the transient ‘burst’ of mRNA induction after stress serves to accelerate change in the corresponding protein levels. We identified several classes of post‐transcriptional regulation, but show that most of the variance in protein changes is explained by mRNA. Our results present a picture of the coordinated physiological responses at the levels of mRNA, protein, protein‐synthetic capacity, and cellular growth.
By characterizing dynamic changes in yeast protein abundance following osmotic shock, this study shows that the correlation between protein and mRNA differs for transcripts that increase versus decrease in abundance, and reveals physiological reasons for these differences.
Natural microenvironments change rapidly, and living creatures must respond quickly and efficiently to thrive within this flux. At all cellular levels—signaling, transcription, translation, metabolism, cell growth, and division—the response is dynamic and coordinated. Some aspects of this response, such as dynamic changes of the transcriptome, are well understood. But other aspects, like the response of the proteome, have remained obscured primarily because of previous limitations in technology. Without coordinated time‐course data, it has remained impossible to correctly characterize the correlations and dependencies between these two essential levels of cell biology.
This work presents an extended picture of the coordinated response of the transcriptome and proteome as cells respond to an abrupt environmental change. To assay proteomic dynamics, we developed a strategy for large‐scale, multiplexed quantitation using isobaric tags and high mass accuracy mass spectrometry. This sensitive yet efficient platform allows for the expedient collection of quantitative time‐course proteomic data at six time points, sufficiently reproducible to permit meaningful interpretation of variation across biological replicates. Time‐course transcriptome data were generated from paired biological samples, allowing us to examine the relationships between changes in mRNA and protein for each gene in terms of direction and intensity, as well as the characteristics of the temporal profiles for each gene.
It was immediately obvious that a single measure of correlation across the entire data set was a meaningless metric. We therefore analyzed relationships between mRNA and protein for different subsets of data. In response to osmotic shock, hundreds of transcripts are highly induced, and their temporal pattern reveals a transient peak of maximal induction, which resolves into a new elevated level as cells acclimate (Figure 2). For this group of genes, there is extremely high correlation between peak mRNA change and protein change (R2∼0.8). But the dynamics of the molecules differ: while mRNA levels transiently overshoot their final levels, proteins gradually rise in abundance toward their new, elevated state. We observed, however, that a measure of efficiency connects the two profiles. The time it takes for a protein to acclimate to its new state correlates with the magnitude of the excess mRNA induction. Thus, the cell imparts an urgency to protein induction by transiently producing excess transcript.
The most surprising result, however, involves transcripts that decrease in abundance. In response to osmotic shock, the cell transiently reduces over 600 transcripts, many of which are among the most highly expressed in unstressed cells. But protein levels for these genes remain, for the most part, almost completely unchanged. The stark absence of protein repression is independent of basal protein abundance, independent of reported protein half‐lives, reproducible across biological replicates, and validated by quantitative western blots. Furthermore, since we do detect a handful of proteins whose abundance is significantly reduced, our technology is capable of identifying protein loss. Thus, we conclude that transcript reduction serves another purpose besides reducing protein levels.
To explore alternate interpretations of the consequence of transcriptional repression, we devised a mass‐action kinetic model, which describes protein changes based on mRNA dynamics in the context of transient changes in the rates of cell division. The model successfully recapitulated the observed data, allowing us to alter modeling parameters to test various hypotheses.
In response to osmotic shock, overall rates of translation temporarily decrease and cell growth transiently arrests before resuming at a slower rate. We reasoned that mRNA reduction might lower the rate of new protein synthesis, but that retarded production is balanced by reduced cell division. We explored both aspects of this logic with our model.
As expected, removing cell division from our model led to a calculated decrease of protein levels, indicating that reduced growth is necessary for maintaining protein levels. However, when we computationally held mRNA levels stable and calculated protein levels in the absence of mRNA repression, we did not find the expected increase in protein abundance.
We then considered the possibility that one function of the regulated repression of these highly abundant transcripts was to liberate proteins essential for translation, such as ribosomes or translation initiation factors. To explore this, we examined a mutant lacking the Dot6p/Tod6p transcriptional repressors, which fails to properly repress ∼250 genes in response to osmotic shock. In the wild type, the mRNA for a Dot6p/Tod6p target (ARX1) decreased seven‐fold, and the remaining transcript was generally unassociated with poly‐ribosomes. In the mutant, however, the mRNA levels were reduced only two‐fold, while the remaining transcript continued to bind ribosomes. Therefore, failure to reduce transcript levels led to a persistent association with poly‐ribosomes, thereby consuming translational machinery.
Our hypothesis is, therefore, that widespread changes in the transcriptome promote efficient translation of new proteins. Transcript increase serves to increase abundance of the encoded proteins, while reduction of some of the most abundant and highly translated mRNAs supports this project by liberating translational capacity. While it is not clear what factors are the limiting elements, it is clear that a full picture of cellular biology requires exploring the dynamics of the cellular response.
The correlation between protein and mRNA change is very high at transcripts that increase in abundance, but negligible at reduced transcripts following NaCl shock.
Modeling and experimental data suggest that reducing levels of high‐abundance transcripts helps to direct translational machinery to newly made transcripts.
The transient burst of transcript increase serves to accelerate changes in protein abundance.
Post‐transcriptional regulation of protein abundance is pervasive, although most of the variance in protein change is explained by changes in mRNA abundance.
Yeast cells remodel a large fraction of their transcriptome to acclimate to the new conditions following environmental shock. Such stressful conditions trigger a common expression program called the environmental stress response (ESR), which includes increased expression of stress‐defense genes and reduced expression of protein‐synthesis and growth‐related messages (Gasch et al, 2000; Causton et al, 2001) that are synthesized at high levels during active growth (Lipson et al, 2009). Increased abundance of stress‐defense genes is important for surviving subsequent stressful insults (Berry and Gasch, 2008), while the role of transcript reduction in the ESR is less clear. It has been observed that the abundance of stress‐reduced transcripts correlates with growth rate under some (Jorgensen et al, 2002, 2004; Regenberg et al, 2006; Brauer et al, 2008) but not all (Gasch et al, 2000, 2001) conditions, leading to the hypothesis that transcript reduction contributes to resource conservation during stress defense. However, the precise function of mRNA reduction during environmental transition remains enigmatic.
Stress responses have been extensively studied at the transcript level because mRNA measurement is robust and broadly accessible, while protein analysis is considerably less well developed. This technology gap has prompted the assumption that changing transcripts mediate proportional protein alterations; however, most recent proteomic studies report poor correlation. For example, many studies comparing absolute abundance (Gygi et al, 1999; Ideker et al, 2001a; Ghaemmaghami et al, 2003; Greenbaum et al, 2003; Washburn et al, 2003) or abundance changes (Ideker et al, 2001a; Griffin et al, 2002; Li et al, 2003; Washburn et al, 2003; de Godoy et al, 2008; Soufi et al, 2009; Fournier et al, 2010) of protein versus mRNA reported only modest correlation between the two. Several smaller investigations cited higher agreement (R2∼0.7), but these were limited to a few hundred mRNA–protein pairs (Futcher et al, 1999; Lu et al, 2007). Many of these studies did not collect mRNA from the same cells from which proteins were measured, and all but a few (Li et al, 2003; Picotti et al, 2009; Soufi et al, 2009; Fournier et al, 2010) neglected temporal changes or did not perform biological replicates. Hence, the true relationship between mRNA and protein levels remains an open question.
To understand the dynamic relationship between transcript and protein, we developed a strategy for large‐scale, multiplexed quantitation by way of isobaric tags and high mass accuracy mass spectrometry (<5 p.p.m.) (Figure 1A). This platform allowed for the expedient collection of time‐course data (response to 0.7 M sodium chloride (NaCl) at six time points over 4 h) at the protein level (two technical replicates of each biological sample in ∼5 days instrument analysis). In contrast to metabolic labeling approaches (Jiang and English, 2002; Ong et al, 2002), which permit up to three‐way sample comparisons and require specific growth conditions, this multiplexed method enables practical acquisition of biological data over many samples and in any conditions. Here, we reveal the quantitative and temporal relationships between changes in mRNA and protein levels in cells acclimating to a sudden change in osmolarity.
We followed the response of actively growing Saccharomyces cerevisiae to an osmotic shock of 0.7 M NaCl. This dose of salt provides a robust physiological response but results in high viability and eventual resumption of cell growth. Samples were collected before and at 30, 60, 90, 120, and 240 min after NaCl treatment (measuring the peak transcript changes that occurs at or after 30 min (Berry and Gasch, 2008)), in biological triplicate time courses that captured cells acclimated to both environments and their transition between states. After lysing cells harvested from each time point, we digested the proteins with trypsin, generating peptides to be labeled with one of the six isobaric tags. Tagged samples were then pooled and fractionated via strong‐cation exchange (SCX) for LC–MS/MS analysis on an LTQ Orbitrap Velos mass spectrometer (Figure 1). Performing our experiment in biological triplicate generated a total of 454 755 peptide–spectral matches (PSMs), 35 828 unique peptides, and 2965 proteins (1% false discovery rate (FDR); see Materials and methods). We wrote custom software, TagQuant (Wenger et al, 2011), to extract reporter ion intensities and exclude tandem mass spectra containing interference resulting from cofragmentation of multiple precursors. Removal of precursors having ⩾25% interference greatly improved quantitative accuracy, precision, and dynamic range. To obtain maximal noise reduction across all time courses and biological replicates, we used PSMs of unambiguous provenance and required at least two unique peptides per protein. This approach is much more conservative than most proteomic analyses but provides maximal accuracy in peptide quantitation. Following this conservative analysis, we confidently measured the relative abundances of 35 000 unique peptides mapping to 2451 proteins with 60% overlap across biological replicates (Figure 1B). Of the 1814 proteins quantified in at least biological duplicate, 780 (43%) showed statistically significant changes in abundance (5% FDR, modified t‐test (Storey and Tibshirani, 2003; Smyth, 2004)). Quantitative western blotting validated these large‐scale results for several selected proteins (Supplementary Figure S1).
Relationships between changing mRNA and proteins
We next generated time‐course transcriptome data from the same biological samples and compared maximal changes in mRNA with changes in corresponding protein abundance during NaCl acclimation (Figure 2). For transcripts that increase in abundance, the log–log linear correlation between mRNA change and protein change was significantly higher than previously reported (R2=0.77), indicating that on a global scale nearly 80% of the variance in changing protein abundance is explained by increases in mRNA (Figure 2B and C). Variation in this correlation revealed that some groups of functionally related genes showed even higher agreement while others were unrelated (Supplementary Figure S2; Supplementary Table S1). This result expands on previous studies demonstrating that the correlation between mRNA and protein abundance varies by gene functional group (Washburn et al, 2003). Many of the proteins with the largest changes in abundance function in processes known to be important for NaCl survival, including glycerol and trehalose metabolism and stress defense (Supplementary Figure S3).
In stark contrast to the high correlation between protein abundance and increased mRNAs, the poor correlation between protein changes and transcripts that decrease in abundance (R2=0.09; Figure 2B and C) reflects the overall lack of protein reduction over time. Although the decreased abundance of some proteins was statistically significant, the magnitude of change was far less than that of the mRNA (e.g. ∼1.1‐fold reduction of ribosomal proteins (RPs) compared with 1.5‐ to 2‐fold reduction of transcripts). Quantitative western blotting validated the lack of protein reduction (Supplementary Figure S1), and furthermore the proteomic results did not suffer from technical inability to measure protein reduction (Supplementary Figure S4). Importantly, the observation was true irrespective of protein abundance or half‐life as measured under standard conditions (Supplementary Figure S5). The lack of reduction over this time frame is predicted for long‐lived proteins; however, this was not expected for the many proteins with short half‐lives (Belle et al, 2006). While the measured changes in mRNA abundance may be influenced by untranslated mRNAs, prior evidence suggests that the majority of most transcripts are associated with ribosomes (Arava et al, 2003), indicating that our measurements largely reflect the translated pool of mRNAs in the cell. Regardless, we conclude that the NaCl‐activated reduction in transcript abundance serves another role besides mediating changes in protein abundance.
Dynamic modeling suggests alternate roles of transcript reduction
We devised a mass‐action kinetic model to describe dynamic protein changes, in which protein abundance is a function of new synthesis and of disappearance through degradation and cell division. To calculate translation rates for each protein, we used global measurements of basal protein abundances and protein half‐lives as reported by Ghaemmaghami et al (2003) and Belle et al (2006), respectively, along with our own measurements of mRNA and growth rate changes (see Materials and methods). Using this model, we predicted changes in protein abundance based on observed changes in mRNA levels. The model's successful performance in both training and test data sets (Supplementary Appendix) allowed us to use it to test various hypotheses through simulations.
Since cell growth was arrested for ∼45 min after NaCl addition before resuming at roughly half the initial rate (Supplementary Figure S6), we first used the model to test the effect of cell‐division arrest. The model suggested that the transient growth arrest maintained protein levels despite transcript reduction. As expected, when the model assumed a constant growth rate, it predicted an ∼1.5‐fold decrease for virtually all proteins from reduced transcripts (Figure 3, blue trace).
We next tested the consequences of transcript reduction. We reasoned that, if reduced transcript abundance is critical for reducing synthesis of these proteins, then simulating transcripts without a decrease in abundance would lead to increased protein levels during growth arrest. Surprisingly, however, the model predicted little increase in protein abundance in the absence of transcript reduction; most proteins were calculated with only ∼1.2‐fold weaker reduction than that predicted from measured mRNA levels (Figure 3, purple trace). Indeed, polysome analysis showed that the maximal reduction in translation initiation (represented by the increase in monosome versus polysome complexes) peaked at 5 min after NaCl exposure (Figure 4A), whereas transcript reduction did not occur until 30 min after treatment. These results strongly suggest that the translational repression, well known to occur during the NaCl response (Uesono and Toh, 2002; Melamed et al, 2008; Warringer et al, 2010), is independent of transcript reduction and is counteracted by transient cell‐division arrest, such that corresponding protein levels do not change appreciably.
Although the reduced levels of these transcripts did not influence abundance of the encoded proteins, it may be critical for proper protein levels in the cellular system if translational capacity is limited. To explore this, we characterized the aspects of the NaCl response in a mutant lacking the transcriptional repressors Dot6p and Tod6p (Lippman and Broach, 2009; Zhu et al, 2009). In the absence of stress, the dot6Δtod6Δ mutant grew indistinguishably from wild‐type cells (Supplementary Figure S6C) and displayed similar polysome profiles (Figure 4B). However, in response to NaCl treatment, the mutant failed to properly repress ∼250 genes in the yeast ESR (Supplementary Dataset S3). Over 90% of these genes (P=10−66, hypergeometric distribution) contain upstream Dot6p/Tod6p binding elements (GATGAG; Hughes et al, 2008; Zhu et al, 2009), consistent with direct repression by the proteins. The mutant also resumed growth at a slower rate after NaCl treatment (Supplementary Figure S6C), and showed delayed resumption of normal translation profiles (Figure 4B), indicating a specific defect in acclimating to the stress.
We followed the polysome association of transcripts encoded by Dot6pTod6p targets, focusing on ARX1 (Figure 4). In unstressed wild‐type cells, ARX1 transcript was associated with polysomes, as expected (Arava et al, 2003). We observed an ∼7‐fold reduction in ARX1 levels 30 min after NaCl treatment—most of the remaining transcript was found in the monosome peak, suggesting reduced or stalled translation initiation. In contrast, the mutant showed only a 2‐fold reduction in ARX1 levels at 30 min after shock. As in wild‐type cells, there was a substantial increase in monosome‐bound ARX1 mRNA after NaCl exposure, reflecting regulated translation initiation. Surprisingly, however, a large fraction of the remaining mRNA was associated with polysomes at 30 min. We obtained virtually the same result for another Dot6p/Tod6p target, NOP2, which also showed a repression defect in the mutant (data not shown). In contrast, polysome profiles were indistinguishable at induced transcript HSP104 (Figure 4D), although our analysis would miss subtle differences. Unfortunately, we were unable to quantify changes in the corresponding proteins using several available antibodies. Nonetheless, these results show that, although the mutant dramatically reduced global translation initiation immediately after stress, failure to repress high‐abundance transcripts led to their continued polysome association as translation was resuming.
mRNA dynamics affect protein acclimation time
Seventy‐four percent of transcripts with increased abundance after NaCl treatment showed a transient ‘burst’ of change before acclimating to final levels (Figures 2D and 5A), consistent with prior studies (Gasch et al, 2000). The majority peaked at 30 min, coincident with the maximal reduction in reduced transcripts. In contrast, only 15% of proteins showed transient change, while most gradually adjusted to final levels. As expected, there was a delay between mRNA changes and protein adjustments; however, we observed a wide range of protein acclimation times, even for those encoded by transcripts peaking at 30 min.
We found that the degree of transient mRNA induction largely determines the time to protein acclimation. Transcripts with the greatest ‘overshoot,’ compared with their final, steady‐state level, produced proteins that acclimated quicker than average (P=8.5 × 10−7 for top quartile, t‐test), whereas transcripts without the transient burst led to delayed protein acclimation (P=2 × 10−3 for bottom quartile; Figure 5A and B). Our mathematical model confirmed these results—simulated transcripts lacking the transient induction burst produced slower protein changes (Figure 5C) that were predicted to acclimate 53 min later on average than proteins from transiently altered mRNAs (Figure 5D). This confirms our original hypothesis (Gasch et al, 2000) that gene expression dynamics follow a second‐order process in which the transient burst serves as a ‘loading dose’ to accelerate protein change before mRNA adjusts to ‘maintenance’ levels at the new growth state. Transcripts with the largest transient increase were enriched for those in the yeast ESR (P=6 × 10−8); however, we found no relationship between the degree of transience and the importance of the protein for NaCl survival (see Supplementary Figure S3).
Pervasive post‐transcriptional regulation
We next sought to define the impact of post‐transcriptional regulation (PTR). We identified five distinct protein classes, together implicating >40% the proteome as affected by PTR (Supplementary Table S2). The first class contained 85 proteins (11% of 791 significantly altered proteins) whose abundance changed in the absence of significant mRNA changes. Another 70 (9%) changed in the opposite direction of their mRNAs—this group was enriched for protein folding chaperones (FDR=0.035), which are known to be regulated at the level of mRNA export (Saavedra et al, 1996). A further 34 proteins changed with greater magnitude than their underlying mRNAs in at least two replicate time courses, unlike most other proteins whose abundance change was well below the mRNA change. This group was enriched for transcripts that are translationally upregulated in response to stress (P=0.016) (Law et al, 2005; Melamed et al, 2008; Halbeisen and Gerber, 2009) and includes several genes encoding signaling proteins.
Capitalizing on biological replicate analysis and the precision of our measurements, we identified a fourth mode of regulation—post‐transcriptional noise reduction that buffered protein change against variation in mRNA levels (Figure 6A and B). We identified 138 proteins (24%) whose biological variability was below the corresponding mRNAs (95th confidence level, see Materials and methods). These were heavily enriched for several functional categories, including stress‐defense proteins, translation factors, and RPs (Supplementary Table S3). Indeed, several RPs are regulated through protein degradation in proportion to subunit stoichiometry (Warner et al, 1985; Tsay et al, 1988), a mode that likely extends to most if not all RPs. Newman et al (2006) previously showed that RPs and translational regulators display low cell‐to‐cell variation within actively growing yeast cultures, likely due to PTR rather than tight transcriptional control. In contrast to that study, which found high within‐culture ‘noise’ for stress‐defense proteins, our results reveal that these proteins are targeted for noise reduction under the stressful conditions used here (P=4 × 10−5, hypergeometric distribution). Thus, the level of ‘noise’ in protein abundance is condition specific. Our conservative estimates indicate that post‐transcriptional noise reduction affects at least a quarter of yeast proteins. The prevalence of this mode of regulation raises major implications for gene expression studies, as it may dampen the effects of transcript variation observed within and across individuals.
Our results highlight the importance of quantitative, dynamic, and replicated measurements in understanding the true relationship between changing mRNA and protein. The high mass accuracy mass spectrometry used here in combination with isobaric tags allowed us to make a number of observations: (1) The magnitude of transcript induction is significantly more predictive of changes in protein level than previously thought, at least during the response to NaCl. (2) Transcript reduction did not reduce levels of the corresponding proteins; this result has not been reported before and the breadth of the effect came as a surprise. (3) The temporal relationship between transcript and protein changes is not linear but rather determined by a more complicated function. (4) There is evidence for pervasive PTR, even though the magnitude of that PTR is small. The biological implications of these observations are discussed below.
The correlation between increased transcripts and their encoded proteins during NaCl acclimation is significantly higher than previously measured, due not only to increased measurement precision but also to dynamic considerations that allowed comparison of the appropriate time points. The high protein–mRNA correlation for increased transcripts confirms that the purpose of transcript elevation is to modulate protein abundance. Although we find evidence for pervasive PTR, most of the variance in observed protein increase remains explained by mRNA changes (Figure 6C). Notably, however, nearly 20% of all changing proteins do not correlate with underlying mRNA changes. Thus, while increases in mRNA are a good predictor of protein change under the conditions studied here, the absence of mRNA induction is not necessarily conclusive.
In stark contrast, transcript reduction did not alter protein levels under these conditions. This result has not been widely reported previously, although analysis of available data sets generally shows fewer and smaller changes in reduced proteins compared with those increasing in abundance, across several different environmental comparisons (Blomberg, 1995; Ideker et al, 2001b; Griffin et al, 2002; Li et al, 2003; Washburn et al, 2003; Picotti et al, 2009; Soufi et al, 2009; Fournier et al, 2010). However, Fournier et al (2010) observed a decrease in RPs after rapamycin treatment that was significantly delayed (∼6 h) compared with the reduction in RP transcripts. In that case, the late protein decrease may be due to ribosome consumption through ribophagy (Kraft et al, 2008), rather than translational repression. The lack of RP reduction soon after rapamycin‐dependent transcript reduction is therefore consistent with our results, and suggests that the lack of correlation between transcript reduction and protein abundance may be common to other stress conditions. This raises major implications for the interpretation of transcriptomic data, since reduced mRNA abundance is often used to infer protein reduction and dispensability.
Reduced levels of abundant transcripts, particularly those affecting ribosome biogenesis, has been proposed to reduce the costly synthesis of encoded proteins in proportion to growth demand (Waldron and Lacroute, 1975; Kief and Warner, 1981; Nomura, 1999; Jorgensen et al, 2002, 2004; Rudra and Warner, 2004; Regenberg et al, 2006; Brauer et al, 2008). Indeed, our mathematical modeling revealed that changes in growth rate explain the lack of protein change, not just for long‐lived RPs but for all proteins from reduced transcripts. However, our results strongly suggest that transcript reduction is not the driving force behind decreased protein synthesis. Reduced translation initiation (inferred from the increased monosome:polysome ratio) occurs much before transcripts drop in abundance (Figure 4B; Uesono and Toh, 2002). We also observed that final levels of most RP transcripts were at or above unstressed levels, even though division occurred at half the initial rate. This observation underscores the discordance between RP transcript abundance and growth rate during osmotic stress.
Instead, the reduced abundance of these transcripts may serve to redirect translational capacity to newly made mRNAs. The timing of transcript reduction occurs when translation profiles begin to recover as cells resume growth and coincides with maximal levels of increasing transcripts. The temporal relationship between increasing and decreasing transcripts holds across many environmental shocks, even though the kinetics of transcript change is not well correlated with the duration of cell‐division arrest during these transitions (Gasch et al, 2000, 2001). We suggest that the coordination of transcript reduction with transcript induction avoids competition for translational machinery, simply through the temporary removal of high‐abundance transcripts when translation is resuming.
In this model, cells may be limited for translation factors, particularly initiation factors that regulate translation in response to stress (Proud, 2007; Gandin et al, 2008; Park et al, 2011), but another possibility is competition for translating ribosomes. Nearly 90% of ribosomes in growing cells are actively translating proteins (Warner, 1999; Arava et al, 2003; von der Haar, 2008), leaving little capacity to synthesize new proteins during adversity. Under this scenario, the level of transcript reduction might be tuned to ribosome demand at increased mRNAs. A simple calculation (Figure 7) based on the number of ribosome per cellular transcript (Arava et al, 2003), transcript counts per cell (Lipson et al, 2009), and our own measurements of mRNA fold‐change suggests that at 30 min after NaCl treatment, the fraction of ribosomes available solely due to transcript reduction (24±5% of translating ribosomes) is approximately the estimated maximum needed to translate the increased transcripts (32±6% translating ribosomes). Although our estimate does not capture known translational regulation at specific transcripts (Uesono and Toh, 2002; Law et al, 2005; Melamed et al, 2008; Halbeisen and Gerber, 2009; Warringer et al, 2010), it supports the hypothesis that transient transcript reduction is linked to translational capacity. Our polysome analysis is consistent with this model, since failure to repress transcripts in the dot6Δtod6Δ mutant led to continued polysome association of aberrantly abundant mRNAs. These results support a recent theoretical study by Scott et al (2010), suggesting that the relationship between continuous growth rate and nutrient quality in Escherichia coli is due to ribosome allocation.
Our model is also consistent with the observed dynamics of gene expression change. Through direct observations and mathematical modeling, we demonstrate that the transient burst in transcript abundance serves to accelerate protein change. It is likely not a coincidence that this response peaks as resuming translational machinery is available to new transcripts. The mRNA response follows a second‐order dynamic system, resulting in a pulse of gene expression change before cells acclimate. This pulse is likely due in part to the interplay between transient alterations in transcription and mRNA stability during stress response and acclimation (Garcia‐Martinez et al, 2007; Shalem et al, 2008; Molin et al, 2009). In addition, feed‐forward signaling loops can also produce pulse‐like responses (Mangan and Alon, 2003; Alon, 2007; Kaplan et al, 2008), and such signaling motifs likely contribute here as well (Gasch, 2002a). The outcome of transient transcript changes is rapid alteration in protein changes, with larger mRNA pulses producing faster protein acclimation. More broadly, the complex relationships between these processes highlight the importance of considering dynamic observations in the study of living systems.
Materials and methods
Cell growth and RNA preparation
Strain BY4741 was grown >7 generations to an optical density (OD600) ∼0.3 at 30°C in YPD medium. A sample of unstressed cells was removed as the ‘time 0’ reference and pre‐warmed medium was added for a final concentration of 0.7 M NaCl. Samples were removed at 30, 60, 90, 120, and 240 min after addition of NaCl (the culture was diluted 1:1 with fresh medium after the 120‐min collection to maintain log‐phase growth). The majority of transcript changes peak at or after 30 min (Berry and Gasch, 2008; data not shown). At each time point, one sample was recovered for microarray analysis and another was recovered for proteomic analysis by 3 min room temperature spin at 3500 r.p.m. after which time the cells were flash frozen in liquid nitrogen and maintained at −80°C until use.
The dot6Δtod6Δ mutant was constructed by replacing TOD6 with the hygromycin‐MX cassette in the BY4741‐dot6Δ∷KANMX strain (Open Biosytems), through homologous recombination. Both deletions were verified by diagnostic PCRs.
Proteomic lysis and digestion
Cells were lysed by three passages through the French press at 4°C in 3 ml of lysis buffer consisting of 50 mM Tris pH 8, 4 M urea, 75 mM NaCl, 1 mM DTT, complete mini EDTA‐free protease inhibitor (Roche Diagnostics, Indianapolis, IN), and phosSTOP phosphatase inhibitor (Roche Diagnostics). The lysate was centrifuged at 14 000 r.p.m. for 10 min. Cysteine residues were reduced and alkylated by incubating lysate with 3 mM DTT (final concentration) for 45 min at 37°C followed by incubation in 15 mM IAA for 1 h at room temperature in the dark. The alkylation reaction was capped by incubating the reaction with DTT for 15 min at room temperature. Proteins from each time point were digested overnight at 37°C after the addition of 1 mM CaCl2, 50 mM Tris (to decrease urea to 1 M) and adjusting to pH 8 at an enzyme:substrate ratio of 1:50 of trypsin (Promega, Madison, WI). Each digest was quenched by the addition of TFA to a final concentration of 0.5% (pH ⩽2), desalted via solid phase extraction on a 50‐mg tC18 SepPak cartridge (Waters, Milford, MA), and the eluent was lyophilized.
Each TMT (Tandem Mass Tag) was resuspended in 41 μl of pure ACN. Lyophilized peptides were resuspended in 18 MΩ H2O to a final concentration of 5 μl/μl, with the addition of 200 mM TEAB bringing the final volume of 132.5 μl. Resuspended peptides were mixed with each respective tag as follows: 30 min, TMT126; 60 min, TMT127; 90 min, TMT128; 120 min, TMT129; 0 min, TMT130; 240 min, TMT131. Each reaction was allowed to incubate at room temperature with intermittent mixing for 1 h. Tagging was quenched with the addition of 8 μl 5% hydroxylamine for 15 min. The differentially labeled TMT samples were pooled in equal amounts and dried down.
Labeled peptides were fractionated using SCX chromatography as described previously (Swaney et al, 2010). Each TMT 6‐plex peptide mixture was resuspended in 500 μl of SCX buffer A (5 mM KH2PO4, 30% ACN, pH 2.65) and loaded onto a polysulfoethyl apartamide column (9.4 × 200 mm2, Poly LC, Columbia, MD) connected to a Surveyor LC quaternary pump (Thermo Electron, San Jose, CA) running at 3.0 ml/min. Peptides were detected via a PDA detector (Thermo Electron). A total of 16 fractions were collected in 3 min intervals over the course of the SCX gradient. The following gradient was used: 2 min of isocratic buffer A, followed by a linear gradient of 0–10% buffer B (5 mM KH2PO4, 30% ACN, 350 mM KCl, pH 2.65) from 5 to 35 min. Buffer B was ramped up to 100% over the next 6 min. Then, a 7‐min transition from buffer B to 100% C (50 mM KH2PO4, 500 mM KCl, pH 8) and buffer D (18 MΩ water) were used to wash the column. Each fraction was lyophilized and desalted on 50 mg tC18 SepPak cartridges. Desalted eluants were lyophilized and stored at −80°C.
SCX fractions were resuspended in 30 μl of 0.2% formic acid. A Waters nanoAquity HPLC and auto‐sampler was used to load samples onto an 8‐cm, 75 μm i.d. pre‐column packed with 5 μm C18 particles, and separated across a 25‐cm, 50 μm i.d. analytical column packed with 5 μm C18 particles. Samples were loaded onto the pre‐column in 98% A (0.2% formic acid in water), 2% B (0.2% formic acid in CAN). Buffer B was increased to 5% over the first 3 min of the separation followed by an increase in B by 0.25% per min over 120 min, followed by a quick ramp up to 70% B over 10 min and held for 5 min. The gradient was dropped back to 2% over a period of 5 min and allowed to re‐equilibrate for 20 min. Eluent was detected on an LTQ Orbitrap Velos (Thermo Fisher Scientific, Bremen, Germany) via an integrated electrospray emitter operating at 2.2 kV. Experiments consist of MS1 analysis in the obitrap mass analyzer followed by 10 data‐dependent high energy collision cell dissociation MS/MS events with a precursor isolation width of 2 Th and mass analysis in the orbitrap as well. Each fraction was run twice. For all experiments, an AGC target value of 1 000 000 charges was used for MS1 and 45 000 charges for MS2. A mass resolution of 30 000 was used for MS1 and 7500 for MS2. Precursors were dynamically excluded for 45 s, and only peptides with assigned charge states of two or greater were selected for MS/MS interrogation. All proteomic data are available on our website (http://www.chem.wisc.edu/~coon/Downloads/Lee_MSB_2011/) and in the Proteome Commons Tranche repository under accession hash 8P18XwqHFqQk+NZ5TJj4DR0e5qLcOYvmhFfHwlyNaK9SWrkO93CYbKok485Hw/ud8Lz/5 × CjvVwNf5yVCMclx5zNvc0AAAAAAABuKg==
Proteomic data analysis
Spectral reduction was performed by DTA Generator. The OMSSA (Open Mass Spectrometry Search Algorithm) (Geer et al, 2004) was used to search spectra against a concatenated target‐decoy version of the Saccharomyces Genome Database (http://www.yeastgenome.org) with fully tryptic enzyme specificity, allowing up to three missed cleavages. Static modifications of carbamidomethylation (+57 Da) on cysteine residues, TMT 6‐plex on N‐termini, and TMT 6‐plex on lysine residues, and variable modifications of oxidation of methionine residues (+16 Da) and TMT 6‐plex on tyrosine residues were specified. An average mass tolerance of ±4.5 Da was used for precursors, while a monoisotopic mass tolerance of ±0.01 Da was used for fragments with orbitrap detection. Identifications were filtered to a FDR of 1% using in‐house software that iteratively checked combinations of expectation value (e‐value) score and precursor mass error to find thresholds that maximize unique peptide identifications. Proteins were reduced for parsimony and filtered to 1% FDR as well. Protein quantitation was evaluated with custom software that corrects for isotopic impurities, normalizes reporter ion intensities, and coalesces peptide quantitation into protein quantitation (Wenger et al, 2011). Peptide and protein summaries are available in Supplementary Datasets S1 and S2, respectively. Reproducibility across biological replicates is summarized in Supplementary Table S4.
Microarray sample preparation and data acquisition
Total RNA was recovered by hot phenol extraction as previously described (Gasch, 2002b) and purified with a Qiagen RNeasy column. cDNA synthesis and microarray labeling was performed as described (Berry and Gasch, 2008), except that total RNA was labeled with a mixture of oligo‐dT and random hexamer at a 1.7:1 molar ratio. Labeled cDNA from each RNA sample was hybridized against genomic DNA labeled as previously described (Pollack et al, 1999), swapping dye orientations for each time course. Samples were hybridized to custom Nimblegen tiled arrays and incubated on a Maui Hybridization Chamber at 42°C for 16 h. Arrays were scanned and analyzed with a GenePix4000 scanner (Molecular Devices, Sunnyvale, CA), and signal from both channels was extracted with the program NimbleScan. Signal intensity from the RNA‐sample channel was quantile normalized for all 18 arrays using the Bioconductor package quantile.normalize() (Bolstad et al, 2003). The log2 ratio was calculated for each gene from the median intensity of gene probes at each time point compared with that of the unstressed cells.
The response of wild‐type and dot6Δtod6Δ cells to 0.7 M NaCl was conducted in biological duplicate, as described above, comparing transcript abundance at 30 min after NaCl addition to the unstressed sample from each strain. Relative basal transcript abundance in wild‐type and dot6Δtod6Δ cells was performed by extracting the signal intensity corresponding to the unstressed samples and comparing across arrays as described for the NaCl time course. Genes affected by DOT6TOD6 deletion were identified as those whose expression was >1.5‐fold different between strains in both replicates: of the 323 genes identified, 244 displayed a repression defect in the dot6Δtod6Δ mutant. Nearly all of these showed ∼1.5 × higher basal levels in the unstressed mutant compared with wild type. All microarray data are available under GEO accession GSE23798 and in Supplementary Dataset S3.
Quantitative western blotting
Cells were lysed in cracking buffer and incubated at 95°C for 5 min with glass beads. Lysate from equivalent number of cells was analyzed by western analysis, using a mixture of mouse anti‐actin (MAB1501, Millipore, Billerica, MA) and a second rabbit polyclonal antibody (Hsp104 antibody (Novus Biologicals, NB120‐20547), and gifts of the Elizabeth Craig laboratory) to simultaneously detect Act1 and a second protein of interest. IRDDye680‐conjugated anti‐rabbit IgG and IRDye800‐conjugated anti‐mouse IgG (LiCor, Lincoln, NE) were used as secondary antibodies at 1:20 000 dilutions in Odyssey blocking buffer. Signal intensities were analyzed by using the Odyssey infrared image system (LiCor).
Polysome profiling and qPCR
Polysome profiles were collected as described in Arava et al (2003) using a 5–50% continuous sucrose gradients analyzed on a Teledyne Isco Fractionation unit. Fractions of ∼60 μl were collected and pooled to capture peaks indicated in Figure 4. Synthetic control RNA (Agilent, Santa Clara, CA) was added to each pool as a normalization control at 15 × 106 copies and RNA from the pool was purified over an RNeasy column (Qiagen, Valencia, CA). Quantitative PCR was done in at least technical duplicate using iQSYBR Green Supermix (Bio‐Rad, Hercules, CA) on a MyiQ2 Bio‐Rad Cycler. Primers spanned a 3′ 100–200 bp region of each ORF. Cycle numbers were normalized to the doped synthetic RNA (Figure 4) or to SEC21 as controls (Supplementary Figure S7). Relative log2 mRNA abundance in each polysome fraction was normalized to mRNA abundance measured in the trough between polysome peaks 1 and 2, to adjust to baseline levels.
Significant mRNA and protein changes were determined for all molecules measured in at least biological duplicate using an ANOVA test in the limma Bioconductor package (Smyth, 2004) with a Q‐value (Storey and Tibshirani, 2003) FDR <0.05 taken as significant (Supplementary Dataset S3). In all cases, correlations were assessed based on the average values for fully triplicated data to minimize effects of technical noise. Transcripts that increased or decreased in abundance were defined based on the direction of the largest absolute change in abundance—in most cases, increases or decreases in abundance were largely consistent across all time points. Molecules subject to transient abundance changes were determined in limma by comparing each time point to the final 240 min sample, with FDR <0.05 taken as evidence of transience. Timing of peak levels for transient changes and of acclimation for gradual changes were also identified by a t‐test in limma comparing adjacent time points, taking the earliest statistically differentiated time point as the peak/acclimation time. Times were manually inspected and corrected in some cases where subtle but consistent differences were not scored correctly. Transcripts whose change in expression at 15 min was greater than the change at 30 min (Berry and Gasch, 2008; unpublished data) were removed from this analysis (amounting to 13 of 140 genes for Figure 5). The magnitude of the transient burst was measured for RNAs peaking at 30 min by calculating the average difference in changes at 30 min versus 240 min.
Estimates of ribosome redistribution were calculated as follows. Measures of mRNA abundance per cellular mRNA pool from Lipson et al (2009) were converted to copies per cell, based on 60 000 transcripts per actively growing yeast cell (Zenklusen et al, 2008). The number of ribosomes per transcript, as measured by Arava et al (2003), was multiplied by the number of transcripts for each gene; the total number of engaged ribosomes was scaled to 171 000 (90% of 190 000 ribosomes; Warner et al, 2001). To estimate the number of ribosomes released at each time point in each time course, we multiplied the number of ribosomes per gene sequence by the fraction of mRNA that disappeared or appeared relative to starting levels. At each time point, all transcripts for which the average log2 value was <0 were taken as reduced, and the inverse was true for increased transcripts. The average and s.d. of the number of ribosomes released across all three time courses is shown in Figure 7.
Proteins subject to PTR were defined as follows: proteins whose maximal change was in the opposite direction of the maximal mRNA change were defined as opposing mRNA changes. Proteins changing greater than mRNA were identified if their steady‐state differences at 240 min were greater than the corresponding mRNA differences in at least two of the three time courses. This is in contrast to most proteins that reach <50% the fold‐change in mRNA (data not shown). Proteins subject to noise reduction were identified as follows: the Pearson correlation was calculated for each replicate profile for mRNA or protein, for all triplicated pairs where both mRNA and protein changed significantly. The average and s.d. of the three pairwise correlations was identified; mRNAs with an average correlation outside the average protein correlation plus two s.d. (95th confidence level) were identified. Functional enrichment was determined by testing the hypergenomic distribution using FDR <0.05 (Storey and Tibshirani, 2003) as significant in all cases.
Mathematical modeling overview
A detailed discussion of the modeling and results is presented in Supplementary Appendix. We developed a mass‐action model in which:
where kt,ρ represents the translation rate for protein ρ at time t, and kd,ρ represents the protein degradation rate (fixed over time) as measured by Belle et al (2006), μ represents the dilution rate due to cell division, and represents the concentration of the protein at time t. Experimental evidence suggested little change in cell volume during the experiment (Supplementary Figure S6B), allowing us to assume a constant cell volume.
We suspected that the translation rate changed with time, particularly during cell‐division arrest (Uesono and Toh, 2002), requiring us to create a model that allowed the ‘constants’ in the above equation to vary. These considerations changed our formula to the following form:
where Pρ is the amount of protein ρ per cell, ks,ρ(t) is the translation rate for protein ρ at time t, mρ (t) is the mRNA amount for the transcript encoding protein ρ, kd,ρ is the decay rate of protein ρ, and μ(t) is the diffusion of proteins caused by cell growth at time t.
Because we made measurements at specific, discrete time points, we described time as having specific values (although they need not be uniformly distributed, as in many simulations): t0=0; t1=30, t2=60, t3=90, t4=120, and t5=240 min. We assumed for simplicity that translation rate and cell growth rate change instantly and only at the time points at which we have taken measurements (namely, at 30 min after NaCl treatment). We tested three models of translation. The first (model 1) calculated a single ks,ρ based on the t0=0 to t1=30 min transition. The second (model 2) calculated a single ks,ρ based on a regression across all time points. Because we expected that translation rates may vary during and after cell‐division arrest (Uesono and Toh, 2002), we calculated a third model (model 3) in which two ks,ρ values were used: ks,ρ0 estimated translation rates between 0–30 min and ks,ρ1 estimated rates based on a regression across the 60–240 min time points. The model was developed on replicate 1 as the training set and validated on replicates 2 and 3 as test sets. In all cases, model 3 outperformed the other models (Supplementary Table S5; Supplementary Figure S8).
Simulations described in the main text were performed using model 3 on the training set. To test the effects of dilution, protein abundance was calculated with μ=80 min throughout the time course instead of only after t=30 min. To test the effect of increased‐mRNA dynamics, the mRNA level at 30 min was replaced with the average values of the remaining time points to remove the transient burst, for all transcripts whose transcript increase peaked at 30 min. The effect of mRNA reduction was tested by calculating protein abundance when transcript levels were held at the median fold‐change across the entire transcriptome (log2∼0). Mathematical models used in the analysis (Supplementary Dataset S4) have also been made available for analysis with COPASI (Supplementary Dataset S5). COPASI code was not used for the analysis described here, but was validated for several candidate mRNAs and made available for public use.
We thank JL Will, L Hoover, and S Haroon for technical assistance; E Craig for providing antibodies and instrument access; AJ Bureta for graphical support; and R Gourse, E Craig, M Culbertson, J Yin, and the Gasch and Coon laboratories for constructive comments. SET was supported by NIH Genomic Sciences Training Grant T32HG002760. This work was funded by NIH R01 GM083989 to APG and an ARRA collaborative supplement to JJC and APG, and NIH R01GM080148 and PO1GM081629 to JJC.
Author contributions: MVL, SET, JJC, and APG designed the experiments; MVL performed proteomic analysis; SET performed transcriptomic analysis; SLH performed mathematical modeling; JH performed all polysome experiments and quantitative westerns; CDW wrote proteomic analysis software; MVL, SET, and APG performed statistical analysis and data interpretation; JJC and APG conceived of the experiments; MVL, SET, SLH, JJC, and APG wrote the manuscript.
Conflict of Interest
The authors declare that they have no conflict of interest.
Supplementary Tables S3–5, Supplementary Figures S1–8
R2 for mRNA‐protein changes by GO category
List of proteins subject to different types of post‐transcriptional regulation
Summary of metrics and quantitation data for all peptides in each of three biological replicates
Summary of metrics and quantitation data for all proteins in each of three biological replicates
Complete normalized dataset and ANOVA FDR values
EXCEL‐based modeling used in the manuscript
COPASI‐compatible modeling code and instructions
This is an open‐access article distributed under the terms of the Creative Commons Attribution License, which permits distribution, and reproduction in any medium, provided the original author and source are credited. This license does not permit commercial exploitation without specific permission.
- Copyright © 2011 EMBO and Macmillan Publishers Limited