Methylome-proteome integration after late-life voluntary exercise training reveals regulation and target information for improved skeletal muscle health

Abstract

Exercise is a potent stimulus for combatting skeletal muscle ageing. To study the effects of exercise on muscle in a preclinical setting, we developed a combined endurance-resistance training stimulus for mice called progressive weighted wheel running (PoWeR). PoWeR improves molecular, biochemical, cellular and functional characteristics of skeletal muscle and promotes aspects of partial epigenetic reprogramming when performed late in life (22-24 months of age). In this investigation, we leveraged pan-mammalian DNA methylome arrays and tandem mass-spectrometry proteomics in skeletal muscle to provide detailed information on late-life PoWeR adaptations in female mice relative to age-matched sedentary controls (n = 7-10 per group). Differential CpG methylation at conserved promoter sites was related to transcriptional regulation genes as well as Nr4a3, Hes1 and Hox genes after PoWeR. Using a holistic method of -omics integration called binding and expression target analysis (BETA), methylome changes were associated with upregulated proteins related to global and mitochondrial translation after PoWeR (P = 0.03). Specifically, BETA implicated methylation control of ribosomal, mitoribosomal, and mitochondrial complex I protein abundance after training. DNA methylation may also influence LACTB, MIB1 and UBR4 protein induction with exercise - all are mechanistically linked to muscle health. Computational cistrome analysis predicted several transcription factors including MYC as regulators of the exercise trained methylome-proteome landscape, corroborating prior late-life PoWeR transcriptome data. Correlating the proteome to muscle mass and fatigue resistance revealed positive relationships with VPS13A and NPL levels, respectively. Our findings expose differential epigenetic and proteomic adaptations associated with translational regulation after PoWeR that could influence skeletal muscle mass and function in aged mice. KEY POINTS: Late-life combined endurance-resistance exercise training from 22-24 months of age in mice is shown to improve molecular, biochemical, cellular and in vivo functional characteristics of skeletal muscle and promote aspects of partial epigenetic reprogramming and epigenetic age mitigation. Integration of DNA CpG 36k methylation arrays using conserved sites (which also contain methylation ageing clock sites) with exploratory proteomics in skeletal muscle extends our prior work and reveals coordinated and widespread regulation of ribosomal, translation initiation, mitochondrial ribosomal (mitoribosomal) and complex I proteins after combined voluntary exercise training in a sizeable cohort of female mice (n = 7-10 per group and analysis). Multi-omics integration predicted epigenetic regulation of serine β-lactamase-like protein (LACTB - linked to tumour resistance in muscle), mind bomb 1 (MIB1 - linked to satellite cell and type 2 fibre maintenance) and ubiquitin protein ligase E3 component N-recognin 4 (UBR4 - linked to muscle protein quality control) after training. Computational cistrome analysis identified MYC as a regulator of the late-life training proteome, in agreement with prior transcriptional analyses. Vacuolar protein sorting 13 homolog A (VPS13A) was positively correlated to muscle mass, and the glycoprotein/glycolipid associated sialylation enzyme N-acetylneuraminate pyruvate lyase (NPL) was associated to in vivo muscle fatigue resistance.

Keywords: BETA; DNA methylation; ageing; mitoribosome; proteomics.

Figure 1. Methylation of conserved CpGs is modified by late‐life PoWeR

A, molecular function analysis of differentially methylated (DM) CpGs in promotor regions after PoWeR (compared to sedentary, SED). B, fold difference with PoWeR compared SED for transcriptional regulators Hes1 (cg11825649, P = 0.008; cg13065537, P = 0.026; cg07747635, P = 0.033), Hoxc10 (cg08900043.2, P = 0.021; cg15869291, P = 0.027; cg08900043.1, P = 0.038), Hoxc11 (cg21777948, P = 0.0006; cg05604535, P = 0.043; cg208086387, P = 0.046), and Tead1 (cg15462582, P = 0.015; cg00553886, P = 0.030) with multiple DM promotor regions and Hoxa9 (cg13036782, P = 0.024) and Hoxd3 (cg21235151, P = 0.034) with a single DM promoter region. C, Nr4a3 DM CpG promotor regions (cg20991883, P = 0.0002; cg19178879, P = 0.0004; cg11403906, P = 0.002) fold difference from SED. Individual fold change data points presented for PoWeR (n = 10) relative to the average of SED (n = 9).

Figure 2. Upregulated proteins in skeletal muscle after late‐life PoWeR

