Fluoride exposure disrupts creatinine homeostasis and undermines the reliability of urinary metabolomics normalization.

2.1. Study population

2.1.1. NHANES population study

Data were obtained from the 2015–2016 cycle of NHANES, which was approved by the Research Ethics Review Board of the National Center for Health Statistics (NCHS), U.S. Centers for Disease Control and Prevention. All participants provided written informed consent. A total of 9971 individuals participated in NHANES during 2015–2016. We excluded 7653 participants without urine fluoride (UF) measurements, as UF concentration was the primary exposure variable, and missing UF data could lead to exposure misclassification. According to NHANES documentation, UF testing was conducted in a subsample of participants, and missingness occurred at random with respect to age, sex, and major health characteristics; therefore, this exclusion was unlikely to introduce systematic bias. We further excluded 29 participants with abnormal urinary albumin-to-creatinine ratios and one participant with an extreme UF value well beyond the population distribution. The final analytic sample comprised 2408 participants with valid, non-zero sample weights. Demographic information—including age, sex, race/ethnicity, and the income-to-poverty ratio—was collected through household interviews. Height and weight were measured during the mobile examination center visit. Body mass index (BMI) was calculated as weight in kilograms divided by height in meters squared. Missing BMI values (n = 30) were imputed using sex- and age-specific median values. For missing income-to-poverty ratio data, values were imputed based on household income questionnaire responses. For individuals with missing responses (n = 149), the overall median value was used; for the remaining cases (n = 77), the median of the corresponding income group was applied. A flowchart depicting the participant selection process is provided in Fig. 1.

1. Download: Download high-res image (150KB)
2. Download: Download full-size image

Fig. 1. Participant selection flow chart.

2.2. Measurement of UF, CR, urine specific gravity, and urine osmolality

In the NHANES population, UF was measured using an ion-selective electrode (ISE) method, conducted by the Inorganic and Radiation Analytical Toxicology Branch (IRATB) of the Division of Laboratory Sciences, National Center for Environmental Health. The limit of detection (LOD) for UF was 0.144 mg/L; values below the LOD were replaced with LOD/?2. Approximately 94 % of participants had UF levels above the LOD. Creatinine was measured using an enzymatic endpoint method performed by the University of Minnesota, with a reportable limit of 5 mg/dL. In the elderly Chinese population, UF was measured according to the national standard of China, “Determination of Fluoride in Urine by Ion-Selective Electrode Method” (WS/T 89–2015), at Harbin Medical University. The linear detection range for UF was 0.1–10 mg/L. Creatinine was measured using a creatininase oxidase method, also conducted by Harbin Medical University, with a detection limit of 30 µmol/L. Urine specific gravity was measured using reflectance photometry and completed by China Servicebio Company. Urine osmolality was measured using the freezing-point depression method and completed by China Guangce Research Institute.

2.3. Urine sample preparation and metabolomics profiling in the Chinese older adult population

In this study, a 1:1 case–control design was applied based on age and sex to compare participants with and without skeletal fluorosis. A total of 25 cases and 25 matched controls were selected for metabolomic analysis. Untargeted metabolomic analysis was conducted by Calibra Lab at DIAN Diagnostics (Hangzhou, Zhejiang, China) on their CalOmics metabolomics platform.

All urine samples were immediately aliquoted after collection and stored at -80 °C until analysis to prevent metabolite degradation. For sample preparation, each urine sample was extracted with methanol at a ratio of 1:4 (urine: methanol), vortexed for 3 min, and centrifuged at 4000 × g for 10 min at 20 °C. The supernatant was collected, divided into four aliquots of 100 uL each, dried under a stream of nitrogen, and reconstituted in an appropriate solvent prior to instrumental analysis.

