Introduction

Osteoporosis (OP) is defined by a decrease in bone mineral density (BMD), degradation of bone microstructure, and increased bone fragility, which finally leads to reduced bone strength and an increased risk of fragility fractures.¹ Early detection, prevention, and management can minimize the burden of osteoporosis on individuals and society, while measurement of bone density by dual-energy X-ray absorptiometry (DXA) in the same bone site is a crucial diagnostic method for OP.² Among them, the femoral neck, lumbar spine, and forearm are the most commonly used bone sites for detecting BMD.³ Recently, the heel region has also been acknowledged as a reliable substitute for assessing BMD.⁴ The costs of treating OP are significant in many countries and are expected to increase, placing a substantial burden on both individuals and society.⁵

Similar to spontaneous bone loss, cervical spondylosis (CS) is one of the main features of aging that refers to widespread degenerative changes in the components of the cervical spine (ie the vertebral disc, tiny joints, the Luschka joint, ligaments and vertebrates).⁶ Until 2017, around 2.887 million individuals globally experienced neck pain symptoms, with an occurrence rate of up to 3551 cases per 10 million people.⁷ Although non-surgical treatment targeting the affected nerve roots or spinal cord can effectively alleviate the pain, many of patients still suffer neurological damage because of the absence of surgical intervention.⁸

Historically, it has been proposed that reduced BMD could worsen nearby disc degeneration, perhaps contributing to the arrival of associated disorders, but most of these studies have focused on the lumbar spine.^9–11 However, significant anatomical and degenerative differences exist between the cervical and lumbar spine, necessitating caution when extrapolating lumbar spine findings to cervical applications.^12,13 A recent observational study revealed that individuals with higher BMD are more likely to develop CS, which goes against traditional beliefs.¹⁴ This might arise from elevated vertebral BMD, particularly in the sub-endplate region, which obstructs medullary arterial branch vascular channels and nutrient diffusion, accompanied by reduced vascular budding and angiogenesis.^15,16 Conversely, OP can delay disk degeneration by increasing diffusive transport of nutrients through mechanical and vascular pathophysiologic pathways within the intervertebral disks.¹⁷ On the other hand, spinal degenerative changes may drive localized BMD elevation, primarily attributed to stress distribution-induced sclerotic bone formation. For instance, in axial radiographic spondyloarthritis, the initial bone mineral loss caused by chronic inflammation becomes obscured by subsequent osteogenesis.¹⁸ As cervical disc degeneration progresses, bone density in the vertebral body and sub-endplate bone increased, particularly in the caudal region.¹⁹ In an experimental study of aging zebrafish, researchers found osteoporosis-like changes, including disrupted bone morphology, microstructural alterations, mineral heterogeneity, and elevated fragility. These changes were accompanied by degeneration of the intervertebral ligament (equivalent to the annulus fibrosus), as well as dehydration and cell abnormalities in the disc center (nucleus pulposus equivalent).²⁰ Unfortunately, although scholars have established significant associations between bone metabolism, cervical disc degeneration, bone turnover markers, and amino acid levels, the understanding of the relationship between BMD and CS still remains incomplete in cohort studies or cross-sectional studies due to potential confounders and inadequate reverse causality.²¹

As an alternative research method, Mendelian randomization (MR) uses single nucleotide polymorphisms (SNPs) as instrumental variables (IVs) to mimic the randomization commonly found in randomized controlled trials (RCTs).²² With the lack of reliable randomized controlled trials, MR is commonly employed to investigate the causal effect of risk factors on outcomes.²³ Furthermore, as the genotype of the organism is established before the disease begins and is not influenced by how it developed, undesirable defects caused by potential confounders or reverse causal relationships can be reduced to a minimum.²⁴ To the best of our knowledge, there are currently no MR studies to explore the causal relationship between BMD and CS. Thus, we employed publicly available genome-wide association study data as the basis for bidirectional two-sample MR and Multivariable Mendelian randomization (MVMR), which aimed to enhance our comprehension of the pathophysiology of CS and offered evidence for clinical treatment.

Materials and Methods

Study Design

