Simultaneous proteome localization and turnover analysis reveals spatiotemporal features of protein homeostasis disruptions

Abstract

The spatial and temporal distributions of proteins are critical to protein function, but cannot be directly assessed by measuring protein bundance. Here we describe a mass spectrometry-based proteomics strategy, Simultaneous Proteome Localization and Turnover (SPLAT), to measure concurrently protein turnover rates and subcellular localization in the same experiment. Applying the method, we find that unfolded protein response (UPR) has different effects on protein turnover dependent on their subcellular location in human AC16 cells, with proteome-wide slowdown but acceleration among stress response proteins in the ER and Golgi. In parallel, UPR triggers broad differential localization of proteins including RNA-binding proteins and amino acid transporters. Moreover, we observe newly synthesized proteins including EGFR that show a differential localization under stress than the existing protein pools, reminiscent of protein trafficking disruptions. We next applied SPLAT to an induced pluripotent stem cell derived cardiomyocyte (iPSC-CM) model of cancer drug cardiotoxicity upon treatment with the proteasome inhibitor carfilzomib. Paradoxically, carfilzomib has little effect on global average protein half-life, but may instead selectively disrupt sarcomere protein homeostasis. This study provides a view into the interactions of protein spatial and temporal dynamics and demonstrates a method to examine protein homeostasis regulations in stress and drug response.

Fig. 1. Overview of the SPLAT strategy.

a Experimental workflow. Control, thapsigargin-treated, and tunicamycin-treated human AC16 cardiomyocytes were labeled with ¹³C₆¹⁵N₂ L-Lysine and ¹³C₆¹⁵N₄ L-Arginine dynamic SILAC labels. For each condition, 3 biological replicate SPLAT experiments were performed (n = 3). After 16 h, cells were harvested and mechanically disrupted, followed by differential ultracentrifugation steps to pellet proteins across cellular compartments. Proteins from the ultracentrifugation fractions were digested and labeled using tandem mass tag (TMT) followed by mass spectrometry. b Dynamic SILAC labeling allowed differentiation of pre-existing (unlabeled, i.e., SILAC light) and post-labeling (heavy lysine or arginine, i.e., +R[10.0083]) synthesized peptides at 16 h. The light and heavy peptides were isolated for fragmentation separately to allow the protein sedimentation profiles containing spatial information to be discerned from TMT channel intensities. c Computational workflow. Mass spectrometry raw data were converted to mzML format to identify peptides using a database search engine. The mass spectra and identification output were processed using Riana (left) to quantitate the time-dependent change in SILAC labeling intensities and determine the protein half-life, and using pyTMT (right) to extract and correct TMT channel intensities from each light or heavy peptide MS2 spectrum. The TMT data were further processed using pRoloc/Bandle to predict protein subcellular localization via supervised learning. d Temporal information and spatial information are resolved in MS1 and MS2 levels, respectively. SPLAT allows the subcellular spatial information of the heavy (new) and light (old) subpools of thousands of proteins to be quantified simultaneously in normal and perturbed cells. HL: Half-life.

Fig. 2. Simultaneous measurements of spatial and temporal kinetics under UPR.

