Integrated Analysis of Key Pathways and Drug Targets Associated With Vogt-Koyanagi-Harada Disease

Background Vogt-Koyanagi-Harada (VKH) disease is a complex disease associated with multiple molecular immunological mechanisms. As the underlying mechanism for VKH disease is unclear, we hope to utilize an integrated analysis of key pathways and drug targets to develop novel therapeutic strategies. Methods Candidate genes and proteins involved in VKH disease were identified through text-mining in the PubMed database. The GO and KEGG pathway analyses were used to examine the biological functions of the involved pathways associated with this disease. Molecule-related drugs were predicted through Drug-Gene Interaction Database (DGIdb) analysis. Results A total of 48 genes and 54 proteins were associated with VKH disease. Forty-two significantly altered pathways were identified through pathway analysis and were mainly related to immune and inflammatory responses. The top five of significantly altered pathways were termed as “inflammatory bowel disease,” “cytokine-cytokine receptor interaction,” “allograft rejection,” “antigen processing,” and “presentation and Herpes simplex infection” in the KEGG database. IFN-γ and IL-6 were identified as the key genes through network analysis. The DGIdb analysis predicted 48 medicines as possible drugs for VKH disease, among which Interferon Alfa-2B was co-associated both with IFN-γ and IL-6. Conclusions In this study, systematic analyses were utilized to detect key pathways and drug targets in VKH disease via bioinformatics analysis. IFN-γ and IL-6 were identified as the key mediators and possible drug targets in VKH disease. Interferon Alfa-2B was predicted to be a potentially effective drug for VKH disease treatment by targeting IFN-γ and IL-6, which warrants further experimental and clinical investigations.


INTRODUCTION
Vogt-Koyanagi-Harada (VKH) disease is an immune-mediated disorder characterized by chronic, bilateral granulomatous panuveitis, often associated with neurological, audiovestibular and cutaneous manifestations (1). VKH disease is more commonly seen in Asians, Hispanics, Native Americans (2), and rare in Africans (3). Bilateral panuveitis, hearing disorder and meningitis are the main clinical features. Treatment with systemic corticosteroid is the mainstay of VKH disease therapy in the acute uveitic phase (4). However, despite proper treatment with a high-dose of corticosteroid, 79% of patients will experience recurrent attacks and develop chronic disease (5). Moreover, a high-dose of corticosteroid over a prolonged period may lead to side effects such as Cushing syndrome, hyperglycemia, and increased incidence of severe infections (6). Currently, there remain unmet medical needs for novel therapies that can etiologically target the molecules or immune mediators involved in the disease.
Although the exact biological mechanisms are still unclear, an increasing number of candidate genes and proteins have been reported to be involved in the development of VKH disease, which may be possible drug targets for the disease (7,8). However, it is still challenging to prioritize these drug targets among many genes and proteins. Here, we used integrated bioinformatic analysis to summarize the candidate genes and proteins associated with VKH disease and identify the potential key pathways and drug targets, which may help to develop new therapeutic agents.

