Identification of Key Genes Involved in Diabetic Peripheral Neuropathy Progression and Associated with Pancreatic Cancer.

INTRODUCTION: Diabetes mellitus (DM) patients suffer from high morbidity and premature mortality due to various diabetic complications and even cancers. Therefore, this study aimed to identify key genes involved in the pathogenesis of diabetic peripheral neuropathy (DPN) and pancreatic cancer (PC).
METHODS: We analyzed three gene expression profiles (GSE95849, GSE28735 and GSE59953) to obtain differentially expressed genes (DEGs). Then, Gene Ontology (GO) analysis and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis were performed by using the Database for Annotation, Visualization, and Integrated Discovery (DAVID). The Search Tool for the Retrieval of Interacting Genes/Proteins (STRING) database was then used to establish a protein-protein interaction (PPI) network. The MCODE and cytoHubba plug-ins of Cytoscape were used to select hub genes. Finally, survival analysis of the hub genes was performed using the Kaplan-Meier plotter and GEPIA online tool.
RESULTS: We first analyzed GSE95849 to obtain DPN-related genes. DEGs were obtained from three groups in GSE95849. The DEGs were enriched in the Toll-like receptor signaling pathway, hematopoietic cell lineage and chemokine signaling pathway. Importantly, we identified three shared genes as hub genes, including TLR4, CCR2 and MMP9. We then analyzed and integrated GSE95849 and GSE28735 to obtain genes common in DM and PC. A total of 58 mutual DEGs were identified, and these DEGs were enriched in the ECM-receptor interaction, focal adhesion and pathways in cancer. Five hub genes (including PLAU, MET, CLU, APOL1 and MMP9) were associated with the overall survival of PC patients. However, the results from the analysis of GSE59953 showed that hyperglycemia or TGF-β1 treatment did not affect the expression level of these hub genes, but the DEGs based on hyperglycemia or TGF-β1 treatment were mostly enriched in the ECM-receptor interaction, focal adhesion and pathways in cancer. Finally, functional enrichment analysis of MMP9 showed that significant genes correlated with MMP9 were associated with the tumorigenicity of cancers, insulin resistance, development of DM and inflammation.
CONCLUSION: In summary, inflammation and immunity-related pathways may play an important role in DM and DPN, while the ECM-receptor interaction, focal adhesion and pathways in cancer pathways may play significant roles in DM and PC. MMP9 may be used as a prognostic marker for PC and may be helpful for the treatment of DM, DPN and PC.

Introduction

Diabetes mellitus (DM) is one of the largest public health problems worldwide. Diabetic patients suffer from high morbidity and premature mortality due to the development of various diabetic complications, including cardiomyopathy, vasculopathy, nephropathy, retinopathy, and neuropathy.1 In addition, diabetes mellitus was found to be associated with increased pancreatic cancer risk.2

Diabetic peripheral neuropathy (DPN) is the most common diabetic complication, characterized by pain, paraesthesia and sensory loss, afflicting over 50% of people with diabetes.3 Previous studies have showed that the pathogenesis and underlying mechanisms leading to DPN are complex. Diabetes-related metabolic factors such as increased glucose, decreased insulin, and increased lipids produce changes underlying the development of diabetic neuropathy. Injury to neurons, microvascular endothelium, and Schwann cells in DM contributes to the pathogenesis of neuropathy. In addition, various stress pathways, including oxidative stress, inflammation, and apoptosis, are involved in the development of diabetic neuropathy.3–5 Crucial advances have been achieved in understanding the pathogenesis of DPN. However, controversy remains due to its multifactorial etiology.

Diabetes mellitus is associated with increased pancreatic cancer risk.2 A previous study revealed that the presence of diabetes mellitus significantly increased the risk of the subsequent development of pancreatic cancer by 2-fold.6 Interestingly, studies found that metformin use can prolong overall survival and lower mortality in patients with PC and DM.7,8 However, the underlying mechanisms and potential interconnections between diabetes and pancreatic cancer are not fully well known.

Therefore, we performed this bioinformatics analysis by using data from the GEO database. We first identified the key genes and pathways involved in DPN and then analyzed the key genes and pathways involved in both DM and PC. Interestingly, we identified MMP9 as the common hub gene of DM, DPN and PC. The analysis results will provide meaningful clues for the treatment of DM, DPN and PC.

Materials and Methods

Microarray Data Information and DEGs Identification

