Reduced sphingolipid biosynthesis modulates proteostasis networks to enhance longevity

As the elderly population increases, chronic, age-associated diseases are challenging healthcare systems around the world. Nutrient limitation is well known to slow the aging process and improve health. Regrettably, practicing nutrient restriction to improve health is unachievable for most people. Alternatively, pharmacological strategies are being pursued including myriocin which increases lifespan in budding yeast. Myriocin impairs sphingolipid synthesis, resulting in lowered amino acid pools which promote entry into a quiescent, long-lived state. Here we present transcriptomic data during the first 6 hours of drug treatment that improves our mechanistic understanding of the cellular response to myriocin and reveals a new role for ubiquitin in longevity. Previously we found that the methionine transporter Mup1 traffics to the plasma membrane normally in myriocin-treated cells but is not active and undergoes endocytic clearance. We now show that UBI4, a gene encoding stressed-induced ubiquitin, is vital for myriocin-enhanced lifespan. Furthermore, we show that Mup1 fused to a deubiquitinase domain impairs myriocin-enhanced longevity. Broader effects of myriocin treatment on ubiquitination are indicated by our finding of a significant increase in K63-linked ubiquitin polymers following myriocin treatment. Although proteostasis is broadly accepted as a pillar of aging, our finding that ubiquitination of an amino acid transporter promotes longevity in myriocin-treated cells is novel. Addressing the role of ubiquitination/deubiquitination in longevity has the potential to reveal new strategies and targets for promoting healthy aging.


