Introduction

Type 2 diabetes (T2D), a chronic metabolic disorder characterized by elevated blood glucose levels and increased risk of numerous serious and life-threatening complications, constitutes one of the biggest health problems in the world (Diamond, 2003). According to Global Burden of Disease (GBD) data, the age-standardized global prevalence of type 2 diabetes was approximately 6.0% in men and 5.0% in women in 2019 (Tinajero and Malik, 2021). It accounts for more than 100 billion dollars of healthcare costs annually in the United States (Diamond, 2003). The chronic comorbidities of T2D could develop gradually, and could lead to serious damage to heart, blood vessel, kidney, eye and foot (Teck, 2022). Other organ systems such as skeletal health could also be influenced by T2D (Lei and Kindler, 2022).

Osteoporosis is a common chronic disease characterized by low bone mass and disruption of bone microarchitecture. Fragility fracture is the ultimate outcome of poor bone health. Our previous studies have suggested that bone mass and fracture could be influenced by many modifiable or non-modifiable factors (Zhu and Zheng, 2021), such as body weight (Zhu et al., 2022), sleep behavior (Qian et al., 2021), inflammatory disease (Xia et al., 2020), birth weight (Xia et al., 2022) and genetic factors (Zhu et al., 2021). T2D is also considered to be a major factor that could affect bone health, it seems that T2D patients would have higher bone mineral density (BMD) and higher fracture risk than individuals without T2D (Khosla et al., 2021). This is the so called “diabetic bone paradox” (Botella Martinez et al., 2016; Romero-Diaz et al., 2021). For example, in an Italian nationwide study of 59,950 women of whom 5.2% had diabetes (predominantly type 2 diabetes), it noted an association between diabetes and any fracture (OR 1.3, 95% CI 1.1–1.4, and OR 1.3, 95% CI 1.2–1.5, for vertebral or hip fractures and non-vertebral, non-hip fractures, respectively) (Adami et al., 2020). Interestingly, the prevalence of vertebral or hip fracture was higher in participants with diabetes but without obesity (OR 1.9, 95% CI 1.7–2.1) than in patients with obesity and diabetes (OR 1.5, 95% CI 1.3–1.8), suggesting that obesity might be partially protective against vertebral or hip fractures in type 2 diabetes (Adami et al., 2020).

However, a recent comparative cohort analysis using routinely collected UK primary care records data from the Health Improvement Network (including 174,244 individuals with incident type 2 diabetes and 747,290 without diabetes) found no evidence to suggest a higher risk of fracture in type 2 diabetes patients, specifically, the risk of having at least one fracture was estimated to be 6% lower for females and 3% lower for males in the type 2 diabetes cohort than for females and males without diabetes (Davie et al., 2021). Lower fracture risk was also observed in the type 2 diabetes patients compared to those without the disease in the age group great than 85 years (Davie et al., 2021). Another large-scale cohort study showed that type 2 diabetes could only explain less than 0.1% of the fracture risk (Axelsson et al., 2023), and if the T2D patients with risk factors (such as low BMI, long diabetes duration, insulin treatment, and low physical activity) were excluded, T2D patients would have lower fracture risk than their matched controls (Axelsson et al., 2023). In a prospective study to examine the relationship between BMD and fracture in older adults with type 2 diabetes, it was reported that femoral neck BMD T score and FRAX score were both associated with fracture risk in individuals with type 2 diabetes, suggesting that BMD is still a useful clinical predictor for the evaluation of fracture risk in type 2 diabetes patients (Schwartz et al., 2011).

As the pathophysiology of fracture is more complicated than the BMD trait, and while there were some explanations for the “diabetic bone paradox” (Osorio, 2011), the integrated analyses with genetic data for the diseases could provide an alternative approach to alleviate the bias of the unknown confounding factors (Davey Smith and Hemani, 2014). Therefore, in this study, we firstly performed weighted genetic risk score (wGRS) regression analysis to assess relationship between the genetically predicted T2D and fracture with genetic summary data and individual genotype data in UK biobank. The two-sample mendelian randomization (MR) approach was used as an independent validation analysis. We applied the MiXeR method and conditional/conjunctional false discovery rate (ccFDR) approach to identify the shared genetic components between the traits. Finally, within the UK biobank dataset, the stratified cox regression analyses were applied to explore the association between T2D and fracture risk by including different number of the risk factors secondary to T2D. As complement, the relationship between T2D and BMD was also investigated.

Results

The genetically predicted type 2 diabetes and fracture

The overall study design was presented in Supplementary Fig. 1. We first assessed the relationship of genetically predicted T2D and fracture in UK biobank dataset with the weighted genetic risk score (wGRS) analysis. Within the 294,571 UK biobank samples (Supplementary Fig. 1), we constructed the weighted genetic risk score (wGRS) for the individuals in the UK Biobank with the 404 SNPs, which were independently associated with type 2 diabetes (Supplementary Table 1). The wGRS of the 404 SNPs were strongly associated with type 2 diabetes in UK Biobank data (OR=1.6, P<2.0×10-16), suggesting that the instruments were powerful for the MR analysis. When we regressed the observed fracture on the wGRS, we found that the genetically determined type 2 diabetes was associated with lower risk of fracture (OR=0.982, 95%CI=0.975-0.989, P=0.006) (adjusting for reference age, sex, BMI, physical activity, fall history, HbA1c and medication treatments) (Fig. 1A). When we classified the fracture sites into weight-bearing bones (neck, vertebrae, pelvic, femur, tibia) and other bones, it indicated that there was trend of protective association between T2D wGRS and weight-bearing bones fracture (OR=0.9772, 95%CI=0.9552-0.9997, P=0.04737, N of fracture=8,992, N of non-fracture=265,262), and other bones fracture (OR=0.9838, 95%CI=0.9688-0.9991, P=0.0386, N of fracture=20,317, N of non-fracture=265,262) (Fig.1A). We further estimated the effect of sex interaction on fracture risk with T2D wGRS × sex interaction term in regression model, and no significant interactions were identified for fracture risk (P=0.5576). Moreover, we conducted the stratified analysis by sex and identified similar trends of association (Supplementary Fig.2A). Meanwhile, the genetically determined type 2 diabetes was associated with higher BMD in pooled samples (Fig. 1B) and in both male and female (Supplementary Fig.2B).