GSE95849, GSE28735 and GSE59953 gene expression profiles were downloaded from the GEO database. DEGs were obtained using the web tool GEO2R (). The data set GSE958499 consisted of 6 Diabetic peripheral neuropathy (DPN), 6 diabetes mellitus (DM), and 6 healthy controls (CN). To obtained DEGs associated with DPN, the conditions of |log FC| ≥ 2, P-value< 0.05 was used for DM vs CN and DPN vs CN group. While |log FC| ≥1.5, P-value< 0.05 was used for DPN vs DM group. GSE2873510,11 contained 45 pairs of pancreatic cancer (PC) and adjacent non-tumor tissues. To obtained DEGs associated with DM and PC, the cutoff criterion of |log FC| ≥1, P-value< 0.05 was used for DM vs CN group of GSE95849 and PC vs CN group of GSE28735. Then, venn diagram was composed to select mutual DEGs between GSE95849 and GSE28735. GSE5995312 contained 4 samples of pancreatic stellate cells (PSCs) that cultured in normal glucose concentration (control), 4 samples of PSCs that exposed to 21 d hyperglycemia (21 d treatment), 4 samples of PSCs that treated with 48 h TGF-β1 treatment (48 h treatment) and 4 samples of PSCs that exposed to 21 d hyperglycemia + subsequent TGF-β1 treatment (21 d + 48 h treatment). To obtained DEGs, the conditions of |log FC| ≥ 0.5, P-value< 0.05 was used for control vs 21 d treatment group and control vs 48 h treatment group. While |log FC| ≥1.5, P-value< 0.05 was used for control vs 21 d + 48 h treatment group.

Genetic Ontology and Pathway Enrichment Analysis of DEGs

In order to explore the biological functional roles of the above DEGs, the DAVID was used to GO and KEGG pathway enrichment analysis. P<0.05 was used as the cutoff criterion.

Establishment of PPI Network and Modules Selection

The STRING database was used to construct PPI network. The Cytoscape software () was used to visualize and analyze the PPI network. The Molecular Complex Detection (MCODE) plug-in Cytoscape was used to select the significant modules.

Identification of Hub Genes

The degrees method of cytoHubba or the MCODE scores of the Cytoscape software was used to select important hub genes among DEGs. For DEGs about DPN, the top 30 genes with high degrees or MCODE scores ≥12 were taken as the cutoff criterion. A Venn diagram was composed to select mutual genes from two methods. The mutual genes from the three groups were selected as the hub genes. For DEGs between DM and PC, the top 15 genes with high degrees were taken as the hub genes. For DEGs from GSE59953, the top 30 genes with high degrees were showed.

Survival Analysis

The Kaplan–Meier plotter ()13 and the Gene Expression Profiling Interactive Analysis (GEPIA) ()14 databases were used for survival analysis of the hub genes. GEPIA was also used for expression levels analysis of the hub genes.

LinkedOmics Analysis

The LinkedOmics database ()15 contains multi-omics data and clinical data for 32 cancer types and a total of 11,158 patients from The Cancer Genome Atlas (TCGA) project. The LinkFinder module of LinkedOmics was used to analyze genes differentially expressed in correlation with MMP9 in the TCGA PC cohort (n=178). Results were analyzed using Spearman correlation coefficient. The Link-Interpreter module of LinkedOmics was used to perform pathway and network analysis of differentially expressed genes. Gene set enrichment analysis (GSEA) was used to perform GO and KEGG pathway analysis of the data from the LinkFinder results. The rank criterion was a false discovery rate (FDR) < 0.05.

Results

Identification of DEGs Among Different Groups of GEO Data Sets

In the present study, identification of key genes and pathways involved in DPN, DM and PC were performed by integrated bioinformatics analysis (Figure 1). We first identified DEGs among different groups of GSE95849, which contained 6 DPN, 6 DM, and 6 CN samples. For DM vs CN group, 497 DEGs were screened. While a total of 1362 DEGs were screened from DPN vs CN group and 842 DEGs from DPN vs DM group (Figure 1). We second analyzed and integrated GSE95849 and GSE28735 to identify genes between DM and PC. A total of 3739 DEGs were selected from GSE95849 and 413 DEGs from GSE28735 (Figure 1).

Figure 1

Workflow of the present study.

Abbreviations: PC, pancreatic cancer; DEG, differentially expressed gene; GEO, Gene Expression Omnibus; GO, Gene Ontology; KEGG, Kyoto Encyclopedia of Genes and Genomes; PPI, protein–protein interaction; DPN, Diabetic peripheral neuropathy; DM, diabetes mellitus; CN, healthy controls; DAVID, Database for Annotation, Visualization, and Integrated Discovery; STRING, Search Tool for the Retrieval of Interacting Genes; MCODE, Molecular Complex Detection; GEPIA, Gene Expression Profiling Interactive Analysis.

GO and KEGG Pathway Analysis of DEGs Among Different Groups in GSE95849