In this study, we employed bidirectional two-sample MR and Multivariable Mendelian randomization (MVMR) approach utilizing summary statistics derived from genome-wide association studies (GWASs) to examine the associations between BMD which across five age groups and four skeletal sites and CS. To ensure the robustness of our MR analysis, we adhered to three fundamental assumptions: (i) the selected genetic variants should exhibit significant associations with the exposure of interest; (ii) the genetic variants employed should not exert any influence on the outcome, except through the chosen exposure; and (iii) the genetic variants should not be correlated with any confounding factors that may impact the relationship between the exposure and the outcome. Previous studies have provided evidence of genetic correlation between multiple traits (body mass index, type 2 diabetes, hyperlipidemia, smoking, etc). and IVs used for exposures.^25–28 Hence, those variables were individually subjected to Multivariable Mendelian randomization (MVMR) with BMDs across five age groups and four skeletal sites. The detailed study design is depicted in Figure 1.

Figure 1 Overview of the MR analysis design. The flowchart illustrates three fundamental MR assumptions: (1) Assumption 1 represents the relevance assumption; (2) Assumption 2 represents the Independence assumption; (3) Assumption 3 represents the exclusion-restriction assumption. Dashed lines indicate irrelevant associations (marked by red cross); Solid lines indicate correlations (marked by green checkmark). Abbreviations: IVs, instrumental variables; BMD, bone mineral density; CS, cervical spondylosis.

Data Sources

The summary statistics for BMD were derived from measurements taken at four common skeletal sites in postmenopausal women and men aged 50 years or older: the femoral neck (FN), lumbar spine (LS), forearm (FA), and heel bone (HB). GWAS data for FN-BMD (N = 32,735), LS-BMD (N = 28,498), and FA-BMD (N = 8143) were obtained from the GEFOS consortium—the largest DXA-based BMD GWAS to date—with genetic values adjusted for sex, age, age squared, and weight.²⁹ Moreover, the GWAS on HB-BMD (N = 426,824), measured via ultrasound, were adjusted for age, sex, genotyping array, and assessment center using a fixed-effects model.³⁰ TB-BMD provides the most reliable and unbiased measure of BMD variation at the same skeletal site across the lifespan. The TB-BMD data (N = 56,284) were obtained from a large GWAS meta-analysis and adjusted for age, weight, height, genomic principal components, and study-specific covariates using a linear regression model.³¹ For unrelated individuals, residuals were calculated separately by sex, while in family-based studies, sex was included as a covariate. Based on TB-BMD, participants were categorized into five age groups: 0–15 (N = 11,807), 15–30 (N = 4180), 30–45 (N = 10,062), 45–60 (N = 18,805), and over 60 (N = 22,504). The characteristics of each contributing study are presented in Table 1.

Table 1 Data Sources Used in This Study

