Network-Based Data Analysis Reveals Ion Channel-Related Gene Features in COVID-19: A Bioinformatic Approach.

Coronavirus disease 2019 (COVID-19) seriously threatens human health and has been disseminated worldwide. Although there are several treatments for COVID-19, its control is currently suboptimal. Therefore, the development of novel strategies to treat COVID-19 is necessary. Ion channels are located on the membranes of all excitable cells and many intracellular organelles and are key components involved in various biological processes. They are a target of interest when searching for drug targets. This study aimed to reveal the relevant molecular features of ion channel genes in COVID-19 based on bioinformatic analyses. The RNA-sequencing data of patients with COVID-19 and healthy subjects (GSE152418 and GSE171110 datasets) were obtained from the Gene Expression Omnibus (GEO) database. Ion channel genes were selected from the Hugo Gene Nomenclature Committee (HGNC) database. The RStudio software was used to process the data based on the corresponding R language package to identify ion channel-associated differentially expressed genes (DEGs). Based on the DEGs, Gene Ontology (GO) functional and pathway enrichment analyses were performed using the Enrichr web tool. The STRING database was used to generate a protein-protein interaction (PPI) network, and the Cytoscape software was used to screen for hub genes in the PPI network based on the cytoHubba plug-in. Transcription factors (TF)-DEG, DEG-microRNA (miRNA) and DEG-disease association networks were constructed using the NetworkAnalyst web tool. Finally, the screened hub genes as drug targets were subjected to enrichment analysis based on the DSigDB using the Enrichr web tool to identify potential therapeutic agents for COVID-19. A total of 29 ion channel-associated DEGs were identified. GO functional analysis showed that the DEGs were integral components of the plasma membrane and were mainly involved in inorganic cation transmembrane transport and ion channel activity functions. Pathway analysis showed that the DEGs were mainly involved in nicotine addiction, calcium regulation in the cardiac cell and neuronal system pathways. The top 10 hub genes screened based on the PPI network included KCNA2, KCNJ4, CACNA1A, CACNA1E, NALCN, KCNA5, CACNA2D1, TRPC1, TRPM3 and KCNN3. The TF-DEG and DEG-miRNA networks revealed significant TFs (FOXC1, GATA2, HINFP, USF2, JUN and NFKB1) and miRNAs (hsa-mir-146a-5p, hsa-mir-27a-3p, hsa-mir-335-5p, hsa-let-7b-5p and hsa-mir-129-2-3p). Gene-disease association network analysis revealed that the DEGs were closely associated with intellectual disability and cerebellar ataxia. Drug-target enrichment analysis showed that the relevant drugs targeting the hub genes CACNA2D1, CACNA1A, CACNA1E, KCNA2 and KCNA5 were gabapentin, gabapentin enacarbil, pregabalin, guanidine hydrochloride and 4-aminopyridine. The results of this study provide a valuable basis for exploring the mechanisms of ion channel genes in COVID-19 and clues for developing therapeutic strategies for COVID-19.

Introduction

The new coronavirus pneumonia disease 2019 (COVID-19) caused by severe acute respiratory syndrome coronavirus type 2 (SARS-CoV-2) is a severe respiratory disease, with the first case identified in December 2019. It has now spread worldwide and poses a serious threat to public health and economic conditions (Wiersinga et al. 2020; Coronaviridae Study Group of the International Committee on Taxonomy of Viruses 2020; Dhama et al. 2020; Bchetnia et al. 2020; Weston and Frieman 2020; Ahsan et al. 2021). After the COVID-19 outbreak, effective diagnostic methods and infectious disease control measures, such as new coronavirus nucleic acid detection, urban blockade and use of masks, were implemented throughout the country, which decreased the spread of the disease in some countries and regions (Lian et al. 2020). Despite the rapid and effective containment of COVID-19 outbreaks in multiple regions, the global spread of COVID-19 has not been effectively controlled in all affected countries to date (Zhang et al. 2020a). According to the data provided by the World Health Organization, as of 22 December 2020, 78,299,811 confirmed cases of COVID-19 and more than 1.7 million deaths were reported worldwide (Ahamad et al. 2020; Aktar et al. 2021; Uddin et al. 2021). In addition, by the end of May 2021, approximately, 169 million COVID-19 cases and more than 3.5 million deaths had been confirmed worldwide (Aktar et al. 2021; Auwul et al. 2021). COVID-19 has emerged as one of the most devastating and long-lasting epidemics affecting human health in the 21st century. Therefore, effectively controlling and treating COVID-19 are grave concerns worldwide. Unlike other severe infectious diseases, COVID-19 lacks a typical clinical presentation. Respiratory symptoms and others reported to be related to COVID-19 include fever, cough, headache, conjunctivitis, diarrhoea and muscle or systemic pain (Rothan and Byrareddy 2020; Pascarella et al. 2020), which may be similar to the symptoms of other respiratory infections. Some patients may not have any symptoms after infection but may spread the disease. Moreover, specific drugs for treating COVID-19 have not yet been developed for clinical use, thus making the rapid diagnosis, effective control and treatment of COVID-19 difficult. Therefore, the identification of novel biomarkers at different omic levels (genomic, transcriptional or proteomic) may facilitate large-scale screening, diagnosis and treatment of COVID-19 (Chen et al. 2021), which is crucial to reveal the underlying pathogenesis of COVID-19, develop novel therapeutic strategies, discover potential therapeutic targets and improve therapeutic efficacy.

Ion channels are specific membrane proteins present in all cell membranes and some organelle membranes (e.g., mitochondria, Golgi apparatus, endoplasmic reticulum and lysosomes) (Wu et al. 2019; Tao et al. 2016). Ion channel genes encode these specific membrane proteins (mainly pore-forming membrane proteins) expressed in each living cell, which are oligomeric protein complexes composed of multiple subunits with ion-selective and voltage-gating properties (Noskov and Roux 2007). These membrane proteins can precisely control the passive influx and efflux of signalling ions into and out of the cell, thereby regulating ion concentrations inside and outside the cell membrane, membrane potential and volume size of the cell (Kondratskyi et al. 2018). Ion channels are involved in various physiological activities (Becchetti 2011), including muscle contraction, hormone secretion, cell proliferation and immune responses (Camerino et al. 2008; Fiske et al. 2006; Roger et al. 2006), and play an important role in maintaining homeostasis of the intracellular environment. Ion channels include sodium (Na +), potassium (K +), calcium (Ca +) and chloride (Cl-) ions and nonspecific cation channels (Lu et al. 2021a). Because ion channels play a key role in diverse biological functions, abnormal expression of ion channel genes plays a crucial role in many diseases (Sun et al. 2020). Ion channel genes are closely associated with tumour initiation and progression {e.g., breast cancer (Nelson et al. 2014), lung cancer (Ko et al. 2014), liver cancer (Lu et al. 2021b) and gastric cancer (Anderson et al. 2019)}, epilepsy (Oyrer et al. 2018), kidney stones, hypertension, insulin secretion deficiency and cardiac arrhythmias (Jentsch et al. 2004). However, limited information is available regarding the molecular characteristics of ion channel genes involved in COVID-19. In addition to impairing the respiratory system, COVID-19 has been reported to present in multiple organs to produce various clinical manifestations, including cardiovascular, urological, musculoskeletal, and neurological symptoms (Chen et al. 2020), whereas ion channels are known to be highly enriched in the nervous system and cardiac organs. It has been suggested that ion channels may be involved in the inflammation, pain, fever, anosmia, ageusia, respiratory, cardiovascular, gastrointestinal and neurological complications caused by COVID-19 infection (Jaffal and Abbas 2021). Epilepsy has been reported in the literature to occur with COVID-19 infection (Nikbakht et al. 2020), and epilepsy is currently considered to be an ion channel-related disorder. Therefore, exploring the relevant features of ion channels in COVID-19 and understanding their biological mechanisms are crucial for the treatment of COVID-19.