The association of genetically predicted type 2 diabetes with fracture and BMD by using different MR methods. (A) the genetically predicted type 2 diabetes and fracture. (B) the genetically predicted type 2 diabetes and BMD. (C) visualized the association of lead SNPs for type 2 diabetes with the risk of fracture. Abbreviations: wGRS, weighted genetic risk score; MR-PRESSO, MR pleiotropy residual sum and outlier.

We also performed the two-sample MR analyses with fracture GWAS summary data (Trajanoska et al., 2018) which is independent of UK Biobank samples. The inverse variance weighting (IVW) method showed a causal effect of genetically predicted type 2 diabetes on low fracture risk (OR=0.965, 95%CI=0.943-0.988, P=0.003) (Fig. 1A, Supplementary Table 2) using 298 SNPs as the instruments (Supplementary Table 3). This causal relationship was also significant in simple median test (OR=0.967, 95%CI=0.936-0.997, P=0.033) (Fig. 1A). There was heterogeneity in IVW results (Q’ p < 0.05), when we excluded pleiotropic variants using restrictive MR pleiotropy residual sum and outlier test (MR-PRESSO) method, the causal association was still detected between type 2 diabetes and fracture (OR=0.967, 95%CI=0.945-0.989, P=0.004) (Fig. 1A). Moreover, the MR-egger regression also suggested an inverse association between type 2 diabetes and fracture (OR=0.9666, 95%CI=0.9497-0.9828, P=0.0002) (Fig. 1C). The individual effect of the SNPs for type 2 diabetes on fracture was corrected by the false discovery rate (<0.05; Benjamini and Hochberg, 1995), two of 298 lead SNPs (including rs4580892 near RSPO3) of type 2 diabetes remained as potential regions which would also have effect on fracture (Fig. 1C and Supplementary Table 3). We also performed multivariable MR analysis to test the effect of T2D on fracture risk adjusted for confounding factors. We found that T2D had a direct effect on decreased fracture risk adjusted for BMI (OR=0.974, 95%CI=0.953-0.995, P=0.017), and BMI mediated 9.03% of the protective effect (Supplementary Table 2). Similarly, with BMD GWAS summary data (Morris et al., 2019), IVW (β=0.041, 95%CI=0.027-0.054, P=8.14×10-9), simple median (β=0.028, 95%CI=0.020-0.036, P=6.59×10-13) and MR-PRESSO (β=0.029, 95%CI=0.022-0.036, P=1.14×10-14) all showed a causal association between type 2 diabetes and BMD (Fig. 1B, Supplementary Table 2 and Supplementary Table 4). The multivariable MR analysis suggested that T2D also showed direct effect on increased BMD after adjusting for BMI (β=0.042, 95%CI=0.026-0.057, P=1.92×10-7) (Supplementary Table 2).

The distinct signal shared by type 2 diabetes and fracture

Leveraging the genetic summary datasets, we first evaluated the genetic correlation among the traits and diseases by using LDSC (Bulik-Sullivan et al., 2015). It is found that the genetic correlation between type 2 diabetes and fracture was not significant, but with inverse direction (rg=-0.0114) (Supplementary Table 5). Instead, we used MiXeR (Frei et al., 2019) to evaluate the polygenic overlap irrespective of genetic correlation between T2D and fracture. As represented in Venn diagrams of shared and unique polygenic components (Fig. 2a), the MiXeR analysis suggested that type 2 diabetes and fracture exhibited polygenic overlap, sharing 428 causal variants, in other words, 18% of variants (428 of 2370) associated with type 2 diabetes might contribute to the risk of fracture (Dice coefficientLJ=25.25%), and genetic correlation was observed (rg=LJ-0.086) (Fig. 2A and Supplementary Table 6). Only 39% of shared variants between type 2 diabetes and fracture showed concordant direction of association, and the correlation of effect sizes within the shared polygenic component was negative (rho_β=-0.336) (Supplementary Table 6).

The venn diagrams, conditional Q-Q plots and genetic variants jointly associated with type 2 diabetes and fracture/BMD at conjunctional false discovery rate (conjFDR) < 0.05. (A) the shared number of variants between type 2 diabetes and fracture. (B) the conditional Q-Q plot of fracture given association with type 2 diabetes at the level of p ≤ 0.1, p ≤ 0.01, p ≤ 0.001. (C) the shared genetic loci between type 2 diabetes and fracture. (D) the shared number of variants between type 2 diabetes and BMD. (E) the conditional Q-Q plot of BMD given association with type 2 diabetes at the level of p ≤ 0.1, p ≤ 0.01, p ≤ 0.001. (F) the shared genetic loci between type 2 diabetes and BMD.