We then performed GO and KEGG pathway analysis of the above DEGs by using DAVID. For DM vs CN group, the DEGs were enriched in defense response, inflammatory response, response to wounding, immune response and chemotaxis in the biological process (BP) category. Moreover, these genes were significantly enriched in the Golgi apparatus, cell fraction, cytoplasmic vesicle, vesicle and membrane fraction in the cellular component (CC) category. They were enriched in cytokine binding, protein domain-specific binding, chemokine receptor activity, chemokine binding and GTPase binding in the molecular function (MF) category (Table 1, left column). For DPN vs CN group, the DEGs were enriched in inflammatory response, response to wounding, defense response, immune response and innate immune response in the BP category. For the CC category, the DEGs were enriched in the Golgi apparatus, cell fraction, cytoplasmic vesicle, vesicle and membrane fraction. For the MF category, the DEGs were enriched in cytokine binding, protein domain-specific binding, chemokine receptor activity, chemokine binding and GTPase binding (Table 1, middle column). For DPN vs DM group, the DEGs were enriched in response to wounding, behavior, inflammatory response, cell activation and defense response in the BP category. With regard to the CC category, the DEGs were enriched in cell fraction, axon part, cell surface, vesicle membrane and soluble fraction. For the MF category, the DEGs were enriched in protein complex binding, chemokine receptor binding, chemokine activity, transcription factor binding and transcription factor activity (Table 1, right column).

Table 1

Functional Enrichment Analysis for the DEGs of Different Groups of GSE95849

	DM vs CN			DPN vs CN			DPN vs DM
Category	Term	Count	P-value	Term	Count	P-value	Term	Count	P-value
BP	GO:0006952~defense response	56	8.29E-16	GO:0006954~inflammatory response	63	2.73E-13	GO:0009611~response to wounding	52	5.2E-08
	GO:0006954~inflammatory response	39	7.92E-15	GO:0009611~response to wounding	85	7.92E-13	GO:0007610~behavior	45	9E-07
	GO:0009611~response to wounding	50	9.89E-15	GO:0006952~defense response	88	1.54E-10	GO:0006954~inflammatory response	35	1.5E-06
	GO:0006955~immune response	53	4.69E-12	GO:0006955~immune response	91	4.89E-09	GO:0001775~cell activation	31	6.2E-06
	GO:0006935~chemotaxis	22	1.21E-09	GO:0045087~innate immune response	27	3.05E-06	GO:0006952~defense response	51	1E-05
CC	GO:0005794~Golgi apparatus	94	2.00E-05	GO:0005794~Golgi apparatus	94	2.00E-05	GO:0000267~cell fraction	63	0.003
	GO:0000267~cell fraction	108	1.16E-04	GO:0000267~cell fraction	108	1.16E-04	GO:0033267~axon part	8	0.00574
	GO:0031410~cytoplasmic vesicle	70	1.96E-04	GO:0031410~cytoplasmic vesicle	70	1.96E-04	GO:0009986~cell surface	25	0.00748
	GO:0031982~vesicle	70	6.59E-04	GO:0031982~vesicle	70	6.59E-04	GO:0012506~vesicle membrane	14	0.00814
	GO:0005624~membrane fraction	81	8.56E-04	GO:0005624~membrane fraction	81	8.56E-04	GO:0005625~soluble fraction	23	0.0082
MF	GO:0019955~cytokine binding	23	4.46E-06	GO:0019955~cytokine binding	23	4.46E-06	GO:0032403~protein complex binding	20	0.00062
	GO:0019904~protein domain specific binding	45	2.16E-05	GO:0019904~protein domain specific binding	45	2.16E-05	GO:0042379~chemokine receptor binding	9	0.00092
	GO:0004950~chemokine receptor activity	9	1.92E-04	GO:0004950~chemokine receptor activity	9	1.92E-04	GO:0008009~chemokine activity	8	0.00291
	GO:0019956~chemokine binding	9	3.47E-04	GO:0019956~chemokine binding	9	3.47E-04	GO:0008134~transcription factor binding	35	0.00606
	GO:0051020~GTPase binding	18	0.001287	GO:0051020~GTPase binding	18	0.001287	GO:0003700~transcription factor activity	58	0.00697

Abbreviations: DM, diabetes mellitus; CN, healthy controls; DPN, diabetic peripheral neuropathy; DEGs, differentially expressed genes; BP, biological process; CC, cellular component; MF, molecular function.

