Screening of potential key ferroptosis-related genes in Chronic Obstructive Pulmonary Disease

Introduction

Chronic obstructive pulmonary disease (COPD) is a common respiratory system disease that often results in distal bronchiectasis, accompanied by tracheal wall degradation, with a high incidence rate and high mortality characteristics.^1–3 Characterized by persistent respiratory symptoms and progressive airflow blockage.³ Cigarette smoke (CS) exposure can cause an elevation in reactive oxygen species (ROS), leading to oxidative stress. This can further disrupt iron metabolism, resulting in inflammation, cellular aging, and cell death.^4–8 Currently, COPD has become the third leading cause of death in the world. Due to the aging population and high smoking rates, the global burden of COPD is increasing. By 2030, the overall prevalence and mortality rate of COPD will increase by 8%.^9,10 At present, there is no effective targeted drug for clinical use, so there is an urgent need to explore promising drug targets.

Ferroptosis, first proposed by Brent R. Stockwell of Columbia University in 2012, is a programmed cell death induced by lipid peroxidation through an iron dependent pathway with unique morphological and biological characteristics.^11,12 Ferroptosis is characterized by increased ROS, iron dependent lipid peroxidation, and plasma membrane damage.¹³ Numerous studies have indicated the crucial role of ferroptosis in the development of acute lung injury,¹⁴ asthma,¹⁵ COPD,^16–18 Alzheimer’s disease¹⁹, and renal failure.¹⁵ Prolonged inhalation of particulate matter found in smoke, which includes iron particles, can contribute to the accumulation of ROS in the epithelial lining of the airways. This accumulation can cause damage to the plasma membrane, leading to the release of damage-related molecular patterns (DAMPs). Consequently, iron homeostasis is altered, resulting in excessive oxidative stress, inflammation, and eventually iron death. These processes can exacerbate the progression of emphysema and airway narrowing in individuals with COPD.^7,20 Glutathione peroxidase 4 (GPX4) is an enzyme that plays a significant role in regulating abnormal amino acid metabolism during iron death. It functions as a crucial mediator in clearing lipid free radicals and acts as a vital protective metabolite in preventing iron death. Upon exposure to smoke, GPX4 plays a negative regulatory role in preventing the accumulation of unstable iron and reducing lipid peroxidation in lung epithelial cells. Additionally, it aids in mitigating non-apoptotic cell death, thereby reducing oxidative stress and lung injury associated with COPD.^21,22 Nuclear Factor E2 Related Factor 2(Nrf2), is an important protein transcription factor that plays a crucial role in intracellular oxidative stress response and antioxidant defense. When cells are subjected to oxidative stress conditions, Nrf2 gradually migrates to the nucleus, binds to the promoter of the target gene, and promotes the production of a series of antioxidant enzymes, such as GPX, which can neutralize harmful oxygen radicals and eliminate toxic compounds, thereby protecting cells from oxidative damage.²³ However, the ferroptosis-related genes of COPD remain largely unknown and require further exploration. Exploring and revealing the potential ferroptosis-related genes of COPD will provide us with potential biomarkers for the treatment of COPD.

We analyzed the dataset GSE148004 in the GEO database to explore the differentially expressed ferroptosis related genes in COPD. Then, differentially expressed ferroptosis related genes were analyzed using protein-protein interaction (PPI), correlation analysis, gene ontology (GO) enrichment analysis, and Kyoto Encyclopedia of Genes and Genomics (KEGG) pathway enrichment analysis.

Materials and Methods

Microarray Dataset Collection and Data Process

