Department of Nephrology, The First Affiliated Hospital of Soochow University, Suzhou, People’s Republic of China

Background: Renal failure related death caused by diabetic kidney disease (DKD) is an inevitable outcome for most patients. This study aimed to identify the critical genes involved in the onset and progression of DKD and to explore potential therapeutic targets of DKD.
Methods: We conducted a batch of protein quantitative trait loci (pQTL) Mendelian randomization analysis to obtain a group of proteins with causal relationships with DKD and then identified key proteins through colocalization analysis to determine correlations between variant proteins and disease outcomes. Subsequently, the specific mechanisms of key regulatory genes involved in disease progression were analyzed through transcriptome and single-cell analysis. Finally, we validated the mRNA expression of five key genes in the DKD mice model using reverse transcription quantitative PCR (RT-qPCR).
Results: Five characteristic genes, known as protein kinase B beta (AKT2), interleukin-2 receptor beta (IL2RB), neurexin 3(NRXN3), slit homolog 3(SLIT3), and TATA box binding protein like protein 1 (TBPL1), demonstrated causal relationships with DKD. These key genes are associated with the infiltration of immune cells, and they are related to the regulatory genes associated with immunity. In addition, we also conducted gene enrichment analysis to explore the complex network of potential signaling pathways that may regulate these key genes. Finally, we identified the effectiveness and reliability of these selected key genes through RT-qPCR in the DKD mice model.
Conclusion: Our results indicated that the AKT2, IL2RB, NRXN3, SLIT3, and TBPL1 genes are closely related to DKD, which may be useful in the diagnosis and therapy of DKD.

Keywords: Mendelian randomization analysis, diabetic kidney disease, clinical correlated genes, biomarker, immune cell infiltration

Introduction

Diabetic kidney disease (DKD) is a major chronic kidney disease (CKD) in China and worldwide, featured by an abnormal glomerular filtration rate (GFR), increased urinary protein, and pathological changes in renal microvascular system. Approximately half of DKD patients are estimated to develop end-stage renal disease (ESRD) requiring long-term dialysis,1,2 which may last for a relatively long time under stable conditions.3 In the majority of patients with DKD, however, death caused by renal failure is an unavoidable outcome. Chronic high blood sugar levels can damage the glomeruli and renal tubules, leading to a breakdown of the glomerular filtration barrier and leakage of large molecular proteins into the urine, resulting in proteinuria. In addition, high blood sugar stimulates the production of inflammatory mediators in glomerular cells, causing an inflammatory response that leads to damage to the glomeruli and renal tubules. Moreover, hyperglycemia can also increase intracellular oxidative stress, resulting in increased production of oxygen free radicals, dysfunction of glomerular endothelial cells, proliferation and hypertrophy of mesangial cells, increased mesangial matrix, disappearance of podocyte fusion, thickening of the glomerular basement membrane, atrophy of renal tubular epithelial cells (TECs) and other pathophysiological changes, further damaging renal cells and structures.4,5 At present, the main means of treating DKD include controlling blood sugar, hypertension, blood lipids, inflammation, and proteinuria. In clinical practice, new therapeutic drugs, such as angiotensin-converting enzyme inhibitors (ACEIs) and sodium-dependent glucose co-transporter 2 inhibitor (SGLT2i), have achieved good results in treating DKD. However, the main treatment for ESRD patients is the expensive renal replacement therapy (RRT).6 Therefore, finding sensitive diagnostic biomarkers and new therapeutic targets to achieve early diagnosis and delay the progression of DKD is highly important.

Recent studies have shown that immune regulation is closely associated with the development of glomerulosclerosis, tubulointerstitial fibrosis, and ESRD in a variety of kidney diseases. Moreover, DKD is increasingly recognized as an inflammatory disease in which immune regulation is involved in its development and progression. Immune cells, including macrophages, T cells, B cells and mast cells, infiltrate the kidneys. Many proinflammatory cytokines and chemokines also play essential roles in the pathogenesis of DKD.7 Therefore, exploring the relationship between key genes related to DKD and immune regulation in the occurrence and development of DKD is crucial for identifying new diagnostic biomarkers and therapeutic targets.

Mendelian randomization (MR) uses genetic variation to determine whether observed associations between risk factors and outcomes are consistent with causal effects.8 MR utilizes genetic variation as an instrumental variable to enhance the ability of causal inference.9 An instrumental variable should satisfy three main conditions: being associated with the exposure (correlation hypothesis), having no common cause with the outcome (independence hypothesis) and being associated with the outcome only through exposure (excluding the restrictive hypothesis). Colocalization analysis requires two traits and considers whether their genetic association at a single locus is explained by overlapping or different variants.10

The field of bioinformatics has been developing for a long time. It aims to use information science and statistical collection methods to understand biological phenomena and solve the proposed research challenges.11 There are now many different technologies available for measuring various characteristics of biological systems. For example, different techniques (such as RNA sequencing and RNA microarrays) can be used to measure the same physical and chemical properties, and different physical and chemical properties can also be measured for the same object (such as protein and RNA content in cells). In recent years, genome-scale technology has led to the systematic construction of very large-scale quantitative datasets containing multiple measurement methods.12 In addition, the emerging technology of single-cell RNA sequencing (scRNA-seq) in recent years allows for comprehensive analysis of individual cells, identifying and analyzing specific cellular components and genes that may play a key role in disease progression. It facilitates the classification of cell populations based on their unique characteristics, helping to understand how different cells interact with each other and with their surrounding environment.13 Compared to RNA-seq, scRNA-seq has advantages, in that it includes heterogeneity within the anatomy of cells and identification of rare disease-associated cells by using a single cell profile in a mixture of cells.14 Dynamic changes in gene expression can be experimentally determined for DKD samples through the detection of glomerular cells via scRNA-seq.15 At present, comprehensive bioinformatics analysis has been extensively applied to analyze vast datasets generated from microarray experiments to identify biomarkers or features of various diseases for early monitoring and prognosis assessment. For example, Mona et al identified biomarker genes from a functional network containing 407 genes differentially expressed between healthy control individuals and lung cancer patients in the Gene Expression Omnibus (GEO) dataset of common gene expression. Low expression of 16 gene markers is closely associated with good lung cancer survival, cellular modulation and DNA repair.16

