A Plasma Proteomic Signature of Skeletal Muscle Mitochondrial Function

Although mitochondrial dysfunction has been implicated in aging, physical function decline, and several age-related diseases, an accessible and affordable measure of mitochondrial health is still lacking. In this study we identified the proteomic signature of muscular mitochondrial oxidative capacity in plasma. In 165 adults, we analyzed the association between concentrations of plasma proteins, measured using the SOMAscan assay, and skeletal muscle maximal oxidative phosphorylation capacity assessed as post-exercise phosphocreatine recovery time constant (τPCr) by phosphorous magnetic resonance spectroscopy. Out of 1301 proteins analyzed, we identified 87 proteins significantly associated with τPCr, adjusting for age, sex, and phosphocreatine depletion. Sixty proteins were positively correlated with better oxidative capacity, while 27 proteins were correlated with poorer capacity. Specific clusters of plasma proteins were enriched in the following pathways: homeostasis of energy metabolism, proteostasis, response to oxidative stress, and inflammation. The generalizability of these findings would benefit from replication in an independent cohort and in longitudinal analyses.


Introduction
Mitochondrial oxidative phosphorylation is the major source of energy production for all cellular functions [1]. Accordingly, impaired mitochondrial function, one of the hypothetical mechanisms that drives the aging process [2], has been associated with the development of phenotypical and functional manifestations of aging and with age-related diseases [3][4][5]. Mitochondrial oxidative capacity measured in vivo in skeletal muscle declines with aging and is associated with lower walking speed, muscle strength, and physical activity independent of age [6][7][8], as well as with chronic inflammation [9]. Consistent with the decline of oxidative capacity with aging, discovery proteomic studies in skeletal muscle from healthy individuals over a wide age-range have shown a substantial decline of mitochondrial proteins with aging, including proteins of electron transport chain complexes (ETC), enzymes of the Krebs cycle as well as structural proteins [10]. In a recent study, aimed to define the proteomic signature of mitochondrial oxidative capacity in skeletal muscle, we identified muscle proteins that were differentially represented in individuals with higher and lower oxidative capacity, measured by phosphorus magnetic resonance spectroscopy ( 31 P-MRS) [11]. As expected, proteins overrepresented in muscle with higher oxidative capacity were enriched for pathways connected with mitochondrial metabolism and translation within mitochondria. Unexpectedly, we also found highly significant enrichment for mRNA processing/alternative splicing pathways, although this finding remains unexplained [11].
Although the identification of specific skeletal muscle proteins that are associated with a direct measure of mitochondrial function is important to gain insight into mechanisms of mitochondrial decline, it has limited clinical use because it is invasive and requires muscle biopsy specimens. In this study, we hypothesized that reduced oxidative capacity in skeletal muscle may be reflected by characteristic changes in circulating proteins and, therefore, we searched for a proteomic signature of muscular mitochondrial oxidative capacity in plasma.

Results
Demographic characteristics of 165 study participants are displayed in Table 1. Participants were 45% female, in the age range of 22-93 years (average 57.7 ± 20), and mostly Caucasian. Out of the 1301 SOMAmers analyzed, we identified 87 proteins significantly (p < 0.01) associated with muscle mitochondrial oxidative capacity (τ PCr ), adjusting for age, sex, and PCr depletion ( Figure 1, Table 2). Sixty proteins were negatively associated with τ PCr , therefore positively associated with a better oxidative capacity, while 27 proteins were associated with a poorer oxidative capacity. The top 10 proteins most strongly associated with τ PCr in the first model were endothelial cell-selective adhesion molecule (ESAM), insulin-like growth factor binding protein 3 (IGFBP-3), contactin 2 (CNTN2), p-selectin (coded by the SELP gene), proto-oncogene tyrosine-protein kinase Fyn (FYN), lactoperoxidase (PERL), dermatopontin (DERM), C-X-C motif chemokine ligand 16 (CXCL16), Ras-related C3 botulinum toxin substrate 3 (RAC3), and tyrosine-protein kinase Lyn (LYN) ( Table 3, Model 1). After adding race and BMI to the multivariable model, 62 proteins were significantly associated with τ PCr . The top 10 significant proteins were similar to those resulting from model 1, with the new appearance of kallikrein 11 (coded by the KLK11 gene), follistatin-related protein 3 (FSTL3), cathepsin F (CATF), and IGFBP-6 ( Table 3, Model 2). pH at baseline 7.07 (0.03, 6.97, 7.14) Figure 1. Volcano plot of the association between protein concentrations (gene ID in the plot) and mitochondrial oxidative capacity (τ PCr ), adjusted for age, sex, and phosphocreatine depletion. Negative beta: proteins were positively associated with a better oxidative capacity; positive beta: proteins were associated with a poorer oxidative capacity. When analyzing the patterns of functional enrichment, many gene ontology (GO) terms were significantly enriched among the 87 proteins significantly correlated with mitochondrial oxidative capacity ( Table 4). The 87 significant proteins were enriched for genes in the oxidative stress, inflammation, metabolism regulation, and proteostasis pathways ( Figure 2). Table 3. Top-10 most significant SOMAmers associated with τ PCr , an inverse measure of mitochondrial oxidative capacity, adjusting for age, sex, and amount of phosphocreatine depletion (Model 1); top-10 most significant SOMAmers associated with τ PCr adjusting for age, sex, phosphocreatine depletion, BMI, and race (Model 2).