We used the conditional/conjunctional false discovery rate (ccFDR) approach (Andreassen et al., 2013) to identify specific shared loci between type 2 diabetes and fracture from the GWAS summary statistics. The stratified conditional QQ plot was utilized to visualize the enrichment of association with fracture across varying significance thresholds for type 2 diabetes. We observed leftward deflected from the expected null line in QQ plot, which suggested the existence of polygenic overlap between type 2 diabetes and fracture (Fig. 2B). The conjunctional false discovery rate (conjFDR) analysis identified 10 genomic loci shared between type 2 diabetes and fracture (Fig. 2C and Supplementary Table 7), with the top SNP rs4580892 in the intron of gene RSPO3 (conjFDR=1.45E-05). The shared loci showed mixed directions of allelic associations, with 7 of 10 shared loci had inverse direction of effect between type 2 diabetes and fracture (Supplementary Table 7). We found that the locus approximately 250kb upstream and downstream of the gene RSPO3 (hg19, chr6: 127189749-127689749) possessed many significant SNPs associated with type 2 diabetes and fracture (Fig. 3A and 3B), with the nearest gene RSPO3. The top SNP rs4580892 had inverse direction of effect between type 2 diabetes and fracture, where rs4580892_T allele was associated with increased type 2 diabetes risk (OR=1.041227, P=8.46×10-9) (Fig. 3a) and decreased fracture risk (OR=0.944, P=3.68×10-10) (Fig. 3B).

The distinct signal (RSPO3) shared by type 2 diabetes and fracture. (A) the regional plot of the association of type 2 diabetes and RSPO3 gene expression (adipose subcutaneous) within hg19: chr6:127189749-127689749 (RSPO3 gene region ± 250kb window). (B) the regional plot of the association of fracture and RSPO3 gene expression (adipose subcutaneous) in the same region. (C) the association of genetically predicted RSPO3 gene expression with type 2 diabetes and fracture risk. (D) the association of genetically predicted RSPO3 gene expression with BMI, waist circumference, waist-hip ratio and MRI-derived visceral adipose. (F) the association of genetically predicted circulating RSPO3 with type 2 diabetes and fracture risk.

Further, after merging the eQTL summary data of RSPO3 in adipose subcutaneous with the summary data of type 2 diabetes and fracture, rs72959041 and rs1936806 were the top cis-eQTL for type 2 diabetes and fracture in GTEx database, respectively. By applying the SMR method (Zhu et al., 2016), we found that the higher expression of RSPO3 (ENSG00000146374) in adipose subcutaneous would be associated with increased the risk of type 2 diabetes (OR=1.102, 95%CI=1.031-1.179, P=0.004), but decreased the risk of fracture (OR=0.695, 95%CI=0.566-0.854, P=0.0005) (Fig. 3C and Supplementary Table 8). Interestingly, in adipose subcutaneous, higher expression of RSPO3 was associated with higher waist circumference (β=0.090, 95%CI= 0.029-0.151, P=0.004) and higher waist-hip ratio (β=0.204, 95%CI=0.095-0.313, P=0.0002) (Fig. 3D and Supplementary Table 8). Meanwhile, higher expression of RSPO3 was associated with higher MRI-derived visceral adipose (β=0.199, 95%CI=0.103-0.294, P=4.36×10-5) (Fig. 3D and Supplementary Table 8). The association between the expression of RSPO3 and BMI was not significant, but the direction is the same as the waist circumference (Fig. 3D and Supplementary Table 8). Moreover, to estimate the impact of RSPO3 protein level on type 2 diabetes and fracture risk, we used the top SNP at rs4580892, a cis-pQTL for circulating RSPO3 (P=2.34LJ×LJ10−11) identified by Sun et al in an independent dataset (Sun et al., 2018), to instrument the circulating protein level of RSPO3. The MR analyses indicated that increased circulating RSPO3 was strongly associated with increased risk of type 2 diabetes (OR=1.24, 95%CI=1.16-1.34, P=7.86LJ×LJ10−9), but reduced fracture risk (OR=0.73, 95%CI=0.66-0.81, P=4.1LJ×LJ10−10) (Fig. 3E).

Not surprisingly, type 2 diabetes showed significant positive genetic correlation with BMD (rg=0.0923, P=2.50×10-6) (Supplementary Table 5). The MiXeR analysis suggested that 29% of variants (691 of 2370) associated with type 2 diabetes might contribute to BMD (Dice coefficient=33.67%) (Fig. 2D and Supplementary Table 6). The leftward deflected from the expected null line in QQ plot suggested the existence of polygenic overlap between type 2 diabetes and BMD (Fig. 2E). and the conjFDR analysis identified 661 genomic loci shared between type 2 diabetes and BMD, and 449 of 661 loci (68%) had concordant associations between type 2 diabetes and BMD (Fig. 2F and Supplementary Table 9).

Observed relationship between type 2 diabetes and fracture

Within the 352,879 UK Biobank participants (Supplementary Fig. 1), 13,817 (3.92%) developed type 2 diabetes during 2006 and 2015, with the mean duration of type 2 diabetes 8.34 years (Supplementary Table 10). Compared to those without diabetes, the participants with type 2 diabetes were older (63.20 vs. 60.55, P <2.2×10−16), and more likely to be men and smokers, and had a higher BMI (32.07 vs. 27.08, P <2.2×10−16) (Supplementary Table 10). We identified 16,147 (4.6%) participants with fracture within the 352,879 UK Biobank participants (Supplementary Table 10).

