This investigation was conducted as a retrospective cohort study, utilizing de-identified, and patient-level data obtained from the TriNetX platform. TriNetX is a federated global health research network that aggregates real-world electronic health records (EHRs) from hospitals and healthcare systems across the United States. The analysis was executed within the US Collaborative Network on October 11, 2025, which encompasses 71 healthcare organizations providing comprehensive longitudinal clinical data. All data extraction and analyses were performed exclusively within the secure TriNetX cloud environment. Given that all data were fully de-identified in strict compliance with HIPAA and GDPR, the requirement for institutional review board approval and informed patient consent was waived by the platform. Nevertheless, the study protocol received independent review and approval from the Taipei Tzu Chi Hospital Institutional Review Board (Approval Number: 14-IRB134) and adhered to the ethical principles of the Declaration of Helsinki.

The initial population comprised adults (aged 18 years or older) with a documented glaucoma diagnosis (ICD-10-CM codes H40–H42, specifically including primary open-angle H40.1 and primary angle-closure H40.2) recorded between 1 January 2005 and 1 January 2020. The initial screen identified 761,232 patients. Cases were subsequently restricted to first-occurrence glaucoma without secondary causes by excluding any secondary glaucoma codes (H40.3–H40.6) recorded within 1 year on or before the index glaucoma diagnosis, resulting in 722,862 eligible patients.

Vitamin D status was determined based on serum 25-hydroxyvitamin D [25(OH)D] concentrations identified by LOINC code 1989-3. Measurements obtained within 1 year prior to the first glaucoma record were included, and for patients with multiple results, the value closest to the index date was selected to best reflect pre-observation physiological status. To ensure stable baseline exposure, individuals with any historical 25(OH)D value < 30 ng/mL were excluded from the VDA cohort, and those with any prior value ≥30 ng/mL were excluded from the VDD cohort, minimizing potential misclassification due to recent status changes. Over 92% of matched participants in both cohorts had complete baseline 25(OH)D data with minimal missing values. Patients were categorized into two exposure groups: VDD (< 30 ng/mL) and VDA (≥30 ng/mL), with 16,602 and 15,638 participants, respectively, prior to matching.

The index date was set as the earliest date a patient met all the criteria for a glaucoma diagnosis. The observational follow-up period commenced 1 day after the index date, with a maximum duration of 1,825 days (5 years), censored by the event of death or the end of available health records, whichever occurred first. Any patient with evidence of the primary outcome prior to the start of the observation window was excluded from the analyses to ensure the capture of incident outcomes.

To mitigate selection bias and minimize the influence of confounding factors, a 1:1 PSM procedure was implemented using a greedy nearest-neighbor algorithm with a caliper width of 0.1. Matching covariates included 46 clinically relevant variables, such as demographic features, comorbid cardiometabolic conditions, medication use, renal and metabolic laboratory indices, and inflammatory biomarkers. Crucially, the exposure-defining variable (serum 25(OH)D levels) was excluded from the covariate list to maintain the primary exposure contrast. First, the standardized mean differences (SMDs) were calculated post-matching to confirm acceptable covariate balance, with an SMD < 0.1 considered indicative of comparability. Second, the overlap of propensity score distributions (common support) was visually inspected to ensure that the matched cohort shared similar characteristics across the entire range of the propensity scores (Supplementary Figure S1). The final pre-matched cohorts comprised 16,602 patients in VDD and 15,638 in VDA. Following PSM, both the VDD and VDA cohorts contained 10,881 individuals.

The primary endpoint was the incidence of neurodegenerative disorders during the 5-year follow-up period, defined as a composite outcome comprising the following ICD-10-CM codes: unspecified dementia (F03), Alzheimer's disease (G30), and Parkinson's disease (G20). Unspecified dementia captured clinically confirmed dementia without subtype specification, Alzheimer's disease included all sub codes representing Alzheimer's pathology—the most common specific dementia subtype—and Parkinson's disease represented primary Parkinsonism. These conditions were selected for their clinical relevance and high prevalence in U.S.-based EHR data, with Alzheimer's disease and unspecified dementia together accounting for most dementia diagnoses and providing a robust basis for evaluating associations with vitamin D status.