Discussion
In this study we characterized the plasma proteomic profile associated with skeletal muscle oxidative capacity. Using the 1.3 k SOMAscan assay, we found that 87 proteins were associated with mitochondrial oxidative capacity, as measured by 31 P-MRS on the quadriceps muscle. Three of these proteins (heat shock protein family D member 1, HSP 60, the serine protease HTRA2, and cyclophilin F, coded by the PPIF gene) are defined as mitochondrial proteins, and many of the other proteins identified are relevant to mitochondrial functions, or to processes that have been linked with impaired mitochondria. This proteomic signature is reflective of the deranged metabolic mechanisms that are either causes or consequences of impaired mitochondrial function, independent of chronological age. Importantly, although a few proteins found in our analysis are specific to either the muscle (e.g., myostatin) or the blood tissues, most of the proteins are ubiquitous and can be detected in different tissues. From the available data, it is difficult to identify which tissues the plasma proteins represent. Whether the signature identified is a direct marker of muscle function or a marker of a more generalized energetic alteration reflected in plasma proteins cannot be determined.
The most relevant proteins to the association with oxidative capacity are summarized based on our bioinformatic analyses of the enriched processes.

Energy Metabolism
Multiple proteins involved in the regulation of metabolism showed a strong association with muscle oxidative capacity. IGFBP-3 and IGFBP-6 are the carriers of insulin-like growth factor-1 (IGF-1), the primary effector of growth hormone. In addition to its insulin-like functions, IGF-1 stimulates cell growth and proliferation in most tissues of the body and inhibits apoptosis [12]. IGFBP-3 and IGFBP-6 prolong the half-life of IGF-1 and regulate the growth promoting effects of IGF-1, altering its interaction with cell surface receptors. IGFBP-3 also exhibits IGF-independent antiproliferative and apoptotic effects. Several other proteins significantly associated with τ PCr are involved in apoptosis, among which caspase 3 and the protease HTRA2. This finding is particularly interesting because increased apoptosis signaling has been implicated in the pathogenesis of age-related sarcopenia [13].
Decreased plasmatic levels of the adipose-derived hormone leptin were found associated with better oxidative capacity. Leptin is one of the main regulators of feeding and energy balance, connecting changes in energy stores to a set of adaptive physiologic responses. Leptin regulates numerous physiologic processes, among which feeding behavior, metabolism, thermogenesis, immune function, and the neuroendocrine axis [14]. Increased leptin has been related to obesity, inflammation, hypoxia, and has been implicated in ROS generation [15,16]. Any of these mechanisms could underlie the association showed by our data.
Plasma concentrations of the AMP-activated protein kinase (AMPK) a2b2g1 complex were positively associated with better mitochondrial oxidative capacity in this analysis. AMPK is the master sensor of energy metabolism and its activity is principally modulated by the AMP/ATP ratio and, to a lesser extent, by the ADP/ATP ratio, which are direct biomarkers of the status of energy availability. AMPK responds to reduced energy availability by downregulating activities that are energy demanding such as protein and lipid synthesis and cell cycle, and improves energy production through increased catabolism [17]. Moreover, AMPK modulates fundamental mitochondrial processes such as biogenesis, fission, and autophagy, and promotes mitochondrial health [17].
Aldolase A and pyruvate kinase (PKM2) are two key glycolytic enzymes that were found positively associated with oxidative capacity. Other than allowing production of ATP through glycolysis, aldolase A also modulates the myocyte's shape and contractility, and its absence has been related to metabolic myopathy [18,19]. Aldolase A is prominently expressed in skeletal muscle.
Nicotinamide phosphoribosyltransferase (PBEF, coded by NAMPT) catalyzes the rate limiting reaction of the mammalian nicotinamide adenine dinucleotide (NAD + ) salvage pathway. NAD + is an essential cofactor regulating several metabolic processes such as glycolysis, fatty acid oxidation, the tricarboxylic acid cycle, and oxidative phosphorylation, but also mitochondrial biogenesis [20]. NAD + levels have been implicated in aging, age-related diseases, and longevity [21]. PBEF is released by multiple cell types including myocytes, and it was associated with lower oxidative capacity in our analysis. Arguably, this could reflect either a spilling from damaged muscle cells, or a compensation for a dietary deficiency of NAD + precursors.
Enrichment for the Kit signaling pathways, represented by proteins such as FYN, VAV, LYN, protein kinase C alpha (PKC-A), and proto-oncogene tyrosine-protein kinase Src (SRC), was found. Although it is best known for its role in hematopoietic stem cell differentiation, Kit plays an important role in the regulation of mitochondrial function and energy expenditure [22]. Kit promotes the expression of the peroxisome proliferator-activated receptor-γ (PPARγ) coactivator-1α (PGC-1α), the master regulator of mitochondrial biogenesis and function [23]. PGC-1α has been found to specifically promote mitochondrial biogenesis in skeletal muscle, and its deficiency has been correlated with metabolic derangements and muscle dysfunction [24].
Another pathway significantly enriched was that of the advanced glycation end products (AGEs). AGEs are a heterogeneous group of bioactive molecules formed by the nonenzymatic glycation of proteins, lipids, and nucleic acids. AGEs accumulate with aging in several tissues and contribute to oxidative stress, chronic inflammation, and are implicated in the pathogenesis of cardiovascular diseases and chronic kidney disease [25].

