Skip to main content

Osteoporosis and osteoarthritis: a bi-directional Mendelian randomization study

Abstract

Objective

To investigate the causal relationship between low bone mineral density (BMD) and osteoarthritis (OA) using Mendelian randomization (MR) design.

Methods

Two-sample bi-directional MR analyses were performed using summary-level information on OA traits from UK Biobank and arcOGEN. Sensitivity analyses including MR-Egger, simple median, weighted median, MR pleiotropy residual sum, and outlier approaches were utilized in conjunction with inverse variance weighting (IVW). Gene ontology (GO) enrichment analyses and expression quantitative trait locus (eQTL) colocalization analyses were used to investigate the potential mechanism and shared genes between osteoporosis (OP) and OA.

Results

The IVW method revealed that genetically predicted low femoral neck BMD was significantly linked with hip (β = 0.105, 95% CI: 0.023–0.188) and knee OA (β = 0.117, 95% CI: 0.049–0.184), but not with other site-specific OA. Genetically predicted low lumber spine BMD was significantly associated with OA at any sites (β = 0.048, 95% CI: 0.011–0.085), knee OA (β = 0.101, 95% CI: 0.045–0.156), and hip OA (β = 0.150, 95% CI: 0.077–0.224). Only hip OA was significantly linked with genetically predicted reduced total bone BMD (β = 0.092, 95% CI: 0.010–0.174). In the reverse MR analyses, no evidence for a causal effect of OA on BMD was found. GO enrichment analysis and eQTL analysis illustrated that DDN and SMAD-3 were the most prominent co-located genes.

Conclusions

These findings suggested that OP may be causally linked to an increased risk of OA, indicating that measures to raise BMD may be effective in preventing OA. More research is required to determine the underlying processes via which OP causes OA.

Introduction

Osteoporosis (OP) is a common metabolic skeletal disease characterized by reduced bone mineral density (BMD), resulting in an increment of bone fragility among the elderly [1]. It causes considerable emotional, physical, and financial stresses to patients, and often leads to disability and poor quality of life. Osteoarthritis (OA) is the most common joint disease in which degenerative changes in joint cartilage cause aseptic inflammation and involve the entire joint tissues [2]. Its symptoms include joint pain, stiffness, and activity limitation. The etiology and pathological changes in OA remain largely unclear, resulting in the notable absence of curable therapies.

Does OA directly affect OP, or vice versa? Numerous observational studies have shown a strong link between OP and OA. Some have shown that patients with OA had a lower BMD, but others have reported that the progression of OA was accompanied by an increase in BMD [3]. There is evidence of greater systemic BMD in OA patients at several joint locations, after adjustment for spinal osteophytes [4]. A longitudinal study, however, has concluded that the development of pre-existing OA may be negatively correlated with BMD. Increased age-specific femoral neck BMD quartiles were linked to a lower risk of knee OA development throughout the course of the 8-year trial, according to research by Zhang et al. in the Framingham population [5]. Additionally, Hart et al. noted a tendency for lower hip BMD in those with progressing knee OA compared to non-progressors in the Chingford Study [6]. Although many studies have attempted to explore this controversy, previous observational data are limited for causal inference due to potential biases introduced by confounders. Mendelian randomization (MR), an emerging approach to evaluate causal links, has been popular in recent years. MR may reduce confounding effects and eliminate the bias of reverse causation since genotypes are independent of postnatal lifestyle and environmental variables and predate the illness process [7].

Pleiotropy may exist between OA-influencing variables such as body mass index (BMI) and genetic variables linked to BMD. Using a Mendelian randomization approach, April Hartley et al. found causal effects of hip and knee OA on BMD independent of BMI, while Liu Lin et al. found that OP reduced the incidence of OA (knee OA and hip OA). However, it is still unknown if there are site-specific causal links between low BMD and OA, and the potential signaling pathways involved are not clearly defined. Therefore, the aims of our study were to assess the causal relationship between OP and OA in the framework of a two-sample bi-directional MR study in which multiple OP measures (FN-BMD, LS-BMD, and TB-BMD) and OA measures (OA at any sites, knee OA and hip OA) will be used, and to further explore the potential signaling pathways and shared genes between OP and OA were identified.

Methods

Study design and data sources

