Identification and mechanistic insights of ubiquitin-proteasome system and pyroptosis-related biomarkers in type 2 diabetes mellitus

Abstract

BACKGROUND

Pyroptosis and ubiquitination have been identified as key processes influencing the pathogenesis of diabetes mellitus (DM).

AIM

To investigate the genes associated with the ubiquitin-proteasome system (UPS) and pyroptosis in type 2 DM (T2DM), and elucidate their mechanisms of action in T2DM.

METHODS

The datasets GSE76894, GSE41762, and GSE86469 were utilized in this study. UPS-related genes (UPSGs) and pyroptosis-related genes (PRGs) were obtained from existing literature. Differential expression analysis was performed to identify differentially expressed genes (DEGs). DEGs were intersected with UPSGs and PRGs to identify differentially expressed UPSGs and PRGs. Ubiquitin-pyroptosis-related biomarkers were determined using Spearman’s correlation, t-tests, and receiver operating characteristic curve analysis. Pathway enrichment of biomarkers was assessed using Gene Set Enrichment Analysis (GSEA). Single sample GSEA (ssGSEA) and Spearman’s correlation were used to analyze the relationship between biomarkers and immune cells. A competitive endogenous RNA network was constructed. Subsequently, drugs related to the biomarkers were identified and a gene-drug network was established. In dataset GSE86469, single-cell sequencing was utilized to determine cell types. Finally, the expression levels of biomarkers were validated through quantitative PCR (qPCR) and western blot analysis.

RESULTS

A total of 581 DEGs were identified in GSE76894. Four genes [ATP binding cassette subfamily C member 8 (ABCC8), retinol binding protein 4 (RBP4), Ras protein-specific guanine nucleotide-releasing factor 1 (RASGRF1), and solute carrier family 34 member 2 (SLC34A2)] were identified as ubiquitin-pyroptosis-related biomarkers in T2DM, based on consistent expression trends and significant differences in GSE76894 and GSE41762. These biomarkers were enriched in oxidative phosphorylation and mitogen-activated protein kinase signaling pathways, which are relevant to DM. ssGSEA revealed significant differences in the enrichment scores of nine immune cell types between groups. A total of 17 microRNAs (miRNAs) and 36 long non-coding RNAs (lncRNAs) were identified, forming numerous miRNA-lncRNA interactions. Additionally, 22 drugs related to the biomarkers, such as gliclazide and tretinoin, were identified. In GSE86469, eight cell types, including alpha and beta cells, were characterized. qPCR and western blot analysis confirmed that the expression trends of RASGRF1 and SLC34A2 were consistent with the findings in GSE76894.

CONCLUSION

This study identified four ubiquitin-pyroptosis-related biomarkers (ABCC8, RBP4, RASGRF1, and SLC34A2) in T2DM through bioinformatics analysis, providing novel insights into the diagnosis and treatment of T2DM.

Key Words: Type 2 diabetes mellitus; Ubiquitin-proteasome system; Pyroptosis; Biomarkers; Intercellular biocommunication; Regulatory networks

Core Tip: This study identified four ubiquitination and pyroptosis-associated biomarkers (ATP binding cassette subfamily C member 8, retinol binding protein 4, Ras protein-specific guanine nucleotide-releasing factor 1, and solute carrier family 34 member 2) in type 2 diabetes mellitus (T2DM) through integrated bioinformatics analysis, establishing a potential mechanistic link between ubiquitin-proteasome system dysregulation and pyroptosis in T2DM pathogenesis. These biomarkers were enriched in oxidative phosphorylation and mitogen-activated protein kinase signaling pathways, and modulated immune cell infiltration. Our research discovered possible drug candidates and a gene interaction network linked to RNA competition, which may aid in detecting and treating type 2 diabetes.

INTRODUCTION

Diabetes mellitus (DM) is a significant global health challenge characterized by chronic hyperglycemia, insulin resistance, and dysregulated glucose homeostasis[1]. The impact of DM extends widely, with profound consequences on individual health and healthcare systems. DM can be classified into various types, with type 2 DM (T2DM) being the most prevalent form. T2DM is characterized by insulin resistance and impaired pancreatic β-cell function, leading to disrupted glucose metabolism[2]. It is often associated with other metabolic syndromes, such as obesity and dyslipidemia, further complicating its clinical management[3].

Pyroptosis is a novelly discovered form of programmed cell death[4,5]. It represents an inflammatory response triggered by pathogens such as bacteria or their endotoxins. Activation of the NLR family pyrin domain containing 3 (NLRP3) inflammasome initiates a cascade of events: It activates caspase-1, which in turn cleaves gasdermin D, releasing its N-terminal domain. This action creates pores in the cell membrane, resulting in membrane rupture, cell osmolysis, DNA lysis, and the release of cellular contents and inflammatory mediators, including interleukin 1 beta (IL-1β) and IL-18, leading to a robust inflammatory response in various diseases[6-8]. Notably, increasing research links pyroptosis to T2DM. In one study, NLRP3^-/- mice exhibited improved glucose tolerance and insulin sensitivity. NLRP3 can promote islet fibrosis and increase insulin resistance in mice[9]. Gao et al[10] indicated that NLRP3-mediated inflammation is activated by high glucose and lipid toxicity, disrupting islet function and supporting mitochondria in releasing oxygen free radicals, exacerbating islet dysfunction.

