# Identification and Verification of Potential Biomarkers in Renal Ischemia-Reperfusion Injury by Integrated Bioinformatic Analysis

Authors: Ziwen Pan; Yuanyuan Yang; Rui Cao; Yang Qiu; Shanglin Li; Yuanyuan Zhao; Sheng Chang; Song Chen; Zhishui Chen; Weijie Zhang; Daqiang Zhao
Journal: BioMed Research International
 (2023)

---

## Abstract

Background. Renal ischemia-reperfusion injury (RIRI) plays an important role in the poor prognosis of patients with renal transplants. However, the potential targets and mechanism of IRI are still unclear. Method. Differential gene expression (DEG) analysis and weighted correlation network analysis (WGCNA) were performed on the GSE27274 dataset. Pathway enrichment analysis on the DEGs was performed. To identify the hub DEGs, we constructed a protein-protein interaction (PPI) network. Finally, the hub genes were verified, and candidate drugs were screened from the DsigDB database. Results. A hundred DEGs and four hub genes (Atf3, Psmb6, Psmb8, and Psmb10) were screened out. Pathway enrichment analysis revealed that 100 DEGs were mainly enriched in apoptosis and the TNF signaling pathway. The four hub genes were verified in animal models and another dataset (GSE148420). Thereafter, a PPI network was used to identify the four hub genes (Atf3, Psmb6, Psmb8, and Psmb10). Finally, eight candidate drugs were identified as potential drugs. Conclusion. Three hub genes (Psmb6, Psmb8, and Psmb10) were associated with RIRI and could be potential novel biomarkers for RIRI.


---

## Body

## 1. Introduction

Several individuals suffer from end-stage kidney disease every year. According to Lv and Zhang [1], 2.618 million individuals underwent renal replacement therapy in 2010. In this population of individuals, 2.05 million individuals underwent dialysis, whereas the remaining patients underwent kidney transplantation. Kidney transplantation is the preferred treatment for patients with end-stage renal disease, providing a better prognosis than dialysis. During kidney transplantation, graft injury begins at the moment the graft is separated from the donor. After a brief warm ischemic period, the grafts are placed in a hypothermal preserving solution for a long cold ischemia time (CIT) [2]. Renal ischemia-reperfusion injury (RIRI) is a process during which tissues and organs undergo an initial interruption of blood flow and subsequent restoration of blood flow. In this process, tissues and organs are subjected to oxidative stress and inflammatory reaction. Early IRI induces later loss by chronic hypoxia, reduced kidney mass, graft vascular injury, and subsequent fibrosis [3]. IRI is an important factor affecting the early functional recovery of a transplanted kidney, usually manifested as acute tubular necrosis and delayed graft function (DGF) [4]. DGF is defined as the requirement for dialysis within 7 days of transplant, and it is distinguished from primary nonfunction (PNF) by the eventual recovery of renal function [5]. Prolonged CIT is an independent risk factor for DGF [6]. DGF is the most common complication after kidney transplantation and a risk factor affecting the survival rate of transplanted kidneys and patients. Clinically, IRI has been a challenge for several years, and it still remains one today. IRI is a multifactorial and intricate process involving the activation of cell death programs, endothelial dysfunction, transcriptional reprogramming, and activation of the innate and adaptive immune systems [7]. Therefore, exploration of the key molecular mechanisms and targets is needed to address the RIRI challenge.With the help of bioinformatics, we gain new insights into the development of disease. Weighted gene coexpression network analysis (WGCNA) is a system biology method for describing the correlation patterns among genes across microarray samples [8]. WGCNA can identify modules of highly correlated genes, and these methods have been widely used in the fields of bioinformatics to identify candidate biomarkers. Gong et al. [9] performed WGCNA and constructed a protein-protein interaction (PPI) network to identify 10 DEGs in human anaplastic thyroid cancer. Zeng et al. [10] performed WGCNA and combined experimental verification to identify two hub genes in nonalcoholic fatty liver disease.Some bioinformatic studies on RIRI have been reported. Guo et al. [11] found 10 hub genes and Hif-1α signaling pathways associated with RIRI. Additionally, Zhu et al. [12] found three genes and two miRNAs that may be potential targets for RIRI. By WGCNA, Lin et al. [13] found that Rplp1 and Lgals1 were associated with the development of acute kidney injury. Building on previous studies, our study further verified the accuracy of our results and found potential drugs to treat RIRI.Above all, we focused on DEGs at both 6 hours and 24 hours after reperfusion and performed WGCNA to find genes linked with the development of IRI. Thereafter, we performed enrichment analysis for DEGs to find potential mechanisms. Finally, we performed a variation of hub genes and found potential drugs. This study sought to offer a novel insight into the pathogenesis of IRI and found potential biomarkers and drugs for treating IRI.


## 2. Materials and Methods

### 2.1. Data Information