INTRODUCTION
Aging is widely accepted as a major risk factor for many chronic diseases and resultant physiological decline leading to mortality [1]. Research on many fronts is revealing potential ways to postpone agerelated decline, maintain normal physiological function longer and improve healthspan. Some of the most promising research seeks to limit nutrient intake or increase daily fasting time as a means to improve healthspan in humans [2][3][4]. Still, these strategies will be difficult for most humans to adhere to in order to gain health benefits. Pharmacological agents offer a potential way to obtain the beneficial effects of nutrient limitation, but such compounds have yet to be identified although progress is encouraging [1,5,6].
We have identified a potential pharmacological agent, the natural product myriocin (Myr, ISP-1), that AGING increases chronological lifespan in budding yeasts (Saccharomyces cerevisiae) by more than two-fold [7]. Myr works, at least in part, by reducing the free pool of most amino acids similar to what amino acid restriction does [5,8,9]. Our interest in Myr stems from its target enzyme serine palmitoyltransferase (SPT), catalyzing the first and rate limiting step in sphingolipid biosynthesis in all eukaryotes [10][11][12]. In addition, Myr was first identified in a search for antibiotics [13] and anti-inflammatory drugs [14]. More recently, it has shown beneficial effects in treating age-associated diseases including diabetes, cancers, neurological and cardiovascular disorders [15][16][17][18][19][20] and other diseases including muscular dystrophies, cystic fibrosis and retinopathy [21][22][23].
Sphingolipids serve as both structural components of cellular membranes and as signal or regulatory molecules influencing many physiological processes, particularly in mammals [15,[24][25][26]. Because most de novo lipid biosynthesis begins in the endoplasmic reticulum and continues in the Golgi Apparatus before the terminal products are distributed to cellular membranes, Myr treatment has the potential to diminish or enhance a variety of processes. Our recent studies identified diminished processes that may foster longer lifespan: newly synthesized Mup1, the major high-affinity methionine (Met) transporter, trafficked normally to the plasma membrane (PM) but was inactive in drug-treated cells resulting in reduced Met uptake, starting after about 2 h of drug treatment, as the fraction of active Mup1 was diluted by cell growth and division. Moreover, Myr promoted endocytic clearance of Mup1, indicating that altered sphingolipid levels trigger remodeling of nutrient transporter composition at the PM [7]. Thus, post-translational effects are vital to Myr-induced down-sizing of amino acid pools.
Previously we found that Myr treatment had large, global effects on transcription after 6-7 cell doublings [27]. In the present work, we examined mRNA levels during the initial stages of Myr treatment to construct an overview of transcriptional changes with the aim of identifying how long it takes cells to respond to drug treatment and to identify novel factors critical for Myr-enhanced lifespan. We find that transcription is strongly up-regulated starting after 4 h of Myr treatment when cells progress through a second cell division cycle. In addition, transcript data suggested a novel role for ubiquitin in lifespan and targeted studies identified ubiquitination of Mup1 as essential for Myr-enhanced longevity. Our transcriptomics data provide a valuable resource to better understand how myriocin treatment enhances longevity.

Four hours of myriocin treatment induce robust transcriptional changes
To examine transcriptional changes induced by Myr treatment, we diluted stationary phase prototrophic BY4741 cells (50-60% of cells synchronized for first two cell cycles over a 6 h time frame [7]) into fresh culture medium (Time 0), with and without Myr treatment. Samples for analysis of mRNA abundance by RNA seq were taken at 1 h intervals over a 6 h time course ( Figure 1A). The normalized RNA seq data for the time course contained transcripts, expressed as transcripts per million (TPM), mapped to 6198 genes, 5169 of which were uniquely annotated and of sufficient signal intensity for subsequent analysis (Supplementary Table 1, Filter tab). Filtered data were analyzed by two-way ANOVA with a statistical cutoff of p = ≤ 0.01 to give a set of 4964 significant genes ( Figure 1B). These were further sorted into drug (D), time (T), both time and drug (T&D) and Interactions (I) ( Figure 1B, Euhler diagram), and assigned a p-value. Genes found to be significant in any of these three categories are referred to as differentially expressed (DEG).
DEGs that change up or down over time were identified by using a 'template assignment tool' (MATERIALS AND METHODS, example in Figure 2A Table 1 -Pattern Graph tab can be used to plot the pattern of any gene). Within each temporal pattern (other than '00') a gene was placed in one of two horizontal reflections ( Figure 1C, columns 'a' or 'b'). Three types of Drug effects are shown in each pattern ( Figure 1C, Y-axis: light blue, black or orange curves) along with the number of genes in each category (see Legend for Figure 1C for more details). These 11 patterns explained the majority (77%) of all DEGs. Results of DEGs in significantly enriched patterns are shown as a heatmap and a graph for Vehicle and Myriocin treatments ( Figure 1D, see Legend for details and pattern labels). A graph of any gene transcript across the 1-6 h timeframe can be plotted using the Graph Reporter tab in Supplementary Table 1. The most significant GO terms enriched in these patterns are shown in Table 1 and a complete list of genes in each GO term is presented in Supplementary Table 2 (Webgestalt patterns combined tab). These GO terms primarily represent major processes, cellular compartments or molecular functions needed to drive cell growth and division such as AGING and myriocin (Myr) treated cultures. RNA was isolated and mRNA was quantified using RNA-seq. (B) Differentially expressed genes (DEGs)-RNA-seq data were filtered to include well-annotated genes with sufficient signal intensity prior to statistical analysis. Two-way ANOVA (Time: 6 time points; Drug: two conditions-Veh and Myr) was applied. The ANOVA test produces 3 p-values for each gene, one for the main effect of time (T-black circle), one for the main effect of drug (D-red circle), and one for the 'Interaction' between Time and Drug (I-green circle; analyzed separately-see Figure 2). DEGs (p ≤ 0.01 on at least one of the three p-values) are shown. Because each DEG could be significant by 1-3 ANOVA p-values, data are displayed in a Euler diagram to indicate overlap. (C) For DEGs significant by Time (T), Drug (D), or both (T&D), post hoc analysis categorized each DEG into one of 11 different temporal patterns (00-no temporal effect; 01 and 05changed at the first or last time point; 02-04-plateaus at intermediate time points; 06-09-spikes at intermediate time points; 10-linear change with time; for each temporal pattern other than '00', a gene could be assigned to one of two horizontal reflections-'a' or 'b'), and within each temporal pattern, one of three drug effects (increased by Myr-↑ orange, decreased by Myr-↓ blue, not changed by Myr-↔ black) was defined. Numbers of genes assigned to each pattern are included on the right side of each pattern diagram ( * significantly more DEGs than expected by chance; p ≤ 0.01, binomial test). (D) Results for each asterisked pattern from C are shown. On the left, heatmaps of standardized average signal at each time point for vehicle and myriocin groups for 5 representative DEGs are shown. On the right, the graphed standardized averages (± SD) for all genes assigned to each pattern are shown (scale bar at bottom). Pattern names are given by the Euler diagram region, then the temporal pattern # and 'a' or 'b' designation, followed by drug effect (in parentheses) and the number of genes assigned. All data for the 5169 RNA seq transcripts that passed filtering metrics are shown in Supplementary Table 1 along with the two-way ANOVA analysis results which can be used to sort for genes in each sector of the Euler diagram. chromatin organization, transcription regulation, mitotic cell cycle and biochemical pathways.
In contrast to the Time x Drug sector, the 1009 genes in the Time sector of the Euler diagram fall primarily into three patterns ( Figure 1C and 1D, T01b(↔), T10A (↔) and T10b(↔)). Genes in these patterns are not significantly changed by Myr treatment but are significantly changed with time and serve as an example of time-related changes. All of these patterns will require further effort to determine their significance, if any, in Myr-induced longevity.
The Drug sector of the Euler diagram ( Figure 1C and 1D, pattern D(↑)) contains 329 up-regulated genes ( Table 1 and Supplementary Table 2). Because of our experimental design, we cannot be sure if some genes in this sector are driving the changes associated with time that we see, but they could be. The most enriched GO terms relate to chromosomes, and transcription (Table 1), indicating a strong response to Myr treatment independent of time in culture.
Next, we analyzed the 1570 genes in the Interaction group for temporal patterns responding to Myr treatment ( Figure 2A). Five patterns contained more than 100 DEGs ( Figure 2A, center panel, numeric results for this grid are available in Supplementary  Table 1, tab called 'I intersections'). The right-hand graph in Figure 2A indicates the number of up-and down-regulated genes in these 5 patterns. Five templates for these patterns ( Figure 2B, left-side graphs) along with heatmaps of representative genes in each pattern for Vehicle and Myriocin treated samples are shown in Figure 2B (center panels). The right-most panels in Figure 2B represent average Z-scored gene expression profiles for the 5 patterns examined. Both the heatmaps and the data plots on the right side of Figure 2A reveal that Myr treatment enhances transcription at the 5 and 6 hour time points. This observation was tested by examining all genes found to be significant by the two-way ANOVA's in the interaction term. For this test we plotted the log 2-fold changes of Myr treatment/vehicle treatment at each time point and color-coded the results if they were significantly down-regulated (blue) or up-regulated (orange) (p ≤ 0.05, pairwise Fisher's LSD). As predicted, the predominant effect of Myr was to upregulate genes at the 5 and 6 h time points ( Figure  2C). Functional overexpression analysis of significant DEGs at these two time points (called: I: 5 and 6(↑)) revealed a more generalized picture of the biological processes and cellular compartments responding to Myr ( Figure 2C, right panel, all genes listed in Supplementary Table 2). As we discuss in detail below, up-regulation of genes starting between the 4th and 5th hour occurs when the majority of cells enter their second cell division cycle [7]. Thus, these data explain much of the effect of time and drug components within the interaction gene set.
We recently reported that Myr treatment has a notable effect on the size of most amino acid pools which remain significantly smaller in drug-treated cells starting before the 1 h time point, and remaining smaller over the six hour time course studied [7]. To determine if transcription has roles in lowering and maintaining smaller amino acid pools, we performed functional enrichment analysis on the 1570 genes in the Interaction group ( Figure 2A) because this group captures a large fraction of genes responding differently over time with Myr treatment. We hypothesized that down-regulation of genes would imply a contribution to lower amino acid pools while up-regulation would imply the opposite effect or an attempt to restore pools to the size found in vehicle-treated cells. The KEGG pathway term 'sce01100:Metabolic pathways', containing 314 genes, is highly enriched in this gene set (p = 1,6E-27, AGING  Table 2). (C) Left: The analysis of interaction-significant DEGs, particularly in 2B, suggested a simplified and more general trend among the genes found to be significant by the two-way ANOVA's Interaction term. That is, Myr treatment appeared to increase expression at time points 5 and 6 regardless of the temporal pattern of expression. To test this, we plotted the log 2-fold changes of Myr treatment/vehicle treatment at each time point, and color-coded those results if they were significantly (p ≤ 0.05, pairwise Fisher's LSD) downregulated (blue) or upregulated (orange). As expected, the predominant effect was upregulation at time points 5 and 6. Right: Functional overrepresentation analysis. Columns-Tot: total number of genes in the dataset assigned to GO term; Sig: number of total genes that were significant and were assigned to the indicated pattern of expression; Exp: number of genes expected to be found in that pathway for that pattern of expression by chance; Ratio: Sig/Exp; p-value: overrepresentation analysis ( [46]); FDR: multiple testing adjusted ORA p-value according the BH procedure ( [47]; Note-complete list of all pathways, and the genes assigned to them, is provided in Supplementary  Columns-Abbreviations: Tot: total number of genes in the dataset assigned to GO term; Sig: number of total genes that were significant and were assigned to the indicated pattern of expression; Exp: number of genes expected to be found in that pathway for that pattern of expression by chance; Ratio: Sig/Exp; p-value: overrepresentation analysis ( [46]); FDR: multiple testing adjusted overrepresentation analysis p-value according the BH procedure ( [47]). Note-complete list of all pathways, and the genes assigned to them, is provided in Supplementary Table 2. Table  2, INTERACTION group and 314 genes Metabolic pathways tabs). Importantly, there is strong enrichment for genes in other KEGG pathways including tryptophan, methionine and branched chain amino acid metabolism. Additionally, within the Interaction group the term 'GO:0003333~amino acid transmembrane transport,' with 16 genes, is enriched. A distinguishing feature of most genes in these three pathways plus the transporters is up-regulation in Myr-treated cells with many genes peaking at 5 and 6 h. These data suggest that Myr-treated cells are attempting to increase amino acid uptake and de novo synthesis as a way to increase amino acid pool levels, however, we consider other interpretations in the Discussion.