Although we found that genetically predicted type 2 diabetes might not be associated with risk of fracture, we observed a higher risk of fracture in the type 2 diabetes patients in the cox proportional hazards regression after adjusted for the reference age, sex, BMI, physical activity, HbA1c, medication treatments and fall history (Model 0) (HR=1.527, 95%CI 1.385-1.685, P<2×10−16) (Fig. 4A and Supplementary Table 11). And the average causal mediation effect (ACME) by BMI was protective with 30.2% of the intermediary effect, respectively (BMI: indirect effect=-0.003, P <2×10-16) (Supplementary Table 12). Similar findings were observed for both male and female (HR=1.587, 95%CI 1.379-1.828, P=1.26×10−10 in male, HR=1.530, 95%CI 1.334-1.756, P=1.27×10−9 in female) (Supplementary Fig.3). When we additionally controlled for BMD, we still observed increased risk of fracture in type 2 diabetes (HR=1.574, 95%CI 1.425-1.739, P<2×10−16) (Model 1) (Fig. 4A and Supplementary Table 11). We also classified the fracture into weight-bearing bones fracture (neck, vertebrae, pelvic, femur, tibia) and other bones fracture. Similar trends of association were observed in model 0 (weight-bearing bones: HR=1.792, 95%CI 1.555-2,065, P=8.25×10−16; other bones: HR=1.337, 95%CI 1.167-1.531, P=2.85×10−5) and model 1 (weight-bearing bones: HR=1.850, 95%CI 1.602-2,136, P<2×10−16; other bones: HR=1.377, 95%CI 1.199-1.580, P=5.54×10−6) (Fig. 4A).

The regression between type 2 diabetes and fracture/BMD in observational study. (A) the relationship between type 2 diabetes and fracture in different sites including all fracture, weight-bearing bones and other bones. Model 0 adjusted for reference age, sex, BMI, physical activity, fall history, HbA1c and medication treatments; Model 1 adjusted for a+ BMD. (B) the stratified analyses between type 2 diabetes, and fracture based on the five risk factors secondary to the diseases adjusted for the age and sex. (C) the relationship between type 2 diabetes and BMD and the stratified analyses based on the five risk factors secondary to the diseases adjusted for the age and sex. The 5 risk factors were: (1) BMI≤25kg/m2; (2) no physical activity; (3) falls in the last year; (4) HbA1c≥47.5mmol/mol; (5) antidiabetic medication treatments.

Inspired by the MR analyses that genetically determined type 2 diabetes might not be a risk factor for fracture, we conducted stratified analyses based on the five risk factors secondary to the diseases, such as BMI≤25kg/m2, no physical activity, falls in the last year, HbA1c≥47.5mmol/mol and antidiabetic medication treatments. Within the 13,817 individuals with type 2 diabetes, 2,303 patients carried none of the above risk factors, 4,128 patients accompanied with one of the risk factors, and 4,252 patients carried at least two risk factors (Supplementary Table 13). We performed stratified cox regression analysis and found that type 2 diabetes with at least two risk factors were associated with an increase of fracture risk (HR=1.39, 95%CI 1.269-1.514, P=1.32×10−7) (Fig. 4B). It is interesting to note that including only one risk factor would attenuate the effect size of type 2 diabetes on fracture risk (HR=1.20, 95%CI=1.060-1.338, P=0.010). Furthermore, the association between type 2 diabetes and fracture was not significant (P=0.452) when analyzing the type 2 diabetes without risk factors (HR=1.08, 95%CI 0.880-1.280, N=2,303) (Fig. 4B). Similar trends of association were observed in both male and female (Supplementary Fig.4A). These results suggested that the risk factors secondary to type 2 diabetes might contribute to the risk of fracture instead of the disease itself.

We also observed that participants with diabetes, despite they were older, had a significantly higher BMD than subjects without diabetes (0.57 vs. 0.54, P<2.2×10−16) (Supplementary Table 10). In the multivariable linear regression analysis, the type 2 diabetes was found to be associated with increased BMD in same model adjusted for age, sex, BMI, physical activity, fall history, HbA1c and medication treatments (β=0.00957, P=1.35×10−10) (Fig 4C, Supplementary Table 11). We examined the relationship between type 2 diabetes and BMD in subgroups with varying numbers of risk factors. We observed that the effect size of type 2 diabetes on BMD decreased when the number of risk factors for type 2 diabetes increased (no risk factors: in pooled β=0.023, P<2×10−16, in male β=0.018, P=2.09×10−6, in female β=0.028, P=2.97×10−11; one risk factor: in pooled β=0.020, P<2×10−16, in male β=0.015, P=1.91×10−7, in female β=0.0245, P=2.46×10−16; at least two risk factors: β=0.017, P<2×10−16, in male β=0.0098, P=7.74×10−4, in female β=0.0250, P<2×10−16) (Fig. 4C, Supplementary Fig.4B).

Discussion

By leveraging the genetic datasets, we found that the genetically predicted type 2 diabetes was associated with higher BMD and lower risk of fracture in both one-sample MR (with 404 IVs) and two-sample MR (with 298 IVs). We also identified ten genomic loci shared between fracture and type 2 diabetes, with the top signal at SNP rs4580892 in the intron of gene RSPO3. And the higher expression of RSPO3 in adipose subcutaneous was associated with increased the risk of type 2 diabetes, but decreased the risk of fracture. Similarly, the increased circulating RSPO3 was strongly associated with increased risk of type 2 diabetes, but reduced fracture risk. In the prospective study, type 2 diabetes was observed to be associated with higher risk of fracture, but BMI mediated 30.2% of the protective effect. However, when stratified by the risk factors secondary to the disease, we observed that the effect size of type 2 diabetes on the risk of fracture decreased when the number of risk factors secondary to type 2 diabetes decreased, and the association became not significant if the type 2 diabetes patients carried none of the risk factors.