KEGG pathway enrichment analysis showed that the DEGs were enriched in the cytokine-cytokine receptor interaction, Toll-like receptor signaling pathway, hematopoietic cell lineage, PPAR signaling pathway and chemokine signaling pathway in the DM vs CN group (Table 2, up column). With regard to the DPN vs CN group, the DEGs were enriched in the Toll-like receptor signaling pathway, apoptosis, hematopoietic cell lineage, citrate cycle (TCA cycle), the NOD-like receptor signaling pathway, the chemokine signaling pathway, the MAPK signaling pathway and the B cell receptor signaling pathway (Table 2, middle column). For the DPN vs DM group, the DEGs were enriched in the NOD-like receptor signaling pathway, the MAPK signaling pathway, the neurotrophin signaling pathway, the B cell receptor signaling pathway, the chemokine signaling pathway, hematopoietic cell lineage, the Toll-like receptor signaling pathway, citrate cycle (TCA cycle) and cytokine-cytokine receptor interaction (Table 2, bottom column). The Toll-like receptor signaling pathway, hematopoietic cell lineage and chemokine signaling pathway were shown in all three groups. Therefore, these three signaling pathways may play an important role in DM and DPN.

Table 2

The KEGG Pathway Enrichment Analysis for the DEGs of Different Groups of GSE95849

Category	Term	Count	%	P-value
DM v CN
KEGG_PATHWAY	hsa04060:Cytokine-cytokine receptor interaction	24	0.57859	1.44E-05
	hsa04620:Toll-like receptor signaling pathway	10	0.24108	0.00606
	hsa04640:Hematopoietic cell lineage	9	0.21697	0.00737
	hsa03320:PPAR signaling pathway	8	0.19286	0.00753
	hsa04062:Chemokine signaling pathway	12	0.2893	0.04589
DPN v CN
KEGG_PATHWAY	hsa04620:Toll-like receptor signaling pathway	18	0.142	0.00431
	hsa04210:Apoptosis	16	0.12622	0.00568
	hsa04640:Hematopoietic cell lineage	15	0.11833	0.01226
	hsa00020:Citrate cycle (TCA cycle)	8	0.06311	0.01298
	hsa04621:NOD-like receptor signaling pathway	12	0.09467	0.01361
	hsa04062:Chemokine signaling pathway	25	0.19722	0.02449
	hsa04010:MAPK signaling pathway	32	0.25245	0.04111
	hsa04662:B cell receptor signaling pathway	12	0.09467	0.04898
DPN v DM
KEGG_PATHWAY	hsa04621:NOD-like receptor signaling pathway	15	0.19214	1.56E-06
	hsa04010:MAPK signaling pathway	31	0.39708	2.19E-05
	hsa04722:Neurotrophin signaling pathway	16	0.20494	0.00127
	hsa04662:B cell receptor signaling pathway	11	0.1409	0.00401
	hsa04062:Chemokine signaling pathway	19	0.24337	0.00573
	hsa04640:Hematopoietic cell lineage	11	0.1409	0.01057
	hsa04620:Toll-like receptor signaling pathway	12	0.15371	0.01205
	hsa00020:Citrate cycle (TCA cycle)	6	0.07685	0.01803
	hsa04060:Cytokine-cytokine receptor interaction	21	0.26899	0.03948

Abbreviation: KEGG, Kyoto Encyclopedia of Genes and Genomes.

PPI Network Analysis and Screening for Hub Genes

We further imported the above DEGs into STRING to obtain a PPI network. The interaction network was then imported into Cytoscape to screen the hub genes using two methods, including the degree method scores (Degrees) and MCODE. First, the cytoHubba plug-in of the Cytoscape software was used to select genes with the highest degrees. The top 30 genes with the highest degrees of connectivity for each group are shown in Table 3. Second, the MCODE plug-in of Cytoscape was used to select the most significant modules (Figure 2A–C). The genes with MCODE score ≥ 12 in each of the three modules are shown in Table 3. The mutual genes from the two methods were selected for each group. In the end, the Venn diagram identified three overlapping genes as representative hub genes from the three groups (Figure 2D), including TLR4, CCR2 and MMP9. Therefore, these hub genes may serve as promising biomarkers of DM and DPN.

Table 3

Identified of the Hub Genes of Different Groups of GSE95849 by Degree and MCODE-Score of Cytoscape Sofeware