Correlation analysis of transcriptomics data identifies clues to amino acid metabolism
As another approach to understand effects of Myr on the transcriptome, we analyzed transcript data using the R program package of Weighted Gene Correlation Network Analysis, WGCNA [28], to identify clusters (modules) of highly correlated genes showing a similar response to Myr over the 1-6 h time course. A gene dendrogram produced by average linkage of hierarchical clustering revealed many modules with the seven most correlated having from 18 to 2490 genes (Supplementary Table 3, column L in Cluster Analysis tab, see also the Module Eigengenes tab). We then examined the relationship of these gene modules to free AGING amino acid pools over the same 1-6 hr time course as an alternative way to determine if transcription played a role in lowering or maintaining smaller amino acid pools in drug-treated cells as we had previously found [7]. Results of this analysis are presented in Figure 3A as a heatmap plot. For amino acids (abbreviation shown on the X-axis), the numbers in a column represent the correlation p-value and the correlation (in parentheses) between a module eigengene (ME, Y-axis) expression value and an amino acid level. Here a ME corresponds to the first principal component of the expression matrix of the module. The ME can be considered the most representative gene in a module and, as used here, it represents the average effect of Myr on genes compared to the untreated drug control (Supplementary Table 3, Module Eigengenes tab). The heatmap shows that the Brown and Green modules with 1126 and 144 genes, respectively, have the greatest negative correlation to Myr treatment. That is, Myr maintains low pool levels but enhances the level of a transcript at one or more time points. The values in the 'InQ' column of the heatmap represent the correlation p-value and correlation value (parentheses) between the module eigengene expression and incubation time course (1-6 h). For these association calculations, the control (no drug) samples were set as 0 and Myr-treated samples were set as 1. Thus, the MEturquoise row has a strong correlation with time ( Figure 3A, column inQ) as shown graphically in Figure 3F. The values in the 'myriocin' column represent the association between the module eigengene expression and the absence/presence of Myr and they indicate that the Green and Brown MEs are the most negatively correlated MEs.
As another way to represent the effects of Myr treatment on MEs and to detect interactions between MEs, we performed network analysis and visualized the results by using Cytoscape. One such network, representing the relationship between genes in the Brown and Green modules, includes the Turquoise module because genes interact with multiple genes in both the Brown and Green modules ( Figure 3B). We further analyzed these three MEs by plotting the ME for each time point with and without Myr treatment ( Figure  3C). For example, the ME representing the Brown module reveals Myr induces transcription starting at 4 h ( Figure 3C) whereas there is no increase in amino acid pool size in drug-treated cells [7]. This opposite effect of Myr on transcription and amino acid pools corresponds to the negative correlation represented by the green shading in Figure 3A. The highest-ranking GO term in the Brown module is Regulation of Transcription ( Figure 3D) where drug treatment enhances transcript levels particularly after 4 h, consistent with the ANOVA analysis ( Figure 2B and 2C). A different negative correlation pattern is seen for the ME representing the small Green module where genes are up-regulated by Myr across the 1-6 h time frame compared to vehicle control but gene expression drops in both control and drug-treated samples starting at 2 h ( Figure 3E). The main GO term in the Green module is Translation (Biological Process, p = 1.7E-88) along with Ribosome Assembly and related processes responding to drug-induced slowing of protein synthesis and growth rate [7]. Lastly, the large Turquoise module with 2490 genes captures Myr-induced transcriptional events involving processes or pathways or cell components each with more than 200 genes ( Figure 3F, 3G). Interestingly, transcripts represented in the Turquoise module increase across the 1-6 h time frame with drugtreated cells always having higher transcript levels, substantiating the ANOVA analysis represented in Figure  1D (patterns T10a(↔), T&D04a(↑)., T&D05a(↑)), Figure  2B (patterns I05a(04a(↑)), I10a(04a(↑)) and I10a(10a(↑))) and Figure 2C. Functional enrichment analysis for all MEs in presented in Supplementary Table 3.

