Transcriptional and cytopathological hallmarks of FSHD in chronic DUX4-expressing mice

Facioscapulohumeral muscular dystrophy (FSHD) is a prevalent (1) and slowly progressive myopathy caused by genetic mutations (2, 3), leading to the loss of epigenetic repression (4) of tandem copies of the DUX4 gene (5) embedded within the D4Z4 macrosatellite repeat at 4qter. A sequence downstream of the repeat array, specific to a disease-permissive allele, provides an atypical polyadenylation signal to the terminal copy of DUX4, presumably leading to the accumulation of DUX4 mRNA and expression of DUX4 protein (6). Expression levels of DUX4 are extremely low, as monoclonal antibodies to date have failed to detect DUX4 by immunostaining of FSHD biopsies. However, the presence of DUX4 can be inferred by transcriptional profiling of biopsies, where average levels of a panel of DUX4 target genes are elevated in FSHD compared with those in control biopsies (7).

DUX4 expression has been detected in rare nuclei of FSHD myoblasts; the first such study showed approximately 1:1,000 DUX4⁺ nuclei (8), and a follow-up study showed even lower frequencies (9). After differentiation, DUX4 expression was found to be increased, in part due to fusion allowing a single transcribing nucleus to provide DUX4 protein to neighboring nuclei, but also due to elevated levels of expression from D4Z4, while a majority of nuclei still did not show DUX4 expression (10). The most sensitive readout of DUX4 expression in primary cells from patients with FSHD to date used a DUX4 fluorescent reporter, where most cell lines showed approximately 1%–2% DUX4⁺ cells by FACS and an increase with differentiation of up to 3% (11). These studies on primary cells reinforce the notion that FSHD involves infrequent DUX4 expression; thus, animal models aiming to recapitulate the pathology of FSHD need to express DUX4 in a similar infrequent manner.

At the individual cell level, DUX4 expression at high levels leads to the death of myoblasts (12) and other cell types (13, 14). Transcriptional changes are rapid, occurring within 4 hours of DUX4 expression (13), and include both upregulation of genes with nearby DUX4-binding sites (15) as well as downregulation of large numbers of genes (16). Efficient DNA-binding requires interaction of both N-terminal homeodomains with DNA (17, 18), and, at a large subset of its binding sites, DUX4 acts as a pioneer factor, using its C-terminus to recruit EP300 and CBP, which acetylate local histones and, thereby, open up previously inaccessible domains (19, 20).

There are now 3 mouse models that show muscle pathology in response to muscle-specific DUX4 expression. The first of these used doxycycline (dox) induction via a muscle fiber–restricted rtTA (21), while the 2 more recent models use muscle fiber–specific Cre-lox recombination (22, 23) to express DUX4 in muscle fibers. Although Cre-based models can give low overall expression through use of an inefficient Cre, the per-cell expression level in those cells that do express is relatively high due to the use of a strong promoter upstream of the floxed stop codon. On the other hand, the iDUX4pA model expresses DUX4 inefficiently and only in rare myofiber nuclei, a situation more closely resembling that of patients with FSHD, making the iDUX4pA model particularly attractive. In the current study, we used the dox-inducible iDUX4pA mouse (21) to develop a long-term chronic disease model based on rare, stochastic, long-term DUX4 expression and evaluated the cytolopathological features of this slowly progressive disease. These studies reveal compelling correlations between the histology of diseased murine muscle and muscle tissue of patients with FSHD, and a remarkable similarity between gene expression of diseased muscle of the iDUX4pA-HSA mouse and that of patients, and point to a critical involvement of fibroadiopogenic progenitors (FAPs) in the pathological process driving FSHD.

