Mutation-independent Proteomic Signatures of Pathological Progression in Murine Models of Duchenne Muscular Dystrophy.

The absence of the dystrophin protein in Duchenne muscular dystrophy (DMD) results in myofiber fragility and a plethora of downstream secondary pathologies. Although a variety of experimental therapies are in development, achieving effective treatments for DMD remains exceptionally challenging, not least because the pathological consequences of dystrophin loss are incompletely understood. Here we have performed proteome profiling in tibialis anterior muscles from two murine DMD models (mdx and mdx52) at three ages (8, 16, and 80 weeks of age), all n = 3. High-resolution isoelectric focusing liquid chromatography-tandem MS (HiRIEF-LC-MS/MS) was used to quantify the expression of 4974 proteins across all 27 samples. The two dystrophic models were found to be highly similar, whereas multiple proteins were differentially expressed relative to WT (C57BL/6) controls at each age. Furthermore, 1795 proteins were differentially expressed when samples were pooled across ages and dystrophic strains. These included numerous proteins associated with the extracellular matrix and muscle function that have not been reported previously. Pathway analysis revealed multiple perturbed pathways and predicted upstream regulators, which together are indicative of cross-talk between inflammatory, metabolic, and muscle growth pathways (e.g. TNF, INFγ, NF-κB, SIRT1, AMPK, PGC-1α, PPARs, ILK, and AKT/PI3K). Upregulation of CAV3, MVP and PAK1 protein expression was validated in dystrophic muscle by Western blot. Furthermore, MVP was upregulated during, but not required for, the differentiation of C2C12 myoblasts suggesting that this protein may affect muscle regeneration. This study provides novel insights into mutation-independent proteomic signatures characteristic of the dystrophic phenotype and its progression with aging.

Duchenne muscular dystrophy (DMD) is a severe, X-linked, pediatric neuromuscular disorder characterized by progressive muscle wasting, loss of ambulation around age 10, and cardiorespiratory failure that is ultimately fatal (1–4). The disease is caused by mutations in the DMD gene that disrupt the translation reading frame, leading to loss of dystrophin protein expression (5, 6). Dystrophin is important for mechanical force transduction and is also involved in signaling functions (7–11), in part because of its role as an organizing center for the dystrophin-associated protein complex (DAPC). Absence of dystrophin leads to myofiber fragility, sensitivity to contractile damage, chronic cycles of myonecrosis and regeneration, persistent inflammation, and progressive fibro/fatty muscle degeneration (12–14).

The DMD gene consists of 79 exons and contains at least seven internal promoters (15) that give rise to the various dystrophin isoforms (e.g. Dp427, Dp260, Dp140, Dp116, Dp71, and Dp40) (supplemental Fig. S1). Some isoforms are ubiquitously expressed (e.g. Dp71) (16), whereas others are expressed in a more tissue-restricted pattern (such as Dp260, which is the retinal isoform of dystrophin) (17). As a result, the genomic locations of DMD-causing mutations may differentially affect the expression of these various isoforms, and by extension disease manifestation. For example, loss of the Dp71 isoform has been associated with cognitive impairment (18, 19).

Several dystrophic mouse strains have been developed to investigate DMD pathophysiology, and test novel therapeutics in vivo (20). The most used model is the mdx mouse, which carries a nonsense mutation in exon 23 leading to loss of the major muscle dystrophin isoform Dp427 (supplemental Fig. S1) (21, 22). Although the mdx mouse recapitulates some aspects of DMD pathology (23), it is generally considered to exhibit mild muscular dystrophy with only a small reduction in life-span (24). mdx mice undergo a brief period of degeneration and regeneration between 2 and 12 weeks of age (25–27), with more pronounced muscle pathology and cardiomyopathy manifesting much later in life (28, 29). Importantly, the exon 23 mutation observed in the mdx mouse does not typically occur in boys with DMD (30), which has motivated the development of more patient-relevant dystrophic mouse models. To this end, Araki et al. generated the mdx52 mouse model in which Dmd exon 52 is deleted, leading to the absence of the Dp260 and Dp140 isoforms in addition to Dp427 (supplemental Fig. S1) (31). Deletions in the so-called “hot-spot” region (DMD exons 45-55) are some of the most observed mutations in boys with DMD, thereby making this model more patient-relevant than the more widely used mdx mouse (6). We recently reported differences in the number of dystrophin-positive revertant fibers and regenerating fibers between mdx and mdx52 mice (32). Specifically, mdx mice contained higher numbers of revertant fibers at all ages tested (2–18 months of age), whereas mdx52 mice contained elevated numbers of centrally-nucleated fibers at 2 months of age only (32). Importantly, the mdx52 model also allows for the testing of patient mutation-relevant exon skipping strategies in vivo (e.g. targeting exon 51 or exon 53) (33–35).

Despite significant research effort, many aspects of dystrophic pathology remain unclear. As such, there is a need to understand the complex molecular mechanisms underlying DMD at the level of gene and protein expression. Transcriptomics methodologies have enabled the simultaneous measurement of tens of thousands of genes in the muscles of both dystrophin-deficient animal models and DMD patient biopsy samples (36–40). We have previously compared global transcript and protein expression in mdx versus control muscle (39). Importantly, mRNA expression often does not correlate with protein levels (39, 41–43), meaning that investigations that rely on transcriptomics alone may potentially be misleading. Although the global quantification of nucleic acid transcripts is relatively simple using digital gene expression analysis (i.e. RNA-sequencing) or hybridization methods (e.g. DNA microarrays), proteomics analysis is substantially more challenging. Early proteomics studies in dystrophic muscle used 2D-electrophoresis to identify differentially abundant proteins, but results were limited to only a handful of differential expression calls (44, 45). The use of other methodologies such as iCAT, in vivo SILAC, and label-free approaches resulted in the identification of further differentially expressed proteins in dystrophic muscle, although these studies were still limited by their low proteomic coverage (∼1000 proteins quantified, or less) (46–48). Global proteome profiling in fibrous tissues such as skeletal muscle is complicated by the presence of very high concentrations of a few structural proteins, such as actins and myosins (49). Peptides derived from these proteins mask signals from lowly abundant proteins and thereby limit the depth of proteome coverage that can be achieved. To increase analytical depth, we recently applied high-resolution sample pre-fractionation based on narrow-range isoelectric focusing of peptides to quantify expression of over 3272 proteins in mouse muscle (39).

To date, there have been relatively few high-resolution proteomics studies in dystrophic muscle, and fewer still that have measured global changes in protein expression in muscle throughout the progression of pathology over time (39, 47, 50). Here, we have performed MS-based proteomic profiling in both the mdx and mdx52 DMD mouse models compared with WT controls (WT) at three ages representing different stages of dystrophic pathology. We show that this state-of-the-art proteomic strategy has uncovered previously unidentified pathological pathways in dystrophic mouse models, which are potential therapeutic targets for DMD.

EXPERIMENTAL PROCEDURES

Animals

Mice were housed under 12:12 h light-dark conditions with food and water ad libitum. All experimental protocols in this study were approved by the Experimental Animal Care and Use Committee of the National Institute of Neuroscience, NCNP, Japan.

mdx52 mice were generated at our facility at the NCNP (31) and have been back-crossed with C57BL/6 mice for more than 10 generations. mdx mice on a C57BL/6 background were kindly provided by Dr T. Sasaoka (Brain Research Institute, Niigata University, Niigata, Japan). C57BL/6 mice were used as controls to match the background of the dystrophic strains (i.e. mdx52 and mdx). Serum and tissues from each strain were collected at 4, 8, 16, 24, 48, and 80 weeks of age (n = 3–5 per group). Tibialis anterior (TA) muscle of all strains at 8, 16, and 80 weeks (n = 3) was subsequently cryosectioned, collecting 50 sections of 10 μm for each sample.

Cell Culture