Deubiquitination of Mup1 impairs Myr-enhanced longevity
A drawback of WGCNA and pathway and pattern analyses is a reduced ability to identify significant biological features involving smaller numbers of genes or to identify genes that do not fit a specific pattern across the 1-6 h time-frame. To circumvent these limitations and to discover transcript changes with potential roles in lowering amino acid pools or novel roles in longevity enhancement, we examined genes with the highest possible significance (1E-16) in the Drug, Time and Interactions columns of data from the two-way ANOVA analysis of mRNA levels. This approach identified only 16 genes out of the 4964 genes (Supplementary Table 1, Filter tab, columns M, N and O). The UBI4 gene, which is stress-induced and encodes for ubiquitin, captured our attention because stresses of various types have known roles in aging and longevity and we previously identified increased stress responses in Myr-treated cells [27,29]. A defining feature of UBI4 transcript abundance in our studies is a rapid increase at the 5-  Table 2, Webgestalt patterns combined tab). Proteolysis is a highly enriched GO term in this pattern of gene expression ( Figure 2C, table on the right). The potential biological significance of the UBI4 transcript pattern may relate to our previous analyses of Mup1-pHluorin trafficking. We found that the fluorescent signal decreases in the PM more rapidly during the 4-7 h time-frame in Myr-treated cells, in a dose-dependent manner, which requires ubiquitin conjugation for endocytosis to occur [7].
To determine if ubiquitination has roles in Myrenhanced chronological lifespan (CLS), we examined untreated ubi4∆ cells and found a faster loss of viability compared to untreated BY4741 cells ( Figure 4B, 50% AGING In this experiment, a negative LOG value indicates increased abundance in the myriocin-treated sample (and vice versa). All measurements were normalized to the Heavy:Light ratio for total ubiquitin (unmodified peptides), which did not significantly differ with Myr treatment at either 2 or 4 h post-myriocin treatment. AGING survival day 3 vs. day 5, respectively). Additionally, the CLS of ubi4∆ cells treated with Myr does not increase and remains the same as untreated cells. Untreated wildtype (WT, BY4741) cells show 50% survival at day 5-6 which is extended to day 11 by Myr treatment (AUC, 95% CI: −6.828 to −2.840, p-value = 0.0025). From these data we conclude that UBI4 is required for Myr to enhance CLS.
Since ubiquitin has many functions, we studied its role in ubiquitin-mediated endocytosis of Mup1 using the same strains as used in our published analyses of endocytosis [7]. Specifically, CLS was assayed using cells with chromosomal MUP1 tagged with pHluorin-UL36 (catalytic active or catalytic dead, CD). UL36 is a viral deubiquitinase (DUB) that can reverse localized ubiquitin conjugation activity, thus protecting the fusion protein from ubiquitin-mediated degradation or trafficking events [30]. Our results reveal that Myr-enhanced longevity is prevented by Mup1 deubiquitination, which is a gain-offunction condition since it promotes PM stability. To determine how MUP1 loss of function affects drugenhanced longevity, we compared the CLS of mup1Δ cells to wild-type cells and cells with the mup1Δ allele replaced with wild-type MUP1. Myr treatment significantly enhanced CLS in each of the three strains ( Figure 4D, for WT BY4741, mup1Δ and MUP1 the AUC for Myr-treated vs. untreated cells, p-values ≤ 0.009, 0.026 and 0.0024, respectively). The viability curves for untreated mup1Δ and MUP1 cells do not differ significantly (AUC, 95% CI: −1.682 to 0.860, p-value = 0.42). We conclude from these data that deleting MUP1 does not prevent Myr-enhanced longevity. Given the increased UBI4 transcript level after 4-6 h of Myr treatment ( Figure 4A), we hypothesized that Myr treatment might affect total ubiquitin and its cellular distribution. To test this, we affinity purified FLAGubiquitin (from yeast cells harboring N-terminal FLAG fusions at two endogenous ubiquitin-encoding loci (RPS31 and RPL40B)) from untreated or Myr-treated yeast cells. Immunoblot analysis revealed that Myr treatment for 2 h or 4 h did not significantly alter the amount of FLAG-ubiquitin recovered or the amount of K48-linked ubiquitin polymers recovered ( Figure 4E). In contrast, we observed a significant increase in K63linked ubiquitin polymers after 4 h of Myr treatment ( Figure 4E, right lane in bottom immunoblot), suggesting that increased conjugation of ubiquitin in K63-linked polymers is part of the cellular response to depletion of sphingolipids.
To confirm this result, and to resolve other linkage types, we performed stable isotope labeling with amino acids in cell culture (SILAC) on yeast cells and subjected heavy-labelled and light-labelled cells to mock-treatment and Myr-treatment, respectively. At 2 hours and 4 hours of treatment, we collected cells, prepared lysates, and affinity purified FLAG-ubiquitin, followed by mixing, tryptic digestion, and processing of peptides for analysis by mass spectrometry. This SILAC-MS analysis resolved various peptides corresponding to the different linkage types of ubiquitin polymers, including K6-linked, K11-linked, K48-linked, and K63-linked polymers. Interestingly, we detected only modest changes in linkage types after 2 hours of Myr treatment, while 4 hours of Myr treatment resulted in increased formation of several polymer types, most notably for K63-linked ubiquitin polymers ( Figure 4F). This finding is consistent with the immunoblot results shown in Figure 4E, confirming that yeast cells increase the formation of K63-linked ubiquitin polymers in response to Myr treatment. Taken together, these results suggest that yeast cells respond to Myr treatment by remodeling ubiquitin pools, partly by increasing production of ubiquitin via transcription of UBI4 and partly by deploying existing ubiquitin pools to promote increased formation of K63linked ubiquitin conjugates.