a Bar graphs of expression ratios of known ER stress markers upon 8 and 16 h of thapsigargin (Thps.) (n = 6 normal; n = 3 thapsigargin). *: P < 0.01; **: P < 0.001; ***: P d < 0.0001; limma multiple-testing corrected (FDR adjusted) P (two-sided). Error bars: expression ratios ± s.d. b PC1 and PC2 of proteins spatial map showing the localization of confidently allocated proteins in normal and thapsigargin-treated AC16 cells. Data point: protein; color: subcellular compartment classification. c Distribution of light (unlabeled) protein features in each of the 12 subcellular compartments (n = 3); fill color represents whether the protein is also annotated to the same subcellular compartment in UniProt Gene Ontology Cellular Component terms. d Histograms of determined log10 protein turnover rates in control and thapsigargin cells (n = 3). Text overlay: median half-life. e Boxplot of log2 turnover rate ratios in thapsigargin /normal cells for proteins localized to the ER (blue) (T) or not (F); or the Golgi (GA; green). P values: Mann–Whitney test (two-sided). A Bonferroni corrected threshold of 0.05/13 is considered significant. Center line: median; box limits: interquartile range; whiskers: 1.5x interquartile range; n = 286, 2234, 62, 2458 proteins over 3 independent experiments per group. f First-order protein kinetic curves in normal (gray), and thapsigargin treated (red) AC16 cells of 4 known ER stress markers with elevated turnover (HSPA5, RCN3, HSP90B1, PDIA4). Point: best-fit k; line: first-order kinetic curve of k; bands: fitting s.e. g Turnover rate ratios (thapsigargin vs. normal) of top proteins with elevated temporal kinetics in UPR. Color: compartment; P.adj: Mann–Whitney test (two-sided) with Benjamini-Hochberg multiple-testing correction. Dashed line: 1:1 ratio; point: median ratio; range: median ± MAD/median of ratios; n=number of peptide observations (parenthesis) over 3 biologically independent samples per group. h Gene set enrichment analysis (GSEA) of turnover rate ratios in thapsigargin treatment. Color: multiple-testing corrected (FDR adjusted) P in GSEA (two-sided); x-axis: normalized enrichment score (NES). Size: proteins in the gene set.

Fig. 3. SPLAT captures extensive protein localization changes under UPR.

a Alluvial plot of differential localization (DL) events (> 0.99 BANDLE DL probability; estimated FDR < 1%) following thapsigargin treatment showing a cohort of proteins moving from the Golgi apparatus (GA) and lysosome towards the plasma membrane (PM) (n = 3). b Spatial map for SLC3A2 (open black circle) in normal (left) and thapsigargin-treated (right) AC16 cells, showing its lysosomal assignment in normal cells and PM assignment in stressed cells. Colors: allocated subcellular compartment. c Ultracentrifugation fraction profile of SLC3A2 and other amino acid transporters SLC7A5, SLC1A4, SLC1A5 and ion channel proteins SLC30A1, ATP1B1, ATP1B3, and ATP2B1 with similar migration patterns. X-axis: fraction 1 to 10 of ultracentrifugation. Y-axis: relative channel abundance. Bold lines: protein of interest; light lines: ultracentrifugation profiles of all proteins classified to the compartment. Colors correspond to subcellular localization in panel b and all AC16 data in the manuscript; numbers within boxes: BANDLE allocation probability to compartment; numbers at arrows: BANDLE DL probability. d Immunofluorescence of SLC3A2 (red) against the lysosome marker LAMP2 (green) and DAPI (blue). Numbers in cell boundary: colocalization score per cell. Scale bar: 90 µm. e Colocalization score (Mander’s correlation coefficient) between SLC3A2 and LAMP2 decreases significantly (two-tailed unpaired t-test P value: 3e–8) in thapsigargin, consistent with movement away from lysosomal fraction (n = 205 normal cells, n = 32 thapsigargin treated cells). Center line: median; box limits: interquartile range; whiskers: 1.5x interquartile range; points: outliers. f Alluvial plot showing the migration of ER, GA, and nucleus proteins toward the peroxisome/endosome containing fraction in thapsigargin treated cells. g Ultracentrifugation fraction profile of stress granule proteins UBAP2, USP10, CNOT1, CNOT2, ZC3H7A, and NUFIP2. RNA Granule Score 7 or above is considered a known stress granule protein. Phi: predicted phase separation participation. Circle denotes a prediction of True within the database. RBP: Annotated RNA binding protein on the RNA Granule Database. One circle denotes known RBP in at least one data set; two circles denote known RBP in multiple datasets.

Fig. 4. SPLAT reveals protein-lifetime dependent differential localization.