Multiple deuterated compounds were used as internal standards to correct for signal drift and ensure quantitative accuracy. The internal standard mixture—composed primarily of deuterated amino acids, fatty acids and sugars (e.g., d3-leucine, d35-octadecanoic acid, d7-glucose)—was added in equal amounts to every sample and to the pooled quality control (PQC) samples. The retention time and signal intensity of each internal standard were continuously monitored throughout the analytical run to ensure instrument stability; when abnormal drift was detected, corrective quality-control actions (e.g., re-sampling and re-analysis) were implemented.

Rigorous quality-control procedures were applied: a pooled QC sample, prepared by proportionally mixing aliquots from all study samples, was injected once every ten study samples to monitor system stability; process blanks were included to detect potential contamination. Non-targeted detection was performed on a UPLC–MS/MS platform (Waters ACQUITY 2D UPLC coupled with Thermo Q Exactive) using four chromatographic–mass spectrometric methods (POSa, POSb, NEGa, NEGb) with C18 or HILIC columns under different ionization modes.

2.4. Statistical analysis

3.1. Nonlinear association between fluoride exposure and CR levels in the NHANES population

We applied RCS modeling to explore and visualize the dose–response relationship between UF and CR (Fig. 2A). The results indicated a significant inverted L-shaped nonlinear association (P for nonlinearity < 0.001). Specifically, CR levels rose rapidly with increasing UF concentrations below 0.51 mg/L, while the trend plateaued at concentrations above this threshold. Based on this inflection point, we conducted piecewise weighted linear regression. After adjusting for age, sex, BMI, race/ethnicity, and income-to-poverty ratio, each 1 mg/L increase in UF was associated with a 224.33 mg/dL increase in CR (95 % CI: 175.10, 273.57; P < 0.001) in the lower exposure range (UF < 0.51 mg/L). In contrast, for UF > 0.51 mg/L, the regression coefficient was 41.41 mg/dL (95 % CI: 11.12, 71.70; P = 0.02).

1. Download: Download high-res image (94KB)
2. Download: Download full-size image

Fig. 2. Associations between UF and CR (A) The RCS results of the NHANES database (B) The RCS results of the Chinese Elderly Database.

3.2. Validation of the UF–CR association in the Chinese older adult population

In the Chinese older adult population, the median age was 63 years (IQR: 55.8–69.0), with 35.5 % male and 64.5 % female participants. The median UF concentration was 1.25 mg/L (IQR: 0.88–1.81), and the median CR level was 21.41 mg/dL (IQR: 12.51–28.91).

After adjusting for age, sex, and BMI, RCS modeling revealed a significant nonlinear association between UF and CR (P for nonlinearity < 0.001). As UF increased, CR followed an inverted L-shaped trajectory, with an inflection point observed near UF = 1.02 mg/L, suggesting a potential threshold or physiological adaptation in CR excretion under increasing fluoride exposure (Fig. 2B). Among covariates, both age (P = 0.001) and sex (P < 0.001) were significantly associated with CR, while BMI was not significant (P = 0.53). Subsequently, piecewise linear regression analysis was conducted based on the identified inflection point to evaluate UF’s effect on CR across exposure strata. After adjustment for age, sex, and BMI, UF levels below 1.02 mg/L were associated with a significant increase in CR: for each 1 mg/L increase in UF, CR increased by 20.83 mg/dL (95 % CI: 13.93–27.73; P < 0.001). At UF levels >1.02 mg/L, the increasing trend persisted but with a reduced effect size (B = 3.27 mg/dL; 95 % CI: 1.79–4.75; P < 0.001). In addition, age was negatively associated with CR (B = –0.21 mg/dL; P < 0.001), and male participants had significantly higher CR levels compared with females (B = 6.38 mg/dL; P < 0.001).

3.3. Impact of fluoride exposure on the effectiveness of CR normalization in metabolomic analyses

