key: cord-0773695-4h5nj8w4 authors: Xie, Li-Min; Huang, Yin-Fei; Liu, Ye-Ling; Liang, Jia-Qi; Deng, Wei; Lin, Geng-Ling; Luo, Huan-Min; Guo, Xu-Guang title: Identification of the Hub Genes and the Signaling Pathways in Human iPSC-Cardiomyocytes Infected by SARS-CoV-2 date: 2022-03-02 journal: Biochem Genet DOI: 10.1007/s10528-022-10206-7 sha: 6e4e84e56fee81d34852c2f7444aacc5e98304fc doc_id: 773695 cord_uid: 4h5nj8w4 Severe Acute Respiratory Syndrome Coronavirus Type 2 (SARS-CoV-2) is an enveloped single-stranded RNA virus that can lead to respiratory symptoms and damage many organs such as heart, kidney, intestine, brain and liver. It has not been clearly documented whether myocardial injury is caused by direct infection of cardiomyocytes, lung injury, or other unknown mechanisms. The gene expression profile of GSE150392 was obtained from the Gene Expression Omnibus (GEO) database. The processing of high-throughput sequencing data and the screening of differentially expressed genes (DEGs) were implemented by R software. The R software was employed to analyze the Gene Ontology (GO) analysis and the Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis. The protein–protein interaction (PPI) network of the DEGs was constructed by the STRING website. The Cytoscape software was applied for the visualization of PPI network and the identification of hub genes. The statistical analysis was performed by the GraphPad Prism software to verify the hub genes. A total of 516 up-regulated genes and 191 down-regulated genes were screened out. The top 1 enrichment items of GO in biological process (BP), Cellular Component (CC), and Molecular Function (MF) were type I interferon signaling pathway, sarcomere, and receptor ligand activity, respectively. The top 10 enrichment pathways, including TNF signaling pathway, were identified by KEGG enrichment analysis. A PPI network was established, consisting of 613 nodes and 3,993 edges. The 12 hub genes were confirmed as statistically significant, which was verified by GSE151879 dataset. In conclusion, the hub genes of human iPSC-cardiomyocytes infected with SARS-CoV-2 were identified through bioinformatics analysis, which may be used as biomarkers for further research. Severe Acute Respiratory Syndrome Coronavirus Type 2 (SARS-CoV-2) is an enveloped single-stranded RNA virus which belongs to the family of coronaviruses and the severe acute respiratory syndrome-related coronavirus (Gorbalenya et al. 2020) . As the seventh known coronavirus, it can be transmitted from person to person (Anderson et al. 2020) . According to the latest data from Johns Hopkins University in the United States, as of 6:35 a.m. Beijing time, July 23rd, 2020, the cumulative number of confirmed cases and deaths of new coronavirus disease 2019 worldwide has exceeded 15 million and 620,000, respectively. Currently, the number is still rising rapidly. The respiratory secretions or droplets of infected people, if inhaled, will cause infection in healthy people. SARS-CoV-2 can spread among people through close contact (Anderson et al. 2020; Shereen et al. 2020; Xu et al. 2020) . Fever and cough are frequently seen in patients, but some patients may develop myalgia and fatigue Huang et al. 2020) . Although the patients mainly show the symptoms of the respiratory system, in fact, the virus can damage many organs at the same time, such as heart, kidney, intestine, brain and liver (Dhakal et al. 2020) . It is reported that 23% of severe patients were diagnosed with cardiac injury, a common complication that would exacerbate the disease in its severity and progressing, which includes palpitations or chest pain, myocarditis, acute myocardial infarction (AMI), heart failure (HF), arrhythmia, as well as venous thromboembolic events (VTE) Chen et al. 2020; Long et al. 2020; Rizzo et al. 2020 ). Among them, the AMI and arrhythmia can be seen in patients with COVID-19, and it is more common in ICU patients (Bansal 2020) . The receptor-binding domain (RBD) of SARS-CoV-2 might bind tightly with human ACE2, which was strongly expressed in pericytes of human hearts, indicating the potential susceptibility of cardiomyocytes to SARS-CoV-2 (Chen et al. 2020) . Additionally, the expression of ACE2 was highly increased in the heart failure patient (Chen et al. 2020 ). In the direct mechanism, the myocardial injury can be caused by the infiltration of viruses into myocardial tissues, resulting in the death and inflammation of cardiac muscles. In the indirect mechanism, cardiomyocyte stress secondary to pulmonary failure, hypoxemia, and other cardiac inflammation can lead to severe systemic inflammatory response syndrome (SIRS) (Akhmerov and Marbán 2020) . The preliminary clinical studies have developed several drugs of therapeutic potential, such as antiviral drugs, nucleoside analogues, neuraminidase inhibitors, therapeutic peptide, RNA synthesis inhibitors, anti-inflammatory drugs, chloroquine phosphate, Lianhuaqingwen capsules and ShuFengJieDu capsules as well (Chakraborty et al. 2020; Lu 2020; Gao et al. 2020) . The high-throughput sequencing technology provides the information of genes, pathophysiological mechanism, diagnosis and therapy, as well as the method of antiviral prevention and vaccination development (Chan et al. 2020) . It not only shortens the experimental cycle and reduces the cost of experiments, but also significantly improves the application value of species (such as SARS-CoV-2), the databases of which need to be improved due to the scarcity of their biological information (Pareek et al. 2011) . The exact mechanism of myocardial damage caused by COVID-19 remains unclear. In this study, the gene expression data of human-induced pluripotent stem cell-derived cardiomyocytes (iPSC) infected by SARS-CoV-2 was analyzed through bioinformatics analysis, the result of which is hopefully helpful for studying the mechanism of the infection of myocardium cells by SARS-CoV-2. It is of great significance for the establishment of cardiac antiviral drug screening platform (Bojkova et al. 2020) , and helpful to promote the development of clinical targeted drug therapy. The aim of our study is identifying hub genes and signaling pathways that SARS-CoV-2 infected cardiomyocytes, which may be available to screen appropriate and effective biomarkers for further research. The gene expression profiles of the GSE150392 dataset was downloaded from the Gene Expression Omnibus (GEO, https:// www. ncbi. nlm. nih. gov/ geo/) database of the National Center for Biotechnology Information (NCBI). The database was designed to stockpile primitive sequences, gene expression datasets and platform records. Based on the platform of [GPL18573 Illumina NextSeq 500 (Homo sapiens)], the dataset contains six human iPSC-cardiomyocytes samples, among which the GSM4548303-GSM4548305 served as the infection group while the GSM4548306-GSM4548308 as the control group. The R software was applied to calibrate and standardize the high-throughput sequencing data of the GSE150392 dataset that was downloaded from the GEO database. The limma package of the R software was employed to screen differentially expressed genes (DEGs) with conditions of |log2FC|≥ 2 and P < 0.05, the pheatmap package to construct the heat-map of DEGs, and the ggplot2 package (Ito and Murphy 2013) to construct the volcano map. The function and pathway enrichment analysis was employed to explain the biological functions of DEGs. The GO analysis, consisting of three categories of cellular component (CC), biological process (BP), molecular function (MF), was applied to describe the characteristics of genes and gene products (Fan and Wei 2020) . The KEGG analysis, a public integrated database of biological systems designed for the systematic analysis of gene functions, rendered it possible to understand the pathways which have changed significantly under experimental conditions. The clusterProfiler package in the R software was applied to conduct the GO and KEGG analysis of DEGs, the results of which would be regarded as statistically significant if P < 0.05. The DEGs were uploaded to the STRING website (Szklarczyk et al. 2017 ) (http:// string-db. org/) to establish a PPI network that aimed to analyze the relationship between proteins in DEGs. The visualization and hub-genes screening were performed by Cytoscape software (Su et al. 2014 ), on which the cytoHubba plug-in was installed to calculate the score of each gene. According to Maximal Clique Centrality (MCC), the 12 genes with the highest degree were identified as hub genes. To confirm the findings from the bioinformatics analysis, the single-cell RNA sequencing's data of the peripheral blood mononuclear cells (PBMCs) from the GSE150728 dataset and the high-throughput sequencing's data of human embryonic stem cell-derived cardiomyocytes (hESC-CMs) from the GSE151879 dataset were applied to verify the hub genes. The GraphPad Prism software was used to analyze the data from the GSE150728 dataset and GSE151879 dataset, and P < 0.05 indicated the significant difference. The microRNAs (miRNAs), short non-coding single-stranded RNAs, play a part in pathogenesis of diseases in which they can control a variety of biological functions by regulating gene expression (Chou et al. 2018; Mitash et al. 2020 ; Saçar Demirci and Adan 2020). After uploading the gene list to the Network Analyst website, the miRTarBase v8.0 database was used to perform network analysis (Chou et al. 2018) , which collects information about miRNA-target interactions (MTI) verified by experiment. Totally 707 DEGs have been screened out according to the criteria of |log2FC|≥ 2 and P < 0.05, including 516 up-regulated genes and 191 down-regulated genes, the results of which were presented in the form of volcano map and heat map (Fig. 1A , B) with the help of ggplot2 package and pheatmap package in the R software. The GO enrichment analysis denoted that the parts playing a role of significance in BP section are response to type I interferon, type I interferon signaling pathway, cellular response to type I interferon, etc. ( Fig. 2A) . As for CC section, sarcomere, contractile fiber part, contractile fiber, I-kappa B/NF-kappa B complex, as well as myofibril rank top 5 (Fig. 2B) . Additionally, the MF part presents that the top enrichment items are receptor ligand activity, receptor regulator activity, cytokine receptor-binding, DNA-binding transcription activator activity, RNA polymerase IIspecific, and cytokine activity (Fig. 2C ). More details of the top 5 significant items in BP, CC, and MF of GO enrichment can be seen in Table 1 . The biological functions of DEGs were identified through the KEGG pathway enrichment analysis. The result of the analysis of is shown in Fig. 2D . According to the P-value, the top 10 enrichment pathways were obtained, Table 2 . A PPI network was constructed in the STRING website, which consists of 82 nodes and 213 edges. By visualizing and screening hub genes with Cytoscape software, a PPI network and top 12 genes were obtained and presented in Fig. 3A , B. The top 12 genes were: ISG15, MX1, IFIT1, IFIT2, RSAD2, OAS2, IFIT3, MX2, OAS1, OAS3, IFI6 and DDX58, which were confirmed as statistically significant verified by the GSE151879 dataset (Fig. 4A) . The relative expression levels of 10 hub genes, namely ISG15, MX1, IFIT1, RSAD2, OAS2, IFIT3, MX2, OAS1, OAS3 and IFI6, were also consistent with single-cell RNA sequencing of PBMCs from the GSE150728 dataset, while DDX58 and IFIT2 showed no statistically significant difference (Fig. 4B) . After network visualization, Gene-miRNA interactions were presented in Fig. 5 . According to the criterion of the degree, the miRNA-hsa-mir-335-5p was finally screened out. The PPI network of DEGs and the top 12 hub genes of PPI network Fig. 4 A Verification of hub genes with GSE151879 dataset. P-value < 0.05 is considered to be statistically significant. B Verification of hub genes with GSE150728 dataset. P-value < 0.05 is considered to be statistically significant In this study, the gene expression profile of the GSE150392 dataset was downloaded from the GEO database, and 707 DEGs were identified, including 516 up-regulated DEGs and 191 down-regulated DEGs. In addition, the GO enrichment analysis and KEGG pathway analysis were applied to analyze the functional enrichment of DEGs. A PPI network was applied to identify the hub genes that might play a significant role in the infection of cardiomyocytes by SARS-CoV-2. As presented in the figure that described the results of GO enrichment analysis of DEGs, it is believed that the smaller the P-value is, the more significant the extent of enrichment in the GO analysis is. As for biological processes (BP), DEGs are most significantly enriched in response to type I interferon, type I interferon signaling pathway and cellular response to type I interferon. Type I interferon has a broad-spectrum antiviral feature against RNA viruses, which may lead to innate and adaptive immune response to function against viruses. Besides, IFN-α promotes the expression of ACE2 that is strongly expressed in human hearts. The RBD of SARS-CoV-2 binding ACE2 leads to infection, while SARS-CoV-2 can employ speciesspecific interferon to drive up-regulation of ACE2 to promote infection (Chen et al. 2020; Ziegler et al. 2020; Mantlo et al. 2020) . In terms of cellular components (CC), The direct damage of cardiomyocytes, cardiac inflammation, immune response, myocardial interstitial fibrosis, and hypoxia caused by SARS-CoV-2 might act on sarcomere and contractile fiber probably, all of which may result in functional changes in the sarcomere and contractile fiber, and gradually evolve into myocardial injury (Babapoor-Farrokhran et al. 2020) . In the respect of molecular function (MF), receptor ligand activity and receptor regulator activity are most significantly enriched in DEGs. The ACE2 appears to be the crucial functional receptor of SARS-CoV-2 (Coutard et al. 2020; Walls et al. 2020) . Interaction between viral protein ligand and cell receptor of SARS-CoV-2, a critical step in replication of viruses, is also a crucial process of virus infection, which might be determined by the type of chemical interactions, activity between receptor and ligand, and the regulator activity (Ortega et al. 2020) . The KEGG pathway enrichment analysis discovered that quite a few significantly enriched KEGG pathways rated by P-value, such as TNF signaling pathway, NF-kappa B (NF-κB) signaling pathway, and Legionellosis, are highly relevant to the infection of SARS-CoV-2. As for the TNF signaling pathway, TNF is the crucial inflammatory cytokine resulting in various acute and chronic inflammations (Soy et al. 2020) . TNF-α level was observed much higher in serum of patients with COVID-19, implying that it is positively correlated with the severity of disease (Soy et al. 2020) . With the progress of SARS-CoV-2 infection, monocytes and macrophages increased, resulting in cytokine storm, which meant the release of cytokine and pro-inflammatory cytokines, including IL-1, IL-8, and TNF (Soy et al. 2020; Runfeng et al. 2020) . The inflammatory response triggered by the SARS-CoV-2 in myocardium drives an up-regulation of TNF signaling pathway. In terms of the following pathway of NF-kappa B (NF-κB), it is involved in inflammation, innate and adaptive immune response, pathological development of tumor or cancer, and even the critical regulators of stress responses, apoptosis and differentiation (Baldwin 2001; DiDonato et al. 2012; Oeckinghaus et al. 2011) . SARS-CoV-2 triggered the stimulation of humoral and cellular immune response as well as mitogen-activated protein kinase (MAPK) pathway and NFkappa-B signaling pathway, and gradually caused deregulation of immune system and up-regulation of inflammatory pathways, especially NF-kappa B (Mozafari et al. 2020; Wu and Yang 2020; Li et al. 2020) . The higher expression of NFkappa B conducted the induction of the naïve T cell activation and proliferation, promoted the expression of chemokine expression and other immune response, and finally resulted in inflammatory cytokine storm in myocardium (Mozafari et al. 2020; Hoesel and Schmid 2013) . The 12 hub genes of cardiomyocytes specifically infected by SARS-COV-2 were identified, namely ISG15, MX1, IFIT1, IFIT2, RSAD2, OAS2, IFIT3, MX2, OAS1, OAS3, IFI6 and DDX58, which may become the biomarkers. Among the 12 genes, Anna et al. (Rahnefeld et al. 2014) found that ISG15, as a part of innate immunity of cardiomyocytes, plays a significant role in the inhibition of viral replication in mouse cardiomyocytes infected by Coxsackievirus B3 (CVB3). Besides, SARS-CoV-2 can induce a large number of ISG expressions, either as effectors, regulators or both, which is the antiviral response of the body (Clemente et al. 2020) . It is suggested that ISG15 may play a crucial role in cardiomyocytes infected by SARS-CoV-2. As for another hub gene (MX1), Eshwar et al. have shown that the β 4 subunit of Cav1.2 channels can promote the expression of IFN-related genes (including MX1) in cardiomyocytes, thus alleviating viral infection (Tammineni et al. 2018 ). IFIT1, IFIT2, IFIT3 and IFIT5 constitute the human IFIT gene family, which expressed very low in most cell types but can be significantly increased by interferon therapy, viral infection and pathogen-related molecular pattern (PAMP) (Pidugu et al. 2019 ). SARS-CoV-2 can induce strong antiviral response by up-regulating antiviral factors such as OAS1-3, IFIT1-3 and CXCL9 / 10 / 11 of Th1 chemokine, and reduce ribosomal protein transcription (Lieberman et al. 2020) . Hsa-miR-335-5p regulates significant transcription factors and participates in activating the pathways affiliated with inflammation and angiogenesis. It may activate WNT and TGF-β signaling pathways, both of which exert pleiotropic and multifunctional effects, regulating a considerable number of biological processes, TGF-β and WNT signaling pathways (Kay et al. 2019; Esquinas et al. 2017 ). In the pathogenesis of cardiac remodeling and fibrosis, TGF-β plays an important role (Dobaczewski et al. 2011) . The gene expression profiles of normal human iPSC cardiomyocyte samples and human iPSC cardiomyocyte samples infected by SARS-CoV-2 were studied through bioinformatics analysis, and finally, the hub genes of human iPSC-cardiomyocytes after SARS-CoV-2 infection were identified. Although sequencing data from other datasets were used for verification, further experimental research analysis (such as qRT-PCR) still need to be conducted to acquire abundant data to verify the predicted results obtained from the bioinformatics analysis. There are still some limitations to this study. It failed to understand clearly the mechanism of several dominant genes being involved in the infection of cardiomyocytes by SARS-CoV-2. The other limitation is that the small sample size and data sets of this study would affect the accuracy of the conclusions. Further research is entailed to demonstrate the molecular mechanism of these biological events involved in the infection of cardiomyocytes by the novel coronavirus. In conclusion, the RNA-Seq data of human iPSC-cardiomyocytes infected by SARS-CoV-2 was studied through bioinformatics analysis. Finally, the hub genes associated with SARS-CoV-2 infection were identified, which may reveal the mechanism of heart damage caused by SARS-CoV-2 or suggest biomarkers. COVID-19 and the heart Consideration of the aerosol transmission for COVID-19 and public health Myocardial injury and COVID-19: possible mechanisms Series introduction: the transcription factor NF-kappaB and human disease 2020) Cardiovascular disease and COVID-19 SARS-CoV-2 2 infects and induces cytotoxic effects in human cardiomyocytes 2020 SARS-CoV-2 causing pneumonia-associated respiratory disorder (COVID-19): diagnostic and proposed therapeutic options Genomic characterization of the 2019 novel human-pathogenic coronavirus isolated from a patient with atypical pneumonia after visiting Wuhan The ACE2 expression in human heart indicates new potential mechanism of heart injury among patients infected with SARS-CoV-2 miRTarBase update 2018: a resource for experimentally validated microRNA-target interactions Deubiquitinating enzymes in coronaviruses and possible therapeutic opportunities for COVID-19 The spike glycoprotein of the new coronavirus 2019-nCoV contains a furin-like cleavage site absent in CoV of the same clade SARS-CoV-2 infection and cardiovascular disease: COVID-19 heart NF-κB and the link between inflammation and cancer Transforming growth factor (TGF)-β signaling in cardiac remodeling Gene and miRNA expression profiles in PBMCs from patients with severe and mild emphysema and PiZZ alpha1-antitrypsin deficiency Identification of novel potential biomarkers and therapeutic targets involved in human atrial fibrillation based on bioinformatics analysis Breakthrough: Chloroquine phosphate has shown apparent efficacy in treatment of COVID-19 associated pneumonia in clinical studies Severe acute respiratory syndrome-related coronavirus: The species and its viruses-a statement of the Coronavirus Study Group The complexity of NF-κB signaling in inflammation and cancer Clinical features of patients infected with 2019 novel coronavirus in Wuhan Application of ggplot2 to pharmacometric graphics Hsa-miR-335 regulates cardiac mesoderm and progenitor cell differentiation Coronavirus infections and immune responses In vivo antiviral host response to SARS-CoV-2 by viral load, sex, and age Cardiovascular complications in COVID-19 Drug treatment options for the 2019-new coronavirus (2019-nCoV) Antiviral activities of type I interferons to SARS-CoV-2 infection The role of MicroRNA in the airway surface liquid homeostasis Inflammation: a bridge between diabetes and COVID-19, and possible management with sitagliptin Crosstalk in NF-κB signaling pathways Role of changes in SARS-CoV-2 spike protein in the interaction with the human ACE2 receptor: an in silico analysis Sequencing technologies and genome sequencing Emerging functions of human IFIT proteins in cancer Ubiquitin-like protein ISG15 (interferon-stimulated gene of 15 kDa) in host defense against heart failure in a mouse model of virus-induced cardiomyopathy COVID-19 in the heart and the lungs: could we "Notch" the inflammatory storm? Lianhuaqingwen exerts anti-viral and anti-inflammatory activity against novel coronavirus (SARS-CoV-2) Computational analysis of microRNA-mediated interactions in SARS-CoV-2 infection COVID-19 infection: origin, transmission, and characteristics of human coronaviruses Cytokine storm in COVID-19: pathogenesis and overview of anti-inflammatory agents used in treatment Biological network exploration with Cytoscape 3 The STRING database in 2017: quality-controlled protein-protein association networks, made broadly accessible The β(4) subunit of Ca(v)12 channels is required for an optimal interferon response in cardiac muscle cells Structure, Function, And Antigenicity of the SARS-CoV-2 spike glycoprotein TH17 responses in cytokine storm of COVID-19: an emerging target of JAK2 inhibitor Fedratinib Imaging and clinical features of patients with 2019 novel coronavirus SARS-CoV-2 Clinical course and outcomes of critically ill patients with SARS-CoV-2 pneumonia in Wuhan, China: a single-centered, retrospective, observational study SARS-CoV-2 receptor ACE2 is an interferon-stimulated gene in human airway epithelial cells and is detected in specific cell subsets across tissues Publisher's Note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations Funding No funding was used to support this study.Data availability All data generated or analyzed during this study are included in this published article. Conflict of interest The authors declare that they have no competing interests. Li-Min Xie 1,2 · Yin-Fei Huang 1,2 · Ye-Ling Liu 1,2 · Jia-Qi Liang 1,2 · Wei Deng 1,2 · Geng-Ling Lin 1,2 · Huan-Min Luo 1,2 · Xu-Guang Guo 1,2,3,4 * Xu-Guang Guo gysygxg@gmail.com 1