We previously performed laser microdissection to collect areas of single fibers from muscle biopsies of 18 patients with IBM [14]. These samples were taken from normal healthy fibers, or in the case of IBM-affected muscles, from affected RV-containing fibers and adjacent normal appearing fibers. We then analyzed these samples by mass spectrometry using label-free spectral count-based relative protein quantification (see Methods). For the study presented here, healthy control and IBM proteomic datasets were generated from healthy control myofiber regions (HCs), unaffected myofiber regions from IBM patients (disease controls, DCs), non-vacuole containing sarcoplasmic regions of affected fibers (AFs), and myofiber regions containing rimmed vacuoles (RVs) (Fig. 2a). Fig. 2

Proteins in rimmed vacuoles from protein aggregation myopathies are supersaturated. Representative images of: (a) healthy control myofibers (HC), control unaffected myofibers in diseased samples (DC), surrounding tissues of affected fibers (AF), and rimmed vacuoles (RV) from human subjects with inclusion body myositis, and (c) DC and AF samples from human subjects with myotilin mutations. Outlines represent areas for LMD. In (c), prior to LMD, muscle was immunostained with an antibody directed to myotilin (green) to identify aggregate containing fibers (AF). b, d, e Comparison of the unfolded supersaturation scores (σ_u) of the proteome (Prt) (N = 15,954) and (b) proteins enriched in RVs (RV) (N = 50), amyloid plaques (Plq) (N = 26), neurofibrillary tangles (NFT) (N = 76), proteins found in TDP-43 inclusions (TDP) (N = 32); (d) proteins enriched in affected fibers from any of three protein aggregation myopathies (hPAM) (N = 50); and (e) proteins enriched in affected fibers from individual protein aggregation myopathies involving desmin (Desm) (N = 6), filamin (Fil) (N = 16), and myotilin (Myot) (N = 46) mutations. f, g Comparison of the folded supersaturation scores (σ_f) for the proteome (Prt) (N = 1605) and (f) the proteins enriched in affected fibers from any of three protein aggregation myopathies (hPAM) (N = 46) and (g) the proteins enriched in affected fibers from individual protein aggregation myopathies involving desmin (Desm) (N = 5), filamin (Fil) (N = 15), and myotilin (Myot) (N = 43) mutations. The fold change (Δ) represents the fold difference in the median σ_u or σ_f scores between each inclusion type and the proteome. The median σ_u or σ_f supersaturation score for the proteome is normalized to 0. Boxes range from the 25th percentile to 75th percentile, while whiskers extend to maximum and minimum data points up to 1.5x interquartile range above and below the limits of the boxes. Remaining outliers are plotted as open circles. Statistical significance determined by one-tailed Wilcoxon/Mann-Whitney test with Holm-Bonferroni correction. *p < 0.05, **p < 0.01, ***p < 0.001, ****p < 0.0001

Comparison of these datasets enabled us to identify a set of proteins enriched within RVs, as compared to DCs. This list of 53 RV-enriched proteins includes 17 proteins previously identified to accumulate in IBM tissue (Additional file 1: Dataset S1). We next asked whether these proteins share similar biophysical features despite their different sequences, structures and functions. We had previously estimated supersaturation of a protein as the product of its predicted aggregation propensity (given by the Zyggregator score (Z_agg) which correlates negatively with its solubility) and its expression level, either based on mRNA levels from microarray data or proteomic analysis [3].

We thus compared the supersaturation levels of RV-enriched proteins to those of co-aggregators within amyloid plaques [24], neurofibrillary tangles [38], and TDP-43 inclusions [5] (Additional file 1: Dataset S2). As an approximation of the supersaturation level for a given protein, we used mRNA levels averaged over dozens of different human tissues unaffected by misfolding disease (Additional file 1: Dataset S4) and aggregation propensities predicted from the primary sequences for the unfolded states of proteins (Z_agg) (Additional file 1: Dataset S3) termed the unfolded supersaturation score (σ_u) (Additional file 1: Dataset S2). While this approach does not benefit from tissue specificity, it was previously shown that this average estimate demonstrated elevated supersaturation scores for proteins associated with aggregation and cellular pathways implicated in neurodegenerative disorders and enabled the direct comparison of inclusions from muscle to the central nervous system [3].

