Genetic Modifiers and Phenotype of Duchenne Muscular Dystrophy: A Systematic Review and Meta-Analysis

The transforming growth factor beta (TGFβ) pathway could modulate the Duchenne muscular dystrophy (DMD) phenotype. This meta-analysis aims to estimate the association of genetic variants involved in the TGFβ pathway, including the latent transforming growth factor beta binding protein 4 (LTBP4) and secreted phosphoprotein 1 (SPP1) genes, among others, with age of loss of ambulation (LoA) and cardiac function in patients with DMD. Meta-analyses were conducted for the hazard ratio (HR) of LoA for each genetic variant. A subgroup analysis was performed in patients treated exclusively with glucocorticoids. Eight studies were included in the systematic review and four in the meta-analyses. The systematic review suggests a protective effect of LTBP4 haplotype IAAM (recessive model) for LoA. It is also suggested that the SPP1 rs28357094 genotype G (dominant model) is associated with early LoA in glucocorticoids-treated patients. The meta-analysis of the LTBP4 haplotype IAAM showed a protective association with LoA, with an HR = 0.78 (95% CI: 0.67–0.90). No association with LoA was observed for the SPP1 rs28357094. The LTBP4 haplotype IAAM is associated with a later LoA, especially in the Caucasian population, while the SPP1 rs28357094 genotype G could be associated with a poor response to glucocorticoids. Future research is suggested for SPP1 rs11730582, LTBP4 rs710160, and THBS1 rs2725797.


Introduction
Duchenne muscular dystrophy (DMD) is an X-linked recessive lethal disease that predominantly affects males, with an incidence of one case per 3500-9000 live births [1]. It is caused by mutations, insertions, or deletions of the dystrophin/DMD gene (location Xp21.2-p21.1) resulting in an absence of the cytoskeletal protein dystrophin. Dystrophin is an actinbinding protein that, through the α/β-dystroglycan complex, links the cytoskeleton to the extracellular matrix protein laminin, stabilizing the surface of muscle cells and preventing their rupture during contraction and relaxation cycles. Therefore, this is essential for the strength, stability, and functionality of the myofibers [2,3]. The absence of dystrophin in the muscle tissue causes a destabilization of the dystrophin-associated glycoprotein complex. This leads to, during contraction and relaxation cycles, sarcolemmal instability, and a decrease in force transmission by the sarcomere occurs. In dystrophic muscle there is an infiltration of mononuclear cells, variation in fiber size, centrally located nuclei, and degeneration and replacement of muscle tissue by fibrotic tissue. This leads to a progressive loss of muscle strength, due to the destruction of muscle tissue, along with the difficulty of muscle contraction of the remaining muscle, due to fibrotic tissue and muscle contractures. Ultimately, loss of ambulation (LoA) occurs in late childhood or early adolescence, in were screened. If necessary, the studies' authors were contacted. The specific search is described in Appendix A.

Inclusion/Exclusion Criteria
The inclusion criteria were as follows: (1) participants-males with DMD, with no age restrictions; (2) exposure-(i) LTBP4 haplotypes (i.e., IAAM/IAAM vs. others), occasionally a single SNP may be considered due to linkage disequilibrium; (ii) SPP1 genetic variants; (iii) other genetic variants in genes directly involved in the TGFβ pathway; (3) outcomes-LoA (main outcome), defined as the age at which the person requires continued use of a wheelchair because of an inability to walk independently, and cardiac function (secondary outcome), including age of onset of dilated cardiomyopathy (DCM) and age of onset of myocardial dysfunction. No language limitations were imposed.
The exclusion criteria were as follows: (1) participants-a population that was not mainly affected by DMD and the effect of genotype on participants with DMD could not be estimated; (2) exposure-genetic variants not directly involved in the TGFβ pathway.
The literature search was conducted independently by two reviewers (CP-M and IC-R), and disagreements were solved by consensus or with a third reviewer (VM-V).

Data Extraction
An ad hoc table was conducted with the following data extracted from the selected studies: (1) reference (authors and publication year), (2) exposure (gene and genetic variant or haplotype), (3) country, (4) sample size, (5) design, and (6) outcomes (LoA and cardiac function). The data extraction was conducted individually by two reviewers (CP-M and IC-R).