A total of 995 urinary metabolites were detected, encompassing a wide range of metabolic pathways. According to Super Pathway classification, amino acids accounted for the largest proportion (28.74 %), followed by xenobiotics (24.82 %) and lipids (19.10 %). Other categories included incompletely annotated compounds (6.13 %), nucleotides (5.43 %), cofactors and vitamins (4.72 %), carbohydrates (4.22 %), energy metabolites (4.22 %), peptides (1.31 %), and a small proportion of hormones and secondary metabolites. This distribution suggests that the metabolomic platform used in this study offered broad coverage of amino acid, lipid, and xenobiotic pathways, supporting its capacity to capture fluoride-associated metabolic alterations.

To evaluate the influence of CR normalization on metabolite detection, we compared the data before and after CR adjustment. PCA revealed that PC1 explained 26.73 % of the variance before normalization, with a discernible group separation between cases and controls. After normalization, the explanatory power of PC1 declined to 19.42 %, and group separation was reduced accordingly (Figure S1). Heatmap visualization showed that, prior to adjustment, a broad range of metabolites—particularly amino acids, lipids, and xenobiotics—were upregulated in the case group, indicating pathway activation (Figure S2). After CR normalization, these differences were attenuated, and some metabolites even reversed in direction (Figure S3), suggesting that normalization may introduce systematic bias under certain exposure scenarios, thereby masking true biological signals.

Differential metabolite analysis corroborated these findings: 132 differential metabolites (124 upregulated, 8 downregulated) were identified before adjustment, compared with 148 (1 upregulated, 147 downregulated) after adjustment. Among them, 47.34 % showed a complete reversal in expression direction, and 24.72 % changed in statistical significance. Only 17 metabolites remained significantly different across both conditions, of which 9 exhibited opposite expression trends (Fig. 3). Notably, CR levels themselves showed an inverse trend between groups before and after normalization (Figure S4), further suggesting that conventional CR-based adjustment may obscure or distort true metabolic differences and introduce misinterpretation of exposure effects.

1. Download: Download high-res image (232KB)
2. Download: Download full-size image

Fig. 3. Volcanic map results of metabolomics differences test in Chinese elderly people (A) Raw data (B) Creatinine-adjusted (C) Specific gravity-adjusted (D) Osmolality-adjusted. N (down) and N (up) indicate the number of down-regulated and up-regulated metabolites, respectively.

In contrast, the results obtained after CR normalization showed markedly lower concordance with both the raw and SG/OSM-normalized profiles. Specifically, CR adjustment produced the largest number of metabolite direction reversals, the greatest shifts in statistical significance, and the lowest overlap in enriched pathways with the unadjusted data. Whereas SG and OSM normalization largely maintained the core metabolic signatures related to sulfur metabolism, endocrine regulation, and cardiometabolic signaling, CR normalization generated a pathway structure that deviated substantially from this pattern. This divergence reflected reduced preservation of exposure-related metabolic signals when creatinine was used as the adjustment factor under fluoride exposure conditions. The SG-normalized results are shown in Fig. 3C and Fig. 4C, and the OSM-normalized results are shown in Fig. 3D and Fig. 4D.

We also considered the presence of “exposed but unaffected” individuals in the analysis. Some control participants had urinary fluoride levels > 1.6 mg/L, while some cases had UF < 1.6 mg/L, reflecting the cumulative nature of skeletal fluorosis and short-term fluctuations in urinary fluoride. Their inclusion may attenuate metabolic differences between groups. To address this, we performed a sensitivity analysis excluding control samples with UF > 1.6 mg/L and case samples with UF < 1.6 mg/L, followed by BH correction. The results showed directionality consistent with the original analysis, and enrichment of fluoride-related pathways, supporting the biological plausibility of our findings. The urinary fluoride distribution between cases and controls and the sensitivity analysis results are provided in Figure S5.

Collectively, these results demonstrate that CR normalization exerts the greatest impact on metabolite directionality, statistical significance, and pathway interpretation, whereas SG and OSM normalization maintain higher fidelity to the underlying metabolic structure of the raw data under fluoride exposure.

