|
|
||||||||
| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||

,1














,1
* Institute of Neurology, Catholic University, Rome, Italy;
Neuromuscular Disease Unit, Department of Pediatrics, G. Gaslini Institute, University of Genova, Genova, Italy;
Unit of Molecular Medicine and Pathology, Bambino Gesù Childrens Hospital IRCCS, Rome, Italy;
Institute of Anatomy and Cell Biology, Catholic University, Rome, Italy;
|| Division of Cardiology and ICU, St. Paolo Hospital, Civitavecchia, Italy;

Center for Neuromuscular Diseases, UILDM, Rome, Italy; and

Don Gnocchi Foundation, ONLUS, Rome, Italy
1Correspondence: M.P. and E.R.: Institute of Neurology, Catholic University, L.go A. Gemelli 8, 0018, Rome, Italy. E-mail: m.pescatori{at}rm.unicatt.it; ericci{at}rm.unicatt.it
| ABSTRACT |
|---|
|
|
|---|
Key Words: microarrays Duchenne muscular dystrophy human disease Wnt, Notch, and BMP signaling
| INTRODUCTION |
|---|
|
|
|---|
In mature muscle, dystrophin is localized adjacent to the cytoplasmic face of the sarcolemmal membrane in a cytoskeletal protein assembly, termed costameres, that links the force generating sarcomeric apparatus to the sarcolemmal membrane and the extracellular matrix (ECM). In this location, dystrophin is associated with a large oligomeric complexthe dystrophin-glycoprotein complex (DGC) (2)
which bridges the costameric cytoskeleton and the ECM. In the absence of dystrophin, greater stress is placed on myofibrillar and membrane proteins on muscle contraction. This produces severe muscle damage and generates a restless cycle of necrotic/regenerative dynamics. The progressive inability of DMD skeletal muscle to properly complete tissue regeneration leaves a fatty fibrous "scar," ultimately determining the muscle contractile dysfunction typical of the disease.
Albeit an increased serum creatine kinase level and abnormal muscle histology are always present, boys with DMD are phenotypically indistinguishable from the normal ones at birth and, in their first years of life, acquire early motor milestones at normal times. A clear defect in muscle function becomes generally apparent by the end of the second year. As the disease is typically diagnosed between the ages of 3 and 7, the first 2 years are often considered and referred to as clinically presymptomatic.
Large-scale parallel gene expression analysis of skeletal muscle biopsies from DMD patients and unaffected controls has been successfully used in some recent studies to describe at a transcriptional level shared muscle tissue alterations in DMD children older than 5 years (3
4
5
6
7)
, when the disease first shows its major consequences on muscle function and ambulation. As a defined gene expression signature was shown to characterize these symptomatic patients, we sought to investigate whether and to which extent alterations may also be present in muscle from DMD infants. To this aim, we used Affymetrix technology to compare the individual expression profiles of 19 DMD patients with age at biopsy scattered along the first 2 years of the disease with those of 14 age-matched controls. This approach allowed us to describe with high resolution the altered transcriptional state that characterizes this early, presymptomatic phase of the disease and highlight some molecular pathways as potential critical targets in the pathophysiology of the disease.
| MATERIALS AND METHODS |
|---|
|
|
|---|
All bioptic specimens used in this study were taken, for diagnostic purposes, under institutionally approved protocols. Since all our patients were minors, their parents were asked to sign an informed consent disclosing future use of the bioptic materials for research. This study was reviewed and authorized by our institutional Ethics Committee, according to our institutional regulation and national laws and guidelines.
RNA extraction
Total RNA was extracted from frozen quadriceps muscle biopsies by TriZol (TriZol reagent, Invitrogen, Carlsbad, CA, USA). RNA was further purified using the RNAeasy mini kit following the RNA cleanup protocol as indicated by the manufacturer (Qiagen, Valencia, CA, USA). RNA purity and integrity were assessed by spectophotometric analysis and agarose gel electrophoresis.
Affymetrix GeneChips
We use one type of gene chip in this study: the Affymetrix HG-U133A GeneChip. CRNA synthesis was performed using 5 µg of total RNA as template, as described in the Affymetrix Gene Expression Manual. Gene chips were washed and stained in an Affymetrix fluidic station 430 and analyzed on an Affymetrix G2500 GeneChip scanner. To prevent overcorrelation, samples were processed eight at a time and arranged so that DMD and control samples were both in each experimental session.
Data analysis
Data analysis was performed using BRB ArrayTools developed by Dr. Richard Simon and Amy Peng Lam (8
9
10)
.
Probe level summaries were generated using an empirically motivated statistical approach: the log scale robust multiarray analysis (rma) procedure developed by Irizzarry et al. (11)
. An expression table reporting normalized, background corrected, probe level summaries for all the 22,283 probesets is posted under supplemental material (DMDDATA_rma). The data set was filtered to exclude from the analysis genes showing minimal variation across the set of arrays. Probesets whose expression differed by at least 1.5-fold from the median in at least 20% of the arrays were retained.
All the experiments in this work met specific GeneChip QC criteria (supplemental material, GeneChip QC). A data table (rma), together with the relative cel files and relevant information about the experiment, is available at http://www.ncbi.nlm.nih.gov/geo/ under accession #GSE6011. Whenever needed, association probability (P value) and correlation coefficient (R2) among genes or samples were computed using Excel.
Multidimensional scaling (MDS)
Multidimensional scaling is a group of methods for representing high-dimensional data graphically in low (usually 2 or 3) dimensions. The objective in multidimensional scaling is to preserve pairwise similarities or distances between objects in the low-dimensional graphical representation. Multidimensional scaling analysis is similar to cluster analysis in that one is attempting to examine relationships among samples, but it also provides a graphical representation of the pairwise similarities or distances among samples without forcing the samples into specific clusters. We measured the similarity between gene expression patterns by computing the correlation distances (1-Pearson correlation) for each pair of samples based on standardized log-transformed expression values across all of the genes passing the general filtering (see Data analysis). During the multidimensional scaling procedure (BRB ArrayTools, Multidimensional scaling), samples were positioned in a 3-dimensional space so that the distance between each pair of samples very closely approximated the correlation distance measurements in the matrix for the corresponding sample pair (8
9
10)
. Samples with gene expression profiles more similar to each other will lie closer and form an aggregation (cluster) in 3-dimensional space. We also ran the same analysis using Euclidean distances, and the result was similar, showing separate clustering of DMD and control samples
Class comparison
We computed the probability of genes being differentially expressed between the two classes, using the unequal variance t test, and identified genes differentially expressed between the two classes by a multivariate permutation test (8
9
10)
to provide 95% confidence that the false discovery rate was <10% (BRB ArrayTools, Class comparison). Although t statistics were used, the multivariate permutation test is nonparametric and does not require the assumption of Gaussian distributions. The results of this analysis are displayed as an html file reporting all the differentially expressed genes ranked by P value and tabulated along with descriptive statistics and links to the Entrez/gene database (supplemental material, DMDCLASSCOMP_filt).
Gene ontology
The evaluation of which gene ontology classes were differentially expressed between normal and DMD samples was performed using a functional class scoring analysis as described by Pavlidis et al. (12)
. For each gene in a GO class, the P value for comparing normal vs. DMD samples was computed. The set of P values for a class was summarized by two summary statistics: 1) the LS summary is the average log P values for the genes in that class and 2) the KS summary is the Kolmogorov-Smirnov statistic computed on the P values for the genes in that class. The statistical significance of the GO class containing n genes represented on the array was evaluated by computing the empirical distribution of these summary statistics in random samples of n genes. We considered a GO category differentially regulated if the significance level of either one of the KS or LS statistics was <0.005 (BRB ArrayTools, Class comparison). The results of this analysis are displayed as an html file (Supplemental Material, DMDGO).
Correlation analysis
To search for genes whose expression was significantly related to patient age, we computed the significance level for each gene to test the hypothesis that the Spearmans correlation between gene expression and age was zero. These P values were then used in a multivariate permutation test (8
9
10)
in which ages were randomly permuted among arrays to provide confidence that the median value of false discoveries was <10 (BRB ArrayTools, Quantitative Trait Analysis). The multivariate permutation test is nonparametric and does not require the assumption of Gaussian distributions.
Patients were arbitrarily grouped into 4 age defined classes:
Class 1: 1.55 months, n = 5
Class 2: 612 months, n = 8
Class 3: 1422 months, n = 6
Class 4: 2861 months, n = 3
The results of the analysis are displayed as an html file reporting correlated genes ranked by P value and tabulated along with descriptive statistics and links to the Entrez/gene database (supplemental material, DMDAGECORR). We ran a similar correlation analysis across all the control profiles. None of the genes showing age dependence in patients showed a similar behavior in controls. To exclude the possibility that genes presented as modulated with age in patients and those discussed in this article were regulated similarly in controls, we browsed the data set and analyzed the trends of expression in controls for the 16 genes reported in Table 3
. This was done by plotting the absolute expression values of the control subjects, ordered by increasing age, and visually inspecting the expression trends, highlighted by the resulting curve. Using this procedure, we identified three genes (MAPRE3, C6orf106, and POSTN) whose expression showed age dependence in both patients and controls.
|
We are aware that the particular grouping chosen is arbitrary and that changing the number, size, and composition of the classes would affect results of the analysis. The number of genes reported is small and likely represents a gross underestimation of the real number of those modulated along the progression of the DMD. This approach therefore was not intended to provide a comprehensive description of the progression of the disease; it was a tool that allowed us to identify some genes whose expression trends, in our patients, might suggest a role in progressive pathophysiological processes.
Pathogenic components analysis
We defined four clusters of genes as representative of four major aspects of DMD pathophysiology (muscle regeneration/immaturity, inflammation, ECM homeostasis, and energy metabolism). From the class comparison results, we sorted genes involved in the relevant pathophysiological processes into four lists based on information obtained from both public databases (NCBI Entrez Gene, Jackson Laboratory Mouse Genome Informatics and Weitzman Institute of Science Gene Cardsl) and upon extensive search of the scientific literature in PubMed. The four lists were made nonredundant (one probeset for gene), and gene clustering (BRB ArrayTools, 2-way clustering algorithm, average linkage, "one minus correlation" distance, median-centered expression values) was used to extract from each list a correlated component showing covariation of the expression values across the data set. For each gene, log2 expression values were expressed as the ratio to the mean of the control population. The resulting fold changes were averaged in each cluster to generate four indexes (average log2 fold change in cluster) we consider indicative of the extent of activation of these four aspects of the disease.
Real-time rt-PCR
Two-step, real-time polymerase chain reaction (PCR) reactions were performed on skeletal muscle RNA from 15 DMD patients (DMD1, DMD2, DMD3, DMD4, DMD5, DMD7, DMD9, DMD11, DMD12, DMD13, DMD15, DMD18, DMD19, DMD21, and DMD22) and 8 unaffected controls (C1, C3, C4, C7, C9, C11, C12, and C13). First-strand cDNA was synthesized from 1 µg total RNA using Superscript II (Life Technologies, Carlsbad, CA, USA) and an oligo dT primer, following the manufacturers instructions. Real-time PCR was performed on an ABI Prism 7000 Sequence Detection System using the Applied Biosystem TaqMan universal PCR master mix (with UNG) and the following gene-specific TaqMan primers sets (Applied Biosystems, Foster City, CA, USA): GAPDH, Hs99999905_m1; MCSP/NG2, Hs00426981_m1; MYH3, Hs00159463_m1; MYH8, Hs00267293_m1; ACTC, Hs00606316_m1; CHRNA1, Hs00175578_m1; COL3A1, Hs00164103_m1; COL1A1, Hs00164004_m1; COL1A2, Hs00164099_m1; CALP6, Hs00560073_m1; CHRNG, Hs00183228_m1; CHRNE, Hs00181084_m1; myostatin, Hs00193363_m1; follistatin, Hs00246260_m1; FSTL1, Hs00200053_m1. Standard curves were generated to verify PCR efficiency. After amplification, the difference between threshold cycles
ct (ct GeneX-ct GAPDH) was calculated for all samples and the mean
ct value in the control population was subtracted from all samples. Log2 expression differences were converted to fold change ratios using the equation FC = 2
ct. The unequal variance t test was used to compute the probability of genes being differentially expressed between the two classes (Excel).
| RESULTS |
|---|
|
|
|---|
|
DMD vs. control class comparison: genes and gene ontology (GO) categories modulated in DMD skeletal muscle
Of the 22283 probe sets represented on the HG U133A GeneChip, 1663 (7.5%) met inclusion criteria and were retained for further analysis. As a result of the t test, 127 probe sets were differentially expressed between the two classes with a P value lower than 1E-07, 202 with P < 1E-06, 281 with P < 1E-05, 399 with P < 1E-04, and 561 with P < 1E-03. To control the proportion of false discoveries within differentially expressed genes, we used a multivariate permutation test to provide the 95% confidence that the false discovery rate was <10% (8
9
10)
. By this cutoff statistic, we identified 777 probe sets as differentially expressed; this was equivalent to applying a probability threshold of 9.7E-03. A table reporting a list of these probe sets, ranked by P value and tabulated along with within-class geometric mean, between-classes fold difference, and annotations, is posted under supplemental material (DMDCLASSCOMP_filt). Since many genes are represented on the U133A GeneChip by multiple probes, these 777 probe sets are representative of 618 genes. Among them, more are overexpressed (421/618) than underexpressed (197/618). Table 1
reports a selection of differentially expressed genes grouped by function. We should mention that the fold change estimates generated by our analysis are smaller than those reported earlier by other authors (3
4
5
6
7)
. The difference can be attributed to the use of a different probe set summarization algorithm, rma, which is more precise than MAS5.0 but produces lower fold change estimates (11)
. The possibility that this discrepancy is due to the different ages of the patients studied can be ruled out, as we included in our study two older patients (aged 4 and 5 years) who showed expression values and fold changes comparable to presymptomatic patients.
|
Fourteen genes showed in all patients expression values greater than the greatest value observed in controls; conversely, no genes had in all patients expression values lower than the lowest in controls (Table 2
). The two genes more frequently underexpressed, outside the control range in 21/22 patients, were the muscle-type glycogen phosphorylase (PYGM) and CAPRI, a Ca2+-dependent, RAS GTPase-activating protein (AP) (RASA4) that interacts with members of the RHO and CDC42 family of small G-proteins. Dystrophin mRNA was reduced in 20 of 22 patients.
|
Complementary to the gene-based approach, we provided additional information by evaluating which gene ontology classes were differentially expressed between normal and DMD muscle (10
, 12)
. A table reporting results of the analysis is posted under supplemental material (DMDGO).
Transcriptome alteration in DMD skeletal muscle
Muscle genes
The group of genes showing the greatest changes in expression encodes myofibrillar components induced along with the activation of a bona fide muscle regeneration program. This includes postnatal re-expression of developmental isoforms of sarcomeric proteins and increased transcription of genes encoding postsynaptic NMJ components and cytoskeletal elements supporting chronic remodeling of sarcomeric structures (13)
. By real-time PCR, we also detected induction of the nAChR gamma subunit (CHRNG), whose expression is shut down at birth and replaced by the postnatal-type epsilon subunit (CHRNE) in normal individuals (14
, 15)
(supplemental material, Table ST2). More than 40 transcripts encoding enzymes involved in glycogen metabolism, glycolysis, TCA cycle, lipid transport, and ß-oxidation were down-regulated. Although the magnitude of the change was not large (0.40.6x), the downward modulation of this gene cluster represents a highly coordinated transcriptional change and an invariant character of this disease (3
, 4)
(Table 1
, Muscle genes).
Dystrophin (Dp427) was decreased in 20/22 patients. Two subjects (DMD1 and DMD11) appeared to escape nonsense-mediated decay (NMD) and displayed close to normal dystrophin mRNA expression. However in both cases the protein was absent from the muscle fibers and the two children displayed typical clinical features. Moreover, they did not appear as a correlated group by cluster analysis (data not shown). A small proportion of DMD patients showing reduced dystrophin NMD have been described by Chelly et al. (16)
.
The relative abundance of most of the transcripts encoding other DGC components was normal except for syntrophin
1 (SNTA1), which was repressed, and sarcoglycan
(SGCE), which was induced (supplemental material, DGC Components Expression).
Postnatal myogenesis relies on the activation of resident myogenic cells that proliferate throughout the regenerating muscle and differentiate to fuse into multinucleated muscle fibers (17)
. As cells are not synchronized and can undergo asymmetric mitosis, different phases of the myogenic program are expected to coexist in regenerating muscle.
In our patients, we observed the concomitant induction of genes involved in cell cycle progression (cyclin G2, cyclin D2, and CDK4) and cell cycle withdrawal (CDKN1A, GAS1, and GADD45A). Genes encoding myogenic transcription factors of the bHLH and MEF family (MyoD, myogenin, and MEF2C) were also induced, indicating diffuse activation of myogenesis in DMD muscle (18)
.
In recent years, increasing attention has been placed on the role of the myostatin (19
20
21)
and insulin-like growth factor 1 (IGF1) (22
23
24
25)
pathways in controlling muscle trophism and myogenesis (26)
. In our patients, we observed increased expression of IGF1, IGF2, and IGF binding proteins 4 and 7. PRSS11, an IGFBP protease acting as negative modulator of this pathway (27)
, was also induced. We studied the expression of myostatin (GDF8), follistatin (FST), and follistatin-like 1 (FSTL1) by real-time rt-PCR. Both FST and FSTL1 transcripts were increased in patients whereas expression of GDF8 was often reduced, but did not reach statistical significance (Fig. 2
and supplemental materials, Table ST2).
|
Three major proteolytic systems operating in muscle have been implicated in the pathophysiology of muscle atrophy and dystrophy: lysosomal proteases (CTS) (28)
, calpains (CAPN) (29
, 30)
, and the ubiquitin proteasome system (31
, 32)
. Three cathepsin genes were up-regulated (cathepsins B, C, and K), CAPN3 was repressed, and CAPN6 was induced. Genes encoding proteasome subunits were variably modulated: PSMB8 and PSME2 genes were induced and the PSMD12 repressed. Enzymes involved in protein ubiquitination (UBE2B, UBE2D1, CUL5, and FBXO3) or deubiquitination and ubiquitin recycling (USP13 and USP25) were repressed.
Inflammation
Genes involved in immune response and inflammation were massively induced in our patients, and our analysis provided a detailed description of the inflammatory response in this presymptomatic phase of the disease (Table 1
, Inflammation).
We observed increased expression of class I and II MHC, components of the complement system (C1, C3, and the H factor), a set of IFN-inducible genes, and markers of infiltrating immune cells. Four chemokines (CCL14, CCL2, CXCL12, and CXCL14) were also induced. Among these, CXCL14 has been shown recently to be a potent chemoattractant and activator of dendritic cells (DC) and is suggested to be involved in DC homing (33)
. A robust inflammatory response was recently reported in four patients aged 810 months and ascribed to early activation of the NF-
B pathway by TLR7 (34)
. Accordingly, both TLR7 and class I MHC and VIM reported by the authors as NF-
B target genes were induced in our patients since early postnatal life. In line with the idea that exaggerated NF-
B activation may exert a detrimental effect over the course of the disease, pharmacological inhibition of this pathway was reported to ameliorate the mdx pathology (35
, 36)
.
Fibrosis
The steady-state amount and composition of the ECM that surrounds cells in solid tissues rely on a balance between deposition of structural components and remodeling (37)
.
Our patients showed increased expression of a large number of genes encoding ECM components and enzymes involved in matrix biosynthesis and remodeling (Table 1
, ECM remodeling). Type I and III collagen, lumican, and asporin showed the largest fold changes whereas other collagen types (type IV, V, VI, XIV, XV, and XVIII), fibronectin, laminins, elastin, and proteoglycans (chondroitin sulfate, heparan sulfate, and small leucine-rich type) were induced to a lesser extent. In this framework, we found increased expression of the matrix metalloproteinase MMP2 and two MMPs inhibitors: TIMP1 and TIMP2. Von Moers et al. (38)
, by describing a similar expression pattern in DMD patients aged 3.515 years suggested that an unbalance in MMP1/TIMP1 stoichiometry characterizes DMD muscle and contributes to progressive fibrosis. Since these alterations were already present in our 1.5-month-old patient, increased ECM synthesis and inhibition of fibrolytic activity are early events in DMD. TIMP1 expression can be induced in response to cytokines and hormones and has been linked to the development of pulmonary and liver fibrosis (39
, 40)
. Among such cytokines, the TGF-ß family of multifunctional cytokines controls proliferation, differentiation, and other cellular activities by acting as negative autocrine growth factors (41)
. Three members of the TGF-ß family (namely, TGFB1, TGFB2, and TGFB3) are known to participate in tissue regeneration. Deregulation in their signaling is implicated in the development of organ fibrosis and scar formation. The increased expression of these cytokines has been described in DMD patients and in mdx mice (38
, 42
43
44
45
46)
.
In our young patients, TGFB1 was the predominant TGF-ß isoform induced, being expressed above all controls in all patients but one (Table 2
, Fig. 3
). In agreement with earlier studies reporting its overexpression in advanced DMD (3
, 4
, 34)
, TGFB3 was expressed above control mean in 18 of 22 patients. However, TGFB3 was induced to a lesser extent than TGFB1, and only in patients older than 6 months (Fig. 3)
.
|
Decorin (DCN) and dermatopontin are known modulators of TGF-ß activity and were both induced. These proteins can individually elicit opposite regulatory effects on TGF-ß signaling but, when coexpressed, form a stable complex, with inhibitory activity (47)
. A recent report described the reduced expression of DCN mRNA in DMD patients aged 28 years (48)
. In our data set, four independent probesets showed superimposable results (0.93<R2<0.95, among three probesets), indicating increased DCN expression.
DMD progression does not represent a major factor affecting muscle tissue expression profile
An original goal of this study was to identify gene expression alterations in the skeletal muscle of presymptomatic DMD children, which may represent molecular landmarks for progression of the disease. In our cohort of infant patients (19/22 aged 1.5 to 22 months), clustering analysis did not display patient groupings by age, showing that age is not a major variable affecting gene expression in DMD muscle.
Since function-related gene clusters rather than single genes may represent and describe the progression of the DMD more comprehensively, we generated four lists of transcripts representative of four key aspects of disease pathophysiology: muscle regeneration, inflammation, ECM remodeling, and energy metabolism (supplemental materials, Genelist_Pathogenic Components). Following the expression trends of these clusters in the natural history of the disease, we show that, from a transcriptional point of view, three of these components are already fully activated in the younger patients studied. The muscle regeneration cluster only is progressively induced during the first months of postnatal life (Fig. 4
).
|
We also compared our results with those reported in earlier studies based on the Affymetrix technology.
Besides the overall good agreement on transcripts and pathways modulated in DMD that has already been mentioned, all genes that were considered relevant and singularly discussed by other authors (3
4
5
6
7)
were confirmed in our patients. Remarkably, we were able to compare at the single gene level our results with those reported by Haslet et al. (4)
, whose study design was similar to ours (supplemental material, Table ST3): of the 105 genes scored by the authors as differentially expressed (P<1E-4 and mean fold change>2), 87 (83%) turned out to be differentially expressed within the same probability threshold in our data set, 95 (90%) within P < 1E-3, and 104 (99%) within P < 5E-2 whereas only one showed a probability higher than 0.05 (FADS3, 0.0538). The comparison showed that all of the genes reported by the authors as modulated in advanced DMD with high confidence were similarly modulated in our presymptomatic patients. This result also suggests that cross-experiment analysis of the distribution of P values within differentially expressed gene lists measures and describes the level of congruence between independent experiments better than the widely used degree of overlap between the same lists.
Correlation analysis identifies genes induced or repressed with patient age
Longitudinal analysis of time-series microarray data has been used to describe dynamic expression changes occurring in the progression of the mdx pathology as well as the response to specific acute events like experimentally induced muscle damage or denervation. As unsupervised clustering would not highlight the effect of the covariate "age" if it affected the expression of a small number of genes, we analyzed our patients as a time series searching for genes that are induced or repressed along the natural history of the disease. We made use of correlation analysis to search for genes modulated across 4 age-defined classes of patients: 1.56, 712, 1324, and 2560 months. Using this approach, we identified a small number of genes whose expression is either increased or decreased in patients but not in controls as a function of age (supplemental material, DMDAGECORR).
Table 3
reports 16 genes correlated to patients ages with P values lower than 5E-04. Among these, GPSM2, encoding a cytoplasmic regulator of G-protein signaling involved in mitotic spindle movements and cell cycle progression, showed the highest P value and correlation coefficient (Fig. 5
A). The increased expression of this gene with patient age was further confirmed by a second probe set interrogating the same transcript. Other genes similarly induced include FRZB1, also independently scored by two probe sets, and OXCT1, a key metabolic enzyme that catalyzes a rate-limiting step in the metabolism of ketone bodies (Fig. 5B
, C). The higher expression of this gene in older patients is consistent with the elevation of plasma ketone bodies reported in advanced DMD (49
, 50)
. FRZB1 is a member of the soluble frizzled related proteins family (SFRPs) of Wnt inhibitors. Activation of Wnt signaling is required to sustain proliferation of vessel-derived mesenchimal progenitors (mesoangioblasts) and is also able to induce myogenic specification of muscle-derived SCA1+ SP cells and circulating AC133+ cells (51
52
53
54)
. FRZB1 and other sFRP members have been shown to antagonize Wnt signaling in these cell types. SFRP4, another member of the inhibitory SFRP family, was similarly induced with patient age.
|
Among genes progressively down-regulated, COL4A1, POSTN, and DLK1 (Fig. 5D-F
) showed the highest P values followed by ADAMTSL3 and GREM2. It should be noted that although the expression of these genes decreases with patient age, they are not underexpressed in patients. These three genes are highly induced in early postnatal DMD muscle but, after a decreasing trend, their expression is normalized by the end of the second year.
DLK1, also known as preadipocyte factor 1 (PREF1), is a human homologue of the Drosophila Notch ligand delta. The product of this gene represses adipogenic differentiation in a variety of experimental models (55
, 56)
and can also affect muscle trophism. In callipyge sheep, the increased expression of DLK1 has been causally associated with the muscular phenotype and hypertrophy of type II fibers that characterize these animals (57)
.
GREM2 is a member of the Cerberus and Dan family of BMP antagonists that play a role in regulating organogenesis, body patterning, and tissue differentiation (58
, 59)
. BMP signaling promotes osteoblastic specification in different types of potentially myogenic cells like the C2C12 cell line, mesoangioblasts, and other populations of muscle-derived stem cells (54
, 60
, 61)
. Taken together, the results of our correlation analysis show that expression of individual genes changes as a function of patient age.
Real-time rt-PCR
We validated by real-time PCR some of the gene expression changes identified by microarray analysis as well as a number of genes we considered relevant to the disease (CHRNE, CHRNG, GDF8, FST, and FSTL1). All of the genes further analyzed showed expression changes consistent with those estimated by gene chip analysis (supplemental material, Table ST2). In general, PCR-based expression change estimates were larger than those observed in the microarray study.
| DISCUSSION |
|---|
|
|
|---|
However, as pathological findings can be observed in muscle sections of presymptomatic patients, it is not surprising that our infant patients shared a dystrophic signature with the older ones. Nevertheless, many of the observations reported in our study had not been described in human patients, and the information provided by our analysis complements and integrates the great amount of transcriptional data available on more advanced DMD muscle.
Notably, the timing of induction of the molecular pathology in DMD appears different from the one described in mdx mice, in which minimal transcriptional alterations were shown to characterize the initial prenecrotic phase of the disease, whereas later phases were contrasted by parallel changes in transcriptome composition (62
63
64)
.
Our attention was drawn to three elements of this complex picture because of their potential relevance to the clinical progression of the disease and to the development of therapeutic strategies: 1) extensive induction of inflammation; 2) the presence of a fibrogenic signature; 3) progressive establishment of an unfavorable pattern of morphogenetic signaling.
In a recent study addressing the issue of the evolution of muscle expression profile in the natural history of the DMD, Chen et al. showed that a chronic inflammatory response characterizes presymptomatic DMD muscle (34)
. Although the general traits of the transcriptional alterations in this preclinical phase of the disease are poorly described, many of the observations reported by the authors hold true in our data set. As highlighted by the large number of genes involved in immune function induced in our young patients, our study provides a robust statistical foundation for their findings and extends the observation to patients as young as 1.5 months. It is worth noting that the inflammatory component identified in our patients is larger than previously estimated, questioning the idea that inflammation is induced to a lesser extent in DMD than in mdx pathology. Accordingly,
67% of the inflammatory genes reported as induced in mdx mice (64)
displayed similar behavior in our patients.
We show that other aspects of the molecular pathology of the disease are also induced early in the natural history of DMD. Among these, a fibrogenic signature consisting of a large number of genes whose products participate to the protein makeup of fibrotic formations or function in ECM remodeling is prominent in our patients, even though only negligible signs of fibrosis can be observed at this stage of the disease. This signature does not appear sufficient per se to induce progressive ECM accumulation and fibrosis. A similar pattern of gene expression can be observed in immuno-mediated muscle inflammatory diseases (refs. 65
, 66
and our unpublished data) where, upon timely pharmacological intervention, muscle regeneration can occur efficiently, and abnormal matrix deposition and scar formation can be prevented. This is also consistent with the observation that, in mdx mice, a common fibrogenic signature is shared by the diaphragm and the quadriceps muscles, although only the former undergoes progressive fibrosis (46
, 62
, 63)
.
Infiltration of muscle by phagocytic inflammatory cells and induction of ECM synthesis and remodeling are commonly observed during muscle regeneration (67
68
69)
, and a large part of the transcriptional alterations described may therefore reflect the permanent regenerative state to which DMD muscle is set. Muscle regeneration does not necessarily result in abnormal matrix deposition and scar formation, but in several instances represents an effective means to ensure tissue plasticity and replace damaged fibers. Concurrent mechanisms, involving a persistent imbalance between ECM synthesis and degradation and a decrease in muscle regenerative capacity, contribute to determining the inability of DMD muscle to complete an effective regenerative process. Although progressive depletion of the satellite cells compartment is still a leading hypothesis to explain the inability of DMD muscle to support lifelong the high turnover rate demanded by recurrent myofiber injury, current evidence suggests that tissue microenvironmental factors may contribute to the progressive loss of regenerative capacity that characterizes the advanced stages of the disease (70
, 71)
. In a recent work, Conboy et al. showed that age-dependent attenuation of Notch signaling is the cause of the decreased ability of aged muscle to regenerate (72
, 73)
. Exogenous stimulation of this pathway was able to rejuvenate satellite cells and restore regenerative potential to old muscles (72
, 73)
, supporting the notion that muscle regeneration may be pharmacologically controlled. This may also suggest that abnormalities in muscle regeneration may arise from postnatal defects in morphogenetic signaling.
Four genes identified by our correlation analysis encode extracellular or membrane-bound proteins functioning in ligand/receptor interactions upstream of one of the Wnt, Notch, and BMP pathways. Changes in the activation state of each of these pathways may influence cell cycle progression and fate decisions in some potentially myogenic, mesenchimal stem cell populations (51
52
53
54
55
56
57
58
59
60
61)
. We may speculate that the altered expression of these signaling molecules offers evidence of the establishment of a nonpermissive tissue microenvironment and that abnormal morphogenetic signaling may interfere with the regenerative capacity of DMD muscle by preventing proliferation of myogenic progenitors and/or their commitment to a myogenic fate.
We consider the correlation across patient age a significant aspect of our analysis, as it provides evidence that the expression of genes changes along the natural history of DMD. It is conceivable that applying a similar approach to a large number of patients whose ages scatter over a range representative of the full evolution of the disease may highlight aspects of DMD pathophysiology that, because of their evolutive nature, have so far escaped the common patients vs. controls experimental design.
In conclusion, the demonstration that most of the molecular aspects contributing to the pathophysiology of DMD are already induced in infant patients, in agreement with previous evidence supported by histological and immunohistochemical findings, supports the idea of precocious therapeutic interventions for the disease (34)
. Furthermore, as microarray analysis has been proposed as a tool to evaluate the effects of therapy in DMD and mdx mice (74
75
76)
, our study, by providing the basic reference knowledge required, laid the ground for application of this technology to evaluate the effects of therapy in preclinical DMD children.
SUPPLEMENTAL MATERIAL
Supplemental Tables ST1, ST2, and ST3 and files MIAME Compliance, DMDDATA_rma, GeneChip QC, DMDCLASSCOMP, DMDGO, DMDAGECORR, Genelist_Pathogenic Components, DGC_Components Expression, and Supplemental Information_Materials and Methods are available at www.fasebj.org. The 37 microarray data set is available at http://www.ncbi.nlm.nih.gov/geo/ and under accession #GSE6011.
| ACKNOWLEDGMENTS |
|---|
Received for publication September 15, 2006. Accepted for publication November 21, 2006.
| REFERENCES |
|---|
|
|
|---|