Transcriptome sequencing and Mendelian randomization analysis identified biomarkers related to neutrophil extracellular traps in diabetic retinopathy

Summary: In the development of diabetic retinopathy (DR), neutrophil infiltration hastens the adhesion between neutrophils and endothelial cells, leading to inflammation. Meanwhile, neutrophil extracellular traps (NETs) produced by neutrophils could clear aging blood vessels, setting the stage for retinal vascular regeneration. To explore the mechanism of NETs-related genes in DR, the transcriptome of NETs from normal and DR individuals were analyzed with gene sequencing and mendelian randomization (MR) analysis. Five NETs-related genes were identified as key genes. Among these genes, CLIC3, GBP2, and P2RY12 were found to be risk factors for Proliferative DR(PDR), whereas HOXA1 and PSAP were protective factors. Further verification by qRT-PCR recognized GBP2, P2RY12 and PSAP as NETs-associated biomarkers in PDR.

Purpose: To investigate neutrophil extracellular traps (NETs) related genes as biomarkers in the progression of diabetic retinopathy (DR).

Methods: We collected whole blood samples from 10 individuals with DR and 10 normal controls (NCs) for transcriptome sequencing. Following quality control and preprocessing of the sequencing data, differential expression analysis was conducted to identify differentially expressed genes (DEGs) between the DR and NC groups. Candidate genes were then selected by intersecting these DEGs with key module genes identified through weighted gene co-expression network analysis. These candidate genes were subjected to mendelian randomization (MR) analysis, then least absolute shrinkage and selection operator analysis to pinpoint key genes. The diagnostic utility of these key genes was evaluated using receiver operating characteristic curve analysis, and their expression levels were examined. Additional analysis, including nomogram construction, gene set enrichment analysis, drug prediction and molecular docking, were performed to investigate the functions and molecular mechanisms of the key genes. Finally, the expression of key genes was verified by qRT-PCR and biomarkers were identified.

Results: Intersection of 1,004 DEGs with 1,038 key module genes yielded 291 candidate genes. Five key genes were identified: HOXA1, GBP2, P2RY12, CLIC3 and PSAP. Among them, CLIC3, GBP2, and P2RY12 were identified as risk factors for DR, while HOXA1 and PSAP were protective. These key genes demonstrated strong diagnostic performance for DR. With the exception of P2RY12, all other key genes exhibited down-regulation in the DR group. Furthermore, the nomogram incorporating multiple key genes demonstrated superior predictive capacity for DR compared to a single key genes. The identified key genes are involved in oxidative phosphorylation and ribosome functions. Drug predictions targeting P2RY12 suggested prasugrel, ticagrelor, and ticlopidine as potential options owing to their high binding affinity with this key genes. The qRT-PCR results revealed that the results of GBP2, PSAP and P2RY12 exhibited consistent expression patterns with the dataset.

Conclusion: This study identified GBP2, P2RY12 and PSAP as NETs-associated biomarkers in the development of PDR, offering new insights for clinical diagnosis and potential treatment strategies for DR.

1 Introduction

Diabetic retinopathy (DR) is the top microvascular complication of diabetes mellitus, affecting approximately 30% to 40% of diabetic patients (1, 2). Globally, DR afflicts over 100 million people, making it the leading cause of irreversible blindness in working-age adults (1, 3). The key characteristics of DR include capillary occlusion, retinal ischemia, and vascular leakage (4). Current treatment options for DR include intravitreal injection of anti-vascular endothelial growth factor (VEGF), retinal laser photocoagulation and vitrectomy surgery (5). Although these therapeutic interventions have been deemed safe and effective through clinical research, they only provide temporary control of DR progression. Moreover, retinal capillary damage is irreversible and abovementioned therapies carry the risk of elevated intraocular pressure and retinal ischemia (6–8). Therefore, it is a priority to explore and identify the underlying pathogenesis of DR, which may provide more effective disease management.

Under the stimulation of inflammatory mediators, neutrophils release chromatin, histones and neutrophil granule proteins, forming a network structure outside the cell that binds and kills pathogens. This network structure is called neutrophil extracellular traps (NETs) (9). Recent reports have shown elevated serum NETs levels in patients with proliferative diabetic retinopathy (PDR) (10, 11). Other evidence suggests that neutrophil infiltration accelerates the adhesion of neutrophil to endothelial cells and promotes inflammation. NETs produced by neutrophils play a role in removing aging-damaged blood vessels and promoting the regeneration of retinal vasculature (12). Neutrophil aggregation near neovascularization has been detected in retinal slices of human PDR patients, and typical NETs have been observed in vitreous and retinal tissues. These findings indicate that NETs participate in the process of DR (10, 13). However, it is unclear which NETs-related genes regulate the pathogenesis of DR.

To gain deeper insights into the potential role that NETs-related genes played in DR, we employed a statistical method known as Mendelian randomization (MR). Observational epidemiological studies are prone to reverse causation, confounding and various biases, which limits their ability to provide conclusive evidence (14, 15). MR analysis utilizes genetic variants as exposure indicators that are not influenced by conventional study designs (14). This method effectively eliminates the confounding effects of extraneous factors and enhances causal inference, contributing to a better understanding and prevention of adverse effects on human health caused by modifiable exposures (15, 16). Despite its proven utility, the application of MR in DR research remains limited. In this study, we utilized MR analysis to investigate the causal associations between NETs-related genes and DR.

To this end, we first identified differentially expressed genes (DEGs) between the DR and normal control (NC) groups through transcriptome sequencing and bioinformatics analysis. Subsequently, we examined and certified the associations between these identified markers and PDR. In our study, we revealed the significance of five NETs-related genes in the pathogenesis of PDR and provided novel strategies for DR control.

2 Materials and methods

This study was carried out in the Second Hospital of Shandong University and adhered to principles of the Declaration of Helsinki involving human subjects. The study was approved by the Ethics Committee of Scientific Research of the Second Hospital of Shandong University (Approval number: KYLL2024043). Signed consent was obtained from each participant.

2.1 Sample collection and RNA-sequencing analysis

Whole blood samples were collected from 10 patients with DR and 10 NC patients. The DR group consisted of patients who were diagnosed as PDR with vitreous hemorrhage whereas without retinal detachment or fibrous formation, with a mean age of 61.50 ± 7.13 years. It included 7 males and 3 females, all of whom had type 2 diabetes. All of the patients with diabetes were diagnosed according to the criteria of the American Diabetes Association (17). PDR was defined as the presence of neovascularization or fibrous proliferation of the disc or elsewhere on the retina (18). Patients with the following conditions were excluded: (1) severe systemic diseases such as metabolic syndrome (excluding type 1 diabetes), ongoing infection or autoimmune diseases and malignant tumors; and (2) any other ocular disease, such as glaucoma, high myopia, retinal diseases, or a history of previous ocular surgery (except mild cataract) (19).

The NC group contained patients who underwent routine cataract surgery at the hospital. The patients had a mean age of 59.20 ± 10.83 years and included 6 males and 4 females. The patients would be excluded if they were found to suffer from ocular diseases other than age-related cataract, as well as participants with diabetes, hyperthyroidism, hypertension, immune-related diseases, and cancer. The specific sample information can be found in Supplementary Table S1.

