Causal effects of specific gut microbiota on musculoskeletal diseases: a bidirectional two-sample Mendelian randomization study

Background Recent observational studies and clinical trials demonstrated an association between gut microbiota and musculoskeletal (MSK) diseases. Nonetheless, whether the gut microbiota composition has a causal effect on the risk of MSK diseases remains unclear. Methods Based on large-scale genome-wide association studies (GWAS), we performed a two-sample Mendelian randomization (MR) analysis to investigate the causal relationship between gut microbiota and six MSK diseases, namely osteoporosis (OP), fracture, sarcopenia, low back pain (LBP), rheumatoid arthritis (RA), and ankylosing spondylitis (AS). Instrumental variables for 211 gut microbiota taxa were obtained from the largest available GWAS meta-analysis (n = 18,340) conducted by the MiBioGen consortium. And the summary-level data for six MSK diseases were derived from published GWAS. The inverse-variance weighted (IVW) method was conducted as a primary analysis to estimate the causal effect, and the robustness of the results was tested via sensitivity analyses using multiple methods. The Bonferroni-corrected test was used to determine the strength of the causal relationship between gut microbiota and various MSK diseases. Finally, a reverse MR analysis was applied to evaluate reverse causality. Results According to the IVW method, we found 57 suggestive causal relationships and 3 significant causal relationships between gut microbiota and MSK diseases. Among them, Genus Bifidobacterium (β: 0.035, 95% CI: 0.013–0.058, p = 0.0002) was associated with increased left handgrip strength, Genus Oxalobacter (OR: 1.151, 95% CI: 1.065–1.245, p = 0.0003) was correlated with an increased risk of LBP, and Family Oxalobacteraceae (OR: 0.792, 95% CI: 0.698–0.899, p = 0.0003) was linked with a decreased risk of RA. Subsequently, sensitivity analyses revealed no heterogeneity, directional pleiotropy, or outliers for the causal effect of specific gut microbiota on MSK diseases (p > 0.05). Reverse MR analysis showed fracture may result in a higher abundance of Family Bacteroidales (p = 0.030) and sarcopenia may lead to a higher abundance of Genus Sellimonas (p = 0.032). Conclusion Genetic evidence suggested a causal relationship between specific bacteria taxa and six MSK diseases, which highlights the association of the “gut-bone/muscle” axis. Further exploration of the potential microbiota-related mechanisms of bone and muscle metabolism might provide novel insights into the prevention and treatment of MSK diseases.


Introduction
Musculoskeletal (MSK) diseases encompass a spectrum of pathologies that affect muscles, soft tissues, joints and the spine (Murray et al., 2012).Globally, all MSK diseases combined accounted for 21.3% of the total years lived with disability (YLDs), second only to mental and behavioral disorders (23.2%) (March et al., 2014;Smith et al., 2014).Disability due to MSK diseases has increased by 45% from 1990 to 2010, in particular sarcopenia, osteoporosis (OP) and osteoarthritis (OA), and is expected to continue to rise with an increasingly obese, sedentary, and ageing population (Vos et al., 2012).In addition, MSK diseases also impose a significant economic burden on both patients and the health care system (Sebbag et al., 2019).Recent epidemiological studies have revealed that MSK diseases-related lost productivity accounts for 2% of the gross domestic product (GDP) in the European Union (Bevan, 2015).According to data from the 2012 National Health Interview Survey (NHIS), the annual cost of direct treatment and lost wages for MSK diseases in the United States was estimated at $213 billion (Rosenfeld et al., 2018).Therefore, the prevention and management of MSK diseases has been globally recognized as a crucial public health issue that needs to be solved urgently.
The gut microbiota is a large and complex community of microbial species inhabiting the human gastrointestinal tract, which plays a critical role in human health and diseases (Li et al., 2022).Human gut microbiota can influence host physiology by regulating multiple processes, including inflammation, oxidative stress, immune function, and anabolic balance (Ticinesi et al., 2019).Growing evidence suggests that many MSK diseases, such as OP, OA, sarcopenia, and low back pain (LBP) are closely correlated with altered gut microbiota (Biver et al., 2019).According to a recent study by Xu et al., the gut microbiota compositions of patients with osteoporosis were significantly different from those of healthy controls, especially the enriched Dialister and Faecalibacterium genera (Xu et al., 2020).Another observational study indicated that the abundance of Lactobacilli, Bacteroides, and Prevotella is higher and the abundance of Enterobacteriaceae is lower in normal individuals than in those with sarcopenia (Castro-Mejia et al., 2020).In addition, the gut microbiota can also be involved in the development of rheumatoid arthritis (RA) and ankylosing spondylitis (AS) by altering intestinal barrier permeability, regulating expression levels of host hormones, and mediating immune responses (Asquith et al., 2019).Although the correlations between gut microbiota and Musculoskeletal (MSK) diseases have been reported by previous researches, these findings are mainly based on observational and cross-sectional analyses, and it is still unclear whether there is a causal link between the gut microbiome and MSK diseases.
Mendelian randomization (MR) is a groundbreaking analytical method that provides an unbiased estimation of the causal link between phenotypes (Smith and Ebrahim, 2003).Compared with traditional observational studies, the MR study uses genetic variation as instrumental variables (IVs) to avoid the influence of traditional confounding factors and reverse causality, which provides robust evidence on the mechanisms of disease pathogenesis and the efficacy of treatments (Davey and Hemani, 2014).Given the aforementioned findings and the MR method, our objective was to conduct a bidirectional two-sample MR analysis using extensive genome-wide association study (GWAS) data to explore the causality between specific gut microbiota and six MSK diseases, including OP, fracture, sarcopenia, LBP, RA, and AS.