We found that proteins enriched in RVs have elevated supersaturation scores (σ_u) in control tissues (RV, median fold change (Δ): 2.4x, p = 1.4•10^− 6). This was similar to proteins observed to co-aggregate (co-aggregators) with plaques (median Δ: 6.0x, p = 4.5•10^− 8) and neurofibrillary tangles (median Δ: 5.2x, p = 1.3•10^− 13) in AD, and TDP-43 (median Δ: 2.1x, p = 1.8•10^− 3) in ALS, respectively (Fig. 2b, Additional file 1: Dataset S7). The elevated supersaturation score for RV-enriched proteins was also present when we calculated tissue-specific supersaturation scores (σuts\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\upsigma}_u^{ts} $$\end{document}₎ using the subset of the cross-tissue microarray expression database that included skeletal muscle expression (RV: median Δ: 2.1x, p = 2.2•10^− 6) (Additional file 2: Figure S1; Additional file 1: Dataset S7). Comprehensive statistical results are shown in Additional file 1: Dataset S12.

To determine whether the phenomenon of supersaturation observed for IBM-associated proteins (Fig. 2b), an acquired PAM, is also observed for proteins associated with hereditary PAMs (hPAM), we extended our studies to proteomic datasets from laser microdissected myofibers of muscle biopsies of patients with three different genetically defined hPAMs (10 patients with DES mutations, 7 patients with FLNC mutations and 17 patients with MYOT mutations) [20, 25, 26]. Samples were taken from affected aggregate-containing fibers (AF) or adjacent normal appearing disease control fibers (DC) (Fig. 2c). We then identified proteins that are enriched within the aggregate-containing fibers, as compared to unaffected disease control fibers (Additional file 1: Dataset S2). The σ_u score is similarly elevated for the proteins enriched in hPAM aggregate fibers (AF) (median Δ: 2.2, p = 6.9•10^− 4) (Fig. 2d, Additional file 1: Dataset S7). We note, however, that sample size limitations led to statistically insignificant results for two of the three individual hPAMs (desminopathy median Δ: 8.5x, p = 9.8•10^− 2; filaminopathy median Δ: 1.9x, p = 8.3•10^− 2, myotilinopathy median Δ: 1.8x, p = 6.7•10^− 3) (Fig. 2e). We then calculated σuts\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\upsigma}_u^{ts} $$\end{document} and found the increased supersaturation of proteins in aggregate-containing tissue is significant in this context (hPAM: median Δ: 4.5x, p = 1.2•10^− 8; desminopathy median Δ: 11x, p = 2.5•10^− 2; filaminopathy median Δ: 5.6x, p = 1.4•10^− 3, myotillinopathy median Δ: 3.9x, p = 7.5•10^− 7) (Additional file 2: Figure S1; Additional file 1: Dataset S7). In addition, we estimated the significance of the increase in supersaturation of the σuts\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\upsigma}_u^{ts} $$\end{document} scores relative to the σ_u scores (p < 1•10^− 6).

Particular to our proteomic datasets, we can evaluate the degree of supersaturation within the context of the skeletal muscle proteome rather than abundances from mRNA levels in public databases. Thus, we asked whether aggregate-enriched proteins in hPAMs were supersaturated, based on their abundance in background of a healthy muscle proteome (Additional file 1: Dataset S6). To do this, we combined protein abundances derived from healthy control muscle using a version of the Zyggregator algorithm that weights residue-level aggregation propensities based on predictions of the relative burial of proteins after folding (ZaggSC\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\mathrm{Z}}_{\mathrm{agg}}^{\mathrm{SC}} $$\end{document}), as described previously [32] (Additional file 1: Dataset S3). This estimate is termed folded supersaturation score (σ_f) as compared with the previous estimate of the unfolded supersaturation score (σ_u) (Additional file 1: Dataset S9). To directly compare these two estimates (σ_u and σ_f), we calculated the σ_f score of proteins enriched in hPAM aggregate-containing fibers (median Δ 9.1x, p = 5.3•10^− 5) (Fig. 2f compared with Fig. 2d). The elevated supersaturation score among proteins enriched in hPAM aggregate-containing fibers relative to HC is a result of both abundances and aggregation propensities higher than those of the proteome (Δ: 6.3x, p = 3.9•10^− 7; ZaggSC\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\mathrm{Z}}_{\mathrm{agg}}^{\mathrm{SC}} $$\end{document} Δ: 2.2x, p = 0.049) (Additional file 2: Figure S2A-B, Additional file 1: Datasets S3 and S9). Similarly, we found elevated σ_f scores when considering proteins enriched in aggregate-containing tissue specific to each hPAM (Fig. 2g, Additional file 1: Dataset S2) (desminopathy median Δ: 31x, p = 1.4•10^− 2; filaminopathy median Δ: 33x, p = 8.8•10^− 4; myotillinopathy median Δ: 13x, p = 5.7•10^− 5).