The ubiquitin-proteasome system (UPS) is a tightly regulated cellular machinery responsible for protein degradation, controlling the abundance of key regulators involved in glucose metabolism, insulin signaling, and inflammatory responses[11,12]. Dysfunctional UPS has been linked to insulin resistance, β-cell dysfunction, and other hallmarks of T2DM[13-15]. Interestingly, ubiquitination also plays a pivotal role in modulating the inflammasome components through ubiquitin (Ub) post-translational modifications, which are crucial for regulating inflammasome activation[16]. The ubiquitination of components such as NLRP3, caspase-1/11, IL-1β, and other constituents of the inflammasome finely tunes several critical nodes within the regulatory network[17]. Notably, ubiquitination significantly governs the expression and activation of the NLRP3 inflammasome[18].

Given the intertwined nature of these cellular processes and their relevance to T2DM, this study aimed to bridge the knowledge gap by comprehensively screening for potential biomarkers associated with the UPS and pyroptosis in T2DM patients. We hypothesize that specific proteins within the UPS and key markers of pyroptosis are dysregulated in T2DM patients, potentially serving as diagnostic and therapeutic targets. Furthermore, we attempted to unravel the intricate regulatory mechanisms governing these processes, shedding light on the pathophysiology of T2DM.

MATERIALS AND METHODS

Data source