In this study, we adopted an analytical strategy involving an integrated bioinformatic approach to explore the mechanism of ion channel-related genes in COVID-19. We used biological datasets and several online databases to identify relevant features of ion channel genes in COVID-19 and identified 29 ion channel-related differentially expressed genes (DEGs) in COVID-19. Based on this, several bioinformatic analyses were performed to understand the involvement of these genes in the biological processes of the organism. In addition, we attempted to elucidate the pathogenic molecular mechanisms of these genes in COVID-19 and predict potential therapeutic agents. These differential genes have good application prospects for the diagnosis and treatment of COVID-19 and provide new perspectives for the discovery of potential biomarkers and drug targets of COVID-19.

Methods

Data Sources

To analyse the biological mechanisms and potential therapeutic targets of ion channel-related genes present in COVID-19, we obtained gene expression datasets (GSE152418 and GSE171110) from the Gene Expression Omnibus (GEO) database (http://www.ncbi.nlm.nih.gov/geo/). GSE152418 is based on the Illumina NovaSeq 6000 (Homo sapiens) (GPL24676) platform for RNA-sequencing (RNA-Seq) data on COVID-19. This dataset contains information on 17 patients with COVID-19 and 17 healthy subjects, with samples collected from peripheral blood mononuclear cells (PBMCs), and was derived from the research contributions of Arunachalam et al. (2020). GSE171110 is based on the Illumina HiSeq 2500 (Homo sapiens) (GPL16791) platform for RNA-Seq data on COVID-19. This dataset contains information on 44 patients with COVID-19 and 10 healthy subjects, with samples collected from whole blood tissue, and was derived from the research contributions of Lévy et al. (2021). Table 1 shows the basic information of both datasets. Ion channel genes were downloaded from the HUGO Gene Nomenclature Committee (HGNC) database (http://www.genenames.org/), and a total of 330 ion channel genes were obtained.

Table 1

Two RNA-seq transcriptome datasets and differential gene number information for COVID-19

ID	Country	Platform	Source	Healthy	COVID-19	Up regulated DEGs count	Down regulated DEGs count	Total DEGs count
GSE152418	USA	GPL24676	PBMC	17	17	1905	175	2080
GSE171110	France	GPL16791	Whole Blood	10	44	2620	1366	3986

Screening of Differential Genes

To identify ion channel-related DEGs, we downloaded the RStudio software (version 2022.02.0 + 443) (https://www.rstudio.com/), which is run on the R software (version 4.1.3) (https://www.r-project.org/). The RNA-Seq data of patients with COVID-19 (GSE152418 and GSE171110) were processed using RStudio based on the edgeR R package (version 3.36.0). The cut-off criteria of false discovery rate (FDR) < 0.05 and absolute value of log2-fold change (|logFC|) ≥ 1.0 were used to screen for significant DEGs in the abovementioned datasets. A Venn diagram was created using the VennDiagram R package (version 1.7.3) to show interacting genes by considering the intersection of DEGs obtained for each of the GSE152418 and GSE171110 datasets with the 330 ion channel-related genes. These interacting DEGs were used for subsequent analyses. Volcano plots were drawn using the EnhancedVolcano R package (version 1.12.0) to show the differential genes in the GSE152418 and GSE171110 datasets.

Functional and Pathway Enrichment Analyses

A comprehensive gene set enrichment web tool, Enrichr (https://maayanlab.cloud/Enrichr/), was used for functional annotation and pathway enrichment analysis of DEGs (Chen et al. 2013; Kuleshov et al. 2016; Xie et al. 2021). Gene set enrichment analysis (GSEA) is an important effort that enables the classification and summarisation of common biological insights to help understand the underlying biological mechanisms of target gene sets (Subramanian et al. 2005). We used Gene Ontology (GO) terms for functional enrichment analysis, which included the three aspects of biological process (BP), cellular composition (CC) and molecular function (MF). The Kyoto Encyclopedia of Genes and Genomes (KEGG), WikiPathways and Reactome databases were used as data sources for pathway enrichment analysis. An adjusted p value of < 0.05 was considered statistically significant GO terms and pathways.

Protein–Protein Interaction Networks and Screening of Hub Genes

In cell biology and systems biology, the evaluation and analysis of protein–protein interaction (PPI) networks and their functions are fundamental and key to the interpretation and understanding of cellular activities (Szklarczyk et al. 2019; Ewing et al. 2007; Ben-Hur and Noble 2005). We input the 29 identified DEGs into the STRING database (https://string-db.org/) to generate PPI networks (Szklarczyk et al. 2017). Furthermore, we downloaded the Cytoscape software (version 3.9.1) (https://cytoscape.org/) and imported the constructed PPI networks into this software for further processing and analysis (Shannon et al. 2003; Smoot et al. 2011). The Cytoscape software is an open platform that includes a number of plug-ins with scalable visualisation options and network analysis (Shannon et al. 2003). We used the cytoHubba plug-in in the Cytoscape software (http://apps.cytoscape.org/apps/cytohubba) to screen for hub genes. CytoHubba is a plug-in for ranking and extracting central, potential or targeted elements of a biological network based on various network features and contains 11 methods to score networks based on different perspectives, with the best one being Maximal Clique Centrality (MCC) at present (Chin et al. 2014). We used the MCC method to identify the top 10 hub genes in the PPI network.

Transcriptional and Post-transcriptional Regulatory Networks Analyses

Transcription factors (TFs) are proteins that attach to specific genes and control the rate of transcription of genetic information (Caramori et al. 2013). MicroRNAs (miRNAs) are a class of short, endogenously initiated and non-coding RNAs that repress or degrade messenger RNAs (mRNAs) through translation, thereby controlling gene expression at the post-transcriptional level (Cai et al. 2009). TFs and miRNAs are essential for molecular biology research. We used the online web tool NetworkAnalyst (Zhou et al. 2019; Xia et al. 2015) (https://www.networkanalyst.ca/) to construct TF–DEG and DEG–miRNA regulatory networks to analyse relevant TFs and miRNAs. The TF–DEG network was established using the JASPAR database. JASPAR is a publicly available repository that provides maps of TFs for multiple species in six major taxa (Khan et al. 2018). The DEG–miRNA network was established using the TarBase database. TarBase is the main experimental validity database for miRNAs interacting with target genes (Sethupathy et al. 2006).

Gene–Disease Association Analysis

We analysed DEGs using the DisGeNET database through the online web tool NetworkAnalyst to examine the association of these DEGs with diseases. DisGeNET is a comprehensive database exploring the association of genes and diseases based on various biomedical aspects of diseases, which synchronises relationships from multiple sources (Pinero et al. 2017). It provides and highlights new insights into the study of human genetic diseases (Pinero et al. 2020).

Protein–Drug Interaction Analysis

It is important to predict protein–drug interaction (PDI) based on target genes or identify potential drug molecules. We used the gene set enrichment network tool Enrichr based on the Drug Signature Database (DSigDB) to identify potential drugs that significantly interact with genes through PPI network pairs of the screened 10 hub genes (Yoo et al. 2015). DSigDB is a free web-based resource repository containing relevant information on drugs and their target genes for GSEA (Yoo et al. 2015). DSigDB currently contains a total of 22,527 gene sets, including 17,389 drugs and 19,531 genes (Auwul et al. 2021; Mahmud et al. 2021). An adjusted p value of < 0.05 was set as the statistical criteria for identifying drugs significantly associated with target genes.

Results

Identification of Differentially Expressed Genes

The COVID-19 datasets GSE152418 and GSE171110 from the GEO database were used for analysis. We used the RStudio software to process the data using the edgeR R package. The criteria for screening DEGs were as follows: FDR < 0.05 and |logFC|≥ 1.0. In the GSE152418 dataset, a total of 2080 genes were screened; of which, 1905 were upregulated and 175 were downregulated. In the GSE171110 dataset, a total of 3986 genes were screened; of which, 2620 were upregulated and 1366 were downregulated. Table 1 shows information regarding the number of DEGs in the two datasets. The VennDiagram R package was used for Venn analysis. The 330 ion channel-related genes obtained from the HGNC database were intersected with the DEGs of the GSE152418 and GSE171110 datasets, and 29 ion channel-related DEGs were eventually identified. Table 2 shows the expression information of 29 DEGs in the GSE152418 and GSE171110 datasets. The Venn diagram is shown in Fig. 1A. The visualization of differential genes in the GSE152418 and GSE171110 datasets is shown in Fig. 1B and C.

Table 2

Expression information of 29 DEGs in the GSE152418 and GSE171110 datasets

Gene	GSE152418			GSE171110
	Log2FC	FDR	Type	Log2FC	FDR	Type
KCNN3	3.875565247	5.06E − 22	Up	3.71805478	7.96E − 14	Up
CACNA2D3	− 1.934427624	7.81E − 11	Down	− 3.224177117	7.06E − 17	Down
SCN1B	1.694963987	1.59E − 10	Up	1.925623978	4.24E − 09	Up
KCNA2	1.188179351	4.16E − 06	Up	1.043307154	0.003261774	Up
AQP10	1.853616625	6.08E − 06	Up	1.891566838	0.001042295	Up
ANO2	2.739331713	1.02E − 05	Up	2.601512862	3.36E − 06	Up
CACNA1E	1.951378596	5.47E − 05	Up	1.782979464	4.19E − 05	Up
GJA4	2.811174491	9.77E − 05	Up	2.209293255	0.000160142	Up
CACNA1A	1.03471284	0.000146293	Up	1.332922487	2.75E − 05	Up
CATSPERD	2.776130738	0.000226835	Up	1.504241193	0.026799943	Up
TRPC1	− 1.033107512	0.00024422	Down	− 1.457710412	1.73E − 06	Down
KCNT2	2.646145885	0.000286884	Up	3.745253923	0.000357053	Up
TRPM1	2.259304259	0.000323608	Up	1.657205891	0.001286955	Up
ASIC4	2.064311007	0.000479613	Up	1.69682146	1.39E − 05	Up
NALCN	2.206464254	0.000638513	Up	1.675033515	0.000173107	Up
KCNG2	1.501487903	0.001142459	Up	3.019338629	1.24E − 05	Up
TRPM3	2.010824749	0.001300631	Up	1.281553025	0.00692089	Up
AQP9	1.060161386	0.001965165	Up	1.185764556	8.64E − 06	Up
GRIK5	1.698086032	0.002450449	Up	1.232928314	0.001037696	Up
KCNJ4	2.11639305	0.00385028	Up	3.155367763	0.001483182	Up
CHRNA2	1.659963162	0.005103148	Up	1.82443945	0.000381008	Up
KCNA5	− 1.170337718	0.005329720	Down	− 1.790680806	0.004616677	Down
GABRD	1.656924477	0.007725023	Up	1.481505214	0.030444915	Up
SCNN1B	2.131785521	0.009464007	Up	3.803559267	3.58E − 08	Up
TRPV5	1.442177547	0.010570825	Up	2.045644869	1.65E − 05	Up
CHRNB2	1.154450485	0.022489232	Up	1.760804737	5.78E − 05	Up
CNGB1	1.505139155	0.027513742	Up	1.175020522	0.000297923	Up
GJB5	1.361208469	0.034618091	Up	1.373879663	0.003175905	Up
CACNA2D1	1.767534655	0.041805155	Up	3.068993135	0.003340099	Up

Fig. 1

Identification of ion channel-related genes in the GSE152418 and GSE171110 datasets of differentially expressed genes (DEGs). A Screening of 29 DEGs. B Volcano plots of GSE152418. C Volcano plots of GSE171110

GO and Pathway Enrichment Analyses

GO and pathway enrichment analyses were performed to identify the biological significance and enriched pathways of the 29 DEGs of interest using the online web tool Enrichr. GO analysis included the following three categories: BP, CC and MF. GO enrichment analysis revealed that DEGs in BP were significantly enriched in ‘inorganic cation transmembrane transport’, ‘calcium ion transport’ and ‘calcium ion transmembrane transport’. The bar graph is shown in Fig. 2A. DEGs in CC were significantly enriched in ‘integral component of plasma membrane’, ‘voltage-gated potassium channel complex’, ‘potassium channel complex’ and ‘voltage-gated calcium channel complex’. The bar graph is shown in Fig. 2B. DEGs in MF were significantly enriched in the categories ‘ion channel activity’, ‘ligand-gated cation channel activity’, ‘calcium channel activity’, ‘voltage-gated cation channel activity’ and ‘cation channel activity’. The bar graph is shown in Fig. 2C. Table 3 enlists the top 10 terms enriched by DEGs in BP, CC and MF.

Fig. 2

GO enrichment analysis bar plot of DEGs obtained by Enrichr web tool. A Biological processes; B Cellular composition; C Molecular function

Table 3

GO functional enrichment analysis of DEGs

Category	ID	Term	Adjusted p-value	Number	Genes
Biological Process	GO:0098662	Inorganic cation transmembrane transport	7.47E-17	14	KCNJ4;KCNG2;TRPC1;KCNA2;CACNA1A;KCNA5;NALCN;TRPM1;ASIC4;SCNN1B;TRPV5;KCNN3;SCN1B;TRPM3
	GO:0006816	Calcium ion transport	4.47E − 12	9	TRPM1;CHRNB2;TRPC1;CACNA2D1;TRPV5;CACNA1A;CACNA1E;NALCN;TRPM3
	GO:0070588	Calcium ion transmembrane transport	2.18E − 09	7	TRPM1;TRPC1;CACNA2D1;TRPV5;CACNA1A;NALCN;TRPM3
	GO:0006812	Cation transport	5.05E − 09	8	TRPM1;KCNJ4;KCNG2;SCNN1B;KCNA2;KCNA5;TRPM3;CNGB1
	GO:0071805	Potassium ion transmembrane transport	3.64E − 08	7	KCNJ4;KCNG2;KCNT2;KCNA2;KCNA5;KCNN3;NALCN
	GO:0030001	Metal ion transport	8.77E − 08	6	CHRNB2;TRPC1;CACNA2D1;TRPV5;CACNA1A;CACNA1E
	GO:0006813	Potassium ion transport	5.42E − 07	6	KCNJ4;KCNG2;KCNA2;KCNA5;KCNN3;NALCN
	GO:0051899	Membrane depolarization	3.07E − 06	4	CHRNB2;CACNA1A;CACNA1E;SCN1B
	GO:0006811	Ion transport	1.40E − 05	5	CHRNB2;CHRNA2;KCNN3;TRPM3;CNGB1
	GO:0090276	Regulation of peptide hormone secretion	7.08E − 05	4	KCNG2;CACNA1A;KCNA5;CACNA1E
Cellular Component	GO:0005887	Integral component of plasma membrane	2.99E − 11	17	KCNJ4;CHRNB2;KCNG2;CHRNA2;AQP10;TRPC1;KCNA2;AQP9; KCNA5;TRPM1;ASIC4;SCNN1B;GJA4;TRPV5;GABRD;SCN1B; TRPM3
	GO:0008076	Voltage-gated potassium channel complex	5.75E − 05	4	KCNJ4;KCNG2;KCNA2;KCNA5
	GO:0034705	Potassium channel complex	5.75E − 05	4	KCNJ4;KCNG2;KCNA2;KCNA5
	GO:0005891	Voltage-gated calcium channel complex	9.89E − 05	3	CATSPERD;CACNA2D1;CACNA2D3
	GO:0034702	Ion channel complex	0.000114755	3	CHRNB2;CHRNA2;TRPC1
	GO:0034703	Cation channel complex	0.000861526	3	TRPC1;SCNN1B;KCNA5
	GO:0005892	Acetylcholine-gated channel complex	0.000861526	2	CHRNB2;CHRNA2
	GO:0005922	Connexin complex	0.001728713	2	GJA4;GJB5
	GO:0005921	Gap junction	0.001968584	2	GJA4;GJB5
	GO:0034706	Sodium channel complex	0.001968584	2	SCNN1B;SCN1B
Molecular Function	GO:0005216	Ion channel activity	1.86E − 13	9	TRPM1;CHRNB2;ASIC4;GRIK5;TRPC1;TRPV5;NALCN;TRPM3; CNGB1
	GO:0099094	Ligand-gated cation channel activity	6.55E − 10	7	KCNJ4;CHRNB2;ASIC4;CHRNA2;GRIK5;SCNN1B;KCNN3
	GO:0005262	Calcium channel activity	6.55E − 10	7	TRPM1;TRPC1;CACNA2D1;TRPV5;CACNA2D3;CACNA1A; CACNA1E
	GO:0022843	Voltage-gated cation channel activity	1.18E − 09	7	KCNG2;CACNA2D1;KCNA2;CACNA2D3;CACNA1A;KCNA5; CACNA1E
	GO:0005261	Cation channel activity	1.18E − 09	7	TRPM1;TRPC1;KCNA2;TRPV5;CACNA1E;NALCN;TRPM3
	GO:0005249	Voltage-gated potassium channel activity	2.20E − 08	6	KCNJ4;KCNG2;KCNT2;KCNA2;KCNA5;KCNN3
	GO:0005267	Potassium channel activity	2.06E − 06	5	KCNG2;KCNT2;GRIK5;KCNA2;KCNA5
	GO:0005245	Voltage-gated calcium channel activity	2.15E − 06	4	CACNA2D1;CACNA2D3;CACNA1A;CACNA1E
	GO:0005272	Sodium channel activity	2.15E − 06	4	ASIC4;GRIK5;SCNN1B;SCN1B
	GO:0015267	Channel activity	4.24E − 06	4	ASIC4;AQP10;AQP9;TRPV5

Pathway analysis is a technique to deduce the interaction between various diseases by modelling a molecular or biological process underlying the organism (Wittig and De Beuckelaer 2001), which can reveal the organism’s response to its intrinsic modifications. Three global databases, KEGG, WikiPathways and Reactome, were used as data sources for pathway enrichment analysis. The analysis showed that DEGs in the KEGG database were significantly enriched in ‘nicotine addiction’, ‘cholinergic synapse’, ‘glutamatergic synapse’, ‘MAPK signaling pathway’, ‘adrenergic signaling in cardiomyocytes’, ‘neuroactive ligand-receptor interaction’ and ‘oxytocin signaling pathway’. The bar graph is shown in Fig. 3A. DEGs in the WikiPathways database were significantly enriched in ‘calcium regulation in the cardiac cell’, ‘MAPK signaling pathway’ and ‘arrhythmogenic right ventricular cardiomyopathy’. The bar graph is shown in Fig. 3B. DEGs in the Reactome database were significantly enriched in ‘neuronal system’, ‘stimuli-sensing channels’ and ‘ion channel transport’. The bar graph is shown in Fig. 3C. Table 4 provides information regarding the top 10 pathways enriched by DEGs in the three databases, namely, KEGG, WikiPathways and Reactome.

Fig. 3

Pathway enrichment analysis bar plot of DEGs obtained by Enrichr web tool. A KEGG; B WikiPathways; C Reactome

Table 4

Pathway enrichment analysis of DEGs

Category	Term	Adjusted p-value	Number	Genes
KEGG 2021 Human	Nicotine addiction	0.000966407	3	CHRNB2;CACNA1A;GABRD
	Cholinergic synapse	0.007295400	3	KCNJ4;CHRNB2;CACNA1A
	Glutamatergic synapse	0.007295400	3	GRIK5;TRPC1;CACNA1A
	MAPK signalling pathway	0.007465173	4	CACNA2D1;CACNA2D3;CACNA1A;CACNA1E
	Adrenergic signalling in cardiomyocytes	0.007465173	3	CACNA2D1;CACNA2D3;SCN1B
	Neuroactive ligand-receptor interaction	0.007465173	4	CHRNB2;CHRNA2;GRIK5;GABRD
	Oxytocin signalling pathway	0.007465173	3	KCNJ4;CACNA2D1;CACNA2D3
	Type II diabetes mellitus	0.009340756	2	CACNA1A;CACNA1E
	GnRH secretion	0.015913653	2	TRPC1;KCNN3
	Arrhythmogenic right ventricular cardiomyopathy	0.018931604	2	CACNA2D1;CACNA2D3
WikiPathway 2021 Human	Calcium Regulation in the Cardiac Cell WP536	0.001559627	4	GJA4;GJB5;CACNA1A;CACNA1E
	MAPK Signalling Pathway WP382	0.005206294	4	CACNA2D1;CACNA2D3;CACNA1A;CACNA1E
	Arrhythmogenic Right Ventricular Cardiomyopathy WP2118	0.042832475	2	CACNA2D1;CACNA2D3
	Cell-type Dependent Selectivity of CCK2R Signalling WP3679	0.087319584	1	TRPC1
	Osteoclast Signalling WP12	0.087319584	1	TRPV5
	15q13.3 copy number variation syndrome WP4942	0.087319584	1	TRPM1
	Nicotine Activity on Dopaminergic Neurons WP1602	0.087319584	1	CHRNB2
	PKC-gamma calcium signalling pathway in ataxia WP4760	0.087319584	1	CACNA1A
	Type II diabetes mellitus WP1584	0.087319584	1	CACNA1A
	GABA receptor Signalling WP4159	0.095884738	1	GABRD
Reactome 2016	Neuronal System Homo sapiens R-HSA-112316	3.05E − 13	12	KCNJ4;CHRNB2;KCNG2;CHRNA2;GRIK5;CACNA2D1;KCNA2;CACNA2D3;CACNA1A;KCNA5;KCNN3;CACNA1E
	Stimuli-sensing channels Homo sapiens R-HSA-2672351	3.93E − 11	8	TRPM1;ASIC4;TRPC1;SCNN1B;TRPV5;NALCN;TRPM3;ANO2
	Ion channel transport Homo sapiens R-HSA-983712	7.95E − 09	8	TRPM1;ASIC4;TRPC1;SCNN1B;TRPV5;NALCN;TRPM3;ANO2
	Transmission across Chemical Synapses Homo sapiens R-HSA-112315	8.10E − 09	8	KCNJ4;CHRNB2;CHRNA2;GRIK5;CACNA2D1;CACNA2D3;CACNA1A;CACNA1E
	Depolarization of the Presynaptic Terminal Triggers the Opening of Calcium Channels Homo sapiens R-HSA-112308	3.43E − 08	4	CACNA2D1;CACNA2D3;CACNA1A;CACNA1E
	Transmembrane transport of small molecules Homo sapiens R-HSA-382551	6.76E − 08	10	TRPM1;ASIC4;AQP10;TRPC1;SCNN1B;AQP9;TRPV5;NALCN;TRPM3;ANO2
	TRP channels Homo sapiens R-HSA-3295583	4.29E − 07	4	TRPM1;TRPC1;TRPV5;TRPM3
	Potassium Channels Homo sapiens R-HSA-1296071	2.47E − 06	5	KCNJ4;KCNG2;KCNA2;KCNA5;KCNN3
	Cardiac conduction Homo sapiens R-HSA-5576891	1.02E − 05	5	KCNJ4;TRPC1;CACNA2D1;CACNA2D3;SCN1B
	Muscle contraction Homo sapiens R-HSA-397014	5.73E − 05	5	KCNJ4;TRPC1;CACNA2D1;CACNA2D3;SCN1B

PPI Networks and Hub Genes

PPI analysis was performed using the STRING database to identify key molecules. We input the 29 DEGs into the STRING database to generate a PPI network, setting the minimum required interaction score (MIS) to 0.150 and hiding disconnected nodes in the network. The analysis revealed that this PPI network contained 29 nodes and 161 edges (Fig. 4A). In the PPI network, the most interconnected nodes were considered hub genes. We imported the PPI network results into the Cytoscape software for network visualisation, and hub genes were screened using the MCC method in the cytoHubba plug-in. The top 10 hub genes screened were KCNA2, KCNJ4, CACNA1A, CACNA1E, NALCN, KCNA5, CACNA2D1, TRPC1, TRPM3 and KCNN3. The hub gene network is presented in Fig. 4B, and relevant information regarding these 10 hub genes is provided in Table 5. These hub genes may be potential biomarkers for COVID-19.

Fig. 4

PPI network analysis of DEGs. A PPI network obtained from STRING database. B 10 hub genes screened by cytoscape software

Table 5

Information of the 10 hub genes screened by PPI network analysis

Rank	Gene symbol	Gene name	Chromosome	Group name	Score
1	KCNA2	Potassium voltage-gated channel subfamily A member 2	1p13.3	Potassium voltage-gated channels	180,284
2	KCNJ4	Potassium inwardly rectifying channel subfamily J member 4	22q13.1	Potassium inwardly rectifying channel subfamily J	177,840
3	CACNA1A	Calcium voltage-gated channel subunit alpha1 A	19p13.13	Calcium voltage-gated channel alpha1 subunits	173,700
4	CACNA1E	Calcium voltage-gated channel subunit alpha1 E	1q25.3	Calcium voltage-gated channel alpha1 subunits	173,592
5	NALCN	Sodium leak channel, non-selective	13q32.3-q33.1	Sodium leak channels, non-selective	138,376
6	KCNA5	Potassium voltage-gated channel subfamily A member 5	12p13.32	Potassium voltage-gated channels	124,560
7	CACNA2D1	Calcium voltage-gated channel auxiliary subunit alpha2delta 1	7q21.11	Calcium voltage-gated channel auxiliary alpha2delta subunits	105,192
8	TRPC1	Transient receptor potential cation channel subfamily C member 1	3q23	Transient receptor potential cation channels	82,254
9	TRPM3	Transient receptor potential cation channel subfamily M member 3	9q21.12-q21.13	Transient receptor potential cation channels	81,414
10	KCNN3	Potassium calcium-activated channel subfamily N member 3	1q21.3	Potassium calcium-activated channels	63,386

Correlation Analysis of DEGs with TFs and miRNAs

To understand the transcriptional and post-transcriptional regulatory characteristics of DEGs, TF–DEG and DEG–miRNA interaction networks were constructed using the network database. The TF–DEG network contained 86 nodes and 193 edges, and TFs closely associated with DEGs were FOXC1, GATA2, HINFP, USF2, JUN and NFKB1. The TF–DEG network is presented in Fig. 5. The miRNA–DEG network contained 161 nodes and 246 edges, and miRNAs closely associated with DEGs were hsa-mir-146a-5p, hsa-mir-27a-3p, hsa-mir-335-5p, hsa-let-7b-5p and hsa-mir-129–2-3p. The DEG–miRNA network is shown in Fig. 6. These TFs and miRNAs may play an important regulatory role on DEGs.

Fig. 5

TF-DEG network obtained by NetworkAnalyst web tool. Red indicates DEGs, and blue indicates TF. The boxes marked in red indicate TF that are more closely associated with DEGs

Fig. 6

DEG-miRNA network obtained by NetworkAnalyst web tool. Red indicates DEGs, and blue indicates miRNA. The boxes marked in red indicate miRNA that are more closely associated with DEGs

Analysis of the Association between DEGs and Diseases

Different diseases may have one or more genes in common, making it possible for diseases to be interrelated (Al-Mustanjid et al. 2020). Discovering associations between genes and diseases can help in the development and design of disease treatment strategies (Moni and Lio 2014). The online web tool NetworkAnalyst was used to obtain a DEG–disease association network, which contained 248 nodes and 280 edges. We found that intellectual disability and cerebellar ataxia were more closely associated with the DEGs. This finding suggests that COVID-19 may be associated with intellectual disability and cerebellar ataxia. The DEG–disease association network is presented in Fig. 7.

Fig. 7

DEG–disease association network obtained by NetworkAnalyst web tool. Red indicates DEGs, and blue indicates disease. The boxes marked in red indicate disease that are more closely associated with DEGs

Drug Prediction Analysis

Assessment and analysis of protein–drug interactions are essential to understand the structural features required for receptor sensitivity (Mahmud et al. 2020; Mosharaf et al. 2020). In addition to investigating protein–drug interactions, we screened for candidates that could affect COVID-19. We considered the screened 10 hub genes as drug targets and performed drug-target enrichment analysis using the online web tool Enrichr based on DSigDB. The results showed that the drugs gabapentin, gabapentin enacarbil, pregabalin, guanidine hydrochloride and 4-aminopyridine, which act on five pivotal genes, namely, CACNA2D1, CACNA1A, CACNA1E, KCNA2 and KCNA5, respectively, may be potential drugs for the treatment of patients with COVID-19. Table 6 provides relevant information regarding these drugs.

Table 6

List of the predicted drugs identified from protein–drug interaction enrichment analysis

Drug Name	Adjusted p-value	Chemical formula	Structure	Genes
Gabapentin	5.13E − 06	C9H17NO2		CACNA2D1; CACNA1A; CACNA1E
Gabapentin enacarbil	5.13E − 06	C16H27NO6		CACNA2D1; CACNA1A; CACNA1E
Pregabalin	5.13E − 06	C8H17NO2		CACNA2D1; CACNA1A; CACNA1E
Guanidine hydrochloride	0.001705731	CH6ClN3		KCNA2; KCNA5
4-aminopyridine	0.001705731	C5H6N2		KCNA2; KCNA5

Discussion

COVID-19 is an emerging and rapidly growing pandemic with increasing infection and mortality rates worldwide (Team CC-R, 2020). COVID-19 has spread worldwide and poses a significant threat to humans. The current situation has prompted researchers to discover effective treatments against COVID-19 (Kumar et al. 2021). The causative agent of COVID-19, SARS-CoV-2, is highly pathogenic to humans. However, it is currently poorly understood, and specific treatments for COVID-19 remain unexplored. Hence, it is difficult to overcome this life-threatening prevalent disease (Li et al. 2021). Therefore, it is important to use bioinformatic methods to analyse the characteristics and pathogenesis of COVID-19 and discover novel therapeutic targets for the development of effective drugs and vaccines, thus providing a basis for public health decision-making (Ma et al. 2021). Ion channels are pore-forming membrane proteins that allow the passage of ions through the channel pore. Their functions include establishing the resting membrane potential (Abdul Kadir et al. 2018), shaping action potentials and other electrical signals by controlling ion flow across cell membranes, controlling ion flow in secretory and epithelial cells and regulating cell volume. Ion channels play an important role in various biological functions (Sun et al. 2020). However, the role of ion channel genes in COVID-19 remains unclear.

In this study, we used an integrated bioinformatic approach to gain insights into the associated features of ion channel-related genes in COVID-19. A total of 29 ion channel-related DEGs were identified in two RNA-Seq datasets (GSE152418 and GSE171110) containing data derived from the blood tissues of 61 patients with COVID-19 and 27 healthy subjects and including 330 ion channel gene sets. To examine the biological significance of these DEGs in the pathogenesis of COVID-19, we performed GO and pathway analyses on the DEGs. GO is a general theoretical model in gene regulation that outlines the functions of genes and their interrelationships (Al-Mustanjid et al. 2020). It develops progressively through the acquisition of biological knowledge regarding gene function and its regulation based on linguistic relationships among various ontological categories (Rana et al. 2019). The GO database was used as an annotation source for ontology to analyse the three categories, namely, BP, CC and MF, of the target genes. According to GO term interpretation, BP is the molecular activity, CC is the cellular structure in which genes regulate their function and MF is a description of activity at the molecular level (Moni and Lio 2015). Pathway analysis is a modern scientific strategy that helps to understand and reveal how biologically or molecularly complex diseases are connected and is the best way to obtain an organism’s response triggered by internal changes (Rana et al. 2020). In this study, GO enrichment analysis revealed that the DEGs were integral components of the plasma membrane (CC) and were significantly enriched in relevant functions such as inorganic cation transmembrane transport (BP) and ion channel activity (MF), and these processes mainly involved calcium and potassium channels. Furthermore, pathway enrichment analysis revealed that the DEGs were significantly enriched in pathways related to nicotine addiction (KEGG), calcium regulation in the cardiac cell (WikiPathways) and the neuronal system (Reactome). Calcium channels are activated upon membrane depolarisation to conduct calcium ions into the cell and organelles while initiating many physiological responses, including secretion, contraction and gene transcription (Zamponi et al. 2015). Calcium channel mutations and their dysfunctions have been associated with several diseases, such as disorders of the cardiovascular system {e.g., hypertension, arrhythmias and heart failure (Liao and Soong 2010; Venetucci et al. 2012)}, periodic skeletal muscle paralysis (Jurkat-Rott and Lehmann-Horn 2006), impaired insulin release and islets β-cell apoptosis in patients with diabetes (Yang et al. 2014), chronic pain and migraine (Bourinet et al. 2014; Kowalska et al. 2021) and numerous brain disorders (Heyes et al. 2015; Ortner and Striessnig 2016). However, the mechanism of action of calcium channels in COVID-19 remains unclear. Several studies (Neuraz et al. 2020; Peng et al. 2021; Kow et al. 2022) have suggested that the use of calcium channel blockers (CCBs), which reduce mortality in patients with COVID-19, has a therapeutic effect on COVID-19. However, other studies have reported (Mancia et al. 2020) that CCBs have no significant therapeutic effect on COVID-19 but increase the risk of tracheal intubation and death in patients with COVID-19 (Mendez et al. 2021). Potassium channels are located on the cell membrane and control the efflux and influx of potassium ions out of and into the cell (Kuang et al. 2015). They play a crucial role in both excitable and non-excitable cells. They are found in almost all species except some parasites (Kuo et al. 2005). The role of potassium channels in COVID-19 remains unknown; however, several studies have highlighted that multiple anti-COVID-19 drugs and inflammatory cytokines can interfere with cardiac potassium channels, such as the use of antibiotics (azithromycin and fluoroquinolones), antimalarials (hydroxychloroquine and chloroquine) and antivirals (lopinavir/ritonavir and atazanavir). In addition, some tyrosine kinase inhibitors (vandetanib) can inhibit hERG potassium channels and/or impair channel transport, thereby causing prolongation of the QT interval and increasing the risk of ventricular arrhythmias (Carpenter et al. 2020; Cubeddu et al. 2022). The smoke of inhaled cigarettes contains nicotine. Smoke particles carry nicotine to the pulmonary organs and are rapidly absorbed into the pulmonary venous circulation and subsequently into the arterial circulation, from where they move rapidly to the brain and bind to nicotinic cholinergic receptors (ligand-gated ion channels that normally bind to acetylcholine), producing and maintaining tobacco addiction (nicotine addiction) by acting on nicotinic cholinergic receptors in the brain and triggering the release of dopamine and other neurotransmitters, which is a major cause of disability and premature death in patients (Benowitz 2010). Although the association of smoking with the morbidity and mortality of a wide range of respiratory infections is well recognised, it remains unclear in COVID-19. Studies have suggested that active smokers do not have a high prevalence of COVID-19, which may be related to the ability of smoking to modulate angiotensin-converting enzyme-2 (ACE2) expression; however, the exact effects remain unclear (Usman et al. 2021). A recent study reported that smoking and nicotine may upregulate ACE2 (Brake et al. 2020). If smoking can upregulate ACE2, it may be a protective factor for COVID-19 (Verdecchia et al. 2020). However, studies published before the COVID-19 pandemic have reported that smoking and nicotine contribute to the downregulation of ACE2 (Oakes et al. 2018), which may promote increased expression of ACE2 receptors and viral receptors in smokers, thus increasing the opportunity for SARS-CoV-2 to invade the body (Berlin et al. 2020). However, the role of nicotine in COVID-19 requires further investigation. The regulatory role of calcium ions as intracellular second messengers (Bers 2008) in the heart is self-evident. It is well known that myocardial contraction is controlled by intracellular calcium ion concentration changes. The concentration of calcium ions in cardiomyocytes should be high enough to activate contractile proteins to pump blood out of the heart. During diastole, the concentration of calcium ions in cardiomyocytes should decrease to a sufficiently low level, which in turn relaxes the heart muscles so that the heart chamber becomes congested (Eisner 2018). This process relies on the regulation of calcium ion concentration, in which calcium channels play an important role. Studies have suggested a potential susceptibility of cardiomyocytes to COVID-19 (Yang et al. 2021). Cardiomyocytes contain abundant calcium ion channels; therefore, calcium regulation in cardiomyocytes may be one of the mechanisms of myocardial injury in patients with COVID-19. Ion channel-related genes, initially considered to be associated with inherited excitability disorders in the muscle and heart, play an important role in the molecular diagnosis of central nervous system diseases (Noebels 2017). Ion channels underlie the genesis of nerve impulses and are therefore an important component of the nervous system. Related studies have reported that SARS-CoV-2 can invade the nervous system (Liu et al. 2021; Mukerji and Solomon 2021) and heart (Van Cleemput et al. 2021) in humans. In this study, enrichment analysis suggested that ion channel-related genes play an important role. Overall, the results of GO and pathway analyses in this study partly explained the molecular basis and mechanism of action of the identified DEGs in the pathogenesis of COVID-19.

PPI networks are used to decode the key signalling molecules in molecular networks (Rahman et al. 2019). In this study, PPI network analysis revealed the most important hub proteins. We built a PPI network based on the 29 DEGs and screened 10 hub genes from them, which may be key drug targets or biomarkers for COVID-19. The KCNA2 gene encodes potassium voltage-gated channel subfamily A member 2, which is a member of the oscillator-like delayed rectifier potassium channel family (Corbett et al. 2016). It is mainly expressed in axons and presynaptic terminals in the central nervous system (Gu et al. 2003; Lorincz and Nusser 2008). It is now known that KCNA2 mutations can cause various neurological disorders, such as epileptic encephalopathy, mental retardation and motor disorders caused by cerebellar dysfunction (Doring et al. 2021). The KCNJ4 gene encodes potassium voltage-gated channel subfamily J member 4, which is an inward rectifier potassium channel family member. Studies have shown that KCNJ4 is associated with the progression and poor prognosis of lung adenocarcinoma (Wu and Yu 2019), dilated cardiomyopathy (Szuts et al. 2013) and prostate cancer (Kim et al. 2016). The CACNA1A gene encodes a subunit of the voltage-dependent P/Q-type calcium channel α-1A (Zhang et al. 2020b), and the CACNA1E gene encodes a subunit of the voltage-dependent R-type calcium channel α-1E (Helbig et al. 2018). These genes are widely expressed throughout the central nervous system and are strongly associated with epilepsy and intellectual developmental disorders (Hommersom et al. 2021; Royer-Bertrand et al. 2021). In addition, CACNA1E is of potential therapeutic value in non-small cell lung cancer (Gao et al. 2022). The NALCN gene encodes a non-selective cation channel that conducts a permanent sodium leak current and regulates the resting membrane potential and neuronal excitability associated with respiration, locomotion and circadian rhythms (Bramswig et al. 2018; Lutas et al. 2016; Shi et al. 2016). NALCN is essential for mammalian survival; however, the gating, ion selectivity and pharmacological properties of NALCN remain unclear (Chua et al. 2020; Kschonsak et al. 2020). The KCNA5 gene encodes potassium voltage-gated channel subfamily A member 5, which is involved in the regulation of several functions including cardiac action potential, vascular smooth muscle cell activity, insulin release and tumour cell proliferation (Bossini-Castillo et al. 2012; Ahmed et al. 2016). The CACNA2D1 gene encodes a calcium voltage-gated channel α2δ-1 subunit, which enhances channel transport, increases the expression of functional calcium channels at the plasma membrane and affects the biophysical properties of the channel (Dolphin 2012). It has been widely implicated in the regulation of neuronal excitability, action potential firing patterns and neurotransmission in nociceptive pathways (Gribkoff 2006). The TRPC1 gene encodes the transient C-potential subfamily channel 1 (Zeng et al. 2021), which is involved in the regulation of intracellular calcium ion concentration and plays an important role in cell proliferation, differentiation, apoptosis and migration and is expressed in almost all normal tissues and many tumours (Berridge et al. 2000; Zeng et al. 2016). Studies have shown that TRPC1 is a therapeutic target against herpes simplex virus type 1 (He et al. 2020). The TRPM3 gene encodes the transient receptor potential channel melastatin subfamily member 3, which is a non-selective calcium ion-permeable cation channel that can be activated by diverse stimuli including heat, osmotic pressure and chemically related activators (Held and Toth 2021). It activates/modulates calcium channels or transporters by guiding calcium ions through their pores or regulating membrane potential, which in turn increases the levels of intracellular calcium (Wu et al. 2010). TRPM3 is a recognised temperature receptor in the peripheral sensory neurons of the dorsal root ganglion, and mutations in TRPM3 in humans have recently been reported to be associated with epilepsy and intellectual disability (Zhao et al. 2020). KCNN3 encodes a neuronal small-conductance calcium-activated potassium channel containing two polyglutamine chains, which plays a key role in determining neuronal firing patterns and regulating intracellular calcium channels (Curtain et al. 2005). In this study, we screened 10 hub genes via the PPI network using COVID-19 datasets; however, the associated features are not reported. These genes may be potential therapeutic targets for COVID-19.

Because gene expression regulation is controlled by TFs and miRNAs at the transcriptional and post-transcriptional levels, changes in these molecules may provide critical information regarding the dysregulated expression of the 29 DEGs. Among the identified TFs, FOXC1, GATA2 and NFKB1 have been identified as important regulators of COVID-19 (Islam et al. 2020). HINFP, a histone cell cycle regulator, is a unique zinc-finger TF (Medina et al. 2008). Upstream stimulatory factor 2 (USF2), a TF involved in various cellular processes that is essential for maintaining reactive oxygen radical levels and mitochondrial morphology and function, is particularly prominent in tumour development (Chi et al. 2020). JUN is a subunit of activator protein 1 (an inducible transcription factor composed of multiple protein complexes), which plays a role in many types of cellular differentiation and inflammatory processes (Chang et al. 2013) and plays a specific role in the regulation of angiogenesis and endothelial cell proliferation (Yoshitomi et al. 2021). miRNAs are emerging as attractive biomarkers in numerous diseases. Among the identified miRNAs, hsa-mir-146a-5p is thought to regulate natural killer cells (innate lymphocytes with cytotoxic properties), thereby allowing the host to play an important role in early defence against infectious pathogens and surveillance against tumours (Pesce et al. 2018). The gene encoding hsa-mir-27a-3p is located on human chromosome 19, and its downregulation promotes tumour cell proliferation and the development of neurological diseases and is an important suppressor in diseases (Wu et al. 2015; Sala Frigerio et al. 2013). hsa-mir-335-5p plays an important role in regulating cardiac differentiation (Kay et al. 2019). hsa-let-7b-5p is a potential therapeutic target for tuberculosis (Tripathi et al. 2018). However, these miRNAs have not been reported to be associated with COVID-19 thus far. In conclusion, the TFs and miRNAs identified based on bioinformatic analyses in this study may serve as important regulators in the pathogenesis of COVID-19. In addition, we analysed the gene–disease relationship to predict the association of DEGs with diseases and found that intellectual disability and cerebellar ataxia were strongly associated with DEGs and may share common pathogenic features with COVID-19. It has been reported that the intellectually challenged population is at a greater risk of developing COVID-19 (Turk et al. 2020). Neurological symptoms have been suggested as potential complications of COVID-19, and cerebellar ataxia is a rare post-infectious or post-parainfection immune-mediated phenomenon associated with COVID-19 (Chan et al. 2021). The DEGs identified in this study may act as a ‘bridge’ in these diseases. Finally, we used DSigDB to identify drugs that may target the five screened pivotal genes, namely, CACNA2D1, CACNA1A, CACNA1E, KCNA2 and KCNA5. These drugs included gabapentin, gabapentin enacarbil, pregabalin, guanidine hydrochloride and 4-aminopyridine. Among these drugs, gabapentin can treat cough induced by acute and chronic COVID-19 infection, thereby reducing clinical symptoms (Song et al. 2021). Pregabalin may play a role in reducing the mortality of COVID-19 (Oddy et al. 2021). Guanidine alkaloids may have strong antiviral activity, thus providing a basis for the study of anti-COVID-19 (El-Demerdash et al. 2021). Further biological and clinical studies of these candidates are recommended to evaluate their potential therapeutic significance in patients with COVID-19.

In this study, we used bioinformatic analyses to investigate the features of ion channel-related genes associated with COVID-19 to identify key candidate genes and their regulatory molecules, examine the gene–disease association and discover potential therapeutic agents. To improve the reliability of the results, we used two datasets with data derived from blood tissues for analysis to avoid the influence of sample size and different tissue samples on the results. However, this study has some limitations owing to the lack of clinical validation of the identified molecules. Therefore, further validation is required to interpret the results.

Conclusions

In this study, we examined the molecular features of ion channel genes associated with COVID-19. On analysing the RNA-Seq transcriptome datasets (GSE152418 and GSE171110) and 330 ion channel-related genes downloaded from the HGNC database, we identified 29 DEGs. GO analysis revealed that these DEGs were integral components of the plasma membrane (CC) and were enriched in inorganic cation transmembrane transport (BP) and ion channel activity (MF). Pathway analysis revealed that the DEGs were enriched in pathways related to nicotine addiction (KEGG), calcium regulation in the cardiac cell (WikiPathways) and the neuronal system (Reactome). PPI networks were constructed using 29 DEGs, and 10 important hub genes (KCNA2, KCNJ4, CACNA1A, CACNA1E, NALCN, KCNA5, CACNA2D1, TRPC1, TRPM3 and KCNN3) were identified. Significant TFs (FOXC1, GATA2, HINFP, USF2, JUN and NFKB1) and miRNAs (hsa-mir-146a-5p, hsa-mir-27a-3p, hsa-mir-335-5p, hsa-let-7b-5p and hsa-mir-129–2-3p) were identified through the TF–DEG and DEG–miRNA networks. The DEG–disease association network revealed that intellectual disability and cerebellar ataxia were highly associated with these DEGs. Drug–target enrichment analysis based on DSigDB identified relevant drugs (gabapentin, gabapentin enacarbil, pregabalin, guanidine hydrochloride and 4-aminopyridine) targeting five hub genes (CACNA2D1, CACNA1A, CACNA1E, KCNA2 and KCNA5, respectively), which may have potential value for the treatment of COVID-19. Because the present study was based on bioinformatic analyses, further clinical studies should be performed to validate the identified molecular features. We hope that the results of this study will be helpful for the rapid control of COVID-19.