A cellular and molecular spatial atlas of dystrophic muscle

Abstract

Asynchronous skeletal muscle degeneration/regeneration is a hallmark feature of Duchenne muscular dystrophy (DMD); however, traditional -omics technologies that lack spatial context make it difficult to study the biological mechanisms of how asynchronous regeneration contributes to disease progression. Here, using the severely dystrophic D2-mdx mouse model, we generated a high-resolution cellular and molecular spatial atlas of dystrophic muscle by integrating spatial transcriptomics and single-cell RNAseq datasets. Unbiased clustering revealed nonuniform distribution of unique cell populations throughout D2-mdx muscle that were associated with multiple regenerative timepoints, demonstrating that this model faithfully recapitulates the asynchronous regeneration observed in human DMD muscle. By probing spatiotemporal gene expression signatures, we found that propagation of inflammatory and fibrotic signals from locally damaged areas contributes to widespread pathology and that querying expression signatures within discrete microenvironments can identify targetable pathways for DMD therapy. Overall, this spatial atlas of dystrophic muscle provides a valuable resource for studying DMD disease biology and therapeutic target discovery.

Keywords: Duchenne muscular dystrophy; asynchronous regeneration; skeletal muscle; spatial transcriptomics.

Fig. 1.

Spatial transcriptomics profiling of skeletal muscle identifies multiple cell clusters within distinct locations. (A) Workflow of spatial gene expression profiling of mouse gastrocnemius/plantaris muscles. (B) H&E-stained WT hindlimb muscle from one representative sample (Left) with Visium spot transcriptome clusters visualized on tissue-covered area (Center) and as a UMAP plot (Right). The 11 unique clusters were identified from n = 3 WT muscle sections. (C) Heatmap showing top canonical markers defining each cell cluster. (D) Myofiber cluster identity displays high overlap with (E) spatial patterning of corresponding myosin heavy chain gene expression and (F) myosin heavy chain isoform protein localization in serial sections. (G) Nerve and NMJ clusters express high levels of (H) myelinating genes (e.g., Mpz) and postsynaptic NMJ genes (e.g., Chrna1), respectively. (I) Serial section immunostaining depicting enrichment for neurofilament-stained nerve bundles and α-bungarotoxin-stained NMJs. Arrows indicate NMJs within NMJ cluster areas as defined by spatial transcriptomics. (J) Fibro./Teno. clusters are enriched at the epimysium surrounding and separating the gastrocnemius/plantaris, and BV + PVAT clusters are enriched at the superficial gastrocnemius region associated with the saphenous vein. (K) Immunostaining of serial sections depicting large CD31+ blood vessels (arrow) at the muscle periphery with surrounding adipose tissue identified by perilipin staining and enrichment for adipocyte gene expression (e.g., Fasn) in BV + PVAT cluster areas. (L) Fibro./Teno. cluster regions show enrichment for ECM proteins, as assessed by wheat germ agglutinin staining of serial sections and by high expression of ECM genes (e.g., Col1a1).

Fig. 2.

Mapping asynchronous regeneration of dystrophic muscle by spatial gene profiling. (A) Unbiased clustering of gene expression profiling identified 13 unique clusters in MDX hindlimb muscle (n = 5) that were mapped to H&E-stained tissue sections. (B) Distribution of immune and regenerative clusters among each MDX sample, with the total spot proportion of these four clusters denoted in parentheses. (C) Heatmap showing top canonical markers defining each cell cluster. (D) Top five Gene Ontology pathways (biological process) enriched in Immune 1 and Immune 2 clusters. Gene ratio is the total number of DEGs divided by the total number of genes in each geneset, and the adjusted P value is based on gene set enrichment analysis. (E) Volcano plot highlighting differential expression of immune (red, blue), myofilament (green), and ECM genes (purple) between Immune 1 and Immune 2 clusters. (F) Cell cluster localization of Immune 1, Immune 2, and Regn. fiber clusters with serial immunostaining for markers of damaged myofibers (mIgG) and macrophages (F4/80). (G) Immunostaining for pro- (CD68) and anti- (CD206) inflammatory macrophages in Immune 1 and Immune 2 cluster areas. (H) Schematic depicting the state of regeneration of spatially defined MDX clusters, with granulocytes depicted in purple and macrophages depicted in orange.

Fig. 3.

Gene expression signatures from unique dystrophic muscle areas align with distinct timepoints of regeneration. (A) Schematic depicting focal cardiotoxin injury of WT mice and integration of MDX cluster gene signatures with spatial transcriptomics datasets obtained at multiple timepoints following injury. (B) Immunostaining images depicting focal areas of cardiotoxin-induced damage (dotted lines) at different timepoints and expression scores of MDX Immune 1, Immune 2, and Regn. fiber clusters plotted on injured tissue at each timepoint.

Fig. 4.

Single-cell deconvolution of spatial transcriptomics data in dystrophic muscle. (A) Reclustering of skeletal muscle single-cell RNAseq from De Micheli et al. (20) and integration with spatial transcriptomics of dystrophic muscle. (B) Chord diagram displaying co-occurrence of inferred cell types in dystrophic muscle. (C) Heatmap showing the proportion of inferred cell types present within spatial clusters (excluding myofibers). Representative spatial localization and abundance of inferred cell types within Immune 1 (D), Immune 2 (E), and Regn. fiber (F) areas of dystrophic muscle.

Fig. 5.

Fibrosis and inflammation within dystrophic muscle is propagated from areas of localized damage. (A) Fiber-type immunostaining and Type IIb fiber cluster localization depicting undamaged Type IIb fiber areas in WT and MDX tissue that were selected for differential gene expression analysis. (B) Top five Gene Ontology (biological process) pathways up-regulated in MDX Type IIb fiber clusters compared to WT. Gene ratio is the total number of DEGs divided by the total number of genes each geneset, and the adjusted P value is based on gene set enrichment analysis. (C) Schematic depicting classification of spatial transcriptomics node spots and neighborhood set spots for Squidpy nearest neighbor analysis. (D) Representative neighborhood enrichment analysis between clusters within the MDX3 sample. (E) Cluster pairs with positive z-score neighborhood enrichment among all samples containing indicated clusters (n = 4 to 5) (F) Localization of Immune 2 and MPC clusters and quantification of the proportion of peripheral Immune 2 spots directly adjacent to MPC spots among all 5 MDX samples (mean ± SEM). (G) Immunostaining for collagen I in serial sections corresponding to Immune 2 and MPC cluster areas. Circled areas align with Immune 2 clusters. (H) Quantification of collagen I immunostaining in Immune 2 areas compared with the remainder of the tissue section (n = 5; bars indicate mean ± SEM; ***P < 0.001, unpaired t test). (I) Quantification of the number of MPC cluster spots in WT and MDX muscle (n = 5; bars indicate mean ± SEM; ***P < 0.001, unpaired t test). (J) Immunostaining for localization of Pdgfra⁺ cells within Immune 2/MPC areas. (K) Localization of Immune 1 and Fiber + inflamm. clusters and quantification of the proportion of peripheral Immune 1 spots directly adjacent to Fiber + inflamm. spots among MDX samples (n = 4; mean ± SEM). (L) RNAscope imaging for Ccl2, Ccl7, and Spp1 mRNA within Immune 1/Fiber + inflamm. areas.