The datasets GSE76894, GSE41762, and GSE86469 were obtained from the gene expression omnibus database (http://www.ncbi.nlm.gov/geo/)[19]. GSE76894 contained 19 T2DM samples and 84 normal control (NC) samples, as a training set. GSE41762 included 20 T2DM samples and 57 NC samples, as a validation set. GSE86469 included three T2DM samples and five NC samples for single-cell RNA sequencing (scRNA-seq). Information related to these datasets is shown in Supplementary Table 1. A total of 797 UPS-related genes (UPSGs)[20] (Supplementary Table 2) and 33 pyroptosis-related genes (PRGs)[21] (Supplementary Table 3) were obtained from previous research.

Differential expression analysis and functional enrichment analysis

In the training cohort, differential gene expression analysis was conducted using the limma package in R to identify differentially expressed genes (DEGs) between individuals with T2DM and healthy controls. Following linear model implementation via the limma package, unadjusted P values were calculated for all genes. To control the false discovery rate in multiple hypothesis testing, the Benjamini-Hochberg method was employed to adjust raw P values. To address this, we calculated the adjusted P values and identified DEGs based on the criteria of |log2 fold change| ≥ 0.5 and adjusted P < 0.05. The results were visualized using volcano plots (ggplot2 package v3.3.6)[22] and heat map (“heatmap3” package (version 1.1.9)[23]. To understand the related biological functions and pathways of DEGs, Gene Ontology (GO), and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis were performed using the “clusterProfiler” package in R (adjusted P < 0.05)[24].

Screening for Ub-pyroptosis-related biomarkers

DEGs were intersected with UPSGs and PRGs to obtain differentially expressed UPSGs (DEUPSGs) and PRGs (DEPRGs). Candidate genes were selected based on significant associations with two or more DEUPSGs (|r| > 0.5, P < 0.05) and one or more DEPRGs (|r| > 0.5, P < 0.05) using Pearson’s correlation. Receiver operating characteristic (ROC) curve analysis (area under the curve [AUC] > 0.7) using the pROC package (version 1.18.0)[25] was employed to obtain key candidate genes in the training and validation sets. Moreover, the significance of the differential expression of key candidate genes between groups (T2DM and control) was analyzed by the t-test. Genes with consistent expression trends and significant differences between groups (P < 0.05) in the training set and validation set were used as biomarkers in this study. Afterward, in the training set, the correlation between biomarkers was calculated based on Pearson’s correlation analysis[26], and the box chart was used to show results through “ggplot2” in R package[27]. Finally, ROC analysis was used to verify the biomarkers.

Gene set enrichment analysis and ingenuity pathway analysis

In order to understand the differences in biomarker-related biological functions and signaling pathways involved in, we implemented Gene Set Enrichment Analysis (GSEA). First, the correlations between biomarkers and other genes were calculated and sorted in the training set. The “msigdbr” package of R language was used to download the C2: KEGG gene set of Homo sapiens as the background set, and GSEA function of R language[28] was used for enrichment analysis of the sequenced genes in the background gene set (adjusted P < 0.05). In turn, all of the differential genes were analyzed by ingenuity pathway analysis (IPA)[29] to obtain the related pathways involved, and the pathways involved in biomarkers were screened.

Immune infiltration analysis and Spearman’s correlation analysis

The Gene Set Variation Analysis (GSVA) package in R language[30] was used for immune infiltration analysis of all samples in the training and validation sets, with immune-related genes[31] as the background gene set to obtain enrichment scores for 28 immune cell types. The Wilcoxon test analyzed differences in immune cell enrichment scores between T2DM and control groups. Spearman’s correlation analysis evaluated the correlation between biomarkers and immune cell enrichment scores with intergroup differences[32].

Construction of competitive endogenous RNA and gene-drug networks

Biomarker-related microRNAs (miRNAs) were screened with clipExpNum ≥ 5, and miRNA-related long non-coding RNAs (lncRNAs) with clipExpNum ≥ 20 in Starbase (http://starbase.sysu.edu.cn/in). A competitive endogenous RNA (ceRNA) network was constructed based on the screened miRNAs and lncRNAs. Additionally, biomarkers-related drugs were identified using the DrugBank database (http://go.drugbank.com), and a gene-drug network was constructed.

Filtering and normalization of scRNA-seq data

scRNA-seq is a transformative technology that provides insights into cell-to-cell variation, high-resolution profiling of disease mechanisms, and potential clinical utility[33]. Cells with fewer than 3500 expressed genes were excluded according to quality control criteria, leaving 638 cells and 26616 genes. The Seurat package in R language[34] was further used to remove genes detected in fewer than three cells. After filtration, 638 cells and 17467 genes remained. Logarithmic standardization and the variance stabilizing transformation method were used to extract genes with large intercell variation coefficients, focusing on the top 2000 high-variable genes for subsequent analysis.

Dimension reduction, unsupervised clustering, and visualization

Cell-to-cell communication was measured by quantifying ligand-receptor pairs between different cell types. Understanding differences in cell-cell interactions is crucial for exploring potential therapeutic targets. The CellChat package in R language[35] evaluates cell-cell communication based on the CellChatDB database (www.cellchat.org), using cell expression data combined with ligand-receptor interactions to simulate cell communication.

Analysis of cell communication between healthy and T2DM samples

Cell-to-cell communication was quantified using ligand-receptor pairs between different cell types, evaluated by the CellChat package based on CellChatDB.

Animal husbandry and tissue collection

This study received approval from the Ethics Committee of The First Affiliated Hospital of USTC (Hefei, China) and the Medical Institution Animal Care and Research Advisory Committee, following Animal Research: Reporting of In Vivo Experiments guidelines (No. 2024-N(A)-80]) Ten 8-week-old C57BL/6 mice were housed under specific pathogen-free conditions. Of these, five mice were assigned to the NC group and fed a standard chow diet, while the remaining five were assigned to the T2DM group and fed a high-fat diet. After a 4-week acclimatization and dietary induction period, mice in the T2DM group received intraperitoneal injections of streptozotocin (25 mg·kg^-1) once daily for five consecutive days[36], whereas the NC group received equivalent volumes of citric acid-sodium citrate buffer and continued on the standard diet. Mice with fasting blood glucose levels ≥ 11.1 mmol/L were considered to have developed T2DM. All mice in the T2DM group successfully exhibited the diabetic phenotype, with an average fasting blood glucose level of 13.9 ± 3.2 mmol/L, compared to 7.5 ± 1.9 mmol/L in the NC group. Following model establishment, all animals were maintained on their respective diets for an additional 8 weeks. At the end of the 12-week experimental period, mice were fasted for 12 hours, anesthetized with 0.7% sodium pentobarbital, and blood and pancreatic tissues were collected post-mortem.

Quantitative PCR validation

Total RNA from pancreatic tissues was extracted with TRIzol reagent (Yeasen Biotech, Shanghai, China). cDNA synthesis was performed using the PrimeScript RT Master Mix Kit (Takara Bio, Tokyo, Japan). The qPCR assay used the SYBR^® Premix Ex Taq^™ II Kit (Takara Bio). The relative quantification of mRNAs was calculated using the 2^-ΔΔCt method. Primer sequences are listed in Supplementary Table 4.

Western blot validation

Pancreas samples were homogenized in RIPA buffer with protease and phosphatase inhibitors. After centrifugation at 14000 g for 15 minutes, the supernatant was collected for protein quantification using the BCA Protein Assay (Beyotime, Shanghai, China). Proteins were denatured in Laemmli buffer, separated by sodium dodecyl sulfate polyacrylamide gel electrophoresis, and transferred to polyvinylidene difluoride membranes. Membranes were blocked in 5% bovine serum albumin in Tris-buffered saline with Tween 20, and incubated with primary antibodies (1:1500, anti-ATP binding cassette subfamily C member 8 [ABCC8], Abcam, Waltham, MA, United States; 1:1000, anti-retinol binding protein 4 [RBP4], Abcam; 1:1500, anti-solute carrier family 34 member 2 [SLC34A2], Bioss Inc., Woburn, MA, United States; 1:1500, anti-Ras protein-specific guanine nucleotide-releasing factor 1 [RasGRF1], Bioss] and horseradish peroxidase-conjugated secondary antibodies. Protein bands were visualized using chemiluminescence (Thermo Fisher Scientific, Waltham, MA, United States) and quantified by image analysis (version 1.40; National Institutes of Health, Bethesda, MD, United States).

RESULTS

Identification and functional enrichment of DEGs

In the training set, a total of 581 DEGs were identified, including 157 upregulated genes and 424 downregulated genes (Figure 1A and B, Supplementary Table 5). Subsequently, functional enrichment analysis revealed that among the 551 GO terms, there were 65 genes categorized into cellular components, molecular functions, and biological processes. Additionally, 22 KEGG pathways were identified (Figure 1C and D; Supplementary Tables 6 and 7). Notably, some of these pathways were related to DM, such as hsa04911 (insulin secretion), hsa04950 (maturity-onset diabetes of the young [MODY]), hsa04930 (T2DM), and GO: 0005267 (potassium channel activity).

Figure 1 Differential expression analysis and functional enrichment analysis in the GSE76894 dataset. A: Volcano plot of differentially expressed genes (DEGs) between type 2 diabetes mellitus (T2DM) and normal control (NC) samples. The horizontal axis, log2 (fold change), represented the logarithm to the base 2 of the fold change in gene expression between T2DM samples and NC samples. The vertical axis, -log10 (P value), was used to measure the significance of the differences. The smaller the P value, the larger the -log10 (P value), and the more significant the differential gene expression was. Red represented genes with up-egulated expression, while green represented genes with downregulated expression; B: Heatmap of the distribution of DEGs between T2DM and NC samples. In the upper area, the horizontal axis represented different samples, the vertical axis represented gene expression levels, the waveform plot showed the distribution density of gene expression, and the black dashed lines corresponded, from bottom to top, to the gene expression levels of the 0%, 25%, 50% (mean), 75%, and 100% quantiles; C: G of Gene Ontology (GO) terms enriched in DEGs. The regions of different colors respectively represented GO entries related to biological process, cellular component, and molecular function. Each rectangular block represented a specific GO term, and generally, the size of its area was associated with the number of genes enriched in that term; D: Graph of Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways enriched by DEGs. Regions in different colors represented different KEGG pathways. Each rectangular block represented a specific KEGG pathway, and the size of its area was usually related to the number of genes enriched in that pathway.

Identification of biomarkers

A Venn diagram showed that a total of 13 DEPRGs (UBE2QL1, UCHL1, RAG1, KCNA6, KCNG3, FBXL16, TRIM9, KBTBD12, KCNB2, MEX3D, ASB9, KCNA5, and KCNC1) and one differentially expressed UPS-related gene (DEUPSG) (IL1B) were identified in the training set (Figure 2A; Supplementary Table 8). Based on the correlation of DEGs, DEUPSGs, and DEPRGs, a total of 15 candidate genes [ABCC8, RBP4, RASGRF1, guanylate-binding protein 1 (GBP1), GBP2, SLC34A2, lymphotoxin beta, endothelin 3 (EDN3), high mobility group nucleosomal binding domain 2 pseudogene 46 (HMGN2P46), SLC6A14, pyridine nucleotide-disulphide oxidoreductase domain 2 (PYROXD2), neuronal vesicle trafficking associated 2 (HMP19), arrestin domain containing 3 (ARRDC3), phospholipase C epsilon-1 (PLCE1), and synuclein alpha interacting protein (SNCAIP)] were identified (Supplementary Table 9). Correlation analysis revealed that ABCC8, RBP4, RASGRF1, and DEUPSGs had higher correlations. Additionally, GBP2 had the highest correlation coefficient with the PRGs IL1B (Figure 2B; Supplementary Tables 10 and 11). After ROC curve analysis of the 15 candidate genes and further screening for genes with AUC > 0.7, a total of 8 key candidate genes (ABCC8, RBP4, RASGRF1, GBP1, GBP2, SLC34A2, ARRDC3, and SNCAIP) were identified (Figure 2C). In both the training and validation sets, t-tests revealed that ABCC8, RBP4, RASGRF1, and SLC34A2 showed consistent expression trends. It is noteworthy that the expression levels of ABCC8, RBP4, and RASGRF1 were significantly lower in the T2DM group compared to the control group (P < 0.05), while the expression levels of SLC34A2 were significantly higher in the T2DM group compared to the control group (P < 0.05). Intergroup differences suggested that these biomarkers are closely related to the pathophysiological processes of T2DM and are capable of distinguishing T2DM patients from the normal population. Therefore, ABCC8, RBP4, RASGRF1, and SLC34A2 were recognized as biomarkers (Figure 2D and E; Supplementary Tables 12 and 13). Moreover, ABCC8, RBP4, and RASGRF1 exhibited AUC values exceeding 0.7 in both datasets, demonstrating their good diagnostic potential (Figure 2F and G). Pearson’s correlation analysis indicated the highest correlation between ABCC8 and RBP4 (r = 0.82) (Figure 2H).

Figure 2 Identification and validation of four biomarkers. A: Differentially expressed genes (DEGs) were intersected with ubiquitin-proteasome system related genes (UPSGs) and pyroptosis-related genes (PRGs) to recognize differentially expressed UPSGs (DEUPSGs) and differentially expressed PRGs (DEPRGs), respectively; B: Correlation analysis of DEGs to DEPRGs and DEUPSGs; C: Receiver operating characteristic (ROC) curves of candidate genes obtained by correlation analysis; D and E: Expression analysis of biomarkers in the type 2 diabetes mellitus and normal control samples in the GSE76894 and GSE41762 datasets; F and G: ROC curves of four biomarkers in the GSE76894 and GSE41762 datasets; H: Correlation analysis of biomarkers. ^aP < 0.05; ^bP < 0.01; ^cP < 0.001; ^dP < 0.0001. AUC: Area under the curve.

Relevant pathways and functions of biomarkers

In the training set, kegg_oxidative_phosphorylation (oxidative phosphorylation) and kegg_mapk_signaling_pathway (mitogen-activated protein kinase [MAPK] signaling pathway) in the TOP pathway were related to DM (Figure 3A-D; Supplementary Tables 14-17). Afterward, the related pathways of differential genes were obtained by IPA analysis and the pathways involved in biomarkers were screened. A total of seven pathways, among which the synaptogenesis signaling pathway with the highest Z-score was involved in RASGRF1 (Figure 3E; Supplementary Table 18). The regulatory relationship between them is shown in the Figure 3F (Supplementary Table 19).

Figure 3 Investigation of potential function for biomarkers. A-D: Gene set enrichment analysis of biomarkers (A: ATP binding cassette subfamily C member 8 [ABCC8]; B: Retinol binding protein 4 [RBP4]; C: Ras protein-specific guanine nucleotide-releasing factor 1 [RASGRF1]; D: Solute carrier family 34 member 2 [SLC34A2]); E: Results of ingenuity pathway analysis for biomarkers; F: Regulatory relationship of biomarkers and their enrichment pathways. ^aP < 0.05; ^bP < 0.01; ^cP < 0.001.

Assessing the association between biomarkers and immunogenicity

The Wilcoxon test was used to analyze differences in immune cell enrichment scores between T2DM and control groups, revealing significant differences in nine immune cell types (P < 0.05) in the training set, such as activated B cells, activated CD4 T cells, CD56 bright natural killer cells, and macrophages, all of which showed higher levels of infiltration in the disease group (Figure 4A, Supplementary Table 20). Similarly, in the validation set, activated B cells, central memory CD4 T cells, effector memory CD8 T cells, myeloid-derived suppressor cells, and neutrophils were notably higher in the disease group (Figure 4B, Supplementary Table 21). In both datasets, ABCC8, RBP4, and RASGRF1 exhibited a negative correlation with neutrophils, whereas SLC34A2 exhibited a positive correlation (Figure 4C and D, Supplementary Tables 22-25). These findings suggest that the identified biomarkers may influence the onset and progression of T2DM by acting on immune cells.

Figure 4 Immune infiltration analysis. A and B: Comparison of immune cell scores between type 2 diabetes mellitus (T2DM) and normal control groups in training set and validation set; C and D: Correlation of biomarkers to differential immune cells in training set and validation set. The color and size of the irregular shapes in the figure represent the size of the correlation. ^aP < 0.05; ^bP < 0.01; ^cP < 0.001; ^dP < 0.0001. MDSC: Myeloid-derived suppressor cell.

ceRNA network regulation and drug prediction of biomarkers

In Starbase, a total of 17 miRNAs related to the biomarkers were screened with clipExpNum ≥ 5. Predicted miRNAs for RBP4 and RASGRF1 included hsa-miR-5195-3p, hsa-miR-382-5p, hsa-miR-4306, hsa-miR-4644, and hsa-miR-506-3p. Additionally, 36 LncRNAs related to 17 miRNAs were screened with clipExpNum ≥ 20, including DANCR, SNHG29, GNAS-AS1, FGD5-AS1, and MALAT1 (Figure 5A; Supplementary Table 26). These formed numerous interactions, such as hsa-miR-361-3P-SNHG1, hsa-miR-577-OLMALING, and hsa-miR-873-5p-AC023509.1. Moreover, 22 drugs related to the biomarkers were identified, excluding RASGRF1 (Figure 5B; Supplementary Table 27). For instance, drugs associated with ABCC8 included ATP, glipizide, chlorpropamide, mitiglinide, and glimepiride. Drugs related to RBP4 included tretinoin, beta carotene, fenretinide, and vitamin A, whereas drugs related to SLC34A2 included calcium phosphate, calcium phosphate dihydrate, sodium phosphate dibasic, sodium phosphate monobasic, and unspecified forms of sodium phosphate monobasic.

Figure 5 Exploration of biomarker regulatory mechanisms and drug prediction. A: Construction of mRNA-microRNA (miRNA)-long non-coding RNA (lncRNA) network for biomarkers. The red graphics represent biomarkers, yellow represent miRNAs, and green represent lncRNAs; B: Regulatory network of biomarkers and targeted drugs. The red graphics represent biomarkers and green represent drugs.

Distribution of cell types in the T2DM and control groups

Based on t-distributed Stochastic Neighbor Embedding clustering, both the T2DM and control groups contained cells from different sources, totaling eight cell types: Acinar, alpha, beta, delta, ductal, gamma/PP, none/other, and stellate (Figure 6A). Visualization of biomarker distribution in cell groups showed that islet tissues mainly contained alpha, beta, delta, and bamma/PP cells. ABCC8 was expressed in all five cell types, RBP4 was mainly expressed in beta and delta cells, RASGRF1 in alpha and beta cells, and SLC34A2 predominantly in ductal cells (Figure 6B).

Figure 6 Expression levels of biomarkers in cells based on the single-cell RNA sequencing data. A: Results of t-distributed Stochastic Neighbor Embedding (t-SNE) clustering results of t-SNE clustering of cells from type 2 diabetes mellitus and normal control samples; B: Expression levels of four biomarkers in eight cell types.

Cell-cell interaction

Based on the CellChatDB database, it was found that in the control group, except for stellate and none/other cells, other cells communicated, with ductal cells having the most ligand-receptor pairs and the highest interaction intensity between alpha and beta cells. In the T2DM group, only alpha, beta, ductal, and none/other cells exchanged signals, with the highest interaction intensity between alpha and beta cells. Communication among delta, gamma/PP, and acinar cells with other cells was reduced in the T2DM compared to the control group (Figure 7; Supplementary Tables 28-31).

Figure 7 Results of cell-to-cell communication containing number of interactions and interaction weights/strength in type 2 diabetes mellitus and normal control samples. A: Normal control; B: Type 2 diabetes mellitus.

Functional comparison between cells

To further understand functional differences among cells, GSVA analysis based on the MSigDB database revealed that alpha, beta, delta, and gamma/PP cells had similar functions, while acinar, ductal, and stellate cells shared similar functions (Figure 8A; Supplementary Table 32). Subsequently, functional differences between T2DM and control groups were screened (P < 0.05), identifying 16 different functions: 9 upregulated (e.g., myogenesis, transforming growth factor beta signaling, and inflammatory response) and 7 downregulated (e.g., adipogenesis, mammalian target of rapamycin complex 1 signaling, and fatty acid metabolism) (Figure 8B; Supplementary Table 33).

Figure 8 Gene set variation analysis. A: Heatmap of pathways obtained by gene set variation analysis in eight cell types; B: Score of pathways enriched in type 2 diabetes mellitus and normal control groups.

Validation of biomarkers

Both western blot and quantitative PCR analyses demonstrated that SLC34A2 expression was significantly upregulated in T2DM samples compared to controls (P < 0.05), whereas RASGRF1 expression was significantly downregulated in T2DM samples (P < 0.01) (Figures 9 and 10). Experimental validation findings demonstrated concordance with bioinformatically derived hypotheses. Experimental validation identified divergent expression patterns in two biomarkers: RBP4 Levels were significantly elevated in T2DM cohorts (P < 0.05), whereas ABCC8 demonstrated no significant differential expression between the T2DM and control groups (P > 0.05). These findings contrasted with bioinformatics predictions. The observed divergence between experimental and bioinformatics findings for RBP4 and ABCC8 may reflect potential confounding factors including cohort sample size constraints, methodological heterogeneity, or underlying biological variability, highlighting the necessity of expanded validation studies in ethnically diverse populations with adequate statistical power. The reproducible experimental validation of SLC34A2 and RASGRF1 confirms their biological plausibility and establishes these molecules as credible candidate biomarkers for T2DM pathogenesis.

Figure 9 Expression level of biomarkers in the type 2 diabetes mellitus and normal control samples by quantitative PCR. ^aP < 0.05. NC: Normal control; T2DM: Type 2 diabetes mellitus.

Figure 10  Expression level of biomarkers in the type 2 diabetes mellitus and normal control samples by western blot analysis. ^aP < 0.05. T2DM: NC: Normal control; Type 2 diabetes mellitus.

DISCUSSION

T2DM, characterized by insulin resistance and inadequate insulin production, represents a significant global health challenge. Understanding the molecular mechanisms underlying the pathophysiology of T2DM is crucial for improving current therapies. Pyroptosis, a highly inflammatory form of programmed cell death characterized by rapid cell membrane rupture and the release of pro-inflammatory contents, has emerged as a significant player in the context of T2DM[37]. Recent studies have highlighted that pyroptosis can influence insulin sensitivity and secretion by modulating the release of inflammatory factors such as IL-1β and IL-18[38,39]. Furthermore, pyroptosis potentially impacts key molecules in insulin signaling pathways, such as insulin receptor substrate 1 [IRS-1] and IRS-2, consequently affecting normal insulin signaling and exacerbating insulin resistance[40]. The UPS is a prominent pathway responsible for the degradation of numerous proteins through polyubiquitination[41]. Research indicates that the UPS plays a pivotal role in insulin secretion by orchestrating critical proteins within the β-cell secretory cascade[42]. Moreover, it governs vital proteins such as IRS-2 and cAMP response element binding protein, which are essential for β-cell survival. Disruption of UPS function may lead to reduced insulin secretion levels in pancreatic β-cells[12,43]. Therefore, an in-depth exploration of pyroptosis and UPS in T2DM holds promise for identifying novel therapeutic targets and strategies for the treatment and management of T2DM.

Integrative analysis of DEGs, DEUPSGs, and DEPRGs has revealed coordinated associations of ABCC8, RBP4, and RASGRF1 with UPS activity, in addition to a functional link between GBP2 and the pyroptosis effector IL-1b. While existing evidence does not directly implicate RBP4 and RASGRF1 in UPS regulation, our findings propose previously uncharacterized mechanistic associations. ABCC8, a core subunit of ATP-sensitive potassium channels essential for insulin secretion, exhibits ubiquitination-suppressive activity in MODY[44], thereby stabilizing insulin signaling pathway components. Under T2DM-associated metabolic stress, UPS dysregulation may disrupt β-cell proteostasis through pathogenic accumulation of misfolded proteins. Zhang et al[45] demonstrated RBP4-driven NLRP3 inflammasome activation promoting pyroptosis in acute myocardial infarction. By contrast, T2DM-associated RBP4, previously linked to oxidative stress pathway enrichment[46], may potentiate UPS activation via reactive oxygen species (ROS)-mediated proteotoxic stress, thereby inducing pyroptotic cell death in pancreatic β-cells and insulin-responsive tissues. Wang et al[47] demonstrated a link between GBP2 and pyroptotic inflammation in diabetes, showing that it thereby facilitates T2DM progression through IL-1b signaling pathways. We hypothesize that ABCC8, RBP4, RASGRF1, and GBP2 play a role in diabetes through distinct pathways involving UPS and pyroptosis mechanisms. Additionally, we report novel T2DM correlations with EDN3, HMGN2P46, SLC6A14, PYROXD2, and HMP19.

ROC analysis and t-tests differentiated ABCC8, RBP4, RASGRF1, and SLC34A2 as biomarkers associated with the UPS and pyroptosis. PCR, western blotting, and bioinformatics analysis demonstrated significant elevation of SLC34A2 expression and differential expression profile of RASGRF1 in T2DM compared to non-T2DM individuals. Le et al[48] reported increased SLC34A2 Levels in the liver and colon of diet-induced obese mice, which likely findings of enhanced oxidative phosphorylation and ROS production. These changes may promote UPS degradation and drive NLRP3 pyroptosis via inflammatory processes. By contrast, RASGRF1 exhibits significant downregulation in T2DM and db/db islets[49], which has been shown to modulate the MAPK pathway[50]. This modulation likely results in the phosphorylation of proteins to promote UPS degradation. Reduction of RASGRF1 may impair UPS efficiency, leading to an increase in pro-inflammatory proteins that subsequently enhance NLRP3 pyroptosis, further contributing to insulin resistance.

Previous studies have demonstrated that oxidative stress suppresses insulin signaling through ROS-mediated mechanisms[47], while the MAPK pathway modulates T2DM development by altering insulin signaling pathways[50]. GSEA analysis associated these pathways with the UPS and pyroptosis pathways, with UPS clearing misfolded proteins and pyroptosis amplifying inflammation, both contributing to T2DM pathogenesis. Bioinformatic analysis via IPA indicates RASGRF1-mediated synaptogenesis signaling pathway dysregulation converges on pro-inflammatory axes in T2DM. ABCC8 modulates the UPS in insulin secretion, while RBP4 associates oxidative stress with both the UPS and pyroptosis pathways. RASGRF1 Links the MAPK pathway to protein degradation and inflammation, and SLC34A2 exacerbates metabolic stress, collectively providing a mechanistic framework for T2DM progression.

T2DM is a chronic low-grade inflammatory disease, and CD4 T cells play a crucial role in inducing activation of pancreatic islet β cells, producing various cytokines, thus promoting the occurrence of T2DM[51]. Moreover, neutrophils are identified to be the first-line defense against the invasion of pathogens in innate immunity. NETosis, the release of neutrophil extracellular traps (NETs) by neutrophils to capture pathogens and contribute to inflammation, is implicated in diabetes and its associated complications[52]. In this study, we found that RASGRF1 is associated with neutrophils and SLC34A2 is linked to activated CD4 T cells, suggesting their potential involvement in immune processes related to diabetes onset, identified through single sample GSEA analysis and the Wilcoxon test. Based on the PCR and western blot validation results, we hypothesize that SLC34A2 promotes the expression of CD4 T cells, which subsequently influences the function of pancreatic islet β cells and stimulates cytokines production, thereby contributing to the development of T2DM. Meanwhile, RASGRF1, by reducing neutrophil, inhibits NET release, thereby delaying the inflammatory response in diabetes.

According to data from the DrugBank database, drug sensitivity analysis suggests that specific drugs, such as ATP, glipizide, chlorpropamide, mitiglinide, and glimepiride, are closely associated with ABCC8. ABCC8 is primarily linked to insulin-secretagogue antidiabetic drugs, leading us to speculate that it might serve as a target in the process of reducing blood sugar with antidiabetic medications. Insulin secretagogues have the potential to enhance insulin secretion by upregulating ABCC8 expression, thereby lowering blood sugar levels[53]. In the case of RBP4, the relevant drugs include tretinoin, beta carotene, fenretinide, and vitamin A. Retinoids, derivatives of vitamin A (retinol), act as intricate regulators of adipogenesis by activating specific nuclear receptors, including the retinoic acid receptor and retinoid X receptor. The interaction between circulating RBP4 and its membrane receptor, retinoic acid 6, coordinates cellular retinol uptake, ultimately leading to the activation of the c-Jun N-terminal kinase pathways and the Janus kinase 2/signal transducer and activator of transcription 5 cascade. These events contribute to insulin resistance and the development of diabetes[54]. Previous studies have extensively explored the use of retinoic acid, the active metabolite of vitamin A, for the treatment of diabetic complications[55]. Based on the findings of this study, we hypothesize that RBP4 could be a target for vitamin A and its metabolites in the treatment of diabetes. Lastly, for SLC34A2, the associated drugs encompass calcium Phosphate, calcium phosphate dihydrate, sodium phosphate dibasic, sodium phosphate monobasic, and others. As a sodium-phosphate co-transporter, SLC34A2 is closely related to phosphate levels. Variations in SLC34A2 can lead to impaired phosphate transport in patients with pulmonary alveolar microlithiasis[56]. Modulating SLC34A2 to mitigate ROS-driven inflammation or enhancing the MAPK signaling axis of RASGRF1 could preserve β-cell functionality. Additionally, immune-targeted therapies targeting SLC34A2/CD4 T cell activation or RASGRF1/NETosis could potentially mitigate inflammation and insulin resistance, providing a potential avenue for personalized management of T2DM.

MiRNAs predicted to target RBP4 and RASGRF1, such as hsa-miR-5195-3p and hsa-miR-382-5p among others, were identified in ceRNA. Additionally, a group of 36 LncRNAs associated with 17 miRNAs was screened, which includes notable examples like FGD5-AS1 and SNHG1. Hsa-miR-5195-3p and hsa-miR-873-5p are miRNAs that are known to influence insulin secretion and sensitivity, two critical factors in the development and progression of diabetes. They exert their effects by modulating the expression of genes involved in glucose metabolism and insulin signaling pathways[53-55]. FGD5-AS1 has emerged as a key regulator in diabetes. It can affect the availability and function of miRNAs like hsa-miR-361-3p and hsa-miR-873-5p, thereby indirectly influencing insulin secretion and action[56]. SNHG1 has been implicated in the regulation of cellular processes central to diabetes pathophysiology, such as pancreatic beta-cell function and insulin resistance. It can interact with various miRNAs and proteins, modulating their activity and thus impacting glucose homeostasis[57]. We speculate that FGD5-AS1 regulates RBP4 expression by acting as a competing endogenous RNA, thereby sponging hsa-miR-873-5p, which in turn potentially influences the development of T2DM. The SNHG1/hsa-miR-5195-3p/RASGRF1 axis is believed to function in a similar manner.

The complex interactions among these miRNAs and lncRNAs in diabetes suggest potential therapeutic targets. Understanding this network could lead to the development of novel strategies for the treatment and management of diabetes, potentially improving patient outcomes.

scRNA-seq revealed that RASGRF1 is predominantly expressed in α and β cells, while SLC34A2 is enriched in ductal cells. Cell-cell communication in the T2DM group was significantly reduced compared to the control group. Recent studies underscore the importance of communication between pancreatic endocrine and exocrine cells in diabetes pathogenesis[58]. Islet cell cooperation is highly synergistic, relying on intricate intercellular signaling to monitor neighboring cell activity[59]. For instance, alpha and delta cells coordinate insulin secretion through cell-to-cell communication[60,61]. Reduction in T2DM likely impairs islet cell cooperation, impairing insulin secretion and glucose regulation[62]. RASGRF1 downregulation in alpha and beta cells impairs Ras-mediated signaling, which is essential for β-cell responsiveness and α-cell glucagon regulation, thereby diminishing cross-talk essential for islet function. Similarly, elevated SLC34A2 in ductal cells may disrupt exocrine-endocrine axis communication, leading to increased local inflammation and diminishing cross-talk essential for islet function. Recent studies have demonstrated that impaired β-cell communication in T2DM fosters insulin resistance through inflammatory mechanisms and modulates islet inflammation and dysfunction[63]. We propose that these biomarkers represent key molecular indices across key cell types, contributing to the progression of T2DM by disrupting intercellular bio-communication pathways, thereby fostering β-cell dysfunction and insulin resistance.

CONCLUSION

In conclusion, this study identified biomarkers related to the UPS and pyroptosis, which play important roles in T2DM. Through our analysis, we explored the functional pathways, immune characteristics, regulatory networks, and potential drug interactions, gaining valuable insights into the mechanisms of T2DM and identifying potential therapeutic strategies. However, the study has certain limitations. Among the limitations, the small sample size and dependence on bioinformatics analysis may limit the robustness and generalizability of our findings. Despite the experimental validation of biomarker expression, the specific functions and mechanisms of these biomarkers require further investigation through additional in vitro and in vivo studies. Despite these limitations, clinical validation processes for the biomarkers and their potential drug interactions are essential, and further investigation of these biomarkers will be conducted in future research.