All statistical computations, including the estimation of absolute risks, risk ratios, and odds ratios, were executed using the integrated Risk Analysis features of the TriNetX platform. The cumulative incidence of time-to-event outcomes was visualized using Kaplan–Meier survival curves, and differences between the VDD and VDA groups were tested using the Log-Rank test. Cox proportional hazards models were subsequently employed to calculate the HRs and their associated 95% Confidence Intervals (CIs) for incident outcomes. Statistical significance was predefined by a two-sided p-value less than 0.05. Primary statistical analyses, including propensity score matching (PSM) and Cox proportional hazards regression, were performed within the TriNetX (Cambridge, MA, USA) analytics platform. For data visualization and figure generation, we utilized GraphPad Prism (version 8.0.1) and Microsoft Excel.

A comprehensive set of sensitivity analyses and one exploratory analysis were performed. We first conducted a Landmark Analysis to assess the time-varying effect of VDD on cumulative incidence and HRs at 1, 3, and 5 years post-index. For Confounding Assessment, we compared the unadjusted HRs with the primary PSM-matched HRs, which quantified the impact of controlling for confounders. We then performed a Severe Deficiency Adjustment by re-analyzing the data using a stricter classification for VDD, defined as 25(OH)D levels < 20 ng/mL, to investigate a potential dose-response relationship. To evaluate the potential mitigating effect of systemic anti-inflammation, we performed a non-steroidal anti-Inflammatory drugs (NSAID) confounding assessment by excluding individuals with recorded long-term NSAID use. Finally, an exploratory healthcare utilization analysis was performed by comparing the mean number of hospital visits between the matched VDD and VDA cohorts to assess baseline differences in health-seeking behavior or comorbidity burden.

To evaluate the potential influence of detection bias and systemic health disparities on study outcomes, we performed supplementary analyses of healthcare utilization over the 5-year follow-up period. Specifically, we compared the VDD and VDA cohorts in terms of hospital visit counts (Supplementary Table S5), frequency of all-cause hospitalizations (Supplementary Table S6), and utilization of ophthalmology-specific services and procedures (Supplementary Table S7), using independent t-tests and Chi-square tests. Finally, to address potential selection bias from missing data, we conducted a Restricted PSM sensitivity analysis by limiting matching covariates to 26 high-completeness variables (>85% data availability), ensuring that the observed association was not sensitive to less complete laboratory parameters (Supplementary Table S8).

As shown in Figure 1, the derivation of the final study cohort commenced with the identification of 761,232 adult patients (aged 18 years or older) diagnosed with glaucoma (ICD-10 codes H40–H42). After applying exclusion criteria for secondary glaucoma, 722,862 individuals remained. Vitamin D status was determined by serum 25-hydroxyvitamin D (25(OH)D) levels: 16,602 individuals met the criteria for VDD (25(OH)D < 30 ng/mL), and 15,638 for VDA (25(OH)D ≥ 30 ng/mL).

To ensure unbiased comparisons, a comprehensive 1:1 PSM was performed using 46 covariates (excluding the exposure-defining variable, serum 25(OH)D), successfully creating two well-balanced cohorts, each containing 10,881 patients. As shown in Table 1, the matching process successfully resolved marked baseline imbalances, including the age difference (VDD 61.5 ± 11.3 years vs. VDA 61.6 ± 11.5 years; Standardized Mean Difference, SMD, SMD = 0.002). Furthermore, the propensity score density plots demonstrated near-perfect overlap between the matched cohorts, confirming robust common support (Supplementary Figure S1). These finalized cohorts were then followed for up to 5 years to assess the association of neurodegenerative outcomes.

Kaplan–Meier survival analyses demonstrated a significant difference in the primary outcome, unspecified dementia, between glaucoma patients with VDD and those with VDA over the 5-year follow-up period, as shown in Figure 2A.

VDD was associated with a higher risk of unspecified dementia, with a significantly lower survival probability throughout the observation window (Log-Rank chi-squared = 7.742; p = 0.005) and an increased hazard of disease (HR = 1.241; 95% CI: 1.066–1.446) (Figure 2A). The absolute 5-year risk was 3.2% in the VDD cohort vs. 2.6% in the VDA cohort (Supplementary Table S1).

By contrast, there was no statistically significant difference in the incidence of AD (Figure 2B) or PD (Figure 2C). For AD, the association was non-significant (HR 1.224; 95% CI: 0.965–1.552; p = 0.095). Similarly, PD incidence was nearly identical between the two groups (HR 1.066; 95% CI: 0.805–1.413; p = 0.655).

Landmark Kaplan–Meier analyses comparing glaucoma patients with VDD and VDA were conducted at 1, 3, and 5 years following the index diagnosis (Table 2).

