Immune Interaction Map of Human SARS-CoV-2 Target Genes: Implications for Therapeutic Avenues

There exists increasing evidence that people with preceding medical conditions, such as diabetes and cancer, have a higher risk of infection with SARS-CoV-2 and are more vulnerable to severe disease. To get insights into the possible role of the immune system upon COVID-19 infection, 2811 genes of the gene ontology term “immune system process GO: 0002376” were selected for coexpression analysis of the human targets of SARS-CoV-2 (HT-SARS-CoV-2) ACE2, TMPRSS2, and FURIN in tissue samples from patients with cancer and diabetes mellitus. The network between HT-SARS-CoV-2 and immune system process genes was analyzed based on functional protein associations using STRING. In addition, STITCH was employed to determine druggable targets. DPP4 was the only immune system process gene, which was coexpressed with the three HT-SARS-CoV-2 genes, while eight other immune genes were at least coexpressed with two HT-SARS-CoV-2 genes. STRING analysis between immune and HT-SARS-CoV-2 genes plotted 19 associations of which there were eight common networking genes in mixed healthy (323) and pan-cancer (11003) tissues in addition to normal (87), cancer (90), and diabetic (128) pancreatic tissues. Using this approach, three commonly applicable druggable connections between HT-SARS-CoV-2 and immune system process genes were identified. These include positive associations of ACE2—DPP4 and TMPRSS2—SRC as well as a negative association of FURIN with ADAM17. Furthermore, 16 drugs were extracted from STITCH (score <0.8) with 32 target genes. Thus, an immunological network associated with HT-SARS-CoV-2 using bioinformatics tools was identified leading to novel therapeutic opportunities for COVID-19.