Long-term DUX4 induction promotes a progressive muscular dystrophy. Previous work to induce the iDUX4pA transgene via the HSA-rtTA driver used daily intraperitoneal injection of dox. To facilitate long-term experiments and minimize stress associated with daily injections, we tested the provision of dox chow. Female iDUX4pA-HSA mice, controls lacking the HSA-rtTA transgene, and WT sibling controls were fed a nutrient chow containing dox for a period of 6 months, starting from 4 weeks of age. This method of induction was well-tolerated; however, iDUX4pA-HSA mice showed a slow progressive muscle atrophy and a deterioration in muscle function, while their WT siblings did not. At the end of 6 months, induced mice displayed a markedly kyphotic posture, and X-ray imaging revealed a prominent outward spinal curvature (Figure 1A). While body weight was slightly reduced compared with controls (Figure 1B), muscles showed pronounced atrophy (Figure 1, C and D). We evaluated muscle function both in live mice as well as on isolated muscles. Muscle strength was evaluated in vivo by measuring torque generated by ankle dorsiflexor muscles following subcutaneous electrical stimulation via the common peroneal nerve. iDUX4pA-HSA mice displayed much lower absolute torque, as expected based on their extensive muscle wasting. Interestingly, even after normalizing to tibialis anterior (TA) muscle weight, iDUX4pA-HSA muscles displayed a significantly reduced relative torque, indicating a clear decline in the quality of force generated by the TA muscle (Figure 1E and Supplemental Figure 1, A and B; supplemental material available online with this article; https://doi.org/10.1172/JCI133303DS1).

Figure 1

Muscle phenotype in iDUX4pA-HSA mice induced with doxycycline for 6 months. (A) X-ray image of WT and iDUX4pA-HSA mice fed with doxycycline chow for 6 months. (B) Body weight (g) of the female mice at the time of analyses at 6 months. (C) Gross appearance of quadriceps from WT (left) and iDUX4pA-HSA (right) mice at the time of dissection. (D) Mass of different muscles normalized to the body weight at 6 months. (E) Isometric in vivo tetanic torque and torque-frequency relationship of the anterior crural muscles expressed relative to the mass of the TA muscle. n = 8. (F) Ex vivo tetanic forces expressed relative to the physiological cross-sectional area of the EDL muscle, the isometric tetanic force with representative force-time tracings, and passive stiffness. n = 8. (G) Representative H&E staining of sections from quadriceps from WT and iDUX4pA-HSA mice. Scale bar: 100 μm. (H) Immunofluorescence staining for DUX4 (red) and laminin (green) in quadriceps. Tissues sections were counterstained with DAPI (nuclei, blue). Scale bar: 40 μm (inset magnification, 4:1). Arrows highlight rare DUX4⁺ cells. (I) Representative FACS profile for CD45⁺/CD68⁺, CD45⁺/CD206⁺, and CD45⁺/Gr1⁺ cells from the pool of muscle at 6 months after induction. (J) Graphical summary analyses of the inflammatory cells presented in A. Data are shown as mean ± SEM; *P < 0.05; **P < 0.01; ***P < 0.001; ****P < 0.0001 by t test; n = 8.

We then performed ex vivo physiological evaluation of isolated extensor digitorum longus (EDL) muscles. This revealed marked differences in parameters of force generation between iDUX4pA-HSA mice induced for 6 months and dox-treated WT sibling controls. As expected based on the degree of muscle wasting, maximal tetanic forces were reduced by more than half (Figure 1F). When normalized to physiological cross-sectional area, isometric, concentric, and eccentric specific forces were all significantly reduced (Figure 1F and Supplemental Figure 1, C and D), again indicating that muscles suffered not only from a decrease in total mass, but also from a decrease in the intrinsic force-generating capacity of what muscle remained. Further supporting the gross physiological differences driven by long-term low-level DUX4 expression, the rate of contraction by the EDL was reduced by approximately half, and passive stiffness was elevated by approximately one-third (Figure 1F). We then subjected EDL muscles to a series of eccentric contractions and measured the force generated at each contraction. WT muscle displayed no significant decline in force generation over 5 contractions. iDUX4pA-HSA muscles demonstrated a modest decline in eccentric force with each contraction as well as reduced isometric force (Supplemental Figure 1E). Using the same in vivo and ex vivo approaches, we compared muscle strength of iDUX4pA mice lacking the HSA-rtTA transgene and did not find any significant difference compared with WT mice (Supplemental Figure 2).

Histologically, all muscles evaluated were markedly dystrophic, with signs of cellular infiltrate, fibers of varying sizes, occasional fat deposition, and plentiful deposition of ECM (Figure 1G). Immunostaining revealed rare DUX4⁺ myofibers in the affected areas of the muscle (Figure 1H, arrows). When exposure was increased to the point where background staining arose, additional spots of very dim positivity that colocalized with DAPI became apparent in the DUX4 samples (Figure 1H). The difficulty of distinguishing such spots from background precludes their quantification; however, a few clearly DUX4⁺ nuclei were present in most muscle cross sections.

To understand the cellular infiltrate, we performed FACS on hind limb muscle cell suspensions of animals induced for 6 months and quantified the presence and abundance of various nonmyogenic cell types, focusing initially on inflammatory markers. Muscle induced for 6 months displayed an approximately 2-fold elevation in frequency of Gr1⁺, CD68⁺, and CD206⁺ hematopoietic cells (Figure 1, I and J). Thus, low level stochastic chronic DUX4 induction in iDUX4pA-HSA promoted progressive muscular dystrophy characterized by loss of muscle mass, decreased contractile ability, and infiltration of connective tissue, fat, and inflammatory cells.

Dystrophic muscles from iDUX4pA-HSA mice recapitulate the gene expression signature of FSHD. In order to better understand disease in the iDUX4pA-HSA model, we performed transcriptional profiling by RNA-seq on whole muscle from mice on dox for a short term (2 weeks) or a long term (16 weeks). This revealed a large set of differentially expressed genes (DEGs) (Figure 2A), with principle component analysis (PCA) demonstrating primary axes correlating with DUX4 status (x axis) and age (y axis), and diseased muscles clustering closely independent of age (Figure 2B). As expected based on the PCA, a large number of DEGs (807 genes) were common to the short-term and long-term DUX4-induced muscle, predominantly genes from cell death–related pathways. However, clear evidence of progression of disease can be discerned in the changes specific to earlier or later time points, in which pathway analysis pointed to greater dysregulation over time of pathways related to both immune response and extracellular matrix (ECM) biology (Figure 2C).

Figure 2

Differentially expressed genes in muscle from iDUX4pA-HSA mice. (A) Heatmap of differentially expressed genes (DEGs) in muscles from control and iDUX4pA-HSA mice induced for 2 and 16 weeks. (B) PC analyses on genes expressed in muscle from control and iDUX4pA-HSA mice induced for 2 and 16 weeks. Note that samples from iDUX4pA-HSA mice cluster together despite age-specific gene expression profile. (C) Venn diagram shows specific and common differentially expressed genes compared with controls at 2 and 16 weeks of induction. KEGG pathway analyses revealed uniquely and commonly affected pathways.

We also found that known DUX4 target genes were upregulated. Using the top 50 upregulated and 50 downregulated DUX4 target genes in mice (13), we derived a composite DUX4 score and found it to be significantly upregulated at both 2 and 16 weeks in the iDUX4pA-HSA samples (Figure 3A). To relate molecular changes in this mouse model to FSHD in humans, we compared the DEGs in our data set to human data obtained from 2 independent studies of MRI-guided biopsies of patients with FSHD (24, 25). Because these data are from muscle tissue identified by MRI as T2-STIR2⁺ (short TI inversion recovery), they represent muscles currently experiencing active disease. We identified the sets of genes differentially expressed at 2 weeks and 16 weeks of DUX4 expression, identified the human homolog of each, and tested for enrichment of these gene sets in each human study. This revealed that the human homologs of DEGs in the diseased mouse muscle showed extremely significant enrichment in the 2 independent human studies (Figure 3, B and C). Notably, enrichment was greater in mice that were induced for an extended period of time (16 weeks) compared with those induced for 2 weeks (Supplemental Figure 3). The genes driving these enrichments were mostly shared between the Wang (25) and Tasca (24) data sets and were involved in cellular processes related to ECM, adhesion, and migration, including many specific cancer pathways (Figure 3, D and E, and Supplemental Figure 3A). This remarkable result indicates that active disease processes mechanistically represented at the RNA level in human FSHD are largely recapitulated in the iDUX4-HSA mouse model.

Figure 3

FSHD signature gene expression in muscle from iDUX4pA-HSA mice. (A) The DUX4 score shows enrichment of DUX4 early target genes in muscles from iDUX4pA-HSA mice induced for 2 and 16 weeks. A list of DUX4 early affected genes was identified in iC2C12-DUX4 cells after 4 hours of DUX4 induction (13). (B) Gene set enrichment analysis (GSEA) of differentially expressed genes (DEG) at 2 and 16 weeks in induced mice and FSHD biopsies identified by Wang et al. (25). (C) GSEA of differentially expressed genes at 2 and 16 weeks in induced mice and FSHD biopsies identified by Tasca et al. (24). (D) Venn diagram for common differentially expressed genes in patients with FSHD from the Wang et al. data set (25) and iDUX4pA-HSA mice at 2 and 16 weeks and top enriched KEGG pathways. (E) Venn diagram shows the number of specific and commonly expressed genes in patients with FSHD from the Tasca data set (24) and iDUX4pA-HSA mouse, and top KEGG pathways.

Vasculature defects in iDUX4pA-HSA muscle. Preliminary FACS analysis of diseased iDUX4pA-HSA muscle revealed markedly reduced Lin⁺ cell numbers. We therefore tested the Lin components (CD45 and CD31) individually and found that CD31⁺ endothelial cells and CD31⁺CD146⁺ pericytes were significantly reduced in frequency in iDUX4pA-HSA muscle (Figure 4, A and B). The reduction was dose dependent and not observed at the lowest dose of dox chow tested (Figure 4B). To determine what aspect of vasculature was affected, we performed immunostaining of TA muscle sections with CD31. This revealed a reduction of capillary density to approximately 60% of controls (Figure 4, C and D). The pattern of capillary loss was not uniform throughout the muscle but rather was localized to foci, within which many fibers were apparently devoid of any adjacent capillary. This result is intriguing, as histological analysis of FSHD muscle biopsies has similarly shown a reduction in capillary number (26).

Figure 4

Compromised muscle vasculature in iDUX4pA-HSA mice. (A) Representative FACS profile for CD45^–CD31⁺ and CD45^–CD146⁺CD31⁺ cells in the pool of muscles (gastrocnemius, quadriceps, pectoralis) of WT and iDUX4pA-HSA mice at 3 months after induction with a 10× serial dilution of doxycycline. (B) FACS analyses for CD45^–CD31⁺, CD45^–CD146⁺, and CD45^–CD31⁺CD146⁺ cells presented in A. Data are shown as mean ± SEM; *P < 0.05; ***P < 0.001; ****P < 0.0001 by 1-way ANOVA, n = 4. (C) Immunofluorescence staining for CD31 (green), laminin (red), and nuclei (DAPI, blue) in quadriceps at 6 months after induction with chow containing 625 mg/kg doxycycline. Scale bar: 100 μm. (D) Quantification of CD31⁺ myofibers in sections presented in C. ****P < 0.0001. (E) 3D render and maximum intensity projections (MIP) of TA muscle from WT and iDUX4pA-HSA mice at 6 months. Disorganized vasculature with nonperfused segments of the vasculature (indicated by arrows) in iDUX4pA-HSA muscle. Scale bar: 100 μm. (F) Heatmap of differentially expressed genes in CD45^–Pdgfrα^–CD31⁺ FACS-sorted cells from control (n = 3) and iDUX4pA-HSA (n = 4) mice after 6 months induction with doxycycline. (G) The DUX4 score shows no enrichment of DUX4 early target genes in CD31⁺ cells from iDUX4pA-HSA mice. (H) MA plot showing the fold change of gene expression on samples presented in F. (I) KEGG pathway analyses showing the most enriched affected processes in CD31 cells from iDUX4pA-HSA mice.

To comprehensively visualize the morphology of capillary beds in iDUX4pA-HSA muscle, we perfused mice with fluorescent lectin in vivo, isolated TA muscles, and subjected them to histological clearing and confocal microscopy. As shown in Figure 4E, the capillary beds of iDUX4pA-HSA muscles showed severe perturbation, with many capillaries breaking from their typical parallel orientation. The extensive muscle damage likely allows lectin to extravase, allowing staining outside of capillaries. Therefore, this approach is suitable for gross visualization but not for quantification; but in spite of this, focal zones devoid of perfusion were clearly evident (Figure 4E, arrows).

To investigate the cell autonomous changes in endothelial cells from the affected muscle and to determine whether these might be caused by leakage of the inducible system into endothelial cells, we performed RNA-seq on FACS-sorted CD45^–CD31⁺ cells from dystrophic iDUX4pA-HSA muscle induced for 6 months and similar muscle from control siblings, also treated for 6 months with dox chow. Upregulated genes predominated among the DEGs (721 of 988) in CD31⁺ cells from dystrophic muscle (Figure 4F and Supplemental Table 1). These did not include DUX4 targets, demonstrating that alterations in endothelial cells are in response to DUX4 expression in fibers, not to leaky DUX4 expression in endothelial cells (Figure 4G). The most affected pathway alterations included the not unexpected pathways PI3K/AKT and cell-cell interaction as well as less expected processes related to protein degradation and phagosome function (Figure 4, H and I). Thus, we conclude that DUX4 expression in fibers leads to vasculature defects, including abnormal capillary conformation, reduced capillaries, and an overall reduction in the number of endothelial cells and pericytes.

FAP expansion is a hallmark of iDUX4pA-HSA dystrophic muscle. Emerging evidence suggests a crucial role of FAPs in the regeneration of the muscle after acute injury as well as for fibrosis during chronic muscle damage (27–30). We previously observed that FAPs were increased after 2 weeks of high-level DUX4 induction, a regimen that produced acute severe muscle damage (21). We therefore evaluated FAP infiltration in iDUX4pA-HSA muscle after long-term induction. Mice were fed with dox chow for 3 months with various concentration of dox, and FAPs were quantified by FACS as CD45^–/CD31^–/Itga7^– and Sca1⁺ or Pdgfrα⁺ cells. Although Pdgfrα is the more specific marker, as Sca1 is also expressed on other cell populations, and indeed there may be heterogeneity in the FAP compartment, both markers are commonly used to identify FAPs (28, 29). This revealed a remarkable dose-dependent increase in Lin^–Sca1⁺ or Lin^–Pdgfrα⁺ cells and decrease of Lin^–Itga7⁺VCAM⁺ myogenic progenitors (Figure 5, A and B), correlating with the degree of muscle atrophy (Supplemental Figure 4A). Increased numbers of FAPs were also associated with an increase in inflammatory cells (Supplemental Figure 4B). To address the physiological effect of this increase in number of FAPs, we quantified the extent of fibrosis by Sirius Red/Fast Green staining. This revealed a dramatic elevation in collagen deposition between the affected myofibers (Figure 5, C and D). The dose-dependent DUX4-provoked inflammation and fibrosis were also confirmed by detection of specific markers by quantitative real-time RT-PCR (RTqPCR) (Supplemental Figure 4C).

Figure 5

Fibrosis and FAP infiltration in DUX4 affected muscle. (A) Representative FACS profiles for CD45^–CD31^–Itga7^–Pdgfrα⁺, CD45^–CD31^–Itga7^–Sca1⁺, and CD45^–CD31^–Itga7⁺VCAM⁺ cells in skeletal muscle from iDUX4pA-HSA mice induced with various concentration of doxycycline for 3 months. Doxycycline was continuously administrated to the mice through custom-made food containing dilutions of doxycycline. (B) Summary of FACS analyses on the cells presented in A. Data are shown as mean ± SEM; **P < 0.01; ***P < 0.001; ****P < 0.0001 by 1-way ANOVA, n = 4. (C) Representative images of Sirius red/Fast green staining of quadriceps from WT and iDUX4pA-HSA mice at 6 months. Scale bar: 100 μm. (D) Quantification of deposition of fibrous tissue in the muscle presented in C. SRFG, Sirius Red/Fast Green. Data are shown as mean ± SEM; ***P < 0.001 by t test, n = 4.

We wished to compare the behavior of FAPs in muscle injury caused by DUX4 to that caused by cardiotoxin and glycerol. Intramuscular injection of cardiotoxin causes rapid myofiber death, rapid regeneration within 10 days, and normal muscle histology (other than centrally located nuclei) by 30 days. On the other hand, glycerol induces injury that is characterized by some degree of fibrosis and fat deposition. We compared TA muscles from the iDUX4pA-HSA mice that were fed for 10 days with dox to those of WT mice injected with cardiotoxin or glycerol (Supplemental Figure 5A). At 10 days after injury, we found a divergent effect on muscle mass; TAs injured by glycerol were hypertrophic while DUX4 induced TAs were atrophic (Supplemental Figure 5B). However, regarding FAP infiltration, the glycerol- and DUX4-induced muscle injuries were similar, with both increasing numbers of Pdgfrα⁺ and Sca1⁺ cells at the 10 day time point, while FAPs were not significantly different at this time point in CTX-injured muscle (Supplemental Figure 5C). It was previously shown that FAP numbers increase but then return to normal within 5 days of myotoxin injury (28); thus, the sustained increase in FAPs is suggestive of an altered profibrotic/adipogenic physiological state.

Chronic low-level DUX4 expression in fibers drives FAPs into an altered pathological state. To more deeply understand the changes occurring in the FAP compartment, we performed transcriptional profiling on sorted Lin^–Pdgfrα⁺ cells from iDUX4pA-HSA mice and littermate controls fed dox chow for 10 days or 6 months. We detected around 750 differently expressed genes in both groups; however, only 169 of these DEGs were shared between the 2 time points (Figure 6, A and B), indicating that FAPs were evolving significantly in the context of chronic dystrophy. In the chronic context, a large number of DEGs were identified; about twice as many genes were downregulated compared with those upregulated with disease (Figure 6A and Supplemental Tables 2 and 3). Pathway analysis suggested major changes in ECM interaction and motility factors, and downregulation of antiinflammatory secreted factors, consistent with these cells acquiring a proinflammatory state (Figure 6B). We compared the profile of DEGs in FAPs from iDUX4pA-HSA chronic disease to those identified in a recent study evaluating FAPs in the context of lytic muscle injury or atrophy driven by denervation (31). Although the majority of DEGs were unique to the iDUX4pA-HSA chronic disease, a large subset were shared. These shared DEGs pointed to changes in pathways involved in focal adhesion and ECM interactions as well as motility (specifically in common with the cardiotoxin-respondent FAPs) (Supplemental Figure 6A). These results demonstrate that FAPs, in addition to being increased in number in dystrophic muscles of chronically induced iDUX4pA-HSA mice, have acquired an altered, presumably pathological, state, that may further contribute to muscle degradation. It was surprising that the dysregulation of FAPs in iDUX4pA-HSA affected muscle showed similarities to both FAPs promoting regeneration following acute injury and FAPs from a chronic atrophy model. This suggests that both types of FAPs may be simultaneously present and that rather than the appropriate sequence of activity that restores muscle to its preinjured state, a temporally confused expression profile locks FAPs into a state that impairs regeneration and promotes fibrosis, fat accumulation, and atrophy.

Figure 6

Genes affected in FAPs from the iDUX4pA-HSA mouse are enriched in DEG sets from patients with FSHD. (A) Heatmap represents a clustering of differentially expressed genes in FAP (CD45^–CD31^–Pdgfrα⁺) FACS-sorted cells from control and iDUX4pA-HSA mice after 10 days and 6 months induction with doxycycline (625 mg/kg). (B) Venn diagram and KEGG pathway analyses on the samples presented at A. (C) Gene set enrichment analysis (GSEA) of differentially expressed genes (DEGs) in PDGFRα⁺ cells isolated from mice induced for 10 days and DEGs identified in FSHD biopsies in 2 different studies (24, 25). (D) GSEA in DEGs identified in PDGFRα⁺ cells from mice induced for 6 months and FSHD samples. (E) Venn diagrams for core-enriched genes identified in GSEA analyses presented in C and D. (F) Venn diagram presents the overlap of core-enriched genes identified in FSHD biopsies (Tasca, ref. 24, and Wang, ref. 25, data sets) and DEGs from whole muscle from mice treated with doxycycline for 2 and 16 weeks and FAPs isolated from acute (10 days) and chronic (6 months) induced mice. Genes that were identified that were differentially expressed in FAPs from the FSHD animal model are highly enriched in DEG sets from muscle of patients with FSHD. (G) The DUX4 score shows no enrichment of DUX4 early target genes in PDGFRα cells from iDUX4pA-HSA mice.

Genes overexpressed in FAPs during chronic disease are enriched in human FSHD biopsies. To investigate the relevance of these gene expression changes in FAPs to FSHD in humans, we again selected the set of human orthologs of DEGs identified in FAPs from iDUX4pA-HSA acute and chronic disease and evaluated this set of genes in transcriptional profiling data from MRI-guided human FSHD biopsies. Remarkably, this set of genes showed highly significant enrichment in both independent MRI-guided FSHD biopsy transcriptional profiling studies published to date (refs. 24, 25 and Figure 6, C and D). When we evaluate the number of genes driving the gene set enrichment analyses (GSEA) scores of whole muscle versus FAPs, i.e., those genes that are commonly regulated between human and mouse whole muscle or human and mouse FAPs, we found that while the signal was stronger in whole muscle, the signal was respectable in FAPs, within a factor of 2 (Figure 6, E and F; Supplemental Figure 6B; and Supplemental Figure 7). Importantly, DUX4 target genes were not elevated in FAPs from diseased muscle (Figure 6G), demonstrating that, similar to what was observed in endothelial cells, alterations in FAPs are the indirect consequence of DUX4 expression in fibers, not of leaky DUX4 expression in FAPs. These data strongly suggest that gene expression changes in FAPs in the iDUX4pA-HSA are relevant to FSHD in humans.

FSHD is a disease that is typically noticed in the second decade of life and involves a very long period of progressive muscle wasting. This study establishes a disease model that mirrors the long time frames of FSHD in humans. By treating with dox chow, which facilitates long-term induction of DUX4, mice can be exposed to disease for periods exceeding half the normal life span of a mouse (the longest running cohort is currently at 1.2 years on dox chow). Interestingly, a second feature of the human disease is resembled in these animals: the difficulty of detecting DUX4 at the protein level by immunostaining of muscle. DUX4 expression is both infrequent and low level to the point of being almost undetectable. Nevertheless, the muscular dystrophy is profound and persistent, with many histological similarities to FSHD, including increased variance in fiber size, infiltration of mononuclear cells, focal damage, and regenerated fibers, all of which are nonspecific and infrequent but detectable in FSHD, as well as reduced capillary density, fibrosis, and sporadic fat deposition, all of which are prevalent in FSHD.

In addition to temporal and histological relevance, the long-term iDUX4pA-HSA model shows impaired contractile function, using both in vitro and in vivo assays. As a consequence of this muscle weakness, the iDUX4pA-HSA mice develop kyphosis, a common feature of the most severe cases of FSHD, particularly in cases with respiratory issues (32, 33). Furthermore, the decline in function is mechanistically distinct from that seen in muscle disease due to mutation in dystrophin-associated glycoprotein complex proteins, as observed by eccentric contraction-coupled force loss. Rapid loss of force-generating capacity following lengthening contractions is a signature characteristic of dystrophin-deficient (mdx) mouse muscle (34, 35). However, the 20% decline measured here in iDUX4pA-HSA EDL muscle is very modest in comparison to that in mdx EDL muscle, which shows decline of approximately 80% or more in the same assay (36), and marks the muscle pathophysiology of this FSHD mouse model as distinct from that seen in DMD models.

The transcriptional profile of long-term induced muscle shared features with that of short-term induced muscle, including in direct DUX4 target genes and in pathways associated with cell death. However, the long-term induced muscle clearly demonstrated signs of an altered inflammatory component, driving Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways associated with parasitic and bacterial pathogens and autoimmunity. Consistent with the notion of a pathological inflammatory response, we observed by FACS the simultaneous presence of both M1 and M2 macrophages, suggesting that the orderly sequence of inflammatory cell–mediated changes that coordinates regeneration after injury cannot occur in this diseased muscle, with obvious attendant problems for the ability of the muscle to undergo regeneration. The question of why FSHD muscle, which does not experience nearly as much muscle fiber lysis as does DMD muscle, is nevertheless unable to regenerate, has not been addressed to date. The data presented here suggest the possibility that this may be mediated by a temporally disordered macrophage component, highlighting the importance of carefully evaluating the regenerative inflammatory compartment in FSHD.

The most remarkable feature of the gene expression changes evoked in diseased muscle of long-term induced iDUX4pA mice is its similarity to human FSHD at the gene expression level. The human cognate of the iDUX4pA-HSA DEG set is strongly enriched in RNA-seq data from MRI-guided biopsies of FSHD muscle. Although mouse models for FSHD based on DUX4 expression have been questioned due to the apparent distinctions in many DUX4 target genes in mice versus humans, the data presented here demonstrate that this criticism is unfounded. The similar transcriptional profiles arise from common pathological processes, pointing to a high level of relevance of the iDUX4pA model to FSHD at the pathophysiological level.

In addition to changes in the inflammatory compartment, another remarkable feature of dystrophy in this model is a significant increase in the number of FAPs and decrease in the number of endothelial cells in long-term diseased muscle. Both changes occur relatively rapidly after induction of DUX4 expression in fibers and appear to be maintained throughout the course of disease. The reduction in endothelial cells correlated with a patchwork loss of capillary beds, evident by histology and whole-mount imaging of cleared muscle. These capillary changes are relevant to FSHD, as significantly reduced capillary density is a known histological feature of FSHD muscle (26). A common nonmuscle manifestation of FSHD is retinal vascular abnormalities, seen as telangiectasia, capillary closure and leakage, and formation of microaneurysms (37, 38). As the underling mechanism of the capillary pathology in the retina is unknown, a common mechanism might be responsible for vascular abnormalities in muscle and retina. Since in muscle, the loss of endothelial cells is most likely a secondary effect of DUX4 expression in fibers, exploring the direct and indirect effects of DUX4 on retinal vascular endothelial cells may be worthwhile.

Although FAPs have not yet been directly evaluated in FSHD biopsies, the increase in FAPs we see here is also likely to be highly relevant to FSHD. Given the dramatic fibrotic component to iDUX4pA-HSA muscle pathology and occasional intramuscular fat, it seems highly likely that FAPs are central players in pathological evolution, and the fact that this pathological evolution maintains transcriptional similarity to human FSHD at the whole-muscle biopsy level suggests that the pathological processes are similar. Furthermore, transcriptional analysis of the isolated FAPs demonstrated that the human cognates of their DEGs are also significantly enriched in MRI-guided FSHD muscle biopsies. In contrast, the endothelial DEGs did not show significant GSEA enrichment, suggesting again that the FAP changes are the more relevant to human FSHD. This remarkable finding strongly suggests that FAPs in human FSHD may be altered in a similar way to FAPs in dystrophic iDUX4pA-HSA muscle. FAPs will be most interesting to investigate in human muscle specimens

Mice. Four-week-old female mice carrying both the iDUX4pA and HSA-rtTA transgenes were fed dox chow that contains 625, 62 and 6.25 mg/kg dox (ENVIGO). For DUX4 burst induction, a single dose (100 mg/kg) of dox dissolved in PBS was injected intraperitoneally.

Muscle histology. H&E and Sirius red/fast green staining were performed on 10-μm OCT-frozen muscle sections (TA and quadriceps) as previously described (21). For immunofluorescence, tissue sections were fixed in 4% paraformaldehyde for 10 minutes, permeabilized with 0.3% Triton-X for 30 minutes, and stained with primary antibody diluted in 3% BSA overnight at 4°C. After 3 rounds of washing with PBS, the appropriate conjugated antibody was applied for 60 minutes at room temperature. Nuclei were visualized with DAPI (1:5000, MilliporeSigma). The following antibodies were used: rabbit anti-DUX4 (MAB95351, also called RD2-47c, 1:50, R&D Systems) (21), anti-mouse laminin (L8271, MilliporeSigma), anti-rabbit laminin (L9393, MilliporeSigma), anti-embryonic myosin heavy chain (F1.652, Developmental Studies Hybridoma Bank), and CD31 (clone 390, BD Biosciences).

FACS analyses. Bilateral TA, gastrocnemius, soleus, quadriceps, or pectoralis muscles were digested together to obtain mononuclear cells. Muscles were minced with a razor blade and digested with collagenase type II and dispase for 90 minutes. Mononuclear cells were stained with the following antibodies (all rat): PE-Cy7-CD45 (clone 30-F11, BD Biosciences), PE-Cy7-CD31 (clone 390, BioLegend); APC-CD31 (clone 390, eBioscience), PE-Pdgfrα (CD140A, clone APA5, BD Biosciences); APC integrin α7 (67-0010-05, AbLab); PE-CD146 (clone ME-9F1, BD Biosciences); biotinylated-VCAM (CD106, clone 429, BD Biosciences); PE-Sca1 (clone D7, eBioscience); PE-Gr1 (clone RB6-8C5, eBioscience); APC-CD11b (M1/70, eBioscience); Pe-CD206 (clone C068C2, BioLegend); and PE-CD68 (clone FA-11, BioLegend) resuspended in PBS/1% FBS. Samples were run on a BD FACSAria instrument, and data were analyzed using FlowJo (BD Biosciences). All experiments were performed on at least 4 biological replicates.

In vivo analysis of muscular torque. Contractile function of the left anterior crural muscles (TA, EDL, and extensor hallucis muscles) was measured in vivo as previously described (39, 40). Mice were anesthetized with isoflurane (1.5% isoflurane and 125 mL O₂ per minute) and placed on a temperature-controlled platform to maintain core body temperature between 35°C and 37°C. The left knee was clamped and the left foot was secured to an aluminum “shoe” attached to the shaft of an Aurora Scientific 300B servomotor (Aurora Scientific). Sterilized platinum needle electrodes were inserted through the skin for stimulation of the left common peroneal nerve. Stimulation voltage and needle electrode placement were optimized with 5–15 isometric contractions (150-ms train of 0.1-ms pulses at 250 Hz). Following optimization, contractile function of the anterior crural muscles was assessed by measuring isometric torque as a function of stimulation frequency (20–300 Hz), with the highest recorded torque defined as peak isometric torque.

Ex vivo muscle preparation and physiology. Mice were anesthetized with 75 mg/kg pentobarbital and analyzed for EDL muscle contractile function and susceptibility to eccentric contractions. EDL muscles were anchored to a force transducer and incubated in Krebs–Ringer bicarbonate buffer. Baseline contractile force measurements and passive stiffness were conducted as previously described (41, 42). For eccentric contractions, EDL muscles were passively shortened to 90% resting length and then maximally stimulated for 200 ms while the muscle was simultaneously lengthened to 110% resting length at 0.5 lengths/s.

Tissue clearing and whole-mount imaging. Seventy-five μL of Dylight 488–conjugated Tomato Lectin (Vector Labs, DL-1174) was injected through the retro-orbital sinus and allowed to circulate for 10 minutes. Avertin-anesthetized mice were perfused with 60 mL of 4% paraformaldehyde with 0.05% glutaraldehyde (MilliporeSigma, G5882) in PBS. Whole TA muscles were isolated and post-fixed overnight in fresh fixative. The samples were washed in PBS several time and then quenched in 4% (w/v) glycine (MilliporeSigma, G7126) and 4% acetamide (w/v) (MilliporeSigma, A0500) in PBS at 37°C overnight. Tissue clearing was performed as previously described (43). Briefly, samples were equilibrated in cold polymerization buffer (4% acrylamide, MilliporeSigma, A9099) and 0.25% VA-044 (Wako Chemical, NC0632395). Nitrogen-purged polymerization buffer was exchanged and samples were incubated for 3 hours at 37°C. Samples were washed in 200 mM SDS solution (MilliporeSigma, L6026) buffered with 20 mM boric acid buffer (pH 8.4, MilliporeSigma, B6768) and then incubated for at 37°C in a solution containing 5% Quadrol (MilliporeSigma, H2383), 10 mM SDS, an d2% Triton X-100 (MilliporeSigma, T8787) in 20 mM boric acid buffer (pH 8.4). The cleared samples were washed in PBS with 0.5% Triton X-100 and then equilibrated in PROTOS solution (44).

Samples were mounted in a custom imaging chamber and imaged on a Nikon A1R confocal system using the CFI Plan Apochromat 10XC Glyc lens (Nikon). Advance denoising was performed in NIS-Elements-Ar (Nikon), and images were visualized in green pseudocolor (Bitplane).

RNA isolation, RTqPCR, and RNA-seq. RNA was extracted using Trizol (Invitrogen) and an RNA extraction kit (Zymo) following manufacture’s protocol. RNA samples were treated with DNAse on the column (Zymo). cDNA was made using 1.0 μg total RNA with oligo-dT primer and the cDNA Synthesis Kit (Applied Biosystems) following the manufacturer’s instructions. qPCR was performed by using Premix Ex Taq (Probe or SybrGreen qPCR, Takara) and commercially available probes from Applied Biosystems (Gapdh, Mm99999915_g1; vwf, Mm00550376_m1; Pdgfra, Mm00440685_g1) or custom-designed primers (Col1a1 F, 5′ GAGCGGAGAGTACTGGATCG and R, 5′ TACTCGAACGGGAATCCATC; Col3a1 F, 5′ TGGTCCTCAGGGTGTAAAGG and R, 5′ GTCCAGCATCACCTTTTGGT; Tgfb1 F, 5′ CTCCCGTGGCTTCTAGTGC and R, 5′ GCCTTAGTTTGGACAGGATCTG; MMP-2 F, 5′ CAAGTTCCCCGGCGATGTC and R, 5′ TTCTGGTCAAGGTCACCTGTC). Gene expression levels were normalized to that of Gapdh and analyzed with 7500 System Software using the ΔCT method (Applied Biosystems). RNA-seq libraries from iDUX4pA-HSA mice were prepared using 500 ng total RNA from hind limb muscles and the KAPA mRNA Hyper Prep kit following manufacturer’s instructions. For the generation of RNA-seq libraries from endothelial progenitors and FAPs, i.e., CD31⁺ and Pdgfrα⁺ cells, 75,000 cells were FACS sorted from hind limb muscle, and total RNA was extracted using the RNeasy Mini Kit (Qiagen). cDNA was synthesized using the NEBNext Single Cell/Low Input cDNA Synthesis & Amplification kit (New England Biolabs) and libraries we generated with the Nextera XT DNA Library Preparation Kit (Illumina). 50-base paired-end sequencing was performed on an Illumina Novaseq instrument at the University of Minnesota Genomics Center.

Bioinformatics analysis. The raw FASTQ files were analyzed using a customized pipeline (gopher-pipelines; https://bitbucket.org/jgarbe/gopher-pipelines/commits/8c73e0e82e51caf7de4c2ee58a949cab67e364b9) developed and maintained by the Minnesota Supercomputing Institute. Specifically, the sequencing quality was checked using FastQC v0.11.7 (http://www.bioinformatics.babraham.ac.uk/projects/fastqc/). Low-quality reads were then filtered using Trimmomatic v0.33 (http://www.usadellab.org/cms/index.php?page=trimmomatic) (45). The remaining reads were then aligned to GRCh38/hg38 or mm10 reference genome, and transcript abundances were quantified using Salmon v0.6.0 alignment algorithm (46). The reads were filtered to include genes with counts per million value over 5 in at least 20% of the samples for subsequent analysis. DEGs analysis was performed in R v3.4.3 using edgeR package (https://bioconductor.org/packages/release/bioc/html/edgeR.html) (47). The DEGs were visualized using the ggplot2 package. Data have been deposited into the GEO database under accession number GSE141562.

PCA was performed using prcomp and visualized using ggplot2 in R. R was also used for generating Venn diagrams.

The DUX4 scores (Figure 3A, Figure 4G, and Figure 6G) were calculated using a previously validated method (48). Briefly, the score was calculated using t statistics by comparing the top 50 upregulated and bottom 50 downregulated DUX4 target genes identified previously in iC2C12-DUX4 cells after 4 hours of induction (13).

KEGG pathway analysis was performed in R using the enrich KEGG function from the clusterProfiler package (49). A FDR-adjusted P value of less than 0.05 was used as the cutoff for statistical significance.

Previously published gene expression analyses were used to compare with our data. Raw fastq files of 43 samples (9 control and 34 MRI-guided FSHD) were downloaded from the GEO database (GSE115650) (25). Raw fastq files of 12 samples were downloaded from the SRA database (SRP148291) (31). All of the fastq files were analyzed using the same method described above. The log-normalized microarray data of 15 samples (7 control and 8 MRI-guided FSHD) were also downloaded from GEO database (GSE26852) (24). GSEA were performed using default parameters, except that the ratio of classes was used as the metric for ranking genes.

Statistics. GraphPad Prism software was used for statistical analyses of the data, except where indicated. Differences between groups were evaluated by 1-way or 2-way ANOVA followed by Tukey’s post hoc tests. Differences were considered significant at P values of 0.05 or lower.

Study approval. Mice were maintained, and in vivo experiments were conducted, under a protocol approved by the University of Minnesota IACUC.