Exploring the Causal Relationship Between Osteoporosis and Rheumatoid Arthritis: A Bidirectional Mendelian Randomization Study

Introduction

Rheumatoid arthritis (RA) is an autoimmune disease with a prevalence of 0.5%–1.0%, more common among the aging population, and affecting women more than men.^1–3 RA primarily affects the synovial membrane and cartilage, leading to bone erosion, which can result in severe disability and increased mortality.^4,5 Osteoporosis (OP) is one of the serious complications in advanced RA, increasing the risk of low-trauma fractures and significantly impacting patients’ health and quality of life.^6,7

Osteoporosis is a chronic metabolic bone disease, commonly diagnosed by dual-energy X-ray absorptiometry (DXA), characterized by decreased bone mineral density (BMD).^8–10 OP is globally prevalent, affecting approximately 200 million people, with incidence rising with age and being more common in women. Genetic factors, lifestyle, and medical history are also important contributors to OP prevalence.^11,12

Epidemiological studies indicate that 60%-80% of RA patients also suffer from OP. Numerous observational studies show a strong relationship between OP and RA.^13–15 RA patients often experience focal or extensive bone involvement.^16–19 Guler-Yuksel et al reported that in 381 newly diagnosed active RA patients, 11% and 25% had osteoporosis in the spine and hip, respectively.²⁰ Synovial inflammation leads to cortical bone loss and marginal bone erosion around the joints, with RA duration and severity being independent risk factors for vertebral fractures. However, the exact mechanism behind the coexistence of RA and OP remains unclear.^21–23

Understanding the causal relationship between OP and RA is crucial for prevention and treatment. However, current observational studies with various confounding factors often fail to provide convincing evidence. Mendelian randomization (MR) is a method used to assess causal relationships between exposures and outcomes, as it employs genetic variables (IVs) that are equally, randomly, and independently distributed. Unlike traditional clinical randomized controlled trials, MR avoids potential confounding factors.

This study aims to explore the potential bidirectional causal effects between RA and OP. Using a two-sample bidirectional Mendelian randomization approach and data from genome-wide association studies (GWAS), this study seeks to clarify whether a direct causal relationship exists between these two diseases. The findings may enhance our understanding of the interaction between RA and OP and provide better management and preventive strategies for high-risk patients, ultimately improving patient care and prognosis.

Methods

Study Design and Data Sources

In this two-sample MR study, we used single nucleotide polymorphisms (SNPs) as instrumental variables (IVs) and GWAS data to determine the causal relationship between osteoporosis and RA. The study design and MR assumptions are summarized in Figure 1. Genetic data for osteoporosis and RA were obtained from published GWAS, with details provided in Table 1.As this study was based on previously published GWAS summary data, approval from an institutional review board was not necessary, and all participants provided informed consent beforehand.

Table 1 Data Sources Used in This Study

Figure 1 The MR design, guided by strict assumptions. Assumption 1: Genetic variants (instrumental variables) must be strongly associated with the risk factor. Assumption 2: These genetic variants should not be linked to confounding factors. Assumption 3: The variants should influence the outcome solely through the risk factor, not via other pathways. MR reduces confounding and reverse causation, enhancing causal inference because genetic variants are randomly allocated at conception, minimizing bias. Abbreviations:MR, Mendelian randomization;BMD,Bone mineral density; PRESSO, Pleiotropy Residual Sum and Outlier.

Genetic Instrument Selection

In this study, genetic instruments for rheumatoid arthritis (RA) and bone mineral density (BMD) were carefully selected based on stringent criteria to ensure the robustness of the Mendelian randomization (MR) analysis. Single nucleotide polymorphisms (SNPs) were identified using a genome-wide significance threshold (P<5E-8) and a minor allele frequency (MAF) greater than 0.01.²⁴ This approach was employed to avoid the inclusion of weak instruments that could introduce bias into the results, thereby enhancing the validity of the causal inferences made in the study.

Linkage Disequilibrium (LD) Management: To manage linkage disequilibrium (LD), a clumping procedure was implemented. This procedure used a stringent R^2 threshold of less than 0.001 within a 10,000 kb window to address significant LD.²⁵ In cases where SNPs were in high LD, the SNP with the lower P-value was retained. This careful selection process ensured that the SNPs used as instrumental variables were independent of each other, reducing the risk of confounding due to LD and improving the precision of the MR analysis.