In unaffected diseased control myofibers (DC), we found that the σ_f scores of proteins enriched in aggregate-containing myofibers (AF) were elevated relative to the proteome (desminopathy median Δ: 16x, p = 1.2•10^− 2; filaminopathy median Δ: 4.6x, p = 5.7•10^− 2; myotilinopathy median Δ: 4.4x, p = 2.7•10^− 3) (Fig. 3a, d, g). In comparison, the σ_f scores of these proteins in affected fibers (AF) were higher (desminopathy median Δ: 1700x, p = 4.3•10^− 5; filaminopathy median Δ: 48x, p = 4.2•10^− 5; myotilinopathy median Δ: 33x, p = 3.6•10^− 8) (Fig. 3b, e, h). To confirm the robustness of these results, we introduced varying amounts of random noise into our data, and found that the results are robust even when noise of at least 5 times the magnitude of the signal, and in many cases as much as 20x the magnitude of the signal, is introduced (Additional file 2: Figures S3-S4). Fig. 3

Protein supersaturation in hereditary protein aggregate myopathies. Comparison of folded supersaturation scores (σ_f) between the proteome and proteins enriched in aggregate-containing myofibers (AF) for: (a, b, c) desminopathy (DC: Prt N = 611, Agg N = 6; AF: Prt N = 1387, Agg N = 6), (d, e, f) filaminopathy (DC: Prt N = 333, Agg N = 16; AF: Prt N = 507, Agg N = 16), and (g, h, i) myotillinopathy (DC: Prt N = 680, Agg N = 46; AF: Prt N = 742, Agg N = 46). Box plots and statistical tests as in Fig. 1. *p < 0.05, **p < 0.01, ***p < 0.001, ****p < 0.0001

We used our IBM proteomic datasets, to segment the data starting from healthy controls (HC) and continuing to unaffected fibers in affected patients (DC), areas from affected fibers surrounding the RV (AF), and the RV itself (RV) (Fig. 2a, Fig. 4a-e, Additional file 1: Datasets S6 and S9). By this approach, we were able to determine how the σ_f scores of the proteins that are enriched in RVs transition from healthy fibers to aggregate-containing fibers. We calculated σ_f based on protein abundances for each of these contexts. Even in healthy controls, the σ_f scores of RV-enriched proteins are higher in the muscle context than what we found in the cross-tissue transcriptional analysis for unfolded supersaturation (σ_u) (median Δ: 7.3x, p = 9.6•10^− 4) (comparison p < 2.0•10^− 5) or with the skeletal muscle unfolded supersaturation score (σuts\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\upsigma}_{\mathrm{u}}^{\mathrm{ts}} $$\end{document}) (median Δ: 2.1x, p = 2.2•10^− 6) (Fig. 4a compared with Fig. 2b and Additional file 2: Figure S1). In healthy control (HC) fibers, this result is driven by the higher median aggregation propensity of RV-enriched proteins rather than an increase in abundance (abundance Δ: 1.9x, p = 0.12; ZaggSC\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\mathrm{Z}}_{\mathrm{agg}}^{\mathrm{SC}} $$\end{document} Δ: 7.4x, p = 3.7•10^− 3) (Additional file 2: Figure S2C-D). Fig. 4

Escalating supersaturation in inclusion body myositis. Comparison of folded supersaturation scores (σ_f) for the proteome (Prt) and proteins enriched in rimmed vacuoles (RV) relative to diseased control myofibers. a Healthy control myofiber (HC) (Prt N = 1605, RV N = 47), b control myofibers unaffected in diseased samples (DC) (Prt N = 1988, RV N = 52), c aggregate-containing affected myofibers (AF) (Prt N = 2396, RV N = 52), and (d) rimmed vacuoles (RV) (Prt N = 2104, RV N = 52). e Comparison of the fold difference in median σ_f between RV and Prt. Box plots and statistical tests as in Fig. 1. *p < 0.05, ***p < 0.001, ****p < 0.0001