4.1. Fluoride exposure substantially disrupts the stability of CR excretion

In both the NHANES adolescent and Chinese older adult populations, we consistently observed a significant inverted L-shaped nonlinear dose–response relationship between UF and CR levels, with CR rising rapidly at low UF concentrations and reaching a plateau at higher concentrations. In NHANES adolescents, the inflection point was approximately 0.51 mg/L UF, whereas in the Chinese older adults it was around 1.02 mg/L UF. This nonlinear pattern, replicated across two independent populations, provides robust evidence that fluoride exposure can directly disrupt CR excretion. The difference in inflection points between the populations may reflect variations in age structure, muscle mass, renal function, skeletal fluoride accumulation, and physiological adaptation. These observations suggest that the stability of CR as a urinary dilution correction factor may be compromised under conditions of fluoride exposure.

Traditionally, CR has been widely used in metabolomics, environmental epidemiology, and toxicology to adjust for urine dilution, based on the assumption that CR production is stable and not significantly influenced by exogenous exposures (Khamis et al., 2017, Li et al., 2025). However, both epidemiological and animal studies indicate that fluoride can impair proximal tubular function, disrupt ion channel activity, and alter bone and muscle metabolism, potentially leading to dynamic changes in CR excretion (Santoyo-Sanchez et al., 2013, Çetin and Yur, 2016, Luo et al., 2017). These mechanisms may lead to dynamic changes in CR excretion, undermining the assumption of its stable elimination. Previous studies have indicated that CR excretion is susceptible to various endogenous and exogenous factors, such as urine flow rate, hydration status, age, sex, body size, and diet. This variability is particularly pronounced in individuals with unstable renal function or under environmental exposure conditions (Goldstein, 2010, Gaines et al., 2010, O’Brien et al., 2016). Therefore, using CR for normalization under fluoride exposure may not only fail to correct for dilution bias but also introduce new sources of systematic error.

4.2. Risks of signal attenuation and underestimation of causal effects induced by CR normalization

Our results further demonstrate that CR normalization substantially interferes with the extraction of metabolomic signals. In the unnormalized data, skeletal fluorosis cases and controls were clearly separated in the PCA model, with differential metabolites predominantly enriched in amino acid, lipid, and xenobiotic pathways. After CR normalization, PC1 variance decreased from 26.73 % to 19.42 %, group separation was markedly reduced, and some differential metabolites reversed direction. In differential metabolite analysis, 47.34 % of metabolites showed complete direction reversal, 24.72 % changed in statistical significance, and only 17 metabolites remained significant under both conditions, among which 9 exhibited opposite trends. These results indicate that CR normalization can obscure true biological signals and potentially mislead exposure-phenotype interpretations. These observations are consistent with previous studies showing that when a normalization factor covaries with the exposure, CR adjustment can attenuate or distort exposure–phenotype associations (Chen et al., 2013, Bulka et al., 2017). Our study provides empirical evidence of this phenomenon using concrete metabolomic data, highlighting the need for caution when using CR as the sole normalization factor under fluoride exposure conditions.

4.3. Impact of CR normalization on mechanistic interpretation through pathway identification

Beyond statistical features, we further observed substantial shifts in the direction of enriched metabolic pathways before and after CR normalization. In the unnormalized data, enriched pathways were primarily associated with energy and endocrine-related mechanisms, such as cysteine and methionine metabolism, the sulfur cycle, circadian rhythm, and longevity regulation—pathways potentially influenced by fluoride exposure. After CR correction, however, the significantly enriched pathways shifted toward dopaminergic synapse activity, protein digestion, and pyrimidine metabolism, indicating a “drift” in mechanistic interpretation. This finding suggests that CR normalization not only affects the identification of differential metabolites but may also distort the interpretation of exposure-related mechanisms at the pathway level. This issue of miscalibrated mechanistic inference has also been highlighted in a systematic NMR-based metabolomics study by Li et al., who recommended against using CR as the sole normalization factor when it covaries with the exposure variable (Li et al., 2022b).

