We performed the gene set enrichment analysis as an additional prioritization method. We first collected three differential gene expression data sets of hosts infected by HCoVs from the NCBI Gene Expression Omnibus (GEO). Among them, two transcriptome data sets were SARS-CoV-infected samples from patient’s peripheral blood94 (GSE1739) and Calu-3 cells95 (GSE33267), respectively. One transcriptome data set was MERS-CoV-infected Calu-3 cells96 (GSE122876). Adjusted P value less than 0.01 was defined as differentially expressed genes. These data sets were used as HCoV–host signatures to evaluate the treatment effects of drugs. Differential gene expression in cells treated with various drugs were retrieved from the Connectivity Map (CMAP) database36, and were used as gene profiles for the drugs. For each drug that was in both the CMAP data set and our drug–target network, we calculated an enrichment score (ES) for each HCoV signature data set based on previously described methods97 as follows:4 \documentclass[12pt]{minimal} 				\usepackage{amsmath} 				\usepackage{wasysym} 				\usepackage{amsfonts} 				\usepackage{amssymb} 				\usepackage{amsbsy} 				\usepackage{mathrsfs} 				\usepackage{upgreek} 				\setlength{\oddsidemargin}{-69pt} 				\begin{document}$${\mathrm {ES}} = \left\{ {\begin{array}{*{20}{c}} {{\mathrm{ES}_{\mathrm{up}}} - {\mathrm {ES}_{\mathrm{down}}},{\mathrm {sgn}}\left( {{\mathrm {ES}_{\mathrm{up}}}} \right)\, \ne \,{\mathrm {sgn}}\left( {{\mathrm {ES}_{\mathrm{down}}}} \right)} \\ {0,\mathrm {else}} \end{array}} \right.$$\end{document}ES=ESup−ESdown,sgnESup≠sgnESdown0,elseESup and ESdown were calculated separately for the up- and down-regulated genes from the HCoV signature data set using the same method. We first computed aup/down and bup/down as5 \documentclass[12pt]{minimal} 				\usepackage{amsmath} 				\usepackage{wasysym} 				\usepackage{amsfonts} 				\usepackage{amssymb} 				\usepackage{amsbsy} 				\usepackage{mathrsfs} 				\usepackage{upgreek} 				\setlength{\oddsidemargin}{-69pt} 				\begin{document}$$a = \mathop{\max}\limits_{1 \le j \le s}\left( {\frac{j}{s} - \frac{{V\left( j \right)}}{r}} \right),$$\end{document}a=max1≤j≤sjs−Vjr,6 \documentclass[12pt]{minimal} 				\usepackage{amsmath} 				\usepackage{wasysym} 				\usepackage{amsfonts} 				\usepackage{amssymb} 				\usepackage{amsbsy} 				\usepackage{mathrsfs} 				\usepackage{upgreek} 				\setlength{\oddsidemargin}{-69pt} 				\begin{document}$$b = \mathop{\max}\limits_{1 \le j \le s}\left( {\frac{{V\left( j \right)}}{r} - \frac{{j - 1}}{s}} \right),$$\end{document}b=max1≤j≤sVjr−j−1s,where j = 1, 2, …, s were the genes of HCoV signature data set sorted in ascending order by their rank in the gene profiles of the drug being evaluated. The rank of gene j is denoted by V(j), where 1 ≤ V(j) ≤ r, with r being the number of genes (12,849) from the drug profile. Then, ESup/down was set to aup/down if aup/down > bup/down, and was set to −bup/down if bup/down > aup/down. Permutation tests repeated 100 times using randomly generated gene lists with the same number of up- and down-regulated genes as the HCoV signature data set were performed to measure the significance of the ES scores. Drugs were considered to have potential treatment effect if ES > 0 and P < 0.05, and the number of such HCoV signature data sets were used as the final GSEA score that ranges from 0 to 3.