C2C12 myoblasts were cultured at 37 °C with 5% CO2 in Dulbecco's modified Eagle's medium (DMEM) containing 20% fetal bovine serum (FBS) and 1% antibiotics/antimycotics (growth medium: GM) (all Invitrogen, Carlsbad, CA). DMEM supplemented with 2% horse serum (HS) and 1% antibiotics/antimycotics (differentiation medium: DM) was used to differentiate C2C12 myoblasts for 3–6 days to form multinucleated myotubes. C2C12 myoblasts were seeded in 24-well and 6-well plates at 100,000 cells and 400,000 cells per well respectively. For transfections, cells were incubated in GM before addition of 50 nm siRNA complexes (either targeting Mvp or a control siRNA, ON-TARGETplus siRNA: Dharmacon, Cambridge, UK). Complex formation was performed in the absence of serum, and cells collected after 3 days in DM.

Serum RNA Extraction and miRNA RT-qPCR

Serum (50 µl) from each sample (n = 2–4) was mixed with TRIzol LS (ThermoFisher Scientific), supplemented with 3 µl of a spike-in of 5 nm cel-miR-39 (5ʹ-UCACCGGGUGUAAAUCAGCUUG-3ʹ) per sample as an exogenous reference control. RNA was then extracted according to manufacturer's instructions, with modifications previously described in (51). Samples were stored at −80 °C. TaqMan MicroRNA Reverse Transcription Kit (Applied Biosystems, ThermoFisher Scientific, Waltham, MA) was used for cDNA synthesis of miRNAs. Primers for miR-1a-3p (5ʹ-UGGAAUGUAAAGAAGUAUGUAU-3ʹ), miR-133a-3p (5ʹ-UUUGGUCCCCUUCAACCAGCUG-3ʹ), miR-206-3p (5ʹ-UGGAAUGUAAGGAAGUGUGUGG-3ʹ), miR-223-3p (5ʹ-UGUCAGUUUGUCAAAUACCCCA-3ʹ) and cel-miR-39 were used, following manufacturer's instructions. TaqMan Gene Expression Master Mix (Applied Biosystems) was used for qPCR, following manufacturer's instructions. Primers for miR-1a-3p, miR-133a-3p, miR-206-3p, miR-223-3p and cel-miR-39 were obtained from ThermoFisher and data were normalized as previously described (52).

HiRIEF-nanoLC-MS/MS-Based Proteomics

Sample Preparation for Mass Spectrometry

Sections of tibialis anterior from C57BL/6, mdx and mdx52, in biological triplicates, at 8, 16, and 80 weeks were lysed with 4% SDS, 25 mm HEPES, 1 mm DTT. Samples were prepared using a modified version of the spin filter aided sample preparation protocol (53). Lysates were heated to 95 °C for 5 min followed by sonication for 1 min and centrifugation at 14,000 × g for 15 min. The supernatant was mixed with 1 mm DTT, 8 M urea, 25 mm HEPES, pH 7.6 and transferred to a 10 kDa cutoff Nanosep centrifugation filtering unit (Pall, Port Washington, NY) for centrifugation at 14,000 × g for 15 min. Proteins were alkylated by treatment with 50 mm iodoacetamide in 8 M urea, 25 mm HEPES for 10 min. The proteins were then centrifuged at 14,000 × g for 15 min followed by 2 more additions and centrifugations with 8 M urea, 25 mm HEPES. Trypsin (Promega, Southamption, UK) in 250 mm urea, 50 mm HEPES was added to the cell lysate at a ratio of 1:50 trypsin/protein and incubated overnight at 37 °C with gentle shaking. The filter units were centrifuged at 14,000 × g for 15 min followed by another centrifugation with Milli-Q water and the flow-through was collected. Peptides were labeled using the TMT10plex Isobaric Label Reagent Set according to the manufacturer's protocol (Thermo Fisher Scientific) and cleaned-up using strata-X-C-cartridges (Phenomenex, Torrance, CA).

IPG-IEF of Peptides

TMT labeled peptides were separated by immobilized pH gradient - isoelectric focusing (IPG-IEF) on pH 3.7–4.9 strips (250 µg peptides per strip) as described previously (39, 53). Peptides were extracted from the strips using a prototype liquid handling robot, supplied by GE Healthcare Bio-Sciences AB. A plastic device with 72 wells was put onto each strip and 50 µl of Milli-Q water was added to each well. After incubation for 30 min, the liquid was transferred to a 96-well plate and the extraction was repeated 2 more times. The extracted peptides were dried in speed vac for storage and dissolved in 3% acetonitrile, 0.1% formic acid before MS analysis.

LC–MS/MS Analysis

Online LC–MS was performed using a hybrid Q Exactive - HF mass spectrometer (Thermo Fisher Scientific). Before analysis on the Q Exactive (Thermo Fisher Scientific), peptides were separated using an Ultimate 3000 RSLCnano system. Samples were trapped on an Acclaim PepMap nanotrap column (C18, 3 μm, 100 Å, 75 μm x 20 mm), and separated on an Acclaim PepMap RSLC column (C18, 2 μm, 100 Å, 75 μm x 50 cm) (Thermo Fisher Scientific). Peptides were separated using a gradient of A (5% DMSO, 0.1% formic acid) and B (90% acetonitrile, 5% DMSO, 0.1% formic acid), ranging from 6% to 37% B in 50 min with a flow of 0.25 µl/min. The Q Exactive was operated in a data-dependent manner, selecting top 5 precursors for fragmentation by higher-energy collisional dissociation (HCD). The survey scan was performed at 70,000 resolution from 400–1600 m/z, with a max injection time of 100 ms and a target of 1 × 106 ions. For generation of HCD fragmentation spectra, a max ion injection time of 150 ms and automatic gain control of 1 × 105 were used before fragmentation at 30% normalized collision energy, 35,000 resolution. Precursors were isolated with a width of 2 m/z and put on the exclusion list for 70 s. Single and unassigned charge states were rejected from precursor selection.

Peptide and Protein Identification

All Orbitrap data were searched using SequestHT under the software platform Proteome Discoverer 1.4 (Thermo Fisher Scientific) against the Ensembl 78 mouse protein database (53,838 protein entries) and filtered to a 1% false discovery rate (FDR). A precursor mass tolerance of 10 ppm, and product mass tolerances of 0.02 Da for HCD-FTMS were used. Further settings used were: trypsin with 2 missed cleavage; iodoacetamide on cysteine and TMT on lysine and N-terminal as fixed modifications; and oxidation of methionine as variable modification. Quantification of TMT-10plex reporter ions was performed using Proteome Discoverer on HCD-FTMS tandem mass spectra using an integration window tolerance of 10 ppm. Only peptides unique to a protein group were used for quantitation. A pool of all samples was used in one TMT tag as a linker (denominator) between TMT sets. For data analysis, proteins that could not be identified in all samples were filtered out from the final list and values were log2 transformed.

Western Blot

Tissue samples were lysed in 75 mm Tris-HCL (pH 6.8), supplemented with 10% SDS, 5% 2-mercapto-ethanol and 3% protease inhibitors (Sigma-Aldrich, Dorset, UK). These were then heated at 100 °C for 3 min and centrifuged at 13,000 × g for 10 min. Cells were lysed in RIPA lysis buffer (ThermoFisher Scientific) supplemented with protease and phosphatase inhibitors (Sigma-Aldrich), incubated at 4 °C for 5 min, then centrifuged at 10,000 × g for 5 min. Protein was stored at −80 °C. Protein lysates were resolved on polyacrylamide gels and transferred onto 0.45 μm polyvinylidene fluoride (PVDF) membranes. Membranes were blocked in TBST (50 nm Tris pH 7.5, 150 mm NaCl and 0.1% Tween 20), supplemented with 5% nonfat dry milk for 1 h at room temperature. These were subsequently probed with the primary antibody in TBST/blocking solution overnight at 4 °C. Membranes were then incubated with the secondary antibody for one hour at room temperature. Membranes were imaged using the ECL Prime Western blotting system (Sigma-Aldrich) and imaged using an Odyssey Fc Imager (LI-COR Biotechnology, Lincoln, NE). Antibodies and their dilutions are detailed in supplemental Table S1. Equal protein loading was determined by Coomassie Brilliant Blue (CBB) staining for tissue samples and GAPDH immunoblot for cell culture samples was used to determine protein loading.