4.4. Caution in the use of CR normalization in epidemiological studies

Our findings indicate that CR normalization should be applied with caution in environmental exposure studies, particularly when the exposure factor may influence renal function or CR metabolism. Under fluoride exposure, CR excretion exhibits nonlinear behavior with population-specific inflection points (approximately 0.51 mg/L UF in NHANES adolescents and 1.02 mg/L UF in Chinese older adults), suggesting that CR may no longer satisfy the assumption of stable, exposure-independent production.

Metabolomic analyses show that CR normalization substantially alters the identification of differential metabolites, reverses their expression trends, and shifts pathway enrichment, potentially leading to underestimation or misdirection of exposure–phenotype relationships. Similar issues have been documented in previous studies. For example, a 2022 systematic comparison using Finnish cohort data evaluated seven urine metabolomics normalization strategies and found that, while CR normalization could partially adjust for urine dilution, it significantly altered the estimated associations between clinical phenotypes (e.g., BMI, mean arterial pressure) and metabolites in regression models, suggesting potential underestimation or directional shifts in exposure–phenotype relationships (Li et al., 2022b).

Therefore, unless the independence and stability of CR can be clearly validated in a given research context—particularly under conditions of exogenous exposure or fluctuating renal function—CR should not be used as the sole normalization method. Researchers should carefully assess the mechanistic validity of CR normalization in the context of exposure pathways, individual physiological status, and potential covariation with the exposure variable, to avoid systematic bias and misinterpretation of exposure effects in both statistical and mechanistic analyses.

4.5. Comparison of normalization strategies and recommendations

Given the limitations of CR normalization, we evaluated the performance of multiple alternative normalization strategies. Our results indicate that, compared with CR adjustment, specific gravity (SG) and osmolality (OSM) normalization better preserved the directionality of metabolite changes, the significance of differential metabolites, and core metabolic pathways, including sulfur metabolism, endocrine regulation, and cardiometabolic signaling. The metabolite reversals, shifts in statistical significance, and pathway alterations observed with CR normalization were not seen with SG or OSM, suggesting that these methods are more robust when the normalization factor covaries with the exposure variable.

Furthermore, combining CR correction with MS total useful signal (MSTUS) normalization has been shown to improve the stability and consistency of metabolomic data (Chen et al., 2013). Other potential strategies include total ion current (TIC) normalization and reference metabolite-based approaches (e.g., IS-PSEURID), which have demonstrated higher robustness in populations with environmental exposures or fluctuating renal function (Chen et al., 2013, Li et al., 2022b).

Overall, environmental epidemiology studies should tailor normalization strategies to the specific exposure context, individual physiological status, and characteristics of the normalization factor. By accounting for dilution effects and exposure-related variability, integrating multiple normalization methods can enhance both the reliability and biological interpretability of metabolomic findings. Our results suggest that, under fluoride exposure, relying solely on CR normalization may introduce systematic bias, whereas SG, OSM, or combined normalization strategies provide a more robust approach for metabolomic data correction.

4.6. Strengths and limitations

This study has several limitations. First, both population datasets were cross-sectional, restricting causal inference. Nevertheless, the reproducibility of the inverted L-shaped UF–CR association across two independent populations reinforces the robustness of the observed nonlinear pattern. Second, the metabolomics analysis was conducted in a relatively small sample (25 cases and 25 controls). However, rigorous quality control procedures, convergent results across sensitivity analyses, and coherent pathway-level signatures reduce concerns about statistical instability. Finally, although multiple normalization strategies were compared, non-targeted metabolomics inherently carries analytical variability; future studies incorporating targeted platforms and longitudinal exposure assessments will be essential to refine mechanistic interpretations.