Active components of Ephedra and Asarum
Through database searches, 653 and 365 components of Ephedra and Asarum, respectively, were selected as highly active compounds with DL ≥ 0.18, and HL ≥ 4 as standards. After we deleted compounds without a target, 24 potential active compounds with good ADME properties, including quercetin, kaempferol, luteolin, naringenin, beta-sitosterol, and stigmasterol, were identified for further analysis (Table 1).
Putative targets of herbs used in AA and COVID-19-related target genes
A total of 647 potential human targets were identified for the above 24 compounds (Asarum: 171; Ephedra: 476). After deleting the duplicated targets, 252 targets were finally obtained. Additionally, after deleting the duplicated genes among 736 disease-related genes obtained from the three disease databases, 713 COVID-19-related target genes were identified .
PPI network of putative AA and COVID-19-related target genes
The PPI network of putative AA-related target genes (Figure 2A) contained 236 nodes and 1889 edges and that of COVID-19-related target genes (Figure 2B) contained 352 nodes and 2588 edges. We then used Cytoscape software and GraphPad Prism software to construct the PPI barplot, which showed the target genes significant to both AA and COVID-19 (Figure 2C and 2D).
To obtain AA targets that interact with the disease, we used the Bioinformatics & Evolutionary Genomics platform to match potential AA target genes with disease-related genes and drew a Venn diagram (Figure 3). Finally, 56 putative AA targets that interacted with disease-related targets were identified.
PPI network analysis of AA action against COVID-19
We imported the 56 key targets into the STRING database to generate a PPI network (Figure 4A) (settings: species Homo sapiens, confidence level > 0.7, STRING database confidence level high > 0.7; medium > 0.4; and low > 0.15).
The PPI results were imported into Cytoscape software, and clustermaker2 and Network analyzer plugins were used to generate a subnetwork (Figure 4B) containing 53 nodes and 388 edges. The network radius was two, and the average number of neighbors was 14.642. The size of the circle and the font indicate the degree value, and different colors indicate the number of adjacent points connected. The top10 degrees were IL6 (degree = 34), TNF (degree = 34), MAPK1 (degree = 31), TP53 (degree = 30), MAPK8 (degree = 30), IL1B (degree = 28), RELA (degree = 27), CASP3 (degree = 25), STAT1 (degree = 23), MAPK3 (degree = 23), and IL10 (degree = 22).
Subsequently, the MCODE plug-in in Cytoscape was used for further cluster analysis, and two module clusters in the network were generated (Table 2 and Figure 5). Module 1 had the highest node, edge number, and score values. This may be the main model by which the targets exert their effects; this module included IL6, TNF, PTGS2, IL1B, MAPK1, MAPK3, MAPK8, MAPK14, among others. This cluster consisted of 14 component nodes and 74 edges, the clustering confidence was 0.887, and the average number of neighbors was 11. STAT1 (marked in yellow in Figure 5) was the seed of the cluster; this target might thus play an important role in connecting other nodes in the PPI network.
Component-target network analysis
We used Cytoscape 3.2.1 software to construct a network and analyze the 24 compounds and 252 target genes identified (Figure 6). The network contained 275 nodes (nodes) and 582 edges (edges); the red square represents the compound, the blue circle represents the target protein, and the edges represent the relationship. The network density was 0.015, the distance was 5, and the median of 24 candidate components was 10 degrees, indicating multiple target genes. The top-three effective components were quercetin, kaempferol, and luteolin, which had 154, 65, and 58 target genes, respectively. This network verified that Ephedra and Asarum had a therapeutic effect through multi-component, multi-target, integrated regulation.
GO function enrichment analysis
We imported the 56 target genes into the DAVID database and used R 3.5.3 software to perform GO functional enrichment analysis. The 10 most significant BPs in Figure 7A included an inflammatory response, immune response, positive regulation of transcription from RNA polymerase II promoter, response to drug, positive regulation of transcription, DNA-template, positive and negative regulation of apoptotic process, apoptotic process, and cellular response to lipopolysaccharide. The top-10 CCs in Figure 7B included the plasma membrane, nucleus, cytoplasm, cytosol, nucleoplasm, extracellular space, extracellular exosome, integral component of plasma membrane, extracellular region, and mitochondrion. The top-10 MFs in Figure 7C included protein binding, identical protein binding, enzyme binding, protein homodimerization activity, ATP binding, protein heterodimerization activity, DNA binding, transcription factor activity, sequence-specific DNA binding, protein kinase binding, and zinc ion binding.
KEGG pathway enrichment analysis
As shown in Figure 8A, the herbs affected COVID-19 through multiple signaling pathways; we selected 20 related pathways with corrected P values ≤ 0.01 and gene frequency > 0.5%. In this network, the color represents the p value, and the length represents the number of genes. The top pathways, indicating the key pharmacological mechanisms of the AA herbs against COVID-19, included the Immune System, Cytokine Signaling in the Immune System, Adaptive Immune System, Signaling by Interleukins, Innate Immune System, Signaling by GPCR, Hemostasis, HTLV-I infection, Signal Transduction,Pathways in cancer, PI3K-Akt signaling pathway.
To understand the relationship between protein targets and signal pathways, we selected 10 pathways with higher degree values to establish a target‒pathway network diagram in Figure 8B. The network had 59 nodes and 217 edges. The green squares represent the pathway and the yellow circles represent the targets.
In the subsequent pathway verification, we found that the IL-17, Influenza A, and TNF signaling pathway were consistent with our previous pathway enrichment analysis. The red marks in Figure 8C‒8E represent potential target genes related to AA.
Molecular docking
Quercetin and kaempferol, selected from the above methods were considered to be the main candidate compounds, and IL-6, TNF, and STAT1 were considered to be the core targets. Molecular docking revealed that two TCM compounds bound to TNF with the lowest binding energy (i.e., -7.8 kcal/mol, and the others were -7.6 [quercetin and IL-6], -7.3 [quercetin and STAT1], -7.2 [kaempferol and STAT1], and -6.9 [kaempferol and IL-6]) (Table 3 and Figure 9).