Our research applied comprehensive bioinformatics technology to determine the key genes that have a causal relationship with the progression of DKD. Least Absolute Shrinkage and Selection Operator (LASSO) logistic regression was used for feature selection of diagnostic markers for the disease. The relationships between key genes and immune infiltration and regulatory signaling pathways were further analyzed. Finally, these key genes were visualized in kidney cells through single-cell analysis. Our findings may provide promising targets for the diagnosis and treatment of DKD.

Methods: This study used MR analysis to identify the key genes that are causally associated with the development of DKD. We used Lasso regression algorithm to screen these key genes. The core of this algorithm lies in its ability to determine insignificant variables in the model by minimizing prediction error, thereby producing the most accurate results on a given dataset. Then, Gene Set Enrichment Analysis (GSEA) and Gene Set Variation Analysis (GSVA) were used to further determine the relationships between key genes related to DKD and immune infiltration and further identify prospective molecular mechanisms affecting the progression of DKD. Finally, these key genes were annotated into 12 cell clusters in the kidney by single-cell analysis and visualized. We show the flowchart of the entire study in Figure 1.

Figure 1 Flowchart of the entire study.

Materials and Methods Data Download The GEO Database

We downloaded the Series Matrix File (GSE185011) from NCBI’s public GEO database (https://www.ncbi.nlm.nih.gov/geo/info/datasets.html), along with its annotation file (GPL24676). This dataset includes a total of 10 expression profiles consisting of 5 normal samples and 5 disease samples. We downloaded the single-cell data file of GSE131882 from the NCBI GEO public database. This file included sample data comprising complete single-cell expression profiles from 3 normal individuals and 3 patients, specifically for single-cell analysis.

Exposed pQTL Data

pQTL data were obtained from the deCODE (https://www.decode.com/summarydata/) database. The data utilized in our study were sourced from the prestigious 2021 version of the deCODE database’s pQTL data.17 This comprehensive dataset encompasses genome-wide association studies (GWASs) of 35,559 European individuals, employing 4907 aptamers to meticulously measure plasma protein levels.

Outcome Data

The participants in the GWAS chosen for this study, which are relevant to the outcomes, were predominantly European in origin. The summary outcome data was derived from data contained within the FinnGene R9 database (finngen_R9_DM_NEPHROPATHY_EXMORE). The GWAS catalog encompasses a vast array of publications, highlighted by the most significant associations and complete summary statistics. Additionally, comprehensive summary statistics provide detailed and organized insights into the findings. In total, there were 4111 patients with DKD and 308,539 control individuals.

pQTL MR Analysis

The deCODE (https://www.decode.com) database is recognized as a pacesetter in global efforts to analyze and comprehend the human genome, containing over 150 million sequence variants from different populations around the world. It plays an important role in advancing genomic research, elucidating the connections between genes and traits, and providing new insights into human health and disease. The outcome IDs filtered by the FinnGen biobank database were extracted from the risk (outcome) summary data to obtain the relevant causal relationships in pQTLs. Subsequently, at the selected gene loci, we selected the corresponding single nucleotide polymorphisms (SNPs) as potential instrumental variables (IVs) based on the significance threshold within the site range (P<1e-8). This method not only improves the accuracy of analysis, but also provides more refined and precise genetic information, facilitating a better understanding of the biological mechanisms behind complex diseases. We calculated linkage disequilibrium (LD) between SNPs with R2 < 0.001 (clumping window size=10,000 kb) and only SNPs with a p-value less than 5e-5 were retained for further analysis. The inverse variance weighted (IVW), MR Egger and the weighted median method were used for pQTL MR analysis. The reliability of causality was evaluated by using a statistical method (Wald ratio only when there is only one SNP statistical method in causality) to obtain an overall estimation of the impact of DKD on all cis and certain cross-region genes in whole blood.

Sensitivity Analysis

MR leave-one-out sensitivity analysis was utilized to assess how specific genetic variants influence the risk of DKD. This approach pinpoints and removes variants that affect the overall estimate by excluding each SNP one by one and recalculating the total impact of the remainder. Removing each SNP generates a new point estimate and 95% confidence interval to evaluate the unique impact of that SNP and the robustness of the overall results. We summarize the estimates after removing each SNP individually, as well as the overall estimate when all SNPs are considered. Through comparison of these estimates, we can assess the impact of excluding any single SNP on the results, thereby evaluating the robustness of our analysis.

Colocalization Analysis

In order to further clarify whether the identified disease-related proteins share similar DKD pathogenic variants within genomic regions and to exclude interference from linkage disequilibrium, colocalization analysis was performed with pQTL and GWAS data from DKD samples. We used the 100-kilobase region surrounding the index SNP to calculate the posterior probability. Our results indicate: H3 shows the probability that two traits (gene expression and DKD) are linked but have distinct causal variants. H4 indicates the probability that the two traits are related and share a single causal variant. We set a colocalization threshold of SNP.PP.H4 > 0.95.

The Feature Selection Process of LASSO Regression

A significant feature of Lasso regression is its tendency towards variable selection. When the predictive factors are fully understood and an indispensable part of the research, Ridge regression is more suitable.18 Lasso regression selectively highlights important predictive factors and effectively removes less critical predictive factors by reducing their coefficients to zero.19 We conducted feature selection of diagnostic biomarkers for DKD using Lasso logistic regression. The “glmnet” software package was utilized to apply the Lasso algorithm for feature selection and modeling. During the feature selection process, Lasso regression adjusts the weight of each feature to achieve sparsity (ie, the coefficients of some features become zero) by minimizing the sum of the loss function and the L1 regularization term. This characteristic makes Lasso regression highly effective for high-dimensional datasets, and it can filter out features that have a significant impact on the target variable. By adjusting the regularization parameter (lambda value) in the Lasso model, we can control the sparsity of the features, thereby selecting an appropriate number of features for modeling and avoiding overfitting.

Immune Cell Infiltration Analysis

The CIBERSORT approach has been widely applied to evaluate the types of immune cell in microenvironment. It comprises 547 biomarkers that can distinguish 22 types of immune cell. We used the CIBERSORT algorithm to analyze the sample data, estimate the relative abundance of 22 types of infiltrating immunocytes, and examine the connection between immune cell content and gene expression.

GSEA

GSEA is designed to analyze the expression patterns of genes within different gene sets. By employing pre-established set of genes, this algorithm systematically ranks each gene based on its differential abundance across two distinct sample types. The process is to check the preset gene set is significantly enriched at the top or bottom of the ranking. GSEA was employed to compare signal differences between the high-expression and low-expression groups, and to identify the molecular mechanisms of the five genes in two groups. We set the substitution number to 1000 and determined phenotype as the substitution type.

GSVA

GSVA is a groundbreaking and innovative approach to the study of transcriptional gene set enrichment. By comprehensively scoring a particular gene set and then determining its biological function, the GSVA transforms gene-level changes to pathway-level changes. In this study, gene sets were downloaded from the Molecular Signatures Database (v7.0 version) and the GSVA algorithm was used, which comprehensively evaluated each individual gene set from an exhaustive perspective to assess underlying biological functional changes in different samples.

Single-Cell Analysis

First, the expression profile was read through the Seurat package, and low-expression genes were filtered out. Standardization, homogenization, principal component analysis (PCA), and analysis of the sequence data were performed. The optimal number of pcs was determined through ElbowPlot, and the positional relationship between each cluster was determined through uniform manifold approximation and projection (UMAP) analysis. The cluster was annotated through known cell markers, and each was annotated to some cells that are important for the occurrence of the disease.

Animal Experiments

In our current research, we used a mouse model of streptozotocin (STZ)-induced DKD. Long-term hyperglycemia was constructed by intraperitoneal injection of 50 mg/kg STZ (pH: 4.2–4.5) into 8-week-old C57BL/6J mice for five consecutive days. An equal volume of sterile sodium citrate was injected into mice of the same age to serve as the control group.20 One week later, mice with blood glucose concentrations exceeding 16.7 mmol/L were included in the experimental group.21,22 During the entire modeling period, a Roche blood glucose meter was used to measure blood sugar every two weeks. Blood and tissue samples were obtained from DKD mice after 12 weeks of persistent hyperglycemia. After 12 weeks of continuous induction of hyperglycemia, we collected 24-hour urine from mice using metabolic cages, euthanized the mice, and collected blood and kidney samples. We fixed mice kidney tissue with 4% paraformaldehyde and performed routine paraffin embedding and sectioning (5 µm), followed by Periodic Acid-Schiff (PAS) staining and Masson staining. Mouse microalbuminuria ELISA Kit (Elabscience, Wuhan, China, E-EL-M0792) was used to monitor the protein content in mice urine and creatinine colorimetric assay kit (Elabscience, Wuhan, China, E-BC-K188-M) was used to monitor the serum creatinine of DKD mice. We also used RT‒PCR to quantitatively detect the expression of five key genes in mice kidney tissue. Following the manufacturer’s instructions, total RNA was extracted from the samples using an RNA-easy reagent. After total RNA was quantified by a NanoDrop2000 (Thermo Fisher, CA, USA), a cDNA synthesis kit (Vazyme, Nanjing, China, R333) and Taq Pro Universal SYBR for qPCR Master Mix (Vazyme, Nanjing, China, Q712) were used for reverse transcription and real-time fluorescence quantification, respectively (RT‒qPCR).

Statistical Analysis

All statistical analyses were conducted in the R-language (version 4.2), with a p-value of less than 0.05 considered statistically significant.

Results Identification of Candidate Genes Using MR Analysis

We downloaded relevant data on DKD from the deCODE database. The outcome id derived from the summary statistics of 312,650 cases (controls: 308,539; cases: 4111) related to DKD was as finngen_R9_DM_NEPHROPATHY_EXMORE. Instruments and extract_outcome_data were sequentially extracted to determine the causal relationships related to the pQTL outcome. We screened the causal relationships of 20 pairs of genes linked to pQTL-positive outcomes using MR analysis (Figure 2, IVW pval < 0.01). Then, we assessed their reliability by conducting a leave-one-out sensitivity analysis. The results showed that it was no significant influence on the overall error range to exclude any single SNP, confirming the robustness of the 20 selected pairs of causal relationships (Figure 3).

Figure 2 MR analysis on screening the causality of 20 pairs of genes corresponding to pQTL-positive outcomes.

Figure 3 Sensitivity analysis of the causality of the 20 genes using the leave-one-out method. MR leave−one−out sensitivity analysis for (A)AKT2,(B)CHGB,(C)CHST11,(D)DNER,(E)EGFL6,(F)EIF3G,(G)ERBB3,(H)FSTL1,(I)GP1BA,(J)IGHG4,(K)IL2RB,(L)NRXN3,(M)RACGAP1,(N)RTP4,(O)SLIT3,(P)TBPL1,(Q)TFRC,(R)TNNI2,(S)UGDH,(T)UNC5D.

Identifying the Reliability of Candidate Genes Using Colocalization Analysis

In addition, we performed colocalization analysis on 20 candidate genes at the pQTL-GWAS level, among which the colocalized SNP.PP.H4 values of the genes AKT2, DNER, EIF3G, IL2RB, NRXN3, RACGAP1, SLIT3, TBPL1, TNNI2, and UNC5D were greater than 0.95 (Figure 4). In order to further elevate the exactitude of the identified key genes of DKD, we used the Lasso regression algorithm to screen again the 10 genes identified by colocalization analysis (Figure 5A and B). The results showed that the Lasso regression algorithm narrowed down the range of selected feature genes to five, which will be key genes for subsequent analysis, namely AKT2, IL2RB, NRXN3, SLIT3, and TBPL1.

Figure 4 Colocalization analysis of 20 candidate genes at the pQTL-GWAS level (A)AKT2,(B)DNER,(C)EIF3G,(D)IL2RB,(E)NRXN3,(F)RACGAP1,(G)SLIT3,(H)TBPL1,(I)TNNI2,(J)UNC5D.

Figure 5 (A and B) Screening key candidate hub genes through LASSO regression analysis.

Analysis of the Regulatory Mechanisms of Key Genes Involved in DKD

The microenvironment, a complex network of cells and tissues that surrounds and interacts with the body in a dynamic fashion, is largely constituted by immune cells, the extracellular matrix, various growth factors, inflammatory factors and has special physical and chemical characteristics. It is essential for maintaining a balance between immune activation and the body’s natural defense mechanisms. Moreover, this microenvironmental realm is characterized by unique physical and chemical properties that can profoundly impact disease diagnosis and patient response to clinical treatment. This study further explored the potential molecular mechanisms by which key genes influence the progression of DKD by analyzing the relationships between key genes and immune infiltration in a DKD dataset. Our study explored the proportions of immune cells in each sample and the correlations between different immune cell types (Figure 6A and B). In addition, the results showed significant differences in naïve B cells and resting Mast cells between the control samples and disease group samples (Figure 6C).

Figure 6 (A and B) The proportions of immune cells in each sample and the correlation among different types of immune cell. (C) Differences in immune cell content between control and disease samples. (D) The relationships between five key genes and immune cells.

Immune Infiltration Analysis of Key Genes

This study further explored the relationships between key genes and immune cells and revealed that key genes are highly correlated with immune cells (Figure 6D). Among them, AKT2 was positively correlated with immune cells such as resting mast cells and negatively correlated with immune cells such as neutrophils. IL2RB was positively correlated with immune cells such as CD8 T cells and negatively correlated with immune cells such as neutrophils. NRXN3 was positively correlated with immune cells such as activated dendritic cells. TBPL1 was positively related to immune cells such as memory B cells. The correlations between these five key genes and different immune factors, including immune regulatory factors, chemokines, and cell receptors, were obtained from the TISIDB database (Figure 7). These analyses suggest that key genes are closely related to the level of immune cell infiltration and play an important role in the immune microenvironment.

Figure 7 Correlations between these five key genes and different immunity factors. (A-E) Correlations between key genes and chemokines, immunoinhibitors, immunostimulators, MHC and receptors.

Analysis of Signaling Regulatory Pathways Involving Key Genes

Next, we investigated the specific signaling pathways associated with the five key genes and delved into how these genes might influence disease progression. Specifically, our GSEA analysis showed that AKT2 is involved in two significant pathways: the interleukin 1-mediated signaling pathway and lipopolysaccharide-mediated signaling pathway (Figure 8A). The pathways enriched by IL2RB included the collagen activated pathway and T cell receptor signaling pathway (Figure 8B). The NRXN3-enriched pathways included the negative regulation of the BMP, the retinoic acid receptor and other signaling pathways (Figure 8C). SLIT3-enriched pathways included the negative regulation of the TORC1 signaling and the neurotrophin TRK receptor signaling pathways (Figure 8D). The pathways enriched by TBPL1 included the hexose catabolic process, monosaccharide catabolic process and other pathways (Figure 8E).

Figure 8 The potential molecular mechanisms related to the five key genes. Panels (A-E) show AKT2, IL2RB, NRXN3, SLIT3 and TBPL1, respectively.

GSVA analysis showed that high expression of AKT2 promotes key signaling pathways, including the IL2 STAT5 and MTORC1 pathways (Figure 9A). Similarly, elevated levels of IL2RB promote the MTORC1 and P53 pathways (Figure 9B). High expression of NRXN3 can promote bile acid metabolism and interferon alpha response (Figure 9C). The elevated levels of SLIT3 can promote signaling pathways known as the apoptosis and allograft rejection pathways (Figure 9D). High expression of TBPL1 can promote the MTORC1, interferon alpha response and other signaling pathways (Figure 9E).

Figure 9 Multiple signaling regulatory pathways involving five key genes related to DKD.(A–E) Key genes were analyzed by GSVA. Blue represents the signaling pathway associated with high-level gene expression, and green represents the signaling pathway associated with low-level gene expression. Figures (A-E) show AKT2, IL2RB, NRXN3, SLIT3 and TBPL1 in sequence.

Correlation Analysis of the Key Genes and Key Regulatory Genes of DKD

We have meticulously gathered and reviewed a diverse array of immunity-related regulatory genes from GeneCards (https://www.genecards.org/), an expansive database that provides comprehensive information on genetic markers and their functional significance. We analyzed the expression levels of genes with the highest correlation scores and found that MYD88 and SOCS1 were the most expressed in DKD. There were differences in control group and disease group (Figure 10A). In addition, correlation analysis on key genes and immune regulatory genes was performed. There was a significant correlation between the expression of key genes and the level of disease regulation genes. Among these genes, IL2RB was significantly positively connected with ATM (r= 0.858), and AKT2 was significantly negatively associated with PLCG2 (r=−0.893) (Figure 10B).

Figure 10 (A) Differential expression of immune related regulatory genes between control group and disease group. (B) Correlation analysis of five key genes and the expression levels of immune regulation-related genes in DKD.

Single-Cell Analysis of Key Genes and Renal Cell Clusters

The GSE131882 single cell data set was downloaded from the NCBI GEO public database. First, using the Seurat package to read the expression profile and the low-expression genes were filtered out. We screened data samples using nFeature RNA, nCount RNA and mitochondrial content (nFeature RNA>50 and percentage. mt 10) (Figure S3A), and then found that the batch effect between samples was not significant through PCA downscaling analysis (Figure S3B). Meanwhile, the optimum PC number: 11 (Figure S3C) was obtained by standardization, normalization, PCA, and homogenization analysis, and the position relation of each cluster was acquired by UMAP analysis (Figure 11A). These 12 clusters annotated through known cell markers (Figure S4) were identified as 12 cell types, including LEUK, MES, ENDO, PODO, CD-ICAI, CD-ICB, CD-PC, DCT/CT, DCT, LOH, CFH, PCT (Figure 11B). We analyzed the expression levels of key genes in these 12 cell clusters (Figure 11C and D). Next, we obtained genes related to DKD progression (ACE and PVT1) from the GeneCards database (https://www.genecards.org/) and visualized the coexpression of these genes with five key genes in 12 types of cells (Figures S1 and S2). In addition, we supplemented the expression of key genes between the control group and DKD group (Figure S5), as well as the ROC curves of the diagnostic efficacy of these five key genes for DKD (Figure S6). Finally, we used AUCell function to perform quantitative analysis of the levels of genes related to immunity, metabolism, signaling pathways and proliferation (Figure 11E). AKT2 is highly expressed in signaling pathways such as oxidative phosphorylation, apoptosis, mtorc1 signaling, and unfolding protein response. IL2RB is significantly upregulated in the allograft-injection related immune pathway. NRXN3 is significantly downregulated in oxidative phosphorylation, unfolded protein response, while significantly upregulated in estrogen-response early, PI3K/AKT/mTOR signaling, and mitotic spike pathways. SLIT3 is significantly downregulated in oxidative phosphorylation signaling and significantly upregulated in complement and epithelial mesenchymal transition signaling. TBPL1 is highly expressed in regulatory processes such as protein secretion and mitotic spindle.

Figure 11 Single-nucleus RNA sequencing. (A) The tSNE algorithm was used to classify the cells into 16 clusters according to the key components of the PCA. (B) The 16 clusters were annotated with 12 cell types: LEUK, MES, ENDO, PODO, CD-ICB, CD-ICA, CD-PC, DCT/CT, DCT, LOH, CFH, PCT. (C) and (D) Expression profiles of key genes in these cells were analyzed. (E) Analysis of the levels of genes related to immunity, metabolism, signaling pathways and proliferation.

Abbreviations: LEUK, leukocytes; MES, mesangial cells; ENDO, endothelial cells; PODO, podocyte; CD-ICA, connecting tubule-type A intercalated cell; CD-ICB, connecting tubule-type B intercalated cell; CD-PC, connecting tubule-principle cell; DCT/CT, DCT, distal convoluted tubule; LOH, loop of Henle; CFH, complement factor H; PCT, proximal convoluted tubule.

Validation in Animal Models

In order to confirm the validity of these five markers in early DKD diagnosis, we constructed a model of DKD induced by STZ. Compared with the control group, the proliferation of renal mesangial cells was increased, the mesangial matrix expanded more significantly, and the abnormal glomerular basement membrane thickening was observed in DKD group of PAS staining. Masson staining demonstrated that the DKD mice had significantly increased the deposition of renal tissue fibers (Figure 12A), which showed the successful establishment of DKD model. Subsequently, we evaluated all glomeruli under the microscope of each mouse. We evaluated the proportion of diseased glomeruli to the entire glomeruli to obtain the mesangial matrix proliferation index. The results showed that the degree of kidney disease in DKD mice was more severe, which was statistically significant compared to the control group of normal mice (Figure 12B). In addition, at the end of modeling, the urinary albumin (Figure 12C) and serum creatinine levels (Figure 12D) of DKD group mice were significantly increased compared to the control group mice. The above experimental results indicate the successful construction of the DKD mouse model. Finally, we detected the mRNA expression levels of AKT2, IL2RB, NRXN3, SLIT3, and TBPL1 molecules in the renal tissue of the DKD mice. The results indicated that in the DKD model, AKT2, IL2RB, and NRXN3 were significantly increased, while TBPL1 expression was significantly reduced. However, in the DKD group, SLIT3 showed a decreasing trend and there was no difference between the two groups (Figure 12E).

Figure 12 Verification of five markers in animal experiments. (A) PAS and Masson staining of mouse kidney. (B) Quantitative analyses of mesangial matrix expansion (n = 6 mice per group). Scale bar, 20μm; original magnification×400. (C) 24-hour albuminuria levels in DKD mice at week 12 of modeling. (D) Blood creatinine levels in DKD mice at week 12 of modeling. (E) mRNA expression levels of five genes in kidney tissue of mice.

Discussion

The occurrence and development of DKD involves multiple factors.23 However, its specific mechanisms remain to be explored. At present, the therapeutic effect of DKD treatment is limited due to individual heterogeneity, so it is necessary to explore and develop novel molecules that can contribute to DKD diagnosis and therapy. Based on the genetic data for 312,650 individuals with DKD (control individuals: 308,539; patients: 4111), our study provides robust evidence that 20 proteins (IVW pval < 0.01) are causally associated with DKD, and five of them showed evidence of genetic colocalization with DKD outcomes. Our research findings highlight potential targets for future treatment of DKD and demonstrate the relevance of genomics and proteomics in identifying drug targets. In addition, we employed comprehensive methodologies to investigate the genetic correlations of these key genes with immune infiltration in DKD patients. We further investigated the related signaling pathways and regulatory mechanisms modulated by these key genes to investigate the possible molecular mechanisms that are involved in DKD progression. Finally, we annotated clusters of known cell markers through single-cell analysis and identified five genes that are important for disease occurrence.

Previous studies have shown that the PI3K/AKT pathway participates in various biological processes involved in the development of DKD, such as inflammation, oxidative stress, and apoptosis.24–26 The AKT kinase family, consisting of AKT1, AKT2, and AKT3, is highly homologous with distinct functional specificity and tissue distribution.27,28 AKT1, which is widely expressed in all kinds of tissues, plays a biological role in regulating cell growth and survival.29 AKT2 has been found to regulate glucose metabolism and is mainly expressed in insulin-responsive tissues, such as skeletal muscle and adipose tissue.30 AKT3, which exhibits high expression levels specifically in nervous tissues, predominantly modulates neuronal development.31 Notably, a study has reported a missense mutation of AKT2 with serious diabetes mellitus.32 Therefore, specifically targeting AKT2 subtypes may be a feasible therapeutic method for the treatment of glucose metabolism disorders.33

Interleukin-2 receptor beta (IL2RB) is an interleukin-2 receptor. At present, some studies have suggested that IL2RB may be involved in T-cell-mediated immune response regulation and plays an important role in maintaining immune homeostasis.34,35 IL2RB has been reported to be a pivotal gene related to the function of T cells in DKD.36 Our research provides new evidence of IL2RB in the immune regulation of DKD.

The human neurexin family is composed of NRXN1, NRXN2, and NRXN3.37 Neurexin plays important roles in synaptic neurotransmitter release and cell adhesion and has been found to be highly expressed in presynaptic nerve terminals.38 Previous studies have shown that NRXN3 is involved in metabolism and obesity regulation and is related to diabetes.39–41 However, whether NRXN3 is involved in the progression of DKD is still unclear. Our research indicated that there is a causal relationship between NRXN3 and DKD and revealed the pathway by which NRXN3 exerts its regulatory function and its relationship with immune regulation.

Slit guidance ligand proteins, referred to as SLITs, include several subtypes, such as SLIT1, SLIT2, and SLIT3, that were discovered to be highly conserved in several species.42 Moreover, these proteins have been shown to significantly modulate multiple cellular biological processes in a variety of tissue types, including kidney tissue.43,44 Notably, SLIT3 is associated with various cellular functions and processes, including axonal guidance, angiogenesis,45–47 inflammatory cell chemotaxis,48–51 and tumor cell metastasis.52,53 Previous studies have shown that SLIT3-deficient mice exhibit reduced expression of type I and III collagen at the transcriptional level, and an overall decrease in collagen content was found in various nonneuronal tissues studied, including the aortic outer membrane, lungs, spleen, kidneys, bones, skin, and heart.54 These findings indicate its significant influence on fibroblast function. SLIT3 deficiency can weaken the production of collagen in the heart and other nonneuronal tissues. It plays an important role in regulating fibroblast activity and fibrous collagen synthesis through autocrine regulation, making it a potential therapeutic target for fibrotic diseases.54Sonia Zambrano et al used SMARTseq2 technology to perform single-cell RNA sequencing (scRNA-seq) analysis of glomeruli in mice with early IgA nephropathy. Several paracrine pathways within the glomerulus, such as mesangial cell-derived Slit3, which may activate Robo receptors in podocytes/endothelial cells, were detected, and the functions of key cell‒cell crosstalk pathways, such as the role of the Slit-Robo signaling axis, were further validated through cell experiments.55 Our research revealed that the Slit3 gene actively participates in the immune regulatory pathway and inflammatory response in DKD and is involved in extensive communication with various immune cells.

Transcription initiation (TI) is a crucial and highly regulated event in the biological cycle of organisms, and tata box binding protein (TBP) plays an important role in TI. TBP interacts specifically with the DNA sequence of promoters of most Class II genes and some Class III genes (such as the TATA box).56 TBPL1 (TRF2) is a metazoan TBP analog with 40% homology to TBP. It plays an important transcriptional role in the embryogenesis of zebrafish, Xenopus, and Drosophila.57 It is crucial for mouse spermatogenesis but has no effect on mouse embryogenesis.58–61 TBPL1 binds to GTF IIA subunit 1 like (GTF2A1L) and is present in the cytoplasm with several heat shock protein complexes in mice. When TBPL1 is recruited as an activating gene promoter in haploid cells, it drives the expression of TAF7L-related genes, thereby establishing specific complexes to activate certain genes. The present study revealed that the C-MYC inhibitor 10,058-F4 upregulated TBPL1 expression in lung tissue-derived cells of mice with idiopathic pulmonary fibrosis (IPF) and inhibited IPF.62 In addition, sequencing of human pulmonary carcinoid tumors revealed significant mutations in TBPL1.63 Two miRNAs, miR-18a and miR-133b, directly target TBPL1 expression and play a protective role in colorectal cancer cells (CRCs).64,65 Our study revealed that TBPL1 may also play an important regulatory role in the occurrence and development of DKD.

Oxidative stress and the inflammatory response have been recognized as important processes in the pathogenesis of DKD;66,67 therefore, manipulating the immune system could present potential therapeutic benefits. SGLT2i is currently a commonly used drug in clinical practice, and the typical changes in ferroptosis include significant lipid peroxidation, impaired antioxidant capacity, and iron overload. SGLT2i treatment reduces overactivation of the HIF1α/HO1 axis, significantly alleviating inflammation and oxidative stress associated with ferroptosis.68 It is also now believed that the beneficial effect of SGLT2i in the treatment of diabetic complications currently depends on its anti-inflammatory and anti-immune effects.69–71 Additionally, ACEIs and angiotensin II receptor blockers (ARBs) have been proved to reduce the infiltration of renal macrophages in animal models of diabetes.72,73 Experimental and clinical evidence shows that the complement cascade can lead to damage in diseases not traditionally considered immune-mediated diseases, including DKD and focal segmental glomerulosclerosis. Many complement inhibition drugs have been approved, and other drugs targeting different components of complement cascade reactions are currently being studied in clinical trials.74 A typical example is the in vivo delivery of oral active C5aR1 inhibitors (PMX53), which reverse the phenotypic changes in DKD and normalize the renal mitochondrial fatty acid profile, cardiac phospholipid remodeling, and citrate cycle intermediates. In addition, the exposure of human renal proximal TECs to C5a in vitro leads to alterations in mitochondrial respiratory function and the production of reactive oxygen species.75 Tissue-resident memory T (RM) cells are a newly discovered group of noncirculating memory T cells that play a crucial role in mediating local immune responses and preventing local reinfection by pathogens. Mechanistically, Sparsentan can reduce the T (RM) cell response by interfering with IL-15 signaling. Therefore, targeting T (RM) cells may be a novel therapeutic method for patients with glomerular disease.76 Recently, it has been widely recognized that inflammatory immunity is essential in the progression of DKD. This regulation process involves infiltration of immune cells into the kidney, elevated chemokine and pro-inflammatory cytokine production. This study may be helpful for identifying new prognostic and therapeutic targets of DKD, and numerous therapeutic approaches for DKD are currently being studied in clinical trials for different aspects of kidney inflammation and immunity. The hope is that these new targets, when combined with existing treatments, will be beneficial in treating diseases.7

Based on genetic data for 312,650 individuals with DKD (control individuals: 308,539; patients: 4111), our study provides robust evidence that 20 proteins (Figure 1A; IVW pval < 0.01) are causally associated with DKD, and five of them showed evidence of genetic colocalization with DKD outcomes. Our research findings highlight potential targets for future treatment of DKD and demonstrate the relevance of proteomics in identifying drug targets. However, further research is needed to evaluate the feasibility of using five identified proteins as drug targets for DKD treatment. With the emergence of more comprehensive proteomics platforms and increasing research on more diverse non-European populations, it is possible to discover more DKD drug targets.

Conclusions

Through our innovative MR analysis combined with transcriptomic and single-cell sequencing analyses, we revealed the relationships between five genes and DKD, as well as their relationships with immune infiltration, and further identified the transcription regulatory factors and intercellular pathways involved. Based on our newly identified DKD mechanism, five new targets have been recognized for future DKD conversion and treatment research. However, the specific mechanisms by which these genes influence DKD remain elusive. The targeting of these key genes in drug development has potential for successful clinical trials, providing novel avenues for DKD treatment.

Data Sharing Statement

All the data are included in the manuscript and supplementary materials.

Ethics Statement

The Ethical Committee (ethics approval ID: 20240301), First Affiliated Hospital, Soochow University approved this study. All animal experiments were conducted in accordance with the guidelines of the National Institutes of Health in the United States and approved by the Animal Ethics Committee of Soochow University (ethics approval ID: 202406A0011).

Author Contributions

All authors made a significant contribution to the work reported, whether that is in the conception, study design, execution, acquisition of data, analysis and interpretation, or in all these areas; took part in drafting, revising or critically reviewing the article; gave final approval of the version to be published; have agreed on the journal to which the article has been submitted; and agree to be accountable for all aspects of the work.

Funding

This work was supported by National Natural Science Foundation of China (81600565 and 81700129), the Natural Science Foundation of Jiangsu Province for Distinguished Young Scholars (BK20190052), the Science and Technology Plan of Suzhou City (SKY2022042) and Science and Technology Plan of Suzhou Health Commission (MSXM2024005).

Disclosure

The authors report no conflicts of interest in this work.

References

1. Guo J, Liu Z, Gong R. Long noncoding RNA: an emerging player in diabetes and diabetic kidney disease. Clin Sci. 2019;133(12):1321–1339. doi:10.1042/CS20190372

2. Tang SCW, Yiu WH. Innate immunity in diabetic kidney disease. Nat Rev Nephrol. 2020;16(4):206–222. doi:10.1038/s41581-019-0234-4

3. Sarnak MJ, Levey AS, Schoolwerth AC, et al. Kidney disease as a risk factor for development of cardiovascular disease: a statement from the American Heart Association Councils on kidney in cardiovascular disease, high blood pressure research, clinical cardiology, and epidemiology and prevention. Circulation. 2003;108(17):2154–2169.

4. Wada J, Makino H. Inflammation and the pathogenesis of diabetic nephropathy. Clin Sci. 2013;124(3):139–152. doi:10.1042/CS20120198

5. Thomas MC, Brownlee M, Susztak K, et al. Diabetic kidney disease. Nature Reviews Disease Primers. 2015;1(1):15018. doi:10.1038/nrdp.2015.18

6. Selby NM, Taal MW. An updated overview of diabetic nephropathy: diagnosis, prognosis, treatment goals and latest guidelines. Diabetes Obesity Metab. 2020;22(Suppl 1):3–15. doi:10.1111/dom.14007

7. Hickey FB, Martin F. Diabetic kidney disease and immune modulation. Curr Opin Pharmacol. 2013;13(4):602–612. doi:10.1016/j.coph.2013.05.002

8. Birney E. Mendelian Randomization. Cold Spring Harbor Perspect Med. 2022;12(4). doi:10.1101/cshperspect.a041302

9. Little M. Mendelian randomization: methods for using genetic variants in causal estimation. J Royal Statis Soc Series A. 2018;181(2):549–550. doi:10.1111/rssa.12343

10. Zuber V, Grinberg NF, Gill D, et al. Combining evidence from Mendelian randomization and colocalization: review and comparison of approaches. Am J Hum Genet. 2022;109(5):767–782. doi:10.1016/j.ajhg.2022.04.001

11. Uesaka K, Oka H, Kato R, et al. Bioinformatics in bioscience and bioengineering: recent advances, applications, and perspectives. J Biosci Bioeng. 2022;134(5):363–373. doi:10.1016/j.jbiosc.2022.08.004

12. Miao Z, Humphreys BD, McMahon AP, Kim J. Multi-omics integration in the age of million single-cell data. Nat Rev Nephrol. 2021;17(11):710–724. doi:10.1038/s41581-021-00463-x

13. Ding S, Chen X, Shen K. Single-cell RNA sequencing in breast cancer: understanding tumor heterogeneity and paving roads to individualized therapy. Cancer Commun. 2020;40(8):329–344. doi:10.1002/cac2.12078

14. Hafemeister C, Satija R. Normalization and variance stabilization of single-cell RNA-seq data using regularized negative binomial regression. Genome Biol. 2019;20(1):296. doi:10.1186/s13059-019-1874-1

15. Fu J, Akat KM, Sun Z, et al. Single-cell RNA profiling of glomerular cells shows dynamic changes in experimental diabetic kidney disease. J Am Soc Nephrol. 2019;30(4):533–545. doi:10.1681/ASN.2018090896

16. Maharjan M, Tanvir RB, Chowdhury K, Duan W, Mondal AM. Computational identification of biomarker genes for lung cancer considering treatment and non-treatment studies. BMC Bioinf. 2020;21(Suppl 9):218. doi:10.1186/s12859-020-3524-8

17. Ferkingstad E, Sulem P, Atlason BA, et al. Large-scale integration of the plasma proteome with genetics and disease. Nature Genet. 2021;53(12):1712–1721. doi:10.1038/s41588-021-00978-w

18. Tibshirani R. Regression shrinkage and selection via the lasso. J Royal Statis Soc. 2018;58(1):267–288. doi:10.1111/j.2517-6161.1996.tb02080.x

19. Ali H, Shahzad M, Sarfraz S, Sewell KB, Alqalyoobi S, Mohan BP. Application and impact of Lasso regression in gastroenterology: a systematic review. Indian J Gastroenterol. 2023;42(6):780–790. doi:10.1007/s12664-023-01426-9

20. Li A, Yi B, Han H, et al. Vitamin D-VDR (vitamin D receptor) regulates defective autophagy in renal tubular epithelial cell in streptozotocin-induced diabetic mice via the AMPK pathway. Autophagy. 2022;18(4):877–890. doi:10.1080/15548627.2021.1962681

21. Shahzad K, Bock F, Dong W, et al. Nlrp3-inflammasome activation in non-myeloid-derived cells aggravates diabetic nephropathy. Kidney Int. 2015;87(1):74–84. doi:10.1038/ki.2014.271

22. Shahzad K, Bock F, Al-Dabet MM, et al. Caspase-1, but not caspase-3, promotes diabetic nephropathy. J Am Soc Nephrol. 2016;27(8):2270–2275. doi:10.1681/ASN.2015060676

23. Alicic RZ, Rooney MT, Tuttle KR. Diabetic kidney disease: challenges, progress, and possibilities. Clin J Am Soc Nephrol. 2017;12(12):2032–2045. doi:10.2215/CJN.11491116

24. Wang H, Gao L, Zhao C, et al. The role of PI3K/Akt signaling pathway in chronic kidney disease. Int Urol Nephrol. 2024;56(8):2623–2633. doi:10.1007/s11255-024-03989-8

25. Wang N, Zhang C. Oxidative stress: a culprit in the progression of diabetic kidney disease. Antioxidants. 2024;13(4). doi:10.3390/antiox13040455

26. Zang L, Saitoh S, Katayama K, Zhou W, Nishimura N, Shimada Y. A zebrafish model of diabetic nephropathy shows hyperglycemia, proteinuria and activation of the PI3K/Akt pathway. Dis Models Mech. 2024;17(5). doi:10.1242/dmm.050438

27. Gonzalez E, McGraw TE. The Akt kinases: isoform specificity in metabolism and cancer. Cell Cycle. 2009;8(16):2502–2508. doi:10.4161/cc.8.16.9335

28. Urasaki Y, Beaumont C, Talbot JN, Hill DK, Le TT. Akt3 regulates the tissue-specific response to copaiba essential oil. Int J mol Sci. 2020;21(8):2851. doi:10.3390/ijms21082851

29. Chen WS, Xu PZ, Gottlob K, et al. Growth retardation and increased apoptosis in mice with homozygous disruption of the Akt1 gene. Genes Dev. 2001;15(17):2203–2208. doi:10.1101/gad.913901

30. Cho H, Mu J, Kim JK, et al. Insulin resistance and a dia