Identification of biomarkers associated with synovitis in rheumatoid arthritis by bioinformatics analyses

Abstract Objectives: Rheumatoid arthritis (RA) is the most common inflammatory arthritis in the world, but its underlying mechanism is still unclear. The present study aims to screen and verify the potential biomarkers of RA. Methods: We searched the Gene Expression Omnibus (GEO) database for synovial expression profiling from different RA microarray studies to perform a systematic analysis. Functional annotation of differentially expressed genes (DEGs) was conducted, including GO enrichment analysis and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis. The protein–protein interaction (PPI) networks of the DEGs were constructed based on data from the STRING database. The expression levels of the hub genes in normal membranes and RA synovium were detected by quantitative real-time polymerase chain reaction (qRT-PCR) and Western blot system. Results: A total of 444 differential expression genes were identified, including 172 up-regulated and 272 down-regulated genes in RA synovium compared with normal controls. The top ten hub genes; protein tyrosine phosphatase receptor type C (PTPRC), LCK proto-oncogene (LCK), cell division cycle 20 (CDC20), Jun proto-oncogene (JUN), cyclin-dependent kinase 1 (CDK1), kinesin family member 11 (KIF11), epidermal growth factor receptor (epidermal growth factor receptor (EGFR), vascular endothelial growth factor A (VEGFA), mitotic arrest deficient 2 like 1 (MAD2L1), and signal transducer and activator of transcription 1 (STAT1) were identified from the PPI network, and the expression level of VEGFA and EGFR was significantly increased in RA membranes (P<0.05). Conclusion: Our results indicate that the hub genes VEGFA and EGFR may have essential effects during the development of RA and can be used as potential biomarkers of RA.


Introduction
Rheumatoid arthritis (RA) is a chronic, aggressive, and progressive autoimmune disease, which can cause bone and cartilage damage as well as disability [1]. The main pathological features of RA are joint involvement, proliferative synovitis, destruction of bone and cartilage. The occurrence and development of RA involve many factors, which are related to various factors, such as genetic factors, infection, and immune dysfunction [2]. Although there is extensive research into the mechanisms and causes of RA formation and progression, the etiology of RA is still not clear [3]. and the experiment was repeated three times. Relative expression levels are calculated using the 2 − C T method, and statistical analysis was performed using SPSS software (version 22.0).

Western blot assay
The expression of vascular endothelial growth factor A (VEGFA) and epidermal growth factor receptor (EGFR) were detected by Western blotting. RA synovitis was lysed in IP cell lysis buffer on ice. The quantitative protein was denatured and transferred to the polyvinylidene difluoride membranes. The membranes were blocked with 5% bovine serum albumin (BSA) for 2 h at room temperature, the polyvinylidene fluoride (PVDF) membrane was incubated with the BSA first antibody overnight at 4 • C. After washing the membrane with Tris-buffered saline Tween (TBST), it was incubated with the secondary antibody for 2 h at room temperature. The β-actin was set as an internal control. The primary antibody of VEGFA (1:1000), EGFR (1:1000), and β-Actin (1:1000) used in the present study were purchased from Signalway Antibody LLC (Maryland, U.S.A.). Protein bands were visualized by Azure Biosystems C600 (Azure Biosystems, U.S.A.).

Case and control groups
Our study was approved by the Ethics Committee of the Second Hospital of Jilin University, Jilin, China. A total of ten patients with RA were included in the present study, including four males and six females aged 49.9 + − 4.21.
There were ten donors in the control group, including four males and six females, aged 47.8 + − 2.99. Informed consent was obtained from all enrolled groups. RA synovial tissue was obtained from patients undergoing total knee arthroplasty at The Second Hospital of Jilin University, and normal synovial tissue was taken from donors undergoing knee arthroscopy at The Second Hospital of Jilin University.

Statistical analysis
Each experiment was repeated at least three times, and the data were presented as the mean + − standard deviation. Comparisons between groups were performed using unpaired Student's t test. P-values <0.05 were considered statistically significant.

Identification of DEGs in RA
A total of 35 RA samples and 29 normal controls were included in the study. The raw data of GSE55457, GSE55235, and GSE12021 dataset were standardized, and the results are shown in Figure 1. The |log2(fold change)| > 1 and P-value <0.05 were used to evaluate significant differences in the expressions of mRNA between the two groups, we extracted 1331, 1817, and 1325 differential expression genes from the expression profile datasets GSE55457, GSE55235, and GSE12021, respectively. Integration of DEGs, a total of 444 differential expression genes were identified, including 172 up-regulated and 272 down-regulated genes in RA synovium compared with normal controls ( Figure 2).

Enrichment analysis of DEGs
To acquire the gene functions of DEGs, gene function enrichment analysis was analyzed by the DAVID database. These DEGs were classified into three groups as follows: BP, CC, and MF. The top five BP, CC, and MF terms are shown in Table 2 and Figure 3A. In the BPs group, the DEGs are mainly enriched in the immune response, B-cell receptor signaling pathway, inflammatory response, regulation of immune response, and chemokine-mediated signaling pathway. In the CC group, the DEGs are mainly enriched in the T-cell receptor complex, immunological synapse, external side of plasma membrane, integral component of plasma membrane, and plasma membrane. Moreover, in the MF group, the DEGs are mainly enriched in chemokine activity, transcriptional activator activity, antigen binding, immunoglobulin receptor binding, and protein homodimerization activity. KEGG analysis shows that DEGs were enriched in pathways such as Cytokine-cytokine receptor interaction Osteoclast differentiation, Chemokine signaling pathway, NF-κB signaling pathway, RA, MAPK signaling pathway, and Cell adhesion molecules. The results were shown in Figures 3B and 4, and Table 3.

PPI network analysis and hub genes screening
Based on the STRING online database, we constructed a PPI network with Cytoscape software, containing 222 nodes and 933 edges ( Figure 5). We used the degree analysis method of plug-in CytoHubba in Cytoscape to screen the hub genes. The results show that the top ten genes, which were mostly connected with other genes or proteins, including protein tyrosine phosphatase receptor type C (PTPRC), LCK proto-oncogene (LCK), cell division cycle 20 (CDC20), Jun proto-oncogene (JUN), cyclin-dependent kinase 1 (CDK1), kinesin family member 11 (KIF11), EGFR, VEGFA, mitotic arrest deficient 2 like 1 (MAD2L1), and signal transducer and activator of transcription 1 (STAT1), were identified as hub genes. GO enrichment analysis was analyzed by the DAVID database.   Abbreviation: CAM, cell adhesion molecule. KEGG enrichment analysis was analyzed by the DAVID database.

Figure 5. The PPI network of DEGs: red squares represent up-regulated DEGs, and blue circles represent down-regulated DEGs
The size of the squares and circles represents the degree value. PPI was analyzed by the STRING database.

Validation of hub genes expression level
To verify the results of bioinformatics analysis, the qRT-PCR system was used to detect the expression levels of the top ten hub genes in RA synovial tissue and control groups. Statistical results indicate that the expression levels of EGFR and VEGFA in RA synovial tissue were significant (P<0.05) ( Figure 6A). Moreover, to further verify whether these genes were highly expressed in RA, the expression of EGFR and VEGFA were examined by Western blot in RA synovial tissue and normal synovial tissues, and the results showed that the expression levels of VEGFA and EGFR in the synovial tissue of RA were significantly increased (P<0.05). The results were shown in Figure 6B,C.

Discussion
RA is the most common inflammatory arthritis in the world, affecting approximately 1% of the world's population [16]. The incidence of RA increases with age, however, most commonly between 50 and 60 years, usually around 30 years [17]. Early treatment of RA has been shown to prevent irreversible dysfunction and tissue damage, so early identification and treatment of RA are critical for best outcomes [18]. Therapeutic antibodies against inflammatory cytokines or immune cells have become an essential strategy for RA treatment. However, the current molecular biology research of RA cannot provide sufficient data for biotherapy [19]. In the past two decades, the biotherapy of RA has made significant progress. Information obtained from human genomics and proteomics projects has led to a better understanding of the etiology and pathogenesis of the disease and has made significant advances in its treatment. The genomics studies in RA mostly focus on a single genetic event or single cohort study. Our study integrated three profile datasets from GEO database and used bioinformatics tools to analyze microarray data. A total of 444 DEGs were detected, including 172 up-DEGs and 272 down-DEGs, and we performed GO function analysis and pathway analysis of the DEGs, respectively. Next, we constructed a PPI network for DEGs using the STRING database, and the top ten hub genes that were further analyzed by Cytoscape software as follows: PTPRC, LCK, CDC20, JUN, CDK1, KIF11, EGFR, VEGFA, MAD2L1, and STAT1. The statistical results validated by the qRT-PCR and Western blot system show that the expression levels of VEGFA and EGFR were significantly increased in RA samples. As potential biomarkers, VEGFA and EGFR are very important for early diagnosis and clinical treatment of RA, because they can help predict the disease development of subjects at risk, provide prognostic information that can be used to make treatment choices and evaluate treatment response and results, and can monitor disease activity and progress.
Angiogenesis is a prominent feature of RA. The formation of new blood vessels provides oxygen and nutrients to the inflammatory cells, which contributes to the persistence of joint diseases [20]. VEGFA is a critical factor in vascular development, and it is a dominant endothelial growth factor, which is up-regulated by proinflammatory cytokines and hypoxia. It has been reported that VEGFA inhibitors block the new blood vessels of RA and inhibit the transmission of nutrient proteins to the site of inflammation [21]. In our study, VEGFA has significantly increased in RA samples. The KEGG analysis showed that VEGFA was enriched in the RA pathway. The level of VEGFA in the joint fluid of patients with RA was significantly higher than the normal control group. Serum VEGFA level was closely related to RA disease activity, especially related to joint swelling [22]. Previous studies showed that the VEGFA rs699947 C/A functional polymorphism might be associated with the risk of RA in elderly patients [23]. The genetic variation of the VEGFA gene is associated with circulating VEGF-A levels in RA patients [24]. In our previous research, the expression levels of VEGFA have significantly increased in knee osteoarthritis synovial membranes [25,26]. Angiogenesis is essential for the expansion of synovial tissue in RA: pre-existing blood vessels promote blood-derived leukocytes to enter the subsynovial layer, producing and enhancing inflammation. Angiogenesis involves multiple steps, and each step is regulated by specific factors [27]. The inflammation in osteoarthritis is different from the inflammation in RA. The pathogenesis of osteoarthritis is characterized by chronic low-grade inflammation in the synovium, whereas patients with RA manifest with persistent, high-grade systemic inflammation [28,29]. Previous studies have found that the expression of VEGF in all layers of the synovial tissue of RA is much higher than that of osteoarthritis, which may be due to the high expression of VEGF in early RA synovial tissue, which promotes the formation of neovascularization. At the same time, a large number of other inflammatory factors are expressed around the blood vessels, which directly leads to the tortuous proliferation of blood vessels [30]. Many anti-rheumatic drugs used clinically affect the vascular structure; cyclosporine A has recently been shown to inhibit endothelial cell migration and angiogenesis induced by VEGF, which may be a new mechanism of cyclosporine A in the treatment of RA [20,31]. Targeting newly formed RA vasculature may alter disease progression, and VEGFA inhibitors that inhibit neovascularization may be a new target for RA therapy. The levels of IL-6 and VEGFA in the serum of patients with RA vary with the progression of the disease, which may be a rheumatoid effective disease marker [32].
EGFR is a receptor tyrosine kinase that plays a vital role in developing various tumors [33]. The pathological process of RA is similar to tumors, involving both neovascularization and inflammatory cell infiltration. EGFR acts as a mediator of synovial hyperplasia, which includes many aspects of angiogenesis and regulation [34]. Therefore, EGFR is also involved in the progression of RA. EGFR plays an essential role in the pathogenesis of RA by promoting angiogenesis, cytokine production of synovial fibroblasts, and proliferation of endothelial cells [35]. The treatment of RA mainly inhibits the immune system and inhibits synovial and vascular proliferation by blocking the EGFR signaling pathway, which can provide synergy with current therapies. In our study, EGFR has significantly increased in RA samples, and the KEGG analysis showed that EGFR was mainly enriched in adherens junction, MAPK signaling pathway, and ErbB signaling pathway. EGFR is highly expressed in the synovial tissue of RA patients and autoimmune arthritis mice. Also, the study found that the EGFR inhibitor erlotinib can reduce autoimmune arthritis in mice, reducing synovial membrane and vascular membrane formation [36]. Previous studies have found a statistically significant difference in the genotype frequency distribution between the RA patients and the control group at the EGFR rs17337023 SNP locus. Moreover, EGFR gene rs17337023 polymorphism is associated with the risk of RA in the Iranian population [37,38]. In conclusion, EGFR and its ligands may be involved in the development of RA. The development of therapeutic drugs targeting EGFR will bring a new treatment plan for the RA [39].
Taken above, using bioinformatics tools and tissue sample verification, we have identified two candidate genes and pathways, and these findings may significantly improve our understanding of the causes and potential molecular events of RA, and the hub genes VEGFA and EGFR maybe have a crucial function in the pathogenesis of RA and could be used as potential biological and therapeutic indicators of RA. There are still many limitations in our study: small tissue sample size was used for the analyses and mRNA to protein synthesis is regulated by many factors and lack of in vivo experiment. To further verify the experimental results, more experimental in vivo with extensive well-designed studies are needed.