Risk of Bias Assessment
To assess the risk of bias, we used the Quality Assessment Tool for Observational Cohort and Cross-Sectional Studies from Study Quality Assessment Tools [29]. This tool has a 14-item checklist. Each item receives a score of good, fair, or poor and includes the following domains: research question, study population, groups recruited from the same population and uniform eligibility criteria, sample size justification, exposure assessed prior to measurement, sufficient timeframe to observe an effect, different levels of the exposure of interest, exposure measures and assessment, repeated exposure assessment, outcome measures, blinding of outcome assessors, follow-up rate, and statistical analysis.
The risk of bias assessment was conducted independently by two reviewers (CP-M and IC-R), and disagreements were solved by consensus or with a third reviewer (VM-V).

Grading the Quality of Evidence
The Grading of Recommendations, Assessment, Development and Evaluation tool (GRADE) was used to assess the strength of the evidence for the main outcome included in the meta-analyses [30]. Depending on the study design, risk of bias, inconsistency, indirect evidence, imprecision, publication bias, larger effect, possible confounding variables, and dose-response gradient, the GRADE tool rates the strength of evidence as high, moderate, low and very low for each intervention and outcome.

Data Synthesis
A narrative synthesis and an ad hoc table were generated for the main outcome. Additionally, forest plots were used to graphically depict the hazard ratio (HR) of LoA with their 95% confidence intervals (95% CIs) for each genetic variant. The forest plots differentiated between patients treated with or without glucocorticoids.
When two or more HRs for the LoA associated with a specific genetic variant and their 95% CI were available, a fixed-effects meta-analysis was performed [31]. If the confidence interval or standard error was not reported, it was calculated from the standard deviation or estimated from the p-value [32]. Heterogeneity, as measured by the I 2 test [28,33], was considered not important <40%, substantial 30-60%, important 50-90%, and considerable >75%. The p-value of heterogeneity was also considered to determine whether heterogeneity was significant. The meta-analyses differentiated whether glucocorticoid-treated patients, patients who were not treated with glucocorticoids, or both types of patients were included.
STATA SE software, version 15 (StataCorp, College Station, TX, USA), was used to conduct the statistical analyses.

Infographics
An infographic of the manuscript was created using Canva software ( Figure A1).

Loss of Ambulation-Cox Regression Analyses
Regarding the glucocorticoid subgroup, the LTBP4 haplotype IAAM was not associated with HR compared with other haplotypes [20,21]. A comparison of the SPP1 rs28357094 TT genotype with the GG/GT genotypes revealed an HR = 0.62 (95% CI = 0.42-0.92) [20]. Finally, regarding SPP1 rs11730582 and rs17524488, only rs11730582 genotype C showed a protective association, with an HR = 0.63 (95% CI: 0.45, 0.89) [21]. no significant effect [20,21]. SPP1 rs28357094 showed no association [20,24]. Regarding the glucocorticoid subgroup, the LTBP4 haplotype IAAM was not associated with HR compared with other haplotypes [20,21]. A comparison of the SPP1 rs28357094 TT genotype with the GG/GT genotypes revealed an HR = 0.62 (95% CI = 0.42-0.92) [20]. Finally, regarding SPP1 rs11730582 and rs17524488, only rs11730582 genotype C showed a protective association, with an HR = 0.63 (95% CI: 0.45, 0.89) [21].    Table 4 shows Kaplan-Meier analyses of the included studies. LTBP4 rs10880 was not significantly associated with the age at onset of DCM or myocardial dysfunction [19,25]. However, there was a trend toward a higher percentage of patients without myocardial dysfunction with the rs10880 TT genotype [25]. In the glucocorticoid subgroup, the rs10880 genotype T (recessive model) was associated with a later age of onset of DCM (log-rank p < 0.05) [19]. SPP1 rs28357094 was not associated with the age at onset of DCM or myocardial dysfunction, although the genotype T tended to be harmful [19].

Assessment of the Risk of Bias
The studies fulfilled between 64.3% and 71.4% of the quality criteria proposed by the Quality Assessment Tool for Observational Cohort and Cross-Sectional Studies from Study Quality Assessment Tools. The objective of the study was not correctly specified in 12.5% of the studies. No study justified the sample size used to achieve a specific statistical power, and the authors used the available sample. An evaluation of different exposure levels or their changes over time was not necessary. Finally, blinding of the analysts was not considered necessary in any study. The complete risk assessment is available in Supplementary Table S2.

Evidence Assessment
Using the GRADE tool, the LTBP4 haplotype IAAM/IAAM, SPP1 rs28357094 TT vs. GG/GT and LTBP4 haplotype IAAM/IAAM vs. other haplotypes in glucocorticoidtreated patients showed low certainty of evidence. The complete assessment is detailed in Supplementary Table S3.