DM v CN				DPN v CN				DPN v DM
Name	Degree	Name	MCODE_Score	Name	Degree	Name	MCODE_Score	Name	Degree	Name	MCODE_Score
TLR2	76	CCR5	15	TP53	181	FN1	14	TNF	91	ICAM1	13
TLR4	76	TLR2	15	TLR4	114	TP53	13	MAPK3	79	CCR2	13
TLR8	67	CXCL2	15	TLR2	110	CASP1	12	JUN	76	NFKBIA	14
IL1B	65	TLR8	15	FN1	104	TLR8	15	CXCL8	69	IL17A	13
MMP9	61	TLR4	15	TLR8	101	SELL	14	FOS	68	CCL19	13
CXCL1	60	CXCR2	15	STAT1	86	CYBB	12	IL1B	62	RELB	13
FPR2	55	CXCL1	15	TLR1	85	CXCL1	14	TLR4	60	CXCL10	13
C3AR1	53	CCR2	15	TLR7	84	TLR6	15	CXCL10	57	TNFRSF1A	14
SELL	45	MMP9	14	CD44	83	FCGR2B	12	MMP9	56	CSF1	13
CCR5	44	IL1B	14	FOS	83	CCR2	14	ICAM1	55	MMP9	13
FPR1	42	C3AR1	14	MMP9	79	MMP9	15	ESR1	53	FOS	12
MYD88	42	FPR2	14	CCR5	76	MAPK14	14	TLR7	52	CD44	13
PTGS2	41	TLR6	14	CYBB	76	FOS	14	CD44	51	JUN	13
CXCR2	41	CCR1	13	C3AR1	73	TLR7	15	IL17A	46	MAPK3	14
HCK	40	TREM1	13	MAPK14	72	MYD88	13	CCR2	42	TLR4	13
TLR6	40	CX3CR1	13	CXCL1	72	STAT1	15	CXCL2	40	IL1B	13
MNDA	39	MYD88	13	HSPA8	71	JAK2	13	BDNF	40	CASP1	14
ELANE	39	ARG1	13	SELL	70	CSF1R	13	NFKBIA	38	CX3CR1	13
TREM1	37	PTGS2	13	ESR1	67	TLR4	15	CAT	37	CXCL2	12
ARG1	37	SELL	12	MYD88	67	PTGS2	14	CSF1	35	TNF	13
CXCL2	37	FPR1	12	PTGS2	66	CD44	14	CX3CR1	35	CCL20	13
CCR1	36	IL1RN	12	CSF1R	65	TLR2	15	EDN1	35	CXCL8	13
TLR5	36			JAK2	65	TLR1	15	CCL20	35	TLR7	13
CCR2	35			FPR2	64	CCR5	14	LPAR1	34
CAMP	35			TLR6	64			CASP1	34
CD53	34			MNDA	63			IRF8	33
TOLLIP	34			CCR2	63			CCL19	32
CD93	33			FCGR2B	62			TNFRSF1A	32
IL1RN	33			CASP1	61			NPY	31
TLR10	31			HCK	60			RELB	30
GPR84	31
CLEC4D	31
CD36	31
CX3CR1	31

Abbreviation: MOCDE, Molecular Complex Detection.

Figure 2

Establishment of modules and identified of hub genes for DPN. (A) PPI network of the module for DM vs CN. (B) PPI network of the module for DPN vs CN. (C) PPI network of the module for DPN vs DM. (D) hub genes identified by Venn diagram from the three modules.

Abbreviations: PPI, protein–protein interaction; DPN, Diabetic peripheral neuropathy; DM, diabetes mellitus; CN, healthy controls.

Identification of Hub Genes Between DM and PC

DM has been identified as a risk factor for PC. We therefore tried to identify hub genes between DM and PC by analyzing the data set GSE95849 and GSE28735 datasets. A total of 3739 DEGs were screened between DM and CN in GSE95849, and 413 DEGs were screened between PC and CN in GSE28735. The 59 mutual DEGs were found via a Venn diagram (Figure 3A) and were imported into DAVID to perform GO and KEGG pathway analysis (). Functional analysis showed that the mutual DEGs were enriched in response to wounding, cell adhesion, biological adhesion, cell migration, immune response, defense response, inflammatory response, wound healing and blood vessel morphogenesis in the BP category. For the CC category, the mutual DEGs were enriched in the extracellular region, membrane fraction, insoluble fraction, cell fraction, extracellular region part, proteinaceous extracellular matrix, extracellular matrix, extracellular space, extracellular matrix part, platelet alpha granule, integral to plasma membrane, intrinsic to plasma membrane and basement membrane. With regard to the MF category, the mutual DEGs were enriched in glycosaminoglycan binding, hyaluronic acid binding, polysaccharide binding, pattern binding, growth factor binding, carbohydrate binding, vitamin binding and amine transmembrane transporter activity. KEGG pathway analysis showed that the mutual DEGs were enriched in ECM-receptor interaction, focal adhesion, pathways in cancer and the adipocytokine signaling pathway (). These results indicated that the three pathways may play an important role in DM and PC.

Figure 3