FinnGen biobank is a large-scale research project in genomics and personalized medicine, which has collected and analyzed genome and health data from 500,000 Finnish biobank donors. Summary statistics for CS were retrieved from the FinnGen biobank, comprising 14,670 cases (7794 female and 6876 male) and 294,770 controls, with an average age of first onset at 48.99 years (https://r10.risteys.finregistry.fi/endpoints/M13_CERVICDISC).³² Data involved were publicly available from the FinnGen biobank (https://r10.finngen.fi/, accessed March 28, 2024) and the trait of this study was labelled “Cervical disc disorders”.

Our study was based on previously published studies or publicly available GWAS summary statistics, with the majority of participants being of European ancestry to reduce potential bias from population heterogeneity. The Research Ethics Committee of the Second Affiliated Hospital of Zhejiang Chinese Medical University deemed this work exempt from the Ethics Committee review.

Selection of Instrumental Variables (IVs)

Initially, for the IV screening phase of our MR investigation, we implemented a genome-wide significance threshold (). Following that, we implemented SNP pruning within a window size of 10,000 kb, ensuring that the linkage disequilibrium (LD) between SNPs remained below the threshold of . The selected SNPs related to BMDs were retrieved from the GWAS data of the outcome, we further screened the extracted SNPs to remove those with non-concordant alleles and palindromic SNPs with ambiguous strands. Subsequently, we carefully assessed the IVs selected through the above process for any indications of weak instrument bias by using the F-statistic. To mitigate bias stemming from weak IVs, we discarded IVs with a F value below 10, thereby ensuring that only those with a F value exceeding 10 were retained. The calculation formula is , , where represents the coefficient of determination between the SNP and the exposure, N is the sample size, and k is the number of SNPs used as instruments.³³

MR Statistical Analysis

The inverse-variance-weighted (IVW) approach was used as primary MR analysis method, which estimates the causal effect of gene-predicted exposure on outcome through a weighted regression of SNP-specific Wald ratios (ie, β outcome/β exposure).³⁴ Cochran’s Q test and statistics were carried out to evaluate the heterogeneity of the instrumental genetic variable, and a indicated significant heterogeneity.³⁵ When significant heterogeneity exists, the MR analysis concludes with the random effects IVW, whereas MR-Egger, weighted median, weighted mode, and simple mode were used as supplementary analysis methods.³⁶ Conversely, fixed-effects IVW is employed when heterogeneity is absent. MR-Egger slopes are relatively effective as estimates of MR in the presence of horizontal pleiotropy (). The weighted median method ensures the stability of causality estimates by eliminating errors in the presence of 50% invalid IVs, and may provide better causality detection than the MR-Egger under certain conditions.

Sensitivity Analysis

To test the robustness of our results, a “leave-one-out” analysis was conducted by deleting one SNP at a time to examine the influence of SNPs with potential horizontal pleiotropy on MR estimates.³⁷ Furthermore, the MR-Egger intercept test was conducted. If the intercept term in the MR-Egger intercept analysis proved statistically significant, it indicated substantial horizontal pleiotropy within the study.³⁸ A statistically significant result from the Cochran’s Q statistic test signified a significant level of heterogeneity within the analysis.³⁹ To further enhance the reliability of our MR analysis, we implemented the MR-PRESSO method to identify and remove any potential outliers in our instrumental variable set. Outliers can distort MR results and affect the validity of causal inference.⁴⁰ Finally, sensitivity analyses were visually illustrated using scatter, funnel, and forest plots to demonstrate the robustness of the MR study.

Statistical Analysis

In summary, our analysis strictly adhered to the following criteria: (i) All selected IVs were initially screened using Phenoscanner to exclude potential confounding effects on genetic variants; (ii) IVW analysis was used as the primary analytical method, with Bonferroni correction applied to determine statistical significance; (iii) Multiple MR analysis methods were integrated to ensure consistency in the direction and magnitude of the estimates; and (iv) A range of MR methodologies was employed to exclude heterogeneity and horizontal pleiotropy. All statistical analyses were conducted using the “TwoSampleMR”, “MRPRESSO”, and “MendelianRandomization” packages in R software (version 4.3.3). Based on Bonferroni correction, the threshold for statistical significance was set at , and only results meeting this threshold were considered statistically significant.⁴¹ All p-values are two-tailed.

Results

Genetic Instrumental Variable Selection

In our study, a total of 658 LD-independent and suitable IVs were chosen from GWASs for CS, and finally extracted 81 SNPs for TB-BMD, 19 SNPs for FN-BMD, 3 SNPs for FA-BMD, and 24 SNPs for LS-BMD, 494 SNPs for HB-BMD, 8 SNPs for TB-BMD (0–15), 1 SNPs for TB-BMD (15–30), 9 SNPs for TB-BMD (30–45), 19 SNPs for TB-BMD (45–60), and 22 SNPs for TB-BMD (over 60), respectively. The MR results of this section were based on the SNPs screened under the genome-wide significance threshold (). These SNPs were derived from GWAS meta-analysis data, excluding linkage disequilibrium (LD), and corresponded to an for every single variant. Thus, the selected IVs are robust, and “weak IVs” are negligible. The list of SNPs remaining was used for the subsequent MR analysis, as shown in Supplementary Tables 1–10.

Univariate Mendelian Randomization

According to the IVW estimator, genetically higher levels of TB-BMD (OR = 1.22, 95% CI: 1.10–1.34, P = 7.33E−05), FN-BMD (OR = 1.13, 95% CI: 1.02–1.25, P = 2.25E−03), LS-BMD (OR = 1.15, 95% CI: 1.05–1.27, P = 2.38E−03), and HB-BMD (OR = 1.09, 95% CI: 1.03–1.15, P = 2.15E−03) were significantly associated with an increased risk of CS. Other MR methods indicated a consistent direction of beta values, including MR-Egger (TB-BMD: OR = 1.13, 95% CI: 0.87–1.46, P = 3.71E−01; FN-BMD: OR = 1.70, 95% CI: 1.00–2.90, P = 6.58E−02; LS-BMD: OR = 0.83, 95% CI: 0.61–1.12, P = 2.31E−01) and weighted median (TB-BMD: OR = 1.15, 95% CI: 0.99–1.33, P = 6.40E−02; FN-BMD: OR = 1.15, 95% CI: 0.99–1.32, P = 6.94E−02; LS-BMD: OR = 1.14, 95% CI: 1.00–1.29, P = 3.98E−02; HB-BMD: OR = 1.04, 95% CI: 0.96–1.13, P = 3.29E−01), which can be interpreted as a positive result.⁴² After the Bonferroni correction, IVW method showed that only in the population of 45–60 years old, TB-BMD (per 1 SD increase) was potentially associated with a 14.8% increase in CS (OR = 1.16, 95% CI: 1.08–1.26; P = 1.84E-04). The weighted median analysis yielded similar significant results (OR = 1.13, 95% CI: 1.01–1.26; P = 3.29E−02), while the MR Egger analysis showed consistent but non-significant results (OR = 1.03, 95% CI: 0.77–1.38; P = 8.46E−01). All of the results are shown in Figure 2.

Figure 2 Estimation of the causal relationship between BMDs and CS using different MR methods (Bonferroni-adjusted p-values highlighted in red indicate significant associations at P < 0.005). Abbreviations: CS, cervical spondylosis; TB-BMD, total body bone mineral density; FN-BMD, femoral neck bone mineral density; LS-BMD, lumbar spine bone mineral density; FA-BMD, forearm bone mineral density; HB-BMD, estimated heel bone mineral density; CI, confidence interval; nSNP, number of single nucleotide polymorphism; OR, odds ratio.

Sensitivity Analysis

As shown in Table 2, the results of MR-PRESSO ensured that there were no outliers in the selected IVs (TB-BMD: P = 1.12E-01; FN-BMD: P = 4.60E-01; TB-BMD (45–60): P = 3.53E-01), when MR-Egger Intercept (TB-BMD: P = 5.37E-01; FN-BMD: P = 1.39E-01; TB-BMD (45–60): P = 4.10E-01) demonstrated the absence of horizontal pleiotropy in both studies. When testing HB-BMD selected IVs, we found P < 1.00E-03 even after eliminating all outliers. Consequently, an additional MR-Egger pleiotropy test was conducted on this study and indicated there was no horizontal pleiotropy. However, there existed a significant horizontal pleiotropy between IVs and outcomes which found by the MR-Egger pleiotropy test (P = 3.40E-02) between LS-BMD and CS. Heterogeneity was hypothesized in the study of TB-BMD (MR Egger: P = 9.10E-02; IVW: P = 9.80E-02), FN-BMD (MR Egger: P = 5.44E-01; IVW: P = 4.47E-01), HB-BMD (MR Egger: P = 6.10E-02; IVW: P = 6.30E-02), TB-BMD (45–60) (MR Egger: P = 3.26E-01; IVW: P = 3.42E-01), which may have effect on the risk of CS. And there was no proof to indicate the heterogeneity existing in the findings. The scatter plots for effect sizes of SNPs for BMDs on CS, include TB-BMD, FN-BMD, HB-BMD, TB-BMD (45–60) are depicted in Figure 3. And, the funnel plot displayed a mostly symmetrical distribution in Supplementary Figure 1, indicating the robustness of MR results. Finally, we present the results of the leave-one-out analysis in Supplementary Figure 2, which demonstrated that all of the identified relationships in the TB-BMD, FN-BMD, HB-BMD and TB-BMD (45–60) were not altered by a single SNP. Moreover, Supplementary Figures 3–9 presented forest plots of individual SNP effect analysis for BMDs by sites and ages on CS.

Table 2 Sensitivity Analysis of Mendelian Randomization Studies of BMD and CS

Figure 3 Scatter plots for the MR analysis, which explore the causal effect of BMDs on CS. (A) The scatter plot for TB-BMD on CS. (B) The scatter plot for FN-BMD on CS. (C) The scatter plot for HB-BMD on CS. (D) The scatter plot for TB-BMD (over 60) on CS. Abbreviations: CS, cervical spondylosis; TB-BMD, total body bone mineral density; FN-BMD, femoral neck bone mineral density; HB-BMD, heel bone mineral density.

Multivariate Mendelian Randomization

To adjust for potential pleiotropic pathways that might confound the relationship between BMDs and CS, a MVMR model was utilized, incorporating adjustments for body mass index (BMI), type 2 diabetes mellitus (T2DM), hyperlipidemia and cigarette per day (CPD). The instrumental variables selected for MVMR implementation are presented in Supplementary Table 11. Multivariate IVW MR methods were used to estimate causality, which revealed that TB-BMD and HB-BMD still had a direct causal influence on CS (Adjusting for TB-BMD: OR = 1.14, 95% CI: 1.01–1.29; P = 3.12E −02; Adjusting for HB-BMD: OR = 1.11, 95% CI: 1.01–1.22; P = 2.45E −02). Figure 4 displayed strong evidence that genetically predicted TB-BMD and HB-BMD had a causal effect on CS. In contrast, the results of the MVMR IVW analysis revealed that TB-BMD (45–60) and FN-BMD did not have a direct causal effect on the risk of CS as previously predicted by UNMR (Adjusting for TB-BMD (45–60): OR = 1.03, 95% CI: 0.92–1.17; P = 5.93E −02; Adjusting for FN-BMD: OR = 0.93, 95% CI: 1.78–1.12; P = 4.62E −02).

Figure 4 Forest plots of MVMR analyses examined causal associations between BMD and CS, adjusted for, T2DM, BMI, CPD, and hyperlipidemia. Abbreviations: BMI, body mass index; T2DM, type 2 diabetes mellitus; CPD, cigarette per day; TB-BMD, total body bone mineral density; FN-BMD, femoral neck bone mineral density; HB-BMD, estimated heel bone mineral density; CI, confidence interval; nSNP, number of single nucleotide polymorphism; OR, odds ratio.

Reverse Mendelian Randomization

In the reverse Mendelian randomization (MR) analysis, we carefully screened 8 IVs with genome-wide significance levels () of LD-independent () from the GWASs for CS (Supplementary Table 13). The F-statistics of the IVs for the selected CS were all greater than 10 to ensure that the selected IVs were sufficiently robust enough to eliminate the potential bias. Meanwhile, no outlier IVs were found in the MR-PRESSO analysis. Details for SNPs excluded were shown in Supplementary Table 12. Our investigation has revealed that the CS is causally associated with a risk of LS-BMD and HB-BMD (LS-BMD: OR = 1.12, 95% CI: 1.03–1.22, P = 1.02E −02; HB-BMD: OR = 0.96, 95% CI: 0.93–0.99, P = 1.13E −02). The results from Cochran’s Q test were as Supplementary Table 14 displayed: LS-BMD (MR Egger: P = 0.386; IVW: P = 4.47E −01); HB-BMD (MR Egger: P = 0.996; IVW: P = 9.81E −01). And MR-Egger regression results did not demonstrate a directional pleiotropy effect amongst all genetic variants (LS-BMD: P = 5.21E −01; HB-BMD: P = 8.77E −01). However, the results of both two BMD levels in IVW did not reach statistical significance after applying Bonferroni correction (significance threshold set at P < 0.005 (0.05/10)), thus we did not perform any other MR pleiotropy tests.

Discussion

The correlation between BMD and CS has long been a focus of medical research, yet definitive conclusions remain absent. To the best of our knowledge, this study represents the first MR investigation to explore the causal relationship between BMD and CS.⁴³ Utilizing large-scale GWAS data, we found that elevated TB-BMD, FN-BMD, HB-BMD, and TB-BMD (45–60 years) exhibit direct causal effects on CS after Bonferroni correction. MVMR analysis further revealed that TB-BMD and HB-BMD maintain independent and significant causal associations with CS even after adjusting for confounding factors.

From a public health perspective, we recommend routine BMD screening for adults aged 45–60 years to assess CS risk. For those identified with high BMD levels, appropriate preventive strategies including postural management improvement and neck core muscle strengthening should be implemented to reduce cervical mechanical load. Notably, although MVMR analysis has demonstrated that BMD independently promotes CS progression regardless of confounders, clinical interventions must still consider systemic factors such as basal metabolic status. For instance, Hong et al found that individuals with higher BMD exhibit elevated BMI, attributing this phenomenon to compensatory forward head posture caused by airway volume reduction due to adipose accumulation. This postural alteration disrupts the correspondence between cervical BMD distribution and load distribution, thereby impacting vertebral strength and health.¹⁴ Thus, health education for at-risk populations should not focus solely on bone metabolism and vertebral quality but should also emphasize comprehensive improvements in systemic health, including appropriate aerobic exercise and strength training to facilitate weight management, enhance circulation, and promote nutrient supply.

BMD serves as a critical indicator of bone strength, yet high bone density does not necessarily equate to optimal bone quality or microstructure.⁴⁴ Our findings align with several previous cross-sectional and retrospective studies. For example, Hong et al reported in a 2021 cross-sectional study of 116 patients that postmenopausal women with high BMD may experience greater neck pain severity and accelerated cervical degeneration accompanied by forward head posture alterations.⁴⁵ Through a micro-computed tomography (μCT) scanning system, Wang et al found that greater vertebral BMD correlated with more severe disc degeneration.⁴⁶ A study of participants under 64 years old yielded consistent results, showing that individuals with progressive intervertebral disc degenerative changes tend to exhibit elevated BMD in the lumbar spine and hips.⁴⁷ To better observe and measure bone morphology, microstructure, and mineral heterogeneity in aging zebrafish spines, Erika et al suppressed bone resorption via cathepsin K regulation, identifying a significant association between increased BMD and disc degeneration.²⁰ Contrary to conventional views, middle-aged adults with high BMD are more likely to impact adjacent discs through increased nucleus pulposus pressure caused by compact vertebral bone.⁴⁸ In patients with OP, having lower BMD actually helps nutrients reach the intervertebral discs more easily. This happens through several pathways: endplate vascularization, reduced endplate resistance and decreased intradiscal strain.¹⁷ This occurs because intervertebral discs are largely avascular, with only minimal capillaries present in the outer annulus. Consequently, their primary nutrient supply relies on diffusion through the upper and lower cartilage endplates. Elevated BMD increases endplate static pressure, thereby impairing the diffusion of essential nutrients (eg, glucose) into discs.⁴⁹ Helen et al further validated this hypothesis by assessing vertebral endplate calcification grades in aged gerbils through BMD analysis. They demonstrated that higher BMD correlates with increased endplate calcification tendency, resulting in diminished vascular supply to discs, whereas lower BMD promotes nutrient diffusion from disc peripheries, slowing degeneration progression.⁵⁰

Given the biological mechanisms linking BMD to CS, future clinical practice should consider direct measurement of cervical BMD as a primary tool for predicting CS progression. Compared to other skeletal regions, this approach maximizes prediction accuracy. However, the complex anatomical structure of cervical spine imposes higher technical demands on BMD measurement. Although DXA remains the gold standard for BMD assessment, it is not entirely suitable for evaluating CS susceptibility due to potential measurement errors caused by extensive degeneration such as posterior longitudinal ligament calcification, ligamentum flavum hypertrophy or calcification, and nodular changes in advanced CS stages.⁵¹ Statistics indicate that preoperative DXA screening rates for spinal fusion surgeries are as low as 12%, with even fewer cervical surgery patients receiving preoperative DXA evaluations.^52,53 Compared to DXA, quantitative computed tomography (QCT) enables more accurate three-dimensional analysis but faces limitations due to radiation exposure and cost.⁵⁴ The lack of effective detection methods partially explains the paucity of research in this field over past decades. Recently, the MRI-based vertebral bone quality (VBQ) score derived from T1-weighted imaging has gained attention for its feasibility in assessing cervical bone quality, demonstrating strong correlation with DXA-measured lumbar BMD.⁵⁵ In a recent retrospective study, Chen et al innovatively analyzed the relationship between cervical VBQ scores and QCT-measured BMD, highlighting VBQ as a quantitative indicator of cervical BMD in CS patients. Additionally, this study revealed that in CS patients over 50 years old, increased BMD in both caudal vertebral bodies and sub-endplate regions significantly correlates with aggravated cervical disc degeneration.¹⁹ Currently, VBQ scoring has been widely adopted for preoperative cervical BMD assessment to predict postoperative outcomes such as cage subsidence and titanium mesh complications.^56,57 This clinically efficient and safe method shows promise as a screening tool for CS susceptibility in postmenopausal women and elderly males, aiding diagnostic and early therapeutic decision-making with broad applications in our aging global population.

Despite substantial evidence supporting our conclusions, conflicting perspectives persist. In a retrospective investigation, Astrid et al used QCT to scan BMD in at least three thoracolumbar segments of 134 patients, concluding that degenerative spinal changes may correlate with reduced regional spinal mineralization.⁵⁸ Liang et al measured Hounsfield unit (HU) values in C3-C7 vertebrae of 324 degenerative CS patients, finding the highest HU values in the upper C4 vertebra with gradual decreases toward C3 and C7. They hypothesized that vertebral BMD reduction might trigger or exacerbate adjacent disc degeneration.⁹ The cause of this alteration might involve insufficient nutrient replenishment of the disc cells, as well as poor bone quality, resulting in thinning of the endplates and an increase in microfractures.⁵⁹ It must be emphasized that these studies are observational in nature, limited by their inability to establish causality. Furthermore, compared to MR analysis, observational studies exhibit weaker control over confounders even with strict standardization protocols, whereas MR effectively addresses these issues.⁶⁰ We found that several prior MR studies exploring the relationship between BMD and intervertebral disc degeneration reached conflicting conclusions, which may stem from heterogeneous data sources involving different spinal segments.^10,11,61 Moreover, extrapolating thoracolumbar findings to the cervical spine risks erroneous conclusions. Jacob et al demonstrated through analysis of anthropometric and regional VBQ measures that thoracic VBQ scores directly reflect lumbar vertebral bone quality, whereas cervical measurements do not.⁶² In our UVMR, LS-BMD also exhibited significant horizontal pleiotropy with CS. Consequently, independent investigation of the causal relationship between BMD and CS remains essential.

In summary, the strengths of this study are as follows: First, we utilized the publicly available GWAS datasets to examine causal relationships between BMD (covering five age groups and four skeletal sites) and CS, representing the largest sample size in related research to date. Simultaneously, we mitigated the impact of confounding factors and reverse causation on our conclusions. Second, by restricting our analysis to individuals of European ancestry, we minimized potential biases from population stratification. Third, our multi-database sampling strategy enhanced geographical diversity while avoiding statistical errors from sample overlap. Fourth, as shared covariates for both OP and CS, age and sex were adjusted for in multivariable statistical models. Our study further stratified BMD data by age and conducted subgroup analyses to assess their effects on CS.

Nonetheless, limitations remain. Initially, as only European populations were analyzed due to the absence of large-scale GWAS data from other ethnicities, generalization of results to other populations requires caution. Secondly, adjusting for confounders in the MVMR analysis reduced the number of available IVs, potentially increasing false-negative findings. The stringent Bonferroni correction may have further limited the detection of true positive associations. These limitations highlight the need for future MR studies with larger sample sizes or RCTs. Finally, more comprehensive mechanistic investigations are required to elucidate the precise biological pathways linking BMD and CS before clinical translation.

Conclusion

This study conclusively demonstrates that elevated BMD levels negatively impact cervical spine health, and provides novel genetic insights into predicting and treating CS. Overall, our findings indicate that individuals with high TB-BMD, FN-BMD, and HB-BMD levels—particularly those aged 45–60 years—should be vigilant for the development of CS. These results lay a foundation for future research, and the positive causal relationships are supported by existing evidence. Further studies are required to explore CS prevention via BMD modulation and to develop personalized therapeutic strategies for CS patients.