Microarray datasets were screened from GEO (http://www.ncbi.nlm.nih.gov/geo/). The search keywords were “COPD” and “macrophage” with the following searching strategies: ((“COPD”[MeSH Terms] OR COPD[All Fields]) AND (“macrophage”[MeSH Terms] OR macrophage[All Fields])) AND “Homo sapiens”[porgn] AND (“gse”[Filter] AND “Expression profiling by array”[Filter]). Finally, GSE148004 (as training set) was obtained. GSE148004 contained 16 whole sputum-derived alveolar macrophage RNA samples from 7 patients with COPD and 9 normal samples.

We transformed the probe into gene symbol in dataset based on the platform’s annotation file, when there were multiple probes mapped to the same gene symbol; the mean value of probes was selected as the gene expression value. Differentially expressed genes (DEGs) between COPD and healthy were analyzed via the “limma package” in R software, with the following cutoff for adjustment: p value < 0.05 and FC (fold changes) > 1.0.

Screening Ferroptosis-Related Differentially Expressed Genes

A total of 369 ferroptosis-related genes were performed from the FerrDd Database (http://wwwzhounan.org/ferrdb).²⁴ The intersection of DEGs and ferroptosis-related genes was visualized by the Venn plot.

Protein–Protein Interaction(PPI) Network Analysis of Ferroptosis-Related DEGs

The protein–protein interaction (PPI) network analysis of ferroptosis-related DEGs was conducted using the STRING online website (https://string-db.org/).²⁵ This creature is limited to “humman”. Set a minimum valid constraint score of 0.4. The results of STRING were then imported into Cytoscape (version 3.9.1),²⁶ and the key subnetworks were extracted using the cytoHubba plugin. Then, correlation analysis of ferroptosis-related DEGs was identified using Spearman correlation in the R software “corrplot” package. Finally, the three genes with the highest scores (maximum correlation criterion, degree algorithm) in the subnetworks were selected as hub genes.

Functional Enrichment Analysis of Ferroptosis-Related DEGs

We used the “heatmap” and “ggplot2” packages of R to the draw the thermal and volcanic maps. And we used the “clusterProfiler” package of R to perform the Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses of ferroptosis-related DEGs. GO analysis included three categories, biological process (BP), cellular component (CC), and molecular function (MF), which was important in the exploration of biological functions. KEGG analysis was used to explore potential pathways.

Construction of Gene-miRNA and TF-Gene Network

In this study, TF-gene regulatory network and gene-miRNA interaction network was established by the NetworkAnalyst platform(https://www.networkanalyst.ca/).²⁷ The parameters used are as follows: Specify organism: H. sapiens (human); Set ID type: Official Gene Symbol; Gene miRNA interaction database: miRTarBase v8.0 database; and TF–gene interaction database: JASPAR database. The analysis results of NetworkAnalyst were then imported into Cytoscape software for further visualization.

Potential Drug Prediction

Potential drugs targeting hub genes were predicted by the Enrichr platform.²⁸ The access of the DSigDB database was acquired using Enrichr(https://amp.pharm.mssm.edu/Enrichr/).²⁹

Statistical Analysis

The statistical analyses were performed using R software (version 4.2.2).

Results

Ferroptosis-Related DEGs in COPD

First, the GSE148004 dataset was downloaded from GEO. We next analyzed the expression of 2238 DEGs in 7 COPD patients and 9 healthy individuals, which were presented using a volcanic map (Figure 1A); Then, we downloaded 369 ferroptosis related genes by using the FerrDd databases. A total of 41 ferroptosis-related genes were presented using the Wayne diagram (Figure 1B), including 19 down-regulated genes and 22 up-regulated genes (Table 1); Following the analysis of the GSE148004 dataset with R software, the 41 differentially expressed ferroptosis-related genes between COPD and normal groups were presented in heatmap (Figure 1C).

Table 1 The 41 Differentially Expressed Ferroptosis-Related Genes in COPD Samples Compared to Healthy Samples

Figure 1 The differentially expressed ferroptosis-related genes (DEGs) in COPD and healthy samples. (A) Volcano plot of GSE148004. The cutoff logFC (fold change) was −1/1, and P-value is 0.05. Red dots indicate upregulated genes; blue dots indicate downregulated genes; gray dots indicate non-statistically significant genes. (B) Venn diagram of the ferroptosis-related DEGs. (C) Heatmap of the 41 differentially expressed ferroptosis-related genes in COPD and healthy samples. Abbreviations: COPD, chronic obstructive pulmonary disease; DEGs, differentially expressed genes.

PPI Network and Correlation Analysis

To determine the interactions among ferroptosis-related DEGs, we performed PPI analysis. The results demonstrated that these ferroptosis-related genes interacted with each other (Figure 2A). and a visual analysis was performed using Cytoscape (Figure 2B). In order to explore the expression correlation of these ferroptosis-related DEGs, a correlation analysis was conducted (Figure 2C).

Figure 2 Protein–protein interactions (PPI) analysis the 41 differentially expressed ferroptosis-related genes. (A) The PPI among 41 differentially expressed ferroptosis-related genes. (B) The PPI network constructed with cytoscape. In the figure, nodes in crimson color represent 3 hub genes and yellow and Orange represent other common genes. (C) Spearman correlation analysis of the 41 differentially expressed ferroptosis-related genes. Abbreviation: PPI, Protein–protein interactions.

GO and KEGG Enrichment Analysis

We used R software to conduct GO and KEGG enrichment analysis to analyze the potential biological functions of these differentially expressed ferroptosis-related genes. The findings revealed that the most significant GO enriched terms in the inflammatory response, iron metabolism, and cellular response to hypoxia after trauma (biological process); mitochondrial outer membrane, organ outer membrane, outer membrane (cellular component); heme binding, tetrapyrrole binding, and L-amino acid transmembrane transporter activities regulate heme binding (molecular function). The outcomes were depicted using bubble, spin, and circle diagrams (Figure 3A–C). In KEGG enrichment analysis, differentially expressed ferroptosis-related genes are mainly involved in the process of central carbon metabolism in cancer, ferroptosis, and HIF-1 signaling pathways, and the results are presented using bubble charts (Figure 3D).

Figure 3 The enrichment of GO terms and KEGG pathways based on 41 ferroptosis-related DEGs in COPD. (A) GO enrichment analysis of the ferroptosis-related DEGs, including biological process(BP), cell composition(CC), and molecular function(MF). (B) GO analysis of the top 5 differentially expressed ferroptosis-related genes. (C) GO analysis of the top 10 differentially expressed ferroptosis-related genes. (D) Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis of 41 differentially expressed ferroptosis-related genes. Abbreviations: GO, Gene Ontology; BP, biological process; CC, cellular component; MF, molecular function; KEGG, Kyoto Encyclopedia of Genes and Genomes.

Predictive Results of miRNA and TF Targeting the Hub Genes

Through the Network analyst, the gene-miRNA networks of hub genes HIF1A, RRARG, and KRAS were predicted and shown in Supplementary Table 1. Based on the prediction results, a 163 node gene and miRNA co-expression network was constructed by Cytoscape and was visualized (Figure 4A), similarly, a 24 node gene and TF co-expression network was constructed by Cytoscape and was visualized (Figure 4B), which showed that has-mir-27a-3p could regulate HIF1A, PPARG, and KRAS expression; and the gene-TF network of hug genes HIF1A, PPARG, and KRAS was constructed, which identified GATA binding protein 2 (GATA2), NFKB transcription factor family member 1 (NFKB1), and transcription factor p65 (RELA) as modulators of PPARG and HIF1A; Histone H4 Transcription Factor (HINFP) and Forkhead box(FOXC1) could regulate PPARG and KRAS.

Figure 4 Predictive networks of microRNAs (miRNAs), transcription factors (TFs). (A) Target gene–miRNA network. (B) Target gene–TF network. Abbreviation: TFs, transcription factors.

Drug Prediction Targeting the Hub Genes

We use the DSigDB database to predict drugs associated with hub genes that may be used to treat COPD by regulating ferroptosis. Finally, a total of 380 targeted drugs were predicted, and the top 10 drugs predicted based on comprehensive scores are listed in Table 2. The top five drugs are simvastatin, N-acetylL-cysteine, actinomycin, diclofenac sodium, and phorbol-1-myristic acid-12-acetate (PMA).

Table 2 Top 10 Drugs with Higher Scores

Discussion

COPD is a major global health problem, and its development is influenced by smoking, inflammation, airway remodelling, and other factors. More and more evidence suggests that ferroptosis may be related to the pathogenesis of COPD. Smoking can lead to oxidative stress in lung cells, leading to ferroptosis in the lung epithelium. By analyzing different ferroptosis-related differential genes in COPD patients and healthy controls, this study has identified 41 ferroptosis-related DEGs in COPD, and explored the diagnostic value of three pivotal genes, HIF1A, PPARG, and KRAS in COPD. Then, the Network analyst database is used to predict the miRNA and TF of the hub gene, and the DSigDB database is used to predict the potential target drugs of the hub gene. These drugs can play an important role in improving the severity of COPD through the ferroptosis mechanism of macrophages.

Hypoxia inducible factor-1 (HIF-1) is an important transcriptional regulator of cellular responses to hypoxia, oxidants, and inflammation.³⁰ HIF1A is the gene of hypoxia inducible factor 1α(HIF1α) that plays important transcriptional regulatory roles in oxygen regulation, iron death, inflammation, and immune responses. Firstly, HIF1α is overexpressed in the lungs of COPD patients which can lead to adaptive changes in lung cells, such as angiogenesis, erythropoiesis, and cell proliferation, to adapt to hypoxic environments.^31–33 Other studies suggest that hydrogen sulfide can inhibit CS induced PHD2/HIF1 α/ MAPK signaling activation reduces the secretion of IL-6, and TNF-α, reduces oxidative stress, ROS loading, cells damage, alters iron homeostasis, thereby alleviating emphysema and improving lung function.^31,34 Secondly, under hypoxic conditions, HIF1α play a key role in iron metabolism by regulating the expression of iron-related proteins, such as DMT1, FPN1, Dcytb, and transferrin receptor (TfR), which may lead to abnormal accumulation of intracellular iron and the production of excessive oxygen free radicals, thereby exacerbating the progression of COPD.^35,36

Peroxidase proliferator-activated receptors (PPARs) are members of the steroid hormone receptor superfamily and are fatty acid sensors that regulate various aspects of lipid metabolism.³⁷ Currently, there are three subtypes of PPARs, including α, β/ δ and γ. PPAR α promotes oxidation of fatty acids in the liver,³⁸ PPAR δ is an important regulator of lipid oxidation in skeletal muscle,³⁹ PPAR γ has been found in lung tissue and alveolar macrophages (AMs),⁴⁰ which could combined with NFκB to inhibit the expression levels of inflammatory mediators, such as TNF-a and CC chemokine ligand(CCL)5, which overall reduced granulocyte aggregation, and lessened the degree of inflammation in the lungs of COPD patients.⁴¹ In addition, hydrogen sulfide can improve CS induced emphysema and airway inflammation by activating the Nrf2 and PPARY signaling pathways, maintaining iron metabolism and lipid metabolism balance.⁴²

Kirsten rat sarcoma virus oncogene (KRAS) is the most frequently mutated member of the Ras family in human tumors, usually associated with pancreatic cancer and lung cancer,⁴³ and one of the most important regulators of cell proliferation, differentiation and survival.⁴⁴ COPD is an independent risk factor for lung cancer. KRAS can also affect COPD and iron death through several pathways, but its specific mechanism is still unclear. Moghadamet used non-classifiable Haemophilus influenzae (NTHi) to induce COPD-like airway inflammation to study the effect of KRAS mutation expression in airway epithelium on induced lung cancer, and found that it has a promoting effect on lung cancer.⁴⁵ Another study showed that the high expression of HIF1α is also sufficient to enhance KRAS-induced lung neoplasia, while promoting tumor growth, increasing tumor vessel survival, and cell proliferation.⁴⁶

Two drugs with the highest scores selected by us through DsigDB are currently in clinical use. Studies have shown that Simvastatin enhances the anti-tumor effects of bevacizumab against lung adenocarcinoma A549 Cells by abating the HIF-1α-Wnt/β-Catenin signaling pathway.^47,48 However, research on simvastatin targeting PPARG and KRAS for the treatment of COPD is relatively scarce. Research has shown that the antioxidant N-acetylcysteine (NAC) has been identified as a novel inhibitor of the hypoxic response pathway. It can serve as a precursor to glutathione (GSH), increasing GSH content, clearing reactive oxygen species, and reducing the ferrous ions required for the oxidation of proline (PHDs), this weakens the stability of HIF-1α under hypoxic conditions.⁴⁹ NAC has been widely used in the treatment of COPD. Additionally, other studies have revealed that the imbalance between T helper (Th)17 cells and regulatory T (Treg) cells plays a vital role in the development and outcome of COPD.^50,51 NAC can downregulate the expression of HIF-1α, thus regulating the balance of Th1/Treg and reducing the release of the inflammatory cytokine IL-17, thereby improving lung function in COPD patients.⁵² Currently, clinical drugs for treating COPD mainly focus on the research direction of HIF1A, but there is relatively little research on PPARG and KRAS. Therefore, future research directions can focus on treatment strategies for PPARG and KRAS to gain a more comprehensive understanding and treatment of COPD.

We also evaluated the biological functions and pathways of ferroptosis-related DEGs, indicating that the HIF-1 signaling pathway is involved in the occurrence and development of COPD. However, this study has some limitations. First, the number of samples chosen is relatively small, and a larger sample size is required. Second, due to time and funding constraints, this study only employs bioinformatics methods. Changes in the mRNA or protein levels of DEGs can be detected as a supplement to this study using RT-PCR or Western blotting, and clinical samples or animal experiments can also be supplemented.

Conclusions

In the present study, 41 potential ferroptosis-related genes of COPD were identified through bioinformatics analysis. The Network analyst and Enrichr platform was used to predict the miRNAs and TF and drugs associated with the hub genes, indicating that these drugs may affect the development of COPD by regulating ferroptosis in macrophages. These results expand our understanding of COPD and may contribute to the treatment of COPD.