Meta-Analysis
Three meta-analyses were conducted that included the following studies: four studies of LTBP4 haplotype IAAM compared with other haplotypes in glucocorticoid-treated and nonglucocorticoid-treated patients, two studies of the involvement of SPP1 rs28357094 in glucocorticoid-treated and nonglucocorticoid-treated patients, and two studies of LTBP4 haplotype IAAM compared with other haplotypes in patients who were exclusively treated with glucocorticoids. Due to the scarcity of studies, we were unable to conduct meta-analyses of the remaining genotypes or of patients who were not treated with glucocorticoids.

Main Findings
This systematic review and meta-analysis provides an overview of the evidence supporting the associations of genetic variants with LoA and cardiac function. Our results suggest a notable effect of genetic variants involved in the TGFβ pathway, especially LTBP4 haplotype IAAM (in the recessive model), but not in patients who were exclusively treated with glucocorticoids, probably due to confounding factors. SPP1 rs28357094 did

Main Findings
This systematic review and meta-analysis provides an overview of the evidence supporting the associations of genetic variants with LoA and cardiac function. Our results suggest a notable effect of genetic variants involved in the TGFβ pathway, especially LTBP4 haplotype IAAM (in the recessive model), but not in patients who were exclusively treated with glucocorticoids, probably due to confounding factors. SPP1 rs28357094 did not display a significant association. However, the use of glucocorticoids by patients carrying the SPP1 rs28357094 genetic variant potentially increased the association in favor of genotype T (recessive model). Further research on THBS1 genetic variants is needed, and the limited evidence available indicates an interaction with LTBP4 that might improve the DMD phenotype. Finally, more research is needed on other proinflammatory and profibrotic pathways that might exert some effect and distort the effect observed in our study.