The manufacturer’s protocol was followed to extract total RNA using TRIzol (Invitrogen, CA, USA). The quantity and quality of total RNA were assessed using NanoDrop ND-1000 (NanoDrop, Wilmington, DE, USA) to ensure its integrity. Subsequently, the assessment of RNA integrity was conducted utilizing a Bioanalyzer 2100 (Agilent, CA, USA) and verified through agarose electrophoresis. The samples were considered acceptable if they met the following criteria: concentration > 50 ng/μL, RNA integrity number (RIN) > 7.0, optical density (OD) 260/280 > 1.8, and total RNA > 1μg. Based on the guidelines provided by the manufacturer, the library was prepared for Illumina sequencing using Hieff NGS Ultima Dual-mode mRNA Library Prep Kit. Importantly, these libraries were sequenced on the Illumina Novaseq 6000 platform in a pair-end (PE) 150 mode. The overview of this study design was shown in Supplementary Figure S1.

2.3 Preprocessing of sequencing data

The quality of the sequencing data was evaluated using FastQC software. The original data were preprocessed using Trimmomatic to filter out low-quality data and remove pollution and joint sequences, resulting in clean data. The filtered clean sequencing data were then compared to the reference genome using the hisat2 comparison tool with default parameters. After the analysis, the gene expression matrix was obtained for subsequent analysis.

2.4 Identification and enrichment analysis of DEGs

The DESeq2 package (v 1.36.0) (20) was utilized to compare gene expression levels between the DR and NC groups. DEGs were identified with screening conditions of p < 0.05 and |log2FoldChange(FC)| > 0.5. Volcano plots were generated using ggplot2 (v 3.4.0) (21) to visualize DEGs, heatmaps were generated using ComplexHeatmap (v 2.12.1) (22) to visualize the expression patterns of the top 10 up- and downregulated DEGs. Enrichment analyses based on the Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) databases were performed using the clusterProfiler package (v 4.4.4) (23) and the org.Hs.eg.db (v 3.15.0) (24) human gene annotation package (p adj < 0.05). The GO system comprised three components: biological process (BP), molecular function (MF) and cellular component (CC). The treemap package (v 2.4-3) was used to generate a tree diagram illustrating the enrichment results.

2.5 Screening of NETs related module genes

A total of 257 NETs-related genes were obtained from published literature (25, 26). NETs scores were calculated for each sample using the GSVA package (v 1.38.2) (27), the differences in NETs scores between DR and NC groups were compared using the Wilcoxon test. The results were visualized using gghalves (v 0.1.4) to plot cloud and rain maps. Based on all gene expression data, NETs scores were used as trait data, a weighted gene co-expression network analysis (WGCNA) was performed using WGCNA package (v 1.72-1) (28) to screen for the modules and their genes that were most relevant to NETs scores.

Cluster analysis was implemented on the samples to determine the necessity of removing outlier samples for accurate subsequent analyses. The optimal soft threshold (β) was determined to maximize the adherence of gene interactions to the scale-free distribution. According to the standard of the hybrid dynamic tree cutting algorithm, a threshold of 200 genes per gene module was established in order to facilitate clustering of genes into distinct modules. The correlation between each module and the NET score was analyzed with the criteria of |cor| > 0.3 and p < 0.05 to identify gene modules associated with the NET score trait. The module exhibiting the highest relevance to the NET scores was defined as the key module and included key module genes.

2.6 The identification and functional investigation of candidate genes

