Identification of Hub Genes and Pathways Associated with Oxidative Stress of Cartilage in Osteonecrosis of Femoral Head Using Bioinformatics Analysis.

OBJECTIVE: This study aimed to identify the hub genes and pathways of genes related to oxidative stress of cartilage in osteonecrosis of femoral head (ONFH), and to predict the transcription factors of the hub genes.
METHODS: The GSE74089 was obtained from the Gene Expression Omnibus (GEO) database, including 4 necrotic tissues and 4 normal tissues, and the differentially expressed genes (DEGs) were identified by limma package in R language. Simultaneously, we searched for the genes related to oxidative stress in the Gene Ontology (GO) database. GO and signaling pathways analysis were performed using DAVID, Metascape, and GSEA. Protein-protein interaction (PPI) network was constructed using the STRING database, and the Degree algorithm of Cytoscape software was used to screen for hub genes. Finally, the NetworkAnalyst web tool was used to find the hub genes' transcriptional factors (TFs).
RESULTS: In total, 440 oxidative stress-related genes were found in GSE74089 and GO database, and 88 of them were significantly differentially expressed. These genes were mainly involved in several signaling pathways, such as MAPK signaling pathway, PI3K-AKT-mTOR signaling pathway, FOXO signaling pathway. The top 10 hub genes were JUN, FOXO3, CASP3, JAK2, RELA, EZH2, ABL1, PTGS2, FBXW7, MCL1. Besides, TFAP2A, GATA2, SP1, and E2F1 may be the key regulatory factors of hub genes.
CONCLUSIONS: We identified some hub genes and signaling pathways associated with oxidative stress in ONFH through a series of bioinformatics analyses.

Introduction

It is estimated that about 8.12 million people over the age of 15 in China suffer from osteonecrosis of femoral head (ONFH). ONFH is a disease with high morbidity and low cure rate, which poses a severe threat to patients’ quality of life and brings a heavy economic burden. Surgery is the primary therapeutic method. The pathological mechanisms of ONFH are very diverse, including abnormal lipid metabolism, microcirculatory ischemia, apoptosis, and so on. It is generally accepted that the decrease of blood perfusion leads to the weakening of the normal repair function of cells, resulting in irreversible tissue destruction eventually. Most research has focused on bone tissue; however, cartilage is also an essential factor in disease development. In the early stages of the disease, radiographic changes can be seen, including roughness of the cartilage surface. The degeneration of cartilage is progressively worse. Changes in mechanical stress caused by cartilage degeneration may also further aggravate the progression of the disease. Therefore, understanding the pathological changes on cartilage is vital to understand the condition. The mechanisms of cartilage degeneration require further research.

Many studies revealed that oxidative stress resulting from the imbalance of reactive oxygen species (ROS) and antioxidant capacity is related to pathologic processes such as cardiovascular diseases, diabetes, aging, and obesity. ROS are involved in the apoptosis of chondrocytes, activation of metalloproteinases, and reduced synthesis of matrix components in cartilage. ROS has been widely studied in osteoarthritis, and regulation of oxidative stress can reduce autophagy in cartilage.[8,9] However, there are few studies on oxidative stress of cartilage in ONFH.

In the genomic era, gene chip has been widely used to explore the mechanisms of diseases, which brings some new enlightenment to the pathogenesis from the gene level. Therefore, this study aims to find out the oxidative stress of cartilage in ONFH from the point of view of bioinformatics analysis and provide a reference for further research of ONFH.

Materials and Methods

Data Source