We found that σ_f increases with the physical proximity to the RVs (DC median Δ: 9.1x, p = 1.2•10^− 4; AF median Δ: 37x, p = 1.2•10^− 9; RV median Δ: 56x, p = 4.3•10^− 10) (Fig. 4b-e). These results are robust against high levels of noise (Additional file 2: Figure S5-S6). In order to determine whether this rise in supersaturation was statistically significant, we performed a simulation of 1,000,000 trials in which we randomly selected 53 proteins to determine how frequently we could achieve the following pattern by chance: 1) elevated supersaturation for RV-associated proteins relative to the proteome in each proteome context (HC, DC, AF, and RV) and 2) a rising median Δ for the supersaturation of RV-associated proteins from HC to DC to AF to RV. In this way, we calculated a significance for achieving this pattern of escalating supersaturation p = 0.011 (see Methods). Given that this analysis included proteins expressed in some contexts but not in others (e.g. present in disease fibers but not healthy control fibers), we confirmed that the results were qualitatively similar when considering only a limited set of proteins expressed in all four contexts and had associated Zyggregator scores available (Additional file 2: Figure S7, Additional file 1: Dataset S9). To further validate these results, we performed these analyses utilizing aggregation predictions from the unfolded state with Z_agg and TANGO [7], which similarly demonstrated a significant escalation in supersaturation (Additional file 2: Figure S8-S9, Additional file 1: Dataset S10).

Like proteins enriched in RVs, proteins enriched in hPAM aggregate-containing fibers also exhibit an escalating σ_f in the sporadic disease context (Additional file 2: Figure S10). The escalation in σ_f is specific for proteins that accumulate in PAMs since proteins that co-aggregate with amyloid plaques (Additional file 2: Figure S11A-E) and neurofibrillary tangles (Additional file 2: Figure S11F-J) in AD do not exhibit escalating σ_f in IBM muscle tissues.

We recently reported a transcriptional suppression of supersaturated proteins and pathways in Alzheimer’s disease [4]. We therefore asked whether a similar phenomenon takes place at the transcriptional and translational levels in IBM. To do so, we determined the proteins differentially expressed in affected fibers (AF) relative to healthy controls in IBM (HC) (Additional file 1: Dataset S11). We found, across independent patient datasets, that 52 proteins are decreased and only one protein, desmin, is increased in affected fibers. Those proteins that are decreased in the surrounding fibers tend to have higher σ_f in healthy controls relative to the rest of the proteome (median Δ: 3.8x, p = 9.8•10^− 5) (Fig. 5a). There are 830 proteins (Prt) in our dataset for which we had ZaggSC\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\mathrm{Z}}_{\mathrm{agg}}^{\mathrm{SC}} $$\end{document} in HC context and abundance values across all four contexts, and of these, only 48 (5.8%) are decreased in abundance in affected fibers. By contrast, of the top 5% most supersaturated proteins in this subset (N = 41) (Top σ_f), seven (17%) are decreased in affected fibers (enrichment p = 0.013) (Fig. 5b). As further validation of this phenomenon, we used RNA sequencing data from healthy muscle and IBM muscle to identify the transcripts of proteins that were downregulated in IBM tissue [15]. The downregulated transcripts correspond to proteins whose supersaturation scores tend to be elevated in healthy controls (median Δ: 2.6x, p = 3.3•10^− 3) (Additional file 2: Figure S12A). There are 1366 transcripts for proteins in this dataset for which we were able to calculate σ_f in HC context. Of the top 5% most supersaturated proteins in this subset (N = 68) (Top σ_f), 15 (22%) are decreased in expression in affected fibers versus 157 (11%) for the proteome (Prt) as a whole (enrichment p = 0.016) (Additional file 2: Figure S12B). Fig. 5