Proxy SNPs and Harmonization: When target SNPs were missing in the outcome genome-wide association studies (GWAS), proxy SNPs with high correlation (R^2 > 0.8) were used as substitutes.²⁶ Additionally, to ensure consistency across datasets, SNPs with discordant alleles and palindromic SNPs (those that could be read forward or backward without change) were excluded. This harmonization process was crucial for maintaining data integrity, allowing for accurate and reliable MR analysis by ensuring that the genetic data were correctly aligned between the exposure and outcome datasets.

Validation of Instrumental Variables: The reliability of the instrumental variables was further validated by calculating F-statistics for each SNP. SNPs with F-statistics less than 10 were excluded from the MR analysis, as these would indicate weak instruments that could compromise the study’s conclusions.²⁷ By setting this threshold, the study ensured that only strong and reliable genetic instruments were used, thereby increasing the robustness of the causal estimates derived from the MR analysis.

MR Analysis

Primary Analysis Method: The inverse variance weighted (IVW) method was employed as the primary analysis tool for the MR study. The IVW method is particularly effective in the presence of heterogeneity among genetic variants and can provide unbiased estimates as long as all genetic instruments are valid.²⁷ This method served as the foundation for the causal inferences drawn from the genetic data.

Additional MR Methods: To enhance the reliability of the findings, the study also applied MR-Egger regression and median-based estimators. MR-Egger is less susceptible to directional pleiotropy, where genetic variants may influence the outcome through pathways other than the exposure of interest. Median-based estimators offer robustness to outliers, further strengthening the study’s conclusions by providing alternative causal estimates that can confirm the IVW results.²⁸

Pleiotropy and Sensitivity Tests: The study employed MR-PRESSO and Cochrane’s Q statistic to detect and correct for horizontal pleiotropy, where genetic variants might affect the outcome through multiple pathways, and to assess heterogeneity among the SNPs used.²⁹ Additionally, a leave-one-out sensitivity analysis was performed, where each SNP was removed one at a time to check the stability of the results.³⁰ These tests ensured that the MR findings were not driven by any single genetic variant or by pleiotropic effects, thereby validating the robustness of the results.³¹