Gene expression profiles of GSE27274 (GPL6101) and GSE148420 (GPL14746) were obtained from the Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/geo/) [14]. GSE27274 contains six sham group samples and 18 kidney IRI samples (six samples at 6, 24, and 120 hours after reperfusion). However, in our analysis, we used six samples at 6 and 24 hours after reperfusion. GSE148420 contains four normal samples and eight kidney IRI samples.


### 2.2. Screening Differentially Expressed Genes

All data were log-transformed, and expression matrices were normalized. The R package of “Limma” was used to identify DEGs between IRI samples and sham group samples. DEGs withlog2FoldChange>1 and falsediscoveryrateFDR<0.05 were selected.


### 2.3. Function Enrichment Analysis of DEGs

To further investigate the molecular mechanisms of potential differentially expressed IRGs, we used the David database (https://david.ncifcrf.gov) [15] to conduct Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analyses. GO terms and KEGG pathway with P value < 0.05 were considered significantly enriched.


### 2.4. Identification of IRI-Associated Modules by WGCNA and Matching Common Genes

In this study, we constructed a coexpression network for normalized gene expression data of the GSE27274 dataset using the “WGCNA” package. The soft threshold was set to 20, whereas scale-freeR2 was 0.8. The minimum number of genes in modules was set as 30 genes. The R package of “VennDiagram” was used to identify the intersection of the differential genes obtained from the GEO database and significant module genes.


### 2.5. PPI Network Analysis for the Identification of Hub Genes and ROC Curve Analysis in Another Dataset

To further identify hub genes, we established a PPI network on the Search Tool for the Retrieval on Interacting Genes (STRING) (https://string-db.org/) [16], which was used to predict associations among proteins. Additionally, we set the minimum required interaction score at 0.9. Cytoscape plug-in, CytoHubba, was used to identify hub genes. The diagnostic accuracy of hub genes was tested in another dataset.


### 2.6. Animal Model

Eighteen male Wistar rats weighing 200–300 g were obtained from Beijing Vital River Laboratory Animal Technology Co., Ltd (Beijing, China). Wistar rats were subjected to ischemia/reperfusion, as described in our previous study [17]. Rats were allocated to three groups (n=18) as follows: (1) the six hours after reperfusion group (n=6), (2) the 24 hours after reperfusion group (n=6), and (3) the sham group (n=6). Rats in the 6 and 24 hours after reperfusion groups were anesthetized. Thereafter, their right kidneys were removed, and their left renal pedicles were isolated from connective tissues and clamped for 48 minutes using a vascular clamp. For those in the sham group, they underwent the same procedure but without vessel clipping. The rats were sacrificed after 24 hours, and blood and kidneys were collected for further analysis.


### 2.7. Creatinine Analysis and Hematoxylin and Eosin (H&E) Staining

An automated biochemical analyzer BS-200 (Mindray, Shenzhen, China) was used to measure the creatinine levels of the collected blood samples. The kidney specimens were embedded in paraffin, sectioned, and stained with H&E. Tubular damage was scored, as described previously [18].


### 2.8. Western Blot Analysis

Total proteins from frozen tissues were extracted using lysis buffer (P0013J, Biyuntian, China). The protein concentration of samples was measured (BCA Protein Assay Kit (Biyuntian, China)), and equal amounts of protein were loaded for electrophoresis and subsequent transfer to polyvinylidene fluoride membranes (MilliporeSigma, Burlington, MA, USA). Using 5% BSA for 1 hour, membranes were blocked, and the primary antibody was incubated overnight at 4°C. On the second day, a secondary antibody was incubated at room temperature for 1 hour. Thereafter, proteins were visualized by ECL. A list of the primary antibodies has been provided in Supplementary Table1.


### 2.9. Real-Time PCR

Total RNA from frozen tissues was extracted using RNA-Easy™ Isolation Reagent Vazyme Cat (Vazyme, Nanjing, China). cDNA was synthesized using a reverse transcription kit (Vazyme). The mRNA levels of these genes were assessed using the StepOne software (ThermoFisher Scientific, China).β-Actin was used as a reference gene. Primer sequences are listed in Supplementary Table 2.


### 2.10. Prediction of Candidate Drugs for Key Genes

Approved drugs associated with RIRI were obtained from the DsigDB database [19], which was obtained from Enrichr (https://maayanlab.cloud/Enrichr/). Candidate drugs with adjusted P<0.01 were screened.


### 2.11. Statistical Analysis

DEGs were analyzed using the R package of “Limma” and were defined withlog2FoldChange>1 and falsediscoveryrateFDR<0.05. Differences between the two groups were identified using the t-test. P<0.05 was considered statistically significant. Statistical analysis was performed using GraphPad Prism 9 software (GraphPad Software, San Diego, CA, USA). Figure 1 depicts the workflow of the study.Figure 1
The workflow of the study.


## 2.1. Data Information

Gene expression profiles of GSE27274 (GPL6101) and GSE148420 (GPL14746) were obtained from the Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/geo/) [14]. GSE27274 contains six sham group samples and 18 kidney IRI samples (six samples at 6, 24, and 120 hours after reperfusion). However, in our analysis, we used six samples at 6 and 24 hours after reperfusion. GSE148420 contains four normal samples and eight kidney IRI samples.


## 2.2. Screening Differentially Expressed Genes

All data were log-transformed, and expression matrices were normalized. The R package of “Limma” was used to identify DEGs between IRI samples and sham group samples. DEGs withlog2FoldChange>1 and falsediscoveryrateFDR<0.05 were selected.


## 2.3. Function Enrichment Analysis of DEGs

To further investigate the molecular mechanisms of potential differentially expressed IRGs, we used the David database (https://david.ncifcrf.gov) [15] to conduct Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analyses. GO terms and KEGG pathway with P value < 0.05 were considered significantly enriched.


## 2.4. Identification of IRI-Associated Modules by WGCNA and Matching Common Genes

In this study, we constructed a coexpression network for normalized gene expression data of the GSE27274 dataset using the “WGCNA” package. The soft threshold was set to 20, whereas scale-freeR2 was 0.8. The minimum number of genes in modules was set as 30 genes. The R package of “VennDiagram” was used to identify the intersection of the differential genes obtained from the GEO database and significant module genes.


## 2.5. PPI Network Analysis for the Identification of Hub Genes and ROC Curve Analysis in Another Dataset

To further identify hub genes, we established a PPI network on the Search Tool for the Retrieval on Interacting Genes (STRING) (https://string-db.org/) [16], which was used to predict associations among proteins. Additionally, we set the minimum required interaction score at 0.9. Cytoscape plug-in, CytoHubba, was used to identify hub genes. The diagnostic accuracy of hub genes was tested in another dataset.


## 2.6. Animal Model

Eighteen male Wistar rats weighing 200–300 g were obtained from Beijing Vital River Laboratory Animal Technology Co., Ltd (Beijing, China). Wistar rats were subjected to ischemia/reperfusion, as described in our previous study [17]. Rats were allocated to three groups (n=18) as follows: (1) the six hours after reperfusion group (n=6), (2) the 24 hours after reperfusion group (n=6), and (3) the sham group (n=6). Rats in the 6 and 24 hours after reperfusion groups were anesthetized. Thereafter, their right kidneys were removed, and their left renal pedicles were isolated from connective tissues and clamped for 48 minutes using a vascular clamp. For those in the sham group, they underwent the same procedure but without vessel clipping. The rats were sacrificed after 24 hours, and blood and kidneys were collected for further analysis.


## 2.7. Creatinine Analysis and Hematoxylin and Eosin (H&E) Staining

An automated biochemical analyzer BS-200 (Mindray, Shenzhen, China) was used to measure the creatinine levels of the collected blood samples. The kidney specimens were embedded in paraffin, sectioned, and stained with H&E. Tubular damage was scored, as described previously [18].


## 2.8. Western Blot Analysis

Total proteins from frozen tissues were extracted using lysis buffer (P0013J, Biyuntian, China). The protein concentration of samples was measured (BCA Protein Assay Kit (Biyuntian, China)), and equal amounts of protein were loaded for electrophoresis and subsequent transfer to polyvinylidene fluoride membranes (MilliporeSigma, Burlington, MA, USA). Using 5% BSA for 1 hour, membranes were blocked, and the primary antibody was incubated overnight at 4°C. On the second day, a secondary antibody was incubated at room temperature for 1 hour. Thereafter, proteins were visualized by ECL. A list of the primary antibodies has been provided in Supplementary Table1.


## 2.9. Real-Time PCR

Total RNA from frozen tissues was extracted using RNA-Easy™ Isolation Reagent Vazyme Cat (Vazyme, Nanjing, China). cDNA was synthesized using a reverse transcription kit (Vazyme). The mRNA levels of these genes were assessed using the StepOne software (ThermoFisher Scientific, China).β-Actin was used as a reference gene. Primer sequences are listed in Supplementary Table 2.


## 2.10. Prediction of Candidate Drugs for Key Genes

Approved drugs associated with RIRI were obtained from the DsigDB database [19], which was obtained from Enrichr (https://maayanlab.cloud/Enrichr/). Candidate drugs with adjusted P<0.01 were screened.


## 2.11. Statistical Analysis

DEGs were analyzed using the R package of “Limma” and were defined withlog2FoldChange>1 and falsediscoveryrateFDR<0.05. Differences between the two groups were identified using the t-test. P<0.05 was considered statistically significant. Statistical analysis was performed using GraphPad Prism 9 software (GraphPad Software, San Diego, CA, USA). Figure 1 depicts the workflow of the study.Figure 1
The workflow of the study.


## 3. Results

### 3.1. Identification of Renal IRI-Related DEGs and Associated Modules by WGCNA

Differences in gene expression between 6-hour and 24-hour reperfusion samples and normal samples were analyzed (Figures2(a)–2(d)). The WGCNA package was used to construct coexpression networks (Figures 2(e)–2(h)) for the GSE27274 dataset. A soft threshold β of 14 (Figure 2(f)) was highly suitable for constructed gene modules. Among all modules, the turquoise modules (Figure 2(h)) with 1311 genes were the most relevant for the development of renal IRI. To further identify the DEGs associated with the development of renal IRI, we drew Venn diagrams to identify intersecting genes (Figure 2(i), Table 1).Figure 2
Analysis of gene expression difference and WGCNA of samples: (a, b) and (c, d) are volcano plots and heat maps for the analysis of differences between 6-hour and 24-hour reperfusion samples and sham samples, respectively; (e) sample cluster map; (f) soft thresholds for scale independence and mean connectivity; (g) the cluster dendrogram of genes; (h) module-trait relationship plotters. Each cell has a correlation coefficient andP value. (i) Venn plot of DEGs and turquoise module genes between 6-hour and 24-hour reperfusion samples and normal samples.
(a)(b)(c)(d)(e)(f)(g)(h)(i)Table 1
Overlap of DEGs originating from 6 hours and 24 hours of reperfusion with genes from turquoise modules.
A3GALT2CLDN5GSNNR1I2SLC25A11ACOX1CLIC1HNF4ANT5ESLC34A2ADAMTS1CLUHPNPCOLCESNRPAADHFE1CRY1ICAM1PEX11ASORDALDH9A1CTSDIER2PGPEP1STX8ANXA1CTSSIFITM3PLSCR1SULT1A1ANXA2CXCL10IRF1PMM1TAX1BP3ANXA5CXCL11JUNBPSMB10TIMP1ARPC1BDMGDHKIF22PSMB6TINAGASLEIF4EBP1KLF6PSMB8TMBIM1ATF3EMP3LAMC2PVRTMEM37ATP6V0E2EPN1LCN2RASD1TNFRSF12ABTG2FBLLGALS3RND1TNFRSF14C1QAFBLIM1LRPAP1RPP21TPBGC1QBFBXW5MAP3K1RPS9TPRKBCATFGBMETRNLRUNX1TUBB6CD14GADD45GMFGE8SFXN1UPP1CD68GJA4MYO7ASH2D4AVARS2CHPT1GLYCAM1NFS1SLC17A5WFDC2CHRDGPRC5CNQO1SLC1A5WSB1


### 3.2. Functional Enrichment Analysis of DEGs

GO (Figure3(a), Table 2) and KEGG (Figure 3(b), Table 3) enrichment analyses of these DEGs showed aging, cell surface, and heparin binding as the most common biological terms for biological process, cellular component, and molecular function. Moreover, pertussis was the most enriched term in the KEGG pathway analysis.Figure 3
Functional enrichment analysis of DEGs: (a) GO enrichment analysis of DEGs; (b) KEGG enrichment analysis of DEGs. BP: biological process; CC: cellular component; MF: molecular function.
(a)(b)Table 2
GO enrichment analysis.
CategoryDescriptionGO IDP valueBPAgingGO:00075681.50E−05BPPositive regulation of apoptotic processGO:00430650.001308075BPResponse to ethanolGO:00454710.00140547BPResponse to drugGO:00424930.001686102BPResponse to xenobiotic stimulusGO:00094100.002022073BPCellular response to lipopolysaccharideGO:00712220.004595786BPResponse to thyroid hormoneGO:00970660.004957165BPCellular response to nutrient levelsGO:00316690.004957165BPInnate immune responseGO:00450870.004957181BPRegulation of cell proliferationGO:00421270.005093297CCCell surfaceGO:00099861.56E−05CCExtracellular spaceGO:00056154.72E−05CCLysosomal membraneGO:00057658.70E−05CCBasement membraneGO:00056042.91E−04CCExtracellular exosomeGO:00700626.57E−04CCCytoplasmGO:00057370.001345772CCVesicleGO:00319820.001931725CCProteasome core complex, beta-subunit complexGO:00197740.00195821CCExternal side of plasma membraneGO:00098970.003278416CCPerinuclear region of cytoplasmGO:00484710.005507469MFHeparin bindingGO:00082012.32E−04MFThreonine-type endopeptidase activityGO:00042980.002408191MFTranscriptional activator activity, RNA polymerase II transcription regulatory region sequence-specific bindingGO:00012280.017925491MFEndopeptidase activityGO:00041750.019055149MFCXCR3 chemokine receptor bindingGO:00482480.026128457MFPhospholipase A2 inhibitor activityGO:00198340.026128457MFL-serine transmembrane transporter activityGO:00151940.031272543MFCalcium-dependent phospholipid bindingGO:00055440.033045687MFPeptidase activityGO:00082330.034029365MFTranscription factor activity, sequence-specific DNA bindingGO:00037000.035206985Table 3
KEGG enrichment analysis.
CategoryTermP valueKEGG_pathwayPertussis1.50E−02KEGG_pathwayAlcoholic liver disease1.60E−02KEGG_pathwayComplement and coagulation cascades2.30E−02KEGG_pathwayProteasome4.40E−02KEGG_pathwayTNF signaling pathway4.50E−02


### 3.3. Constructed PPI Network for Hub Gene Identification

To construct a PPI network, we uploaded 100 DEGs to the STRING database (version: 11.5;https://cn.string-db.org/) (Figure 4(a)). Thereafter, we analyzed the top six hub genes using the Cytoscape plug-in. Based on five algorithms of the Cytohubba plug-in (Figure 4(b)), four genes (Psmb6, Psmb8, Psmb10, and Atf3) were identified as key genes of RIRI. These genes were selected for further validation in the RIRI model.Figure 4
Building a PPI network of the DEGs and screening out hub genes: (a) PPI network exported from the STRING database; (b) Venn diagram for identification of hub genes by combining 5 algorithms of Cytohubba plug-in.
(a)(b)


### 3.4. Constructed Rat RIRI Model

As shown in Figure5(a), H&E staining of rat kidneys helped in finding the RIRI group (6 and 24 hours after reperfusion groups), which had significant proximal tubular injury. Concurrently, serum creatinine and BUN levels were significantly higher in the RIRI group than in the sham group (Figure 5(b)). Additionally, the expression of kidney injury molecule-1 (Kim-1) was positively correlated with the severity of RIRI. Kim-1 was more highly expressed in the RIRI group than in the sham group (Figures 5(c) and 5(d)). Altogether, these results indicated that the RIRI model was successfully established.Figure 5
Results of animal model: (a) hematoxylin and eosin (H&E) staining (200x, left) and tubular injury scores (right). The arrow indicates that the injury was serious. (b) The levels of serum creatinine and blood urea nitrogen (BUN) in different groups. (c, d) Western blot analysis of Kim-1 expression in different groups (∗∗∗∗P<0.0001;∗∗∗P<0.001;∗∗P<0.01).
(a)(b)(c)(d)


### 3.5. External Validation and Verification of Hub Genes

We performed ROC curve analysis to predict the values of four hub genes. The AUC value forPsmb6, Psmb8, Psmb10, and Atf3 was 0.719, 1, 1, and 0.906, respectively, in the GSE148420 dataset (Figures 6(a)–6(d)). Consistent with our analysis of GSE27274, these four hub genes were highly expressed in the RIRI group (6 and 24 hours after reperfusion groups). However, due to lack of statistical significance, they were placed in Figures S1(a)–1(h). Additionally, we performed western blot analysis and found that Psmb6, Psmb8, and Psmb10 were highly expressed in the RIRI group (Figures 6(f)–6(k)), whereas Atf3 was highly expressed in the 6 hours after reperfusion group (Figures 6(e) and 6(l)).Figure 6
ROC curve to verify gene expression in the GSE148420 dataset and western blot in the animal model: (a–d) ROC curve of four genes in the GSE148420 dataset; (e–h) the expression of four proteins in different groups (∗P<0.05;∗∗P<0.01;∗∗∗∗P<0.0001).


### 3.6. Identification of Candidate Drugs Associated with Renal IRI

Eight candidate drugs (Table4) with adjusted P value < 0.01 were identified in the DSigDB database. Psmb6 and Psmb8 were the most common drug targets, indicating these two genes might play a key role in the development of RIRI.Table 4
Prediction of candidate drugs for hub genes.
TermAdjustedP valueGenesVelcade (bortezomib) BOSS1.63E−06PSMB6; ATF3; PSMB8; PSMB10Carfilzomib CTD 000047871.63E−06PSMB6; PSMB8; PSMB10Carfilzomib BOSS1.63E−06PSMB6; PSMB8; PSMB10Velcade (bortezomib) TTD 000117782.16E−06PSMB6; PSMB8; PSMB10Carfilzomib6.66E−04PSMB6; PSMB8Bortezomib0.001156865PSMB6; PSMB8Clopamide HL60 down0.001156865PSMB6; PSMB8; PSMB10Azacyclonol HL60 down0.008485756PSMB6; PSMB8; PSMB10


## 3.1. Identification of Renal IRI-Related DEGs and Associated Modules by WGCNA

Differences in gene expression between 6-hour and 24-hour reperfusion samples and normal samples were analyzed (Figures2(a)–2(d)). The WGCNA package was used to construct coexpression networks (Figures 2(e)–2(h)) for the GSE27274 dataset. A soft threshold β of 14 (Figure 2(f)) was highly suitable for constructed gene modules. Among all modules, the turquoise modules (Figure 2(h)) with 1311 genes were the most relevant for the development of renal IRI. To further identify the DEGs associated with the development of renal IRI, we drew Venn diagrams to identify intersecting genes (Figure 2(i), Table 1).Figure 2
Analysis of gene expression difference and WGCNA of samples: (a, b) and (c, d) are volcano plots and heat maps for the analysis of differences between 6-hour and 24-hour reperfusion samples and sham samples, respectively; (e) sample cluster map; (f) soft thresholds for scale independence and mean connectivity; (g) the cluster dendrogram of genes; (h) module-trait relationship plotters. Each cell has a correlation coefficient andP value. (i) Venn plot of DEGs and turquoise module genes between 6-hour and 24-hour reperfusion samples and normal samples.
(a)(b)(c)(d)(e)(f)(g)(h)(i)Table 1
Overlap of DEGs originating from 6 hours and 24 hours of reperfusion with genes from turquoise modules.
A3GALT2CLDN5GSNNR1I2SLC25A11ACOX1CLIC1HNF4ANT5ESLC34A2ADAMTS1CLUHPNPCOLCESNRPAADHFE1CRY1ICAM1PEX11ASORDALDH9A1CTSDIER2PGPEP1STX8ANXA1CTSSIFITM3PLSCR1SULT1A1ANXA2CXCL10IRF1PMM1TAX1BP3ANXA5CXCL11JUNBPSMB10TIMP1ARPC1BDMGDHKIF22PSMB6TINAGASLEIF4EBP1KLF6PSMB8TMBIM1ATF3EMP3LAMC2PVRTMEM37ATP6V0E2EPN1LCN2RASD1TNFRSF12ABTG2FBLLGALS3RND1TNFRSF14C1QAFBLIM1LRPAP1RPP21TPBGC1QBFBXW5MAP3K1RPS9TPRKBCATFGBMETRNLRUNX1TUBB6CD14GADD45GMFGE8SFXN1UPP1CD68GJA4MYO7ASH2D4AVARS2CHPT1GLYCAM1NFS1SLC17A5WFDC2CHRDGPRC5CNQO1SLC1A5WSB1


## 3.2. Functional Enrichment Analysis of DEGs

GO (Figure3(a), Table 2) and KEGG (Figure 3(b), Table 3) enrichment analyses of these DEGs showed aging, cell surface, and heparin binding as the most common biological terms for biological process, cellular component, and molecular function. Moreover, pertussis was the most enriched term in the KEGG pathway analysis.Figure 3
Functional enrichment analysis of DEGs: (a) GO enrichment analysis of DEGs; (b) KEGG enrichment analysis of DEGs. BP: biological process; CC: cellular component; MF: molecular function.
(a)(b)Table 2
GO enrichment analysis.
CategoryDescriptionGO IDP valueBPAgingGO:00075681.50E−05BPPositive regulation of apoptotic processGO:00430650.001308075BPResponse to ethanolGO:00454710.00140547BPResponse to drugGO:00424930.001686102BPResponse to xenobiotic stimulusGO:00094100.002022073BPCellular response to lipopolysaccharideGO:00712220.004595786BPResponse to thyroid hormoneGO:00970660.004957165BPCellular response to nutrient levelsGO:00316690.004957165BPInnate immune responseGO:00450870.004957181BPRegulation of cell proliferationGO:00421270.005093297CCCell surfaceGO:00099861.56E−05CCExtracellular spaceGO:00056154.72E−05CCLysosomal membraneGO:00057658.70E−05CCBasement membraneGO:00056042.91E−04CCExtracellular exosomeGO:00700626.57E−04CCCytoplasmGO:00057370.001345772CCVesicleGO:00319820.001931725CCProteasome core complex, beta-subunit complexGO:00197740.00195821CCExternal side of plasma membraneGO:00098970.003278416CCPerinuclear region of cytoplasmGO:00484710.005507469MFHeparin bindingGO:00082012.32E−04MFThreonine-type endopeptidase activityGO:00042980.002408191MFTranscriptional activator activity, RNA polymerase II transcription regulatory region sequence-specific bindingGO:00012280.017925491MFEndopeptidase activityGO:00041750.019055149MFCXCR3 chemokine receptor bindingGO:00482480.026128457MFPhospholipase A2 inhibitor activityGO:00198340.026128457MFL-serine transmembrane transporter activityGO:00151940.031272543MFCalcium-dependent phospholipid bindingGO:00055440.033045687MFPeptidase activityGO:00082330.034029365MFTranscription factor activity, sequence-specific DNA bindingGO:00037000.035206985Table 3
KEGG enrichment analysis.
CategoryTermP valueKEGG_pathwayPertussis1.50E−02KEGG_pathwayAlcoholic liver disease1.60E−02KEGG_pathwayComplement and coagulation cascades2.30E−02KEGG_pathwayProteasome4.40E−02KEGG_pathwayTNF signaling pathway4.50E−02


## 3.3. Constructed PPI Network for Hub Gene Identification

To construct a PPI network, we uploaded 100 DEGs to the STRING database (version: 11.5;https://cn.string-db.org/) (Figure 4(a)). Thereafter, we analyzed the top six hub genes using the Cytoscape plug-in. Based on five algorithms of the Cytohubba plug-in (Figure 4(b)), four genes (Psmb6, Psmb8, Psmb10, and Atf3) were identified as key genes of RIRI. These genes were selected for further validation in the RIRI model.Figure 4
Building a PPI network of the DEGs and screening out hub genes: (a) PPI network exported from the STRING database; (b) Venn diagram for identification of hub genes by combining 5 algorithms of Cytohubba plug-in.
(a)(b)


## 3.4. Constructed Rat RIRI Model

As shown in Figure5(a), H&E staining of rat kidneys helped in finding the RIRI group (6 and 24 hours after reperfusion groups), which had significant proximal tubular injury. Concurrently, serum creatinine and BUN levels were significantly higher in the RIRI group than in the sham group (Figure 5(b)). Additionally, the expression of kidney injury molecule-1 (Kim-1) was positively correlated with the severity of RIRI. Kim-1 was more highly expressed in the RIRI group than in the sham group (Figures 5(c) and 5(d)). Altogether, these results indicated that the RIRI model was successfully established.Figure 5
Results of animal model: (a) hematoxylin and eosin (H&E) staining (200x, left) and tubular injury scores (right). The arrow indicates that the injury was serious. (b) The levels of serum creatinine and blood urea nitrogen (BUN) in different groups. (c, d) Western blot analysis of Kim-1 expression in different groups (∗∗∗∗P<0.0001;∗∗∗P<0.001;∗∗P<0.01).
(a)(b)(c)(d)


## 3.5. External Validation and Verification of Hub Genes

We performed ROC curve analysis to predict the values of four hub genes. The AUC value forPsmb6, Psmb8, Psmb10, and Atf3 was 0.719, 1, 1, and 0.906, respectively, in the GSE148420 dataset (Figures 6(a)–6(d)). Consistent with our analysis of GSE27274, these four hub genes were highly expressed in the RIRI group (6 and 24 hours after reperfusion groups). However, due to lack of statistical significance, they were placed in Figures S1(a)–1(h). Additionally, we performed western blot analysis and found that Psmb6, Psmb8, and Psmb10 were highly expressed in the RIRI group (Figures 6(f)–6(k)), whereas Atf3 was highly expressed in the 6 hours after reperfusion group (Figures 6(e) and 6(l)).Figure 6
ROC curve to verify gene expression in the GSE148420 dataset and western blot in the animal model: (a–d) ROC curve of four genes in the GSE148420 dataset; (e–h) the expression of four proteins in different groups (∗P<0.05;∗∗P<0.01;∗∗∗∗P<0.0001).


## 3.6. Identification of Candidate Drugs Associated with Renal IRI

Eight candidate drugs (Table4) with adjusted P value < 0.01 were identified in the DSigDB database. Psmb6 and Psmb8 were the most common drug targets, indicating these two genes might play a key role in the development of RIRI.Table 4
Prediction of candidate drugs for hub genes.
TermAdjustedP valueGenesVelcade (bortezomib) BOSS1.63E−06PSMB6; ATF3; PSMB8; PSMB10Carfilzomib CTD 000047871.63E−06PSMB6; PSMB8; PSMB10Carfilzomib BOSS1.63E−06PSMB6; PSMB8; PSMB10Velcade (bortezomib) TTD 000117782.16E−06PSMB6; PSMB8; PSMB10Carfilzomib6.66E−04PSMB6; PSMB8Bortezomib0.001156865PSMB6; PSMB8Clopamide HL60 down0.001156865PSMB6; PSMB8; PSMB10Azacyclonol HL60 down0.008485756PSMB6; PSMB8; PSMB10


## 4. Discussion

IRI is an inevitable procedure in kidney transplantation, and the severity of kidney injury has a close correlation with graft function. Despite several studies on IRI, efforts to improve long-term kidney transplant outcomes remain challenging. Only a better understanding of RIRI pathogenesis and identification of possible targets can improve the progress of grafts.A total of 100 DEGs were identified, and pathway enrichment analyses were performed. From the GO enrichment analysis, aging, positive regulation of apoptotic process, response to ethanol, response to drug, response to xenobiotic stimulus, and cellular response to lipopolysaccharides (LPS) were the top biological processes. Cell surface, extracellular space, and lysosomal membrane were the top cellular components. Heparin binding, threonine-type endopeptidase activity, transcriptional activator activity, and RNA polymerase II transcription regulatory region sequence-specific binding were the top molecular functions. As a special form of cell death, apoptosis, characterized by cell shrinkage, DNA fragmentation, and activation of caspases, was first proposed by Kerr et al. [20]. During RIRI, tubular epithelial cells are mostly severely damaged, and kidney injury can lead to a sharp decline in renal function. LPS is the major cell wall component of Gram-negative bacteria [21]. In the classical way, activation of TLR4 by LPS requires the interaction between LPS-binding proteins, followed by a series of reactions. LPS-TLR4 binding activates downstream signaling of the TLR4 pathway, which further leads to the activation of NF-κB and inflammation.KEGG enrichment analysis was mainly enriched in the TNF signaling pathway, proteasome, complement, and coagulation cascades. The inflammatory cytokine, TNF-α, has been established to play an important role in IRI. Several studies [22, 23] have indicated that TNF-α mRNA levels were raised significantly in the RIRI group, compared with the sham group. Furthermore, Nagata et al. [24] reported that RIRI can be ameliorated by the anti-TNF-α agent, infliximab. The complement system, an important component of the innate immunity system, is an important participant in the immune-inflammatory reaction, which comprises >40 blood-circulating, membrane-associated, and intracellular proteins. Many reports [25, 26] have confirmed that complement is the key mediating factor of chronic renal tubulointerstitial fibrosis after IRI. Activation of the complement system involves three activation pathways: the classical pathway, alternative pathway, and lectin pathway. During RIRI, strong inflammatory complement substances with biological activity are produced by the activation of complement; they include anaphylatoxins (C3a and C5a) and the cytolytic membrane attack complex (MAC) [27].Additionally, we constructed a PPI network and identified four hub genes:Atf3, Psmb6, Psmb8, and Psmb10. According to our western blot analysis results, the proteins of Psmb6, Psmb8, and Psmb10 were consistently highly expressed in the 6-hour and 24-hour reperfusion groups, whereas Atf3 was only highly expressed in the 6 hours after reperfusion group. Our rt-PCR results indicated a higher tendency in the RIRI group compared to the sham group. However, this tendency was not statistically significant. The lack of statistical significance was probably due to the small sample size and reduced mRNA level changes of these genes compared to protein level changes.Psmb6, Psmb8, and Psmb10 are allβ-subunits of the 20S proteasome core components. The proteasome complex 26S [28] consists of the 20S catalytic core and 19S regulatory particle. The 20S core catalyzes the chamber, which contains β-subunits (β1, β2, and β5), and the 19S regulatory particle binds the polyubiquitinated substrates. Proteasome 20S has been reported to be a probable biomarker of various diseases, such as renal cell carcinoma [29], psoriasis [30], and myocardial ischemia-reperfusion injury [31].Currently, a study [32] on hypoxia-induced pulmonary vascular in psmb6 has been reported. The expression of psmb6 increased after chronic hypoxia. During RIRI, the production of a large amount of reactive oxygen species (ROS) promoted the activation of Nrf2, which in turn plays a critical role in the protection against ROS-mediated injury. The Nrf2 agonist can alleviate tunicamycin-induced endoplasmic reticulum (ER) stress and increase the expression of Psmb6 [33]. These findings were consistent with our results. Basler et al. [34] found that psmb8 inhibitor ameliorated the pathological symptoms of dextran sulfate sodium-induced colitis by reducing inflammation. Li et al. [35] found that the deletion of psmb10 attenuates ang II-induced atrial inflammation and oxidative stress. Hypoxia, inflammation, and oxidative stress, as widely known, are closely related to IRI [24, 36]. However, these three proteasome subunits have not been reported to be associated with IRI. Our results demonstrated for the first time that the expression levels of these proteins were closely related to the development of RIRI.Based on the DSigDB database, we found some drugs that were related to these targets by screening for drugs with adjustedP<0.01. These drugs are listed in Table 4. Clopamide [37] is an oral diuretic used to treat some diseases, such as cardiac failure, nephrosis, chronic kidney failure, and cirrhosis. Azacyclonol [38], a metabolite of terfenadine, is a central nervous system inhibitor. Among these drugs, carfilzomib and borterzomib [39, 40] have been reported to be associated with IRI. Both drugs are proteasome inhibitors, and their potentially therapeutic targets are Psmb6 and Psmb8. Wu et al. [39] found that carfilzomib could reverse BNIP3L degeneration and restore mitophagy to alleviate brain IRI. However, borterzomib has different effects on the outcome of IRI-indifferent organs. Liu et al. [40] found that borterzomib alleviated myocardial IRI by activating the Nrf2/HO-1 signaling pathway. Nevertheless, borterzomib exacerbates RIRI, despite decreased secretion of infiltrating T cells and proinflammatory factors [41].This study had some limitations. First, these findings were obtained fromin vivo experiments but not from in vitro experiments and clinical samples. Second, these proteins had higher expression levels in the IRI group. However, their role in RIRI is unknown. Hence, further studies on the role of hub genes in RIRI should be conducted.


## 5. Conclusions

By bioinformatic analysis, we found some novel genes associated with the development of RIRI, and these hub genes are potential therapeutic targets for RIRI. Drugs based on target prediction may improve RIRI outcomes.



---

*Source: 7629782-2023-02-02.xml*