Establishment of modules and identified of hub genes between DM and PC. (A) Venn diagram of mutual hub genes based on GSE95849 and GSE28735. (B) Entire PPI network. (C) PPI network of the most significant module. (D) Top 15 hub genes with the high degrees.

The 59 mutual DEGs were then imported into STRING to construct a PPI network, which had 31 nodes and 64 edges (Figure 3B). The interaction network was then imported into Cytoscape to screen the hub genes using the Degree algorithm. The interaction network with the top 15 genes with the highest degrees was created and is shown in Figure 3C. The degrees of the top 15 genes (including FN1, MMP9, PLAU, VCAN, LCN2, MET, CCL20, ANGPT2, CLU, LYVE1, APOL1, CEACAM1, TCN1, SLP1 and LTBP1) are shown in Figure 3D. Interestingly, the significant module generated by MCODE contained four genes, including MMP9, PLAU, MET and ANGTP2 (Figure 3C).

We then used the GEPIA and Kaplan–Meier plotter online databases to analyze the expression level and prognostic value of the 15 hub genes. The analysis showed that only five hub genes (including PLAU, MET, CLU, APOL1 and MMP9) were associated with the overall survival of PC patients ( and Figure 4B and C). The expression of PLAU, MET, APOL1 and MMP9 was significantly higher in PC tissues than in paired normal tissues ( and Figure 4A). However, the expression level of CLU was not significantly different between PC tissues and paired normal tissues (). These results indicated that the five hub genes may serve as promising biomarkers of DM and PC.

Figure 4

Expression level and prognostic value analysis of MMP9 for PC patients. (A) Expression level of MMP9 by GEPIA. (B) Overall survival analysis of MMP9 by GEPIA. (C) Overall survival analysis of MMP9 by Kaplan–Meier plotter.

Note: *Means the p-value < 0.01.

Diabetes mellitus was associated with increased pancreatic cancer risk. In addition, a study found that elevated fasting plasma glucose (FPG) levels are risk factors for PC.16 A study found that activated human pancreatic stellate cells (PSCs) are the major cellular constituents of the pancreatic ductal adenocarcinoma stroma.17 Katalin Kiss et al kept PSCs under chronic hyperglycemic conditions; the PSCs were subsequently either treated or untreated with TGFβ1 to model the gene expression changes, and the GSE59953 dataset was generated with this data.12 In vivo, they found that chronic hyperglycemia induced the transdifferentiation of PSCs and enhanced the malignant molecular communication with PSCs.12 Therefore, we used GSE59953 dataset to analyze the effect of chronic hyperglycemia (CHG) on gene expression. For the control vs 21 d treatment group, a total of 854 DEGs were identified. For the control vs 48 h treatment group, a total of 205 DEGs were identified. For the control vs 21 d + 48 h treatment group, a total of 1167 DEGs were identified. KEGG pathway analysis showed that the DEGs were also mostly enriched in the ECM-receptor interaction, focal adhesion and pathways in cancer (). We also constructed a PPI network and selected the significant module to identify hub genes (– and ). In the end, the Venn diagram identified six overlapping genes (CXCL8, EGR1, FN1, FOS, PPARG and SPP1) as representative hub genes from the three groups (). Unfortunately, the hub genes identified from GSE95849 and GSE28735 were not significantly altered under hyperglycemic conditions. However, they may share similar pathways, such as the ECM-receptor interaction, focal adhesion and pathways in cancer.

Enrichment Analysis of MMP9 Functional Networks in PC

The comprehensive results revealed that MMP9 was found to be both involved in DPN progression and was associated with PC. We used the function module of LinkedOmics to analyze MMP9 mRNA sequencing data from 178 PC patients in the TCGA. As shown in the volcano plot in Figure 5A, MMP9 was positively correlated with 1550 genes (dark red dots) but negatively correlated with 729 genes (dark green dots) (FDR < 0.01). The 50 significant gene sets positively and negatively associated with MMP9 are shown in the heat map in Figure 5B and C. The top 20 significant genes positively and negatively associated with MMP9 are shown more detail in . MMP9 expression showed a strong positive or negative association with the expression of C1QTNF6 (Spearman correlation = 0.547, p-value= 1.00e–15), TWIST1 (Spearman correlation = 0.542, p-value = 1.00e–15), TGFB1(Spearman correlation =0.523, p-value =1.00E-15), IL4I1 (Spearman correlation =0.513, p-value =1.00E-15), OSCAR (Spearman correlation =0.498, p-value =1.00E-15), CASC2 (Spearman correlation = −0.449, p-value =4.98E-10), PTPRN2 (Spearman correlation =−0.422, p-value =6.12E-09), KIF3B (Spearman correlation =−0.417, p-value =7.05E-09) and DCLRE1A (Spearman correlation =−0.407, p-value =1.67E-08), which were associated with the tumorigenicity of cancers, insulin resistance, development of diabetes mellitus and inflammation. KEGG pathway analysis by GSEA showed that the genes correlated with MMP9 were significantly enriched in cytokine-cytokine receptor interaction (FDR=0, p-value=0), ECM-receptor interaction (FDR=0, p-value=0), hematopoietic cell lineage (FDR=0, p-value=0), B cell receptor signaling pathway (FDR=0, p-value=0), citrate cycle (TCA cycle) (FDR=0, p-value=0) and focal adhesion (FDR=0, p-value=0) (Figure 5D–J). Therefore, MMP9 has the potential to be used as a target for PC diagnosis and treatment.