a Histogram of distances in light and heavy proteins in normalized fraction abundance profiles in normal and thapsigargin-treated (Thps.) AC16 cells. X-axis: the spatial distribution distance of two proteins is measured as the average euclidean distance of all TMT channel relative abundance in the ultracentrifugation fraction profiles across 3 replicates; y-axis: count. Blue: distance for 1,614 quantified light-heavy protein pairs (e.g., unlabeled vs. SILAC-labeled EGFR). Grey: distribution of each corresponding light protein with another random light protein. P value: Mann–Whitney test (two-sided). b Proportion of heavy-light protein pairs with confidently assigned localization to the same location (purple) in normal (left; 93%) and thapsigargin (right; 89%) cells. c Ranked changes in heavy-light pair distance in thapsigargin treatment. The positions of EGFR and ITGAV are highlighted. Inset: Z score distribution of all changes. d Spatial map of the light and heavy EGFR in normal and thapsigargin-treated cells. Each data point is a light or heavy protein species. Color: assigned compartment. Numbers: euclidean distance in fraction profiles over 3 replicates. e Corresponding fraction profiles; x-axis: ultracentrifugation fraction; y-axis: fractional abundance. Post-labeling synthesized EGFR is differentially distributed in thapsigargin and shows ER retention (blue), whereas the preexisting EGFR pool remains to show a likely cell surface localization (pink) after thapsigargin. f, g As above, for ITGAV. h Confocal imaging of EGFR immunofluorescence supports a partial relocalization of EGFR from the cell surface toward internal membranes following thapsigargin. Numbers: ratio of EGFR at the plasma membrane vs. whole cells. Blue: DAPI; green: EGFR; scale bar: 90 µm. i Cell areas; Mann-Whitney two-sided P: 0.72. Center line: median; box limits: interquartile range; whiskers: 1.5x interquartile range. j Total EGFR intensity per cell; Mann-Whitney two-sided P: 2.2e-05. Center line: median; box limits: interquartile range; whiskers: 1.5x interquartile range. k Edge/total intensity ratios in normal and thapsigargin-treated AC16 cells (n = 71 normal cells; n = 93 thapsigargin cells; Mann–Whitney two-sided P: 1.7e–4). Center line: median; box limits: interquartile range; whiskers: 1.5x interquartile range.

Fig. 5. Comparison of ER stress induction methods.

a Bar charts showing activation of known ER stress markers upon tunicamycin (Tunic.) treatment for 8 h and 16 h. X-axis: ER stress markers; y-axis: expression ratio (n = 6 normal AC16; n = 3 tunicamycin). *: P < 0.01; **: P < 0.001; ***: P < 0.0001; multiple-testing corrected (FDR adjusted) P (two-sided). Error bars: expression ratios ± s.d. b Histograms of the determined log10 protein turnover rates in control and tunicamycin-treated cells (n = 3). c Boxplot showing the log2 turnover rate ratios in tunicamycin over normal AC16 cells for proteins that are localized to the ER (T) or not (F); Golgi apparatus (GA), or the lysosome. P values: two-tailed t-test. Center line: median; box limits: interquartile range; whiskers: 1.5x interquartile range; n = 286, 2234, 62, 2458, 139, 2381 proteins over 3 biological independent samples. A Bonferroni corrected threshold of 0.05/13 is considered significant. d Gene set enrichment analysis (GSEA) of turnover rate ratios in tunicamycin treatment. Color: multiple-testing corrected (FDR adjusted) P in GSEA (two-sided) x-axis: normalized enrichment score (NES). Size: proteins in the gene set. e Example of best fit curves in the first-order kinetic model at the protein level between normal (gray), tunicamycin (blue) and thapsigargin (red) treated AC16 cells showing known ER stress markers with elevated turnover in both ER stress inducers (HSPA5, HSP90B1, and PDIA4) as well as stress response proteins with elevated turnover only in tunicamycin (PDIA3, DNAJB11, NIBAN1). Point: best-fit k; line: first-order kinetic curve of k; bands: fitting s.e. f Alluvial plot showing the migration of ER, GA, and peroxisome/endosome proteins toward the lysosome (left). On the right, the ultracentrifugation fraction profiles of the differentially localized proteins RRBP1, FKBP11, GANAB, MANF, IKBIP, and P3H1 are shown that are targeted toward the lysosome in both tunicamycin and thapsigargin treatment (BANDLE differential localization probability > 0.95). Numbers in boxes are the BANDLE allocation probability in each condition (n = 3).