RT-qPCR

RNA was extracted using TRIzol reagent according to manufacturer's instructions. Briefly, samples were homogenized in TRIzol (Thermo Fisher Scientific) after which chloroform was added. Samples were centrifuged resulting in phase separation whereby the supernatant was removed and RNA precipitated with 500 µl of isopropanol, followed by a 10 min incubation at room temperature and centrifugation. Pellets were subsequently washed with 75% ethanol, air-dried and resuspended in nuclease-free water. RNA was stored at −80 °C.

High-capacity cDNA Reverse Transcription Kit (Thermo Fisher Scientific) was used for cDNA synthesis according to manufacturer's instructions. Power SYBR Green master mix was used for gene expression analysis according to manufacturer's instructions, using 1:10 diluted cDNA (primer sequences can be found in supplemental Table S2). Values were normalized to the geometric average of two stable reference genes, Rplp0 and Tbp for tissue samples and Rplp0 for cells. RT-qPCR data were analyzed using the Pfaffl method and PCR efficiencies determined by LinRegPCR (Amsterdam Medical Center, Amsterdam, the Netherlands).

Immunofluorescence

Cells were fixed in 4% paraformaldehyde for 10 min and permeabilized for 15 min in 0.2% Triton X-100. Cells were then blocked for 30 min in 5% Bovine Serum Albumin and incubated with anti-MHC antibody at a 1:20 dilution (MF20, obtained from the Developmental Studies Hybridoma Bank (DSHB), created by the NICHD of the NIH and maintained at The University of Iowa, Department of Biology, Iowa City, IA) overnight at 4 °C. Cells were subsequently probed with secondary antibody (1:500 dilution, goat anti-mouse AlexaFluor 488) and Hoechst33342 (1:5,000 dilution, Thermo Fisher Scientific) for 1 h at room temperature, then stored at 4 °C until imaged. Myogenic Index and Fusion Index were defined as the percentage of nuclei within MHC-positive cells, and the percentage of nuclei within MHC-positive myotubes containing at least 3 nuclei, respectively.

Statistics and Bioinformatics

Statistical analysis was performed using R i386 3.2.0 (The R Project) and R Studio (Boston, MA). TMT ratios were log2 transformed and significance was tested by t test for two-sample comparisons and one-way ANOVA for comparisons between more than two groups using the rowttests (genefilter library) and row.oneway.anova (HybridMTest library) R functions respectively. P-values were corrected for multiple comparisons using the Benjamini-Hochberg method with the p.adjust function from the base R. A stringent adjusted P-value threshold of 0.01 was used to limit the number of false discoveries. Further statistical analysis on individual proteins was performed with GraphPad Prism 7 (GraphPad Software Inc, La Jolla, CA), whereby one-way ANOVA comparisons were made at each age with Tukey post hoc correction. PCA analysis was performed using the prcomp R function. GraphPad Prism 7 was used for visualizing PCA, volcano, and other plots. TMeV (Institute for Genomic Research, Rockville, MD) (54) was used to generate heatmaps and identify k-means protein clusters. The optimum number of k-clusters was determined using Figures of Merit in TMeV and the elbow method using the wssplot R function. k-means plots were generated in R.

Gene list enrichment analysis was performed using ToppFun (55). Canonical pathway analysis and upstream regulator analysis was performed using Ingenuity Pathway Analysis (IPA, Qiagen Bioinformatics, Redwood City, CA). P < 0.01 and Z≥× 2 × were considered significant. Venn diagrams were generated using Venny 2.1.0 (BioinfoGP Service, Centro Nacional de Biotechnología (CNB-CSIC), Madrid, Spain).

RESULTS

Experimental Design

To investigate differential protein expression in dystrophic muscle we collected tibialis anterior (TA) muscles from two dystrophic mouse models (mdx and mdx52). Muscles were collected from age- and sex-matched C57BL/6 WT controls (all mice in this study were male). The choice of ages was based on our previous observations of differences in muscle regeneration and the number of revertant fibers between the mdx and mdx52 strains (32), and on the analysis of serum miRNA levels (supplemental Fig. S2). Our group and others have shown that the muscle-enriched microRNAs (the myomiRs: miR-1a-3p, miR-133a-3p, and miR-206-3p) are highly elevated in serum of dystrophic animal models and DMD patients, indicative of increased muscle turnover (56–62). We therefore analyzed serum myomiR expression at different ages throughout the course of dystrophic pathology to characterize myomiR abundance patterns in the two dystrophic mouse models with age. Serum from animals at 4, 8, 16, 24, 48, and 80 week-old mice was collected, RNA extracted, and myomiRs measured by small RNA TaqMan RT-qPCR. miR-223-3p was included as a nonmyomiR endogenous control that was expected to be relatively stable based on previous studies (52, 56, 63). MyomiR expression was elevated in mdx and mdx52 relative to WT at all ages (supplemental Fig. S2). Interestingly, myomiR levels in mdx52 serum peaked earlier than for the mdx mice (i.e. 8 weeks versus 16 weeks), consistent with the earlier peak of muscle regeneration observed in this model (32). At the later ages, myomiR abundance remained significantly increased relative to controls. Based on these findings, we selected the 8, 16, and 80-week ages in which to perform proteomic profiling across the mouse strains. The 80-week samples were included as a representative of “aged” muscle in which dystrophic pathology is more advanced (12, 58).

Using these samples, we aimed to identify differences between the mdx and mdx52 models that may account for their distinct phenotypes, establish a mutation-independent proteomic signature of the dystrophic condition, and to investigate alterations in the dystrophic proteome associated with pathological progression.

Proteomics Analysis in Dystrophic Muscle

TA protein lysates from male C57BL/6, mdx, and mdx52 mice at 8, 16, and 80 weeks of age (n = 3) were analyzed by HiRIEF-LC–MS/MS (53) and samples labeled using TMT chemistry, allowing for multiplexing of up to 10 samples during a single LC–MS/MS run. Three separate runs were performed whereby a common bridge sample (generated by mixing protein lysates from each sample in equal amounts) was included in each set of 10 samples. TMT ratios were calculated relative to the bridge sample, thereby allowing for direct comparison of samples analyzed on different runs, and effectively increasing the level of multiplexing to 27-plex.

A total of 7111 proteins were identified, of which 4974 (70%) were quantified in all samples and used for proteomics analysis (Fig. 1). Details of peptides analyzed for protein identification and quantification are shown in supplemental Fig. S3. Notably, methodological performance was highly similar when comparing between the three separate LC–MS/MS runs (supplemental Fig. S3). This performance is a substantial improvement on our previous data using similar methodology (iTRAQ labeling with HiRIEF-LC–MS/MS, where 3057 proteins were quantifiable (39)). Normalized TMT ratios for the full data set are provided in supplemental File S1.

Fig. 1.

Protein proportions and dynamic range in murine skeletal muscle tissue. Peptide-spectrum-match (PSM) data were pooled for all experimental groups for proteins detected in all samples. Proteins were ranked by the number of PSMs/protein and the percentage of total protein signal estimated for each protein by dividing the PSMs/protein by the total number of PSMs, multiplied by 100%. The resulting data are shown as a cumulative frequency plot. Dystrophin-associated protein complex proteins are highlighted as black boxes and the top 25 proteins accounting for 50% of detected protein mass are indicated.

Analysis of the protein cumulative distribution frequency showed that only 25 proteins accounted for 50% of the detected protein mass, which included myosin heavy chain proteins (MYH1, MYH2, MYH3, MYH4, MYH6, MYH7, MYH8, and MYH13), actins (ACTA1, ACTA2), actinin (ACTN3), titin (TTN), nebulin (NEB), tropomyosin (TPM1), and creatine kinase (CKM) (Fig. 1). In contrast, the least abundant 1042 proteins comprised only 1% of the matched spectra (Fig. 1).