Proteostasis
In this analysis, we found that the plasma levels of several proteins implicated in proteostasis were significantly different across mitochondrial function. Heat shock protein 27 (HSP 27), a low-weight molecular chaperone that maintains denatured proteins in a folding-competent state, in skeletal muscle plays an important role in stress resistance and actin organization [26]. HSP 60 and cyclophilin F participate in mitochondrial import and correct folding of proteins. In addition, cyclophilin F is a major component of the mitochondrial permeability transition pore (MPTP) highly involved in connecting mitochondrial metabolism and apoptosis [27]. Cyclophilin D, coded by PPID, is another enzyme assisting and accelerating the correct folding of proteins. HSP 27, HSP 60, cyclophilin F, and cyclophilin D plasma concentrations were positively associated with better oxidative capacity.
DnaJ homolog subfamily B member 1 (DNJB1) interacts with HSP 70 and can stimulate its ATPase activity, facilitating ATP hydrolysis and protein folding. Previous studies found that decreased HSP 70 response was associated with age-related functional impairments in skeletal muscle [28], and that overexpression of HSP 70 in transgenic mice conveyed protection against age-related dysfunction [29], supporting the concept of chaperones as essential molecules to restore the normal cell function after an insult [30]. DNJB1 levels were associated with better mitochondrial oxidative capacity.
Interestingly, levels of two proteins involved in protein degradation and turnover were associated with poorer oxidative capacity: cathepsin F, a major component of the lysosomal system, and ubiquitin-conjugating enzyme E2 G2 (UB2G2), which targets abnormal proteins and catalyzes the attachment of ubiquitin.