The excess risk of dementia associated with VDD strengthened over the follow-up period. At 1 year, the difference was not yet statistically significant (HR 1.323, 95% CI 0.963–1.819; p = 0.083). However, the elevated risk became significant at the 3-year landmark (HR 1.396; p = 0.001) and remained significant at 5 years (HR 1.241; 95% CI 1.066–1.446; p = 0.005).

Across all landmark analyses, the incidence of AD and PD did not show statistically significant survival differences, indicating that the impact of VDD is primarily associated with the broader dementia diagnosis and the effect is most pronounced in the mid-to-long term.

Due to limited case counts for specific outcomes, the subgroup analysis (Figure 3) and effect modification assessment were performed exclusively for unspecified dementia. Analysis across patient subgroups revealed that the 5-year dementia risk associated with VDD was markedly influenced by the presence of specific baseline characteristics.

Advanced age (65 years or older vs. 18–64 years) emerged as the most significant predictor (HR 6.84, 95% CI 5.15–9.08; p = 0.0001). Other major risk factors included hypertension (HR 2.70, 1.64–4.44; p < 0.0001), hypoalbuminemia (albumin < 4 vs. ≥4 g/dL: HR 2.64, 2.08–3.36; p < 0.0001), elevated inflammatory burden (CRP level 10 mg/L or greater vs. less than 10 mg/L; HR 1.38, 1.04–1.84; p = 0.026), and diabetes (HR 1.29, 1.03–1.92; p = 0.027). Conversely, higher BMI (30 kg/m² or greater vs. less than 30 kg/m²) was associated with a paradoxical lower observed risk (HR 0.59, 0.44–0.77; p < 0.001). Glaucoma subtype was not a significantly associated with dementia risk (HR = 0.814; 95% CI: 0.63–1.05; p = 0.115).

Sensitivity analyses confirmed the robustness of the primary findings under varied conditions. Before propensity score matching (PSM), VDD was associated with a lower hazard of unspecified dementia (HR = 0.871; 95% CI: 0.762–0.995; P = 0.042), which reversed after matching, showing a higher hazard (HR = 1.241; 95% CI: 1.066–1.446; P = 0.005) (Supplementary Table S2). Applying a stricter deficiency threshold (25[OH]D < 20 ng/mL) revealed a stronger association with dementia risk (HR = 1.493; 95% CI: 1.179–1.890; P < 0.001) than the primary threshold of < 30 ng/mL (Supplementary Table S3). Excluding long-term NSAID users attenuated the association, rendering it non-significant (HR = 1.178; 95% CI: 0.998–1.391; P = 0.053) (Supplementary Table S4). Finally, healthcare utilization during follow-up was examined (Supplementary Tables S5–S7). After matching, the VDA cohort had a higher baseline visit frequency (120.49 vs. 110.86 visits; P < 0.0001; Supplementary Table S5), while 5-year hospitalization rates were similar between groups (1.93 vs. 1.89 visits; P = 0.641; Supplementary Table S6). The VDA group also exhibited more ophthalmology-related encounters and procedures than the VDD group (3.17 vs. 2.87 visits; P < 0.0001; Supplementary Table S7).

The restricted sensitivity analysis confirmed that the hazard ratio remained stable when matching was limited to high-completeness variables (HR 1.239, p = 0.008; Supplementary Table S8).

The association between VDD and increased dementia risk is highly consistent with the known pleiotropic, non-skeletal actions of vitamin D. 1,25-dihydroxyvitamin D functions as a crucial steroid hormone by binding the VDR, which is widely expressed throughout the CNS, including the hippocampus, cortex, and key ocular neuroglial elements like retinal ganglion cells and optic nerve head astrocytes (19, 28, 31).

Vitamin D supports neural protection by aiding amyloid clearance and reducing neuronal toxicity, while deficiency may promote vascular dysfunction and brain atrophy associated with cognitive decline (27, 32) In neurological tissues, active vitamin D signaling regulates processes critical to neurodegeneration by down-modulating chronic inflammation and oxidative stress, such as suppressing the NF-kappa-B signaling pathway (18, 33). Since both glaucoma and dementia involve chronic low-grade neuroinflammation (34–37), VDD may act as an inflammatory sensitizer, compounding neuronal vulnerability and accelerating CNS decline (18, 38).