Protein supersaturation is associated with downregulation. In this analysis, only proteins that are detected in HC, DC, AF, and RV, and for which there are defined σ_f scores in HC are included. a Folded supersaturation scores (σ_f) for the proteome (Prt) (N = 830) and proteins downregulated from HC to AF (N = 50). Box plots and statistical tests as in Fig. 1. b Percentage of proteins downregulated in the proteome (Prt) (48/830), proteins enriched in rimmed vacuoles (RV) (0/47), and top 5% most supersaturated proteins (based on HC context) (Top σ_f) (7/41). Significance determined by the Fisher Exact test, with Holm-Bonferroni correction. c Protein abundances in HC, DC, AF, and RV are plotted for the 47 proteins enriched in RVs included in the subset analyzed in this figure. Desmin is highlighted in yellow, the only RV-enriched protein whose abundance is increased between HC and AF. d Protein abundances in HC, DC, AF, and RV are plotted for the 47 proteins with the highest supersaturation scores (top 5%). Desmin again is highlighted in yellow, the only highly supersaturated protein with rising abundance. Those proteins that are significantly downregulated between HC and AF are highlighted in blue. In (c) and (d), the background black line and grey bar represent median and 25th–75th percentile range for the 830 proteins in the proteome in this subset. *p < 0.05, ****p < 0.0001

The RV-enriched proteins (RV) are an exception, as none of them is downregulated at the protein (Fig. 5b) or transcript (Additional file 2: Figure S12B) level. The individual abundance trajectories of these RV-enriched proteins trend towards rising abundances, although only one of the proteins, desmin, significantly rises in abundance between HC and AF (Fig. 5c). By contrast, among the top 5% most supersaturated proteins, there is a trend towards declining abundances, with a disproportionate number of proteins decreasing in abundance significantly in this group (Fig. 5d). These results suggest that supersaturated proteins are typically downregulated to control their abundance, but when they fail to be downregulated, they escape regulation and are more likely to deposit into inclusions.

As mentioned above, only one protein (desmin) is among the top 5% most supersaturated proteins, is an RV-enriched protein, and increases significantly in abundance. Mutations of desmin are associated with desminopathy, and this protein is found to be enriched in aggregate-associated tissues in myotilinopathy, filaminopathy and IBM. Thus, desmin represents the clearest example of an escape protein and is also highly associated with protein misfolding in muscle tissue.

The datasets used in this work and the proteins in each of them are described in Additional file 1: Datasets S1 and S2.

Raw files were converted into the Mascot generic format (MGF) format using Proteome Discoverer 1.4 (Thermo Fischer Scientific, Germany). MGF files were searched against a combined database containing the Swiss-Prot part of the UniProt Knowledgebase (UniProtKB) [35] or Homo sapiens (release 2014/05/28, 20,265 curated entries). For the generation of shuffled decoy entries DecoyDatabaseBuilder was used [29]. Identifications were performed by Mascot 2.5 (Matrixscience Ltd., [28]) with a peptide mass tolerance of 10 ppm, fragment mass tolerance of 0.5 Da, one allowed missed cleavage and carbamidomethylation (C), oxidation (M) as well as phosphorylation (S, T, Y) as variable modifications. Label free relative quantification by spectral counting was performed as described in [14].

We previously reported abundances as spectral counts normalized by the total number of spectral counts in a given sample [14, 20, 25, 26]. Here, we performed an additional normalization step to account for the fact that longer proteins will generate more peptides in mass spectrometry than smaller proteins of the same abundance [17]. Akin to the normalized spectral abundance factor, we divided normalized spectral counts in our data sets by the protein length. We then divided these values by the sum of all such normalized values in a given sample. We then averaged these normalized protein abundances across replicates and log₁₀-transformed these values to arrive at a final abundance value.

Microarray data was obtained from BioGPS pre-processed using gcrma as previously described. For cross-tissue analysis, cell ine and malignant tissue expression levels were excluded. Transcript identifiers were converted to UniProt IDs, with cases of ambiguous conversion or absence of reviewed UniProt IDs excluded from analysis. Values ≤0 were excluded. Expression levels were then log₁₀-transformed then averaged across all values for a given UniProt ID. A similar procedure was done for the skeletal muscle analysis, but limiting it to the two arrays of skeletal muscle data.

Processed RNA sequencing data was obtained, with expression levels reported in FKPM (GEO Datasets GSE102138) [15]. Any values ≤0 were excluded. Identifiers were converted to reviewed UniProt IDs, with ambiguous conversions excluded from further analysis. In cases in which one multiple identifiers mapped to a single UniProt IDs, these FKPM values were averaged. The values were then log₁₀-transformed. Significantly upregulated and downregulated transcripts were identified based on the reported q-values. In cases in which there were multiple q-values associated with a given UniProt ID, the largest q-value was used. The q-values reported were two-tailed, which we converted to one-tailed q-values for the purpose of our analysis. We used a threshold of significance of p < 0.05.