The “diabetic bone paradox” suggested that T2D patients would have higher areal bone mineral density (BMD) but higher fracture risk than individuals without T2D (Botella Martinez et al., 2016; Romero-Diaz et al., 2021). Other measurements, such as trabecular bone score (Fazullina et al., 2022; Ho-Pham and Nguyen, 2019) and chest CT texture analysis (Kim et al., 2023), could provide additional valuable information in the evaluation of fracture risk, especially in type 2 diabetes patients. As reviewed previously, heterogeneity could exist from study to study, and conflicting observational findings were reported (Khosla et al., 2021). Mendelian randomization (MR) could be an alternative approach to infer the relationship between exposure and outcome, as this method exploits the idea that genotypes are distributed randomly at conception, facilitating their use as instrumental variables (IV) to alleviate the bias of the unknown confounding factors (Davey Smith and Hemani, 2014; Zhao et al., 2019). Trajanoska K et al (Trajanoska et al., 2018) assessed the effect of 15 selected clinical risk factors on the risk of fracture by using two-sample MR analysis, and reported non-significant relationship between type 2 diabetes and fracture risk, but the direction of effect was negative (OR=0.99). However, only 38 SNPs were extracted as instruments from a GWAS published in 2012 (Morris et al., 2012). In the present study, we extracted 298 T2D-associated independent SNPs from Mahajan A et al (Mahajan et al., 2022), which is the largest-scale GWAS meta-analysis to date published in 2022, as the IVs in two-sample MR analysis. We reported that genetically determined type 2 diabetes was associated with lower risk of fracture, even in multivariable MR analysis adjusted for BMI. In addition, we also calculated the wGRS with 404 T2D-associated independent SNPs in the UK Biobank dataset, and performed regression analysis of wGRS of type 2 diabetes on the fracture risk (one-sample MR). Again, we found that the genetically predicted type 2 diabetes was associated with lower risk of fracture in one-sample MR analysis. To be note, two-sample MR results could be served as an independent replication to the one-sample MR results, because the effects of the outcome (fracture risk) for two-sample MR were derived Trajanoska K et al (Trajanoska et al., 2018), while the one-sample MR used the UK Biobank dataset, the study samples had no overlap. Further, consistent with previous studies (Ahmad et al., 2017; Mitchell et al., 2021), the MR analysis in the present study suggested that the genetically predicted type 2 diabetes was associated with higher BMD. That is to say, by alleviating the bias of the unknown confounding factors through MR analysis, the genetically predicted type 2 diabetes did not show this bone paradox.

The genetic correlation between type 2 diabetes and fracture estimated by LDSC (Bulik-Sullivan et al., 2015) was not significant. It is hard for LDSC to identify the genetic pleiotropy with mixed-effect directions, which is what usually happens between two complex traits. Therefore, in this study, we employed the MiXeR method (Frei et al., 2019), which could identify the unique and shared polygenic SNPs between two traits regardless of genetic correlation. The MiXeR analysis suggested that type 2 diabetes and fracture exhibited polygenic overlap, and 61% of shared variants between type 2 diabetes and fracture showed discordant direction of association, and the correlation of effect sizes within the shared polygenic component was negative, suggesting an inverse genetic relationship between type 2 diabetes and fracture. Additionally, the conditional/conjunctional false discovery rate analysis (Andreassen et al., 2013) suggested a top locus at 6q22 (SNP rs4580892) jointly associated with type 2 diabetes and fracture. This SNP is an intronic variant in gene RSPO3 (ENSG00000146374). We found that higher expression of RSPO3 in adipose subcutaneous would be associated with increased risk of type 2 diabetes, but decreased risk of fracture. RSPO3 is a known WNT-signaling modulator (Baron and Kneissel, 2013; Lerner and Ohlsson, 2015), which could bind to LRP5/6 to enhance the activity of osteoblast (Richards et al., 2012). Karin et al demonstrated that RSPO3 is expressed in osteoprogenitor cells and osteoblasts, and that osteoblast-derived RSPO3 is the principal source of RSPO3 in bone and is an important regulator of vertebral trabecular bone mass and bone strength in adult mice (Nilsson et al., 2021). Interestingly, we also found that higher expression of RSPO3 was associated with higher waist circumference and higher waist-hip ratio. It was reported that RSPO3 could impact body fat distribution (Loh et al., 2020), and fat distribution is an independent predictor of type 2 diabetes (Manolopoulos et al., 2010). Therefore, we speculated that the different role of the shared genetic components between bone metabolism and type 2 diabetes might provide one possible explanation to the inverse association pattern, and obese tendency might mediate this pattern. In fact, in UK Biobank, the participants with type 2 diabetes had a higher BMI compared to those without diabetes (31.74 vs. 27.05, P<2.2×10−16), and BMI mediated 30.2% of the intermediary effect between type 2 diabetes and fracture in our study.