The gene expression profile of GSE74089 was downloaded from the Gene Expression Omnibus (GEO) database (http://www.ncbi.nlm.nih.gov/geo), which is a database repository of high-throughput gene expression data, hybridization arrays, chips, and microarrays. The data set GSE74089 was based on the GPL13497 platform (Agilent-026652 Whole Human Genome Microarray 4 x 44K v2) and consisted of 8 femoral cartilage samples, including 4 ONFH samples (3 males and 1 female) and 4 normal samples (3 males and 1 female). Due to the restricted sample size, the authors analyzed the male samples and the female samples together.

Data Preprocessing and Integration

Differentially expressed genes (DEGs) were screened using the limma package in R. The cutoff criteria for statistical significance were absolute log2 fold change (logFC) >1.5 and P value <.05. The DEGs’ heat map and volcano plot were drawn using the ggplot package on the R platform. Then, search for “oxidative stress” in the Gene Ontology (GO) database (http://geneontology.org/), and use “Homo sapiens” as a screening condition to collect genes related to oxidative stress. The intersection of the genes selected from the GO database and the genes from the GSE74089 was considered as oxidative stress-related genes (OS-genes) of cartilage in ONFH, where the differential genes were regarded as OS-DEGs.

GSEA

To better understand the biological mechanism associated with ONFH, we performed GSEA (Gene Set Enrichment Analysis) software (version 4.1.0) to explore the candidate molecular pathways of OS-genes. Hallmark and Canonical Pathways gene sets were downloaded from the Molecular Signatures Database (http://www.gsea-msigdb.org/gsea/msigdb/index.jsp). False discovery rate (FDR) <0.25 and P < .05 were regarded as the cutoff criteria.

GO and KEGG Analysis

DAVID (https://david.ncifcrf.gov/home.jsp) (version 6.8) is an online analysis tool suite with Integrated Discovery and Annotation function, which provides typical batch annotation and gene-GO term enrichment analysis to highlight the most relevant GO terms associated with related genes. In the GO analysis, the identified OS-DEGs were classified into 3 classes: Molecular Function (MF), Biological Process (BP), and Cellular Component (CC). Then, Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis was conducted online using Metascape (http://metascape.org) to predict the signaling pathways in which OS-DEGs may participate. Only terms with P values <.01 and several enriched genes ≥3 were considered as significant.

Protein-Protein Interaction Network Analysis

The protein-protein interaction (PPI) network was derived based on the Search Tool for the Retrieval of Interacting Genes (STRING) database (https://string-db.org), which covered almost all functional interactions between the expressed proteins, and interaction with a combined score >0.4 was considered statistically significant. The results of this analysis were visualized with Cytoscape (version 3.8.0) software. Through the Molecular Complex Detection (MCODE) plug-in of Cytoscape software, the most closely connected modules were selected from the PPI network for further analysis (setting parameters as degree cutoff = 2, node score = 0.2, k-core = 2, maximum depth = 100). Similarly, the key modules were analyzed by KEGG analysis. We used cytoHubba, a plug-in of Cytoscape software, to filter the hub genes from the whole PPI network and calculated it by the Degree method.

Prediction Transcription Factor

To further understand the regulation of these hub genes, we found out the transcriptional factors (TFs) of the hub genes and constructed the genes-TFs network map through NetworkAnalyst (https://www.networkanalyst.ca) using the JASPAR database.

The flow chart of the study is depicted in .

Figure 1.

The flow chart of this study. ONFH = osteonecrosis of femoral head; DEGs = differentially expressed genes; GSEA = Gene Set Enrichment Analysis; GO = Gene Ontology; KEGG = Kyoto Encyclopedia of Genes and Genomes; PPI = Protein-protein interaction; TFs = transcriptional factors.

Result

Screening of Candidate Genes

After analysis of the data set GSE74089, a total of 20,844 genes were found, including 3,626 DEGs. A total of 488 genes related to oxidative stress were selected from the GO database, of which 440 could be matched in the expression profile of GSE74089. The Venn diagram showed that 88 genes crossed between DEGs and oxidative stress-related genes, including 71 upregulated genes and 17 downregulated genes. The volcanic map of GSE74089 and the heat map of 88 OS-DEGs are shown in . The Venn diagram is shown in .

Figure 2.

(A) Volcano plot of differentially expressed genes. The red nodes represent upregulated DEGs with P value <.05 and logFC >1.5; the blue nodes represent downregulated DEGs with P value <.05 and logFC <−1.5. (B) Heat map and clustering analysis of 88 OS-DEGs. The marked genes are the hub genes obtained by subsequent analysis. DEGs = differentially expressed genes; OS = oxidative stress; ONFH = osteonecrosis of femoral head.

Figure 3.

The Venn diagram. DEGs: Differentially expressed genes screened from GSE74089. OS: Oxidative stress-related genes downloaded from Gene Ontology database.

The expression information of 440 OS-gene was uploaded to GSEA software, and the Hallmark and the KEGG gene set database were used to analyze genes at the overall level of the expression profile of OS-genes. The results showed that most of the upregulated genes were involved in oxidative phosphorylation, hypoxia, TGF-β signaling, PI3K-AKT-mTOR signaling, and ubiquitin-mediated proteolysis, as shown in .

Figure 4.

Gene set enrichment analysis of 440 OS-genes by the Hallmark gene set database (left) and the KEGG gene set database (right). FDR <0.25 and P value <.05 were regarded as the cutoff criteria. Most of the upregulated genes were involved in oxidative phosphorylation, hypoxia, TGF-β signaling, PI3K-AKT-mTOR signaling, and ubiquitin-mediated proteolysis. OS = oxidative stress; FDR = false discovery rate; NES = normalized enrichment score; NOFH = necrosis of the femoral head; KEGG = Kyoto Encyclopedia of Genes and Genomes; TGF = transforming growth factor; AKT = protein kinase B; mTOR = mammalian target of rapamycin.

To analyze the role of OS-DEGs in ONFH, the differential genes were analyzed by GO and KEGG. GO BP showed that they were mainly enriched in cellular response to hydrogen peroxide, cellular response to oxidative stress, response to oxidative stress. CC involved specifically includes cytoplasm, neuronal cell body. MF terms mainly contained protein binding, protein homodimerization activity, identical protein binding, and so on. KEGG pathway annotation showed that the targeted genes were highly enrolled in several pathways, such as pathways in cancer, MAPK signaling pathway, PI3K-AKT signaling pathway, and apoptosis ( ).

Figure 5.

(A) TOP10 KEGG pathway enrichment map of OS-DEGs (DEGs related to OS). (B) GO enrichment map with OS-DEGs. Only shown the top 5 enriched terms in BP, CC, MF. KEGG = Kyoto Encyclopedia of Genes and Genomes; OS = oxidative stress; DEGs = differentially expressed genes; GO = Gene Ontology; BP = Biological Process; CC = Cellular Component; MF = Molecular Function; GO = Gene Ontology; MAPK = mitogen-activated protein kinase; TNF = tumor necrosis factor.

PPI Network Analysis

The selected OS-DEGs were imported into the String database to construct the protein interaction network, and then Cytoscape was used to create the PPI network map ( ). The network diagram contains 68 nodes and 198 edges. The MCODE plug-in analyzes the critical modules of the PPI network diagram, and there are 2 modules with scores >4. KEGG analysis showed that the first module was mainly enriched in Alzheimer disease, mitophagy-animal, MAPK signaling pathway, IL-17 signaling pathway, FOXO signaling pathway, and the second module was mostly clustered in mitophagy-animal, PI3K-Akt signaling pathway, apoptosis, microRNAs in cancer pathway ( ). Using the cytoHubba plug-in to find hub nodes in the network graph, the top 10 genes were JUN, FOXO3, CASP3, JAK2, RELA, EZH2, ABL1, PTGS2, FBXW7, MCL1.

Figure 6.

PPI network diagram of OS-DEGs. The darker the green color, the higher the degree of nodes. The network diagram contains 68 nodes and 198 edges. The first 10 genes are selected as JUN, FOXO3, CASP3, JAK2, RELA, EZH2, ABL1, PTGS2, FBXW7, MCL1. PPI = protein-protein interaction; OS = oxidative stress; DEGs = differentially expressed genes.

Figure 7.

(A) Key modules 1 of PPI Network. (B) The chord diagram of modules 1 KEGG analysis. (C) Key modules 2 of PPI Network. (D) The chord diagram of modules 2 KEGG analysis. PPI = protein-protein interaction; KEGG = Kyoto Encyclopedia of Genes and Genomes; MAPK = mitogen-activated protein kinase.

Transcription Factors-Genes Pairs

We found 4 transcription factors with degree ≥5 and constructed the TFs-genes network diagram ( ). They were TFAP2A (degree: 6; betweenness: 93), GATA2 (degree: 5; betweenness: 96), SP1 (degree: 5; betweenness: 85), and E2F1 (degree: 5; betweenness: 85).

Figure 8.

Transcriptional factors-genes network diagram. The red circle represents the hub genes, and the blue quadrilateral represents the transcription factor. The 4 transcription factors are TFAP2A (degree: 6; betweenness: 93), GATA2 (degree: 5; betweenness: 96), SP1 (degree: 5; betweenness: 85), and E2F1 (degree: 5; betweenness: 85).

Discussion

Cartilage degeneration plays an essential role in the development of ONFH. ROS participate in chondrocyte apoptosis, activation of metalloproteinases, and reduced synthesis of matrix components. In this study, ONFH-related and oxidative stress-related genes were found from the GEO and GO databases, respectively. After analysis, 440 OS-genes were obtained, including 88 OS-DEGs.

GO results showed that most of these genes act in the cytoplasm and nucleoplasm. The identified proteins were mainly involved in response to oxidative stress and response to hydrogen peroxide. Hydrogen peroxide is a by-product of oxygen metabolism in organisms and plays an integral part in oxidative stress. This suggests that hydrogen peroxide is produced in cartilage. GSEA and KEGG analysis results showed that the target proteins were mainly involved in MAPK signaling pathway, PI3K-AKT-mTOR signaling pathway, FOXO signaling pathway, mitophagy-animal, TGF-β signaling pathway. Most of the genes participating in these pathways were upregulated. MAPK signaling pathways regulate a variety of biological processes through multiple cellular mechanisms.[11,12] Activation of MAPK signaling pathway increases the expression of MMP1 and MMP3 in chondrocytes. In other cartilage diseases, such as osteoarthritis, inhibition of the MAPK signaling pathway also has a protective effect on cartilage. PI3K-AKT-mTOR is an essential and complex signaling pathway. Some experiments have shown that inhibition of PI3K-AKT-mTOR is beneficial for cartilage regeneration.[15,16] FOXO transcription factors have been shown to protect chondrocytes by regulating autophagy and defending against oxidative stress. Knockdown of FOXO3 increases dexamethasone-induced apoptosis and ROS level in chondrocytes. Activation of the FOXO signaling pathway may be the result of anti-oxidative stress in chondrocytes. ROS are mainly produced in mitochondria. Damaged mitochondria can be cleared by mitochondrial phagocytosis, thus limiting the damage of ROS.

The top 10 hub genes were CASP3, JAK2, RELA, EZH2, ABL1, PTGS2, FBXW7, MCL1, JUN, and FOXO3. All of these genes were significantly upregulated except for RELA. To verify the expression of these genes in ONFH, relevant literature was searched and summarized in . These genes need to be further verified by experiments, especially in the cartilage in ONFH.

Table 1.

The Expression of Related Genes at RNA Level or Protein Level in ONFH by Literature Search.

	Express	RNA/Protein	Sample	Conclusion
PTGS2	Upregulation	RNA	Cartilage	PTGS2 is gradually upregulated as ONFH advances. 20
CASP3	Upregulation	RNA	Cartilage	GC induces gene expression related to apoptosis in chondrocytes. 21
RELA	Upregulation	Protein	Cartilage	The level of phospho-RELA expression in ONFH is higher than in femoral neck fracture. 5
	Upregulation	RNA	Bone tissue	Activation of NF-κB by TLR4 signaling plays an important role in the pathogenesis of necroptosis. 22
EZH2	Upregulation	RNA and protein	Bone tissue	Reducing the level of EZH2 expression could inhibit the development of ONFH. 23
FOXO3	Upregulation	RNA and protein	Bone marrow mesenchymal stem cell	Resveratrol can reverse the inhibition of Wnt signaling by FOXO. 24
JUN	Upregulation	Protein	Osteoblast	JNK/c-Jun signaling is involved in autophagy and apoptosis in ONFH. 25

GC = glucocorticoid; JNK = c-Jun N-terminal kinase; ONFH = osteonecrosis of femoral head; NF-κB = nuclear factor kappa B.

The protein encoded by CASP3 is a cysteine-aspartic acid protease that acts in the execution phase of cell apoptosis. GC induces the expression of apoptosis-related genes in chondrocytes in ONFH, such as Caspase-3, Caspase-9, ASK-1, and JNK-1. RELA is the essential subunit that mediates NF-κB signal transduction and participates in the metabolism of chondrocytes. RELA can induce anti-apoptotic genes to protect chondrocytes from apoptosis. Previous reports have shown that RELA is upregulated in bone tissue and cartilage in femoral head necrosis, and cartilage degeneration is related to the NF-κB signaling pathway, which requires further experimental verification.[5,22] Overexpression of EZH2 can accelerate cartilage destruction by enhancing the expression of inflammatory cytokines and MMP. Reducing the level of EZH2 expression could inhibit the development of ONFH. PTGS2 (also known as COX-2) significantly upregulates the expression of matrix metalloproteinase-1 (MMP1) in articular cartilage. Cox-2 is upregulated in inflamed joint tissue and leads to increased production of PGE2, which is associated with inflammation, apoptosis, and angiogenesis. Cox-2 inhibitors have a potential protective effect on cartilage. FBXW7 is the F-box protein subunit of an Skp1-Cul1-F-box (SCF)-type ubiquitin ligase complex that plays a central part in the degradation of oncoproteins such as c-Myc, c-Jun, Notch, and cyclin E. FBXW7 protects chondrocytes by inhibiting the HIF-1α/VEGF pathway and promoting the integration of collagen II, proteoglycans, and SOX-9. Glucocorticoid induces apoptosis and autophagy of osteoblasts in ONFH via the ROS/JNK/c-Jun signaling pathway. Overexpression of C-Jun inhibits the expression of type II collagen in articular chondrocytes. JAK2 encodes a non-receptor tyrosine kinase that acts in cytokine and growth factor signaling. The JAK2/STAT pathway is activated in chondrocytes in osteoarthritis, and blocking the JAK2/STAT3 pathway can inhibit MMP overexpression and collagen II degradation.[33,34] JAK2 may play a similar role in ONFH, which needs further experimental verification. ABL1 is a proto-oncogene that encodes a protein tyrosine kinase engaged in various cellular processes, including cell division, adhesion, differentiation, and stress response. ROS accelerate the development of disease by activating ABL1. The action of ABL1 in ONFH is unclear.

The TFs-genes network showed that TFAP2A, GATA2, SP1, and E2F1 play important roles in ONFH. Experiments have shown that hydrogen peroxide increases the methylation of TFAP2A (transcription factor AP-2 alpha) and activation of TFAP2A could increase the expression of genes related to oxidative stress.[36,37] The encoded protein of GATA2 acts in regulating the transcription of genes involved in the development and proliferation of hematopoietic and endocrine cell lineages. H2O2 increases GATA-2 binding activity and thus inhibits erythropoietin gene expression. The functions of SP1 are complex, including proliferation, differentiation, DNA damage response, apoptosis, and angiogenesis. The proliferation of chondrocytes can be regulated by inhibiting the activity of Sp1. Overexpression of E2F1 sensitizes cells to apoptosis. Meanwhile, E2F1 can increase the proliferation and migration of chondrocytes. The role of these genes in cartilage needs further study.

Conclusion

In summary, several genes and pathways diagnosed in the current study may be associated with the molecular mechanism of ONFH. These results contribute to our understanding of oxidative stress in ONFH. Unfortunately, this data set does not describe the etiology of ONFH. Most important of all, the results need to be verified by further large-scale research.