Study design
GWAS summary data focusing on the relationships between genetics and disease have been used to detect a causal relationship between the gut microbiome and six MSK diseases (OP, fracture, sarcopenia, LBP, RA, and AS) through bidirectional two-sample MR analysis (Visscher et al., 2012;Jia et al., 2019).Our first step was to determine whether the gut microbiome contributes to the prevention or promotion of six MSK diseases, by selecting the gut microbiome as the exposure and the six MSK diseases as the outcome.Furthermore, we examined changes in gut microbiota following the development of six MSK diseases.According to Bownden et al., the two-sample MR needs to satisfy the following three assumptions (Bowden et al., 2015): (Murray et al., 2012) The instruments of genetic variations should be robustly associated with gut microbiota; (b) The genetic variations should not be associated with any confounders of gut microbiota and six MSK diseases; (c) The genetic variations should affect six MSK diseases solely through gut microbiota, not via other pathways (Sanderson et al., 2019) (Figure 1A).All studies included in the cited GWASs were approved by relevant review committees and all participants provided informed consent.The flow chart of the MR analysis was presented in Figure 1B.

Genome-wide association studies sources
Our gut microbiota summary data came from a large multi-ethnic GWAS meta-analysis that included 18,340 European individuals and individuals from 24 cohorts (MiBioGen Consortium).Three different variable regions of the 16S rRNA gene were targeted in order to profile the microbial composition.And only the taxa found in >10% of the samples were included in the quantitative microbiome trait loci (mbQTL) mapping study for each cohort.Furthermore, after adjustment for age, sex, technical covariates, and genetic principal components, Spearman's correlation analysis was conducted to identify genetic loci that affected the covariate-adjusted abundance of bacterial taxa (Kurilshikov et al., 2021).
In the study, we extracted GWAS summary statistics for six MSK diseases from publicly GWAS analyses.The GWAS summary data on OP were mainly extracted from UK Biobank, including 7,547 cases and 455,386 controls of European ancestry (Sudlow et al., 2015).Summary statistics for fracture were derived from a GWAS metaanalysis, including 44,502 cases and 415,887 controls of European ancestry (Morris et al., 2019).Summary data for ALM and handgrip strength were obtained from the UK Biobank study.An analysis of 450,243 UK Biobank cohort participants was conducted to quantify ALM-related values by summing fat-free mass (Pei et al., 2020).For handgrip strength, UK Biobank provided GWAS summary statistics on right and left handgrip strength based on 461,089 and 461,026 United Kingdom people, respectively (Sudlow et al., 2015).The data on LBP came from a meta-analysis including 13,178 cases and 164,682 controls of European ancestry (https://gwas.mrcieu.ac.uk/datasets/ ukb-d-M13_LOWBACKPAIN).For RA, we collected summary statistics from a GWAS meta-analysis that includes 14,361 RA cases and 43,923 controls of European ancestry (Okada et al., 2014).The genetic data of AS risk were derived from a large GWAS involving 9,069 cases and 1,550 controls of European (Cortes et al., 2013).Detailed information on the demographic characteristics of selected summary-level GWASs applied in this study was shown in Supplementary Table S1.

Selection of genetic instrumental variables
During the selection of the optimal IVs, the following quality control steps were taken to ensure the validity and accuracy of the conclusions on the causal relationship between the gut microbiome and the six MSK diseases (Murray et al., 2012) The instrumental variables selected for analysis need to be highly related to the corresponding exposures (We choose significant SNPs based on a loose cutoff of p < 1 × 10 −5 to ensure sufficient instrumental variables for screening).(March et al., 2014) The instrumental variables are mutually independent and avoid the offset caused by linkage disequilibrium (LD) between the SNPs (r 2 < 0.01, LD distance = 500 kb).(Smith et al., 2014) We eliminated instrumental variables with an F-statistic less than 10 to minimize potential weak instrument bias F = R 2 (n-k-1)/k(1-R 2 ) (n is the sample size, k is the number of included instrumental variables, and R 2 is the exposure variance explained by the selected SNPs).

Statistical analysis
The inverse variance weighted (IVW) method was employed as the main analysis, to obtain an unbiased estimate of the causal relationship between gut microbiota and MSK diseases.Furthermore, the weighted median, MR Egger, simple mode, and weighted mode methods were applied as additional methods to estimate causal effects under different conditions.The weighted median method could combine data on multiple genetic variants into a single causal estimate and provide a consistent estimate if at least half of the weight is derived from valid IVs (Bowden et al., 2016).The MR-Egger method could evaluate whether genetic variants have directional pleiotropy and provide a consistent estimate of the causal effect (Bowden et al., 2015).The intercept of MR-Egger regression was calculated to assess horizontal pleiotropy, and p value >0.05 indicated that the possibility of pleiotropy effect in the causal analysis is weak (Verbanck et al., 2018).Cochran's Q test was derived from IVW estimation and used to detect heterogeneity among instrumental variables (Pierce and Burgess, 2013).In addition, we applied the Mendelian randomization pleiotropy residual sum and outlier (MR-PRESSO) method to determine horizontal pleiotropy and correct potential outliers.The leave-one-out method was used for sensitivity analysis, which sequentially removed one of the SNPs and used the remaining SNPs as instrumental variables for two-sample MR analysis to judge the degree of influence of causal association effect by a single SNP.In addition, we also performed a reverse-direction MR analysis to determine whether there was a reverse-direction causal relationship.
To obtain a more accurate assessment of the causal relationship, we applied Bonferroni correction to determine the significance of multiple-testing at each feature level (p < 0.05/n, where n is the number of bacterial taxa included in each feature level).Hence, the multipletesting significance was 5.56 × 10 −3 , 3.13 × 10 −3 , 2.50 × 10 −3 , 1.56 × 10 −3 , and 4.20 × 10 −4 , respectively, for phylum, class, order, family, and genus.The 'TwoSampleMR' package and the 'MRPRESSO' package in R software (version 4.1.3)were used for all MR analyses.

Selection of instrumental variables
To begin with, 14,587 SNPs correlated with the gut microbiota were identified as IVs from the MiBioGen Consortium at a relatively loose significance level (p < 1 × 10 −5 ).It contained 211 bacterial traits, including 131 genera, 35 families, 20 orders, 16 classes, and 9 phyla.After a series of quality control steps, 2,236 SNPs were finally included in the analysis.In addition, the F-statistics of all IVs were >10, indicating no evidence of weak instrument bias.

Sensitivity analysis, Bonferroni-corrected test, and reverse analysis
We performed a series of sensitivity analyses to test the heterogeneity and horizontal pleiotropy of the selected IVs.Based on Cochran's Q test, we observed no significant heterogeneity (p > 0.05) (Table 9).All p values of the MR-Egger intercept tests were > 0.05, which indicated no horizontal pleiotropy.Furthermore, we also did not discover any outliers through the MR-PRESSO global test (Table 9).Detailed scatter plots for each MR method analysis were shown in Supplementary Figures S2-S9.And results from a leaveone-out analysis demonstrated that no SNP was an influential outlier (Supplementary Figures S9-S16).According to the results of the Bonferroni-corrected test, Genus Bifidobacterium (β: 0.035, 95% CI: 0.013-0.058,p = 0.0002) was significantly correlated with left handgrip strength.In addition, higher level of Genus Oxalobacter retains a strong causal association with an increased risk of LBP (OR: 1.151, 95% CI: 1.065-1.245,p = 0.0003) and higher level of Family Oxalobacteraceae (OR: 0.792, 95% CI: 0.698-0.899,p = 0.0003) retains a strong causal association with a decreased risk of RA.In addition, reverse analysis demonstrates that fracture may result in a higher abundance of Family Bacteroidales (p = 0.030) and sarcopenia may lead to a higher abundance of Genus Sellimonas (p = 0.032) (Supplementary Table S10).

Discussion
To our knowledge, this was the first comprehensive and in-depth investigation of causal associations between gut microbiota and six common MSK diseases (OP, fracture, sarcopenia, LBP, RA, and AS) based on publicly available GWAS data.According to the findings of our study, a total of 57 gut microbiota were potentially causally associated with the progression of MSK diseases.After multiple testing corrections, we observed a significant causal relationship between Genus Bifidobacterium, Genus Oxalobacter, and Family Oxalobacteraceae with the risk of sarcopenia, LBP and RA, respectively.These findings might provide new ideas for the future treatment of MSK diseases by targeting the specific gut microbiota.
The gut microbiota consists of trillions of bacteria in the gastrointestinal tract, which have functions such as improving intestinal permeability, attenuating the inflammatory response, and participating in the immune regulation of the skeletal system (Rizzoli, 2019).In recent years, the "gut-bone/muscle" axis has received increasing attention in the field of bone health and orthopedic diseases and multiple studies indicated that gut microbiota composition is involved in the MSK disease pathogenesis through multiple pathways (Tu et al., 2021;Zhang et al., 2021).In our study, we discovered that genetically proxied higher abundance of Genus Bifidobacterium was significantly correlated with increased grip strength.And consistent directional effects for all analyses were observed in both MR Egger and weighted median methods, which suggests that Bifidobacterium might be a promising target for sarcopenia prevention.According to recent studies, Bifidobacterium has been identified as a critical taxon for frailty and sarcopenia in the elderly (Ticinesi et al., 2019).In a case-control study of elderly Chinese women, metagenomic sequencing of the gut microbiota revealed that the abundance of Prevotella and Bifidobacterium in the healthy control group was higher than those in the sarcopenia group (Wang et al., 2023).An animal experiment reported that Bifidobacterium longum could improve exercise-associated peripheral fatigue indicators (forelimb grip strength) and oxidative stress-related damage indicators in mice, which indicated Bifidobacterium may be involved in the occurrence and progression of age-related physical frailty and sarcopenia (Huang et al., 2020).In addition, a double-blind randomized clinical trial of older adults 65 aged 65 and over demonstrated that long-term intake of a prebiotic supplement containing Bifidobacterium significantly improved grip strength and reduced the level of fatigue (Buigues et al., 2016).Bifidobacterium might affect sarcopenia through its metabolites, such as short-chain fatty acids (acetate, propionate, and butyrate), which are important for immune regulation and metabolic   homeostasis (Wang et al., 2015).Bifidobacterium also contributed to the absorption and utilization of vitamin D and minerals such as calcium, phosphorus, and iron, which are essential for muscle metabolism (Montazeri-Najafabady et al., 2019).Furthermore, Bifidobacterium was rich in genes encoding proteins related to the carbohydrate and amino acid transport and metabolism, influencing protein synthesis and nutrient utilization, which might also be one of the potential mechanisms by which Bifidobacterium prevents sarcopenia (Tu et al., 2023).
Our MR analysis also identified a causal relationship between Genus Oxalobacter, Family Oxalobacteraceae and the prevalence of LBP and RA, respectively.Oxalobacter, is a unique anaerobic bacterium that plays an important role in degrading dietary oxalate and stimulating oxalate secretion by the gut mucosa (Liu et al., 2017;  Daniel et al., 2021).Recently, intensive studies have reported that Oxalobacter is closely related to the pathogenesis of kidney stones, hyperoxaluria, and chronic kidney disease (CKD) (Sidhu et al., 1999;Turkmen and Erdur, 2015;Falony, 2018).In addition, Li et al. observed that the diversity and relative abundance of gut microbiota in individuals with RA were significantly different from those in healthy controls.Specifically, the relative abundance of Pelagibacterium was higher and Oxalobacter and Blautia were lower in patients with RA (Li et al., 2021), which was consistent with our findings.Numerous studies emphasized the key role of the gut microbiota in RA pathogenesis, through a variety of mechanisms, including the production of proinflammatory metabolites, impairment of the intestinal mucosal barrier, and regulation of hormones throughout the body (Zhao et al., 2022).However, studies on the role of Oxalobacter in LBP are currently limited.Thus, it is necessary to further study the possible role of Oxalobacter in LBP and RA development.
In addition, we also observed some promising gut microbiota were associated with the risk of OP and AS.Several mechanisms seem  to offer preliminary explanation for the relationship between gut flora and OP.On the one hand, the gut microbiota could regulate the Treg/ Th17 cells balance or relevant cytokines through the immune system, affecting the intestinal and systemic immune states, thereby establishing a dynamic balance between osteoblasts (OB) and osteoclasts (OC), which is important for normal bone mass maintenance (Locantore et al., 2020).On the other hand, microbial metabolites such as SFAs, secondary bile acid (SBA), and indole derivatives participate in the reconstruction of bone resorption, metabolism and fracture healing by providing energy to gut epithelium cells and promoting calcium and phosphorus absorption (Tu et al., 2021).In addition, intestinal flora imbalance could promote the systemic chronic inflammatory response, which is closely related to the occurrence of AS (Guggino et al., 2021).Stebbings et al. demonstrated that the abundance of Lachnospiraceae, Prevotellaceae, and Bacteroidaceae was higher and that of Ruminococcaceae and Rikenellaceae was lower in AS patients compared with the healthy controls (Stebbings et al., 2002).According to an animal study in mice, regulating the gut microbiota and improving gut barrier function could delay AS progression (Liu et al., 2019).On the basis of the above evidence, we can further investigate how specific gut microbiota affect OP and AS.
Reverse MR analysis showed that fracture may lead to a higher abundance of Family Bacteroidales and sarcopenia may result in a lower abundance of Genus Sellimonas.According to a case-control study of 50 elderly patients, in patients with fracture, the relative abundance of Bacteroidales was higher and Lachnospirales was lower compared to healthy controls, suggesting that altering the microbiome may be an effective way to reduce fracture risk (Rosello-Anon et al., 2023).Another study of postmenopausal Japanese women reported that Bacteroidels and Rikenellaceae may be involved in bone metabolism and fracture risk (Ozaki et al., 2021).Furthermore, an animal study has also revealed that the proportion of Firmicutes and Bacteroidetes increased significantly by 16S RNA sequencing in postmenopausal osteoporosis mice (Wen et al., 2020).In addition, Sellimonas has been shown to be involved in restoring the balance of the intestinal flora (Munoz et al., 2020).Despite this, there has been a very limited amount of literature about Sellimonas.Thus, the underlying mechanism warrants future investigation.Our study had several strengths.First of all, we used MR analysis to infer that the relationship between the gut microbiome and MSK diseases should be less susceptible to confounding and reverse causation than traditional observational analyses.Additionally, we analyzed the causal effect of each taxon on MSK diseases from the genus to the phylum level, which provides guidance for the prevention and treatment of MSK diseases by targeting specific gut microbiota in clinical practice.Third, the combination of stringent quality control procedures and multiple sensitivity analysis approaches enhances the credibility and robustness in the causal relationships inferred from the MR study.However, the present study also had some limitations that should be noted.Firstly, due to the lack of demographic data in the original study (e.g., gender, age, and race), we could not perform further subgroup analyses to obtain more specific effect relationships.Secondly, the GWAS data on gut microbiota used in this study are based on the population cohort from the largest macro-genome sequencing study to date.In the future, summary data of other gut microbiota needs to be obtained to provide a more comprehensive assessment of the causal relationship between gut microbes and the risk of MSK diseases.Last but not least, this study was confined to individuals of European origin and other populations require further MR studies since causal relationships may vary from race to race.
Overall, our study comprehensively evaluated the potential causal association between gut microbiota and six common MSK diseases.These findings provided new insights into the prevention, progression, and treatment of MSK diseases through targeting specific gut microbiota.In the future, more clinical trials and mechanism studies are needed to explore the exact mechanisms underlying the interactions between the gut microbiota and the prevalence of MSK diseases.Forest plot of the causality between gut microbiota with the risk of rheumatoid arthritis.The estimates: Inverse variance weighted (IVW) results of gut microbiota and rheumatoid arthritis risk; p-value: p-value of the estimate.OR, odds ratio; SNP, single-nucleotide polymorphism.

Conclusion
In summary, our study provided evidence of a potential causal relationship between specific gut microbiota and MSK diseases through the bi-directional MR analysis.These findings opened up new avenues for exploring the underlying mechanisms by which gut microbiota influence MSK diseases and contributed to the development of targeted interventions and personalized treatments.Forest plot of the causality between gut microbiota with the risk of ankylosing spondylitis.The estimates: Inverse variance weighted (IVW) results of gut microbiota and ankylosing spondylitis risk; p-value: p-value of the estimate.OR, odds ratio; SNP, single-nucleotide polymorphism.
0.011-0.060,p = 0.005) were also positively associated with left grip (Figures4A,B).And the MR estimates of MR Egger and weighted median indicated that genetically predicted Family Bifidobacteriaceae, Order Bifidobacteriales were positively related to right and left handgrip strength (P <0.05), suggesting Family Bifidobacteriaceae, Order Bifidobacteriales were protective factors for sarcopenia (Table

FIGURE 2
FIGURE 2Forest plot of the causality between gut microbiota with the risk of osteoporosis.The estimates: Inverse variance weighted (IVW) results of gut microbiota and osteoporosis risk; p-value: p-value of the estimate.OR, odds ratio; SNP, single-nucleotide polymorphism.

FIGURE 3
FIGURE 3Forest plot of the causality between gut microbiota with the risk of fracture.The estimates: Inverse variance weighted (IVW) results of gut microbiota and fracture risk; p-value: p-value of the estimate.OR, odds ratio; SNP, single-nucleotide polymorphism.

FIGURE 4
FIGURE 4 Forest plot of the causality between gut microbiota with sarcopenia-related traits.The estimates: Inverse variance weighted (IVW) results of gut microbiota and sarcopenia-related traits; p-value: p-value of the estimate.(A) Hand grip strength (Right); (B) Hand grip strength (Left); (C) Appendicular lean mass.SNP, single-nucleotide polymorphism.

TABLE 1
Mendelian randomisation (MR) results of causal effects between gut microbiome and the risk of osteoporosis.

TABLE 2
Mendelian randomisation (MR) results of causal effects between gut microbiome and the risk of fracture.

TABLE 3
Mendelian randomisation (MR) results of causal effects between gut microbiome and hand grip strength (Right).

TABLE 4
Mendelian randomisation (MR) results of causal effects between gut microbiome and hand grip strength (Left).

TABLE 5
Mendelian randomisation (MR) results of causal effects between gut microbiome and appendicular lean mass.

TABLE 6
Mendelian randomisation (MR) results of causal effects between gut microbiome and the risk of low back pain.
SNP, single nucleotide polymorphism; MR, mendelian randomization; SE, standard error; 95% CI, 95% confidence interval.FIGURE 5Forest plot of the causality between gut microbiota with the risk of low back pain.The estimates: Inverse variance weighted (IVW) results of gut microbiota and low back pain risk; p-value: p-value of the estimate.OR, odds ratio; SNP, single-nucleotide polymorphism.

TABLE 7
Mendelian randomisation (MR) results of causal effects between gut microbiome and the risk of rheumatoid arthritis.

TABLE 8 Mendelian
randomisation (MR) results of causal effects between gut microbiome and the risk of ankylosing spondylitis.

TABLE 9
Sensitivity analysis of the Mendelian randomization (MR) analysis results of gut microbiota and musculoskeletal diseases.