For the human proteome set, we calculated the Z_agg,ZaggSC\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {Z}_{agg}^{SC} $$\end{document} and TANGO scores as previously described [7, 34]. For TANGO, we set the parameters at pH = 7.4, T = 310 K, and ionic strength = 0.1 M. The supersaturation score σ is calculated as the sum

1σ=C+Z\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ \sigma =C+Z $$\end{document}where C is the log₁₀ of the concentration and Z is aggregation propensity score; the concentrations are derived from the protein abundance levels. In each dataset, values were recentered such that the median σ score for each database was 0.

In order to determine vacuole-enriched proteins in the IBM data set, we compared abundance values in the RV dataset to those in the DC dataset. For this analysis, we only included proteins that had a non-zero abundance in both the DC and RV datasets, which constituted a total of 1302 proteins. For these proteins, we performed a one-tailed paired t-test. We then used the Benjamini-Hochberg method to calculate q-values for each of these proteins, using as a threshold of significant q < 0.05, for a False Discovery Rate of 5%.

We performed noise testing to evaluate the robustness of our results for the comparison of supersaturation scores among the IBM data sets, as well as the hPAM data sets. We defined one hundred noise levels on the basis of the standard deviation of a series of Gaussian distributions with mean of 0. The range of standard deviations was log₁₀(1.1) to log₁₀(10.1). At each noise level l, we performed 100 trials t, in which we drew a random number n_{l, t, p} from that the noise level distribution for each of the p proteins in the database. The noise-introduced supersaturation score σ_{p, l, t} was defined as

2σp,l,t=σp+nl,t,p\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\sigma}_{p,l,t}={\sigma}_p+{n}_{l,t,p} $$\end{document}

For trial t of noise level l, the set S_{l, t} of noise values is

3Sl,t=nl,t,1nl,t,2…nl,t,p\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {S}_{l,t}=\left\{{n}_{l,t,1},{n}_{l,t,2},\dots, {n}_{l,t,p}\right\} $$\end{document}

The set m_{l, t} of linear magnitudes of noise for trial t of noise level l is

4ml,t=10lnnl,t,110lnnl,t,2…10lnnl,t,p\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {m}_{l,t}=\left\{{10}^{\ln \left\lceil {n}_{l,t,1}\right\rceil },{10}^{\ln \left\lceil {n}_{l,t,2}\right\rceil },\dots, {10}^{\ln \left\lceil {n}_{l,t,p}\right\rceil}\right\} $$\end{document}

For noise level l, the set M_l of median noise values for its constituent trials is

5Ml=medianml,1medianml,2…medianml,100\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {M}_l=\left\{ median\left({m}_{l,1}\right), median\left({m}_{l,2}\right),\dots, median\left({m}_{l,100}\right)\right\} $$\end{document}

In each Gaussian noise plot, the values plotted on the x-axis were the median of M_l with error bars representing the standard error of the mean as calculated using default settings in the Python package SciPy.

For each trial at each noise level, we determined the sets of noise-modified σ scores for the data sets under consideration. A one-tailed Wilcoxon/Mann-Whitney U test was performed for each of these trials, with multiple hypothesis correction performed based on the same families used for the original analysis, with one difference. At each noise level, the median of the p-values for the 100 trials was plotted with error bars representing the standard error of the mean as calculated using default settings in the Python package SciPy. We performed a one-sided one-sample t-test using the distribution of p-values for a given trial to test the null hypothesis that the mean of these p-values is not significantly less than 0.05. For those cases in which we could not reject the null hypothesis, we plotted the points in grey; otherwise, we plotted the points in color.

For each trial at each noise level, we determined the sets of noise-modified σ scores for the data sets under consideration. The linear difference d_{l, t} between the medians of the supersaturation scores of the control set C_{l, t} and experiment set E_{l, t} being tested at noise level l and trial t is

6dl,t=10medianEl,t−medianCl,t\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {d}_{l,t}={10}^{median\left({E}_{l,t}\right)- median\left({C}_{l,t}\right)} $$\end{document}