The causal link between low BMD at each location and OA with various symptoms was investigated in a two-sample bi-directional MR investigation (Fig. 1). GWAS summary statistics of low BMD traits were extracted from publicly available data [GWAS catalog https://www.ebi.ac.uk/gwas/publications/26367794 for femoral neck (FN) [8] and lumber spine (LS) [8], and GWAS catalog https://www.ebi.ac.uk/gwas/publications/29304378 for total body (TB) [9]]. The summary statistics of OA were obtained from the most recent version from GWAS of UK Biobank data [10] (GWAS catalog https://www.ebi.ac.uk/gwas/publications/30664745 and https://msk.hugeamp.org/downloads.html) [11]. The latter summary statistics were used for sensitivity analysis. All data sources were publicly available, and thus, no ethical approval was required.

Fig. 1
figure 1

Workflow of bi-directional MR analysis. A, the fundamental idea of MR analysis: If we cannot randomize the exposure, we can find a randomized instrumental variable to disentangle. Instrumental variables were highly correlated with exposure factors. Instrumental variables were independent of the outcome. Instrumental variables were not correlated with confounding factors. B, Workflow of our bi-directional MR analysis. MR: Mendelian randomization; BMD: bone mineral density; OA: osteoarthritis; FN: femoral neck; LS: lumbar spine; TB: total body

Instrumental variables for BMD and OA

OP is defined clinically through the measurement of BMD (T score <  − 2.5), which remains the single best predictor of fracture [12]. In this research, BMD was used to characterize three phenotypes of OP measurements at various locations, including FN (n = 32,735), LS (n = 28,498), and systemic TB (n = 66,628). For the three BMD phenotypes, genetic differences in GWAS for poor BMD (− 2.5 T score as cutoff) were employed as instrumental factors. Based on two recent GWAS, we selected instrumental variables for FN-BMD, LS-BMD, and TB-BMD in our main analysis as independent genetic variants at the study-specific genome-wide significance level (P = 5 × 10−8) [8, 9]. The Genetic Factors for Osteoporosis (GEFOS) that published a meta-analysis of 32,735 FN-BMD and 28,498 LS-BMD in the European population in 2012 were included [8]. The other was the largest meta-analysis to date on FN- and LS-BMD, including 17 GWAS and 32,961 individuals of European descent [13]. Three BMD (per SD) phenotypes were adjusted for age, sex, weight, and height in the previous GWAS studies. Summary statistics for a GWAS meta-analysis research comprising 66,628 European participants were included in the GWAS dataset for TB-BMD [9]. From the genomic-wide meta-analysis spanning the UK Biobank and arcOGEN with 384,838 Europeans, summary-level data for three OA phenotypes, including OA at any sites, knee OA, and hip OA, were collected [10]. Sensitivity analysis was performed using the latest OA databases released in 2021, which conducted a GWAS meta-analysis containing 826,690 individuals (177,517 with OA). They discovered 100 risk mutations across 11 OA phenotypes that were independently linked with the illness, 52 of which were previously unrelated to the condition [11].

Genetic variations linked to low BMD were considered as instrumental SNPs in the forward stage (FN-BMD, LS-MD [8], and TB-BMD [9]). We used linkage disequilibrium (LD) [14] clustering based on r2 > 0.001 and removed variations within 1 Mb distance from other SNPs with a greater connection to fulfill the independence across instrumental SNPs for each exposure. Additionally, we harmonized the impact of these instrumental SNPs where possible and eliminated those that were absent from the GWAS of the outcomes to make sure that all associated risk factors and result alleles were on the same strand. Then, SNPs associated with low BMD were selected (FN-BMD 47, LS-BMD 47, TB-BMD 43) (Table S1). In the reverse stage, genetic variants associated with OA were used as instrumental SNPs [10]. Similarly, we used LD clumping based on r2 > 0.001 and removed variations within 1 Mb distance from other SNPs with a strong correlation to guarantee independence across instrumental SNPs for each exposure. To ensure that all corresponding risk factors and outcome alleles were on the same strand, we harmonized the effect of these instrumental SNPs where possible, and those not present in the GWAS of the outcomes were removed. At the same time, SNPs that might have horizontal pleiotropy were removed and the frequency of efficient alleles was used to ensure that palindrome instruments were correctly aligned where possible. Finally, SNPs associated with OA were selected (OA at any sites 27, knee OA 10, hip OA 26) (Table S2).

Statistical analysis

MR analysis

The causal relationship between each exposure and each outcome was assessed using the inverse variance weighted (IVW) approach with a fix effect model. We excluded instrumental variables (IVs) that were substantially linked with outcome. The IVW approach was often regarded as the most trustworthy indicator in MR analysis when there was no sign of directional pleiotropy (P for MR-Egger intercept > 0.05) [15]. The IVW method used the log (OR)/β coefficient of the disease SNP divided by the log (OR)/β coefficient of the exposed SNP to obtain the Wald estimates for each SNP, and then combined these Wald estimates using a method similar to meta-analysis [16]. When each genetic variation satisfied the IV hypothesis, the IVW method could provide a consistent estimate of the causal effect of exposure on the outcome. Cochran’s Q statistics were used to evaluate the IV heterogeneity. In order to further confirm the reliability of MR estimations, we used the MR-Pleiotropy Residual Sum and Outlier techniques (MR-PRESSO), which can identify and eliminate probable pleiotropic IVs and offer the outlier-adjusted estimates, to the IVs for IVW analysis (P-value < 0.05). In addition, we used weighted median, simple median, and MR-Egger as complementary analysis methods to test the robustness of the IVW method using random-effect model estimation. Utilizing aggregate data, the weighted median estimate was produced to give protection against ineffective instruments and reliable estimations of causation if at least 50% of the weight originates from IVs [17].

Subgroup analysis of age

Considering IVs of TB-BMD are from GWAS analyses of all age groups which may lead to a decrease in the accuracy of MR analyses [9], we did a detailed analysis by age groups. We used the IVW method of the random-effects model to evaluate the causal effects of each exposure and each outcome. TB-BMD data on the basis of age was divided into four stages, which were 15 to 30 (n = 4,180), 30 to 45 (n = 10,062), 45 to 60 years old (n = 18,805), and over 60 years old (n = 22,504), respectively [9]. MR analysis method was used to estimate the correlation between OA and the corresponding phase TB-BMD.

Gene ontology enrichment analysis

The gene ontology (GO) [18] which involves biological process, molecular function, and cellular component were analyzed to further validate whether the potential targets were indeed related to OA. We connected the causative BMD lead SNPs identified in several OA symptoms to the adjacent genes. After using the Benjamini–Hochberg method for adjustment [19], P < 0.05 was used as the significance cutoff in our study.

Expression quantitative trait locus analysis

We tested the relationships between OA and OP features using a gene-based strategy and extrapolated GTEx V746 expression levels for the musculoskeletal tissues. The elastic net model is used by MetaXcan [20] to impute the cis-genetic component of expression into a much bigger set for a collection of reference participants for whom both gene expression and genetic variation have been assessed. It then did a transcriptome-wide association study to identify meaningful expression-trait connections and corroborated the imputed gene expression to the trait of interest. With 20,000 genes spread over 48 tissues, we employed a conservative Bonferroni correction to account for the gene-tissue pairings, which resulted in a significance threshold of 5.20 × 10−8. We calculated the likelihood of each GWAS and expression quantitative trait locus (eQTL) signal colocalizing in each significant MetaXcan results using Coloc47 to lessen the impact of LD confounding on the MetaXcan results when various causal SNPs were affecting expression levels and phenotypic traits in a GWAS [21].

R 3.5.3 and STATA 14.0 were used to conduct all statistical analyses. Unless otherwise stated, statistical significance was defined as P < 0.05. Because the objective of our study was the causal relationship between different types of BMD and OA, we did not apply a multiple testing correction to the reported P. GO enrichment analysis was performed by websites tool “FUMA” at http://fuma.ctglab.nl.

Results

Stage 1: forward MR analysis of the effects of BMD on OA

By applying a two-sample MR approach, we first investigated the causal associations of low BMD on OA. We identified 47 (FN-BMD), 47 (LS-BMD), and 43 (TB-BMD) IVs from GWAS that met genome-wide significance level (P < 5 × 10−8) based on the elimination of certain missing data (r2 < 0.001). The heterogeneity test showed no significant heterogeneity among selected IVs (Q_P > 0.05, Table 1), indicating our MR results were not biased by heterogeneity or horizontal pleiotropy. Table S1 shows the effect sizes of a few IVs.

Table 1 Mendelian randomization estimates for bone mineral density on osteoarthritis

The IVW method showed genetically predicted low FN-BMD were positively associated with knee OA (β = 0.117, 95% CI: 0.049–0.184) and hip OA (β = 0.105, 95% CI: 0.023–0.188), but not with OA at any sites (β = 0.029, 95% CI: − 0.014–0.071). Genetically predicted low LS-BMD was positively associated with OA at any sites (β = 0.048, 95% CI: 0.011–0.085), knee OA (β = 0.101, 95% CI: 0.045–0.156) and hip OA (β = 0.150, 95% CI: 0.077–0.224). Genetically predicted low TB-BMD was positively associated with hip OA (β = 0.092, 95% CI: 0.010–0.174), but not with OA at any sites (β = 0.025, 95% CI: − 0.020–0.070), or knee OA (β = 0.043, 95% CI: − 0.023–0.109). Causal estimates from the weighted median and simple median revealed similar findings of FN-BMD on hip OA and OA at any sites: LS-BMD on knee OA, hip OA, and OA at any sites (β ranged from 0.065 to 0.174). The weighted median showed causal association between low FN-BMD and knee OA (β = 0.100, 95% CI: 0.003–0.196). MR-Egger showed no causal effect of BMD on OA phenotypes. Detailed data are shown in Table 1. To test the authenticity of these analyses, we used the summary data on OA released in 2021 [11] for a secondary calculation, and the results were largely consistent (Table S2).

Stage 2: reverse MR analysis of the effects of OA on BMD

Using OA-associated SNPs as IVs, we further investigated the causative relationships between OA and low BMD in the reversal direction. We obtained 26 (OA at any sites), 10 (knee OA), and 26 (hip OA) IVs without effects of LD-independent (r2 < 0.001), reaching P < 1 × 10−5 from GWAS. When examining the causal relationship of OA at any sites on FN-BMD and TB-BMD, we deleted one of the SNPs (rs528981060) because of its significant heterogeneity on outcomes. The heterogeneity test showed no significant heterogeneity (Q_P > 0.05, Table 2) in selected IVs of OA on low BMD, demonstrating that neither horizontal pleiotropy nor heterogeneity affected our MR findings. The negative control analysis showed that OA was not associated with selected IVs, suggesting that the selected IVs for exposures in this study were appropriate (Table S3).

Table 2 Mendelian randomization estimates for osteoarthritis on bone mineral density

The IVW method showed that genetically predicted three OA exposures were not associated with three BMD outcomes (β ranged from 0.018 to 0.061). The estimates from MR-Egger, simple median, and weighted median were largely similar (β ranged from 0.025 to 0.181) (Table 2).

Subgroup analysis by age

Subgroup analyses by age were showed in Table S4. In those with 30–45 years of age, TB-BMD only had a strong causal effect on OA at any sites (P = 0.017) and hip OA (P = 0.021). TB-BMD had a strong causal effect on all the OA phenotypes in participants between 45 and 60 years of age (OA at any sites, P = 0.001; knee OA, P = 0.001; hip OA, P = 0.009). In those above 60 years of age, TB-BMD also had a strong causal effect on OA at any sites, knee OA and hip OA. However, we did not find valid data for MR estimates (data not shown) in participants below 30 years of age. Due to the absence of FN-BMD and LS-BMD subgroup data in different age groups, we were unable to complete the relevant MR analysis.

GO enrichment analysis

The findings of the GO enrichment study are shown in Fig. 2. BMD on different OA sites revealed high enrichment in a number of important regulatory pathways according to GO enrichment analyses. For knee OA, 10 GO biological processes were observed to be involved, such as ossification, cell signaling by WNT, and bone mineralization (Fig. 2A). For hip OA, 68 GO biological processes were observed to be involved and the most significant ones were presented in Fig. 2, such as cell signaling by WNT, cell surface receptor signaling pathway, ossification, skeletal system development and neurogenesis (Fig. 2B). For OA at any sites, 15 GO biological processes were observed to be involved, such as regulation of gene expression, RNA biosynthetic process, tissue morphogenesis and osteoblast (Fig. 2C).

Fig. 2
figure 2

Gene Ontology enrichment analysis of the casual effect of Bone Mineral Density on Osteoarthritis. A, Gene Ontology enrichment analysis of the casual effect of Bone Mineral Density on knee Osteoarthritis. B, Gene Ontology enrichment analysis of the casual effect of Bone Mineral Density on Hip Osteoarthritis. C, Gene Ontology enrichment analysis of the casual effect of Bone Mineral Density on Osteoarthritis at any sites

eQTL co-location analysis

The eQTL colocation analysis showed that four BMD susceptibility genes (ANAPC1, GALNT3, LIN7C, and DDN) had variants co-located on musculoskeletal tissues (Table 3). We found that gene GALNT3 was widely presented in all three types of BMD, whose variants co-located on LS-BMD were the most prominent (rs1346004, rs11680288). Among all susceptibility genes with variants co-located discovered, the effect of gene DDN was the largest (P < 0.001). We also found four hip OA susceptibility genes (LTBP3, MLXIP, SMAD3, and MAPT) that had variants co-located on musculoskeletal tissues. For hip OA co-located variants, the effect of gene SMAD3 was the largest (P < 0.001), followed by MAPT, LTBP3, and MLXIP (all P < 0.001). No genes had variants co-located on knee OA or OA at any sites.

Table 3 Significant eQTL signals for osteoporosis and osteoarthritis from musculoskeletal tissues based on the GTEx data

Discussion

Using the bi-directional two-sample MR method, we found evidence for a causal association of low BMD on OA, and the association was more prominent in people over 45 years old. No clear evidence of a causal relationship from OA to BMD was found. The eQTL co-location analysis revealed that four BMD susceptibility genes (ANAPC1, GALNT3, LIN7C, and DDN) and four hip OA susceptibility genes (LTBP3, MLXIP, SMAD3, and MAPT) had variants co-located on musculoskeletal tissues, suggesting these genes may play roles in regulating the development of OP and OA.

Previous studies reported that OP may be one of the etiologies to promote OA development. A meta-analysis showed that the worldwide incidence of OP increased with aging, and could promote OA through a number of mechanisms [22]. For example, a prospective study showed older participants with radiographic hip and knee OA had higher total hip bone loss over 2.6 years [23]. Developing OP has been shown to aggravate cartilage lesions in an experimental model of OA in rabbits [24]. Similar findings were made by Bellido M. et al., who discovered that the fragility and poor quality of subchondral bone may also increase cartilage injury in addition to promoting bone remodeling [25]. An increased rate of bone turnover has been shown in animal models of early OA [26]. In cadaver specimens of people with early OA, the elastic modulus of trabecular subchondral bone from proximal tibiae decreased despite an increase in bone volume [27]. This local bone softening may accompany a drop in BMD, which may be the result of inadequate mineralization brought on by increased bone remodeling [27]. In summary, our findings are generally consistent with the above observational studies that low BMD may be associated with an increased risk of OA. In contrast, there were some observational studies showing that higher femoral neck and total body BMD were associated with an increased risk of incident OA [28]. A cross-sectional and longitudinal study reported a positive correlation between subchondral BMD and OA severity in patients [29]. In general, these studies indicated the detrimental effect of higher BMD on OA development.

Regarding the possible causal relationship from OA to OP, results varied among studies. The MR Study by Hartley et al. suggested that a reduction in BMD would increase the incidence of OA, while Lin et al. suggested that a reduction in bone mineral density had a protective effect on OA [30, 31]. The possible explanations for the discrepancies include (1) different sources of exposure and differences in the final selection of SNPs, (2) although we tried to minimize phenotypic heterogeneity, it could not disappear, (3) different methods for assessing OA and/or OP [27], and (4) lack of adjustment for essential confounders such as age, gender, weight, and/or BMI [32]. The protective effect of high BMD on OA may be through its effect on reducing the risk of joint space loss [5], and the mechanism of OP causing OA could be the bone loss in the subchondral bone resulting in the collapse of the articular surface and leading to uneven stress on the articular cartilage, which leads to secondary osteophyte proliferation and cartilage damage [33]. It is worth noting that Lin et al. [31] only analyzed the causal relationship from OP to OA, and the summary data of exposure and outcome in their study were both from UK Biobank, which may cause bias in the results. Based on two separate OA summary datasets, our findings largely agreed with Hartley et al.’s [30], indicating that improving BMD might help to prevent OA from a genetic standpoint. Kindly noted the previously study assessed BMD by heel ultrasound as the only exposure, while our study used dual-energy X-ray assessed FN-BMD, LS-BMD, and TB-BMD as the exposures.

In our study, we found that the causal effect of FN-BMD on knee OA was stronger than that of hip OA, while LS-BMD had a strong causal effect on OA at any sites, knee OA and hip OA. Hackinger et al. [34] discovered a genetic link between LS-BMD (but not FN) and OA, which is consistent with common biological pathways causing both BMD and OA. On the other hand, TB-BMD only had a weak causal effect on hip OA, possibly because of the special trabecular structure of the hip [35]. Moreover, although we used two databases for verification, we cannot rule out the possibility of false positives. In the reverse MR analyses, no evidence for a causal effect of OA on BMD was found. This could be due to the influence of the genetic variation on the result is not entirely mediated by the binary exposure (OA definition), power estimations are probably conservative in this case. However, estimating approaches for a binary exposure call for strong assumptions that are unlikely to be physiologically tenable in typical MR situations, making it difficult to perform tests for causal effects without employing exposure information [36].

GO enrichment analysis found that a large number of GO biologic processes play key roles in the potential causal relationship from low BMD to OA, such as ossification [37], cell surface receptor signaling pathway [38], and cell signaling by WNT [39]. Using OP mice mode, Wu et al. discovered that nine cellular elements, including the cytosol, nucleus, cytoplasm, neuronal cell body, protein complex, caveola, and endoplasmic reticulum, were involved in the anti-OP effects on OA [37]. In addition, our GO analysis revealed that the primary ways in which OA affected OP were via inflammatory response, aging, responses to estradiol, glucocorticoids, and hydrogen peroxide [37]. The study by Guan et al. suggested that XianLing GuBao Capsule (XLGB) might treat OP through the PI3K/AKT and MAPK pathways [38]. Similarly, Xiao et al. found that PI3K/AKT/NF-κB and MAPK pathways would be the potential mechanism for XLGB in the treatment of OA [39]. Additionally, aberrant bone remodeling and reduced mineralization have been linked to abnormal OPG and RANKL expression in OA osteoblasts [40]. WNT protein is considered as a risk factor for the onset and development of OA. There are two primary types of WNT signaling pathways: canonical and non-canonical. Of these, the canonical WNT signaling pathway promotes osteogenesis. Osteocytes generate the inhibitor of this process, sclerostin, which prevents osteogenesis [41]. Since we substituted BMD SNPs for OP, it is reasonable to assume that OP also has these biological processes and WNT signaling could be potential pathways to be targeted for treating both diseases.

Through eQTL colocation analysis [21], 3 SNPs (rs1346004, rs11680288, and rs7586085) were co-located at gene GALNT3. An observational study found that the GALNT3 gene significantly correlated with OP and the low expression of the GALNT3 gene can promote the occurrence and deterioration of OP [42]. Our data showed that three newly discovered genes had eQTL effects with OP (ANAPC1, LIN7C, and DDN). In addition, previous reports deemed that the LIN7C gene significantly correlates with knee OA [43]. Our study also found four genes with eQTL effect on hip OA (LTBP3, MLXIP, SMAD3, and MAPT). The LTBP3 gene may be involved in the occurrence and development of OA by regulating the TGF-β signaling pathway [44]. SMAD3 gene is the target gene of circ-RNA (Circ0083429), which regulates the mRNA level of SMAD3 through the sponging of miRNA-346, thus affecting OA. Studies have found that the SMAD3 gene plays an important role in cartilage bone reconstruction and maintenance, and have confirmed the differential expression of SMAD3 in intact and degraded knee and hip cartilage [34]. Top-level differentially expressed genes in OA bone, such as those in the WNT signaling pathway (TWIST1, IBSP, S100A4, MMP25, RUNX2, and CD14) and TGF-β/SMAD3 signaling pathway, are known to function in osteoblasts, osteocytes, and osteoclasts (ADAMTS4, ADM, MEPE, GADD45B, COL4A1, and FST) [45]. MLXIP and MAPT have not been reported to be related to the development of hip OA. Therefore, these newly discovered genes need further confirmation.

The strengths of the current study lie on: First, to minimize possible false-positive results, summary data from two different organizations were used for outcomes, and the latest summary data of OA were also used for verification [10, 11], proving the credibility of our conclusions. Second, four methods of MR analysis on 18 groups of BMD GWAS and OA GWAS summary data were performed, and subgroup analyses by age were performed for TB-BMD. Third, we used GO Enrichment analysis and eQTL to search for potential enrichment genes and to determine whether there are possible co-acting pathways.

Limitations of our study include: First, the majority of the GWAS summary data stem from people of European origin, therefore, we should be cautious when applying our findings to people with other racial and ethnic backgrounds. Second, we evaluated a linear correlation between OP and OA in our MR analysis, but we did not account for the potential of different shapes of association [15]. Third, because summary statistics rather than raw data were used in the analysis, it was not possible to perform subgroup analyses, such as stratified by sex or ethnicity. Last, the results of GO and eQTL were based only on bioinformatic analyses, in vitro and in vivo experiments are needed in the future to validate our findings.

Conclusions

These findings suggested that OP may be causally linked to an increased risk of OA, indicating that measures to raise BMD may be effective in preventing OA. More research is required to determine the underlying processes via which OP causes OA.

Availability of data and materials

All data comes from public databases.

References

  1. Sözen T, Özışık L, Başaran N. An overview and management of osteoporosis. Eur J Rheumatol. 2017;4(1):46–56. https://doi.org/10.5152/eurjrheum.2016.048.

    Article  PubMed  Google Scholar 

  2. Glyn-Jones S, Palmer AJ, Agricola R, et al. Osteoarthritis. Lancet. 2015;386(9991):376–87. https://doi.org/10.1016/s0140-6736(14)60802-3.

    Article  CAS  PubMed  Google Scholar 

  3. Teichtahl AJ, Wang Y, Wluka AE, et al. Associations between systemic bone mineral density and early knee cartilage changes in middle-aged adults without clinical knee disease: a prospective cohort study. Arthritis Res Ther. 2017;19(1):98. https://doi.org/10.1186/s13075-017-1314-0.

    Article  PubMed  PubMed Central  Google Scholar 

  4. Nevitt MC, Lane NE, Scott JC, et al. Radiographic osteoarthritis of the hip and bone mineral density. The Study of Osteoporotic Fractures Research Group. Arthritis Rheum. 1995;38(7):907–16. https://doi.org/10.1002/art.1780380706.

  5. Zhang Y, Hannan MT, Chaisson CE, et al. Bone mineral density and risk of incident and progressive radiographic knee osteoarthritis in women: the Framingham Study. J Rheumatol. 2000;27(4):1032–7.

    CAS  PubMed  Google Scholar 

  6. Hart DJ, Cronin C, Daniels M, Worthy T, Doyle DV, Spector TD. The relationship of bone density and fracture to incident and progressive radiographic osteoarthritis of the knee: the Chingford Study. Arthritis Rheum. 2002;46(1):92–9. https://doi.org/10.1002/1529-0131(200201)46:1%3c92::Aid-art10057%3e3.0.Co;2-#.

    Article  PubMed  Google Scholar 

  7. Fan J, Zhu J, Sun L, Li Y, Wang T, Li Y. Causal association of adipokines with osteoarthritis: a Mendelian randomization study. Rheumatology (Oxford). 2021;60(6):2808–15. https://doi.org/10.1093/rheumatology/keaa719.

    Article  CAS  PubMed  Google Scholar 

  8. Zheng HF, Forgetta V, Hsu YH, et al. Whole-genome sequencing identifies EN1 as a determinant of bone density and fracture. Nature. 2015;526(7571):112–7. https://doi.org/10.1038/nature14878.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  9. Medina-Gomez C, Kemp JP, Trajanoska K, et al. Life-Course Genome-wide Association Study Meta-analysis of Total Body BMD and Assessment of Age-Specific Effects. Am J Hum Genet. 2018;102(1):88–102. https://doi.org/10.1016/j.ajhg.2017.12.005.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  10. Zengini E, Hatzikotoulas K, Tachmazidou I, et al. Genome-wide analyses using UK Biobank data provide insights into the genetic architecture of osteoarthritis. Nat Genet. 2018;50(4):549–58. https://doi.org/10.1038/s41588-018-0079-y.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  11. Boer CG, Hatzikotoulas K, Southam L, et al. Deciphering osteoarthritis genetics across 826,690 individuals from 9 populations. Cell. 2021;184(24):6003–5. https://doi.org/10.1016/j.cell.2021.11.003.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  12. Johnell O, Kanis JA, Oden A, et al. Predictive value of BMD for hip and other fractures. J Bone Miner Res. 2005;20(7):1185–94. https://doi.org/10.1359/jbmr.050304.

    Article  PubMed  Google Scholar 

  13. Estrada K, Styrkarsdottir U, Evangelou E, et al. Genome-wide meta-analysis identifies 56 bone mineral density loci and reveals 14 loci associated with risk of fracture. Nat Genet. 2012;44(5):491–501. https://doi.org/10.1038/ng.2249.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  14. Wang Y, Li T, Fu L, Yang S, Hu YQ. A Novel Method for Mendelian Randomization Analyses With Pleiotropy and Linkage Disequilibrium in Genetic Variants From Individual Data. Front Genet. 2021;12: 634394. https://doi.org/10.3389/fgene.2021.634394.

    Article  PubMed  PubMed Central  Google Scholar 

  15. Wu F, Huang Y, Hu J, Shao Z. Mendelian randomization study of telomere length and bone mineral density. Aging (Albany NY). 2020;13(2):2015–30. https://doi.org/10.18632/aging.202197.

    Article  PubMed  Google Scholar 

  16. Zheng J, Baird D, Borges MC, et al. Recent Developments in Mendelian Randomization Studies. Curr Epidemiol Rep. 2017;4(4):330–45. https://doi.org/10.1007/s40471-017-0128-6.

    Article  PubMed  PubMed Central  Google Scholar 

  17. Bowden J, Davey Smith G, Haycock PC, Burgess S. Consistent Estimation in Mendelian Randomization with Some Invalid Instruments Using a Weighted Median Estimator. Genet Epidemiol. 2016;40(4):304–14.

  18. Yu XH, Yang YQ, Cao RR, Bo L, Lei SF. The causal role of gut microbiota in development of osteoarthritis. Osteoarthritis Cartilage. 2021;29(12):1741–50. https://doi.org/10.1016/j.joca.2021.08.003.

    Article  PubMed  Google Scholar 

  19. Kraus VB, Feng S, Wang S, et al. Trabecular morphometry by fractal signature analysis is a novel marker of osteoarthritis progression. Arthritis Rheum. 2009;60(12):3711–22. https://doi.org/10.1002/art.25012.

    Article  PubMed  PubMed Central  Google Scholar 

  20. Gillies CE, Putler R, Menon R, et al. An eQTL Landscape of Kidney Tissue in Human Nephrotic Syndrome. Am J Hum Genet. 2018;103(2):232–44. https://doi.org/10.1016/j.ajhg.2018.07.004.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  21. Zuber V, Grinberg NF, Gill D, et al. Combining evidence from Mendelian randomization and colocalization: Review and comparison of approaches. Am J Hum Genet. 2022;109(5):767–82. https://doi.org/10.1016/j.ajhg.2022.04.001.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  22. Salari N, Ghasemi H, Mohammadi L, et al. The global prevalence of osteoporosis in the world: a comprehensive systematic review and meta-analysis. J Orthop Surg Res. 2021;16(1):609. https://doi.org/10.1186/s13018-021-02772-0.

  23. Ding C, Cicuttini F, Boon C, et al. Knee and hip radiographic osteoarthritis predict total hip bone loss in older adults: a prospective study. J Bone Miner Res. 2010;25(4):858–65. https://doi.org/10.1359/jbmr.091012.

    Article  PubMed  Google Scholar 

  24. Calvo E, Castañeda S, Largo R, Fernández-Valle ME, Rodríguez-Salvanés F, Herrero-Beaumont G. Osteoporosis increases the severity of cartilage damage in an experimental model of osteoarthritis in rabbits. Osteoarthritis Cartilage. 2007;15(1):69–77. https://doi.org/10.1016/j.joca.2006.06.006.

    Article  CAS  PubMed  Google Scholar 

  25. Bellido M, Lugo L, Roman-Blas JA, et al. Subchondral bone microstructural damage by increased remodelling aggravates experimental osteoarthritis preceded by osteoporosis. Arthritis Res Ther. 2010;12(4):R152. https://doi.org/10.1186/ar3103.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  26. Pastoureau PC, Chomel AC, Bonnet J. Evidence of early subchondral bone changes in the meniscectomized guinea pig. A densitometric study using dual-energy X-ray absorptiometry subregional analysis. Osteoarthritis Cartilage. 1999;7(5):466–73. https://doi.org/10.1053/joca.1999.0241.

  27. Day JS, Ding M, van der Linden JC, Hvid I, Sumner DR, Weinans H. A decreased subchondral trabecular bone tissue elastic modulus is associated with pre-arthritic cartilage damage. J Orthop Res. 2001;19(5):914–8. https://doi.org/10.1016/s0736-0266(01)00012-2.

    Article  CAS  PubMed  Google Scholar 

  28. Nevitt MC, Zhang Y, Javaid MK, et al. High systemic bone mineral density increases the risk of incident knee OA and joint space narrowing, but not radiographic progression of existing knee OA: the MOST study. Ann Rheum Dis. 2010;69(1):163–8. https://doi.org/10.1136/ard.2008.099531.

    Article  CAS  PubMed  Google Scholar 

  29. Cao Y, Stannus OP, Aitken D, et al. Cross-sectional and longitudinal associations between systemic, subchondral bone mineral density and knee cartilage thickness in older adults with or without radiographic osteoarthritis. Ann Rheum Dis. 2014;73(11):2003–9. https://doi.org/10.1136/annrheumdis-2013-203691.

    Article  PubMed  Google Scholar 

  30. Hartley A, Sanderson E, Granell R, et al. Using multivariable Mendelian randomization to estimate the causal effect of bone mineral density on osteoarthritis risk, independently of body mass index. Int J Epidemiol. 2021. https://doi.org/10.1093/ije/dyab251.

  31. Lin L, Luo P, Yang M, Wang J, Hou W, Xu P. Causal relationship between osteoporosis and osteoarthritis: A two-sample Mendelian randomized study. Front Endocrinol. 2022;13:1011246. https://doi.org/10.3389/fendo.2022.1011246.

    Article  Google Scholar 

  32. Emdin CA, Khera AV, Kathiresan S. Mendelian Randomization. JAMA. 2017;318(19):1925–6. https://doi.org/10.1001/jama.2017.17219.

    Article  PubMed  Google Scholar 

  33. Lajeunesse D, Reboul P. Subchondral bone in osteoarthritis: a biologic link with articular cartilage leading to abnormal remodeling. Curr Opin Rheumatol. 2003;15(5):628–33. https://doi.org/10.1097/00002281-200309000-00018.

    Article  PubMed  Google Scholar 

  34. Hackinger S, Trajanoska K, Styrkarsdottir U, et al. Evaluation of shared genetic aetiology between osteoarthritis and bone mineral density identifies SMAD3 as a novel osteoarthritis risk locus. Hum Mol Genet. 2017;26(19):3850–8. https://doi.org/10.1093/hmg/ddx285.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  35. Tamimi I, Cortes ARG, Sánchez-Siles JM, et al. Composition and characteristics of trabecular bone in osteoporosis and osteoarthritis. Bone. 2020;140: 115558. https://doi.org/10.1016/j.bone.2020.115558.

    Article  CAS  PubMed  Google Scholar 

  36. Burgess S, Labrecque JA. Mendelian randomization with a binary exposure variable: interpretation and presentation of causal estimates. Eur J Epidemiol. 2018;33(10):947–52. https://doi.org/10.1007/s10654-018-0424-6.

    Article  PubMed  PubMed Central  Google Scholar 

  37. Wu Y, Gao LJ, Fan YS, Chen Y, Li Q. Network Pharmacology-Based Analysis on the Action Mechanism of Oleanolic Acid to Alleviate Osteoporosis. ACS Omega. 2021;6(42):28410–20. https://doi.org/10.1021/acsomega.1c04825.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  38. Liu Y, Liu Q, Yin C, et al. Uncovering Hidden Mechanisms of Different Prescriptions Treatment for Osteoporosis via Novel Bioinformatics Model and Experiment Validation. Front Cell Dev Biol. 2022;10: 831894. https://doi.org/10.3389/fcell.2022.831894.

    Article  PubMed  PubMed Central  Google Scholar 

  39. Xiao J, Zhang G, Mai J, et al. Bioinformatics analysis combined with experimental validation to explore the mechanism of XianLing GuBao capsule against osteoarthritis. J Ethnopharmacol. 2022;294:115292. https://doi.org/10.1016/j.jep.2022.115292.

  40. Cafferata EA, Monasterio G, Castillo F, et al. Overexpression of MMPs, cytokines, and RANKL/OPG in temporomandibular joint osteoarthritis and their association with joint pain, mouth opening, and bone degeneration: A preliminary report. Oral Dis. 2021;27(4):970–80. https://doi.org/10.1111/odi.13623.

    Article  PubMed  Google Scholar 

  41. Maeda K, Kobayashi Y, Koide M, et al. The Regulation of Bone Metabolism and Disorders by Wnt Signaling. Int J Mol Sci. 2019;20(22). https://doi.org/10.3390/ijms20225525.

  42. Wang XY, Yang B, Liu CS, et al. Research on correlation between GALNT3 gene and osteoporosis. Eur Rev Med Pharmacol Sci. 2018;22(1 Suppl):69–75. https://doi.org/10.26355/eurrev_201807_15366.

    Article  CAS  PubMed  Google Scholar 

  43. Yerges-Armstrong LM, Yau MS, Liu Y, et al. Association analysis of BMD-associated SNPs with knee osteoarthritis. J Bone Miner Res. 2014;29(6):1373–9. https://doi.org/10.1002/jbmr.2160.

    Article  CAS  PubMed  Google Scholar 

  44. Zhao T, Zhao J, Ma C, Wei J, Wei B, Liu J. Common variants in LTBP3 gene contributed to the risk of hip osteoarthritis in Han Chinese population. Biosci Rep. 2020;40(6). https://doi.org/10.1042/bsr20192999.

  45. Hopwood B, Tsykin A, Findlay DM, Fazzalari NL. Microarray gene expression profiling of osteoarthritic bone suggests altered bone remodelling, WNT and transforming growth factor-beta/bone morphogenic protein signalling. Arthritis Res Ther. 2007;9(5):R100. https://doi.org/10.1186/ar2301.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

The authors thank all the participants and staff who made this study possible, especially participants who have contributed data to the UK Biobank and arcOGEN.

Patient and public involvement

Patients and/or the public were not involved in the design, conduct, reporting, or dissemination of this study.

Transparency

The lead authors (YQ/ZZ) affirm that the manuscript is an honest, accurate, and transparent account of the study being reported; that no important aspects of the study have been omitted; and that any discrepancies from the study as planned have been explained.

Funding

The present study was supported by the National Natural Science Foundation of China (32000925), Guangzhou Science and Technology Program (202002030481), Guangdong Basic and Applied Basic Research Foundation (2019A1515111169), and Wu Jieping Medical Foundation Program (3206750.2020–03-12).

Author information

Authors and Affiliations

Authors

Contributions

YQ, SC, MH, GR, YZ, ZG, YZ, CD, and ZZ designed the study. YQ, SC, MH, TF, MZ, PC, YZ, and ZZ managed and analyzed the data. YQ drafted the initial study protocol. YQ and MH wrote the first draft of the article. All authors contributed to the data interpretation and preparation of the report. All authors provided intellectual content to the manuscript and approved the submission. The corresponding author attests that all listed authors meet authorship criteria and that no others meeting the criteria have been omitted.

Corresponding authors

Correspondence to Changhai Ding, Yan Zhang or Zhaohua Zhu.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare no competing financial interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: Table S1.

Results of two-sample Mendelian randomization analyses of bone mineral density on osteoarthritis. Table S2. Mendelian randomization estimates for bone mineral density on osteoarthritis in alternative summary data. Table S3. Results of two-sample Mendelian randomization analyses of osteoarthritis on bone mineral density. Table S4. Subgroup analysis using Mendelian randomization estimates for bone mineral density on osteoarthritis by age.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Qu, Y., Chen, S., Han, M. et al. Osteoporosis and osteoarthritis: a bi-directional Mendelian randomization study. Arthritis Res Ther 25, 242 (2023). https://doi.org/10.1186/s13075-023-03213-5

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s13075-023-03213-5

Keywords