A, molecular function analysis of upregulated proteins after PoWeR relative to SED. B, ASB6 (UniProt ID: Q91ZU1, adj. P = 0.0000005), RPL38 (Q9JJI8, adj. P = 0.0000009), RP2 (Q9EPK2, adj. P = 0.0000009), NDUFB5 (Q9CQH3, adj. P = 0.000004) and WASF2 (Q8BH43, adj. P = 0.000004) upregulation fold difference from SED after PoWeR. C, Xpo1 (cg09455096, P = 0.003 and cg08330824, P = 0.009) and Lpp (cg14306812, P = 0.01) promoter CpG methylation and corresponding protein upregulation (XPO1 adj. P = 0.016; LPP adj. P = 0.023) fold difference from SED after PoWeR. Individual fold change data points presented for PoWeR (n = 7) relative to the average of SED (n = 7). D, the relationship of CpG percentage methylation (any genomic context) after PoWeR to corresponding protein fold‐change for significant proteins (r = 0.123, r ² = 0.015, P = 0.06). E, the relationship of CpG percentage difference from sedentary after PoWeR (any genomic context) to its corresponding protein fold‐change for significant proteins (r = −0.017, r ² = 0.0003, P = 0.80).

Figure 3. BETA integration of the DNA methylome and proteome predicts epigenetic regulation of proteins after late‐life PoWeR

A, BETA analysis of the relationship between CpG methylation status and upregulated (P = 0.03) and downregulated (P = 1.0) proteins after PoWeR compared to background. Key ribosomal proteins (RPL6, adj. P = 0.002; RPL18, adj. P = 0.006; RPL29, adj. P = 0.0006; RPL31, adj. P = 0.005; RPLP0, adj. P = 0.001; RPS14, adj. P = 0.01), translational initiation factors (EIF2B3, adj. P = 0.004; EIF2B4, adj. P = 0.010; EIF2B5, adj. P = 0.010; EIF3D, adj. P = 0.006), mitoribosomal proteins (MRPL1, adj. P = 0.0006; MRPL2, adj. P = 0.001; MRPL4, adj. P = 0.002; MRPL10, adj. P = 0.0001; MRPL11, adj. P = 0.003; MRPL14, adj. P = 0.003, MRPL28, adj. P = 0.0004; MRPL37, adj. P = 0.00004; MRPL43, adj. P = 0.003; MRPL45, adj. P = 0.0002; MRPL46, adj. P = 0.001; MRPL, 53 adj. P = 0.010; MRPS11, adj. P = 0.0002; MRPS22, adj. P = 0.001; MRPS25, adj. P = 0.0001), mitochondrial complex I proteins (NDUFA2, adj. P = 0.0009; NDUFA3, adj. P = 0.002; NDUFA8, adj. P = 0.002; NDUFB3, adj. P = 0.0003; NDUFB9, adj. P = 0.002;NDUFB11, adj. P = 0.001; NDUFC2, adj. P = 0.0008; NDUFS1, adj. P = 0.002; NDUFS2, adj. P = 0.002; NDUFS6, adj. P = 0.010; NDUFS7, adj. P = 0.0009; TIMMDC1, adj. P = 0.004; TMEM186, adj. P = 0.0007), and proteins critical toskeletal muscle health (LACTB, adj. P = 0.000009; MIB1, adj. P = 0.0009; UBR4, adj. P = 0.001). B–F, histograms plotting the fold‐change difference with PoWeR relative to SED. Individual fold change data points are presented for PoWeR (n = 7) relative to the mean of SED (n = 8). G, the relationship of the BETA score and the methylation distance from the transcriptional start site (in base pairs) for upregulated proteins. H, Lisa cistrome results for the top 10 transcription factors (TFs) using the upregulated proteins that BETA predicted to be significantly influenced by their methylation status as determined by −log10 P‐value. I, Lisa results for the top 10 TFs using all upregulated proteins as determined by −log10 P‐value.

Figure 4. Relationship of the skeletal muscle proteome to skeletal muscle mass and in vivo function

A, the relationship of VPS13A protein abundance to absolute gastrocnemius mass (Pearson's coefficient r = 0.84, r ² = 0.71, P = 0.00008, Kendall's τ = 0.66, P = 0.0006) and normalized mass (Pearson's coefficient r = 0.82, r ² = 0.68, P = 0.0002, Kendall's τ = 0.58, P = 0.0025) (SED, n = 8; PoWeR, n = 7). B, the relationship of MYL12B (Pearson's coefficient r = 0.87, r ² = 0.75, P = 0.0001, Kendall's τ = 0.77, P = 0.0003) and NPL (Pearson coefficient r = 0.77, r ² = 0.60, P = 0.002, Kendall's τ = 0.72, P = 0.0006) protein abundance to in vivo plantar flexor fatigue resistance at 60 Hz stimulation (SED, n = 6; PoWeR, n = 7). C, the relationship of MYL12B (Pearson's coefficient r = 0.91, r ² = 0.83, P = 0.0003, Kendall's τ = 0.64, P = 0.0095) and NPL (Pearson's coefficient r = 0.80, r ² = 0.64, P = 0.006, Kendall's τ = 0.56, P = 0.025) (SED, n = 3; PoWeR, n = 7) protein abundance to in vivo total work. Individual data points presented for SED (black) and PoWeR (red). Data for fatigue resistance and total work calculations were not recorded for three SED mice and one PoWeR mouse since these mice exhibited cumulative workloads that suggested an inability to perform the task.