DISCUSSION
We performed transcriptomics analysis to determine if Myr treatment caused a major transcriptional shift during a specific time-frame and to search for novel processes or pathways required for Myr-enhanced longevity. From our analysis, we are able to draw several conclusions. First, the major effect of Myr treatment is to up-regulate gene transcription ( Figures  1D and 2B and 2C -heatmaps and graphs). Second, transcription is robustly up-regulated after 5 to 6 h of Myr treatment when the majority of cells advance through the second cell division cycle [7] (Figure 2B and 2C). Third, transcriptional changes are not the major statistically significant force promoting initial amino acid pool lowering since there is no enrichment for genes involved in amino acid metabolism in transcript patterns ( Figure 1C and 1D). Additionally, enrichment analysis of genes in the Interaction group ( Figure 2C) found enrichment for several amino acid metabolic pathways and transmembrane transporters (Supplementary Table 2, Tabs labelled INTERACTION group and 314 Metabolic pathwayssee yellow highlighting), but these genes are mostly up-regulated by Myr treatment, suggesting that cells are attempting to increase amino acid pool size. Furthermore, correlation analysis found an anticorrelation or negative relationship between gene modules and amino acid pool sizes ( Figure 3A). For example, the large Turquoise module containing genes up-regulated across the 1-6 h time-frame ( Figure 3F) reveals enrichment for metabolic pathways including ones for amino acid metabolism ( Figure 3G and Supplementary Table 3, tab labelled GO terms Turquoise module). Notably, we showed previously [7] that amino acid pool lowering, starting after about 2 h of cell growth, is at least partly due to inactivation of amino acid transporters at the PM following Myr treatment. Fourth, ubiquitination of the methionine transporter Mup1 at the PM is vital for Myrenhanced longevity ( Figure 4B, 4C). However, we cannot altogether exclude a role for transcription of a small number of genes as playing roles in reducing or maintaining smaller amino acid pools. Lastly, Myr treatment strongly up-regulates K63-linked ubiquitination of proteins ( Figure 4E, 4F). Importantly, we have also recently reported that reducing sphingolipid synthesis inhibits the endocytosis of many nutrient transporters, while specifically promoting endocytic clearance of the methionine transporter Mup1 [31]. Taken together, our analysis reveals a complex cellular response to Myr that involves modulation of several biological processes including endocytic trafficking, proteostasis, and amino acid homeostasis.
The transcriptomics data presented here provide a valuable resource not only for constructing a more mechanistic understanding of how Myr treatment reprograms yeast cells to live longer, but it also provides a strategy for discovering new mechanisms to promote longevity. For instance, GO terms involving transcription are highly enriched in the D(↑) temporal pattern (Table 1, genes in this pattern are listed in Supplementary Table 2) where they represent 4 of the 6 top GO terms. Many of these encode transcription regulators including GCN4 whose transcription is upregulated by amino acid starvation in order to promote de novo amino acid biosynthesis [32]. Thus, the failure of cells to increase amino acid pools in response to increased GCN4 transcription implies activation of nontranscriptional mechanisms by Myr treatment to lower and maintain smaller amino acid pools. The GO term Metabolic Pathways (314 genes) was enriched in the Interaction group during a part or all of the 1-6 h timeframe ( Figure 2C, Supplementary Table 2,  INTERACTION group and 314 genes Metabolic  pathways tabs). Functional enrichment analysis of these 314 genes identified KEGG pathways for carbon metabolism, secondary metabolites, amino acid metabolism (tryptophan, methionine, branched-chain amino acids, lysine, arginine, proline, histidine) and other types of metabolism as being significantly enriched. Several of these pathways have known roles in longevity (e.g., methionine, branched-chain amino acids, glycogen, trehalose metabolism, etc.), suggesting that the Myr-sensitive pathways defined in this work are novel but include elements of pathways defined in prior work as having roles in longevity. Endocytosis is another GO term found to be highly enriched in the Interaction analysis of transcripts across the 1-6 h timeframe (Supplementary Table 2, INTERACTION group tab). GO term analysis of this group of 54 genes found enrichment for genes involved in ubiquitin-mediated endocytosis which provide clues to identify the proteins, lipids and cellular machinery controlling the Myrinduced endocytosis of Mup1 that we previously observed [7]. Detailed analysis of the gene interaction network shown in Figure 3B should aid in identifying the most important Node/Hub genes along with their edges potentially involved in controlling Mup1 endocytosis. Lastly, a novel feature of the Interaction group is enrichment of GO terms representing nearly every cellular component, many of which are membrane-bound (Supplementary  Table  2, INTERACTION group tab). We hypothesize that Myr treatment starts early (~ 5 h or the second cell division) to widely and coordinately rewire cellular physiology towards enhanced longevity well before many hallmarks of aging and longevity are detectible. For instance, we previously showed that autophagy begins after 3-5 cell divisions [27] whereas many autophagy genes in the INTERACTION group are strongly up-regulated in the second cell division (4-5 h) following Myr treatment. Our transcriptomics dataset will be a valuable resource for determining if coordinated rewiring is promoted by cell cycle mechanisms or cross-talk between cellular components or by both mechanisms.
In addition to the limitations inherent in WGCNA, pattern and pathway analyses such as genes not falling into any category, there are other potential limitations in our analyses to consider. For example, improved detection and resolution of drug effects might be obtained with increased cell cycle synchronization, more time points and more sample replicates but these are not likely to alter the overall results of this work. Additionally, analysis of transcripts before 1 hr might identify potential mechanisms for maintaining smaller amino acid pools in drug-treated cells. Adding Myr after cells have reentered the cell division cycle would be an alternative approach for studying early effects of Myr on transcription, but integrating such data with results presented here and in our earlier studies could be difficult [7,27,29].
Healthy proteostasis relies on proper functioning of protein degradation networks to maintain proteome quality control and prevent accumulation of misfolded and damaged proteins. Thus, a decline in proteostasis is a hallmark of aging and age-related diseases [33]. Analysis of the yeast transcriptional response to Myr treatment revealed a significant induction of the mRNA transcript encoded by UBI4 ( Figure 4A)one of four genes in the yeast genome encoding for ubiquitin. The other three [34] ubiquitin genes (RPS31, RPL40A, and RPL40B) encode ubiquitin fusions to ribosomal subunits and generate most cellular ubiquitin in normal growth conditions, which inherently couples the translational and degradative capacities in the cell. By comparison, UBI4 encodes a linear (heat-to-tail) ubiquitin pentamer that does not contribute significantly to the total pool of ubiquitin in normal growth conditions. However, in conditions of heat stress, oxidative stress, or starvation UBI4 is transcriptionally induced which serves to increase the degradation capacity of the cell and uncouple it from ribosome biogenesis and translational capacity [34,35]. Thus, our observation that Myr induces transcription of UBI4 but not the other genes encoding ubiquitin (transcript levels of other ubiquitin genes can be displayed by using the 'Graph Reporter' tab in Supplementary Table 1) is suggestive of a proteotoxic and/or starvation stress response involving activation of protein degradation networks and sizeable proteome remodeling. The importance of UBI4 in Myr-enhanced longevity is underscored by our finding that the CLS of untreated ubi4Δ yeast cells decreased compared to WT cells and was not significantly enhanced by Myr treatment ( Figure 4B). In addition, our finding that cells with chromosomal MUP1 tagged with pHluorin-UL36 fail to show Myr-enhanced CLS ( Figure 4C) suggests a novel role for ubiquitination in longevity. Importantly, we previously reported that UL36 fusion to Mup1 impedes its Myr-induced endocytosis [7], underscoring the linkage of ubiquitin-mediated endocytosis of Mup1 to Myr-mediated longevity. One limitation of this analysis is the potential for proximal effects on other proteins in the vicinity of the Mup1-pHluorin-UL36 protein [30,36]. Thus, we cannot exclude the possibility that ubiquitination of Mup1-associated factors is critical for Myr-mediated longevity. Future studies will be required to evaluate if such proximal effects are involved in the impairment of Myr-enhanced lifespan. The outcomes of these studies will potentially provide new targets or strategies for improving human healthspan.