At noise level l, we plotted the median of set {d_{l, 1}, d_{l, 2}, …, d_{l, 100}} with error bars representing the standard error of the mean as calculated using default settings in the Python package SciPy. We performed a one-sided one-sample t-test using the distribution of fold change values for a given trial to test the null hypothesis that the mean of these fold changes is not significantly greater than 1. For those cases in which we could not reject the null hypothesis, we plotted the points in grey; otherwise, we plotted the points in color.

In Fig. 5b and Additional file 2: Figure S12B, the Fisher exact test is used to calculate enrichment of data sets for particular categories of proteins.

To test the significance of our observations of rising supersaturation (Fig. 4, Additional file 2: Figures S7–11) we used a simulation. The null hypothesis was that it would arise by chance that 1) the median Δ > 0 for a set of proteins of interest in each context and 2) median Δ of those proteins would rise successively from HC to DC to AF to RV contexts. To test this, we performed the following procedure K times, where K =1,000,000. For each trial k, we randomly selected N proteins from the proteome (where N is equal to the number of proteins of interest, for instance 53 in the case of RV-enriched proteins or 51 in the case of hPAM-enriched proteins). When selecting N, we used the total number of proteins meeting a particular criterion, even if a smaller number of those proteins was actually present in the original dataset. For these N proteins, D is the set of median Δ compared to the proteome for each of the four contexts: 7D≡med∆HCmed∆DCmed∆AFmed∆RV\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ D\equiv \left\{ med{\Delta }_{HC}, med{\Delta }_{DC}, med{\Delta }_{AF}, med{\Delta }_{RV}\right\} $$\end{document}

If the supersaturation rose successively at each from HC to DC to AF to RV, and median Δ > 0 in each context, we assigned a score E_k of one; otherwise, we assigned a score E_k of zero. We then summed this score over the 1,000,000 trials. 8D=med∆HCmed∆DCmed∆AFmed∆RV\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ D=\left\{ med{\Delta }_{HC}, med{\Delta }_{DC}, med{\Delta }_{AF}, med{\Delta }_{RV}\right\} $$\end{document}9Ek=1,ifminD>0andmed∆RV>med∆AF>med∆DC>med∆HC0,otherwise\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {E}_k=\left\{\begin{array}{c}1,\kern0.5em if\min (D)>0\ and\ med{\Delta }_{RV}> med{\Delta }_{AF}> med{\Delta }_{DC}> med{\Delta }_{HC}\\ {}0,\kern26em otherwise\end{array}\right.\kern0.5em $$\end{document}

We estimated the significance of the escalation in supersaturation as follows: 10E=E1,…,EK\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ E=E1,\dots, EK $$\end{document}11p=∑k=1KEkK\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ p=\sum \limits_{k=1}^K\frac{E_k}{K} $$\end{document}12p=∑k=1KEkK\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ p=\sum \limits_{k=1}^K\frac{E_k}{K} $$\end{document}

In order to test the isolated contribution of escalating median Δ, we removed the constraint of median Δ > 0, and calculated a score Ekr\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {E}_k^r $$\end{document}: 13Ekr=1,ifmed∆RV>med∆AF>med∆DC>med∆HC0,otherwise\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {E}_k^r=\left\{\begin{array}{c}1,\kern0.5em ifmed{\Delta }_{RV}> med{\Delta }_{AF}> med{\Delta }_{DC}> med{\Delta }_{HC}\\ {}0,\kern17.75em otherwise\end{array}\right. $$\end{document}14p=∑k=1KEkK\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ p=\sum \limits_{k=1}^K\frac{E_k}{K} $$\end{document}

We considered all cases analyzed by our original constraints on family for the purpose of multiple hypothesis correction and all cases analyzed by the relaxed criteria a separate family. Multiple hypothesis correction was performed using the Holm-Bonferroni method. P-values for both constraints are reported in Additional file 1: Dataset S12.

To test the significance of differences in median Δ between different contexts (Figs. 2, 3 and 4), we used a simulation. The null hypothesis was that the difference in median Δ (∆_∆), of at least the magnitude reported would arise by chance. The reported difference in median Δ we refer to as ∆∆0\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\Delta }_{\Delta }^0 $$\end{document}. To test this, we performed the following procedure K times, where K =1,000,000. For each trial k, we randomly selected N proteins from the proteome by the same procedure as above for escalating supersaturation. For these N proteins, we calculated the median Δ in contexts C₁ and C₂. Note that we performed this analysis in a one-tailed fashion. 15SΔΔ=ΔΔ1…ΔΔKΔΔ=ΔΔ1…ΔΔK\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {S}_{\varDelta_{\varDelta }}=\left\{{\varDelta}_{\varDelta}^1,\dots, {\varDelta}_{\varDelta}^K\right\}{\varDelta}_{\varDelta }=\left\{{\varDelta}_{\varDelta}^1,\dots, {\varDelta}_{\varDelta}^K\right\} $$\end{document}where 16∆∆k=med∆2k−med∆1k\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {\Delta }_{\Delta }^k= med{\Delta }_2^k- med{\Delta }_1^k $$\end{document}

We assigned a score E_k to each trial and from all the trials together derived a p-value, as follows: 17Ek=1,∆∆k>∆∆00,otherwise\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ {E}_k=\left\{\begin{array}{c}1,\kern0.5em {\Delta }_{\Delta }^k>{\Delta }_{\Delta }^0\\ {}0,\kern0.75em otherwise\end{array}\right. $$\end{document}

18p=∑k=1KEkK\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ p=\sum \limits_{k=1}^K\frac{E_k}{K} $$\end{document}

We considered all cases analyzed in this fashion as a single family. Multiple hypothesis correction was performed using the Holm-Bonferroni method. P-values are reported in Additional file 1: Dataset S12.

In order to perform adequate multiple hypothesis correction while avoiding increasing Type II error by overcorrecting, it was necessary to group our results into a series of families on which multiple hypothesis correction would be performed meaningfully. We used the following principles to help divide the analyses in these studies into a set of coherent families. Except when they were being compared directly, hPAM and IBM data sets were considered part of separate families. IBM families were organized cross data subsets (that is, HC, DC, AF, and RV included in the same family). hPAM families were organized in three families: 1) HC, 2) DC, and 3) AF. This was organized in this way because there were multiple individual hPAMs, but analyses for the composite group of hPAM aggregate-enriched proteins could only be performed logically on the HC dataset as the other data sets were disease-specific. Analyses using σ_u were considered distinct from analyses using σ_f. All σ_u analyses were considered as part of a single family. Among IBM data sets, we performed a series of analyses in which we compared σ_f levels between the proteome and particular subsets of proteins (RV-enriched, hPAM-enriched, plaque-enriched, NFT-enriched) across the four IBM data sets (HC, DC, AF, RV). We considered analyses involving each of these subsets as separate families. Additional file 1: Dataset S12 shows a summary of all statistical tests performed in this analysis, and groups those tests by their respective families.

