key: cord-0786399-tgvocrwh authors: Dai, Wenjuan; Cao, Dawei; Zhang, Wei; Wei, Yangyang; Ding, Daqing; Li, Bei; Gao, Yan; Zhao, Lixuan; Jiang, Yi; Kong, Xiaomei title: Integrated Bioinformatics Analysis Reveals Key Candidate Genes and Cytokine Pathways Involved in COVID-19 After Rhinovirus Infection in Asthma Patients date: 2020-12-14 journal: Med Sci Monit DOI: 10.12659/msm.928861 sha: bf77362e5c445f58a2dec4ae22999bd74ef21320 doc_id: 786399 cord_uid: tgvocrwh BACKGROUND: Rhinovirus (RV) is the most common pathogen involved in asthma, and COVID-19, caused by SARS-COV-2, may be more severe in asthma patients. Here, we applied integrated bioinformatics to identify potential key genes and cytokine pathways after RV infection in asthma, and analyzed changes in angiotensin-converting enzyme 2 (ACE2), the cellular receptor of SARS-COV-2. MATERIAL/METHODS: The gene expression profile dataset GSE149273 was downloaded from NCBI-GEO, which included 90 samples of non-infected, RVA, and RVC. Differentially expressed genes (DEGs) were identified using t tests in the limma R package, and subsequently investigated by GO, KEGG, and DO analysis. Moreover, the expression of ACE2 and the proportion of immune cells were further analyzed to determine the effects of RV on cytokines. RESULTS: A total of 555 DEGs of RVA and 421 of RVC were identified. There were 415 DEGs in RVA and RVC, of which 406 were upregulated and 9 were downregulated. The functional enrichment analysis showed that most DEGs were obviously enriched in cytokines, and were mainly enriched in “influenza” and “hepatitis C, chronic”. In addition, the expression of ACE2 increased significantly and the proportion of immune cytokines significantly changed after RV infection. Our results suggest that RV can activate the cytokine pathway associated with COVID-19 by increasing ACE2. CONCLUSIONS: The DEGs and related cytokine pathways after asthma RV infection identified using integrated bioinformatics in this study elucidate the potential link between RV and COVID-19. Respiratory virus infections are closely associated with the exacerbation of asthma, in which rhinovirus (RV) is the most common pathogen [1] . Coronavirus disease 2019 , caused by severe acute respiratory syndrome coronavirus 2 (SARS-COV-2), may be more serious in patients with chronic lung diseases, including asthma [2] . In asthma patients, biological factors may affect the susceptibility to SARS-COV-2 or the severity of COVID-19. Kim et al. found that RV infection was also the most common additional respiratory pathogen in 116 SARS-COV-2-positive specimens [3] . It is well known that there is a great diversity in RV, in which RV-B infection usually causes mild or asymptomatic symptoms, and RV-A or RV-C are usually associated with more severe illness or exacerbation of asthma [4] . RV infection induces the expression of interferon (IFN)-stimulated genes and subsequent production of cytokines [5] . Angiotensin-converting enzyme 2 (ACE2), a transmembrane protein with an extracellular carboxypeptidase domain located on the cell membrane of a variety of epithelial cells, was originally described as playing a role in the renin-angiotensin system [6] . The various roles of ACE2 have gradually been discovered [7] . It participates in the regulation of cardiovascular physiology, dietary amino acid homeostasis, innate immunity, and intestinal microbial ecology, and plays an important role in heart failure, hypertension, myocardial infarction, and cardiovascular complications of diabetes [8] . It is also a cellular receptor for SARS-COVs, including SARS-COV-2 which caused the COVID-19 epidemic [9] . The spike protein of virus binds to the ACE2 of the host cell, mediating viral invasion of cell membranes and transmission to other cells [10] . Therefore, we speculated that there are differences in ACE2 gene expression in RV-infected asthma patients, and it may be associated with severe COVID19 infections. Gene expression microarray, as an efficient and large-scale technique to obtain genetic data, is widely used to study the gene expression profile of diseases and plays an important role in screening differentially expressed genes (DEGs) [11] . With the wider application of this technology, a large quantity of microarray and sequencing data has been collated and stored in the public database platforms, which are available for free download. Further analysis and integration of the data through relevant statistical software or online analysis tools can lead to more in-depth study of molecular mechanisms of disease, providing a new method for the study of disease-related genes. However, the results of different studies on the identification of DEGs are not completely consistent, and there are still some limitations in single-cohort studies. Many studies have identified the DEGs of different diseases, including ovarian cancer [12] , leukemia [13] , and gastric cancer [14] . The identification of DEGs has been widely used in the analysis of various types of diseases and can provide effective therapeutic targets for the disease. In this study, the original microarray dataset GSE149273 [15] was downloaded from the NCBI-Gene Expression Omnibus (GEO), which contained a total of 90 samples. DEGs were identified in RV-A16 infections (RVA) or RV-C15 infections (RVC) compared to non-infection (control) in patients with asthma. Additionally, gene ontology (GO), Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis, and Disease ontology (DO) analysis were performed to identify the intersecting DEGs of RVA and RVC. The gene expression of ACE2 and the proportion of 22 kinds of immune cells in each sample were further analyzed to determine the effects of FVA or FVC infection on cytokines. Overall, the findings of the present study highlighted key candidate genes and cytokine pathways in RV infection in people with asthma. These may provide a reliable molecular marker to better understand RV infection and analyze its relationship with COVID-19. The GSE149273 profile dataset was downloaded from the NCBI-GEO database (https://www.ncbi.nlm.nih.gov/geo/). The platform for GSE149273 is GPL21290 Illumina HiSeq 3000 (Homo sapiens), which includes 90 samples from 30 asthma patients with non-infection, RVA, and RVC. Platform and series matrix file(s) were downloaded as TXT files. The original data were converted into an expression matrix, and then the mRNA expressions of 90 samples were normalized using the R software package with the robust multi-array average algorithm. The t test method in the limma R package was subsequently used to identify DEGs between RVA and RVC infected samples and non-infected samples, respectively. The log fold change >2 and P<0.05 were considered as the screening criteria for DEGs identification. The intersection of the analysis results of DEGs of RVA and RVC was demonstrated using VennDiagram in the R package, and the heatmap.2 of R package was used to draw the heatmaps of RVC, RVA, and control groups. To further explore the biological processes and signal pathways in which the intersecting DEGs of RVA and RVC might be involved, functional analyses of the DEGs at the intersection e928861-2 of RVA and RVC were performed. GO analysis was performed to demonstrate the dominant function of the intersecting DEGs from molecular function and biological process using the clusterProfiler R package. KEGG pathway analysis was performed to find the correlation between the intersecting DEGs and signal pathways using the clusterprofiler R package. DO analysis was performed to explore the relationship between the intersecting DEGs and diseases using DisGeNET data (https://www.disgenet.org/). All the results were visualized, and P<0.05 and gene counts >5 were considered to show a statistically significant difference in this functional enrichment analysis. The t test was used to analyze the differences in gene expression of ACEs between RVA infection or RVC infection and control. P<0.05 was statistically significant. Moreover, the proportion of 22 immune cells in each sample was estimated to analyze the infiltration of immune cells after infection using cibersortx (https://cibersortx.stanford.edu/index.php), and a heatmap was drawn. The vioplot in the R package was used to visualize the change in the proportion of 22 immune cells relative to the control after RVA or RVC infection. P<0.05 was considered statistically significant. The GSE149273 profile dataset consisted of 30 patients with a median age of 35 years and a balanced sex ratio, 60% of whom were non-Hispanic white individuals. Thirty individuals were treated with non-infection, RVA, and RVC infection respectively, and 90 samples were obtained. The mRNA expression levels of the 90 samples in the microarray dataset GSE149273 were standardized, and the results are shown in Figure 1 . The GSE149273 dataset was screened using the limma R package (logFC >2 and P<0.05), and 555 DEGs were obtained between RVA and control and 421 DEGs between RVC and control. Among them, 529 upregulated genes and 26 downregulated genes were identified in RVA, while 408 upregulated genes and 13 downregulated genes were classified in RVC. In addition, there were a total of 415 DEGs both in RVA and RVC, among which 406 were upregulated and 9 were downregulated ( Table 1 ). The differential expression of multiple genes from FVA or FVC and control in the microarray is shown in Figure 2A and 2B, and the intersection of RVA and RVC DEGs is shown in Figure 2C . The cluster heatmaps of FVA and FVC are shown in Figure 3 . To identify which biological functions and pathways were most significantly associated with the DEGs, a total of 415 genes in the intersection of RVA and RVC were analyzed for GO and KEGG pathways. GO analysis of DEGs in this study was divided into 2 functional groups: a biological processes group and a molecular function group ( Figure 4A ). In the biological processes group, the DEGs were mainly enriched in "response to virus", "defense response to virus", and "response to interferon-gamma". In the molecular function group, the DEGs were mainly enriched in "receptor ligand activity", "cytokine activity", and "cytokine receptor binding". The top 10 GO terms for each category are listed in Table 2 . KEGG pathway analysis demonstrated that DEGs were significantly enriched in "Cytokine-cytokine receptor interaction", "Influenza A", and "NOD-like receptor signaling pathway" ( Figure 4B ). All KEGG pathway terms are listed in Table 3 . These results indicate that most DEGs were significantly enriched in virus and cytokine. To further clarify the relationship between DEGs and various diseases, 415 DEGs were subjected to DO analysis. The specific conditions of diseases enriched by each gene are shown in Figure 5 . After enrichment analysis, it was found that the DEGs were mainly enriched in "Influenza", "Hepatitis C, Chronic", and "Respiratory Syncytial Virus Infections". The enrichment e928861-8 ICAM1 IDO1 IFI27 IFI35 IFI44 IFI44L IFI6 IFIH1 IFIT1 IFIT2 IFIT3 IFIT5 IFITM3 IFNL1 IFNL3 IL6 IL7 IRF7 IRF9 ISG15 LAMP3 LY6E MX1 MX2 MYD88 NOS2 OAS1 OAS2 OAS3 OASL PARP12 PLSCR1 PTGS2 RSAD2 RTP4 SAMD9 SCO2 SERPING1 SOCS1 SOCS3 STAT1 STS TDRD7 TICAM1 TLR2 TLR3 TNF TNFSF10 TOR1B TREX1 TRIM22 TRIM25 UBE2L6 XAF1 ZCCHC2 CCL22 CSF3 CXCL13 CXCL3 IFITM1 IL15RA IL19 IL1RN IL2RG IRF1 JAK2 results are shown in Figure 6A , and the interaction network between diseases and DEGs is shown in Figure 6B . The ACE2 gene is one of the 415 DGEs, and this study further confirmed that the expression of ACE was significantly increased after RVA or RVC infection (P<0.01), along with upregulated DGEs (Figure 7) . Functional enrichment analysis showed that most DEGs were significantly related to cytokines, so we analyzed the content and proportion of 22 kinds of immune cells in 90 samples. The content of immune cells is shown in Figure 8A and the heatmap of cell proportions is shown in Figure 8B . The results showed that the proportion of "NK cells activated", "Macrophages M1", and "Dendritic cells activated" was increased after infection with RVA or RVC, while that of "T cells CD8" and "Mast cells resting" was decreased. Based on the correlation analysis of the immune cells, we found that "T cells follicular helper", "Dendritic cells activated", "Neutrophils", and "NK cells activated" had a positive correlation with "Macrophages M1", and "T cells CD8" was also positively correlated with "Dendritic cells resting" ( Figure 8C ). To identify changes in the expression of individual immune cells after infection, we visualized the proportion of 22 kinds of immune cells relative to control after RVA or RVC infection. The results showed that the proportions of "B cells naïve", "T cells CD8", "NK cells resting", "Macrophages M0", "Macrophages M2", "Dendritic cells resting", and "Mast cells resting" were significantly decreased, while that of "T cells follicular helper", "NK cells activated", "Macrophages M1", "Dendritic cells activated" and "Neutrophils" increased significantly after RVA infection ( Figure 9A ). And the proportion of "B cells naïve", "T cells CD8", "NK cells resting", "Macrophages M0", "Dendritic cells resting", and "Mast cells resting" was significantly decreased, while "T cells gamma delta", "NK cells activated", "Macrophages M1", "Dendritic cells activated", and "Neutrophils" was significantly increased after RVC infection ( Figure 9B ). In general, RV infection significantly activated the expression of NK cells, dendritic cells, and neutrophils, and reduced the expression of CD8 T cells. SARS-COV-2, causing COVID-19, is a novel virus first detected in December 2019 in Wuhan, China and the disease is thought to be more severe in patients with chronic lung disease [16] . RV airway infection is the most common virus causing the common cold and is a major contributor to the development and exacerbation of asthma [17] . RV infection induces the expression of IFN-stimulated genes and subsequent cytokines [18] , and ACE2, the receptor of SARS-COV-2, has recently been identified as an IFN-stimulated gene [19] . Therefore, RV infection may be associated with severe COVID-19 infection through key candidate genes and cytokine pathways. Studies have shown that there are large clinical variations due to differences in autoantibodies among individuals in the course of SARS-COV-2 infection, and type I IFNs in protective immunity against SARS-Cov-2 are crucial [20] . There are many immunological tests for SARS-COV-2 infection, but the protective immune genes for SARS-CoV- 2 have not yet been fully clarified, and the study of genes and cytokines can provide a unique opportunity to control the infection and disease [21] . Microarrays and high-throughput sequencing have been widely used to detect the expressions of various genes in the human genome and to predict potential targets for related diseases, which can help identify target genes for diagnosis or treatment of diseases [22] . In the present study, the genomic dataset was obtained from GEO, including 90 samples with non-infection, RVA infection, and RVC infection. Bioinformatics analysis was conducted, identifying 555 DEGs between RVA and control and 421 DEGs between RVC and control. A total of 415 DEGs were involved both in RVA and RVC, among which 406 were upregulated and 9 were downregulated. Functional enrichment analysis showed that DEGs were mainly enriched in "response to e928861-9 C virus", "response to interferon-gamma", "receptor ligand activity", "cytokine activity", "cytokines-cytokine receptor interaction", "NOD-like receptor signaling pathway", "Influenza A", and "respiratory syncytial virus infection". These results indicate that most DEGs were significantly enriched in virus and cytokine, which is consistent with normal expression after virus infection. SARS-COV-2 infection is mediated by a Spike virus protein and ACE2, and fuses with target cells through a type I mechanism [23] . ACE2 is expressed in the alveoli, heart, esophagus, kidney, bladder, and ileum, suggesting that all of these organs may be affected by SARS-COV-2 [24, 25] . The expression of ACE2 mediates SARS-COV-2 infection in host lung cells [10] , so it can be reasonably concluded that increased expression of ACE2 in lung cells will increase susceptibility to SARS-COV-2 infection or lead to more severe COVID19 disease. Generally, after cellular detection of viral entry into a host cell, IFN-induced IFNstimulated genes are essential for host antiviral defense, and studies have shown that ACE2 is an INF-stimulated gene in human airway epithelial cells [19] . Inducers of the IFN pathway, including influenza infection, can lead to increased expression of ACE2 in cells, thereby increasing susceptibility to COVID-19. Our study found that ACE2 was one of the 415 DEGs in RV and control, and there was a significant increase in the expression of ACE2 in asthma patients infected with RV, indicating that RV infection can increase susceptibility to SARS-COV-2 or lead to more severe COVID19. Immune cells are the cells involved in or related to the immune response, and they play an important role in the human body. Studies have shown that NK cells are activated in e928861-12 COVID-19 [26] , and can also be activated in humans after influenza virus infection, while dendritic cells are necessary for NK cell activation [27] . Dendritic cells play an important role in defense against lung infections [28] . In human dendritic cells derived from monocytes infected by SARS-COV, the chemokines are upregulated [29] . Furthermore, lymphocytosis may be due to the direct effects of cytokines on T cell populations and/or the indirect effects on other cell types such as dendritic cells [30] and neutrophils [31, 32] in severe disease. Our research found that NK cells, dendritic cells, and neutrophils were activated after RV infection in asthma, which was consistent with the previous description. In addition, we found a e928861-13 significant decrease in CD8 T cells of people with asthma infected with RV, which was similar to earlier observations about SARS-CoV-1 infection. Some studies reported that lymphopenia was caused by a sharp decrease in the number of T cells CD4 and CD8 in COVID-19 patients [33] . T cells play an important role in viral infections, in which CD4 T cells help to produce antibodies and coordinate the response of other immune cells, while CD8 kills infected cells to reduce the viral burden, and lower numbers of CD8 cells have been shown to be correlated with severity of COVID-19 [34] . The hyperactivation of monocytes, macrophages, and dendritic cells is the driving factor of human immunopathology. Our study found a significant correlation between CD8 and dendritic cells after RV infection, indicating that asthma patients infected with RV would be more susceptible to COVID-19, which was also consistent with the higher incidence of COVID-19 among young asthmatics, accounting for 27% of the hospitalized 18-49 age group in the United States [35] . All results showed that viral infections associated with asthma exacerbations and SARS-COV-2 infections showed synergistic biomolecular interactions. RV infection can aggravate COVID-19 infection by increasing key gene and cytokine pathways associated with COVID-19. Rhinovirus and asthma exacerbations Risk factors associated with acute respiratory distress syndrome and death in patients with coronavirus disease 2019 pneumonia in Wuhan, China Rates of co-infection between SARS-CoV-2 and other respiratory pathogens Human rhinovirus species and season of infection determine illness severity A comprehensive evaluation of nasal and bronchial cytokines and chemokines following experimental rhinovirus infection in allergic asthma: increased interferons (IFN-g and IFN-l) and type 2 inflammation (IL-5 and IL-13) A novel angiotensin-converting enzyme-related carboxypeptidase (ACE2) converts angiotensin I to angiotensin 1-9 SARS-CoV-2 receptor and regulator of the renin-angiotensin system: Celebrating the 20 th anniversary of the discovery of ACE2 Role of the ACE2/angiotensin 1-7 axis of the renin-angiotensin system in heart failure Structural basis for the recognition of SARS-CoV-2 by full-length human ACE2 SARS-CoV-2 cell entry depends on ACE2 and TMPRSS2 and is blocked by a clinically proven protease inhibitor Expression Atlas update -a database of gene and transcript expression from microarray-and sequencing-based functional genomics experiments Identification of differentially expressed genes and signaling pathways in ovarian cancer by integrated bioinformatics analysis These DEGs were mainly related to cytokines. Our results contribute to a better understanding of the molecular mechanism of RV infection. We found that there was significant activation of NK cells, dendritic cell, and neutrophils, and a significant decrease in CD8 T cells after RV infection in people with asthma. Moreover, the expression of ACE2, which is closely related to COVID-19, was significantly increased after RV infection. The evidence suggests that RV infection activates the cytokine pathways and increases the expression of ACE2, thus exacerbating COVID-19. These results enhance our understanding of RV infection and its correlation with COVID-19. However Integrated computational biology analysis to evaluate target genes for chronic myelogenous leukemia Investigating potential molecular mechanisms of carcinogenesis and genes as biomarkers for prognosis of gastric cancer based on integrated bioinformatics analysis Rhinovirus infections in individuals with asthma increase ACE2 expression and cytokine pathways implicated in COVID-19 Clinical characteristics of coronavirus disease 2019 in China Role of viral respiratory infections in asthma and asthma exacerbations Human rhinovirus-induced epithelial production of CXCL10 is dependent upon IFN regulatory factor-1 SARS-CoV-2 receptor ACE2 is an interferon-stimulated gene in human airway epithelial cells and is detected in specific cell subsets across tissues Autoantibodies against type I IFNs in patients with life-threatening COVID-19 A global effort to define the human genetics of protective immunity to SARS-CoV-2 infection Molecular mechanisms of breast cancer metastasis by gene expression profile analysis ACE2, the receptor that enables infection by SARS-CoV-2: Biochemistry, structure, allostery and evaluation of the potential development of ACE2 modulators Single-cell RNA-seq data analysis on the receptor ACE2 expression reveals the potential risk of different human organs vulnerable to 2019-nCoV infection Single-cell RNA expression profiling of ACE2, the receptor of SARS-CoV-2. bioRxiv Natural killer cell activation related to clinical outcome of COVID-19. medRxiv NKp46 and NKG2D recognition of infected dendritic cells is necessary for NK cell activation in the human response to influenza infection Division of labor between lung dendritic cells and macrophages in the defense against pulmonary infections Chemokine up-regulation in SARScoronavirus-infected, monocyte-derived human dendritic cells Interleukin 10 (IL-10) and viral IL-10 strongly reduce antigen-specific human T cell proliferation by diminishing the antigen-presenting capacity of monocytes via downregulation of class II major histocompatibility complex expression Neutrophil-to-lymphocyte ratio as an independent risk factor for mortality in hospitalized patients with COVID-19 Predictive value of the neutrophil-to-lymphocyte ratio (NLR) for diagnosis and worse clinical course of the COVID-19: Findings from ten provinces in China Immunology of COVID-19: Current state of the science Metabolic disturbances and inflammatory dysfunction predict severity of coronavirus disease 2019 (COVID-19): A retrospective study. medRxiv Hospitalization rates and characteristics of patients hospitalized with laboratory-confirmed coronavirus disease 2019 -COVID-NET, 14 states