INTRODUCTION
Coronaviruses (CoV) are a large family of viruses that cause illnesses ranging from a common cold to more severe diseases. A novel coronavirus (nCoV named SARS-CoV-2) belonging to this family was identified in 2019 in China and in March 2020 the WHO declared a worldwide pandemic of SARS-CoV-2 infection. As of February 13, 2021, more than 108 million cases have been confirmed, with more than 2.28 million deaths attributed to  Even though SARS-CoV and MERS-CoV have previously caused serious outbreaks, the intensity of these was much lower than the ongoing SARS-CoV-2 pandemic. Recent studies have demonstrated that SARS-CoV-2 utilizes the angiotensinconverting enzyme II (ACE-2) receptor (1,2), MERS-CoV, the dipeptidyl peptidase-4 (DPP4/CD26) as an entry (3). The transmembrane serine protease 2 (TMPRSS2) and the proprotein convertase FURIN cleave and activate the viral glycoprotein spike (S), which facilitates virus-cell membrane fusion. SARS-CoV-2 requires both, the ACE2 receptor and TMPRSS2 for protein priming to enter the cell (4), while FURIN facilitates the active binding of SARS-CoV-2 through the ACE2 receptor (5), which is a risk factor for a more severe cause of COVID-19 (4,6).
However, there exists limited information on why SARS-CoV-2 is more lethal and spreads faster than its two ancestors and how the immune system responds to SARS-CoV-2 infection. Recent publications have reported differences in the genome structure of the CoV family members (7) and immunological responses next to inflammation (5,8). The viral entry has been shown to cause a cytokine storm, which may attack lung cells (9). This results in fluid filling the air sacks and alveolar cell apoptosis causing respiratory insufficiency, ultimately leading to death (10). Some strategies have been developed for the treatment of SARS-CoV-2 infection, such as potential vaccines and drugs targeting SARS-CoV-2 or host cell components, e.g., host cell receptors, necessary for viral replication. The complete knowledge of interactions between proteins in SARS-CoV-2 infected cells represents an essential milestone toward understanding the cellular mechanisms and functions of proteins including the genes involved in immune responses against virus infection (11). The gene ontology term (GO: 0002376) comprises genes involved in the development and/or function of the immune system that are necessary for the calibration of immune responses to potential internal or invasive threats.
The demographic data relating to patients and their past medical history are crucial for the development of effective treatment opportunities. Therapeutic targets are urgently needed to manage COVID-19. One possibility is targeting the expression, the transcriptional regulation, or activity of host receptors and associated proteins known to play a critical role in the pathogenicity of CoV infections (4,(15)(16)(17). Indeed, studies using TMPRSS2 transgenic knockout mice have shown that the loss of TMPRSS2 reduces CoV replication in the lungs, and elicits a weaker proinflammatory response, and results in a milder lung pathology (15). In addition, SARS-CoV-2 entry into cells is also decreased upon functional inhibition of TMPRSS2 by the serine protease inhibitor camostat (4). Likewise, ACE2 antibodies or soluble recombinant ACE2 can attenuate viral entry and infection by SARS-CoV-2 (4,16). Thus, a better understanding of the regulatory mechanisms that control expression levels of HT-SARS-CoV-2 might be key for the development of effective novel treatments for SARS-CoV-2 infection.
DPP4/CD26 is a ubiquitous membrane-bound aminopeptidase and has multiple physiological functions, such as the T cell receptor-mediated activation and proliferation of T cells (18) as well as the regulation of glucose homeostasis. Its importance has been highlighted by the approval of DPP4 inhibitors as an established glucose-lowering therapy in type 2 diabetes (19). Modeling the structure of the SARS-CoV-2 spike protein predicts an interaction of DPP4 in addition to ACE2 (20,21). Interestingly, a correlation between DPP4 and ACE2 was found to suggest that both membrane proteins are relevant for virus entry (22). Indeed, DPP4 acted as a CoV co-receptor, suggesting a similar mechanism for the entry of SARS-CoV-2 (23). The coexpression of ACE2 and DPP4 as receptors of the spike glycoprotein postulates that different human CoVs target similar cell types across different human tissues and explains the presence of comparable clinical features in patients infected with different CoVs. Evaluation of SARS-CoV infections revealed that the virus was not only present in the tissues of the lung, liver, kidney, and intestine, but also of the pancreas indicating the pancreas as a potential CoV target (24,25).
In this study, the coexpression of HT-SARS-CoV-2 and immune system process genes was evaluated followed by STRING analysis to determine the functional interactions between HT-SARS-CoV-2 and immune system genes. The STRING database network integrates direct and indirect/functional protein-protein interactions, such as stable physical associations, transient binding, substrate chaining, information relay, and others (26). Despite many limitations due to the low number of patients, the retrospective nature of evidence, and limited patient follow-up, these data provide early insights into how the management of patients with cancer and/or diabetes mellitus might be affected by the SARS-CoV-2 pandemic. This is an important issue, particularly since cancer and diabetes mellitus are risk factors for disease progression and such patients were shown to present more severe symptoms and unfavorable outcomes upon SARS-CoV-2 infection (27)(28)(29).

Samples and Patients
The R2 Genomics Platform is a web-based genomics analysis and visualization platform and allows biomedical researchers to integrate, analyze and visualize clinical and genomics data. The transcriptome data of healthy and cancer patients were sourced from The Cancer Genome Atlas (TCGA) dataset (http://cancergenome.nih.gov/) and microarray data of the NCBI Gene Expression Omnibus (NCBI GEO) (30). The datasets were available on R2 Genomics Analysis and Visualization Platform (http://r2. amc.nl) and ImmuCo (31). Human samples and immune cell, subpopulations were then analyzed and samples included:

Gene Ontology Enrichment Analyses
The Gene Ontology source (GO; http://geneontology.org) provides structured, computable knowledge regarding the functions and products of genes. GO enrichment analyses included three categories, biological process (BP), molecular function (MF), and cellular components (CC). For our studies, the GO term "immune system process (GO: 0002376)" was selected, which comprises 2811 genes.

Networking Analyses and Drug Screening
The protein-protein network was constructed using the STRING Database (String v10) (36) as a protein-protein interaction network, which also analyzed functional enrichments. STRING consisted of 24.6 million proteins with more than 2,000 million interactions. The coexpression scores in STRING v10 were computed by employing a revised and improved pipeline. Using high throughput microarray gene expression data from the NCBI GEO, we mined databases and literature as well as predictions based on genomic context analysis. The gene expression values are then correlated gene-by-gene (Pearson correlation) and the resulting values were calibrated against common members in KEGG pathway maps to compute STRING scores. Text mining score was computed by parsed search of statistically relevant cooccurrences of gene names in the scientific texts (SGD, OMIM, FlyBase, PubMed). The parameters included in the analysis for screening targets were a medium confidence score of 0.4 and a medium FDR stringency of 5%.
In the search for significant drug molecules and determine potential drug-associated mechanisms upon SARS-CoV-2 infections, we employed the Search Tool for Interactions of Chemicals (STITCH, stitch.embl.de), a chemical protein interaction database that contains 0.5 million compounds and 9.6 million proteins. STITCH ("search tool for interactions of chemicals") integrates information about interactions from metabolic pathways, crystal structures, binding experiments, and drug-target relationships. The inferred information from phenotypic effects, text mining, and chemical structure similarity is used to predict relations between chemicals. STITCH further allowed exploration of the network of chemical relations, also in the context of associated binding proteins. Each proposed interaction can be traced back to the source of that data.

Correlation Analysis
Gene correlation analyses were performed with the R2: Genomics Analysis and Visualization Platform and ImmuCo. Correlation statistics were reflected by Pearson correlation coefficients and Pearson's r is a measure of linear correlation between two sets of data. It is the covariance of two variables; it is essentially a normalized measurement of the covariance, such that the result always has a value between −1 and 1. A p-value < 0.05 was considered statistically significant.

Survival Analysis
PrognoScan (http://www.prognoscan.org/) is a large collection of publicly available cancer microarray datasets with clinical annotation and a tool for assessing the biological relationship between gene expression and prognosis. In the PrognoScan database, the association of gene expression with the survival of patients were evaluated by the minimum p-value approach. Briefly, patients were first arranged by expression levels of a given gene. They were then divided into high-and lowexpression groups, and the risk differences between any 2 groups were estimated by the log-rank test. Meta-survival analysis of pancreatic cancers was computed and derived from the GENT2 (http://gent2.appex.kr/gent2/).

Statistical Analysis
Statistical computations and drawings of figures were performed with several packages (ggplot2 and circlize) in the statistical software environment R, version 3.3.2 (http://www.r-project.org).

Identification of Protein Associations
The three human proteins ACE2, TMPRSS2, and FURIN used by SARS-CoV-2 to efficiently enter into human cells were studied for their functional protein associations using the STRING database. TMPRSS2 was networked as a center point associated with both ACE2 and FURIN, while ACE2 and FURIN showed no interactions between each other (Supplementary Figure 1). To study the immune system interactions of the human targets, the GO term immune system processes (GO: 0002376) in the Gene Ontology resource was analyzed.
In total, 2811 genes (Supplementary Table 1) were studied for functional protein association networks using STRING  Frontiers in Immunology | www.frontiersin.org The association of HT-SARS-Cov-2 with immune relevant genes was determined by STRING analyses and represented as coexpression, text mining, and score.
The STRING network analysis screened 102 of the HT-SARS-CoV-2association genes from the 2811 immune system process' genes. Correlation of HT-SARS-Cov2 with immune process relevant genes was performed by RZ2 genomics using mixed healthy and mixed cancer data sets and presented as correlation values (R) including the statistical significance (P).

Determination of Genes Correlating With HT-SARS-CoV-2 in Cancer
The ACE2, TMPRSS2, and FURIN expression were correlated to every single gene in the five different datasets described in the Materials and Methods section above. The top ∼ 1500 genes of the respective association, i.e., having a highly significant correlation (p < 0.01) in the corresponding dataset, are summarized in Supplementary Table 3. From the common STRING associations between the human target genes, 19 associations of eight commonly networking genes were plotted to the three HT-SARS-CoV-2 ( This collected data enabled classification into four distinct categories: (i) associations between human target genes; (ii) triple association of DPP4; and, (iii) double associations of FURIN with either ACE2; or (iv) with TMPRSS2, respectively. The top segment of Table 2 Table 2). Normal and cancer samples mostly followed the same pattern. However, ACE2/GP2 (R 0.108 in normal and R −0.045 in cancer) and FURIN/GAPDH (R −0.198 in normal and R 0.023) present differential correlations between normal and cancer samples.
To explore gene association patterns in the wide range of datasets, ACE2, TMPRSS2, and FURIN were correlated to DPP4, SRC, and ADAM17 expression regardless of the tissue and the disease. The analysis demonstrated that ACE2 and TMPRSS2 displayed a positive correlation, whereas FURIN showed a negative correlation with their counterparts ( Table 3 and  Frontiers in Immunology | www.frontiersin.org Correlation of HT-SARS-CoV-2 with immune process relevant genes was performed by R2 genomics using data sets obtained from normal pancreatic, pancreatic cancer, and diabetic pancreatic tissues. Data are presented as correlation values (R), including the statistical significance (P).

Correlation Pattern of HT-SARS-CoV-2 in Pancreatic Tissues
ACE2 expression in the pancreas was reported to cause pancreatic damage after SARS-CoV-2 infection (37). Interestingly, ACE2 and DPP4 were connected irrespective of the pancreatic nature i.e., healthy, cancer and diabetes mellitus. The healthy pancreatic tissues exhibited distinct associations when compared to pancreatic cancer tissues or tissue samples from patients suffering from diabetes. DPP4 was differentially expressed with TMPRSS2 and FURIN, while TMPRSS2 was highly expressed in healthy and only scanty expressed in the diseased pancreatic tissue. The HT-SARS-CoV-2 ACE2 and TMPRSS2 were positively associated with healthy pancreatic tissue, while a non-significant positive association was noted in the context of cancer and diabetes mellitus. A positive correlation between TMPRSS2 and FURIN was found in the diseased, but not in the healthy pancreas ( Table 3). FURIN expression was not correlated with DPP4 except for pancreatic cancer, where a negative correlation was found. The FURIN/ADAM17 axis remained negatively associated with pancreatic tissues. Interestingly, TMPRSS2 positively correlated with SRC in the pancreatic tissues in mixed healthy and cancer datasets. Other double associated genes lacked any specific correlation pattern. FURIN and MYC were positively associated with healthy tissues and pancreatic cancer. However, FURIN positively acted with MYC in pancreatic cancer, but not in the mixed cancer dataset. AKT1 with TMPRSS2 and FURIN were positively correlated in all contexts except for pancreatic cancer and healthy pancreas, respectively ( Table 3). Interestingly, many networking genes were not correlated in pancreatic tissues. These included ACE2/ADAM17, FURIN/GAPDH, ACE2/FURIN, ACE2/GP2, and ACE2/TFRC, respectively.
ACE2 and TMPRSS2 were positively correlated in immune cells like B cells, T cells, dendritic cells (DC), macrophages, monocytes, and NK cells. In contrast, TMPRSS2 and FURIN exhibited positive associations in B cells and monocytes, while other immune cell subpopulations lacked any correlations. Except for B cells, immune cells displayed no correlations with the central association markers, DPP4/ACE2 or TMPRSS2. FURIN and DPP4 revealed positive correlations in CD4 + and CD8 + T cells and negatively correlated in DC and monocytes. FURIN and ADAM17 showed a general significant association with a positive association in B cells, T cells, monocytes, and NK cells, but negative correlations in DC and macrophages ( Table 4 and Supplementary Table 4). Irrespective of any immune cell type, ACE2 and GP2 were significantly positively correlated to each other. It is noteworthy that many networking genes were not associated with the immune cells under study.

Identification of Drug Candidates for HT-SARS-CoV-2 and Co-expressed Genes
Altogether ∼0.5 million compounds were studied on the STITCH database. Among these compounds, 16 target drugs were found according to the screening conditions ( Table 5). The active components involved 32 drug-associated genes (Supplementary Table 4), which were acquired by STITCH analysis. As shown in Figure 4, ACE2/DPP4 ( Figure 4A) and FURIN/ADAM17 ( Figure 4C) were predicted to have drug connections, whereas TMPRSS2/SRC ( Figure 4B) had no predicted drug connection. Several inhibitors of the angiotensin receptor (ACE) and DPP4, in particular gliptin members, such as sitagliptin, teneligliptin, linagliptin, vildagliptin, and saxagliptin, were identified. For FURIN and ADAM17, a therapeutic involvement via calcium ions and Zn (II) and for ADAM17, an involvement via IK-682 were identified ( Table 5). Among these three interactions, two druggable connections were suggested for a putative therapeutic option.

Survival Analysis on HT-SARS-CoV-2 and Druggable Targets
HT-SARS-CoV-2 and their drug connections were tested for their association with overall survival (OS) rates in 11 different cancer types (Supplementary Table 6). Although many cancer types did not show a significant survival ratio, ACE2 was significant in brain cancer (HR 0.44; p 0.003), breast cancer (HR 1.23; p 0.002), lung cancer (HR 0.7; p 0.009), ovarian cancer (HR 0.63; p 0.048) and renal cell carcinoma (HR 0.17; p 0.021). ADAM17 showed high hazard ratios in bladder (HR 4.65; p 0.04), brain (HR 3.65; p 0.0004) and colorectal (HR 2.0; p 0.024) cancers. The meta-survival analysis on pancreatic cancer datasets only showed significance in ADAM17 (HR 1.21; p 0.08) (Supplementary Figure 4). OS did not follow any relations with the correlation and networking of HT-SARS-CoV-2.

DISCUSSION
In this comprehensive study, the critical immune-related associations with the HT-SARS-CoV-2 ACE2, TMPRSS2, and FURIN were identified using bioinformatics analyses. SARS-CoV-2 requires the proteins ACE2, TMPRSS2, and FURIN for cell entry (4), which are the same proteins used by SARS-CoV (38), the cause of the 2002-2004 SARS epidemic. The SARS-CoV-2 viral genome encodes for 29 different proteins. A recent study by Gordon and et al. expressed and purified 26 of these proteins and analyzed them for protein-protein interactions (39). Using this approach, 332 high confident protein-protein interactions between SARS-CoV-2 and human proteins and 66 druggable human proteins were identified. Interestingly, two sets of pharmacological drugs had antiviral activities. These included inhibitors of mRNA translation of two sigma receptors.
FURIN and ADAM17 are both intertwined with Notch, an evolutionarily conserved communication system between adjacent cells, and thus targeting of Notch might represent an alternative approach to inhibit FURIN and upregulate ADAM17 (40). In this context, it is noteworthy that ADAM17 is involved in the shedding of ACE2 (41). Such hypothesis-driven studies based on the knowledge of the molecular details of virus-cell interaction are still crucial for the identification of therapeutic targets to treat COVID-19.
SARS-CoV-2 infection of pancreatic endocrine cells via ACE2 appears an unlikely central pathogenic feature of SARS-CoV-2 as it relates to diabetes (42). However, we found that ACE2 was highly expressed in pancreas islets. This could be an explanation for the finding that SARS-CoV-2 infections can cause damage of the islets and subsequent acute diabetes (25,43). In our studies, we found that ACE2 and DPP4 were positively connected irrespective of the pancreatic nature, i.e., normal, cancer and diabetes. A positive correlation of ACE2 and DPP4 with TMPRSS2 was lacking in pancreatic cancer and diabetes.
Studies have shown a strong correlation between the severity of COVID-19 and concentrations of different cytokines, such as IL2, IL7, IL10, GCSF, MCP1, and TNF-alpha (44). Notch and IL-6 cooperate to activate the immune system and perpetuate the cytokine storm. In macrophages, Delta Like Canonical Notch Ligand 4 (Dll4)/Notch signaling promotes the production of inflammatory cytokines, including IL-6 (40). IL-6, in turn, increases the expression of the Notch ligands (Dll1 and DII4), thus amplifying the Notch signaling, thereby establishing a feedback loop that promotes the further production of IL-6. In T cells, Notch signaling triggered by Dll1/Dll4 ligands promotes inflammatory Th1/Th17 cytokines, while Jagged1 ligands dampen the IL-6-induced Th17 activation (40).
Currently, there are few reports on the specific mechanisms of interaction networks and associations of TMPRSS2, ACE2, Correlation of HT-SARS-CoV-2 with immune process relevant genes was performed by ImmuCo using data sets obtained from PBMC, B cells, and T cells. Data are presented as correlation values (R), including the statistical significance (P). Drug candidates of HT-SARS-CoV-2, and immune process relevant genes were determined by STRING analysis and experimental interaction, database annotation, automatic text mining presented. Besides, the data were summarized in a combined score.
and FURIN with their targets. Using bioinformatics tools, we here show that ACE2 may participate in Notch signaling to exert immune effects. However, the relationship between these molecules has to be analyzed and further experimentally validated concerning the possible targeting of ACE2/DPP4 and FURIN/ADAM17-related factors in diseased COVID-19 patients. The delicate balance between ACE2, ADAM17, and TMPRSS2 interactions could be decisive for the clinical outcome of COVID-19 (45). Angiotensin receptor blockers (ARBs) similar to angiotensinconverting enzyme (ACE) inhibitors, act by preventing the formation of angiotensin II rather than by blocking the binding of angiotensin II to muscle cells on blood vessels. The therapeutic indication for ARBs prescription comprises control of high blood pressure, treating heart failure, and preventing kidney failure, especially in patients suffering from diabetes mellitus. Since ACE inhibitors are associated with cough, accumulation of bradykinin, and angioedema, ARBs, like losartan and captopril, might be a more favorable therapeutic option by blocking the binding and attachment of the SARS-CoV-2 receptor binding domain to ACE2-expressing cells, thereby inhibiting the infection of host cells (46). Ferrario and colleagues reported that the administration of either ACE inhibitors or ARBs increased the levels of ACE2 mRNA in Lewis rats, compared to rats receiving placebo (47). In particular, the cardiac levels of ACE2 mRNA expression increased by 4.7-fold or 2.8-fold in rats treated with losartan, which was accompanied by increased ACE2 activity. Treatment with the ACE inhibitor captopril can significantly increase ACE2 protein expression in rats with acute lung injury (48). Furthermore, Mourad and Levy suggested aliskiren treatment as an option in the context of SARS-CoV-2 infection, since it can reduce the expression of ACE2 (49).
Early trials in the 1990s showed that DPP4 inhibitors improve glycaemia in animals. Subsequent clinical studies during the 2000s showed a glucose-lowering authority of DPP4 inhibitors in humans with type 2 diabetes as monotherapy or in combination with other therapies, i.e., metformin, sulfonylureas, tiazolidinediones, or exogenous insulin. Five DPP-4 inhibitors (sitagliptin, vildagliptin, alogliptin, saxagliptin, and linagliptin) were approved by regulatory authorities for the treatment of type 2 diabetes and entered the market between 2006 and 2013 (50).
The inhibition of DPP4 led to the suppression of breast cancer tumor growth (51) and displayed a positive trend in OS of colorectal cancer (52) and prostate cancer (53). Analysis of colorectal cancer and lung cancer patients showed a similar trend toward beneficial effects associated with DPP4 inhibition (54); however, further studies are required on the clinical impact of DPP4 inhibition on tumor patients. Based on our results it is suggested that the modulation of DPP4 might interfere with the infection and/or progression of COVID-19 and therefore might represent a therapeutic strategy. Since DPP4 inhibitors have been shown to modulate inflammation and exert anti-fibrotic activity, they might have the potential to affect the hyperinflammatory state associated with severe COVID-19 (55).

CONCLUSIONS
COVID-19 is a devastating novel disease resulting in thousands of deaths worldwide. It is essential to find causal therapies for those patients with a severe course, among them patients with comorbidities, as quickly as possible. The results of the present study suggest novel drug candidates for COVID-19 and underscore that new treatment options can be postulated by using different bioinformatics tools, which now have to be experimentally proven. Screening based on drug targets and their mechanisms is still key to developing new drugs. The present study provides comprehensive bioinformatics information on genes that are potential targets and candidate drugs for SARS-CoV-2. This information may help us to better understand the relevant molecular targeting mechanisms of this disease. However, further explorations are required to evaluate the biological functions of these biomarkers and drugs in the pathogenesis of the disease.

DATA AVAILABILITY STATEMENT
The original contributions generated in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.

AUTHOR CONTRIBUTIONS
KS and BS designed the study and wrote the article. KU and KS performed the data analyses. MB and CW discussed the data and revised the manuscript. All authors contributed to the article and approved the submitted version.