Identifying Candidate Genes and Proteins Associated With VKH Disease
We manually collected candidate genes and proteins associated with VKH disease by a thorough review of the literature in any language published from May 1981 to November 2019, using a similar approach used earlier by others (9,10). We used the following terms to search the PubMed database: "idiopathic uveoencephalitis" OR uveoencephalitis OR "uveomeningitic syndrome" OR Vogt-Koyanagi-Harada. Studies were eligible if they included a comparison of the expression levels of a gene or protein between VKH patients and controls. Key exclusion criteria included: (i) studies only conducted in animal models, (ii) meta-analysis of published results, and (iii) review and comment papers. We also matched eligible genes and proteins with our previously published database (UVEOGENE, http:// www.uvogene.com) (11). Two researchers were trained in each step with pilot tests before collecting and managing data independently. Regular meetings were held to clear up any misunderstandings or disagreements.

Functional Analysis
Functional analysis was performed based on the DAVID online tools (version 6.8, http://david.ncifcrf.gov). In the gene ontology (GO) database, the analysis of candidate genes and proteins was divided into three categories, termed as biological processes (BP), molecular functions (MF), and cellular components (CC) (12). The GO database provides annotations to describe the properties of genes and gene products of different organisms and shows enriched genes' potential function. BP is an ordered combination of molecular functions to describe a wide range of biological processes. The MF is used to describe the function of a gene or gene product, and the CC is designed for describing subcellular structures, locations, and macromolecular complexes of genes. We submitted the list of identified candidate genes to the DAVID online tools and obtained the significant enrichment of the above categories. A significant threshold of P < 0.05 was used.

KEGG Pathway Enrichment Analysis
Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis was performed on candidate genes to enrich significantly altered pathways, using the DAVID online tools regarding the KEGG database. KEGG is an encyclopedia of genes and genomes for understanding high-level biological functions and utilities (13). The KEGG database provides a collection of graphical diagrams (pathway maps), in which some of the known metabolic pathways and regulatory pathways are displayed to describe the linkage, interaction and function of enriched candidate genes across different cells and organisms. In this study, enriched pathways were considered significant if including at least two genes and reaching a significance level of P < 0.05.

Protein-Protein Interaction Analysis
Protein-protein interaction (PPI) analysis was performed based on candidate genes and proteins using the STRING database (https://stringdb.org/). The STRING database aims to collect and integrate the information by integrating known and predicted protein-protein association data from many organisms (14). Default parameters in STRING were used. Cytoscape software (version 3.7.2, California, USA) was used to construct and visualize the PPI network. Cytoscape is a graphical display tool used for visualizing complex biological networks (15). Furthermore, we used the Cytoscape plugin molecular complex detection (MCODE) to explore significant modules in the PPI network using the following scores and parameters: k score= 2, degree cutoff= 2, node score cutoff= 0.2, and maximum depth= 100 (16).

Identification of Hub Genes
CytoHubba plugin in Cytoscape was used to identify hub genes in the PPI network by calculating and analyzing the network structure. Three algorithms were used to generate intersecting genes, including the maximum neighborhood component (MNC), the degree, and the closeness, as described previously (17). Briefly, the degree is the number of the edges of a gene in the network, representing the interaction pairs with others. The closeness could be used to evaluate the genes in the network according to the calculated centrality. The MNC represents the number of nodes in the maximum connected subgraph. The genes generated by the above algorithms' intersection were more likely to locate in a core position and were considered as the hub genes.

Drug-Gene Interactions
The Drug-Gene Interaction Database (DGIdb, http://dgidb. genome.wustl.edu/) is a drug prediction database, which can be used to screen drugs that potentially target certain genes of interest (18). The DGIdb provides gene-drug interactions and potential druggability according to their gene category. We utilized the DGIdb analysis to obtain all possible gene-drug interactions for the top-two ranked molecules of hub genes. Cytoscape was used to visualize the acquired interactions.

Genes and Proteins Acquisition
After screening a total of 1,323 articles from PubMed, we obtained 128 eligible articles. Based on the inclusion and exclusion criteria, we identified 48 genes and 54 proteins associated with VKH disease ( Table 1). The expression of these candidate genes was reported to be significantly changed in VKH disease as compared with healthy individuals. Besides, among these candidate genes, 36 genes overlapped with genes recorded in the Uveogene database for VKH disease (11). These 36 genes had been reported to have susceptible single nucleotide polymorphisms and loci for VKH disease in the Uveogene database. Figure 1 illustrates the flow chart of this study.

GO Analysis of Genes and Proteins
GO enrichment analysis was performed with the DAVID online tools to examine the identified genes and proteins' biological characteristics. The analysis of BP (biological processes) showed a total of 251 functions, 200 of which were significantly enriched (P < 0.05). The top-ranked functions included the categories "immune response," "inflammatory response," "positive regulation of T cell proliferation," "positive regulation of tyrosine phosphorylation of Stat3 protein," and "interferongamma-mediated signaling pathway." The CC (cellular components) analysis included 29 functions, of which 26 were significantly enriched (P < 0.05). For the CC analysis, the identified genes and proteins were mostly enriched in the "external side of plasma membrane," "extracellular space," "integral component of luminal side of endoplasmic reticulum membrane," "extracellular region," and "MHC class II protein complex." The MF (molecular functions) analysis included 37 functions, 26 of which were significantly enriched (P < 0.05). Changes in MF were significantly enriched in "cytokine activity," "peptide antigen binding," "MHC class II receptor activity," "growth factor activity," and "chemokine activity." The top ten functional enrichment analyses of GO (BP, MF, CC) are shown in Figure 2, and the significant GO (BP, MF, CC) are provided in Supplementary Table S1. The corresponding genes enriched in GO analysis ( Figure 2) are listed in Supplementary Table S2.

KEGG Analysis of Genes and Proteins
The DAVID online tools were utilized for the KEGG enrichment pathway analysis to show the potential involvement of pathways related to identified candidate genes and proteins. The analysis of KEGG identified 42 significantly altered pathways (P < 0.05) (Supplementary Table S3). The top-ranked pathways were   mainly involved in categories termed as "inflammatory bowel disease (IBD)," "cytokine-cytokine receptor interaction," "allograft rejection," "antigen processing," and "presentation and Herpes simplex infection" (Figure 3). The corresponding genes enriched in the KEGG pathway ( Figure 3) are listed in Supplementary   Table S4. Additionally, the IBD pathway was the most significant in the enrichment analysis with 19 genes involved ( Figure 4). The cytokine-cytokine receptor interaction had the largest number of genes, 27 genes enriched in the pathway ( Figure 5). Moreover, several genes, including IFN-g, IL6, IL12, IL4, IL23R and IL21, were shared by the IBD pathway and the cytokine-cytokine receptor interaction pathway and were related to the signaling transductions by members of the interleukin-family indicating that these members of the interleukin-family might play an important role in the pathogenesis of VKH disease.

PPI Network Analysis and Related Gene Modules
The PPI network consisted of STRING and Cytoscape, which were used to determine the most important genes and proteins clusters. All the 87 nodes and 754 edges in the PPI network are shown in Figure 6. Besides, the MCODE plugin in Cytoscape identified two main modules. Cluster 1 (score = 23.28) consisted of 26 nodes and 291 edges, and cluster 2 (score = 6.889) consisted of 10 nodes and 31 edges ( Figure 6). In cluster 1, several genes, including IFN-g, IL6, IL4, IL23R and IL21, were significantly enriched. These genes were also related to the enriched pathways, including the IBD

Hub Genes Recognition
CytoHubba was utilized to search for the key genes and calculated all the molecular nodes and edges. Consequently, 102 target genes were involved in the PPI network complex, forming 87 nodes and  754 edges (Figure 7). To search for the important nodes in the PPI network, all nodes were ranked by the three algorithms, including the degree, closeness and MNC provided by cytoHubba. The cytoHubba plugin Cytoscape was used to analyze the hub genes in the PPI network, and the following genes with the top ten grades were identified as hub genes: IL6, IFN-g, IL4, CTLA4, IL1B, STAT3, CCL2, CD40, FOXP3 and IL2RA. Among these, IFN-g and IL-6 were the top two of the ten grades and considered to be the key genes in this model, given the fact that the products of genes were at the core of the PPI network (Figure 7). The descriptions of gene symbols and the details shown in Figure 7 are provided in Supplementary Table S5.

Drug-Gene Interactions
IFN-g and IL-6, as the identified key genes, were entered into the DGIdb to obtain potential drugs. 48 drugs were identified in the DGIdb analysis; the information about the source, scores, and interaction type of target drugs is provided in the Supplementary  Table S6. Most of the target drugs were inhibitors, monoclonal antibodies and immunomodulatory agents. Among these identified drugs, Interferon Alfa-2B was predicted to act on both IFN-g and IL-6 ( Figure 8).

DISCUSSION
In the present study, we used systematic analyses to show key pathways and potential VKH disease drugs. Two significant modules and a top ten hub genes were detected with IFN-g and IL-6 as key genes. 48 target drugs were potentially useful drugs for the treatment of VKH disease. Interferon Alfa-2B targeted both IFN-g and IL-6 and predicted a potentially useful drug to treat VKH, and warrants further experimental and clinical investigations.
The results of GO enrichment analysis indicated that immune response and inflammatory response play a significant role in VKH disease, which confirms earlier studies in this field (19)(20)(21)(22)(23). There still are many mediators related to these pathways, which have not been reported to be associated with VKH disease. Our Frontiers in Immunology | www.frontiersin.org December 2020 | Volume 11 | Article 587443 enrichment analysis suggests that these molecules along with their involved pathways are closely linked with the development and pathogenic processes of VKH disease, which warrants further experimental investigations. Various mediators, including HOXB3, GH1, KIR2DL4 have not been reported in VKH disease, but these were predicted to have a high relevance for VKH disease as shown by our PPI network analysis. Previous studies have suggested that HOXB3, GH1 and KIR2DL4 are involved in autoimmune disease such as Thyroid-associated orbitopathy (TAO), Type 1 diabetes (T1DM), and Systemic lupus erythematosus (SLE) (24)(25)(26). Their functional role in VKH disease requires further studies. We also identified a variety of pathways related to certain autoimmune diseases such as the IBD pathway, the Toll-Like receptor pathway, the CD27-CD70 pathway and the CD40-CD40L pathway which have been shown to be associated with systemic lupus erythematosus (SLE), multiple sclerosis (MS), rheumatoid arthritis (RA), systemic sclerosis (SSc), Sjögren's syndrome (SS), psoriasis, uveitis and other autoimmune diseases (11,(27)(28)(29). These findings suggest a certain degree of shared pathogenic pathways between VKH disease and other autoimmune diseases.
Two key genes, IFN-g and IL-6, were identified in our study by cytoHubba, a plug-in Cytoscape. Lymphocytes produce IFN-g in response to various immune stimuli. The essential role of IFN-g in human anti-viral immunity has been illustrated earlier (30,31). Several studies have indicated that VKH is an autoimmune disease mediated by Th1/IFN-g and Th17/IL-17 pathways (32). Examples include studies that showed that the expression level of IFN-g is significantly higher in peripheral blood mononuclear cell (PBMC), aqueous or serum of VKH patients as compared with those in control subjects (33)(34)(35)(36). Besides IFN-g we also identified the cytokine IL-6 as a key player in VKH disease. IL-6 is a four-helix cytokine composed of 184 amino acids (37) with various physiological functions, including regulating the proliferation and differentiation of immune cells. IL-6 modulates almost all aspects of the innate immune system. It has been shown that IL-6 plays a significant role in regulating the balance between IL-17 producing Th17 cells and regulatory T cells (Treg) (38,39). Both T-cell subsets play an important role in the pathogenesis of VKH disease. Several studies have shown that the concentration of IL-6 in PBMC, monocyte-derived macrophages (MDMs), or aqueous humor from VKH patients is significantly higher than that observed in controls (40)(41)(42). This evidence support that IFN-g and IL-6 are key mediators related to VKH disease, which suggests that they might be an attractive drug target for this disease. According to the analysis of the DGIdb, Interferon Alfa-2B is specific for both IFNgamma and IL-6. The Drug-Gene Interaction database (DGIdb) mines available resources and predicts potentially effective FIGURE 6 | The PPI network was constructed by STRING and Cytoscape. There were 87 nodes and 754 edges in the PPI network. Besides, two modules were recognized by MCODE in Cytoscape. Cluster 1 (score = 23.28) was the most powerful module containing 26 nodes and 291 edges. Cluster 2 (score = 6.889) consisted of 10 nodes and 31 edges. Molecules represented as a pink diamond in the circle with a blue background represents genes belonging to cluster 1, pink triangles in the circle with a purple background represent genes belonging to cluster 2, and genes in the white background do not belong to either cluster. HOXB3, KIR2DL4, GH1, EBI3 and GPR29 in pink graphics were predicted genes using interactions between submitted genes in the PPI network by Cytoscape, and other nodes were of submitted genes in the PPI network and showed interactions between multiple genes. FIGURE 7 | Using Cytohubba, a plugin Cytoscape, the top 10 genes were identified as hub genes. This analysis revealed two key genes: IFN-g and IL-6 in the cluster. All the gene nodes and edges were calculated. Molecules represented as pink diamonds in the circle with a blue background represented genes belonging to the cluster showing the top ten genes in the network by the Cytoscape plugin Cytohubba, and pink circles without a colored background represent genes not belonging to this cluster. HOXB3, KIR2DL4, GH1, EBI3 and GPR29 in pink graphics were predicted genes using interactions between submitted genes in the PPI network by Cytoscape, and other nodes were of submitted genes in the PPI network and showed interactions between multiple genes. FIGURE 8 | Drug targeting the two genes, IFN-g and IL-6, shown as pink circles. 48 drugs as the main therapeutic agents shown in triangles and expressed in different groups according to scores. Firstly, the triangle in green represented GLUCOSAMINE obtained five scores through the DGIdb analysis, and it was the most relevant drug. Secondly, triangles in purple, including FONTOLIZUMAB and SILTUXIMAB acquired four scores, and these were the second most relevant drugs. Thirdly, drugs with three scores were OLSALAZINE and GINSENG ASIAN as shown by dark blue triangles, and INTERFERON ALFA-2B shown as dark blue diamonds since Interferon Alfa-2B was able to target both of the two key genes. Moreover, most drugs obtained two scores and are shown as light blue triangles (the target drugs with scores can be seen in Supplementary Table S6). Lastly, eight drugs, including FUMARIC ACID, APREMILAST, VX-702, CDP-6038, SIRUKUMAB, OLOKIZUMAB, ELSILIMOMAB and PF-04236921, obtained one score and are shown as blue-grey triangles. The scores of 48 drugs identified by the DGIdb show the importance of drugs, with a higher score indicating greater importance. therapeutic targets or prioritized drug development based on specific genes (18,43,44). We performed drug-gene interaction networks through bioinformatics analysis to identify target drugs that may act on both IFN-g and IL-6. Of the 48 target drugs obtained from the DGIdb, most were inhibitors, monoclonal antibodies or immunomodulators. Among these potential drugs, Interferon Alfa-2B was found to be the drug that could target both IFN-g as well as IL-6. The DGIdb provides evidence showing that Interferon Alfa-2B affects the expression of both IFN-g and IL-6 (45,46). Interferon, a class of cytokines, can interfere with virus replication, reduce cell proliferation, and alter immunity (47,48). In recent studies, the role of IFN-Alfa in the pathogenesis of autoimmune diseases has been recognized (49,50). Interferon Alfa-2B is an effective drug for treating autoimmune diseases such as idiopathic thrombocytopenic purpura (ITP) and uveitis (51,52). The use of Interferon Alfa-2B has also been reported in the treatment of severe chronic uveitis in patients with Behcȩt's disease (53) and uveitic cystoid macular edema (54). VKH disease, along with Behcȩt's disease and uveitic cystoid macular edema, are anatomically classified as noninfectious posterior or pan-uveitis and can be treated with the same class of immunomodulatory drugs such as cyclosporin A (55,56). The efficacy of Interferon Alfa-2B in the treatment of VKH disease has not yet been reported, but our analyses highlight the potential of Interferon Alfa-2B for the treatment of this disease, which necessitates further studies.
This study has some major limitations. It should be noted that we were dependent on existing data by integrative bioinformatic analysis, and our analyses were based on currently available information obtained from existing research surveys, suggesting that new information from future studies may influence the results presented here. Due to the lack of available data, dynamic networks' development is not yet based on a genetic-epigenetic association. Moreover, our analyses are largely exploratory, and these results need to be further confirmed by experimental in vitro and in vivo studies.

CONCLUSION
In this study, systematic analyses were performed to identify key pathways and drug targets in VKH disease via bioinformatics analysis. Two significant modules and a top ten hub genes in VKH disease were detected with IFN-g and IL-6 as the top two genes. The study furthermore predicted Interferon Alfa-2B as a potentially useful drug for the treatment of VKH disease.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.