Figure 5

Functional enrichment analysis of MMP9 networks in PC by LinkedOmics. (A) Volcano plot of genes differentially expressed in correlation with MMP9. (B, C) Heat maps of genes positively and negatively correlated with RBM8A (TOP 50). (D) KEGG pathway analysis of MMP9 by GSEA. (E) Cytokine-cytokine receptor interaction (FDR=0, p-value=0), (F) ECM-receptor interaction (FDR=0, p-value=0), (G) Hematopoietic cell lineage (FDR=0, p-value=0), (H) B cell receptor signaling pathway (FDR=0, p-value=0), (I) Citrate cycle (TCA cycle) (FDR=0, p-value=0) and (J) Focal adhesion (FDR=0, p-value=0) pathways showed by GSEA. The rank criterion was a false discovery rate (FDR) < 0.05.

Discussion

In the present study, we performed bioinformatics analysis for two purposes. One was to identify the key genes and biological pathways responsible for the progression of DPN. The other purpose was to identify the key genes and biological pathways between DM and PC. By using the GSE95849 dataset, functional enrichment analysis showed that most of the DEGs from the three compared groups (DM vs CN, DPN vs CN and DPN vs DM) were enriched in defense response, inflammatory response, immune response and response to wounding. KEGG pathway analysis found that the DEGs from the three compared groups shared three pathways, including the Toll-like receptor signaling pathway, hematopoietic cell lineage and chemokine signaling pathway. Moreover, in comparison with DM, the DPN-specific DEGs were also significantly enriched in the NOD-like receptor signaling pathway, MAPK signaling pathway, neurotrophin signaling pathway, B cell receptor signaling pathway, citrate cycle (TCA cycle) and cytokine-cytokine receptor interaction. Interestingly, most of these pathways were associated with inflammation and immune-related functions. Therefore, inflammation and immune-related pathways may play an important role in DM and DPN.

Importantly, we identified three hub genes, including TLR4, CCR2 and MMP9. TLR4 belongs to the Toll-like receptor (TLR) family, which plays an important role in the innate immune response.18 In addition, TLR4 was found to play a critical role in the pathophysiological process of diabetes19 and was identified as a potential diagnostic biomarker for DPN.20 CC chemokine receptor 2 (CCR2) plays an important role in the recruitment and activation of macrophages and in turn contributes to both inflammatory and neuropathic pain states.21 CCR2 plays a critical role in a range of inflammatory diseases in the central nervous system (CNS).21 In addition, chemokine C-C motif ligand 2 (CCL2)-CCR2 signaling was involved in the maintenance of orofacial neuropathic pain.22 Blockade of the CCL2/CCR2 pathway improves diabetic nephropathy.23 Pharmacological modulation of the CCL2/CCR2 pathway results in the attenuation of neuropathic pain.24 Evidence has shown that matrix metalloproteinases (MMPs) can cleave multiple substrates of the CNS and contribute to promoting and dampening diseases of the CNS.25,26 Moreover, MMPs are involved in diabetes-induced retinal neuropathy.27 In addition, MMPs can degrade extracellular matrix (ECM) proteins, which in turn results in tissue remodeling and can contribute to cell migration.28,29 Moreover, MMP9 was found to be involved in the pathogenesis of DM and diabetic complications such as diabetic retinopathy.30 Therefore, these hub genes may serve as promising biomarkers of DM and DPN.