Patients provided informed consent. Study protocols were approved by the local ethics committee (reg. Number 3882–10) at Ruhr-University Bochum, Bochum, Germany. For each patient 250,000 μm² of HC, DC, AF or RV tissue was collected by LMD (LMD 6500, Leica Microsystems, Wetzlar, Germany). Sample lysis and digestion were carried out as previously described [25]. Briefly, samples were lysed with formic acid (98–100%) for 30 min at room temperature (RT), followed by a sonication step for 5 min (RK31, BANDELIN electronic, Berlin, Germany). Samples were kept frozen at − 80 °C until digestion.

Prior to digestion the formic acid was removed and the collected samples were digested in 50 mM ammonium bicarbonate at pH 7.8. Samples were reduced and alkylated by adding dithiothreitol and iodoacetamide. Trypsin (Serva) was added to a final concentration of 1 μg. Digestion was carried out overnight at 37 °C and stopped by adding TFA to acidify the sample. Samples were purified using OMIX C18 Tips (Varian, Agilent Technologies, Böblingen, Germany) completely dried vacuum and again solved in 63 μl 0.1% TFA, as described in [25].

Sixteen microliter per sample were analysed by nano-liquid chromatography tandem mass spectrometry (nanoLC-ESI-MS/MS). The nano high performance liquid chromatography (HPLC) analysis was performed on an UltiMate 3000 RSLC nano LC system (Dionex, Idstein, Germany) as described in [26]. Peptides were separated with a flow rate of 400 nl/min using a solvent gradient from 4 to 40% B-solvent for 95 min. Washing of the column was performed for 5 min with 95% B-solvent and was then returned to 4% B-solvent. The HPLC system was online-coupled to the nano electrospray ionization (ESI) source of an Orbitrap elite mass spectrometer (Thermo Fisher Scientific, Germany). Mass spectrometric measurements were performed as previously described [14].