Inflammation and Response to Reactive Oxygen Species
A strong connection between mitochondrial impairment and chronic inflammation has been gaining increasing attention [31]. Dysfunctional mitochondria produce an excessive amount of reactive oxygen species (ROS), which trigger inflammation both directly and through oxidative damage to proteins, lipids, and nucleic acids [32]. Furthermore, products of damaged mitochondria released in the extracellular space act as damage-associated molecular pattern (DAMP) agents, activating the immune response [32]. Among the proteins associated with poorer mitochondrial function in our analysis, many were markers of inflammation. Proteins involved in the signaling of cytokines, in chemotaxis, and in the response to oxidative stress, showed a marked prevalence in the relevant clusters identified. Several proteins represented either components or activators of the mitogen-activated protein kinase (MAPK) signal transduction pathway, such as MAPK14, MAPK2, tyrosine-protein kinase Lyn (LYN), sphingosine kinase 1 (Q9NYA1), PKC-A, PKC-B, SLAM family member 5 (SLAF5), and fibroblast growth factor receptor 4 (FGFR4). The MAPK signaling cascade regulates survival and death, proliferation, and differentiation of cells, and is an important activator of the inflammatory response. A strict control is therefore crucial, and its disruption has been linked to the development of many diseases [33].
Many proinflammatory molecules have been identified among the members of the senescence -associated secretory phenotype (SASP). SASP is the secretome of senescent cells that contains hundreds of compounds, some of which have not been yet identified [34]. There is some evidence that the SASP can induce senescence in surrounding cells and cause damage accumulation both in cells and intercellular matrix, events that may contribute to the phenotypes of aging as well as to many chronic diseases [35,36]. It is well known that the energic crisis caused by mitochondrial dysfunction can induce cellular senescence, and the finding that common plasma biomarkers of senescence are dysregulated according to mitochondrial function is not fully surprising [37].

Growth Differentiation Factor 15
Previous studies have proposed growth differentiation factor 15 (GDF15) as a biomarker of mitochondrial dysfunction in aging and several age-related diseases [37]. Increased blood levels of GDF15 have been observed in aging and in mitochondrial disease, and this protein has been related to cardiovascular and brain disease [38,39]. GDF15 has been identified among the molecules expressed by senescent cells which constitute the SASP [34]. Although GDF15 did not appear among the most significant proteins in our analysis, it showed an association (p = 0.026) with τ PCr (Model 1); as expected, poor oxidative capacity was associated with higher levels of GDF15 (β coefficient = 0.011).

Participants
This study was conducted in 165 community-dwelling volunteers participating in the Baltimore Longitudinal Study of Aging (BLSA, N = 76) and the Genetic and Epigenetic Signatures of Translational Aging Laboratory Testing (GESTALT, N = 89) studies.
The BLSA is a prospective open cohort study that has continuously enrolled participants aged 20 and older since 1958. The GESTALT study started in April 2015, aimed at discovering new and sensitive molecular biomarkers of aging in different cell types. Volunteers are eligible to participate in BLSA and GESTALT if they meet strict healthy inclusion criteria, where participants are free of major pathologies (with the exception of controlled hypertension) as well as functional and cognitive impairments at enrollment, and are followed for life regardless of changes in health and functional status (Supplementary Appendix A). This analysis was performed on samples collected during visits in which participants met the inclusion criteria.
All assessments, which took place at the Clinical Research Unit of the Intramural Research Program of the National Institute on Aging, National Institutes of Health (NIH) during a 2.5-3.5-day visit, were performed by certified nurse practitioners and certified technicians according to standardized procedures. The protocol for both studies was approved by the NIH Intramural Institutional Review Board (BLSA (03AG0325) and GESTALT (15AG0063) were approved on 12 May 2020). After receiving detailed descriptions of the procedures at every visit, all subjects provided written informed consent.
Demographic and health characteristics were assessed either through self-report questionnaires or using standard criteria and algorithms [40]. Body weight was measured in kilograms using a calibrated scale to the nearest 0.1 kg. Body height was measured in centimeters by a stadiometer to the nearest 0.1 cm [41]. Body mass index (BMI) was calculated by dividing body weight by the square of height in meters.