On the contrary, in observational study, we found that type 2 diabetes was associated with higher risk of fracture patients even adjusted for a bunch of confounding factors such as the age, sex, BMI, physical activity, HbA1c, medication treatments and fall history. When fractures were categorized into different sites (weight-bearing bones and other bones), the association between type 2 diabetes and fracture remained evident. Inspired by the MR analyses that genetically determined type 2 diabetes might not be a risk factor for fracture, we started to perform the stratified analyses based on the risk factors secondary to the diseases. There are many secondary factors associated with type 2 diabetes might contribute to fracture risk. Our previous study suggested that keeping moderate-high BMI (overweight) might be of benefit to old people in terms of fracture risk (Zhu et al., 2022), and an intensive lifestyle intervention, such as weight loss, in T2D patients might increase fracture risk (Johnson et al., 2017). The hyperglycemia could cause osteocyte senescence and premature programmed cell death, leading to decreased ability to sense and respond to mechanical stimuli such as oscillatory shear stress, ultimately contributing to skeletal fragility (Eckhardt et al., 2020). Another major and complicated factor that might influence the risk of fracture in T2D patients is the use of diabetes mellitus medications. For example, the use of insulin (Napoli et al., 2014) or thiazolidinediones (Zhu et al., 2014) was reported to be associated with an increased risk of fracture. In addition, the risk of falls, which might be triggered by some diabetic complications such as visual impairment and peripheral neuropathy, was suggested to increase in the patients with T2D (Schwartz et al., 2002). Besides, low physical activity was identified as one of the most important independent diabetes-related risk factors for fracture through Gradient Boosting Machines (Axelsson et al., 2023). Therefore, in this study, we stratified the T2D patients with five secondary risk factors (BMI≤25kg/m2, no physical activity, falls in the last year, HbA1c≥47.5mmol/mol and antidiabetic medication treatment), and found that the effect size of type 2 diabetes on the risk of fracture decreased when the risk factors secondary to type 2 diabetes decreased, and the association became not significant if the type 2 diabetes patients carried none of the risk factors. That is to say, the diabetic bone paradox might not exist if the secondary risk factors of type 2 diabetes were eliminated. In a recent large-scale cohort study, four factors (duration of T2D, low physical activity, BMI, and insulin treatment) were identified as the important risk factors for fracture among the T2D patients, and the patients without the risk factors had lower fracture risk than their matched controls (Axelsson et al., 2023). One previous study using the UK primary care data found no evidence to suggest a higher risk of fracture in type 2 diabetes patients (Davie et al., 2021) and reported that significantly lower fracture risk was observed for overweight individuals (BMI 25-30kg/m2) in type 2 diabetes than their counterparts without type 2 diabetes (Davie et al., 2021).

In summary, we found that the genetically determined type 2 diabetes might not be associated with higher risk of fracture. And the shared genetic architecture between type 2 diabetes and fracture suggested a top signal near RSPO3 gene. In addition, the stratified prospective regression analysis suggested that the effect size of type 2 diabetes on the risk of fracture decreased if risk factors secondary to type 2 diabetes could be eliminated. These genetic and observational evidences prevailed us to hypothesize that the risk factors secondary to type 2 diabetes might contribute more to the risk of fracture than the disease per se. And it is important to manage the complications of type 2 diabetes to prevent the risk of fracture.

Methods

Study participants and weighted genetic risk score (wGRS) analysis

The UK Biobank data, the application 41376 as we used before (Bai et al., 2020), was applied in this study under a prospective design. We identified the individuals with T2D and fracture using the ICD codes and self-report status. The detail information on the field ID and codes for data extraction from UK Biobank was listed in Supplementary Table 14. We excluded participants if they were identified as follows: 1) ethnically identified as non-European (n =30,481); 2) diagnosed as type 1 diabetes (n=4,455); 3) diagnosed with diseases associated with bone loss (n=21,560); 4) diagnosed as fracture with known primary diseases (n=7,222) (Supplementary Table 15). For the remaining 439,982 samples, we further excluded 145,411 participants with relatedness (kinship) with others in the wGRS analysis (294,571 participants left) (Supplementary Fig. 1). In addition, we classified the fracture into weight-bearing bones (neck, vertebrae, pelvic, femur, tibia) (N of fracture=8,992, N of non-fracture=285,579) and other bones (skull and facial, ribs, sternum, forearm, wrist and hand, foot and other unspecified body regions) (N of fracture=20,317, N of non-fracture=274,254) using the ICD codes and self-report status (Supplementary Table 16).

The summary-statistic data for the type 2 diabetes were obtained from a very recent GWAS consisting of 80,154 individuals with type 2 diabetes and 853,816 controls in European population (with 10,454,875 SNPs) (Mahajan et al., 2022). We drew a set of independent genetic variants with genome-wide significance (P<5×10−8) from the type 2 diabetes summary-statistic data by LD clumping based on r2 < 0.1 in 500 kb window to serve as instrumental variables (N=404) (Supplementary Table 1). We constructed the weighted genetic risk score (wGRS) for the individuals in the UK biobank (294,571 samples with genotypes) as a linear combination of the selected SNPs weighted by their β coefficients on type 2 diabetes: wGRS = β1 ×SNP1 + β2 ×SNP2 + … + βn ×SNPn. n is the number of instrumental variables (here N=404 after LD clumping based on r2 < 0.1 in 500 kb window). Next, cox proportional hazards regression and linear regression analyses were performed to analyze the association between the wGRS and fracture/BMD adjusted age, sex, BMI, physical activity, fall history, HbA1c and medication treatments. Besides, regression modeling was used to estimate the effect of gene-environment interaction (T2D wGRS × sex) on fracture risk and BMD. In addition to the T2D wGRS × sex interaction term, the model was adjusted for covariates: age, sex, BMI, physical activity, fall history, HbA1c and medication treatments.

Two-sample Mendelian Randomization (MR) analyses

