Abstract
Objectives:
To identify the core hypoxic injury pattern of DKD, construct a DKD risk model based on hypoxic injury-related (HIR) score, and explore the potential therapeutic targets of DKD.
Methods:
DKD-related microarray-based transcriptomic analyses, single-nucleus RNA sequencing (snRNA-seq) and spatial transcriptomics were retrieved from the Gene Expression Omnibus (GEO) database. Seven HIR gene sets were obtained from various public databases. Core hypoxic genes were identified using different machine-learning algorithm. LASSO and nomogram were applied to construct a HIR risk score for cellular hypoxic damage. The detailed expression of hub gene would be showed in the single cell and kidney region. The prognostic value of the HIR score was externally validated using plasma proteomics from the United Kingdom Biobank.
Results:
Five core hypoxic injury pathways in DKD were identified: Hypoxia, Autophagy, Ferroptosis, Endoplasmic Reticulum (ER) Stress, and Apoptosis. The HIR risk score was constructed based on three hub genes: CASP3, DUSP1, and ZFP36. The HIR score demonstrated high diagnostic efficiency for DKD patients. Higher HIR scores were associated with significantly infiltrated immune cells and poorer kidney function. In United Kingdom Biobank validation, the HIR score significantly improved the prediction of kidney outcomes, renal death, and secondary endpoints beyond demographic and metabolic variables (AUC increments 0.04–0.06), and correlated negatively with eGFR and positively with lipoprotein(a). The calculated tissue-level HIR scores also showed a highly significant and robust increase in the renal microenvironment of BTBR ob/ob mice.
Conclusion:
These results provided a predictive model for clinical evaluation in patients with DKD and also a new insight into the role of HIR genes in the pathogenesis of DKD.
Graphical Abstract
1 Intoduction
Diabetic Kidney Disease (DKD) is a renal disease secondary to diabetes, primarily characterized by increased urinary albumin excretion and decreased glomerular filtration rate. It is a common chronic complication of diabetes. Reports indicate that 27%–54% of diabetic patients develop DKD (; ), leading to substantial healthcare costs. Both the Kidney Disease Outcomes Quality Initiative (KDOQI) of the National Kidney Foundation and the Chinese guidelines for the prevention and treatment of diabetic kidney disease stated that eGFR and/or the urine albumin-to-creatinine ratio (UACR) were the primary clinical diagnostic criteria for DKD (; ), while kidney biopsy was considered the gold standard for DKD diagnosis. However, eGFR and albuminuria have limitations in sensitivity and specificity, and kidney biopsy’s invasive nature limits its widespread use. The identification of new biomarkers could assist in diagnosing DKD and predicting disease progression.
Chronic hypoxia is an early indicator of renal pathological changes in DKD (Nangaku, 2006). Hyperglycemia induces renal cell hypoxia through various pathological mechanisms. Firstly, hyperglycemia and hyperperfusion lead to glomerular hyperfiltration and increased tubular reabsorption workload, raising oxygen consumption and causing functional hypoxia in the tubules (Singh et al., 2008). Secondly, persistent hyperglycemia triggers oxidative stress, damaging mitochondrial function, disrupting the cellular respiratory chain, and impairing oxygen synthesis (Perco and Mayer, 2018). Thirdly, dysregulation of glucose and lipid metabolism causes endothelial dysfunction, impaired vasodilation, and increased vasoconstriction, leading to increased vascular resistance, impaired blood flow regulation, and reduced renal blood flow and perfusion (). Hypoxic injury plays a crucial role in DKD development, causing decreased intracellular ATP levels, excessive reactive oxygen species production, and activation of apoptotic pathways. This results in the upregulation of fibrotic factors and interstitial fibrosis, further inducing typical glomerular pathological changes such as thickening of the glomerular basement membrane, mesangial cell proliferation, and endothelial cell damage ().
Hypoxia-induced biological processes play important roles in DKD progression, including ferroptosis, mitochondrial reactive oxygen species release, mitochondrial autophagy, endoplasmic reticulum-associated stress, apoptosis, and cell cycle arrest (; Zheng et al., 2023). Currently, our understanding of the relationship between hypoxic injury and DKD remains incomplete, including the relative significance of hypoxia-injury pathways involved in DKD progression and their potential synergistic interactions. Given that hypoxia-injury involved multiple cellular pathways and biological factors, bioinformatics methods were essential to identify key molecules and construct disease models. This study aimed to develop a hypoxic injury-related (HIR) risk score—using core HIR molecules involved in DKD progression to predict renal function and immune infiltration in DKD patients. This would reveal the heterogeneity of DKD patients based on hypoxia injury risk characteristics and explore the specific cell types and pathways mediating this heterogeneity.
2 Materials and methods
2.1 Data acquisition
All datasets were downloaded from the Gene Expression Omnibus (GEO) database (http://www.ncbi.nlm.nih.gov/geo/). We used microarray-based mRNA expression profiles (GSE96804, GSE30528, GSE30529, GSE104954) (Pan et al., 2018; Woronieck et al., 2011; ) for initial analysis, including 63 control samples and 67 DKD samples, all from microdissected human kidney samples covering both glomerular and tubular tissues. The single-nucleus RNA sequencing (snRNA-seq) dataset included 6 control samples and 6 DKD samples, merged from GSE131882 and GSE195460 (Wilson et al., 2019; ; Wilson et al., 2022), resulting in 54,112 snRNA-seq. The spatial transcriptomics data originated from GSE261545. Details of the collected datasets can be found in Supplementary Table S1. Key regulatory genes for the 7 hypoxic injury modes were obtained from various sources, including the KEGG database (KEGG.db v4.0), GSEV-Molecular Signatures Database, Genecard database, DisGeNet database, Reactome database, and Comparative Toxicogenomics Database. This included 330 genes related to hypoxia response (hypoxia), 330 genes related to apoptosis (apoptosis), 265 genes related to mitochondrial autophagy (autophagy), 65 genes related to ferroptosis (ferroptosis), 263 genes related to mitochondrial reactive oxygen species release (ROS), 170 genes related to endoplasmic reticulum stress (ER-stress), and 157 genes related to cell cycle arrest (cell cycle). Clinical data for DKD patients were obtained from the Nephroseq v5 online database.
2.2 Data processing
We combined the microarray data of mRNA expression profiles using the R packages “limma” and “sva”, and removed batch effects using the ComBat function (Swindell et al., 2019). The distribution patterns of the samples were displayed through box plots, principal component analysis (PCA) plots, and 3D projection plots, as detailed in Supplementary Material 1. Genes with log2FC > 0.5 and P value <0.05 were considered differentially expressed genes (DEGs) (Wang et al., 2023), and volcano plots and heatmaps were generated.
The snRNA-seq and spatial transcriptome data were processed using the “Seurat” package (Xu et al., 2021). High-variable genes expressed in at least 3 cells and marker genes expressed in at least 50 features were extracted. Cells with gene expression <200 or >5000, mitochondrial gene expression >20%, ribosomal gene expression >20%, or hemoglobin gene expression >5% were excluded. The SCTransform function and PCA were used to normalize and scale the raw counts. Unsupervised clustering analysis and uniform manifold approximation and projection (UMAP) were employed to identify discrete cell clusters in the snRNA-seq dataset. Cell annotation information for snRNA was obtained from the K.I.T database (http://humphreyslab.com/SingleCell/) and Cellmarker 2.0 database (), with all clusters and subclusters manually annotated. The annotation for spatial pathological section was referred to the original article ().
2.3 Functional enrichment analysis and pathway activity calculation
We used the R package “clusterProfiler” for KEGG and Gene Ontology (GO) enrichment analysis (). To identify potential biological pathways, a significance threshold of 0.05 was set, and gene set enrichment analysis (GSEA) was conducted using the “GSEABase” package (Yin et al., 2022). The GSVA package was used to explore potential pathway-level changes in gene expression to determine pathway activity differences between samples.
2.4 Machine learning and construction of receiver operating characteristic curves
Based on pathway enrichment scores obtained from GSVA, we optimized the selection of gene sets for HIR modules using the weighted gene co-expression network analysis (WGCNA) algorithm (). The top 25% of genes by variance were used as clustering templates. We evaluated the predictive ability of seven hypoxic injury modes using the “pROC” package to calculate AUC values (Robin et al., 2011). Three machine learning algorithms were used to screen core hypoxic injury genes: Least Absolute Shrinkage and Selection Operator (LASSO) (), the Support Vector Machine-Recursive Feature Elimination (SVM-RFE) (), and the Regularized Random Forest (RRF) (Pfeifer et al., 2022). The intersection of the results was used for predictive model construction, with coefficients generated by LASSO. The discriminatory ability of the HIR score for distinguishing DKD patients from controls was evaluated using ROC curve analysis with the pROC R package, and model robustness was rigorously assessed through multiple complementary validation strategies: bootstrap resampling (1,000 iterations) for 95% confidence interval estimation, nested cross-validation (5-fold outer, 3-fold inner, repeated 50 times), 10-fold cross-validation (repeated 100 times), leave-one-out cross-validation (LOOCV), Lasso and Ridge regularized logistic regression using the glmnet R package, and permutation testing (1,000 random permutations of group labels), with the overfitting gap defined as the difference between the apparent AUC and the nested cross-validation AUC. A nomogram was drawn to evaluate the model’s predictive ability (Wu et al., 2020).
2.5 Immune infiltration analysis
To decipher the cell-type-specific immune landscape within the heterogeneous renal parenchyma, single-sample Gene Set Enrichment Analysis (ssGSEA) was implemented using the GSVA R package framework. The baseline immunogenomic marker sets deployed for calculating the relative infiltration scores were retrieved from the authoritative multi-tissue immune atlas validated by Charoentong et al. (; ). Crucially, to justify the organ-specific suitability of these signatures for the rigid, fibrotic microenvironment of diabetic kidney disease (DKD), the marker matrices were cross-referenced and benchmarked against established human kidney single-cell and single-nucleus transcriptomic datasets. This tissue-specific alignment successfully accounted for resident myeloid lineages and interstitial stroma-associated immunomodulatory cells while filtering out blood-borne cellular biases. Hierarchical clustering of the infiltration profiles was subsequently executed using the pheatmap architecture via the Ward. D2 method paired with squared Euclidean distance metrics. Statistical variations in differentially infiltrated immune cell lineages between groups were analyzed and visualized using standard exploratory data analysis pipelines.
2.6 Single-cell downstream analysis
The “FeaturePlot” function was used to display gene expression locations and quantitative differences. The cell pseudotime analysis was conducted using the BiocGenerics and monocle packages to explore differentiation status and dynamic gene expression. Cell communication analysis was performed using the “cellchat” package to reveal pathway connections between primary cells and other cells ().
2.7 Spatial transcriptome analysis
The conditional autoregressive-based deconvolution (CARD) was applied for spatial division and inter-regional difference analysis (). The hub genes would be mapped to specific locations in pathological sections to reflect its actual expression. Regional differential genes would be used for KEGG enrichment analysis to verify the importance of hypoxia response.
2.8 HIR score for the prognosis of DKD based on United Kingdom biobank proteomics
We externally validated three hypoxia-related genes (CASP3, DUSP1, ZFP36) using the United Kingdom Biobank, a prospective population-based cohort of approximately 500,000 adults aged 39–70 years recruited between 2006 and 2010, with linked national mortality, hospitalisation and cancer registries. High-throughput plasma proteomics were measured using the Olink Explore 1536 panel in a probability subsample. For our specific implementation, proteins demonstrating missingness rates >20% were excluded from downstream analysis. The circulating levels of the available target proteins matching our core risk components (e.g., CASP3) were extracted to compute the integrated clinical HIR score. Longitudinal hard endpoints—including functional kidney decline, renal death, and secondary macrovascular outcomes—were cross-referenced using ICD-10 registry codes linked to electronic health records to benchmark the fully adjusted incremental predictive metrics. Diabetic kidney disease (DKD) cases were identified through an integrative approach combining ICD-10 codes (E10.2, E11.2, etc.), established clinical guidelines and laboratory criteria (eGFR <60 mL/min/1.73 m2 and/or urine albumin-to-creatinine ratio ≥30 mg/g sustained for ≥3 months). After excluding participants with pre-existing kidney failure, missing proteomic data or failing quality control, 918 DKD patients were included. The composite kidney outcome comprised renal death (primary outcome), end-stage kidney disease, sustained eGFR decline ≥40% from baseline, and sustained UACR increase ≥30%. Receiver operating characteristic (ROC) curves assessed predictive performance of individual protein markers alone and in combination with demographic factors (age, sex, etc.) and metabolic traits (HbA1c, glucose, lipids, blood pressure, renal markers). DeLong tests compared AUC differences across models, and 2000-iteration bootstrap analyses estimated AUC confidence intervals. Correlative mapping between the protein abundance encoded by the HIR gene and routine clinical indices employed Spearman rank correlation, with coefficients (β) reported and significance defined as P < 0.05.
2.9 The validation of HIR score in animal model
To experimentally validate the calculated hypoxia injury risk (HIR) signature, 8-week-old male black and tan brachyury (BTBR) ob/ob diabetic mice (n = 5) and age-matched male BTBR wild-type (WT) mice (n = 5) were maintained under standard 12-h light/dark cycles with ad libitum food and water until ethical euthanasia at 24 weeks of age. Harvested renal tissues were fixed in 4% paraformaldehyde, embedded in paraffin, sectioned at 4 μm and subjected to Periodic Acid-Schiff (PAS) staining to track structural lesions. Concurrently, total RNA was isolated from cryopreserved, mechanically pulverized kidney samples using TRIzol reagent, purified, and reverse-transcribed into cDNA. Quantitative real-time PCR (qRT-PCR) amplification was subsequently executed using SYBR Green Master Mix to evaluate the transcription profiles of the core target genes (ZFP36, DUSP1, and CASP3), with target fold-changes calculated using the comparative 2−ΔΔCT method relative to the internal house-keeping control GAPDH.
3 Results
3.1 Data merging
The four DKD datasets were merged and normalized after batch effects were removed using the ComBat algorithm. The distribution patterns between normal and DKD samples were visualized using principal component analysis (PCA), box plots, and 3D projection plots (Supplementary Figure S1). These results confirmed that batch effects were successfully removed.
3.2 Identification and enrichment analysis of DEGs between control and DKD groups
The significant differences were observed between normal and DKD samples, identifying 798 DEGs. The volcano and heatmaps (Figures 1A,B) indicated that 387 genes were upregulated and 411 genes were downregulated in the DKD group. The pathway enriched were most associated with metabolism such as “Tryptophan metabolism,” “Fatty acid degradation,” and “Tyrosine metabolism” (Figure 1C). GSEA has revealed the potential pathways such as “ECM-receptor interaction” and “Cell adhesion molecules” (Figure 1D). The top activated pathway was “ECM-receptor interaction,” while the top three suppressed pathways were “Carbon metabolism,” “Oxidative phosphorylation,” and “Peroxisome” (Figure 1E).
FIGURE 1
3.3 Identification and enrichment analysis of HIR DEGs in control and DKD groups
We then focused the differences in HIR genes between normal controls and DKD patients. The heatmap and volcano plot (Figures 2A,B) showed the key HIR genes which expressed deferentially in DKD patients. GO enrichment indicated that HIR DEGs were enriched in biological processes associated with hypoxia regulation (Figure 2C). KEGG analysis revealed that DEGs participated in downstream regulation of hypoxia like “Apoptosis” and the “HIF-1 signaling pathway” (Figure 2D). GSVA showed that 2 pathways, apoptosis and cell cycle, were upregulated in DKD patients, while 4 HIR pathways showed decreased activity in DKD (hypoxia, autophagy, ROS, ER stress) (Figure 2E). Furthermore, a heatmap of pathway activation scores clearly demonstrated distinct clusters between the two groups (Figure 2F). In agreement with the GSVA results, pathways including hypoxia, autophagy, ROS, and ER stress exhibited consistently reduced activity patterns in the DKD group compared to normal controls, whereas apoptosis and cell cycle pathways were visibly upregulated in DKD samples.
FIGURE 2
After data pre-processing (Supplementary Figure S2), 12 cell types were manually identified: mesenchymal cells (MES), glomerular parietal epithelial cells (PEC), proximal convoluted tubular cells (PCT), Loop of Henle cells (LOH), distal convoluted tubular cells (DCT), convoluted tubular cells (CT), collecting duct-principal cells (CD-PC), collecting duct-intercalated cell type A (CD-ICA), collecting duct-intercalated cell type B (CD-ICB), podocytes (PODO), endothelial cells (ENDO), and monocytes (MONO). The AUCell analysis showed the activity of HIR geneset decreased in the following order: CD-ICA, CD-ICB, PEC, CD-PC, ENDO, PODO, LOH, DCT, CT, MES, MONO, PCT (Supplementary Figure S3).
The proportions of several immune cells were found to be increased in DKD patients in Supplementary Figure S4. The Hypoxia, Autophagy, Ferroptosis, and ROS geneset, had a significant negative correlation with the infiltration of many immune cells like effector memory CD8 T cells, activated CD4 T cells and activated dendritic cells, while the Apoptosis geneset showed a positive correlation. The biomarkers APOA4, IGFBP3, CD5L, and FGF23 were independent risk factors for composite endpoints of DKD microalbuminuria, doubling of serum creatinine, ESRD, and cardiovascular events (Titan et al., 2011; ). Supplementary Figure S5 showed that Apoptosis was significantly positively correlated with three biomarkers (IGFBP3, CD5L, FGF23).
3.4 Analysis of the module closely related to HIR geneset in DKD
The WGCNA analysis revealed that the yellow module exhibited the strongest correlation with Autophagy and cell cycle, while the brown module showed the strongest correlation with Apoptosis, the red module with ER stress, the turquoise module with ROS and Ferroptosis, and the grey module with Hypoxia (Supplementary Figure S6). The corresponding intersections were depicted in Supplementary Figure S7. Summarizing the intersections across the seven group genes, KEGG analysis showed significant expression of “AGE-RAGE signaling pathway in diabetic complications”, “Apoptosis”, and “reactive oxygen species”, consistent with the enrichment results of HIR geneset.
In total of 32 genes were subjected to LASSO regression, resulting in the identification of 11 genes, namely, DUSP1, NFIL3, PIK3R3, S100A4,ATF3, IGF1, CASP3, CLU, IGFBP6, ZFP36, HIF1A as hub genes in the DKD dataset (Supplementary Figures S8A–C). By SVM-RFE analysis, 13 genes have been calculated as the best model (Supplementary Figures S8D,E) and ranked by importance. By RRF analysis, 5 genes has been filtered out (Supplementary Figures S8G,H) and ranked by accuracy.
Using the intersection of three machine learning algorithms, we identified three core HIR genes: ZFP36, CASP3, and DUSP1 (Figure 3A). LASSO regression analysis determined the factors for constructing the HIR risk score. The formula was established as follows: ZFP36 expression * (−0.90653) + CASP3 expression * (2.75490) + DUSP1 expression * (−3.16704) (Figures 3B,C). Compared to normal controls, ZFP36 and DUSP1 were significantly downregulated, while CASP3 was upregulated in the kidney tissues of DKD patients (Figures 3D–F). The HIR score was significantly increased in the DKD group compared to the control group (Figure 3G). The higher the HIR score was, the higher the prediction efficiency for DKD would be. In the test set, the average HIR score of the control group was −9.80, with a median score of −9.10 and a standard deviation of 2.94. The average HIR score of the DKD group was −1.84, with a median score of −1.49 and a standard deviation of 2.16. The ROC curve demonstrated that the HIR score exhibited good discriminatory performance for diagnosing DKD, with an apparent AUC of 0.9893 (95% bootstrap CI: 0.9730–0.9995). At the optimal threshold of −5.9704 (determined by Youden index), the model achieved a sensitivity of 95.52% and specificity of 98.41% (Figure 3H). Comprehensive cross-validation confirmed model robustness, with nested cross-validation AUC of 0.9897, 10-fold cross-validation AUC of 0.9892, and leave-one-out cross-validation AUC of 0.9893, yielding a negligible overfitting gap of −0.0004 (Supplementary Figure S9). The Nomogram-based model also demonstrated high weight and predictive capability for the three hub genes (Figures 3I,J).
FIGURE 3
3.5 Identification of DEGs grouped by HIR risk score
According to the median HIR risk score, DKD patients were divided into two groups. We identified 202 DEGs between the low-risk and high-risk groups, including 71 upregulated and 131 downregulated genes (Supplementary Figures S10A,B). GO analysis revealed that the DEGs were enriched in typical process of kidney disease formation like “renal system development” “leukocyte migration” and “glycosaminoglycan binding” (Supplementary Figure S10C). KEGG analysis indicated some inflammatory process like “AGE-RAGE” and “PI3K-Akt signaling pathway” (Supplementary Figure S10D). The GSEA results showed that the MAPK and Ras signaling pathways were activated (Supplementary Figures S10E,F).
The high-risk group exhibited higher levels of immune cell infiltration compared to the low-risk group (Supplementary Figure S11A). Additionally, we found significant positive correlations between the HIR risk score and three core hypoxia genes with certain immune cell types, particularly Type 17 helper T cells and plasmacytoid dendritic cells (Supplementary Figure S11B).
3.6 Validation of HIR risk score in public database
The GSE99339, GSE162830 and GSE142025 were selected as the validation cohort. The expression of three hub core genes (CASP3, ZFP36, DUSP1) showed significant differences in the DKD group (Figures 4A–F) in two validation dataset. The HIR risk score of the DKD group was significantly higher than that of the control group (Figures 4G–I). We generated an ROC curve to evaluate the diagnostic performance of the HIR risk score for DKD, with an AUC of 0.961 and 0.883 (Figures 4J,K). The clinical data showed a significant negative correlation between the risk score and GFR levels (Figure 4L).
FIGURE 4
3.7 Exploration of hub gene expression at single-cell level and spatial transcriptome
The expression characteristics of 3 hub genes were shown in a single-cell level. CASP3 was found to be dispersed among kidney cells (Figure 5A) and differentially upregulated in CD-ICB, CT, and PEC cells in the DKD group compared to the control group (Figure 5B). It expressed more in high HIR score cell types (Figure 5C) and was most enriched in CD-ICA cells (Figure 5D). DUSP1 was concentrated in CD-PC, CT, and ENDO cells (Figure 5E), but differentially downregulated in ENDO and LOH cells (Figure 5F). It also expressed low in high HIR score cell types (Figure 5G) and was least enriched in LOH cells (Figure 5H). ZFP36 was dispersed among kidney cells (Figure 5I) and differentially downregulated in CD-ICB, CD-PC, CT, ENDO, and LOH cells (Figure 5J). It showed low in HIR score cell types (Figure 5K) and was least enriched in LOH cells (Figure 5L).
FIGURE 5
Differential analysis showed that both DUSP1 and ZFP36 were significantly downregulated in LOH cells (Figure 6A). So the LOH cells were segregated into high-risk and low-risk groups based on the median HIR risk score (Figure 6B). The Comparison based on HIR risk scores revealed that the high-risk subgroup of LOH cells ranked by the front (Figure 6C). The pseudotime analysis of LOH unveiled a developmental trajectory progressing from shallow to deep through four pivotal branching points (Figure 6D). The expression of DUSP1 and ZFP36 peaked in mature LOH cells and was lowest in intermediate PCT cells, indicating an enrichment pattern associated with cell differentiation (Figure 6E).
FIGURE 6
The cell communication revealed that the high-HIR risk subgroup of LOH cells communicated with other cell types through a complex network (Supplementary Figure S12A). Among them, EGF served as a key intermediate pathway for both incoming signaling patterns (Supplementary Figure S12B) and outgoing signaling patterns (Supplementary Figure S12C).
For spatial transcriptome, 18,085 genes,2802 sampling points of a kidney tissue section were included in the clustering. The data quality control was shown in Supplementary Figures S13A,B. The section could be divided into eight regions (Figures 7A,B), which were as follows: glomeruli, injured tubules (Inj−T), arteries within the capsule (Artery−C), loop of Henle and collecting ducts (LH−CD), Proximal tubules (PT), tubules with endoluminal protein casts (Cast-T), arteries within the renal parenchyma (Artery-K) and tumor (Figures 7C–F). DUSP1 and ZFP36 have low expression in LH−CD and are widely distributed in all regions, while CASP3 has a low distribution, which were consistent with results of snRNA-seq (Figures 7G–I). Top 10 representative markers for each region was shown in Supplementary Figure S13C. KEGG enrichment of those genes revealed the high expression of HIF-1 signalling pathway in injured tubules and proximal tubules (Supplementary Figure S13D).
FIGURE 7
3.8 HIR score for the prognosis of DKD based on United Kingdom biobank proteomics
Among 918 DKD patients, 63 (6.9%) developed renal composite endpoint, as shown in Supplementary Table 2. The flowchart of patients’ screening was shown in Supplementary Figure S14. The composite kidney outcome group (N = 139) was significantly older, had longer diabetes duration, lower eGFR, higher ACR and urea, and greater cystatin C than controls (N = 779). When stratified, the primary outcome subgroup (renal death, N = 127) mirrored these differences, with larger waist circumference and higher CRP. The secondary outcome subgroup (N = 16) exhibited the most severe renal impairment with lowest eGFR (47.7 mL/min/1.73 m2), highest urea (13.2 mmol/L), and highest cystatin C (1.78 mg/L).
The multi-dimensional validation of HIR score in the proteinomics of United Kingdom Biobank for the prognosis of DKD was shown in Figure 8. In DKD patients, with the kidney composite outcome (kidney outcome), primary outcome (renal death, primary outcome), and secondary outcome (ESRD/KRT/eGFR decline ≥40%, secondary outcome) as endpoints, the predictive efficiency of 6 different models was tested. In the baseline model that only contained demographic and metabolic indicators, after adding the HIR score, the AUC of the three outcomes increased by approximately 0.043 (0.663–0.706), 0.053 (0.676–0.729), and 0.059 (0.789–0.848), respectively, and the 95% confidence intervals all shifted significantly upward. When the model already included eGFR, Urea, ACR, and other renal function indicators, the incremental contribution of adding the HIR score decreased (kidney: +0.008; primary: +0.006; secondary: +0.015) (Figure 8A). The AUC distribution obtained from 10 repeated 5-fold cross-validations (n = 50) showed that in the kidney and primary outcomes, the comparison between “demographic + metabolic” and “demographic + metabolic + HIR” was statistically significant (p = 0.005 and p = 0.001), and the secondary outcome also reached p = 0.044 (Figure 8B). The correlation heatmap showed that the HIR score has the strongest negative correlation with eGFR (the blue bars point to the left), and has a positive correlation with Lipoprotein(a) and Urea (Figure 8C).
FIGURE 8
3.9 In Vivo experimental validation of the HIR_score and core signatures
To assess the biological validity and translational relevance of the computed computational signatures, we established an in vivo diabetic kidney disease (DKD) model using BTBR ob/ob mice alongside age-matched BTBR WT controls (Figure 9A). Phenotypically, 24-week-old BTBR ob/ob mice displayed profound diabetic nephropathy traits, characterized by significantly elevated body weight, blood glucose levels, kidney weight, and urinary albumin-to-creatinine ratios (ACR) compared to the WT cohort (Figure 9B). Histopathologically, PAS and Masson’s trichrome stainings demonstrated pronounced structural abnormalities, featuring prominent glomerular basement membrane thickening, expansion of the mesangial matrix, and extracellular matrix deposition in the diabetic group (Figure 9C). At the molecular transcript level, qRT-PCR analysis from pulverized renal tissues demonstrated that the expression profiles of the core prognostic candidates matched our algorithmic predictions. Specifically, the expression of the protective components ZFP36 and DUSP1 was significantly downregulated, whereas the pro-apoptotic executive vector CASP3 was markedly elevated in the DKD group. The calculated tissue-level HIR_score showed a highly significant and robust increase in the renal microenvironment of BTBR ob/ob mice (Figure 9D).
FIGURE 9
4 Discussion
Our study assessed the activity of 7 different types of hypoxic injury in DKD and further explored their activity in various kidney cells using snRNA-seq data. We also investigated the hypoxic expression patterns of hub genes in core cells and kidney regions. Additionally, we validated the relationships between these 7 types of hypoxic injury and immune infiltration, genes associated with microalbuminuria in DKD, and GFR. We constructed a HIR score based on 3 core hypoxic injury genes (CASP3, DUSP1, and ZFP36) using observational microarray transcriptome analysis. Both internal and external datasets showed that the HIR risk score was a powerful diagnostic indicator and a marker of disease severity in DKD patients. These results provided new insights into the role of HIR genes in the pathogenesis of DKD.
We found significant changes in the expression of Hypoxia, Autophagy, Ferroptosis, endoplasmic reticulum (ER) stress, and Apoptosis in the DKD group. The response to hypoxia is the initial and central aspect of cellular hypoxic injury. Hypoxia transcriptionally induces a robust set of genes controlled by HIFs but also a range of other transcription factors, including nuclear factor-κB (NF-κB) (Yang et al., 2025). The vast majority of O2 sensitive genes are in fact direct HIF targets. These genes at the cellular and organism level help adapt to diminishing levels of O2. However, persistent activation of hypoxia-induced genes can result in pathologies, including metabolic and kidney disease (). Notably, our bulk tissue-level GSVA indicated a decreased overall activity of the hypoxia pathway in DKD, which seemingly conflicts with the paradigms of hypoxia-driven progression. This phenomenon might be rationalized by several factors. First, chronic adaptation and long-term evolutionary pressure in advanced DKD can prompt a feedback downregulation or exhaustion of classical hypoxia responsive machinery. Second, the prominent fibrotic remodeling and cellular dropout (e.g., specialized tubular segments) in bulk samples could cause a signal dilution effect, obscuring localized hypoxic stress, a notion supported by our single-cell and spatial evidence showing distinct cell-type specific hypoxic vulnerabilities.
In our study, we identified ZFP36, DUSP1, and CASP3 as new indicators of hypoxic injury. CASP3, or caspase-3, is positively correlated with DKD progression. It is involved in both hypoxia tolerance regulation and downstream apoptotic responses, which acts as the primary executioner caspase in apoptosis and is a key component of the cytotoxic T lymphocyte (CTL) killing mechanism (). Under hypoxia, caspase-3 served as the main executioner caspase and lead to the morphological and biochemical changes (Jiang et al., 2025). Extensive research has demonstrated that CASP3 is a risk factor for the onset and progression of DKD. In animal experiments, Sekiko et al. found that caspase-3 was significantly upregulated in the renal tubular cells of diabetic mice (Taneda et al., 2010). Wen et al. observed that injecting caspase-3 inhibitors into streptozotocin-induced diabetic mice significantly alleviated albuminuria, improved renal function, and mitigated pathological changes (). Clinical observations showed that in patients with type 1 diabetes and early diabetic nephropathy, increased caspase-3 levels were associated with a higher risk of DKD (Wen et al., 2020).
This study showed that ZFP36 and DUSP1 were significantly negatively correlated with DKD progression. The zinc finger protein 36 homolog (ZFP36), also known as tristetraprolin (TTP), is a protein that plays a regulatory role in hypoxic injury by stabilizing HIF-1α mRNA. The ZFP36 protein directly binds to the 3′-UTR of HIF-1α mRNA, thereby reducing the stability of HIF-1α mRNA. It plays a physiological role in regulating cellular adaptation and apoptosis under prolonged hypoxia (; ). A decrease in ZFP36 has been shown to be a risk factor for the progression of DKD (; Tang et al., 2024; ). DUSP1, or Dual Specificity Phosphatase 1, is an important downstream regulator of the hypoxia tolerance response. In hypoxic conditions, DUSP1 is upregulated through an HIF-1α-dependent mechanism. It helps cells adapt to low oxygen environments by negatively regulating the MAP kinase signaling pathway and dephosphorylating tyrosine and threonine/serine residues on MAP kinases such as ERK, JNK, and p38 MAPK. This lowers the apoptotic threshold after hypoxic injury, enhancing cell survival (Rininger et al., 2012). A decrease in DUSP1 has been shown to be a risk factor for the progression of DKD (Sheng et al., 2019; ).
Single-cell transcriptomics revealed that ZFP36 and DUSP1 were significantly lower expressed in Loop of Henle cells of DKD, which have a high HIR risk score simultaneously. These cells communicated with most other cell types via the EGF pathway. This aligned with previous research findings. The Loop of Henle (LOH) is a critical segment of the kidney involved in reabsorption. Its ascending limb, particularly the thick ascending limb (TAL), efficiently reabsorbs Na+, K+, and Cl−, creating a hyperosmotic environment in the renal medulla, facilitating water reabsorption in the collecting duct (). Due to its deep medullary location, LOH cells are especially sensitive to hypoxic conditions. Hypoxia can lead to metabolic disturbances and dysfunction in these cells, triggering a cascade of hypoxic injury in DKD (). Necrosis of the TAL in isolated perfused kidneys has been considered an important model of renal hypoxia due to its pathological sensitivity to hypoxic injury (Shanley and Johnson, 1989). P. J. Ratcliffe and colleagues established a controlled hypoxia model in perfused rat kidneys to compare the extent of damage between glomerular and tubular cells. They found that the medullary TAL was more sensitive to hypoxic injury compared to other renal segments (Ratcliffe et al., 1988). This characteristic likely explained why these cells uniquely express hub genes associated with hypoxic injury. Notably, the spatial map revealed a striking concordance between regional hypoxia and the spatial distribution of the HIR signature. Within the LH-CD region-anatomically confined to the hyperosmotic renal medulla-marked suppression of DUSP1 and ZFP36 coincides with strong induction of the HIF-1 signalling pathway, confirming that these tubulo-medullary segments are the epicentre of hypoxic stress in DKD.
To bridge the gap between tissue-level hypoxic signatures and clinical applicability, we further validated the prognostic value of the HIR score using large-scale proteomic data from the United Kingdom Biobank. The HIR score significantly improved the discrimination of kidney outcomes, renal death, and secondary endpoints over baseline demographic and metabolic models, with AUC increments of approximately 0.04–0.06 (P < 0.05 via DeLong test, robust across 10-times repeated 5-fold cross-validation). However, in alignment with the realistic constraints of modern renal diagnostics, the incremental predictive gain of the HIR score diminished substantially (approx. +0.006 to +0.015) after the fully adjusted inclusion of conventional clinical pillars (eGFR, urea, and ACR). This marked attenuation honestly indicates a significant informational overlap, suggesting that our proteomic signature shares common pathological covariance with standard functional markers of renal decline. Consequently, the core utility of the HIR score may not lie in completely outperforming or displacing existing non-invasive laboratory standards. Instead, its independent clinical value resides heavily in providing granular, patient-specific mechanistic insights, serving as a molecular lens that reads out cellular-level hypoxic stress patterns and tissue-level microenvironment remodeling, which cannot be captured by routine creatinine or albuminuria surveillance. Phenotypically, the HIR score correlated negatively with eGFR and positively with lipoprotein(a), reinforcing the known interplay between hypoxia-driven tubulointerstitial damage and dyslipidemia in DKD progression. These findings present the HIR signature not as a mere algorithmic replacement for standard tests, but as a pathologically informative risk-stratification tool for developing non-invasive blood- or urine-based mechanistic assays in routine clinical practice.
In the future, the HIR score will be integrated into precision clinical workflows as a three-tiered management pipeline: first, serving as an early-stage “molecular biopsy” to identify high-risk fast-progressors before irreversible eGFR decline occurs, enabling timely, targeted interventions (e.g., SGLT2 inhibitors or HIF stabilizers); second, transitioning into non-invasive liquid biopsy formats, specifically plasma multiplex digital-ELISAs and urine extracellular vesicle mass spectrometry (PRM) assays, to circumvent the procedural constraints of repeat tissue biopsies; and third, functioning as a dynamic companion diagnostic tool to track real-time therapeutic efficacy and guide proactive, individualized treatment adjustments during routine outpatient follow-ups.
This study has certain limitations. First, while our study successfully cross-validated the diagnostic accuracy and directional expression of ZFP36, DUSP1, and CASP3 within the heterogeneous microenvironment of BTBR ob/ob diabetic kidneys, the direct intracellular causalities governed by these vectors remain to be fully characterized. Future mechanistic inquiries are warranted to bridge this gap. We propose to establish a high-glucose combined with hypoxia microenvironment in human renal proximal tubular epithelial cells (HK-2). By deploying stable gain-of-function (plasmid overexpression) and loss-of-function (siRNA/shRNA knockdown) assays for ZFP36, DUSP1, and CASP3 respectively, their precise causal contributions can be delineated.
Second, although we have successfully projected our tissue-derived HIR score onto circulation profiles and verified its prognostic robustness using the large-scale retrospective United Kingdom Biobank proteomic cohort, certain clinical translation barriers remain. The current validation relies on existing public proteomic platforms rather than a customized, clinic-ready diagnostic assay. To further bridge this gap, true non-invasive deployment requires transitioning from retrospective database mining to targeted clinical validation. Third, risk prediction for DKD involved not only gene expression but also methylation modifications, proteomic features, and post-transcriptional regulation. To enhance the reliability of the data, integrating multi-omics approaches would be beneficial for refining the model.
In summary, we identified five key hypoxic injury pathways in DKD development and progression. The HIR score, based on hypoxic injury genes, effectively distinguished between normal controls and DKD patients and predicts renal function changes and kidney events. This provided valuable insights for clinical evaluation and management of DKD.
Statements
Data availability statement
Publicly available datasets were analyzed in this study. This data can be found here: GEO at https://www.ncbi.nlm.nih.gov/gds: GSE96804, GSE30528, GSE30529, GSE104954, GSE99339, GSE131882, GSE195460, GSE261545.
Ethics statement
The studies involving humans were approved by all participants provided electronic informed consent before enrollment, and the protocol received approval from the North West Multi-Centre Research Ethics Committee (reference 11/NW/0382). The United Kingdom Biobank granted data access under application number 564216. The studies were conducted in accordance with the local legislation and institutional requirements. The participants provided their written informed consent to participate in this study.
Author contributions
LJ: Formal Analysis, Data curation, Writing – original draft, Resources, Conceptualization. CC: Data curation, Methodology, Conceptualization, Writing – original draft. HZ: Visualization, Conceptualization, Validation, Project administration, Methodology, Supervision, Writing – original draft. TZ: Software, Visualization, Validation, Funding acquisition, Writing – review and editing, Investigation, Supervision. XW: Funding acquisition, Validation, Resources, Supervision, Methodology, Conceptualization, Writing – original draft.
Funding
The author(s) declared that financial support was received for this work and/or its publication. This work was supported by the the National Natural Science Foundation of China (Nos. 81873140, Nos. 82505310), National High Level Hospital Clinical Research Funding-Excellence & Innovation Initiative of China-Japan Friendship Hospital (No. ZRZC2025-ZYC01) and Elite Medical Professionals Initiative of China-Japan Friendship Hospital (No. ZRJY2025-QM10).
Conflict of interest
The author(s) declared that this work was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Generative AI statement
The author(s) declared that generative AI was not used in the creation of this manuscript.
Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fbinf.2026.1894439/full#supplementary-material
References
1
A/L B Vasanth Rao, V. R.TanS. H.CandasamyM.BhattamisraS. K. (2019). Diabetic nephropathy: an update on pathogenesis and drug development. Diabetes Metab. Syndr.13 (1), 754–762. 10.1016/j.dsx.2018.11.054
2
AfkarianM.ZelnickL. R.HallY. N.HeagertyP. J.TuttleK.WeissN. S.et al (2016). Clinical manifestations of kidney disease among US adults with diabetes, 1988-2014. JAMA316 (6), 602–610. 10.1001/jama.2016.10924
3
AgarwalRD. K. L.DuffinK. L.LaskaD. A.VoelkerJ. R.BreyerM. D.MitchellP. G. (2014). A prospective study of multiple protein biomarkers to predict progression in diabetic chronic kidney disease. Nephrol. Dial. Transpl.29 (12), 2293–2302. 10.1093/ndt/gfu255
4
AtkinsonE. A.BarryM.DarmonA. J.ShostakI.TurnerP. C.MoyerR. W.et al (1998). Cytotoxic T lymphocyte-assisted suicide. Caspase 3 activation is primarily the result of the direct action of granzyme B. J. Biol. Chem.273 (33), 21261–21266. 10.1074/jbc.273.33.21261
5
BechtE.GiraldoN. A.LacroixL.ButtardB.ElarouciN.PetitprezF.et al (2016). Estimating the population abundance of tissue-infiltrating immune and stromal cell populations using gene expression. Genome Biol.17 (1), 218. 10.1186/s13059-016-1070-5
6
BrezisM.ShinaA.KidroniG.EpsteinF. H.RosenS. (1988). Calcium and hypoxic injury in the renal medulla of the perfused rat kidney. Kidney Int.34 (2), 186–194. 10.1038/ki.1988.164
7
BurlakaI. (2022). Analysis of apoptotic, clinical, and laboratory parameters in type 1 diabetes and early diabetic nephropathy: clustering and potential groups evaluation for additional therapeutic interventions. J. Clin. Res. Pediatr. Endocrinol.14 (3), 313–323. 10.4274/jcrpe.galenos.2022.2022-1-21
8
CaiG.ChenX. M. (2016). Prevention, diagnosis, and treatment of chronic kidney diseases in older adults: current status and prospective. Integr. Med. Nephrol. Androl.3 (3), 71–73. 10.4103/2394-2916.187785
9
CaoY.TangW.TangW. (2019). Immune cell infiltration characteristics and related core genes in lupus nephritis: results from bioinformatic analysis. BMC Immunol.20 (1), 37. 10.1186/s12865-019-0316-x
10
ChamboredonS.CiaisD.Desroches-CastanA.SaviP.BonoF.FeigeJ. J.et al (2011). Hypoxia-inducible factor-1α mRNA: a new target for destabilization by tristetraprolin in endothelial cells. Mol. Biol. Cell22 (18), 3366–3378. 10.1091/mbc.E10-07-0617
11
CharoentongP.FinotelloF.AngelovaM.MayerC.EfremovaM.RiederD.et al (2017). Pan-cancer immunogenomic analyses reveal genotype-immunophenotype relationships and predictors of response to checkpoint blockade. Cell Rep.18 (1), 248–262. 10.1016/j.celrep.2016.12.019
12
D'AngeloG. M.RaoD.GuC. C. (2009). Combining least absolute shrinkage and selection operator (LASSO) and principal-components analysis for detection of gene-gene interactions in genome-wide association studies. BMC Proc.3 (Suppl. 7), S62. 10.1186/1753-6561-3-s7-s62
13
FangZ.LiJ.CaoF.LiF. (2022). Integration of scRNA-Seq and bulk RNA-seq reveals molecular characterization of the immune microenvironment in acute pancreatitis. Biomolecules13 (1), 78. 10.3390/biom13010078
14
GraysonP. C.EddyS.TaroniJ. N.LightfootY. L.MarianiL.ParikhH.et al (2018). Metabolic pathways and immunometabolism in rare kidney diseases. Ann. Rheum. Dis.77 (8), 1226–1233. 10.1136/annrheumdis-2017-212935
15
HebertS. C.AndreoliT. E. (1984). Control of NaCl transport in the thick ascending limb. Am. J. Physiol.246 (6 Pt 2), F745–F756. 10.1152/ajprenal.1984.246.6.F745
16
HuC.LiT.XuY.ZhangX.LiF.BaiJ.et al (2023). CellMarker 2.0: an updated database of manually curated cell markers in human/mouse and web tools based on scRNA-seq data. Nucleic Acids Res.51 (D1), D870–D876. 10.1093/nar/gkac947
17
IsnardP.LiD.XuanyuanQ.WuH.HumphreysB. D. (2024). Histopathologic analysis of human kidney spatial transcriptomics data: toward precision pathology. Am. J. Pathol.195 (1), 69–88. 10.1016/j.ajpath.2024.06.011
18
JiangL.ZhangX.WuX. (2025). Hypoxia-induced metabolic reprogramming in the renal tubules of patients with diabetic kidney disease and the associated traditional chinese medicine intervention strategies. IMNA12 (4), e25–00019. 10.1097/IMNA-D-25-00019
19
KaelinW. G. JrRatcliffeP. J. (2008). Oxygen sensing by metazoans: the central role of the HIF hydroxylase pathway. Mol. Cell30 (4), 393–402. 10.1016/j.molcel.2008.04.009
20
KimT. W.YimS.ChoiB. J.JangY.LeeJ. J.SohnB. H.et al (2010). Tristetraprolin regulates the stability of HIF-1alpha mRNA during prolonged hypoxia. Biochem. Biophys. Res. Commun.391 (1), 963–968. 10.1016/j.bbrc.2009.11.174
21
LangfelderP.HorvathS. (2008). WGCNA: an R package for weighted correlation network analysis. BMC Bioinforma.9, 559. 10.1186/1471-2105-9-559
22
LeeP.ChandelN. S.SimonM. C. (2020). Cellular adaptation to hypoxia through hypoxia inducible factors and beyond. Nat. Rev. Mol. Cell Biol.21 (5), 268–283. 10.1038/s41580-020-0227-y
23
LiuF.GuoJ.ZhangQ.LiuD.WenL.YangY.et al (2015). The expression of tristetraprolin and its relationship with urinary proteins in patients with diabetic nephropathy. PLoS One10 (10), e0141471. 10.1371/journal.pone.0141471
24
LiuD.ChenX.HeW.LuM.LiQ.ZhangS.et al (2024). Update on the pathogenesis, diagnosis, and treatment of diabetic tubulopathy. Integr. Med. Nephrol. Androl.11 (4), e23–e29. 10.1097/imna-d-23-00029
25
LuC.WuB.LiaoZ.XueM.ZouZ.FengJ.et al (2021). DUSP1 overexpression attenuates renal tubular mitochondrial dysfunction by restoring Parkin-mediated mitophagy in diabetic nephropathy. Biochem. Biophys. Res. Commun.559, 141–147. 10.1016/j.bbrc.2021.04.032
26
LuoY.ZhangL.ZhaoT. (2023). Identification and analysis of cellular senescence-associated signatures in diabetic kidney disease by integrated bioinformatics analysis and machine learning. Front. Endocrinol. (Lausanne)14, 1193228. 10.3389/fendo.2023.1193228
27
MaY.ZhouX. (2022). Spatially informed cell-type deconvolution for spatial transcriptomics. Nat. Biotechnol.40 (9), 1349–1359. 10.1038/s41587-022-01273-7
28
MutoY.WilsonP. C.LedruN.WuH.DimkeH.WaikarS. S.et al (2021). Single cell transcriptional and chromatin accessibility profiling redefine cellular heterogeneity in the adult human kidney. Nat. Commun.12 (1), 2190. 10.1038/s41467-021-22368-w
29
NangakuM. (2006). Chronic hypoxia and tubulointerstitial injury: a final common pathway to end-stage renal failure. J. Am. Soc. Nephrol.12 (17), 17–25. 10.1681/ASN.2005070757
30
NavaneethanS. D.BansalN.CavanaughK. L.ChangA.CrowleyS.DelgadoC.et al (2024). KDOQI US commentary on the KDIGO 2024 clinical practice guideline for the evaluation and management of CKD. Am. J. Kidney Dis.85 (2), 135–176. 10.1053/j.ajkd.2024.08.003
31
PanY.JiangS.HouQ.QiuD.ShiJ.WangL.et al (2018). Dissection of glomerular transcriptional profile in patients with diabetic nephropathy: sRGAP2a protects podocyte structure and function. Diabetes67 (4), 717–730. 10.2337/db17-0755
32
PercoP.MayerG. (2018). Molecular, histological, and clinical phenotyping of diabetic nephropathy: valuable complementary information. Kidney Int.93 (2), 308–310. 10.1016/j.kint.2017.10.026
33
PfeiferB.HolzingerA.SchimekM. G. (2022). Robust random forest-based all-relevant feature ranks for trustworthy AI. Stud. Health Technol. Inf.294, 137–138. 10.3233/SHTI220418
34
RatcliffeP. J.EndreZ. H.ScheinmanS. J.TangeJ. D.LedinghamJ. G.RaddaG. K. (1988). 31P nuclear magnetic resonance study of steady-state adenosine 5'-triphosphate levels during graded hypoxia in the isolated perfused rat kidney. Clin. Sci. (Lond).74 (4), 437–448. 10.1042/cs0740437
35
RiningerA.DejesusC.TottenA.WaylandA.HaltermanM. W. (2012). MKP-1 antagonizes C/EBPβ activity and lowers the apoptotic threshold after ischemic injury. Cell Death Differ.19 (10), 1634–1643. 10.1038/cdd.2012.41
36
RobinX.TurckN.HainardA.TibertiN.LisacekF.SanchezJ. C.et al (2011). pROC: an open-source package for R and S+ to analyze and compare ROC curves. BMC Bioinforma.12, 77. 10.1186/1471-2105-12-77
37
ShanleyP. F.JohnsonG. C. (1989). Adenine nucleotides, transport activity and hypoxic necrosis in the thick ascending limb of henle. Kidney Int.36 (5), 823–830. 10.1038/ki.1989.268
38
ShaoB.WeiB. (2025). Glycolysis-driven renal fibrosis in chronic kidney disease: emerging mechanisms and therapeutic opportunities. IMNA12 (4), e25–00038. 10.1097/IMNA-D-25-00038
39
ShengJ.LiH.DaiQ.LuC.XuM.ZhangJ.et al (2019). DUSP1 recuses diabetic nephropathy via repressing JNK-Mff-mitochondrial fission pathways. J. Cell Physiol.234 (3), 3043–3057. 10.1002/jcp.27124
40
SinghD. K.WinocourP.FarringtonK. (2008). Mechanisms of disease: the hypoxic tubular hypothesis of diabetic nephropathy. Nat. Clin. Pract. Nephrol.12 (4), 216–226. 10.1038/ncpneph0757
41
SwindellW. R.KruseC. P. S.ListE. O.BerrymanD. E.KopchickJ. J. (2019). ALS blood expression profiling identifies new biomarkers, patient subgroups, and evidence for neutrophilia and hypoxia. J. Transl. Med.17 (1), 170. 10.1186/s12967-019-1909-0
42
TanedaS.HondaK.TomidokoroK.UtoK.NittaK.OdaH. (2010). Eicosapentaenoic acid restores diabetic tubular injury through regulating oxidative stress and mitochondrial apoptosis. Am. J. Physiol. Ren. Physiol.299 (6), F1451–F1461. 10.1152/ajprenal.00637.2009
43
TangC.YangC.WangP.LiL.LinY.YiQ.et al (2024). Identification and validation of glomeruli cellular senescence-related genes in diabetic nephropathy by multi-omics. Adv. Biol. (Weinh)8 (2), e2300453. 10.1002/adbi.202300453
44
TitanS. M.ZatzR.GraciolliF. G.dos ReisL. M.BarrosR. T.JorgettiV.et al (2011). FGF-23 as a predictor of renal outcome in diabetic nephropathy. Clin. J. Am. Soc. Nephrol.6 (2), 241–247. 10.2215/CJN.04250510
45
WangY.OharaT.ChenY.HamadaY.LiC.FujisawaM.et al (2023). Highly metastatic subpopulation of TNBC cells has limited iron metabolism and is a target of iron chelators. Cancers (Basel)15 (2), 468. 10.3390/cancers15020468
46
WenS.WangZ. H.ZhangC. X.YangY.FanQ. L. (2020). Caspase-3 promotes diabetic kidney disease through gasdermin E-Mediated progression to secondary necrosis during apoptosis. Diabetes Metab. Syndr. Obes.13, 313–323. 10.2147/DMSO.S242136
47
WilsonP. C.WuH.KiritaY.UchimuraK.LedruN.RennkeH. G.et al (2019). The single-cell transcriptomic landscape of early human diabetic nephropathy. Proc. Natl. Acad. Sci. U. S. A.116 (39), 19619–19625. 10.1073/pnas.1908706116
48
WilsonP. C.MutoY.WuH.KarihalooA.WaikarS. S.HumphreysB. D. (2022). Multimodal single cell sequencing implicates chromatin accessibility and genetic background in diabetic kidney disease progression. Nat. Commun.13 (1), 5253. 10.1038/s41467-022-32972-z
49
WoronieckaK. I.ParkA. S.MohtatD.ThomasD. B.PullmanJ. M.SusztakK. (2011). Transcriptome analysis of human diabetic kidney disease. Diabetes60 (9), 2354–2369. 10.2337/db10-1181
50
WuJ.ZhangH.LiL.HuM.ChenL.XuB.et al (2020). A nomogram for predicting overall survival in patients with low-grade endometrial stromal sarcoma: a population-based analysis. Cancer Commun. (Lond).40 (7), 301–312. 10.1002/cac2.12067
51
XiaJ.SunL.XuS.XiangQ.ZhaoJ.XiongW.et al (2020). A model using support vector machines recursive feature elimination (SVM-RFE) algorithm to classify whether COPD patients have been continuously managed according to gold guidelines. Int. J. Chronic Obstr. Pulm. Dis.15, 2779–2786. 10.2147/COPD.S271237
52
XuC.LopezR.MehlmanE.RegierJ.JordanM. I.YosefN. (2021). Probabilistic harmonization and annotation of single-cell transcriptomics data with deep generative models. Mol. Syst. Biol.17 (1), e9620. 10.15252/msb.20209620
53
YangR.RenQ. Y.YuanL. S.LiuW.ShiK.ZhouY.et al (2025). Qi-huang-yi-shen formula alleviate tubular epithelial-mesenchymal transition in diabetic kidney disease through HIF-1a/snail signaling pathway. Integr. Med. Nephrol. Androl.12 (1), e24–e40. 10.1097/imna-d-24-00040
54
YinY.TianY.RenX.WangJ.LiX.ZengX. (2022). Qualification of necroptosis-related lncRNA to forecast the treatment outcome, immune response, and therapeutic effect of kidney renal clear cell carcinoma. J. Oncol.2022, 3283343. 10.1155/2022/3283343
55
ZhengX.LiangY.ZhangC. (2023). Ferroptosis regulated by hypoxia in cells. Cells12 (7), 1050. 10.3390/cells12071050
Summary
Keywords
diabetic kidney disease, hypoxia injury, machine learning, single-cell transcriptome analysis, spatial transcriptomics, United Kingdom biobank
Citation
Jiang L, Chieh C, Zhang H, Zhao T and Wu X (2026) Decoding the hypoxic injury landscape and hypoxic risk model construction in diabetic kidney disease: a multi-omics study. Front. Bioinform. 6:1894439. doi: 10.3389/fbinf.2026.1894439
Received
29 May 2026
Revised
05 July 2026
Accepted
09 July 2026
Published
05 August 2026
Volume
6 - 2026
Edited by
Rahul Kaushik, Artificial Intelligence Center for Health and Biomedical Research, Japan
Updates
Copyright
© 2026 Jiang, Chieh, Zhang, Zhao and Wu.
This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Haojun Zhang, 18811701807@163.com; Tingting Zhao, 13001248964@163.com; Xiai Wu, 20180931313@bucm.edu.cn
Disclaimer
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.