The inflammatory mechanism proposed is strongly supported by the sensitivity analysis focusing on NSAID use. Long-term NSAID use is epidemiologically linked to reduced dementia risk, and the protective effect appears more pronounced for NSAIDs without amyloid-β lowering properties, and cumulative dose does not play a major role (39, 40). The finding that the association between vitamin D deficiency and dementia risk diminishes to non-significance after excluding long-term NSAID users underscores a potentially complex relationship. It suggests that the anti-inflammatory effects of NSAIDs might not effectively reduce or could even influence the excess dementia risk linked to vitamin D deficiency, possibly acting as a systemic inflammation proxy or marker (41, 42). This attenuation indicates that long-term NSAID use may confound or modify the relationship between vitamin D status and dementia, emphasizing the role of systemic inflammation in dementia pathogenesis. Limited data exist on how vitamin D deficiency and NSAID use interact in relation to dementia risk, with no strong evidence for a synergistic or antagonistic effect. Current evidence suggests vitamin D deficiency raises the risk of dementia, while long-term NSAID use may modestly reduce it, likely via anti-inflammatory mechanisms (41, 42). However, NSAIDs are not recommended clinically for dementia prevention due to their potential adverse effects and inconsistent findings across studies, with research often limited to prescription use and possibly missing over-the-counter consumption (39).

Subgroup analysis revealed that the dementia risk associated with VDD is highly dependent on underlying systemic health. Markers of advanced age (HR 6.84), hypertension (HR 2.70), and hypoalbuminemia (HR 2.64) were powerful predictors, suggesting that VDD may function primarily as an indicator of an overall vulnerable systemic profile rather than the sole driver of pathology (28, 43).

Also, the elevated risk was also pronounced in patients with markers of systemic inflammation: elevated CRP (≥10 mg/L) (HR 1.38) and diabetes (HR 1.29), reinforcing that VDD compounds risk in individuals with high inflammatory and metabolic burdens, pointing toward a shared pathway involving microvascular compromise (44–46). The finding that higher BMI (≥30 kg/m²) was associated with a paradoxical lower observed risk (HR 0.59) must be interpreted with caution. Potential explanations include reverse causation, where pre-clinical dementia leads to weight loss, or that higher BMI may reflect better underlying nutritional reserves acting as a buffer against severe systemic decline and frailty (47, 48).

On the contrary, our analysis showed that glaucoma subtype (PACG vs. POAG) did not significantly influence the hazard for dementia (HR = 0.814, p = 0.115). This contrast emphasizes that systemic metabolic and inflammatory factors, rather than local ophthalmic classifications, are the primary drivers of neurodegenerative outcomes in this population.

The decision to focus on “unspecified dementia” is validated by clinical practice: many patients initially present with cognitive decline that cannot be confidently classified, and thus “unspecified” or mixed dementia diagnoses account for the majority of incident cases (49, 50). These early, heterogeneous presentations provide greater statistical power. Conversely, the lower event rates for specific neurodegenerative subtypes (AD and PD) in our cohort resulted in limited power to detect subtype-specific associations.

Although comprehensive PSM accounted for multiple covariates, differences in healthcare utilization persisted. The VDA cohort showed higher frequencies of hospital visits and ophthalmic services (Supplementary Tables S5 and S7), which may reflect an underlying tendency for greater healthcare engagement or a different comorbidity profile—factors often cited as sources of unmeasured confounding in large-scale EHR studies (51–53). However, rather than weakening the results, this pattern strengthens their validity. Greater clinical contact in the VDA group would typically increase the chance of dementia detection, yet the VDD cohort still exhibited a higher dementia risk. This suggests that detection bias is unlikely to explain the association and, if present, would have reduced rather than exaggerated the observed effect. Comparable hospitalization rates between cohorts (Supplementary Table S6) further indicate that the findings are not driven by greater systemic disease burden in the VDD group. Consistent results across sensitivity analyses reinforce the robustness of these findings. When limiting the propensity score model to 26 high-completeness variables (Supplementary Table S8), the estimated hazard ratio for dementia remained virtually unchanged, confirming that the association between vitamin D deficiency and dementia risk is stable and not materially influenced by data completeness or model specification.

This retrospective EHR study has inherent limitations precluding causal inference, including reliance on ICD-10 codes prone to inaccuracies, a single baseline 25(OH)D measure omitting longitudinal changes, and absence of granular ophthalmic staging (e.g., visual field indices, cup-to-disc ratios) despite no glaucoma subtype-dementia interaction (P = 0.115). The “unspecified dementia” endpoint limits etiological specificity; however, supplementary analyses refute detection bias from differential healthcare use, as the VDA group had more ophthalmology encounters (P < 0.0001). Although 1:1 propensity score matching and no multiple imputation were used, excellent covariate balance (all SMDs < 0.1), large sample size, and stable sensitivity analyses support robust, clinically relevant findings despite potential generalizability constraints from the federated EHR platform.