Proteins that differ between any of the nine experimental groups were identified by one-way ANOVA analysis and 3492 proteins (70.2%) were identified as significantly different (Benjamini-Hochberg adjusted P < 0.01). TMT ratios for differentially expressed proteins were visualized by heatmap with unsupervised hierarchical clustering (Fig. 2A) and principal component analysis (PCA, Fig. 2B). Biological replicates were tightly clustered, and a clear separation between the WT animals and the two dystrophic strains was apparent. mdx and mdx52 samples were intermingled, suggestive of a common pattern of proteomic perturbation in these dystrophic models. The greatest source of variation in the data, represented by principal component 1 (PC1, containing 37.5% of the variation) on the PCA biplot (Fig. 2B), primarily reflected the difference between WT and dystrophic animals. Conversely, principal component 2 (PC2, containing 24.4% of the variation) reflected the progression of age, with the low, medium, and high values of PC2 corresponding to 8-, 16-, and 80-week-old mice, respectively (Fig. 2B). Highly similar clustering patterns were observed when these analyses were performed on all proteins (supplemental Fig. S4).

Fig. 2.

Differential protein expression in dystrophic muscle. Quantitative proteomics analysis of TA muscles from WT (C57BL/6), mdx and mdx52 mice at 8, 16 and 80 weeks of age (all n = 3). Significantly different proteins identified by one-way ANOVA (Benjamini-Hochberg adjusted P < 0.01) are visualized by (A) hierarchical clustering and heatmap, and (B) principal component analysis (PCA). Equivalent analyses for all proteins are shown in supplemental Fig. S4.

To assess the performance of the proteomics analysis we examined the expression of proteins that are known to be differentially expressed in dystrophic muscle. Dystrophin (DMD) was observed to be highly downregulated in the dystrophic animals at all ages, consistent with its genetic disruption. Accordingly, expression of the DAPC components: NOS1, SGCA, SGCB, SGCG, DAG1, SNTA1, and DTNA were all significantly downregulated in both mdx and mdx52 mice at all ages (supplemental Fig. S5). Loss of dystrophin leads to the mislocalization and subsequent destabilization of DAPC components (64). As a result, downregulation of DAPC proteins is an established feature of dystrophic muscle (39, 65). In contrast, UTRN (utrophin) was significantly upregulated in dystrophic muscle (supplemental Fig. S5). Utrophin is a paralog of dystrophin that is primarily expressed at the neuromuscular and myotendinous junctions and is known to be upregulated at the dystrophin-deficient sarcolemma (39, 66–69). Similarly, aquaporin-4 (AQP4), periostin (POSTN), myostatin (MSTN), and biglycan (BGN) were also differentially expressed in dystrophic muscle as reported previously (39, 70–73). In summary, consistent differential expression of these proteins in both dystrophic strains across all ages (supplemental Fig. S5) demonstrates the robustness of both our proteomics analysis and multiplexing strategy.

A Mutation-Independent Dystrophic Proteome Signature

No significant changes in protein expression were detected when comparing the mdx and mdx52 dystrophic strains by t test (supplemental Fig. S6), consistent with the close clustering of these samples in the heatmap and PCA analyses described above (Fig. 2). Absolute fold changes for the mdx52 versus mdx comparisons were lower than for other comparisons (supplemental Fig. S7), and the coefficients of variation for the mdx52 samples are similar to, or lower than, other experimental groups (indicating that these samples are not inherently noisy, supplemental Fig. S7). Together these data suggest that the lack of significant differential expression calls in the mdx52 versus mdx comparisons is a consequence of the high similarity of the muscle proteomes in these strains, rather than a deficiency in statistical power.

Separate t-tests were performed to compare mdx versus WT (supplemental Fig. S6) or mdx52 versus WT at each age (supplemental Fig. S6). Relatively few differentially expressed proteins (<35) were observed at the 8 week and 16 week ages. In contrast, several hundred differentially expressed proteins were identified in the 80-week muscle samples, of which 199 were common between the mdx and mdx52 strains (supplemental Fig. S6).

Considering the similarity in proteomes between the two dystrophic strains, and the relatively small number of differences detected when comparing each dystrophic model against the WT controls, we decided to pool the TMT ratios for the mdx and mdx52 samples at each age in order to increase group size (to n = 6) and thereby boost statistical power. Differentially expressed proteins in the pooled dystrophic muscles were assessed by t test at each age and 805, 837, and 1618 differentially expressed (adjusted P < 0.01) proteins identified at 8, 16, and 80 weeks respectively (Fig. 3A–3C). Proteins that were uniquely upregulated at 8 weeks and 16 weeks were enriched for RNA processing and ribosome/translation-associated gene ontology (GO) terms respectively (Fig. 3D, 3E). 909 proteins were differentially expressed at 80 weeks only and were highly enriched for GO terms associated with mitochondria and metabolic processes (Fig. 3F). 228 proteins were commonly differentially expressed at all three ages (Fig. 3G), which were enriched for GO terms associated with muscle function, sarcolemma, extracellular matrix, and the DAPC (Fig. 3H).

Fig. 3.