Interpretation
In dystrophic mice there is a 36-base-pair insertion/deletion site in a proline-rich domain of Ltbp4. The 12 amino acid insertion confers resistance to proteolysis of the TGFβ-LTBP4 complex, reducing TGFβ activity and therefore muscle fibrosis. Although this indel is not present in humans, the haplotypes mentioned in the previous sections (i.e., IAAM and VTTT) represent more than 80% of the total of all possible combinations in the human population. The haplotype IAAM behaves like the 36-base-pair insertion in Ltbp4 of dystrophic mice, giving LTBP4 increased binding avidity to TGFβ, thereby reducing its activity [38]. Our meta-analysis shows that the LTBP4 haplotype IAAM (or the rs10880 genotype T) is a protective factor resulting in prolonged ambulation. However, in patients treated exclusively with glucocorticoids, the haplotype IAAM had no effect, probably due to the inclusion of only two studies, one including patients from China, in which the LTBP4 haplotype IAAM does not appear to be associated with LoA, perhaps due to the genetic/ethnic factors described below. Considering the Kaplan-Meier analyses, Barp  , which might imply differences in the care and management or genetic/ethnicity of these patients. Thus, in the cohort analyzed by Barp A et al., the use of glucocorticoids was still not a universal standard, and a relatively low use of these drugs was reported. Chen M et al. investigated Asian cohorts, in which the effect of the LTBP4 haplotype IAAM was different than that on Caucasian cohorts. Furthermore, previous evidence [38] has identified rs710160 as a genetic variant that may modulate the association obtained for the LTBP4 haplotype IAAM with LoA. Thus, rs710160 genotype C together with the IAAM haplotype leads to less profibrotic signaling, potentially resulting in milder DMD. Additionally, rs710160 has a significant linkage disequilibrium with rs10880 in Caucasian populations but not in other populations, which might explain why the haplotype IAAM has a higher association with LoA in the Caucasian population but not in the Asian population, as described in the study by Chen M et al.
Thrombospondin-1, encoded by the THBS1 gene (location 15q14), is a potent activator of the latent TGFβ complex. To activate TGF-β1, thrombospondin-1 interacts with the N-terminal region of the latency-associated protein, which binds non-covalently to TGF-β1. Thus, a trimolecular complex is formed, a conformational change occurs, and the reactivity of TGF-β1 is altered. TGF-β1 activation by thrombospondin-1 may be essential for the development of the heart, liver, bones, testes, and hematopoietic systems, among other tissues and organs [39]. Moreover, it exerts an anti-angiogenic effect by decreasing the activities of the nitric oxide and vascular endothelial growth factor (VEGF) pathway, and blocking endothelial cell migration, [39,40], which is harmful in animal models [41]. Finally, thrombospondin-1 could be useful in maintaining low but constant levels of active TGF-β1 [39]. In our study, THBS1 rs2725797 genotype T is associated with lower thrombospondin-1 activity, with an additive effect with LTBP4 rs10880 allele C, which may improve the phenotype of the disease, as suggested by the data presented in this review.
SPP1 is a cytokine secreted by macrophages and myoblasts. It belongs to the family of small integrin-binding ligand N-linked glycoprotein secreted phosphoproteins and is expressed in numerous tissues in response to tissue damage/regeneration and inflammatory response. Alternative splicing and post-translational modifications make it a difficult cytokine to study. In mdx models without osteopontin, mice had reduced TGFβ levels, fibrosis, and increased strength, with reduced infiltration of inflammatory cells, neutrophils, and natural-killer T cells and increased numbers of regulatory T cells and M2 macrophages. Although the reduction in osteopontin could delay regeneration in acute damage, in mdx models (and probably in DMD) it decreases fibrosis caused by chronic damage [38,42]. Although in our study SPP1 rs28357094 was not associated with LoA in a population that was not stratified by glucocorticoid use, in the glucocorticoid subgroup, it appears to exert a glucocorticoid-dependent effect on LoA in regression analyses. In fact, in subgroup analyses stratified by glucocorticoid treatment, this genetic variant showed no association with LoA in the glucocorticoid-untreated subgroup or when a low percentage of glucocorticoidtreated patients was included. SPP1 rs28357094 genotype G was detrimental to patients using glucocorticoids, with earlier LoA. Allele G is associated with lower osteopontin transcriptional activity under basal conditions [43]. However, SPP1 likely contains glucocorticoid receptor elements, which modulate osteopontin expression through interactions with NF-kB, estrogens, and glucocorticoids. Thus, the SPP1 rs28357094 genotype G might increase osteopontin expression by 3-fold in the presence of glucocorticoids, increasing profibrotic signaling [44]. SPP1 rs11730582 might also use a similar mechanism. Genotype C, which results in higher initial osteopontin levels [45], is significantly associated with a later LoA in glucocorticoid-treated patients, suggesting that its effect is glucocorticoiddependent, similar to rs28357094. Therefore, SPP1 rs28357094 and perhaps rs11730582 are postulated to function not as disease modifiers but as predictors of a good or poor response to glucocorticoid treatment.
Regarding cardiac complications, no clear association with genetic variants can be observed due to the few included studies. LTBP4 rs10880 genotype T tended to delay the age of DCM onset, with a significant association in patients treated with glucocorticoids [19]. This result is consistent with the data obtained for LoA. This association may be due to an additive effect of the LTBP4 haplotype IAAM and glucocorticoids, but this requires further research. In contrast, another study did not detect an association of LTBP4 with the development of left ventricular dysfunction. Interestingly, it was found that, among patients without left ventricular dysfunction, there was a higher proportion of patients with the rs10880 genotype T [25]. These discrepancies are, perhaps, due to the relatively small sample size, the age of the cohorts, or the categorization of patients into glucocorticoidtreated/untreated groups. Furthermore, the SPP1 rs28357094 genotype T tended to be harmful in the development of DCM, including in the glucocorticoid subgroup [19]. This fact is interesting, since it is in the opposite direction to what happens in skeletal muscle. Conversely, it is in agreement with what has been observed in animal models, in which an overexpression of Spp1 causes myocarditis and myocardial dilation [46]. Future research is needed in this regard. Cardiac complications, and especially heart failure, are one of the leading causes of death in patients with DMD. Both the SPP1 genotype and glucocorticoids use could have some negative impact on the development of dilated cardiomyopathy. However, it should be noted that this does not necessarily imply progression to heart failure, and especially glucocorticoids could delay the onset of heart failure through other mechanisms.
The importance of the TGFβ pathway in the DMD phenotype may have several clinical implications: first, drug development aimed at the downregulation of TGFβ signaling to reduce muscle fibrosis. This approach includes angiotensin 1-7, halofuginone, anti-TGFβ1 antibodies, ixazomib, and angiotensin-II type 1 receptor blockers, which antagonize or downmodulate the TGFβ pathway with promising results in animal models [47][48][49][50][51]. However, some of these drugs have produced undesirable pleiotropic effects, which might be a problem in achieving a clinical benefit [49]. Another pathway is myostatin inhibition using compounds such as follistatin, ACE-031, domagrozumab, and the GDF11 propeptide that inhibit or antagonize the effect of myostatin [52][53][54][55][56][57]. However, the results from animal models and patients with DMD have raised some concerns related to other biological functions of myostatin that affect the metabolism and oxidative capacity of muscle fibers [52,58]. Second, genotyping of patients with DMD for LTBP4 and SPP1 could be considered, especially in the Caucasian population. Although knowledge of the genotype of LTBP4 haplotype/rs10880 and SPP1 rs28357094 would not alter the medical treatment of the patient, it can potentially provide the clinician and the patient and their family with more individualized prognosis as to their possible evolution and expected response to glucocorticoids. Third, genotyping can be considered in clinical trials of new drugs, such as treatments designed to restore dystrophin expression. Perhaps, by subgroup analysis according to genotype, part of the variability found in the results could be explained. It is possible that the patient's genotype can determine the response to small changes in dystrophin expression. Finally, the study of genetic variants in unknown DMD patients could be very interesting. Clinical DMD is well established, with diagnosis in early childhood. However, there are rare cases of dystrophinopathies, including Becker muscular dystrophy (BMD) and DMD, that remain undiagnosed until their debut in adolescence or adulthood, mainly due to cardiac involvement. [59][60][61][62][63]. Moreover, there are phenotypically intermediate forms of diagnosed dystrophinopathies, with slowly evolving DMD or rapidly evolving BMD [38]. These exceptional cases could be an opportunity to study the effect of genetic factors (known and unknown) as well as environmental factors on the progression of DMD and other dystrophinopathies.