Proteomic Assessment
Plasma proteins were measured using overnight fasted plasma that was collected at a resting state and subsequently stored at −80 • C. Discovery proteomics was performed using the 1.3 k SOMAscan Assay (SomaLogic, Inc.; Boulder, CO, USA) at the Trans-NIH Center for Human Immunology and Autoimmunity, and Inflammation (CHI), National Institute of Allergy and Infectious Disease, National Institutes of Health (Bethesda, MD, USA). SOMAmer reagents are individually generated via an iterative process called SELEX (Systematic Evolution of Ligands by EXponential enrichment), which consists of affinity selection cycles aimed at increasing the specificity and avidity of oligomers to a target protein epitope. As a result, SOMAmer reagents are designed to be highly specific and sensitive. For a quantitative assessment of variability in the SOMAscan assay, see [42,43]. A discussion of caveats and limitations is provided in [44].
Of the 1322 SOMAmer reagents included in this version of the kit, 12 hybridization controls, four viral proteins (HPV type 16, HPV type 18, isolate BEN, isolate LW123), and five SOMAmers that were reported to be nonspecific (P05186; ALPL, P09871; C1S, Q14126; DSG2, Q93038; TNFRSF25, Q9NQC3; RTN4) were removed, leaving 1301 SOMAmer reagents for the final analysis. There are 46 SOMAmer reagents that are documented to target multicomplex proteins of two or more unique proteins (UniProt IDs). Conversely, there are 49 UniProt IDs that are measured by more than one SOMAmer reagent. The full list of SOMAmer reagents and their protein targets is provided as Supplementary Information. Thus, the 1301 SOMAmer reagents collectively target 1297 UniProt IDs. Of note, there are four proteins in the final protein panel that are rat homologues (P05413; FABP3, P48788; TINNI2, P19429; TINNI3, P01160; NPPA) of human proteins.
The experimental process for proteomic assessment and data normalization has been previously described [42]. The data reported are SOMAmer reagent abundance in relative fluorescence units (RFU). The abundance of the SOMAmer reagent represents a surrogate of protein concentration in the plasma sample.
Data normalization was conducted in three stages. First, hybridization control normalization removed individual sample variance on the basis of signaling differences between microarray or Agilent scanner. Second, median signal normalization removed intersample differences within a plate due to technical differences such as pipetting variation. Last, calibration normalization removed variance across assay runs. Furthermore, there was an additional interplate normalization process that utilized a CHI calibrator of pooled plasma from healthy subjects that allowed normalization across all experiments conducted at the CHI laboratory [42]. An interactive Shiny web tool was used during the CHI QC process [45].