Strains, culture conditions, lifespan assays and statistical significance
Strains (Table 3), culture conditions, lifespan assays and their statistical significance were similar to ones described previously [7]. Specifically, yeast cells, made prototrophic by transformation with pHLUM (carries HIS3, LEU2, URA3 and MET, Addgene, Watertown, MA, USA [37]) and stored at −80°C, were streaked onto minimal glucose plates and incubated at 30°C for 3 days. One to three colonies were inoculated into 5 ml of SD complete medium (SDC) [38] and incubated on a gyratory water bath shaker for 18-24 h at 30°C. Typical culture densities were 7-9 A600 nm/ml (1.5 × 10 7 cells/A600 nm). Cells were diluted into 25 ml fresh SDC warmed to 30°C in a 125 ml long-neck flask, containing vehicle (EtOH) or Myr (details described in the next paragraph), to give a starting cell density of 0.15 A600 nm. Flasks were incubated at 30°C on a gyratory air-bath shaker (200 rpm). After 72 h of incubation cells were diluted (1:50 × 1:100) and from 10 to 40 μL was spread on a YPD plate [38] and incubated 48-72 h at 30°C followed by colony counting. CLS assays were performed at least twice using triplicate cultures. Concentrations of Myr used in CLS assays ranged from 475-525 ng/ml depending on the sensitive of strains compared to wild-type, prototrophic BY4741. Drug-sensitivity was measured by culture density (A600 nm) after 24 h of growth using cultures started at 0.15 A600 nm units of cells/ml. All BY4741 yeast strains used for lifespan assays and for RNA extraction were made prototrophic by transformation with the pHLUM plasmid. The DUB (UL36) fusion yeast strains used in this study were generated by homologous recombination in BY4741 and SEY6210 background strains using the reagents and strategy previously described (Hepowit et al. 2022

RNA extraction
Culture conditions, medium and prototrophic BY4741 cells were the same as for lifespan and amino acid pool assays [7]. RNA was extracted from 5 A600 nm unit/ml of yeast cells harvested by filtration on a membrane filter at each time point (1-6 h, Figure 1A) [7]. Filtered cells were washed once with 5 ml of ice-cold nanopure water and the filter was quickly transferred to a chilled 1.5 ml microfuge tube containing 0.5 ml cold nanopure water. Tubes were vortexed 10 sec followed by centrifugation for 15 sec. Supernatant fluid (450 μl) was transferred to a new tube and frozen in a dry-ice EtOH bath followed by storage at −80°C. Cold acid-washed glass beads (300 μl, 0.5 mm dia.) were added to a frozen cell pellet followed by addition of 300 μl of RLT buffer (RNAeasy mini kit, Qiagen, Germantown, MD, USA). Tubes were vortexed 5 min at room temperature and placed on ice for 1 min. After 4 cycles 300 μl of ice-cold RLT buffer was added followed by mixing and centrifugation at 13,000 × g for 2 min at room temperature. Supernatant fluid (450 μl) was transferred to a new microfuge tube and then mixed with 1 ml of 70% (v/v) EtOH before transfer to a RNeasy spin column and processed according to the manufacturer's instructions. Processed samples were frozen in a dry-ice EtOH bath and stored at −80°C. Resulting FASTQ files mapped to the yeast genome (R64), resulted in count files and normalized using the transcripts per million (TPM) algorithm [39] using WebMev [40]. Resulting data were downloaded as flat files and loaded into Excel for further analysis. From a total of 6198 mapped genes, 5169 were uniquely annotated with gene symbols and had sufficient nonzero readings for further analysis. The filtered data were analyzed by two-way ANOVA for the main effects of drug and time, as well as for interaction and assigned a p-value. A gene found to be significant in any of these three categories was referred to as differentially expressed (DEG). Significant by the time term or both the drug and time term, data were further analyzed by post-hoc pairwise Fisher's protected Least Significant Difference (pLSD), and log 2-fold change comparison to further isolate the effects of drug over time.

RNA-seq analysis
To identify DEGs that change up or down over time we built a 'template assignment tool' that took each DEG and attempted to assign it to each of 63 patterns. The pattern to which a DEG fit best (by Pearson's correlation, had to be R > 0.85) became the pattern to which it was assigned. This is reported in Figure 1C (99.86% of all DEGs were assigned to a template). To determine which templates had significantly more DEGs than expected by chance, we ran a binomial testthis is the basis for some templates being bolded and having an asterisk in Figure 1C-they had a binomial p < 0.01. The 12 patterns of expression that had more DEGs than expected by chance assigned to them explained the majority (77%) of all DEGs. Data have been deposited in the GEO (GSE199904) [NCBI tracking system #22817261]. GO terms were analyzed by using DAVID software [41] or online WebGestalt.

Gene regulatory network using WGCNA
The WGCNA (v1.70-3) [28] was used to identify gene modules and build unsigned co-expression networks, which include negative and positive correlations. Briefly, WGCNA constructs a gene co-expression matrix, uses hierarchical clustering in combination with the Pearson correlation coefficient to cluster genes into groups of closely co-expressed genes termed modules, and then uses singular value decomposition (SVD) values as module eigengenes to determine the similarity between gene modules or to calculate association with sample traits (for example, incubation time or amino acid levels). The top 2,000 variable genes were used to identify gene modules and network construction. Soft power 8 was chosen by the WGCNA function pick SoftThreshold. Then, the function TOMsimilarityFromExpr was used to calculate the TOM similarity matrix via setting power = 8, networkType = "signed. The distance matrix was generated by subtracting the values from similarity adjacency matrix by one. The function flashClust (v.1.01) was used to cluster genes based on the distance matrix, and the function cutreeDynamic was utilized to identify gene modules by setting deepSplit = 3. Cytoscape (v.3.8.2) was applied for the network visualizations.