Confounder and Outcome Control: To control for potential confounders, the SNPs (P<1×10⁶) in the GWAS Catalog website (https://www.ebi.ac.uk/gwas) were eliminated with potential confounders and outcomes, such as smoking, infections, autoimmune diseases, obesity (Supplementary Tables 3,4). This step was crucial to ensure that the observed associations were not influenced by these external variables. Furthermore, F-statistics were calculated for all SNPs used in the MR analysis to confirm that only strong instruments were included. This careful selection of SNPs helped to eliminate weak instruments, thereby increasing the reliability of the causal inferences drawn from the study.

Statistical Significance

All statistical analyses were performed using the Two-Sample MR package in R statistical software version 4.4.0. (R Foundation). Associations with P-values below 0.005 were regarded as strong evidence for a causal relationship, while those with P-values between 0.005 and 0.05 were considered suggestive of causality. This thresholding approach allowed for a nuanced interpretation of the data, distinguishing between more and less definitive evidence. All statistical analyses were conducted using the TwoSampleMR package in R, a widely used tool in the field of Mendelian randomization.

Results

Effect of Bone Mineral Density at Different Sites on RA

We obtained 82, 3, 23, 21, and 52 SNPs from GWAS corresponding to TB-BMD, FA-BMD, FN-BMD, LS-BMD, and Heel BMD (Figure 2). IVW estimates indicate a genetic causal relationship between lower bone density and increased RA risk (TB-BMD: OR: 1.094, 95% CI=1.023–1.170, P=0.009; FA-BMD: OR: 1.159, 95% CI=1.019–1.320, P=0.025; LS-BMD: OR: 1.175, 95% CI=1.070–1.291, P=0.001; Figure 2, Supplementary Table 1). MR-Egger, weighted median, and weighted mode methods yielded similar results (Figure 2, Supplementary Table 1). All Cochran’s Q tests showed no heterogeneity (p>0.05, Figure 2). All Egger regression tests were negative (p>0.05; Figure 2), indicating that our MR results were not influenced by horizontal pleiotropy. Leave-one-out sensitivity tests showed no single SNP had a potential impact on the final results Supplementary Figure 1-1 depicts scatter plots, funnel plots, and leave-one-out analyses of BMD on RA.

Figure 2 Effect of Bone Mineral Density at Different Sites on RA. Abbreviations: .RA, Rheumatoid arthritis; FN-BMD, femoral neck bone mineral density; LS-BMD, lumbar spine bone mineral density; TB-BMD, total body bone mineral density; FA-BMD, forearm bone mineral density; Heel-BMD, heel bone mineral density; CI, confidence interval, inverse-variance-weighted; PRESSO, Pleiotropy Residual Sum and Outlier; OR, Odds ratios; CI, confidence interval.

To further strengthen our MR assumption, we examined the traits related to our instrumental SNPs. Traits association analysis (see Supplementary Table 3) showed that rs4846580,rs633995 in TB-BMD, rs7524102 in LS-BMD,rs2941584,rs3802858,rs7045925 in Heel-BMD,and rs7524102 in FN-BMD are associated with smoking,BMI,ulcerative colitis, which may have some effect on the risk of RA. Sensitivity analysis by removing the SNPs revealed similar results, though with lower statistical power (see Supplementary Figure 3-1).

Effect of Bone Mineral Density at Different Ages on RA

We obtained 8, 1, 9, 20, and 21 SNPs from GWAS for ages 0–15 years, 15–30 years, 30–45 years, 45–60 years, and over 60 years, respectively (Figure 3, Supplementary Table 1). IVW estimates indicate a genetic causal relationship between lower bone density and increased RA risk (0–15 BMD: OR: 1.216, 95% CI=1.019–1.451, P=0.03; 45–60 BMD: OR: 1.110, 95% CI=1.027–1.200, P=0.009; Figure 3). MR-Egger, weighted median, and weighted mode methods yielded similar results (Figure 3). All Egger regression tests were negative (p>0.05; Figure 3), indicating that our MR results were not influenced by horizontal pleiotropy. Cochran’s Q tests indicated heterogeneity in the 15 or less age group (IVW=17.56, P=0.01; MR-Egger=16.24, P=0.01; Figure 3, Supplementary Table 1). Therefore, we conducted MR-PRESSO to further examine this relationship, but no outliers were identified (Figure 3, Supplementary Table 1). Leave-one-out sensitivity tests showed no single SNP had a potential impact on the final results (Supplementary Figure 1-2) depicts scatter plots, funnel plots, and leave-one-out analyses of BMD on RA.

Figure 3 Effect of Bone Mineral Density at Different Ages on RA. Abbreviations: RA, Rheumatoid arthritis;CI, confidence interval, inverse-variance-weighted; PRESSO, Pleiotropy Residual Sum and Outlier; OR, Odds ratios; CI, confidence interval.

To further strengthen our MR assumption, we examined the traits related to our instrumental SNPs. Traits association analysis (see Supplementary Table 3) showed that rs1936792 in 60 or more-BMD is associated with BMI, which may have some effect on the risk of RA. Sensitivity analysis by removing the SNPs revealed similar results, though with lower statistical power (see Supplementary Figure 3-2).

Effect of RA on Bone Mineral Density at Different Sites

In our study, we used 24 independent SNPs as instrumental variables (IVs) for RA (Figure 4; Supplementary Table 2). The analysis indicated no significant relationship between genetically predicted RA and decreased bone density. The primary IVW results showed no statistical association between RA risk and lower bone density levels (Figure 4). This was consistent with MR-Egger regression and median-based estimators (weighted median and weighted mode) (Figure 4). Tests for horizontal pleiotropy showed no directional pleiotropy (Figure 4). The weighted median analysis indicated no causal effect of RA on bone density at different sites (Figure 4). Heterogeneity tests revealed heterogeneity in TB-BMD, LS-BMD, and Heel-BMD (TB-BMD: IVW=42.07, P=6.13E-03, MR-Egger=42.00, P=4.20E-03; LS-BMD:IVW=38.45, P=0.02,MR-Egger=38.18, P=0.02; Heel-BMD: IVW=50.36,P=1.97E-4,MR-Egger=50.20,P=1.22E-4;Figure 4, Supplementary Table 2). Consequently, we performed MR-PRESSO for further testing, finding no significant outliers for TB-BMD, and excluding rs6679677 for LS-BMD and rs2618444 for Heel-BMD. After removing these outliers, no causal relationship between RA and osteoporosis was found (Figure 4; Supplementary Table 2). Leave-one-out sensitivity tests showed no single SNP had a potential impact on the final results (Supplementary Figure 2-1) depicts scatter plots, funnel plots, and leave-one-out analyses of RA on bone density at different sites.

Figure 4 Effect of RA on Bone Mineral Density at Different Sites. Abbreviations: RA, Rheumatoid arthritis; FN-BMD, femoral neck bone mineral density; LS-BMD, lumbar spine bone mineral density; TB-BMD, total body bone mineral density; FA-BMD, forearm bone mineral density; Heel-BMD, heel bone mineral density; CI, confidence interval, inverse-variance-weighted; PRESSO, Pleiotropy Residual Sum and Outlier; OR, Odds ratios; CI, confidence interval.

To further strengthen our MR assumption, we examined the traits related to our instrumental SNPs. The results of the trait association analysis (see Supplementary Table 4) indicated that the SNPs rs35139284, was associated with certain autoimmune conditions and several potential confounders, such as BMI, which may have some impact on OP. Sensitivity analysis by removing the SNPs revealed similar results (see Supplementary Figures 4-1).

Effect of RA on Bone Mineral Density at Different Ages

We found no significant causal effect of RA on bone density at different ages, whether using IVW, MR-Egger, weighted median, or weighted mode methods (Figure 5; Supplementary Table 2). Horizontal pleiotropy tests showed no directional pleiotropy (p>0.05; Figure 5; Supplementary Table 2). Heterogeneity tests indicated heterogeneity in the 30–45 and 45–60 age groups (age 30–45 BMD: IVW=37.42, P=0.03, MR-Egger=36.87, P=0.02;age 45–60 BMD: IVW=46.37, P=2.70E-3, MR-Egger=42.00, P=6.25E-3; Figure 5; Supplementary Table 2). Therefore, we performed MR-PRESSO, excluding rs3757387 for the 30–45 age group and rs35139284, rs6679677, and rs9494894 for the 45–60 age group. After removing these outliers, no causal relationship between RA and OP was found (P>0.05, Figure 5; Supplementary Table 2). IVW leave-one-out analysis indicated that most identified relationships were not altered by a single SNP associated with BMD. However, rs6679677 (pval.outcome 2.14E-4) may interfere with the causal relationship between BMD and RA risk in those over 60, so it was excluded, and reanalysis still showed no causal relationship between RA and OP (P>0.05, Figure 5;Supplementary Table 2). Leave-one-out sensitivity tests showed no single SNP had a potential impact on the final results (Supplementary Figure 2−2). Supplementary Figure 2–2 depicts scatter plots, funnel plots, and leave-one-out analyses of RA on bone density at different ages.

Figure 5 Effect of RA on Bone Mineral Density at Different Ages. Abbreviations: RA, Rheumatoid arthritis; IVW, inverse-variance-weighted; PRESSO, Pleiotropy Residual Sum and Outlier; OR, Odds ratios; CI, confidence interval.

To further strengthen our MR assumption, we examined the traits related to our instrumental SNPs. The results of the trait association analysis (see Supplementary Table 4) indicated that the SNPs rs35139284, was associated with certain autoimmune conditions and several potential confounders, such as BMI, which may have some impact on OP. Sensitivity analysis by removing the SNPs revealed similar results (see Supplementary Figures 4-2).

Discussion

The relationship between RA and OP has been a topic of interest in the medical community, with various hypotheses proposed to explain their connection. Most current studies focus on the unidirectional causal relationship where RA leads to osteoporosis, with little research on whether OP can lead to RA. In this study, we used bidirectional Mendelian randomization to determine whether genetically predicted osteoporosis has a causal relationship with RA, and vice versa. Using the largest publicly available GWAS summary data, we found a causal relationship between genetically decreased bone density and increased RA risk. However, we did not find a causal relationship between genetically increased RA risk and decreased BMD/OP. Sensitivity analyses confirmed the reliability of our results. Our study presents a novel finding: OP can lead to the development of RA.

Clinical research supports the notion that OP can lead to RA. Arnd Kleyer et al hypothesized that bone loss might begin before the onset of arthritis, and patients with the most severe bone pathology might be more prone to developing RA.³² Early and recent clinical observations consistently show that both localized and systemic bone loss is detected early in the RA disease process.³³ This observation is surprising because synovitis would take some time to destroy bone to a clinically detectable level. Van der Heijde et al has shown that most patients exhibit radiographic bone erosion within the first few years of the disease.³⁴ Our study’s findings are supported by clinical research suggesting that osteoporosis can lead to the development of RA, with bone loss potentially occurring before the onset of arthritis.^35–38 Previous studies have demonstrated that both localized and systemic bone loss can occur early in RA, even before significant joint inflammation is evident.^39,40 These observations align with the study’s conclusion that osteoporosis may increase the risk of RA, highlighting the importance of bone health in preventing or mitigating the onset of arthritis.

At the genetic level, our research indicates that the development of OP can lead to joint damage and destruction. Recently, two BMD-related Wnt16 gene SNPs (rs2707466 and rs2908004) were found to be associated with hip and knee osteoarthritis phenotypes in Caucasian patients in a gender-dependent manner, which may impact RA development.^41,42 Over the past decade, Wnt16, a member of the Wnt family, has been widely studied for its close association with bone mineral density, cortical bone thickness, bone strength, and osteoporosis fracture risk.⁴³ Notably, local Wnt16 treatment has been shown to alleviate osteoarthritis, inhibit bone resorption, and promote new bone formation in bone defect models.^43–45 Wnt16 is now a potential therapeutic target for skeletal diseases and osteoarthritis.⁴⁶ Our study underscores the significance of genetic factors, particularly the Wnt16 gene variants, in the relationship between osteoporosis and RA. Wnt16 gene variants are known to be associated with bone mineral density and osteoarthritis, suggesting that osteoporosis might contribute to the development of RA through genetic mechanisms involving Wnt16. This genetic link provides a potential pathway by which bone health can influence joint health and the progression of rheumatoid arthritis. In recent years, an increasing number of studies have highlighted the RANK/RANKL/OPG axis as an important pathway for regulating bone remodeling. Dysfunction of this pathway leads to an imbalance between bone formation and resorption, resulting in osteolytic lesions and structural abnormalities. The RANK/RANKL/OPG signaling pathway is a key regulator of several critical aspects of osteoarthritis, including cartilage degradation, abnormal bone remodeling, and synovial inflammation, playing an essential role in the pathogenesis and progression of osteoarthritis.⁴⁷ This further supports the possibility that osteoporosis may have an impact on the development of rheumatoid arthritis (RA).

While our study found no direct genetic causality from RA to osteoporosis, it proposes that RA may lead to osteoporosis through secondary effects. RA patients may suffer significant pain, stiffness, and loss of mobility, which undoubtedly leads to disuse osteoporosis.⁴⁷ In addition to mechanical factors, RA patients may develop medication-induced osteoporosis, especially due to long-term use of nonsteroidal anti-inflammatory drugs and glucocorticoids, which cause gastrointestinal side effects leading to reduced intake and impaired nutrient absorption.⁴⁸ This highlights the need for a comprehensive approach to managing RA that considers its potential impact on bone health.

However, this work is not perfect with some limitations to be mentioned here. First, we did not perform a stratified analysis of the causal effect of gender on the association between osteoporosis and RA. Moreover, the study’s focus on a population of European descent may limit the applicability of its findings to other ethnic groups, underscoring the need for further research in more diverse populations. Third, the age-based BMD groups, where the sample size for certain SNPs may be small, potentially reducing statistical power.

In conclusion, Our study suggests a potential causal relationship between osteoporosis and RA, suggesting that osteoporosis may predispose individuals to RA. Further research is needed to understand the mechanisms and to confirm these findings across diverse populations.