Phosphorus Magnetic Resonance Spectroscopy
Using a 3T MR scanner (Achieva, Philips Healthcare, Andover, MA, USA), in vivo 31 P-MRS measurements of the concentrations of the phosphorus-containing metabolites phosphocreatine (PCr), inorganic phosphate (Pi), and ATP were obtained from the vastus lateralis muscle of the left thigh, following a standardized protocol described previously [7,46]. Participants were positioned supine on the bed of the scanner, with a foam wedge placed underneath the knee to induce slight flexion, and with ankles, thighs, and hips secured with straps to reduce movement during exercise. Participants were required to perform a ballistic knee extension exercise inside the magnet with their left leg, while a resistance was added by foam pads placed above the left leg in order to enhance the intensity of the exercise. A series of pulse-acquire 31 P spectra were obtained before, during, and after the exercise, which had an average duration of 30 s, with a repetition time of 1.5 s, using a 10-cm 31 P-tuned surface coil (PulseTeq, Surrey, UK) fastened above the left thigh. An example of the acquired spectra is represented in Figure 3 [7]. Signals were averaged over four successive acquisitions for signal-to-noise ratio enhancement, so that the data consisted of 75 spectra obtained with a temporal resolution of 6 s. The duration of exercise was optimized by consistently requiring a depletion in PCr of 50-67% relative to initial baseline values, in order to standardize the measure of oxidative function across different subjects and to provide sufficient dynamic range to fit the PCr recovery curve. Whenever PCr depletion did not reach the threshold of 33%, data collected were excluded from further analysis. If intramuscular acidosis, defined as intracellular pH lower than 6.8, was detected at the end of the exercise, the test was repeated at a lower intensity after waiting for the participant to return to a resting condition [47].
The pH was determined according to the chemical shift of Pi relative to PCr [48]. Spectra were processed with jMRUI software (version 5.2, MRUI Consortium), and metabolite concentrations were calculated by nonlinear least squares fitting implemented through AMARES [49,50]. Post-exercise PCr recovery rates were calculated by fitting time-dependent changes in PCr peak area to the monoexponential recovery function: where PCr(0) is the end-of-exercise PCr signal area (i.e., the PCr signal area at the beginning of the recovery period), ∆PCr is the decrease in signal area from its pre-exercise baseline value, averaged from the multiple baseline scans, to PCr(0) resulting from in-magnet exercise, and τ PCr is the PCr exponential recovery time constant, measured in seconds [7]. This time constant is inversely proportional to the maximum in vivo oxidative capacity of skeletal muscle, with longer τ PCr reflecting slower recovery and therefore lower oxidative capacity [51]. Since the energy demands during post-exercise PCr resynthesis are minimal, 1/τ PCr reflects the maximum mitochondrial ATP production rate [7,[52][53][54]. ATPmax was finally estimated as [PCr baseline ] * (1/τ PCr ) [55].

Statistical Analysis
Protein RFU values were converted to a z-score after natural log-transformation. Association of each protein with mitochondrial oxidative capacity (τ PCr ) was assessed using linear regression models adjusted for age, sex, and amount of PCr depletion. A second model was examined with further adjustments for race (white, black, other), and BMI. The analyses were performed using RStudio (v. 1.2.1335). A nominal p value of 0.01 was considered statistically significant.

Enrichment Analysis
To evaluate whether among the proteins significantly correlated with τ PCr appeared to be enrichment of specific biological processes or molecular functions, a gene enrichment analysis was run on the 87 plasma proteins significantly associated with τ PCr . For this purpose, the bioinformatic tool ClueGO was used [56]. ClueGO permitted the identification of functional gene ontologies (GO) and pathways associated with the most significant proteins in our analysis. To visualize the enriched pathways, the Cytoscape (v. 3.8.0) plug-in was used. The resulting enriched pathways were displayed as a network, in which each pathway was represented as a node and the edges connected similar pathways; pathways were functionally grouped, and the similarity between pathways was determined by kappa statistics [56]. Enrichment significance was represented by node size, and node colors were used to differentiate pathway clusters. The relevant pathways were filtered for > 4 genes after Bonferroni correction.

Limitations and Conclusions
The SOMAscan platform used in this study assesses 1301 proteins, which represent only a fraction of the proteins that are potentially important for aging or mitochondrial function. Furthermore, the SOMAscan technology does not provide an absolute measure of protein abundance, which makes it difficult to compare associations with mitochondrial function across proteins. In addition, while aptamers are designed to detect proteins in their native conformation, there is also a possibility of cross-reactivity between similar proteins, which can limit the accuracy of SOMAscan for proteins with high sequence homology [57].
Importantly, the list of proteins identified do not completely reflect those of the proteins that have been found to be associated with aging. This could be due to the fact that while mitochondrial oxidative capacity declines with aging, the rate of decline is highly heterogeneous across individuals, probably because of the effect of genetic heterogeneity and subclinical and clinical pathology. Due to the limited sample size and number of proteins studied we could not determine an exhaustive proteomic signature and fully understand what biological pathways are shared between aging and mitochondrial function.
The inclusion criteria of BLSA and GESTALT ensure that study participants are exceptionally healthy, and therefore our findings may not be generalizable to a population affected by substantial morbidity or disability. In addition, the cross-sectional nature of our analysis makes it impossible to discriminate whether the proteins associated with differential mitochondrial function in this study represent causes or consequences. Further longitudinal studies should better disentangle the mechanisms underlying the associations identified.
Finally, although 31 P-MRS has been long considered a procedure able to provide a measure of mitochondrial oxidative capacity generalizable to many tissues, the measure of mitochondrial oxidative capacity used in this study is relative to the quadriceps muscle and it is possible that the energetic status of this muscle may not be representative of the energetic status in other muscle groups or tissues. Hence, the interference of other tissues may have reduced the signal-to-noise ratio weakening the results of this study.
In conclusion, mitochondrial oxidative capacity of skeletal muscle was associated with specific clusters of plasma proteins in this study, mainly representing the following pathways: homeostasis of energy metabolism, protein turnover, and inflammation. These findings need to be replicated in an independent population, possibly with longitudinal data and alternative measurements of mitochondrial function, before the results can be used to develop a clinical tool to assess mitochondrial function using the blood proteome.