Fig. 6. Applicability in human iPSC-derived cardiomyocytes.

a Schematic of human iPSC-CM differentiation, carfilzomib treatment, and SPLAT analysis. b Confocal microscopy images showing sarcomeric disarray in iPSC-CMs upon 48 hrs of 0.5 µM carfilzomib; green: cTNT, red: alpha-actinin; blue: DAPI; scale bar: 20 µm. c–h Cell viability (%), normalized Seahorse oxygen consumption rate (OCR; pmol/min), basal respiration, maximal respiration, proton leak, and ATP production upon 0–5 µM carfilzomib for 24 or 48 hrs; ∙: P < 0.1; *: P < 0.05; **: P < 0.01, ANOVA with Tukey’s HSD post-hoc (two-sided) at 95% confidence level; n = 5 biologically independent experiments. Error bars: mean ± s.d. for bar charts in panels c, e, f, g, h; mean ± s.em. for the OCR graph in panel d. Colors in panel d: dosage, same as panel c. O: Oligomycin; AA/R: Antimycin A/Rotenone. i Spatial map with 13 assigned subcellular localizations in iPSC-CMs at the baseline (top) and upon 0.5 µM carfilzomib treatment (n = 2). Source data are provided as a Source Data file. j Histogram of log10 protein turnover rates (k), with median half-life of 97.4 h and 100.0 h in normal and carfilzomib-treated iPSC-CM. k Proteasome activity in iPSC-CMs treated with 0 (Ctrl) vs. 0.5 µM carfilzomib (Cfzb) for 48 h. P value: two-tailed t-test; n = 3 biologically independent samples; error bar: mean ± s.d. l Autophagy assay for iPSC-CMs treated with 0 (Ctrl) vs. 0.5 µM carfilzomib (Cfzb) for 48 h, and positive control (Pos); data were normalized to DAPI and normal cells. P value: two-tailed t-test; n = 10 biologically independent samples; error bar: mean ± s.d. m log2 Turnover rate ratios between carfilzomib-treated and untreated iPSC-CM from the spatiotemporal proteomics data. P values: Mann–Whitney test (two-sided); with a Bonferroni threshold of 0.05/14 considered significant. Center line: median; box limits: interquartile range; whiskers: 1.5x interquartile range; n = 60, 2572, 175, 2457, 51, 2581 proteins over 2 biologically independent samples.

Fig. 7. Proteostatic pathways and lesions in carfilzomib mediated cardiotoxicity in iPSC-CMs.

a Changes in protein turnover rates between carfilzomib vs. normal iPSC-CMs across selected cellular compartments; ∙: P < 0.1; *: <0.05; ** <0.01; Mann-Whitney test (two-sided) with Benjamini-Hochberg multiple-testing correction. Dashed line: 1:1 ratio; point: median ratio; range: median ± MAD/median of ratios; n=number of peptide observations (parenthesis) over 2 biologically independent samples per gorup. Source data are provided as a Source Data file. b Kinetic curve representations of proteins with accelerated temporal kinetics in carfilzomib, including PSMC2 which corresponds to the ratio in panel A, as well as additional ERAD proteins and chaperones; gray: normal iPSC-CM; green: carfilzomib. Point: best-fit k; line: the first-order kinetic curve of k; bands: fitting s.e. c Kinetic curve representations of slowdown of protein kinetics in DSP, DMD, MYH6, and MYH7, corresponding to the ratios in panel a. Point: best-fit k; line: first-order kinetic curve of k; bands: fitting s.e. d, e Spatial map (PC1 vs. PC2) and ultracentrifugation fraction profiles of d. BAG3 and e. PSME4 in normal and carfilzomib-treated human iPSC-CM, showing a likely differential localisation in conjunction with kinetics changes. White-filled circles: light and heavy BAG3 or PSME4 in each plot. The kinetic curves of BAG3 and PSME4 are in panel b. Numbers at arrows correspond to BANDLE differential localization (DL) probability.