The ggVennDiagram package (v 1.2.2) (29) was utilized to intersect DEGs and key module genes to identify differentially expressed NETs-related genes, which were recorded as candidate genes. To further investigate the functions of the candidate genes, functional enrichment analysis was implemented via the Metascape database (https://metascape.org/gp/index.html#/main/step1). On the other hand, a protein-protein interaction (PPI) network was constructed based on the STRING (https://string-db.org) website so as to explore whether there were interactions between candidate genes (confidence = 0.4). The network was visualized by Cytoscape software (v 3.9.1) (30).

2.7 Screening of instrumental variables in MR

MR studies need to comply with three basic principles: (1) there is a durable and significant correlation between IVs and exposure factors; (2) IVs are not associated with confounding factors; (3) IVs can only influence outcome variables through exposure factors. The DR-related Genome-wide association study (GWAS) dataset (trait ID: finn-b-DM RETINOPATHY EXMORE) was provided by Integrative Epidemiology Unit (IEU) Open GWAS (https://gwas.mrcieu.ac.uk/) database. This dataset included 190,594 European samples (cases: controls = 14,584: 176,010) and 16,380,347 single nucleotide polymorphisms (SNPs).

In MR analysis, each candidate gene was treated as an exposure factor, DR was treated as the result. The exposure factors were examined and IVs were filtered using the extract instruments function in the TwoSampleMR package (version 0.5.6) (31). SNPs significantly related to exposure factors (p < 5×10^-8) were searched, while SNPs showing linkage disequilibrium (LD) (clump=TRUE, r² = 0.001, kb=10000) were excluded. The directionality test function in the TwoSampleMR package was employed to test the directionality of the exposure factors. Furthermore, the F statistic was calculated to evaluate the strength of each IV. If the F statistic was greater than 10, it was considered that there was no weak instrumental bias, indicating a strong predictive potential of the IV for the outcome.

2.8 MR analysis between candidate genes and DR

The SNPs of exposure factors and outcomes were unified by the function harmonise data. MR analysis was then performed using the mr function combined with five algorithms: MR Egger (32), weighted median (33), inverse variance weighting (IVW) (34), simple mode (35) and weighted mode (36). The IVW method was deemed the most crucial among the five methods and its outcomes were conclusive. A p value less than 0.05 indicated a causal relationship between candidate genes and DR. A b value greater than 0 represented a risk factor, while a b value less than 0 indicated a protective factor. The MR results were visualized using scatter plots, funnel plots and forest plots.

Sensitivity analysis was conducted to assess the reliability of the MR results. Initially, the Cochran’s Q test was employed to evaluate the presence of heterogeneity in the tests, with a p value greater than 0.05 indicating no significant heterogeneity. Subsequently, a horizontal pleiotropy test was conducted to determine the presence of confounding factors (p > 0.05). Eventually, we employed a Leave-One-Out (LOO) analysis to systematically remove each SNP and assess whether the IVW method was affected by any particular individual SNPs. The genes with potential causal relationship with DR in MR analysis were defined as candidate key genes for subsequent analysis.

2.9 Access to key genes

The candidate key genes identified from MR analysis were included in the least absolute shrinkage and selection operator (LASSO) regression analysis for further screening. Specially, the LASSO regression analysis was implemented via glmnet package (v 4.1-6) (37). The genes identified as key genes had the lowest cross-validation error rate. Additionally, the diagnostic ability of each key genes was assessed by plotting the receiver operating characteristic (ROC) curve using pROC package (v 1.17.0.1) (38), larger values of area under the curve (AUC) representing more accurate diagnostic ability of the key genes for DR. Furthermore, the expression of key genes was assessed between DR and NC groups, the results were visualized using ggpubr (v 0.4.0) (39) to generate a box plot. Spearman correlation analysis was performed among the key genes and the results were demonstrated using circle plots created with the circlize package (v 0.4.15) (40).

2.10 Construction and evaluation of nomogram

The key genes obtained above were employed to construct a nomogram through rms package (v 6.2-0) (41) to facilitate the clinical judgment of the risk rate of DR. The nomogram assigned a score to each key genes, with each score corresponding to a specific key genes. Subsequently, the risk rate of DR was predicted based on the cumulative score. A higher score indicated an increased risk rate. In addition, ROC curves, decision curves and clinical impact curves were plotted to evaluate the predictive value of nomogram. Specially, the decision curve analysis (DCA) was generated via rmda package (v 1.6).

2.11 Gene set enrichment analysis

To further explore the pathways and functions associated with the key genes, GSEA was implemented based on the KEGG database via clusterProfiler package (v 4.4.4). First, the correlation coefficients between the key genes and the remaining genes in the dataset were determined. The genes were then sorted based on these correlation coefficients. Subsequently, GSEA was conducted according to the sorting results, with a significance threshold of p adj < 0.05). Additionally, the GSEA ridge plot was drawn using GseaVis package (v 0.0.5) to visualize the enrichment results.

2.12 Drug prediction and molecular docking

The Drug-Gene Interaction Database (DGIdb) (https://dgidb.genome.wustl.edu/) was utilized to predict key genes-related drugs with potential therapeutic effects in DR. Drugs reported in the literature were selected for molecular docking analysis. The 3D Conformer structures of the drugs was retrieved from PubChem database (https://pubchem.ncbi.nlm.nih.gov/). Meanwhile, the crystal structure of the key genes with Protein Data Bank (PDB) ID 4ntj was obtained from the PDB database (https://www1.rcsb.org/). Subsequently, molecular docking was performed using AutoDock vina and visualized by PyMol software (v 2.5).

2.13 Quantitative real-time polymerase chain reaction

To further confirm the results of the public database analysis, we collected five paired DR and NC whole blood samples and implemented RNA isolation and quantitative real-time polymerase chain reaction (qRT-PCR). Total RNA of 10 samples was separated by the TRIzol (Ambion, Austin, USA) based on the manufacturer’s guidance. The inverse transcription of total RNA into cDNA was implemented by using the SureScript-First-strand-cDNA-synthesis-kit (Servicebio, Wuhan, China) based on the producer’s indication. Then, qPCR was carried out utilizing the 2xUniversal Blue SYBR Green qPCR Master Mix (Servicebio, Wuhan, China) under the direction of the manual. The primer sequences for PCR were tabulated in Table 1. The expression was uniformized to the internal reference GAPDH and computed employing the 2^−ΔΔCt method. Finally, key genes with qRT-PCR results consistent with the dataset were identified as biomarkers.

Table 1

Table 1. Primer information.

2.14 Statistical analysis

The data were processed and analyzed using R software (version 4.2.1). The nonparametric Wilcoxon test was employed to assess differences between different groups, with a p value less than 0.05 considered to statistical significance.

3 Results

3.1 A total of 1,004 DEGs were associated with inflammation-related signaling pathways in DR

Differential expression analysis generated 1,004 DEGs between DR and NC groups. Among these DEGs, 580 were up-regulated in DR samples, while 424 were downregulated, in addition, heatmaps were drawn to visualise the expression patterns of up- and down-regulated TOP10 differential genes. (Figures 1A, B). Subsequently, enrichment analysis was conducted to probe the signaling pathways involved in these DEGs. The results revealed that these DEGs were enriched in 741 GO entries (including 590 BPs, 72 CCs and 79MFs) and 26 KEGG pathways. Specially, these enriched GO entries included cell activation involved in immune response, regulation of inflammatory response, leukocyte activation involved in immune response, etc (Figure 1C). Meanwhile, KEGG enrichment results showed that these DEGs were enriched mainly in NOD-like receptor signaling pathway, IL-17 signaling pathway and others (Figure 1D).

Figure 1

Figure 1. Differential expression analysis between the DR and NC groups. (A) Volcano plot displaying differentially expressed genes (DEGs) between DR patients and healthy controls for combined expression profiling. Pink nodes indicate upregulated DEGs; blue nodes indicate downregulated DEGs; gray nodes indicate genes that are not significantly expressed. (B) Heatmap plot of the top 10 DEGs. Red indicates DR samples, black indicates normal control samples, pink indicates high gene expression, and blue indicates low gene expression. (C) GO enrichment plot. (D) KEGG enrichment plot.

3.2 The identification of 1,038 key module genes was accomplished

The results of the GSVA demonstrated a significant difference in NETs score between DR and NC groups. The NETs score was significantly lower in the DR group than in the NC group, indicating that NETs have an impact on the occurrence and development of DR (Figure 2A). Therefore, the NETs score could be used as a trait to find the key module genes related to it through WGCNA. Cluster analysis revealed that there were no outliers in the samples (Figure 2B). When the ordinate scale-free R² crossed the threshold of 0.85 (red line), the first soft threshold β was determined to be 5, and the mean connectivity also tended to 0, indicating that the network approximated the scale-free distribution at this time. Therefore, the best soft threshold β was selected as 5 (Figure 2C). Furthermore, 16 modules were obtained (Figure 2D). After merging with a similarity of 0.5, 13 modules were ultimately obtained (Figure 2E). The results of correlation analysis demonstrated that MEgreenyellow (cor = 0.8, p = 3× 10^-5) and MEpink (cor = 0.48, p = 0.03) modules were strongly correlated with NETs scores (Figure 2F). Therefore, these two modules were defined as key modules, containing 1,038 key module genes. The scatter plot of the correlation between key module genes and NETs score traits was shown in Figure 2G.

Figure 2

Figure 2. Screening of NETs related module genes. (A) Cloud and rain map. Red indicates control samples, and blue indicates normal DR samples. The NETs were significantly lower in the DR group. (B) Sample dendrogram and trait heatmap plot. Cluster analysis revealed that there were no outliers in the samples. (C) Selection of the best soft threshold. (D) Clustering tree based on the module eigengenes of the modules. (E) Hierarchical cluster dendrogram of identified genes. Each color represents a module, and the vertical line represents a gene. (F) Heatmap of the correlations between module traits and the NET score. Each color represents a module, and the box contains the correlation coefficient and p value (p value in parentheses). Blue represents a negative correlation, and red represents a positive correlation. (G) The scatter plot of the correlation between key module genes and NETs score traits.

3.3 The 291 candidate genes were mainly involved in inflammation and immune-related signaling pathways in DR

By overlapping 1,004 DEGs and 1,038 key module genes, 291 candidate genes were identified (Figure 3A). Enrichment analysis revealed that these candidate genes were mainly enriched in inflammatory and immune-related signaling pathways, such as inflammatory response, innate immune response, cytokine signaling in immune system, and regulation of leukocyte activation (Figure 3B). Moreover, we explored the links between these pathways and found that each function may have the same genes and functions interact with each other (Figure 3C). Based on a confidence level of 0.4, a PPI network consisting of 222 nodes and 820 edges was constructed. Notably, C5AR1, YROBP, FCGR3A, CYBB and other candidate genes exhibited interactions with multiple genes (Figure 3D).

Figure 3

Figure 3. Identification and functional investigation of candidate genes. (A) Venn diagram. (B) Enrichment bar chart of 291 candidate genes. (C) Metascape enrichment interaction network. (D) Protein-protein interaction network.

3.4 The nomogram based on key genes exhibited superior diagnostic efficacy for DR

The MR analysis results demonstrated a causal relationship between 11 exposure factors and DR (IVW p < 0.05) (Table 2). Among them, CLIC3 (b = 0.327, p = 5.32E-05), GBP2 (b = 0.303, p = 0.013), HLA-A (b = 0.250, p = 0.004), SPI1 (b = 0.109, p = 0.001), P2RY12 (b = 0.085, p = 0.010) and OASL (b = 0.054, p = 0.005) were identified as risk factors for DR, whereas DENND3 (b = -0.065, p = 0.045), HLA-F (b = -0.150, p = 0.015), SLC22A15 (b = -0.156, p = 0.022), HOXA1 (b = -0.172, p = 0.012) and PSAP (b = -0.189, p = 1.35E-07) were found to be protective factors for DR (b < 0, p < 0.05). Supplementary Table S2 provides the F statistics for these 11 candidate key genes. These 11 candidate key genes were subsequently subjected to LASSO analysis for screening purposes. When lambda_min was 0.01331, five genes were screened by LASSO regression analysis, which were HOXA1, GBP2, P2RY12, CLIC3 and PSAP (Figure 4A). These genes were defined as key genes for subsequent analysis. The ROC curve demonstrated that the AUC values of these five key genes exceeded 0.7, indicating their robust diagnostic potential for DR (Figure 4B). Furthermore, a comparison of the expression levels of these five key genes between the DR and NC groups revealed a significant disparity. With the exception of P2RY12, all other key genes exhibited significantly downregulated expression in the DR group (Figure 4C). The correlation analysis indicated a highly significant negative correlation between GBP2 and P2RY12 (cor = -0.641, p = 0.002), while a significant positive correlation was observed between GBP2 and PSAP (cor = 0.660, p = 0.002) (Figure 4D; Supplementary Table S3).

Table 2

Table 2. MR analysis between 11 exposure factors and DR.

Figure 4

Figure 4. Access to biomarkers. (A) LASSO regression analysis. The horizontal axis deviation represents the proportion of residuals explained by the model, showing the relationship between the number of feature genes and the proportion of residuals explained (dev). The vertical axis represents the coefficient of genes (left); the horizontal axis represents log (Lambda), and the vertical axis represents the error of cross-validation (right). In practical analysis, we hope to find the position with the smallest cross-validation error. In the right graph, the dotted line on the left is the position with the smallest cross-validation error. Based on this position (lambda. min), the corresponding horizontal axis log (Lambda) is determined. The number of feature genes is displayed above, and the optimal log (Lambda) value is found. The corresponding gene and its coefficient are shown in the left graph, as well as the proportion of residuals explained by the model. (B) The ROC curve of these five biomarkers. The AUC values of these five biomarkers exceeded 0.7. (C) Boxplot of the expression of these five biomarkers. P2RY12 exhibited upregulation, and all other biomarkers exhibited downregulation in the DR group. (D) Chord diagram of the five identified biomarkers. There was a highly significant negative correlation between GBP2 and P2RY12 (cor = -0.641, p = 0.002), while a significant positive correlation was observed between GBP2 and PSAP (cor = 0.660, p = 0.002). (E) The nomogram of the five identified biomarkers. (F) The nomogram ROC curve. (G) Decision curve analysis. (H) Clinical impact curve analysis.

Based on these five identified key genes, a nomogram was constructed (Figure 4E). The AUC value of the ROC curve was 1, indicating that the nomogram outperformed gene prediction alone (Figure 4F). Meanwhile, the decision curve analysis demonstrated that the nomogram model exhibited significant benefits within the high-risk threshold range of 0-1, surpassing the clinical utility of the HOXA1, GBP2, P2RY12, CLIC3 and PSAP curves (Figure 4G). Furthermore, the clinical impact curve further substantiated the superior predictive capability of the nomogram model (Figure 4H).

3.5 Notable causal relationship between key genes and DR

MR analysis revealed that CLIC3 (b = 0.327, p = 5.32E-05), GBP2 (b = 0.303, p = 0.013), and P2RY12 (b = 0.085, p = 0.010) were risk factors associated with DR, while HOXA1 (b = -0.172, p = 0.012) and PSAP (b = -0.189, p = 1.35E-07) exhibited protective effects against the development of DR. The scatter plot showed that CLIC3, GBP2 and P2RY12 had positive slopes, indicating an increased risk for DR. Conversely, HOXA1 and PSAP exhibited a negative slope indicating their potential protective role against DR (Figure 5A). The forest plot demonstrated that CLIC3, GBP2 and P2RY12 increased the risk of DR, whereas HOXA1 and PSAP decreased the risk of DR (Figure 5B). The funnel plot further reflected that the MR analysis conformed to Mendel ‘s second law (Figure 5C). Sensitivity analysis was employed to assess the robustness of the MR results. The heterogeneity test revealed that CLIC3 and GBP2 exhibited significant heterogeneity with DR (the p-value of Cochran’s Q test in IVW results was less than 0.05). However, it is important to note that this observed heterogeneity did not impact the established causal relationships (Table 3). The p-values of the five key genes all exceeded 0.05, indicating the absence of horizontal pleiotropy between key genes and DR (Table 4). Furthermore, the LOO analysis showed no significant alteration in the results upon exclusion of each SNP, suggesting the absence of any substantial bias point (Figure 5D). In summary, these findings implied a substantial causal association between key genes and DR.

Figure 5

Figure 5. MR analysis between biomarkers and DR. (A) The scatter plot of the MR analysis. CLIC3, GBP2, and P2RY12 had positive slopes, whereas HOXA1 and PSAP had negative slopes. (B) The forest plot of the MR analysis. CLIC3, GBP2, and P2RY12 increased the risk of DR, whereas HOXA1 and PSAP decreased the risk of DR. (C) The funnel plot of the MR analysis. (D) LOO analysis. There was no significant alteration in the results upon exclusion of each SNP.

Table 3

Table 3. Heterogeneity test of the five identified biomarkers.

Table 4

Table 4. Horizontal pleiotropy test of the five identified biomarkers.

3.6 Oxidative phosphorylation and ribosome were functional pathways related to the identified key genes

GSEA was conducted to further explore the pathways associated with the key genes and their functions. The five key genes were associated with oxidative phosphorylation and ribosome functioning pathways. These pathways were also implicated in Parkinson Disease, Diabetic Cardiomyopathy, Herpes Simplex Virus 1 Infection, etc. (Figures 6A-E).

Figure 6

Figure 6. GSEA ridge plot of five identified biomarkers. (A) HOXA1. (B) GBP2. (C) P2RY12. (D) CLIC3. (E) PSAP.

3.7 Identification of biomarkers by qRT-PCR

The qRT-PCR results revealed that among the five key genes, GBP2, PSAP and P2RY12 whose expression results were consistent with the expression pattern of the dataset were identified as biomarkers. Specifically, GBP2 and PSAP were significantly downregulated in the DR group, while P2RY12 was significantly upregulated. The expression trends of HOXA1 and CLIC3 aligned with the dataset, however, no significant difference was observed between the DR and NC groups (Figure 7).

Figure 7

Figure 7. The relative mRNA expression of 5 biomarkers in the DR and control groups determined by qRT-PCR. The mRNA levels of GBP2 (p=0.0342) and PSAP (p=0.0419) were down-regulated, P2RY12 (p=0.0163) was up-regulated *p<0.05.

3.8 Stable binding between biomarkers and drugs

Only the drugs corresponding to P2RY12 were retrieved from DGIdb. Among them, prasugrel (42, 43), ticagrelor (44) and ticlopidine (42, 45) had associated references. Consequently, molecular docking of these three drugs with P2RY12 was conducted (Figure 8). These results indicated that there was a robust binding ability between biomarkers and predicted drugs, which might provide potential drug targets for the treatment of DR.

Figure 8

Figure 8. Molecular docking plot of the three drugs with P2RY12. (A) Molecular docking plot of 4 ntj with PRASUGREL. ARG-122, LYS-303, and SER-304 exhibited hydrogen bond interactions with PRASUGREL molecules. The docking affinity between this active molecule and P2RY12 was -6.8 kcal/mol. (B) The molecular docking plot of 4 ntj with TICAGRELOR. SER-83 and LYS-80 engaged in hydrogen bonding interactions with TICAGRELOR molecules. The docking affinity between this active molecule and P2RY12 was measured to be -8.2 kcal/mol. (C) The molecular docking plot of 4 ntj with TICLOPIDINE. TYR-109, ASN-159, ARG-256, and other residues exhibited hydrogen bond interactions with the TICLOPIDINE molecule. The docking affinity between this active molecule and P2RY12 was -6.6 kcal/mol.

4 Discussion

In the initial phase of this investigation, we identified the five most interesting NETs-related genes through MR and LASSO analysis. The diagnostic performance of these five NETs-related genes was excellent for PDR, as demonstrated by ROC curve analysis. CLIC3, GBP2 and P2RY12 were identified as risk factors for PDR, while HOXA1 and PSAP were found to have protective functions. P2RY12 was upregulated in the DR group, while the other key genes were downregulated. Three of key genes, GBP2, P2RY12 and PSAP, were validated in PDR patients by qPCR analysis. These genes were enriched in pathways such as oxidative phosphorylation and ribosomes, confirming the role of NETs in DR. Upon successful validation, these NETs-associated biomarkers could serve as valuable diagnostic tools and provide guidance for the development of novel therapeutic strategies.

In this study, GBP2 was downregulated in PDR and was identified as a potential risk factor. GBP2 is a member of the guanylate-binding protein (GBP) family. Previous studies have reported that GBP2 can regulate the release of a large number of pro-inflammatory factors such as interleukin (IL) 1β and IL-18 (46). The expression of GBP2 was markedly downregulated in the retina of mice with oxygen-induced retinopathy, which was associated with pathological retinal angiogenesis. Overexpression of GBP2 significantly inhibited neovascularization (47). These findings suggest a potential link between GBP2 levels and DR susceptibility. But the roles of GBPs in DR processes are not completely understood. This research opens a new way to explore GBP2 as a biomarker of PDR and implicated mechanism may be related to oxidative phosphorylation.

In our research, P2RY12 was identified as another potential risk factor for PDR. P2RY12 belongs to the family of G-protein-coupled receptors (48). It is used to regulate platelet activation and aggregation. The over-expression of P2RY12 in the diabetes mellitus lead to platelet aggregation and atherosclerosis development (49). Alterations in platelet that favor thrombosis occur early in the diabetic state and contribute to microvascular disease (50). In recent studies, P2RY12 can also regulate retinal microglia in retinal tissue, which initiate inflammation by releasing proinflammatory cytokines, reactive oxygen species (ROS), and reactive nitrogen species (RNS) (51).

PSAP is a lysosomal regulatory protein. It was downregulated and was identified as a potential protective factor in this study. It was previously reported that silencing of PSAP expression suppressed glycolysis as well as oxidative phosphorylation. Inhibition of PSAP may led to a reduction in atherosclerosis development and in plaque inflammation (52). It was also discovered that knockdown of PSAP strongly inhibited death receptor 6-induced apoptosis (53). But the roles of PSAP in PDR processes are rarely reported.

Our verification test also showed the down-regulated expression of CLIC3 and HOXA1 in DR group though without statistical significance. CLIC3 was a risk factors for PDR, while HOXA1 was found to have protective function. CLIC3 is a member of chloride intracellular channel (CLIC) protein family and may be related to tumor invasion (54). Knockdown of CLIC3 in human gastric cancer cells significantly accelerates cell proliferation. Decreased expression of CLIC3 in gastric cancer may result in poor prognosis of the patients (55). In diabetes mellitus, CLIC3 is significantly downregulated with a result of increased clinical complications (56). HOXA1 is a member of the homeobox (HOX) gene family. The misexpression of HOXA1 in differentiated cells could turn it into an oncogene to participate in cancer development (57). HOXA1 promotes apoptosis, inflammation and phosphorylated NF-κB p65 levels (57, 58). The non-significant results of CLIC3 and HOXA1 may be partly due to the small sample size in the present study.

Several reports have demonstrated that hyperglycemia can trigger sterile inflammation and activate neutrophils (59). The release of NETs by neutrophils can cause endothelial cell damage and vascular remodeling (13). NETs formation can be triggered by various stimuli, including microorganisms, urate crystals, autoantibodies, lipopolysaccharides, ROS, nitric oxide, proinflammatory cytokines, as well as interactions between neutrophils and activated platelets or endothelial cells (9). These findings suggest a reciprocal interaction between sterile inflammation, NETs formation and endothelial dysfunction ultimately leading to a vicious cycle.

There have been reports that the expression of GBP2 and P2RY12 is correlated with neutrophil infiltration (60, 61). The mechanisms that trigger neutrophil infiltration are poorly understood but may result from the release of pro-inflammatory cytokines, such as tumor necrosis factor-α (TNF-α), IL-1β, IL-6, and IL-8 (9). According to previous studies, both GBP2 and P2RY12 can regulate the release of these pro-inflammatory factors (46, 51). PSAP was found to regulate mitogen-activated protein kinase (MAPK), phosphatidylinositol 3’-kinase/AKT serine/threonine kinase 1 (PI3K/Akt) and transforming growth factor-β (TGF-β) pathways (62). These pathways are strongly associated with NETs (63–65). Our research found these biomarkers were enriched in oxidative phosphorylation and ribosomes, suggesting that the upregulation of P2RY2 and the downregulation of GBP2 and PSAP may also regulate NETs by promoting oxidative stress and inflammatory reaction in PDR. All these findings indicate that the NETs related biomarkers play an important role in the pathogenesis of DR.

Furthermore, we predicted drugs targeting P2RY12, among which prasugrel (42, 43), ticagrelor (44) and ticlopidine (42, 45) showed enhanced binding affinity for this biomarker. These newer antiplatelet drugs inhibit platelet aggregation by blocking the P2Y12 platelet receptor. Recent studies have shown that ticagrelor significantly decreases the levels of proinflammatory cytokines, such as TNF-α and IL-6, while increasing the number of circulating endothelial progenitor cells, thereby improving vascular endothelial function (66). Ticlopidine has also been found to increase nitric oxide production in human neutrophils (67). As a result, we predicted that ticagrelor and ticlopidine can be used as a targeted drugs for the treatment of DR. Further studies are needed to elucidate the specific regulatory mechanisms underlying the relationships among P2RY12, the targeted drugs and DR.

There are limitations of the study: Using blood-derived genes may not completely reflect changes in DR-specific tissues like the retina. And the insufficient amount of data and the lack of survival analysis on adequate clinical samples may lead to results bias. Especially, this may be the reason that the expression of HOXA1 and CLIC3 did not reach statistical significance. In the future, we will reach out to more patients and medical institutions to expand the sample size, focusing on an in-depth exploration of specific gene regulation mechanisms in eye tissues to provide a more detailed analysis of their complex regulatory pathways, thereby improving the representativeness and comprehensiveness of the research data. Simultaneously, we are actively exploring collaboration opportunities with external experts and research institutions to gain additional technical support and data analysis resources to further enhance our research capacity.

This study was conducted based on transcriptome self sequencing method, three biomarkers (GBP2, P2RY12 and PSAP) associated with NETs were screened for strong diagnostic value in PDR. In addition, MR analysis revealed that GBP2 and P2RY12 are risk factors for DR and PSAP a protective factor, which may provide a theoretical reference for future studies. However, future longitudinal studies are required to establish whether NETs are formed prior to clinically detectable DR and to determine the diagnostic significance of NETs-related biomarkers for the early detection and stratification of DR patients.

Data availability statement

The data presented in the study are deposited in the NCBI SRA repository, accession number PRJNA1171113.

Ethics statement

The studies involving humans were approved by the Institutional Review Board of the Second Hospital of Shandong University (Approval number: KYLL2024043). The studies were conducted in accordance with the local legislation and institutional requirements. The participants provided their written informed consent to participate in this study.

Author contributions

LH: Formal Analysis, Funding acquisition, Resources, Writing – original draft, Writing – review & editing. SW: Formal Analysis, Investigation, Writing – original draft, Writing – review & editing. LZ: Writing – review & editing. JH: Writing – review & editing. YZ: Writing – review & editing. XQ: Conceptualization, Funding acquisition, Supervision, Writing – original draft, Writing – review & editing.

Funding

The author(s) declare that financial support was received for the research, authorship, and/or publication of this article. This work was supported in part by the Natural Science Foundation of Shandong Province (grant no. ZR2023MH195) and the Horizontal Project of Shandong University (grant no. 26010212002212). The funders had no role in the study design, data collection and analysis, decision to publish, or preparation of the manuscript.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2024.1408974/full#supplementary-material

References

1. Yau JWY, Rogers SL, Kawasaki R, Lamoureux EL, Kowalski JW, Bek T, et al. Global prevalence and major risk factors of diabetic retinopathy. Diabetes Care. (2012) 35:556–64. doi: 10.2337/dc11-1909

PubMed Abstract | Crossref Full Text | Google Scholar

2. Ruta LM, Magliano DJ, LeMesurier R, Taylor HR, Zimmet PZ, Shaw JE. Prevalence of diabetic retinopathy in Type 2 diabetes in developing and developed countries. Diabetic Med. (2013) 30:387–98. doi: 10.1111/dme.2013.30.issue-4

PubMed Abstract | Crossref Full Text | Google Scholar

3. Ting DSW, Cheung GCM, Wong TY. Diabetic retinopathy: global prevalence, major risk factors, screening practices and public health challenges: a review. Clin Exp Ophthalmol. (2016) 44:260–77. doi: 10.1111/ceo.2016.44.issue-4

PubMed Abstract | Crossref Full Text | Google Scholar

4. Clapp C, Thebault S, Jeziorski MC, Martínez de la Escalera G. Peptide hormone regulation of angiogenesis. Physiol Rev. (2009) 89:1177–215. doi: 10.1152/physrev.00024.2009

PubMed Abstract | Crossref Full Text | Google Scholar

5. Cheung CMG, Lai TYY, Ruamviboonsuk P, Chen SJ, Chen Y, Freund KB, et al. Polypoidal choroidal vasculopathy: definition, pathogenesis, diagnosis, and management. Ophthalmology. (2018) 125:708–24. doi: 10.1016/j.ophtha.2017.11.019

PubMed Abstract | Crossref Full Text | Google Scholar

6. Arrigo A, Aragona E, Bandello F. VEGF-targeting drugs for the treatment of retinal neovascularization in diabetic retinopathy. Ann Med. (2022) 54:1089–111. doi: 10.1080/07853890.2022.2064541

PubMed Abstract | Crossref Full Text | Google Scholar

7. Wang W, Lo A. Diabetic retinopathy: pathophysiology and treatments. Int J Mol Sci. (2018) 19:1816. doi: 10.3390/ijms19061816

PubMed Abstract | Crossref Full Text | Google Scholar

8. Distefano LN, Garcia-Arumi J, Martinez-Castillo V, Boixadera A. Combination of anti-VEGF and laser photocoagulation for diabetic macular edema: A review. J Ophthalmol. (2017) 2017:1–7. doi: 10.1155/2017/2407037

Crossref Full Text | Google Scholar

9. Kaplan M, Radic M. Neutrophil extracellular traps: double-edged swords of innate immunity. J Immunol (Baltimore Md: 1950). (2012) 189:2689–95. doi: 10.4049/jimmunol.1201719

Crossref Full Text | Google Scholar

10. Wang L, Zhou X, Yin Y, Mai Y, Wang D, Zhang X. Hyperglycemia induces neutrophil extracellular traps formation through an NADPH oxidase-dependent pathway in diabetic retinopathy. Front Immunol. (2019) 9:3076. doi: 10.3389/fimmu.2018.03076

PubMed Abstract | Crossref Full Text | Google Scholar

11. Venkatesh P, Sheemar A, Soni D, Takkar B, Basu S. Inflammatory mediators in diabetic retinopathy: Deriving clinicopathological correlations for potential targeted therapy. Indian J Ophthalmol. (2021) 69:3035–49. doi: 10.4103/ijo.IJO_1326_21

PubMed Abstract | Crossref Full Text | Google Scholar

12. Gomułka K, Ruta M. The role of inflammation and therapeutic concepts in diabetic retinopathy—A short review. Int J Mol Sci. (2023) 24:1024. doi: 10.3390/ijms24021024

PubMed Abstract | Crossref Full Text | Google Scholar

13. Binet F, Cagnone G, Crespo-Garcia S, Hata M, Neault M, Dejda A, et al. Neutrophil extracellular traps target senescent vasculature for tissue remodeling in retinopathy. Science. (2020) 369:eaay5356. doi: 10.1126/science.aay5356

PubMed Abstract | Crossref Full Text | Google Scholar

14. Davey Smith G, Hemani G. Mendelian randomization: genetic anchors for causal inference in epidemiological studies. Hum Mol Genet. (2014) 23:R89–98. doi: 10.1093/hmg/ddu328

PubMed Abstract | Crossref Full Text | Google Scholar

15. Davey Smith G, Ebrahim S. [amp]]lsquo;Mendelian randomization’: can genetic epidemiology contribute to understanding environmental determinants of disease?*. Int J Epidemiol. (2003) 32:1–22. doi: 10.1093/ije/dyg070

PubMed Abstract | Crossref Full Text | Google Scholar

16. Zou X, Ye S, Tan Y. Potential disease biomarkers for diabetic retinopathy identified through Mendelian randomization analysis. Front Endocrinol. (2024) 14:1339374. doi: 10.3389/fendo.2023.1339374

Crossref Full Text | Google Scholar

17. Chamberlain JJ, Rhinehart AS, Shaefer CF, Neuman A. Diagnosis and management of diabetes: synopsis of the 2016 american diabetes association standards of medical care in diabetes. Ann Internal Med. (2016) 164:542–52. doi: 10.7326/M15-3016

Crossref Full Text | Google Scholar

18. Wilkinson CP, Ferris FL, Klein RE, Lee PP, Agardh CD, Davis M, et al. Proposed international clinical diabetic retinopathy and diabetic macular edema disease severity scales. Ophthalmology. (2003) 110:1677–82. doi: 10.1016/S0161-6420(03)00475-5

PubMed Abstract | Crossref Full Text | Google Scholar

19. Wu Z, Liu B, Ma Y, Chen H, Wu J, Wang J. Discovery and validation of hsa_circ_0001953 as a potential biomarker for proliferative diabetic retinopathy in human blood. Acta Ophthalmol. (2020) 99:306–13. doi: 10.1111/aos.14585

PubMed Abstract | Crossref Full Text | Google Scholar

20. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. (2014) 15:550. doi: 10.1186/s13059-014-0550-8

PubMed Abstract | Crossref Full Text | Google Scholar

21. Gustavsson EK, Zhang D, Reynolds RH, Garcia-Ruiz S, Ryten M, Mathelier A. ggtranscript: an R package for the visualization and interpretation of transcript isoforms usingggplot2. Bioinformatics. (2022) 38:3844–6. doi: 10.1093/bioinformatics/btac409

PubMed Abstract | Crossref Full Text | Google Scholar

22. Gu Z, Eils R, Schlesner M. Complex heatmaps reveal patterns and correlations in multidimensional genomic data. Bioinformatics. (2016) 32:2847–9. doi: 10.1093/bioinformatics/btw313

PubMed Abstract | Crossref Full Text | Google Scholar

23. Yu G, Wang L-G, Han Y, He Q-Y. clusterProfiler: an R package for comparing biological themes among gene clusters. OMICS: A J Integr Biol. (2012) 16:284–7. doi: 10.1089/omi.2011.0118

Crossref Full Text | Google Scholar

24. Wang L, Wang D, Yang L, Zeng X, Zhang Q, Liu G, et al. Cuproptosis related genes associated with Jab1 shapes tumor microenvironment and pharmacological profile in nasopharyngeal carcinoma. Front Immunol. (2022) 13:989286. doi: 10.3389/fimmu.2022.989286

PubMed Abstract | Crossref Full Text | Google Scholar

25. Jiang T, Wang Y, Chen X, Xia W, Xue S, Gu L, et al. Neutrophil extracellular traps (NETs)-related lncRNAs signature for predicting prognosis and the immune microenvironment in breast cancer. Front Cell Dev Biol. (2023) 11:1117637. doi: 10.3389/fcell.2023.1117637

PubMed Abstract | Crossref Full Text | Google Scholar

26. Wu J, Zhang F, Zheng X, Zhang J, Cao P, Sun Z, et al. Identification of renal ischemia reperfusion injury subtypes and predictive strategies for delayed graft function and graft survival based on neutrophil extracellular trap-related genes. Front Immunol. (2022) 13:1047367. doi: 10.3389/fimmu.2022.1047367

PubMed Abstract | Crossref Full Text | Google Scholar

27. Hänzelmann S, Castelo R, Guinney J. GSVA: gene set variation analysis for microarray and RNA-Seq data. BMC Bioinf. (2013) 14:7. doi: 10.1186/1471-2105-14-7

Crossref Full Text | Google Scholar

28. Langfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinf. (2008) 9:559. doi: 10.1186/1471-2105-9-559

Crossref Full Text | Google Scholar

29. Gao C-H, Yu G, Cai P. ggVennDiagram: an intuitive, easy-to-use, and highly customizable R package to generate venn diagram. Front Genet. (2021) 12:706907. doi: 10.3389/fgene.2021.706907

PubMed Abstract | Crossref Full Text | Google Scholar

30. Shannon P, Markiel A, Ozier O, Baliga NS, Wang JT, Ramage D, et al. Cytoscape: A software environment for integrated models of biomolecular interaction networks. Genome Res. (2003) 13:2498–504. doi: 10.1101/gr.1239303

PubMed Abstract | Crossref Full Text | Google Scholar

31. Huang P, Zhang PF, Li Q. Causal relationship between cannabis use and cancer: a genetically informed perspective. J Cancer Res Clin Oncol. (2023) 149:8631–8. doi: 10.1007/s00432-023-04807-x

PubMed Abstract | Crossref Full Text | Google Scholar

32. Bowden J, Davey Smith G, Burgess S. Mendelian randomization with invalid instruments: effect estimation and bias detection through Egger regression. Int J Epidemiol. (2015) 44:512–25. doi: 10.1093/ije/dyv080

PubMed Abstract | Crossref Full Text | Google Scholar

33. Bowden J, Davey Smith G, Haycock PC, Burgess S. Consistent estimation in mendelian randomization with some invalid instruments using a weighted median estimator. Genet Epidemiol. (2016) 40:304–14. doi: 10.1002/gepi.2016.40.issue-4

PubMed Abstract | Crossref Full Text | Google Scholar

34. Burgess S, Scott RA, Timpson NJ, Davey Smith G, Thompson SG. Using published data in Mendelian randomization: a blueprint for efficient identification of causal risk factors. Eur J Epidemiol. (2015) 30:543–52. doi: 10.1007/s10654-015-0011-z

PubMed Abstract | Crossref Full Text | Google Scholar

35. Hemani G, Zheng J, Elsworth B, Wade KH, Haberland V, Baird D, et al. The MR-Base platform supports systematic causal inference across the human phenome. eLife. (2018) 7:e34408. doi: 10.7554/eLife.34408

PubMed Abstract | Crossref Full Text | Google Scholar

36. Hartwig FP, Davey Smith G, Bowden J. Robust inference in summary data Mendelian randomization via the zero modal pleiotropy assumption. Int J Epidemiol. (2017) 46:1985–98. doi: 10.1093/ije/dyx102

PubMed Abstract | Crossref Full Text | Google Scholar

37. Li Y, Lu F, Yin Y. Applying logistic LASSO regression for the diagnosis of atypical Crohn’s disease. Sci Rep. (2022) 12:11340. doi: 10.1038/s41598-022-15609-5

PubMed Abstract | Crossref Full Text | Google Scholar

38. Robin X, Turck N, Hainard A, Tiberti N, Lisacek F, Sanchez J-C, et al. pROC: an open-source package for R and S+ to analyze and compare ROC curves. BMC Bioinf. (2011) 12:77. doi: 10.1186/1471-2105-12-77

Crossref Full Text | Google Scholar

39. Cheng Q, Chen X, Wu H, Du Y. Three hematologic/immune system-specific expressed genes are considered as the potential biomarkers for the diagnosis of early rheumatoid arthritis through bioinformatics analysis. J Trans Med. (2021) 19:18. doi: 10.1186/s12967-020-02689-y

Crossref Full Text | Google Scholar

40. Mehrotra R, Shen G, Wang W-L. Circlize package in R and Analytic Hierarchy Process (AHP): Contribution values of ABCDE and AGL6 genes in the context of floral organ development. PloS One. (2022) 17:e0261232. doi: 10.1371/journal.pone.0261232

PubMed Abstract | Crossref Full Text | Google Scholar

41. Sachs MC. plotROC: A tool for plotting ROC curves. J Stat Software. (2017) 79:2. doi: 10.18637/jss.v079.c02

Crossref Full Text | Google Scholar

42. Sangkuhl K, Shuldiner AR, Klein TE, Altman RB. Platelet aggregation pathway. Pharmacogenet Genomics. (2011) 21:516–21. doi: 10.1097/FPC.0b013e3283406323

PubMed Abstract | Crossref Full Text | Google Scholar

43. Girard-Gagnepain A, Amirache F, Costa C, Lévy C, Frecha C, Fusil F, et al. Baboon envelope pseudotyped LVs outperform VSV-G-LVs for gene transfer into early-cytokine-stimulated and resting HSCs. Blood. (2014) 124:1221–31. doi: 10.1182/blood-2014-02-558163

PubMed Abstract | Crossref Full Text | Google Scholar

44. Chan NC, Eikelboom JW, Ginsberg JS, Lauw MN, Vanassche T, Weitz JI, et al. Role of phenotypic and genetic testing in managing clopidogrel therapy. Blood. (2014) 124:689–99. doi: 10.1182/blood-2014-01-512723

PubMed Abstract | Crossref Full Text | Google Scholar

45. Ravelich SR, Shelling AN, Wells DN, Peterson AJ, Lee RSF, Ramachandran A, et al. Expression of TGF-β1, TGF-β2, TGF-β3 and the receptors TGF-βRI and TGF-βRII in placentomes of artificially inseminated and nuclear transfer derived bovine pregnancies. Placenta. (2006) 27:307–16. doi: 10.1016/j.placenta.2005.03.002

PubMed Abstract | Crossref Full Text | Google Scholar

46. Wang Y, Gao W, Shi X, Ding J, Liu W, He H, et al. Chemotherapy drugs induce pyroptosis through caspase-3 cleavage of a gasdermin. Nature. (2017) 547:99–103. doi: 10.1038/nature22393

PubMed Abstract | Crossref Full Text | Google Scholar

47. Xu X, Ding X, Wang Z, Ye S, Xu J, Liang Z, et al. GBP2 inhibits pathological angiogenesis in the retina via the AKT/mTOR/VEGFA axis. Microvascular Res. (2024) 154:104689. doi: 10.1016/j.mvr.2024.104689

Crossref Full Text | Google Scholar

48. Gurbel PA, Kuliopulos A, Tantry US. G-protein-coupled receptors signaling pathways in new antiplatelet drug development. Arteriosclerosis Thrombosis Vasc Biol. (2015) 35:500–12. doi: 10.1161/ATVBAHA.114.303412

Crossref Full Text | Google Scholar

49. Taghizadeh M, Kargarfard M, Braune S, Jung F, Naderi M. Long-term aerobic exercise training in type two diabetic patients alters the expression of miRNA-223 and its corresponding target, the P2RY12 receptor, attenuating platelet function. Clin Hemorheol Microcirculation. (2022) 80:107–16. doi: 10.3233/CH-211209

Crossref Full Text | Google Scholar

50. Colwell JA, Winocour PD, Halushka PV. Do platelets have anything to do with diabetic microvascular disease? Diabetes. (1983) 32(Suppl 2):14–9. doi: 10.2337/diab.32.2.s14

PubMed Abstract | Crossref Full Text | Google Scholar

51. Wang L, Qian Y, Che X, Jiang J, Suo J, Wang Z. Isolation and characterization of primary retinal microglia from the human post-mortem eyes for future studies of ocular diseases. Front Cell Neurosci. (2022) 15:786020. doi: 10.3389/fncel.2021.786020

PubMed Abstract | Crossref Full Text | Google Scholar

52. van Leent MMT, Beldman TJ, Toner YC, Lameijer MA, Rother N, Bekkering S, et al. Prosaposin mediates inflammation in atherosclerosis. Sci Trans Med. (2021) 13:eabe1433. doi: 10.1126/scitranslmed.abe1433

Crossref Full Text | Google Scholar

53. Zhang J, Zhao ZJ, Fu X, Niu H, Hu C, Dong Y, et al. Proapoptotic mitochondrial carrier homolog protein PSAP mediates death receptor 6 induced apoptosis. J Alzheimer’s Dis. (2020) 74:1097–106. doi: 10.3233/JAD-191086

Crossref Full Text | Google Scholar

54. Norman J, Zanivan S. Chloride intracellular channel 3: A secreted pro-invasive oxidoreductase. Cell Cycle. (2017) 16:1993–4. doi: 10.1080/15384101.2017.1377031

PubMed Abstract | Crossref Full Text | Google Scholar

55. Kawai S, Fujii T, Shimizu T, Sukegawa K, Hashimoto I, Okumura T, et al. Pathophysiological properties of CLIC3 chloride channel in human gastric cancer cells. J Physiol Sci. (2020) 70:15. doi: 10.1186/s12576-020-00740-7

PubMed Abstract | Crossref Full Text | Google Scholar

56. Wang Y, He X, Zheng D, He Q, Sun L, Jin J, et al. Integration of metabolomics and transcriptomics reveals major metabolic pathways and potential biomarkers involved in pulmonary tuberculosis and pulmonary tuberculosis-complicated diabetes. Microbiol Spectr. (2023) 11:e0057723. doi: 10.1128/spectrum.00577-23

PubMed Abstract | Crossref Full Text | Google Scholar

57. Han Z, Hu H, Yin M, Lin Y, Yan Y, Han P, et al. HOXA1 participates in VSMC-to-macrophage-like cell transformation via regulation of NF-κB p65 and KLF4: a potential mechanism of atherosclerosis pathogenesis. Mol Med. (2023) 29:104. doi: 10.1186/s10020-023-00685-8

PubMed Abstract | Crossref Full Text | Google Scholar

58. Zhang X, Tang X, Pan L, Li Y, Li J, Li C. Elevated lncRNA-UCA1 upregulates EZH2 to promote inflammatory response in sepsis-induced pneumonia via inhibiting HOXA1. Carcinogenesis. (2022) 43:371–81. doi: 10.1093/carcin/bgac004

PubMed Abstract | Crossref Full Text | Google Scholar

59. Forrester JV, Kuffova L, Delibegovic M. The role of inflammation in diabetic retinopathy. Front Immunol. (2020) 11:583687. doi: 10.3389/fimmu.2020.583687

PubMed Abstract | Crossref Full Text | Google Scholar

60. Zhang S, Chen K, Zhao Z, Zhang X, Xu L, Liu T, et al. Lower expression of GBP2 associated with less immune cell infiltration and poor prognosis in skin cutaneous melanoma (SKCM). J Immunother. (2022) 45:274–83. doi: 10.1097/CJI.0000000000000421

PubMed Abstract | Crossref Full Text | Google Scholar

61. Liverani E, Rico MC, Tsygankov AY, Kilpatrick LE, Kunapuli SP. P2Y12 receptor modulates sepsis-induced inflammation. Arteriosclerosis Thrombosis Vasc Biol. (2016) 36:961–71. doi: 10.1161/ATVBAHA.116.307401

Crossref Full Text | Google Scholar

62. Jiang Y, Zhou J, Hou D, Luo P, Gao H, Ma Y, et al. Prosaposin is a biomarker of mesenchymal glioblastoma and regulates mesenchymal transition through the TGF-β1/Smad signaling pathway. J Pathol. (2019) 249:26–38. doi: 10.1002/path.v249.1

PubMed Abstract | Crossref Full Text | Google Scholar

63. Chen L, Liu Y, Wang Z, Zhang L, Xu Y, Li Y, et al. Mesenchymal stem cell-derived extracellular vesicles protect against abdominal aortic aneurysm formation by inhibiting NET-induced ferroptosis. Exp Mol Med. (2023) 55:939–51. doi: 10.1038/s12276-023-00986-2

PubMed Abstract | Crossref Full Text | Google Scholar

64. Zhou X, Yang L, Fan X, Zhao X, Chang N, Yang L, et al. Neutrophil chemotaxis and NETosis in murine chronic liver injury via cannabinoid receptor 1/Gαi/o/ROS/p38 MAPK signaling pathway. Cells. (2020) 9:373. doi: 10.3390/cells9020373

PubMed Abstract | Crossref Full Text | Google Scholar

65. Adel RM, Helal H, Ahmed Fouad M, Sobhy Abd-Elhalem S. Regulation of miRNA-155-5p ameliorates NETosis in pulmonary fibrosis rat model via inhibiting its target cytokines IL-1β, TNF-α and TGF-β1. Int Immunopharmacol. (2024) 127:111456. doi: 10.1016/j.intimp.2023.111456

PubMed Abstract | Crossref Full Text | Google Scholar

66. Jeong HS, Hong SJ, Cho SA, Kim JH, Cho JY, Lee SH, et al. Comparison of ticagrelor versus prasugrel for inflammation, vascular function, and circulating endothelial progenitor cells in diabetic patients with non-ST-segment elevation acute coronary syndrome requiring coronary stenting: A prospective, randomized, crossover trial. JACC Cardiovasc Interventions. (2017) 10:1646–58. doi: 10.1016/j.jcin.2017.05.064

Crossref Full Text | Google Scholar

67. De La Cruz P, Arrebola M, González-Correa A, Martinez-Cerdán E, Moreno A, de la Cuesta FS. Effects of clopidogrel and ticlopidine on experimental diabetic ischemic retinopathy in rats. Naunyn-Schmiedeberg’s Arch Pharmacol. (2003) 367:204–10. doi: 10.1007/s00210-002-0657-4

Crossref Full Text | Google Scholar