Limitations
Some limitations should be acknowledged. First, the scarcity of studies included in the meta-analyses limits the statistical power and external validity of the results in larger cohorts of patients with DMD. Second, due to the limited number of studies, we were unable to perform publication bias, metaregression, or sensitivity analyses. This could question the association for LTBP4 and LoA. In systematic reviews, it is not uncommon for the first published studies to show stronger associations than subsequent ones, and they are even less likely to be published unless they question previous findings. Third, in some studies, the LTBP4 haplotype was determined using a single SNP (mainly rs10880). Despite the large linkage disequilibrium, it could slightly underestimate the observed result. Fourth, in the Kaplan-Meier survival analysis, we were not always able to consider ethnic differences, the effects of other genetic variants, or a different proportion of glucocorticoid-treated patients, which would potentially modify the observed association. Fifth, participants who are candidates for exon 8 or 44 skipping tend to have slower disease progression [64,65]. Most studies did not consider this possible confounding genotype, and although it is possible that there are no statistically significant differences in the prevalence of these genotypes in the included studies, it cannot be ruled out, considering the ethnic and geographical diversity of the participants, overestimating or underestimating the effect observed for the genetic variants studied. Sixth, other genetic variants, such as LTBP4 rs710160, and their possible differences in prevalence in different populations, might modify the association expected for the LTBP4 haplotype IAAM. Seventh, some confidence intervals of hazard ratios were estimated from the p-value, which might differ from the true confidence interval. Eighth, in general, the authors did not determine the level of fibrosis of the patients' muscle tissue, which, in theory, should correlate negatively with LoA or cardiac function. This is probably due to the complexity of the procedures to obtain this information.

Conclusions
Some of the variability in the DMD phenotype can be explained by the genetic variants directly involved in the TGFβ pathway. Thus, the LTBP4 haplotype IAAM may predict a subsequent LoA, independent of glucocorticoid use, although glucocorticoids may exert an additive effect on this outcome. This haplotype might also be associated with better cardiac function, although more research is needed to confirm or reject this potential association. SPP1 rs28357094 and possibly rs11730582 exert a glucocorticoid-dependent effect, with the rs28357094 genotype G serving as a predictor of a poor response. Further research is needed to establish the possible interaction of LTBP4 rs710160 and THBS1 rs2725797. Moreover, future studies for the LTBP4 haplotype IAAM in non-Caucasian cohorts are required to confirm the findings in this review, as well as the possibility that the association of the LTBP4 haplotype IAAM and LoA is conditioned by the genetic variant rs710160, being relevant in non-Caucasian populations. Knowledge of the effect of the TGFβ pathway on the DMD phenotype provides a possibility of developing new treatments, and genotyping of the genetic variants involved could provide the clinician and the patient with more information on their possible progression or response to certain pharmacological treatments.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/10 .3390/ph14080798/s1, Table S1: Excluded studies with reasons, Table S2: Risk of bias assessment, Table S3: Grades of Recommendation, Assessment, Development, and Evaluation of HR of LoA. Data Availability Statement: The datasets generated and analyzed are available from the corresponding author.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix A.2 Grey Literature
Not specified.