Another important finding is that we identified 15 hub genes between DM and PC by analyzing the GSE95849 and GSE28735 datasets. Five genes (MMP9, PLAU, MET, CLU and APOL1) were associated with the overall survival of PC. These results indicated that they may be identified as biomarkers and therapeutic targets of PC by further studies. We then focused on MMP9, as it was also identified as the key gene involved in DPN in the present study. We used the LinkedOmics database to conduct the functional enrichment analysis of MMP9 in PC. KEGG pathway analysis showed that the genes correlated with MMP9 were significantly enriched in cytokine-cytokine receptor interaction, ECM-receptor interaction, hematopoietic cell lineage, B cell receptor signaling pathway, citrate cycle (TCA cycle) and focal adhesion. Interestingly, KEGG pathway analysis of the 59 mutual DEGs showed that they were also significantly enriched in the ECM-receptor interaction and focal adhesion. In addition, MMPs can degrade extracellular matrix (ECM) proteins, which can in turn result in tissue remodeling and contribute to cell migration.28,29 MMPs were found to be involved in immunity and inflammation by regulating cytokine and chemokine function.31 Anti–inflammatory macrophages enhance the invasion of pancreatic cancer by inducing MMP9 and ADAM8 expression.32 In addition, MMP-9 expression was associated with poorer prognosis of pancreatic cancer.33 A previous study has found that the inhibition of MMPs in vitro can reduce metastases of pancreatic cancer.34 Interestingly, MMP-9 has been identified as an activator of nanosystems for targeted drug delivery in pancreatic cancer.35 In the preclinical models of pancreatic cancer, the addition of an anti-MMP9 antibody can prolong animal survival.36 These results may suggest that the MMP9 antibody could be exploited to improve clinical PC therapy. Moreover, most of the significant genes that were correlated with MMP9, such as C1QTNF6, TWIST, HK3, CASC2, TGF-β1 and OSCAR, were associated with the tumorigenicity of cancers, insulin resistance, development of diabetes mellitus and inflammation. C1q/tumor necrosis factor-related protein-6 (C1QTNF6, also known as CTRP6) has been found to be associated with insulin resistance and the development of diabetes in the Chinese population.37 CTRP6 has been identified as an important metabolic/immune regulator linking obesity to adipose tissue inflammation and insulin resistance.38 A study found that C1QTNF6 was upregulated in gastric carcinoma and contributed to the migration and proliferation of gastric carcinoma cells.39 Study showed that targeting twist family bHLH transcription factor 1 (TWIST1) through loss of function can inhibit the tumorigenicity of glioblastoma.40 High hexokinase 3 (HK3) expression was associated with epithelial-mesenchymal transition (EMT) in colorectal cancer.41 High levels of TGF-β1 were associated with the susceptibility to type 2 diabetes mellitus and type 2 diabetic nephropathy.42 In addition, stimulation with TGF-beta signaling can enhance stem cell properties in colorectal cancer.43 Upregulation of the lncRNA cancer susceptibility 2 (CASC2) suppressed the cell proliferation and metastasis of breast cancer via the inactivation of the TGF-β signaling pathway.44 lncRNA CASC2 expression was downregulated in patients with type 2 diabetes combined with chronic renal failure compared with healthy controls.45 The osteoclast associated Ig-like receptor (OSCAR)-surfactant protein D interaction can activate TNF-α release from human CCR2+ inflammatory monocytes46 and OSCAR-collagen signaling in monocytes plays a pro-inflammatory role and contributes to the pathogenesis of rheumatoid arthritis.47 These data indicated that MMP9 has the potential to be used as a target for the diagnosis and treatment of DM and PC.

Although bioinformatics technologies have the potential to identify candidate genes in various diseases, some limitations remain in this study. First, the RNA samples in GSE95849 were obtained from peripheral blood mononuclear cells (PBMCs) and not from peripheral nerves, so it may be difficult to make the conclusion that these hub genes contributed to the development of DPN. However, it is very difficult to conclude the collection of nerve biopsies in studies of patients with diabetic neuropathy. Transcriptome studies in PBMCs have emgeted as a new kind of transcriptional analysis and have been widely used in diabetic studies.48,49 Therefore, the results analyzed from this dataset are dependable. Second, the samples in GSE95849 and GSE28735 were obtained from DM and PC, respectively. Therefore, the linked between DM and PC was not directly assessed. Third, the samples in GSE59953 were from PSCs; thus, it is also difficult to fully reflect the real word. Furthermore, further biological experiments are needed to validate our results because our study was performed by bioinformatics analysis based on GEO datasets.

However, our data provide a comprehensive bioinformatics analysis of the DEGs and pathways that may be involved in DPN, DM and PC. We further confirmed that inflammation, immunity and their signaling pathways play an important role in DM and DPN. The ECM-receptor interaction, focal adhesion and pathways in cancer may play significant roles in DM and PC. More importantly, we identified MMP9 as the most important hub gene, which served as an inflammatory/immune regulator in DPN and linked DM to pancreatic cancer. We propose that careful targeting of MMP9 could be helpful for the treatment of DM, DPN and PC, but more studies are needed.