Differentially expressed proteins in pooled dystrophic muscle samples at 8, 16 and 80 weeks of age. TMT ratios for the mdx and mdx52 samples (Dystrophic) were pooled and differential expression relative to WT controls assessed at each age (Student's t tests, Benjamini-Hochberg adjusted P < 0.01). Volcano plots for (A) 8 weeks, (B) 16 weeks, and (C) 80 weeks of age. Gene list enrichment analysis was performed for proteins that were uniquely differentially expressed at (D) 8 weeks only, and (E) 16 weeks only. G, Venn diagram of differentially expressed proteins (dystrophic versus WT) showing the overlap between ages. H, Gene list enrichment analysis in proteins that were commonly differentially expressed at all ages.

TMT ratios for all dystrophic (n = 18) and all WT (n = 9) animals were further pooled, and statistically significant differences between these groups tested by t test. The purpose of this analysis was to identify proteins that are differentially expressed in dystrophic muscle independent of mutation type or age (with the added benefit of increased statistical power). Of the 4974 proteins detected, 1795 (36%) proteins were significantly different between pooled dystrophic and WT samples (adjusted P < 0.01) (Fig. 4A) of which 891 were upregulated and 904 downregulated.

Fig. 4.

Differentially expressed proteins and pathway analysis in dystrophic muscle independent of age or mutation type. All TMT ratios for the mdx and mdx52 samples (Dystrophic, n = 18) were pooled and differentially expressed proteins determined relative to WT controls (n = 9) independent of age (Student's t tests, adjusted P < 0.01). Differentially expressed proteins were (A) visualized by volcano plot and analyzed with IPA in order to identify (B) the top 20 most significant (P < 0.01) enriched canonical pathways, (C) the top 10 up- and downregulated pathways based on activation state (i.e. Z-score), and (D) putative upstream regulators (based on the expression ratios of their known downstream targets). Canonical pathways and upstream regulators are ranked by Z-score (positive values indicate activation and negative values indicate inhibition of the predicted regulator). Solid purple lines indicate P-value. Broken purple lines indicate the P = 0.01 significance threshold.

The top 10 consistently up- and downregulated proteins for all dystrophic versus WT comparisons are shown in supplemental Figs. S8 and S9 (excluding DAPC component proteins that were highly differentially expressed, as in previous analyses). The most highly upregulated proteins in dystrophic muscle included: AHNAK nucleoprotein 2 (AHNAK2), HtrA serine peptidase 3 (HTRA3), cathepsin S (CTSS), tubulin beta 6 class V (TUBB6), and an expressed pseudogene of carboxylesterase 2D (CES2D-PS). The enzymes SAM decarboxylase 2 (AMD2) and 3-hydroxybutyrate dehydrogenase 1 (BDH1) were among the most downregulated proteins in dystrophic muscle. The proteoglycan keratocan (KERA) was downregulated at the 8 and 16 week ages only. EGF like, fibronectin type III and laminin G domains (EGFLAM, also known as pikachurin) was consistently downregulated in dystrophic muscle and has been previously shown to interact with α-dystroglycan (DAG1) (74) suggesting that downregulation of EGFLAM may be a consequence of DAPC disruption.

Protein expression ratios from the pooled analysis were analyzed using ingenuity pathway analysis (IPA) to identify perturbed canonical pathways (Fig. 4B, 4C) and predicted upstream regulators (Fig. 4D) in dystrophic muscle. This analysis identified multiple findings consistent with established features of dystrophic pathology (e.g. calcium signaling, inflammasome pathway, nitric oxide signaling, and TNF signaling). Furthermore, many affected canonical pathways and predicted upstream regulators were identified with no known association with DMD pathology (e.g. NRF2-mediated oxidative stress response and ATM signaling, KDM5A, and RICTOR, Fig. 4B, 4C). Mitochondrial Dysfunction was the most significantly affected canonical pathway, consistent with our previous study (39). Notably, the increased depth of proteome coverage enabled the detection of protein expression changes that were not possible in our previously published analysis (39). For example, multiple protein components of Complex IV were downregulated in dystrophic muscle in addition to the other four complexes of the electron transport chain that we reported previously (39). Metabolic pathways associated with the mitochondria were similarly downregulated (i.e. TCA Cycle II, Oxidative Phosphorylation, Fatty Acid β-oxidation I).

To identify the most important proteins contributing to dystrophic pathology, differentially expressed proteins were filtered to include only those with a fold change greater than 2 in at least one comparison, and DAPC members were excluded. The resulting list of 200 proteins was classified into muscle, extracellular matrix, immune response, and lipid metabolism based on gene list enrichment analysis (Fig. 5). Multiple myosin heavy and light chains were differentially expressed in dystrophic muscle (i.e. Downregulation of MYH2, MYL2, MYL3, MYL10 and upregulation of MYH3, MYL4, MYL6B), consistent with our previous study (39). Other structural components of muscle were perturbed, including downregulation of myozenin-2 (MYOZ2), myoglobin (MB), the Na+/Ca2+ exchanger SLC8A1, and two isoforms of sarcolemma associated protein (SLMAP). Conversely, cardiac muscle troponin T (TNNT2), gelsolin (GSN), and p21 (RAC1) activated kinase 1 (PAK1) were upregulated in dystrophic muscle (Fig. 5A). In general, differential expression of these muscle-associated proteins was similar across ages, consistent with being consequences of dystrophin absence in general, rather than as downstream pathologies, which may develop with disease progression.

Fig. 5.

Categorization of the most differentially expressed proteins in dystrophic muscle. Proteins that were significantly different (adjusted P < 0.01) in at least one inter-sample comparison with a fold change greater than 2 (excluding the DAPC components: DMD, DAG1, DNTA, NOS1, SNTA1, SGCA, SGCB, SGCG, SGCD) were analyzed using gene list enrichment analysis. Protein expression was visualized by heatmap for: A, muscle system process (GO:0003012), B, extracellular matrix (GO:0031012), C, immune effector process (GO:0002252), and D, lipid metabolic process (GO:0006629).

Analysis of differentially expressed extracellular matrix-associated proteins revealed reciprocal regulation of different collagen proteins, with COL1A1, COL1A2, and COL2A1, COL22A1 all downregulated, and collagens COL3A1 and COL5A2 upregulated. (Elevation of collagen-3 has previously been linked to fibrotic pathology in dystrophin-deficient muscle (75)). Other upregulated extracellular matrix-associated proteins included myocilin (MYOC), sushi repeat-containing protein (SRPX), prolargin (PRELP), the serine protease HTRA1, the annexins: ANXA1 and ANXA4, and the S100 proteins: S100A4 and S100A11. Conversely, fibromodulin (FMOD), fibrillin-2 (FBN2), and tenascin C (TNC) were all downregulated in dystrophic muscle. Interestingly, two thrombospondin genes were differentially expressed, although in opposite directions, with THBS3 being upregulated and THBS4 downregulated in dystrophic muscle. Similarly, the serpins: SERPINB1 and SERPINF1, were up and downregulated, respectively. In contrast, with the muscle-associated proteins described above, some extracellular matrix factors were perturbed at the 80-week age only or exhibited a progressive increase in expression with age (i.e. COL3A1, MYOC, SRPX, BGN, THBS3, and HTRA1). Together these findings are consistent with widespread matrix remodeling and later-onset fibrotic damage (Fig. 5B).

The majority of immune response-associated proteins were upregulated in dystrophic muscle and include macrophage-related genes (integrin alpha M (ITGAM) and CD5L), complement factors (C7), and immunoglobulin genes (IGHM and IGKC) (Fig. 5C). Differentially expressed proteins associated with lipid metabolism included fatty acid binding protein 3 (FABP3), phospholipase A2 group VII (PLA2G7), and the long-chain fatty acid transporter SLC27A1, which were downregulated, and apolipoprotein D (APOD), phospholipid transfer protein (PLTP), and legumain (LGMN) that were upregulated in dystrophic muscle (Fig. 5D). As with the extracellular matrix-associated proteins, the immune-response and lipid metabolism categories contained some proteins that were differentially expressed at all ages, and some that were perturbed to a greater extent in the 80-week-old samples.

Proteomic Changes Associated with the Progression of Dystrophic Pathology

We next sought to identify proteins that were differentially expressed between dystrophic and WT samples, and which exhibited age-associated changes in expression (and therefore pathological progression). One-way analysis of variance was performed on the dystrophic samples in order to identify significant differences (adjusted P < 0.01) between 8-, 16-, and 80-week-old pooled dystrophic mice, and this list was then filtered using the list of proteins that were found to be differentially expressed in any of the analyses performed above. The overlapping 1196 proteins were then classified into 5 k-clusters based on their expression patterns. The character of each cluster was then illustrated by line graphs of the mean expression of all proteins (Fig. 6B) and heatmaps (Fig. 6C). Gene list enrichment analysis was performed to assign biological meaning to each cluster.

Fig. 6.

Changes in dystrophy-associated protein expression with pathological progression. A, Venn diagram showing the overlap between significantly different (adjusted P < 0.01) proteins when comparing all the Dystrophic versus WT comparisons and the 8-week versus 16-week versus 80-week comparisons (performed in the Dystrophic samples only). The resulting 1,196 proteins were analyzed using k-means clustering. B, Line plots showing the mean±S.E. expression values across all proteins contained in each k-cluster. C, Expression for individual proteins contained in each cluster is shown by heatmap. A subset of proteins in cluster 4 that are upregulated in 80 week dystrophic muscle only are highlighted with the † symbol. Cluster 5 was further analyzed, and protein expression was visualized by heatmap for proteins associated with the following gene ontology terms: D, extracellular matrix (GO:0031012), E, endoplasmic reticulum membrane (GO:0005789), and F, enzyme regulator activity (GO:0030234). All heatmaps were sorted using hierarchical clustering in the gene dimension only, for clarity.

Cluster 1 contained proteins that progressively declined with age but were generally elevated in dystrophic muscle throughout. This cluster was enriched for GO terms associated with focal adhesion and mRNA splicing. Conversely, Cluster 4 contained proteins that were progressively upregulated with age in both WT and dystrophic muscle. In general, expression of proteins in this cluster was downregulated in dystrophic muscle, although a subset of proteins (marked with † in Fig. 6C) were upregulated. This subset included proteins such as integrin linked kinase (ILK), EH domain containing 3 (EHD3), asporin (ASPN), indolethylamine N-methyltransferase (INMT), osteoglycin (OGN), major vault protein (MVP), and immunoglobulin like and fibronectin type III domain containing 1 (IGFN1). Notably, several proteins associated with other muscular dystrophies were also found in this sub-list including caveolin 3 (CAV3), lamin A/C (LMNA), muscleblind like splicing regulator 1 (MBNL1), and fukutin related protein (FKRP) (linked to distal myopathy, Emery-Dreifuss muscular dystrophy, myotonic dystrophy, and limb-girdle muscular dystrophy respectively). The proteins in Cluster 4 were associated with the GO terms sarcomere, sarcolemma, extracellular matrix, neuromuscular junction, and calcium channel complex.

Cluster 2 proteins were downregulated in dystrophic muscle, but the magnitude of downregulation was less apparent in the 16-week-old samples. Proteins in Cluster 3 were downregulated in dystrophic muscle at the 80 week age only. Gene list enrichment analysis showed Cluster 2 and 3 were primarily enriched for mitochondria-related GO terms (data not shown).

Cluster 5 was particularly interesting, as it contained proteins that were upregulated in 80-week-old dystrophic mice, and so are therefore likely associated with more advanced pathological progression. Proteins in Cluster 5 were enriched for the GO terms: extracellular matrix (Fig. 6D), endoplasmic reticulum membrane (Fig. 6E), and enzyme regulator activity (Fig. 6F). The effect sizes for protein expression changes were highest for the extracellular matrix associated changes, which included several proteins discussed above (i.e. BGN, COL3A1, HTRA1, SRPX, and THBS3). Notable other extracellular matrix-associated changes in dystrophic muscle included upregulation of annexin A6 (ANXA6), COL6A6, decorin (DCN), platelet derived growth factor subunit B (PDGFB), S100A6, transforming growth factor beta induced (TGFBI), and the laminins: LAMA4, LAMA5, and LAMB1 (Fig. 6D). The endoplasmic reticulum membrane-associated changes included upregulation of UBX domain protein 4 (UBX4), gamma-glutamyl carboxylase (GGCX), cytochrome P450 family 39 subfamily A member 1 (CYP39A1), ras homolog family member G (RHOG), and dystrophia myotonica protein kinase (DMPK) (Fig. 6E). The enzyme regulator activity-associated changes include upregulation of phosphoinositide-3-kinase regulatory subunit 1 (PIK3R1) and caspase-1 (CASP1).

In summary, the proteins in Clusters 1 and 4 might be considered factors that are differentially expressed with aging, but to a greater magnitude in dystrophic muscle. Conversely, Clusters 3, 4, and 5 represent proteins that are progressively differentially expressed in dystrophic muscle with more advanced pathology.

Pathway Analysis in Dystrophic Muscle throughout Pathological Progression

Dystrophic versus WT protein expression ratios at each age were analyzed separately using IPA in order to identify changes in canonical pathway and upstream regulator status associated with disease progression. Overall results were like those described above (where statistical power was greater on account of larger sample sizes) (Fig. 4C). Many canonical pathways were commonly perturbed across the different ages such as TNFR1 Signaling, PI3K/AKT Signaling, Inflammasome Pathway, Calcium Handling, Nitric Oxide Signaling in the Cardiovascular System (supplemental Fig. S10). However, other pathways exhibited apparent age-associated changes in perturbation. For example, pathways including TWEAK Signaling, Sirtuin Signaling Pathway, Phospholipase C Signaling, Signaling by Rho Family GTPases, and Integrin Signaling were upregulated at all ages but increased their degree of activation at 80 weeks. Conversely, multiple metabolism-associated pathways became progressively inactivated with disease progression, including: Ketolysis, Ketogenesis, Fatty Acid Activation, Gluconeogenesis I, and Fatty Acid β-oxidation I (among others). Pathways associated with tRNA Charging and amino acid degradation were also inhibited at 80 weeks. Interestingly, two pathways, Synaptogenesis Signaling Pathway and GP6 Signaling Pathway, exhibited marked reversals in their activation status (i.e. inactivated at 8 and 16 weeks, and activated at 80 weeks) (supplemental Fig. S11).

Similar findings were observed when considering predicted upstream regulators at each age (supplemental Fig. S12). Many of the predicted regulators identified were common to the pooled analysis described above (Fig. 4D), including IFNG, TNF, RICTOR, KDM5A, MAP4K4, INSR, TFAM, AMPK, PPARGC1A, PPARA, PPARG, and PPARD. Aging-associated changes in predicted upstream regulators were observed for TGFB1, SEMA7A, PLIN5, IKBKG, DNMT3A and TWIST1 that were activated with age, and MYC, ALDH1A2, VEGFB, TRIM24, and CPT1C that were inhibited with age (among many others) (supplemental Fig. S12). Interestingly, a number of pathways and upstream regulators were differentially modulated to a lesser extent at the 16 week age, relative to the 8 and 80 week ages (supplemental Figs. S11, S12). Such a pattern of expression is also apparent for some proteins in the heatmaps shown in Figs. 5B–6D and 6C, Cluster 2. These patterns are in accordance with the period of stabilization observed in the mdx mouse following the “crisis” period that occurs between 3 and 8 weeks of age (76).

Together these analyses have identified disease progression-associated changes in activation state for canonical pathways and upstream regulators in dystrophic muscle. Findings from the pooled analysis (Fig. 4) were generally also consistent in the analyses performed separately at each age, despite lower statistical power in the latter (supplemental Figs. S10, S11, S12).

CAV3, MVP, and PAK1 Are Upregulated In Dystrophic Muscle

Three proteins were selected for further validation of the proteomics data: (1) caveolin-3 (CAV3), a consistently elevated sarcolemmal protein, (2) major vault protein (MVP), one of 8 proteins with the highest significance in the global ANOVA analysis comparing all nine groups (adjusted P = 6 × 10−11) and (3) the serine/threonine-protein kinase, PAK1, which was one of the most highly upregulated proteins in the dystrophic group at all ages (>3-fold in 8-week mdx). These proteins were analyzed by Western blot and compared with expression levels in the HiRIEF-LC–MS/MS data (Fig. 7A, 7B). Similar results were obtained in three independent blots (one representative blot is shown for each protein, Fig. 7B). Expression patterns were comparable between HiRIEF-LC–MS/MS and Western blot experiments (Fig. 7B), thereby validating the MS data using an orthogonal protein quantification methodology. The expression of the corresponding mRNA transcripts was determined in parallel, which mirrored the protein expression data in the case of Mvp and Pak1, consistent with their transcriptional upregulation in dystrophic muscle, although some changes between WT and the dystrophic groups did not reach statistical significance at the P < 0.05 level (Fig. 7C). In contrast, Cav3 mRNA expression was not significantly increased in mdx and mdx52 samples compared with the WT at 8 and 80 weeks of age.

Fig. 7.

Validation of CAV3, MVP and PAK1 upregulation in dystrophic muscle. Relative protein expression of CAV3, MVP and PAK1 was determined by (A) HiRIEF-LC–MS/MS, and (B) Western blot. Protein loading was assessed by Coomassie Brilliant Blue (CBB) staining. C, Relative mRNA expression of Cav3, Mvp and Pak1 as determined by RT-qPCR. Genes-of-interest expression was normalized to the geometric mean of two stable reference genes: Rpl10 and Tbp. Statistical analyses are one-way ANOVA at each age, with Tukey post hoc correction. All values are mean+S.E., n = 3 for protein analysis and n = 3-4 for RNA analysis, *P < 0.05, **P < 0.01 and ***P < 0.001.

MVP is Upregulated During, but Dispensable for, Myogenic Differentiation

We next sought to further investigate the function of MVP in the commonly used C2C12 murine myoblast cell line. C2C12 cells can be propagated as undifferentiated myoblasts, but upon several days culture in low serum differentiation medium (DM) these cells fuze into syncytial myotubes (Fig. 8A). Mvp mRNA and MVP protein were expressed at low levels in proliferating myoblasts and progressively upregulated during myogenic differentiation (Fig. 8B, 8C). However, depletion of MVP expression by RNAi did not affect myogenic differentiation (Fig. 8D, 8E) despite high levels of knockdown (Fig. 8F, 8G). These data suggest that MVP may play a role in, but is dispensable for, myogenic differentiation/muscle regeneration.

Fig. 8.

MVP is upregulated during, but dispensable for, myoblast differentiation. C2C12 murine myoblasts were cultured in growth media (GM) and then switched to differentiation media (DM) for 3 and 5 days (DM3 and DM5 respectively). A, Myogenic differentiation was confirmed by immunofluorescence (IF) staining for myosin heavy chain (MHC). Expression of Mvp mRNA was determined by (B) RT-qPCR and MVP protein by (C) Western blot. C2C12 cells were transfected with 50 nm control (siCtrl) or Mvp-targeting (siMvp) siRNAs and fixed after 3 days of culture in DM. Myogenic differentiation was assessed by (D) IF staining for MHC and quantified using (E) Myogenic and Fusion Indices. Knockdown was confirmed by (F) RT-qPCR and (G) Western blot. RT-qPCR data were normalized to Rplp0 expression. GAPDH served as a protein loading control. All values are mean+S.E., n = 3 for RT-qPCR and n = 4 for IF quantification, P < 0.001, one-way ANOVA and Tukey post hoc test or student's t test as appropriate. Images were taken at 10 × magnification, scale bars indicate 400 μm.

DISCUSSION

The study of the proteomic alterations that underlie diseases such as DMD remains important for understanding disease pathology, identifying novel therapeutic targets, and discovering biomarkers of disease progression or response to therapy. In this study, we have used MS-based profiling to characterize the proteomic signature of dystrophic muscle in two different mouse models of DMD, and at different stages of the disease. Despite the technical challenges associated with proteomics analysis in skeletal muscle, we report quantification of 4974 proteins across all nine conditions in triplicate (i.e. WT, mdx, and mdx52 at 8 weeks, 16 weeks, and 80 weeks of age) (Fig. 1), the highest number of quantified proteins in dystrophic muscle described to date.

We selected the TA for proteomic analysis for several reasons, (1) we have previously generated a wealth of comparable proteomics/transcriptomics/miRNomics data in this muscle (39, 61), (2) the TA is accessible for intramuscular injection of therapeutics and toxicants allowing for comparison with data derived from such studies, and (3) experimental restoration of dystrophin restoration in this muscle is relatively straight forward (56, 57, 61, 77). The diaphragm muscle of the mdx mouse exhibits more severe degeneration (78) and so high-resolution analyses in this tissue may reveal additional insights into dystrophic pathology. Others have performed proteomics analyses in dystrophic diaphragm, although at lower depth than described here (72, 79).

Initial analyses aimed to identify proteomic changes that might explain the differences in muscle regeneration phenotype observed in the TA when comparing between the mdx and mdx52 strains at 8 and 16 weeks of age (32). However, our data suggest that the skeletal muscle proteomes of these two dystrophic mouse models are very highly similar (Figs. 2, supplemental Fig. S6, S7). Importantly, much of the proteome remains invisible on account of the peptide masking effect (49), failure to consistently detect low abundance protein-derived peptides in all samples, or peptides from some proteins falling outside the narrow pI range used for prefractionation. Based on an update of our previous calculations of protein-coding transcript detection (39), we estimate that ∼35% of the proteome was detected in the analyses described here. As such, putative proteins that are differentially expressed between mdx52 and mdx mice may be invisible to our analysis for reasons independent of statistical power. Alternatively, expression changes in specific subpopulations of cells (e.g. infiltrating immune cells) or subsets of myofibers (e.g. regenerating fibers) may not be detectable in the context of unchanged protein expression in the bulk population. It is therefore possible that the differences between the mdx and mdx52 models may be explained by expression changes in proteins that could not be detected using this methodology. Future proteomic studies should aim to understand dystrophy-associated expression changes in distinct cell/fiber subtypes.

The pooling of expression data from the mdx and mdx52 strains allowed us to identify key proteins associated with dystrophic pathology, irrespective of mutation (Fig. 3, 4). Pathway analysis was used in an effort to understand the biological relevance of differential protein expression. Mitochondrial disruption was identified by both IPA and gene list enrichment across multiple analyses (Figs. 3, 4, supplemental Fig. S10, and Fig. 6C Clusters 2 and 3). A reduction in the number of mitochondria in dystrophic muscle is a possible explanation for these findings and is further supported by the downregulation of the protein components of all five electron transport complexes. Furthermore, the analysis of protein expression across ages suggests that this perturbation of mitochondrial biology worsens with pathological progression (Fig. 3, 6C). Nevertheless, Mitochondrial Dysfunction was still observed in the 8 week-old dystrophic samples (Fig. 3, S10). Percival et al. reported a reduction in the density of the subsarcolemmal mitochondria pool in mdx TA muscle (80), and mdx-derived myoblasts contain fewer mitochondria than their WT counterparts (81), consistent with this hypothesis.

Upstream regulator analysis identified multiple master regulators of mitochondrial biogenesis as being perturbed in dystrophic muscle. These included the peroxisome proliferator activated receptor family members: PPARA, PPARG, and PPARD, their associated transcriptional co-activators: PPARGC1A (PGC-1α), and PPARGC1B (PGC-1β), and the mitochondria-specific transcription factor TFAM (Fig. 4, supplemental Fig. S10, S12). (Notably, downregulation of PGC-1α protein itself has been reported in the vastus lateralis muscle of the dystrophic GRMD dog (46)). Overlap between differentially expressed proteins enriched in PPAR, PGC-1α, PGC-1β, and TFAM was minimal (data not shown), suggesting bona fide perturbation of these distinct regulators as opposed to enrichment because of high numbers of differentially expressed proteins that are common between gene lists. Downregulation of PGC-1α activity, and the resulting loss in mitochondrial numbers/function, may therefore account for the observed inhibition of the canonical pathways associated with mitochondrial dysfunction, oxidative phosphorylation, fatty acid oxidation, and gluconeogenesis (82) (Fig. 4). Multiple factors may contribute to downregulation of PGC-1α. For example, AMPK (5′ AMP-activated protein kinase) and SIRT1 activate PGC-1α by phosphorylation and deacetylation respectively (83), both of which were predicted to be downregulated by upstream regulator analysis (Fig. 4).

Several inflammation-associated canonical pathways and upstream regulators were predicted to be upregulated, consistent with the well-described persistent inflammation that accompanies myonecrosis and regeneration in dystrophic muscle (84). These included tumor necrosis factor (TNF), interferon gamma (INFG, INF-γ) (85) and the pro-fibrotic transforming growth factor beta (TGFB1, TGF-β) (Fig. 4). Pro-inflammatory cytokines are likely derived from the infiltrating immune cells. However, TNF and INFG activate the NF-κB (Nuclear factor kappa-light-chain-enhancer of activated B cells), a transcription factor that acts as a master regulator of inflammation, which in turn results in feed-forward activation of TNF and IFNG expression in the myofibers themselves. Activation of NF-κB is a well-described feature of dystrophic muscle (13, 86, 87) but was not identified as being upregulated in our data set, possibly because of limited proteome coverage. Conversely, the NF-κB inhibitor (NFKBIA, IκBα) was paradoxically predicted to be upregulated. This is consistent with the reported auto-regulatory feedback loop by which NF-κB promotes the expression of its inhibitor (88). Interestingly, NF-κB negatively regulates PGC-1α via a direct binding interaction (89, 90), and TNF has been shown to also suppress PGC-1α in a cardiac cell model (91), suggestive of a link between inflammatory and metabolic pathologies in dystrophic muscle.

PI3K/AKT Signaling was the most highly activated canonical pathway in dystrophic muscle, consistent with previous observations in the mdx diaphragm (92). Activation of PI3K/AKT can result in muscle hypertrophy (93), suggesting that this may explain the increase in fiber size observed in mdx muscle. Similarly, AKT-mediated inactivation of FOXO transcription factors leads to the inhibition of atrophic signals, which may explain the reduced expression of MSTN in dystrophic muscle. (FOXO3 signaling was predicted to be downregulated based the data reported here).

Insulin receptor signaling was found to be the most strongly downregulated upstream regulator, suggesting that insulin signaling is not causative of PI3K/AKT pathway activation in this context. Alternative possibilities include Integrin Signaling, ILK (integrin-linked kinase), and RICTOR (a subunit of mTOR complex 2 (mTORC2)) that were all upregulated in dystrophic muscle (Figs. 4, supplemental Figs. S10, S11, S12) and are possible factors contributing to PI3K/AKT activation (94, 95). ILK signaling was also identified by pathway analysis in our previous study in mdx TA (39), and others have similarly reported ILK activation in the mdx tissues (47, 92). Additionally, there is some evidence that inflammatory cytokines such as TNF and IFNG can stimulate PI3K/AKT (96).

PI3K/AKT itself contributes to the activation of the NF-κB via the phosphorylation of IKK, which in turn inhibits IκBα by phosphorylation (92) and has been shown to directly inhibit PGC-1α via phosphorylation in liver cells (97). Taken together, these data point to interactions between metabolic, mitochondrial biogenesis, inflammatory, and muscle growth pathways in dystrophic muscle (Fig. 9).

Fig. 9.

Cross-talk between metabolic, inflammatory and muscle growth pathways in dystrophic muscle. Upstream regulator and canonical pathway analyses were integrated to generate an explanatory schema of processes occurring in dystrophic muscle. Drugs with the potential to treat DMD are shown in purple. Activation of NF-κB was not predicted based on our data but has been well-described in dystrophic muscle. ILK Signaling was enriched in dystrophic muscle, although the analysis did not determine an activation state.

Notably, many of the regulatory proteins identified above have been explored as potential targets for pharmacological manipulation (Fig. 9). For example, PGC-1α gain-of-function and pharmacological activation of PPARD have been shown to be beneficial in the mdx mouse (98), and reverse mitochondrial defects in mdx myoblasts (81), respectively. Activation of AMPK and SIRT1 (by AICAR (99) and resveratrol (100) respectively) has been shown to attenuate pathology in the mdx mouse, also a likely consequence of PGC-1α upregulation. Treatment with NF-κB inhibitors (such as CAT-1041 and edasalonexent) ameliorated dystrophic pathology in mdx mice and GRMD dogs (101) and has shown promise in a phase 1 clinical trial in DMD patients (102). Treatment with the anti-TNF antibody infliximab (brand name: Remicade) improved dystrophic pathology by reducing myonecrosis (103), and cardiac fibrosis (104) in dystrophic mice, although detrimental changes in heart function were observed attributable to inhibition of the PI3K/AKT pathway (104).

Pathway analysis with IPA is appealing as it takes expression fold changes into account, generates Z-score outputs for pathway activation state that are easily interpretable, allows for exploration of regulatory relationships, and includes the prediction of perturbed upstream regulators. However, any pathway analysis is necessarily dependent on the quality of its gene list annotations. As such, we have found that IPA frequently does not identify muscle associated pathways/regulators (data not shown). In contrast, multiple muscle-associated GO terms were identified using a less sophisticated gene list enrichment approach (55). This was particularly apparent in the list of the most differentially expressed proteins (Fig. 5A). We have therefore performed both types of analysis in this study.

A recent study by Capitanio et al. described a label free MS-based proteomics analysis of vastus lateralis muscles biopsied from DMD and BMD patients, compared with healthy controls (n = 3) (48). To determine whether the findings reported here could be recapitulated in human muscle we compared all proteins found to be differentially expressed in DMD or BMD relative to controls against all proteins found to be significantly altered in any of the dystrophic versus WT comparisons described in this study. 135 proteins were commonly identified between DMD/BMD patients and dystrophic mice, thereby providing independent evidence that many of the findings reported here are relevant to the situation in dystrophic patient muscle (supplemental Fig. S13). Indeed, this level of overlap is perhaps surprising when considering the differences in organism, muscle type, mutation type, disease-stage progression, and the methodologies used. Conversely, 91 proteins uniquely identified in the human study, of these 48 were also detected in our datasets, but were not differentially expressed, with the remaining 43 proteins invisible in our analyses (supplemental Fig. S13). 2458 proteins were uniquely identified as being differentially expressed in the dystrophic mouse samples, which emphasizes the benefit of the high-resolution proteomics approach used in this study, in terms of deep proteome coverage. Differentially expressed proteins that were common to dystrophic human and mouse muscle were enriched for GO terms associated with muscle structure, the cytoskeleton, and energy production (supplemental Fig. S13) consistent with findings described above (Fig. 3). Proteins that were uniquely differentially expressed in the human patient data set were enriched for antioxidant activity, peroxidase activity, exocytosis, and vesicle structure and function (supplemental Fig. S13). Notably the P-values for these enriched GO terms were relatively low. Proteins that were uniquely differentially expressed in our murine datasets were enriched for RNA binding, multiple metabolic processes, and the mitochondrion (supplemental Fig. ) similar to findings described above (Fig. 3). P-values for the enrichment of these GO terms were comparatively very high, which is to be expected when considering the greater number of proteins available for analysis.

A possible limitation of this study is the influence of batch effects, as experimental samples for each age were split between three separate pools of TMT-labeled peptides. This study design ensured that statistical tests performed within each age are not influenced by batch effects. However, comparisons between ages may be subject to this technical source of variation, although there is evidence that such effects are minimal in our data. Firstly, the normalization of TMT ratios to the bridge sample minimizes batch effects to some extent, and the proteomics performance for the bridge sample was highly similar for each TMT10 set (supplemental Fig. S3). Secondly, comparing TMT ratios across all samples for known proteins that are differentially expressed in dystrophic muscle shows minimal differences between ages/TMT10 sets (supplemental Fig. S5), and the majority of proteins are not changed over time (Fig. 3, 4, supplemental Fig. S6, File S1). Thirdly, proteins/pathways that are altered with age are consistent with established features of dystrophic pathology (such as the accumulation of fibrotic damage) (Fig. 6, supplemental Fig. S10, S11, S12).

This study has highlighted multiple proteins of interest, of which we selected MVP for further investigation. MVP is the major constituent of vault particles, ribonucleoprotein structures (the function of which is poorly understood) (105). MVP was upregulated during, but dispensable for, myogenic differentiation (Figs. 7, 8) suggesting that MVP upregulation may be associated with muscle regeneration in dystrophic muscle. Alternatively, MVP may reflect a different facet of dystrophic pathology, as this protein has been shown to be elevated in response to serum-starvation-induced apoptotic stress (106). Notably, the other protein components of vault particles (i.e. TEP1 and PARP4) were similarly upregulated in dystrophic muscle, although to a lesser extent, suggesting an increase in the number of vaults may explain these data. Recent evidence has shown that vault particles can act as scaffolds and are involved in cell survival signaling (107, 108). Moreover, they have been reported to regulate growth and survival of respiratory smooth muscle cells (108). To our knowledge, this is the first study demonstrating an association between MVP/vault particles and dystrophic pathology, which is deserving of further investigation.

Using two dystrophic mouse models and the inclusion of different disease time points, this study provides a wider perspective of, and offers new insights into, the pathological mechanisms involved in Duchenne muscular dystrophy. The data presented here describe a wealth of novel differentially expressed proteins, perturbed biological pathways, and predicted upstream regulators. This data set will constitute a valuable resource for the DMD research community, and for those studying the effects of aging on the WT muscle proteome.

DATA AVAILABILITY

Protein expression data and statistical analysis are summarized in File S1. All other data relating to this manuscript are freely available on request from the authors. The MS raw proteomics data have been deposited to the ProteomeXchange Consortium via the JPOST partner repository with the data set identifier PXD017169 (JPST000734).