To validate the wGRS results, we also performed the two-sample MR analyses that is independent of UK Biobank samples. The summary-statistic data for fracture (a discovery set of 37,857 fracture cases and 227,116 controls with 2,539,800 SNPs) (Trajanoska et al., 2018) and BMD (426,824 samples and 13,753,401 SNPs) (Morris et al., 2019) were extracted from the GEFOS consortium (http://www.gefos.org/), while the summary data for type 2 diabetes (Mahajan et al., 2022) is the same as used in wGRS analysis.

We used the inverse variance weighting (IVW) (Burgess et al., 2013), simple median and MR-PRESSO (Verbanck et al., 2018) approaches in two-sample MR analyses. For the outcome of fracture, we merged the two summary datasets for T2D and fracture (Mahajan et al., 2022; Trajanoska et al., 2018) and got 2,479,475 overlapping SNPs, of which 6,946 SNPs were genome-wide significant for type 2 diabetes. After LD clumping based on r2 < 0.1 in 500 kb window, 298 independent genetic variants were left (Supplementary Table 2). Similarly, we got 9,204,694 overlapping SNPs for type 2 diabetes and BMD, and 389 independent genetic variants were left after LD clumping (Supplementary Table 3). After harmonizing the effects so that they reflect the same allele, 289 (for fracture) and 380 (for BMD) SNPs were finally used in the IVW and simple median MR analysis. Because the presence of horizontal pleiotropy could bias the MR estimates, we additionally used the MR-PRESSO. The two-sample MR analyses were conducted in R version 4.0.2 using TwoSampleMR (Hemani et al., 2018), MendelianRandomization (Yavorska and Burgess, 2017) and MR-PRESSO (Verbanck et al., 2018) packages. Moreover, we regressed the effects of these 298 SNPs of both traits to highlight the overall effect of T2D on fracture with ‘grs.summary’ function in the R package ‘gtx’ (http://www2.uaem.mx/r-mirror/web/packages/gtx/gtx.pdf).

Multivariable Mendelian Randomization (MVMR) analysis

Next, we conducted multivariable MR analysis (Sanderson et al., 2019; Sanderson et al., 2021) to examine the direct effect of T2D on fracture and BMD adjusted for BMI with ‘MVMR’ R package (https://github.com/WSpiller/MVMR). After adjusting for confounders, the effect of exposure on the outcome was considered to be a direct effect. Specifically, we first extracted the overlapping SNPs from the summary data for T2D (Mahajan et al., 2022), BMI (Locke et al., 2015) and fracture (Trajanoska et al., 2018). Then the independent significant SNPs (P<5×10−8 and R2<0.1) for either T2D or BMI were pooled as instruments. Additionally, we performed SNP harmonization to correct the orientation of alleles. The final IVs used in MVMR were presented in Supplementary Table 17.

Infer the shared genetics

With the summary-statistic GWAS data of type 2 diabetes (Mahajan et al., 2022), BMD (Morris et al., 2019) and fracture (Trajanoska et al., 2018), we performed genome-wide genetic correlation analysis between type 2 diabetes and fracture/BMD by using linkage disequilibrium score regression (LDSC) (Bulik-Sullivan et al., 2015; Bulik-Sullivan et al., 2015) which estimates the degree of shared genetic factors between two traits, We used MiXeR (https://github.com/precimed/mixer, version 1.2.0) to quantify polygenic overlap (e.g. how many unique and shared polygenic SNPs for type 2 diabetes and fracture) irrespective of genetic correlation (Frei et al., 2019). MiXeR models additive genetic effects as a mixture of four components, representing null SNPs in both traits (π0); SNPs with a specific effect on the first and on the second trait (π1 and π2, respectively); and SNPs with non-zero effect on both traits (π12). The dice coefficient of two traits was estimated as (Frei et al., 2019). We constructed conditional quantile-quantile (QQ) plots which reveals the distribution of P values for fracture/BMD conditioning on the significance of association with type 2 diabetes at the level of P□<□0.1, P□<0□.01 and P□<0□.001 to visualize polygenic enrichment (Schwartzman and Lin, 2011). We used the conditional/conjunctional false discovery rate (ccFDR) approach (https://github.com/precimed/pleiofdr) to identify the specific shared loci (Andreassen et al., 2013). We used the conditional false discovery rate (condFDR) to detect SNPs associated with fracture given associations with type 2 diabetes. We denoted condFDR for fracture given associations with type 2 diabetes as condFDR(fracture|T2D) and vice versa, and considered the significance cutoff < 0.01. We used conjunctional FDR (conjFDR) to identify SNPs jointly associated with type 2 diabetes and fracture. After repeating the condFDR procedure for both traits, the conjFDR analysis reported the loci that exceed a condFDR significance threshold for two traits simultaneously (the maximum between the condFDRs for both traits), conjFDR < 0.05 was set as the significance.

We employed the Summary-data-based Mendelian Randomization (SMR) method developed by colleagues (Zhu et al., 2016) to test the association of the expression level of gene RSPO3 with BMI (Locke et al., 2015), waist circumference, waist-hip ratio (Shungin et al., 2015) and MRI-derived visceral adipose (Agrawal et al., 2022), type 2 diabetes and fracture using summary-level data from GWAS and expression quantitative trait loci (eQTL) data of the subcutaneous adipose tissue (9,962,255 SNPs included) from the GTEx database (release v8) (https://www.gtexportal.org/home/; Consortium et al., 2017). In the SMR analysis, the top cis-eQTL genetic variants were used as the instrumental variables (IVs) for gene expression. Additionally, we downloaded the cis-pQTL summary data for the circulating RSPO3 reported in the study by Sun et al (Sun et al., 2018), and performed MR analyses to determine the association of circulating RSPO3 with type 2 diabetes and fracture risk using the inverse variance-weighted (IVW) MR approach.

Observational analyses

For the 439,982 UK biobank samples (see foregoing description of study participants), we only focused the participants diagnosed with T2D within the 10-year period from 1 January 2006 to 31 December 2015, leaving 425,772 participants (with 14,860 type 2 diabetes patients) (Supplementary Fig. 1). Here, each type 2 diabetes patient had a diagnosis date, taking this date as the reference date, we first calculated the onset age, then among the participants who were free of T2D, we selected up to 27 participants (whenever possible) whose age at the reference date (± 3 years) could be matching to the onset age as referents (Supplementary Fig. 5). In total, 363,884 non-T2D referents were individually matched with 6-year age band at the reference date. We prospectively followed these type 2 diabetes patients and referents from the reference date until diagnosis of fracture, death, emigration, 19 April 2021 (diagnose a fracture of the last person in the cohort), whichever came first. Survival time was calculated based on whether the patient had a fracture. If individuals had a fracture, the survival time is calculated as the time of the first diagnosis of fracture minus the reference date. If individuals did not have a fracture, it was defined as the minimum time of the reference date to diagnose a fracture of the last person in the cohort (19 April 2021), death, or emigration date. We excluded 25,865 participants with fracture diagnosis date, or death or emigration date before the reference date, leaving 352,879 participants included in the final analysis (13,817 type 2 diabetes patients and 339,062 referents) (Supplementary Fig 1 and Supplementary Table 10). Cox proportional hazards regression, as a statistical method to analyze the effect of risk factor on the time it takes for a specific event to happen, was used to test the relationship between T2D and fracture. Meanwhile, multiple linear regression analyses were performed to test the association between T2D and BMD. Here, the BMD was estimated from quantitative ultrasound measurement at heel. Use of the device generates two variables, including speed of sound (SOS) and BUA (the slope between the attenuation of the sound signal and its frequency as it travels through the bone and soft tissue). Heel BMD was calculated by the following formula: BMD = 0.002592 ×(BUA+SOS)−3.687.

First, we adjusted for clinical risk factors including reference age, sex, BMI, physical activity, fall history, HbA1c and medication treatments to examine the relationship between T2D and fracture/BMD (Model 0). To examine the intermediary effect of risk factors on the relationship between T2D and fracture, the mediation analysis was performed using the R packages of “mediation”. Individuals treated with any glucose-lowering medication including insulin product, metformin, troglitazone, pioglitazone, rosiglitazone, tolbutamide, glibenclamide, gliclazide, glipizide, gliquidone, glimepiride, chlorpropamide, tolbutamide, repaglinide and nateglinide were recorded as having received medical treatment. We also included BMD as an additional confounding factor for fracture analysis as complement to the basic model (Model 1). In addition, we classified the fracture into weight-bearing bones (neck, vertebrae, pelvic, femur, tibia) (N of fracture=8,992, N of non-fracture=285,579) and other bones (skull and facial, ribs, sternum, forearm, wrist and hand, foot and other unspecified body regions) (N of fracture=20,317, N of non-fracture=274,254) using the ICD codes and self-report status (Supplementary Table 16). As we did in wGRS analysis, we classified the fracture into weight-bearing bones and other bones fracture. Briefly, 6,582 (1.92%) participants were identified as weight-bearing bones and 9,586 (2.77%) participants were identified as other bones. Second, we carried out stratified analyses between type 2 diabetes, fracture and BMD based on the risk factors secondary to the disease. We took 5 clinical factors to classify the individuals at risk, for example, if an individual had BMI≤25kg/m2, no physical activity, falls in the last year, HbA1c≥47.5mmol/mol and antidiabetic medication treatment, this individual was identified to have 5 risk factors, and so forth. Based on the number of risk factors, we grouped 13,817 individuals with T2D into subgroup for analysis. This analysis was adjusted for reference age, sex. For the BMD analysis, the age when attended assessment center was included in the analysis instead of reference age.

Data Availability

All data produced in the present work are contained in the manuscript

Article and author information

Authors’ contributions

H-F.Z. conceptualized and designed the study. P-P.Z. conducted the analyses and drafted the paper. Z-F.S., and L.X., helped in clinical interpretation of the results and drafted the manuscript. P.L., W-J.X., C-D.Y., Z-W.X., M-Y.Y., Y.Q., J-D.Z., J-X.G., and D.K. helped in the data acquisition, result interpretation and manuscript editing. The authors read and approved the final manuscript.

Competing interests

The authors declare that they have no competing interests.

Funds

This work was supported by the "Pioneer" and "Leading Goose" R&D Program of Zhejiang (#2023C03164), the National Natural Science Foundation of China (#82370887), the Chinese National Key Technology R&D Program, Ministry of Science and Technology (#2021YFC2501702), and the funds from the Westlake Laboratory of Life Sciences and Biomedicine (#202208014).

Acknowledgements

We thankfully acknowledge the High-performance Computing Center at Westlake.

Ethics

All individuals provided written informed consent. The North West Multi-Centre Research Ethics Committee approved the UK Biobank ethical application (reference number: 16/NW/0274).

Supplementary Figures

Flow chart of the overall study design.

The regression between type 2 diabetes GRS and fracture risk/BMD in male and female. Models were adjusted for reference age, sex, BMI, physical activity, fall history, HbA1c and medication treatments.

The regression between type 2 diabetes and fracture risk/BMD in male and female in observational study. Model 0 adjusted for reference age, sex, BMI, physical activity, fall history, HbA1c and medication treatments; Model 1 adjusted for a+ BMD.

The stratified analyses between type 2 diabetes and fracture/BMD based on the five risk factors (BMI≤25kg/m2, no physical activity, falls in the last year, HbA1c≥47.5mmol/mol and antidiabetic medication treatments) in male and female. Models were adjusted for the age and sex.

The matching information of 14,860 patients with fracture.