# Deciphering Pharmacological Mechanism of Buyang Huanwu Decoction for Spinal Cord Injury by Network Pharmacology Approach

Authors: Zhencheng Xiong; Feng Yang; Wenhao Li; Xiangsheng Tang; Haoni Ma; Ping Yi
Journal: Evidence-Based Complementary and Alternative Medicine
 (2021)

---

## Abstract

Objective. The purpose of this study was to investigate the mechanism of action of the Chinese herbal formula Buyang Huanwu Decoction (BYHWD), which is commonly used to treat nerve injuries, in the treatment of spinal cord injury (SCI) using a network pharmacology method. Methods. BYHWD-related targets were obtained by mining the TCMSP and BATMAN-TCM databases, and SCI-related targets were obtained by mining the DisGeNET, TTD, CTD, GeneCards, and MalaCards databases. The overlapping targets of the abovementioned targets may be potential therapeutic targets for BYHWD anti-SCI. Subsequently, we performed protein-protein interaction (PPI) analysis, screened the hub genes using Cytoscape software, performed Gene Ontology (GO) annotation and KEGG pathway enrichment analysis, and finally achieved molecular docking between the hub proteins and key active compounds. Results. The 189 potential therapeutic targets for BYHWD anti-SCI were overlapping targets of 744 BYHWD-related targets and 923 SCI-related targets. The top 10 genes obtained subsequently included AKT1, IL6, MAPK1, TNF, TP53, VEGFA, CASP3, ALB, MAPK8, and JUN. Fifteen signaling pathways were also screened out after enrichment analysis and literature search. The results of molecular docking of key active compounds and hub target proteins showed a good binding affinity for both. Conclusion. This study shows that BYHWD anti-SCI is characterized by a multicomponent, multitarget, and multipathway synergy and provides new insights to explore the specific mechanisms of BYHWD against SCI.


---

## Body

## 1. Introduction

Spinal cord injury (SCI) is a group of disorders in which the integrity and continuity of the spinal cord are disrupted by trauma, tumor, hypoxia, inflammation, and other factors [1]. According to the World Health Organization (WHO) report, about 250,000 to 500,000 people worldwide suffer from SCI every year, and the population is mainly under 30 years old [2]. The Beijing SCI Epidemiology Survey Group reported in 2002 that the incidence of SCI in Beijing was about 60/106, an increase of nearly 10 times compared to the late 1980s, and that the per capita hospitalization cost was 27819.3 CNY [3]. The high disability rate and high consumption of SCI have imposed a heavy burden on the families of patients and society, while how to protect and restore injured neurons after SCI has been considered a great challenge in clinical and experimental research [4]. There is still a lack of effective treatment for SCI, but it is known that traditional Chinese medicine (TCM) has a long history and rich practical experience, which may be explored for the treatment of SCI [2].Modern Chinese medicine scholars classify SCI as “Ti Duo” and “Wei Zheng” based on the TCM classics “Ling Shu Jing-Han Re Bing” and “Su Wen-Wei Lun” [3]. Based on the analysis of the etiology and pathogenesis of SCI by TCM, as well as the pathological changes of the damaged Governor Vessel and the symptoms corresponding to SCI, the main treatment of TCM is to unblock Governor Vessel and strengthen “Yang Qi,” supplemented by tonifying the kidney and filling the essence, strengthening the tendons and bones, activating blood circulation and removing blood stasis, and warming the meridians [3]. Buyang Huanwu Decoction (BYHWD), recorded in “Yi Lin Gai Cuo” by Wang Qingren in the Qing Dynasty, has been used in China for hundreds of years to treat stroke-related disabilities and improve the recovery of neurological function [5]. The whole formula of BYHWD consists of seven herbs: Chishao (Radix Paeoniae Rubra), Chuanxiong (Rhizoma Chuanxiong), Danggui (Radix Angelicae Sinensis), Dilong (Pheretima), Huangqi (Radix Astragali), Honghua (Flos Carthami), and Taoren (Semen Persicae) [5]. Studies have shown that BYHWD may promote the repair of SCI by upregulating the expression of Notch1 gene [6], and BYHWD combined with bone marrow mesenchymal stem cells (BMSCs) transplantation can promote the recovery after SCI by saving axotomized red nucleus neurons in rats [4]. Although some progress has been made in the study of BYHWD for SCI, the complexity of herbal prescriptions with multiple components, targets, and pathways makes the study of the specific mechanism of action of BYHWD for SCI difficult to clarify.In the last decade, with the development of modern multiomics theories such as genomics, proteomics, and metabolomics, the introduction of a systems biology perspective, and the application of bioinformatics, the concept of network pharmacology has emerged. According to the network basis of drug action, network pharmacology methods, including protein-protein interaction (PPI) network, Gene Ontology (GO) annotation, and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis, are applied to infer the specific pathways of multimolecular drug action in the human body, which is consistent with the multicomponent, multitarget, and multipathway synergy of Chinese herbal medicine and helps to clarify the mechanism of action of BYHWD (Figure1).Figure 1
The overall flowchart of this study.


## 2. Materials and Methods

### 2.1. Acquisition of BYHWD-Related Targets and SCI-Related Targets

The compounds related to BYHWD were obtained by searching “Chishao or CHI SHAO orRadix Paeoniae Rubra,” “Chuanxiong or CHUAN XIONG or Rhizoma Chuanxiong,” “Danggui or DANG GUI or Radix Angelicae Sinensis,” “Dilong or DI LONG or Pheretima,” “Huangqi or HUANG QI or Radix Astragali,” “Honghua or HONG HUA or Flos Carthami,” and “Taoren or TAO REN or Semen Persicae” in Traditional Chinese Medicine Systems Pharmacology (TCMSP, version 2.3, https://tcmspw.com/tcmsp.php) database [7], and then the active compounds that met the criteria (oral bioavailability (OB) ≥ 30 and drug-likeness (DL) ≥ 0.18) were further screened by the absorption, distribution, metabolism, and excretion (ADME) principle [8, 9], and the corresponding target protein information of the active compounds was also obtained. Subsequently, by using the UniProt database (https://www.uniprot.org/) (species should be set to “Homo sapiens”), the target protein name is converted into the corresponding gene name and UniProt ID [10]. The compounds and targets related to BYHWD were also obtained directly by searching the above keywords in the BATMAN-TCM platform (http://bionet.ncpsb.org/batman-tcm/) according to the set criteria (score cutoff >25, adjusted P value < 0.05) [11]. Finally, the results of the two databases are summarized, integrated, and deduplicated to obtain BYHWD-related targets.DisGeNet database (https://www.disgenet.org/, v7.0) (score ≥ 0.01) [12], GeneCards (https://www.genecards.org/) (relevance score ≥ 20) [13], MalaCards (https://www.malacards.org/) [14], Comparative Toxicogenomics Database (CTD) (http://ctdbase.org/, last update by June 2020) (inference score ≥ 20) [15], and Therapeutic Target Database (TTD) (http://db.idrblab.net/ttd/, last update by June 1, 2020) [16] were used to screen out potential targets related to SCI by searching for the keyword “Spinal cord injury or Spinal cord diseases.” The results of the five databases are summarized, integrated, and deduplicated to obtain SCI-related targets.


### 2.2. Acquisition of Potential Therapeutic Targets for BYHWD Anti-SCI

Potential therapeutic targets of BYHWD in the treatment of SCI should be derived from the overlapping targets of BYHWD-related targets and SCI-related targets obtained from the Venn online tool (http://www.bioinformatics.com.cn/). Subsequently, we used Cytoscape software (version 3.7.2) to construct a network diagram of “BYHWD-active compounds-target genes-SCI” [17].


### 2.3. Construction of PPI Network and Acquisition of Hub Genes

PPI is composed of proteins that interact with each other to participate in various aspects of life processes such as biological signaling, regulation of gene expression, energy and material metabolism, and cell cycle regulation [18]. The potential therapeutic targets of BYHWD against SCI were analyzed by a PPI network using the STRING database (http://string-db.org/; version 11) with the species set to “Homo sapiens” and a confidence level greater than 0.4 [19]. Subsequently, the PPI network was visualized by Cytoscape software (version 3.7.2).To obtain the hub genes for BYHWD anti-SCI, we use the CytoHubba plug-in in Cytoscape software, which provides 12 topological analysis methods, including Degree, Maximal Clique Centrality (MCC), BottleNeck, Betweenness, ClusteringCoefficient, Closeness, Density of Maximum Neighborhood Component (DMNC), EcCentricity, Edge Percolated Component (EPC), Maximum Neighborhood Component (MNC), Radiality, and Stress [20, 21].


### 2.4. GO and KEGG Pathway Enrichment Analysis of Potential Therapeutic Genes for BYHWD Anti-SCI

R is commonly used for statistical analysis and plotting software, and its open-source and free features have been favored by many programming enthusiasts, resulting in the birth of many convenient and efficient packages, such as the ClusterProfiler package for GO annotation and KEGG pathway enrichment analysis [22–24]. In this study, we used ClusterProfiler package (adjusted P value < 0.05) to explore the biological processes (BP), cellular components (CC), molecular functions (MF), and signaling pathways associated with BYHWD anti-SCI.


### 2.5. Molecular Docking of Hub Genes and Key Active Compounds

Based on the previously mined BYHWD-related targets and related active compounds, the active compounds and herbs corresponding to the top 10 hub genes need to be sorted out. The interrelationship between “herbs in BYHWD-active compounds-top 10 hub genes” was visualized by a Sankey diagram (http://sankeymatic.com/). In the following, we will select key active compounds as small molecule ligands and corresponding large molecule receptors (proteins encoded by the top 10 hub genes) for molecular docking analysis. The raw files of key active compounds (MOL2 format) can be downloaded from PubChem database (https://pubchem.ncbi.nlm.nih.gov/), and the raw files of the top 10 core genes (PDB format) can be downloaded from the RCSB protein data (http://www.rcsb.org/), and after corresponding processing by AutoDock Tool [25], they are finally converted to PDBQT format for molecular docking in Pymol software (https://pymol.org/2/; version 2.4.1) [26]. The molecular docking score reflects the binding affinity, with smaller values representing higher affinity for the binding of small molecule ligands and large molecule receptors [27].


## 2.1. Acquisition of BYHWD-Related Targets and SCI-Related Targets

The compounds related to BYHWD were obtained by searching “Chishao or CHI SHAO orRadix Paeoniae Rubra,” “Chuanxiong or CHUAN XIONG or Rhizoma Chuanxiong,” “Danggui or DANG GUI or Radix Angelicae Sinensis,” “Dilong or DI LONG or Pheretima,” “Huangqi or HUANG QI or Radix Astragali,” “Honghua or HONG HUA or Flos Carthami,” and “Taoren or TAO REN or Semen Persicae” in Traditional Chinese Medicine Systems Pharmacology (TCMSP, version 2.3, https://tcmspw.com/tcmsp.php) database [7], and then the active compounds that met the criteria (oral bioavailability (OB) ≥ 30 and drug-likeness (DL) ≥ 0.18) were further screened by the absorption, distribution, metabolism, and excretion (ADME) principle [8, 9], and the corresponding target protein information of the active compounds was also obtained. Subsequently, by using the UniProt database (https://www.uniprot.org/) (species should be set to “Homo sapiens”), the target protein name is converted into the corresponding gene name and UniProt ID [10]. The compounds and targets related to BYHWD were also obtained directly by searching the above keywords in the BATMAN-TCM platform (http://bionet.ncpsb.org/batman-tcm/) according to the set criteria (score cutoff >25, adjusted P value < 0.05) [11]. Finally, the results of the two databases are summarized, integrated, and deduplicated to obtain BYHWD-related targets.DisGeNet database (https://www.disgenet.org/, v7.0) (score ≥ 0.01) [12], GeneCards (https://www.genecards.org/) (relevance score ≥ 20) [13], MalaCards (https://www.malacards.org/) [14], Comparative Toxicogenomics Database (CTD) (http://ctdbase.org/, last update by June 2020) (inference score ≥ 20) [15], and Therapeutic Target Database (TTD) (http://db.idrblab.net/ttd/, last update by June 1, 2020) [16] were used to screen out potential targets related to SCI by searching for the keyword “Spinal cord injury or Spinal cord diseases.” The results of the five databases are summarized, integrated, and deduplicated to obtain SCI-related targets.


## 2.2. Acquisition of Potential Therapeutic Targets for BYHWD Anti-SCI

Potential therapeutic targets of BYHWD in the treatment of SCI should be derived from the overlapping targets of BYHWD-related targets and SCI-related targets obtained from the Venn online tool (http://www.bioinformatics.com.cn/). Subsequently, we used Cytoscape software (version 3.7.2) to construct a network diagram of “BYHWD-active compounds-target genes-SCI” [17].


## 2.3. Construction of PPI Network and Acquisition of Hub Genes

PPI is composed of proteins that interact with each other to participate in various aspects of life processes such as biological signaling, regulation of gene expression, energy and material metabolism, and cell cycle regulation [18]. The potential therapeutic targets of BYHWD against SCI were analyzed by a PPI network using the STRING database (http://string-db.org/; version 11) with the species set to “Homo sapiens” and a confidence level greater than 0.4 [19]. Subsequently, the PPI network was visualized by Cytoscape software (version 3.7.2).To obtain the hub genes for BYHWD anti-SCI, we use the CytoHubba plug-in in Cytoscape software, which provides 12 topological analysis methods, including Degree, Maximal Clique Centrality (MCC), BottleNeck, Betweenness, ClusteringCoefficient, Closeness, Density of Maximum Neighborhood Component (DMNC), EcCentricity, Edge Percolated Component (EPC), Maximum Neighborhood Component (MNC), Radiality, and Stress [20, 21].


## 2.4. GO and KEGG Pathway Enrichment Analysis of Potential Therapeutic Genes for BYHWD Anti-SCI

R is commonly used for statistical analysis and plotting software, and its open-source and free features have been favored by many programming enthusiasts, resulting in the birth of many convenient and efficient packages, such as the ClusterProfiler package for GO annotation and KEGG pathway enrichment analysis [22–24]. In this study, we used ClusterProfiler package (adjusted P value < 0.05) to explore the biological processes (BP), cellular components (CC), molecular functions (MF), and signaling pathways associated with BYHWD anti-SCI.


## 2.5. Molecular Docking of Hub Genes and Key Active Compounds

Based on the previously mined BYHWD-related targets and related active compounds, the active compounds and herbs corresponding to the top 10 hub genes need to be sorted out. The interrelationship between “herbs in BYHWD-active compounds-top 10 hub genes” was visualized by a Sankey diagram (http://sankeymatic.com/). In the following, we will select key active compounds as small molecule ligands and corresponding large molecule receptors (proteins encoded by the top 10 hub genes) for molecular docking analysis. The raw files of key active compounds (MOL2 format) can be downloaded from PubChem database (https://pubchem.ncbi.nlm.nih.gov/), and the raw files of the top 10 core genes (PDB format) can be downloaded from the RCSB protein data (http://www.rcsb.org/), and after corresponding processing by AutoDock Tool [25], they are finally converted to PDBQT format for molecular docking in Pymol software (https://pymol.org/2/; version 2.4.1) [26]. The molecular docking score reflects the binding affinity, with smaller values representing higher affinity for the binding of small molecule ligands and large molecule receptors [27].


## 3. Results

### 3.1. Acquisition Results of BYHWD-Related Targets and SCI-Related Targets

A total of 775 compounds (BATMAN-TCM: 219; after deduplication: 189) in BYHWD were preliminarily retrieved from the TCMSP database and BATMAN-TCM platform (score cutoff >25; adjustedP value < 0.05), of which 119 (BATMAN-TCM: 8) were from Chishao, 189 (87) were from Chuanxiong, 125 (77) were from Danggui, 87 (23) were from Huangqi, 189 (22) were from Honghua, and 66 (2) were from Taoren. Relevant compounds and targets’ information for Dilong was not retrieved in either database. After setting the ADME criteria (OB ≥ 30% and DL ≥ 0.18), a total of 103 active compounds (after removing duplicates and nontarget compounds: 60) were screened from the TCMSP database, of which 29 (After removing nontarget compounds: 14) were from Chishao, 7 (6) were from Chuanxiong, 2 (2) were from Danggui, 20 (17) were from Huangqi, 22 (17) were from Honghua, and 23 (19) were from Taoren. The BATMAN-TCM platform provides the gene names directly, while the TCMSP database first provides the target protein names of the active compounds and therefore needs to be converted to gene names via the UniProt database. Finally, after integrating the results of the two databases, we obtained 244 BYHWD-related active compounds and 744 BYHWD-related targets.After deduplication, a total of 923 SCI-related targets were obtained by mining the DisGeNET (number: 84), TTD (number: 5), CTD (number: 553), GeneCards (number: 432), and MalaCards (number: 19) databases.


### 3.2. Acquisition Results of Potential Therapeutic Targets for BYHWD Anti-SCI

The overlapping targets of BYHWD-related targets and SCI-related targets were considered as potential therapeutic targets for BYHWD anti-SCI. As shown in Figure2, we screened the overlapping targets of both by constructing a Venn diagram, and a total of 189 targets were obtained. Through multiple databases mining and construction of the Venn diagram, we obtained BYHWD-related compounds and potential therapeutic targets for BYHWD anti-SCI and finally used Cytoscape software to construct a “BYHWD-active compounds-target genes-SCI” network (Figure 3). In Figure 3, the lines between two nodes represent the existence of mutual relationships, and the larger the node is, the more relationships exist.Figure 2
Venn diagram for BYHWD-related targets and SCI-related targets.Figure 3
BYHWD-active compounds-target genes-SCI” network. The purple diamond represents the herbs contained in BYHWD; the red octagon represents SCI; the cyan rectangle represents the potential therapeutic targets for BYHWD anti-SCI; the light green ellipse represents the active compounds contained in BYHWD; the blue ellipse represents the common active compounds of the herb.


### 3.3. Acquisition Results of Hub Genes for BYHWD Anti-SCI and PPI Network Construction

In the STRING database, we first set the species to “Homo sapiens” and then entered 189 potential therapeutic targets for BYHWD anti-SCI to obtain the PPI network (Figure 4(a)), which involved 189 nodes and 3717 edges. The obtained TSV files were then imported into Cytoscape software (version 3.7.2) for further analysis and visualization (Figure 4(b)). To obtain the hub genes for BYHWD anti-SCI, based on the above the PPI network, we used the CytoHubba plug-in of Cytoscape software, which currently contains 12 topological analysis methods. We targeted the top 10 hub genes for BYHWD anti-SCI, and the results of the 12 algorithms each contained 10 genes. We integrated and analyzed the results to obtain a total of 45 different genes, which were then sorted by the number of algorithms to which they were attributed. Table 1 shows the basic information on hub genes for BYHWD anti-SCI. The top 10 hub genes (AKT1, IL6, MAPK1, TNF, TP53, VEGFA, CASP3, ALB, MAPK8, and JUN) sorted are completely consistent with some algorithms (Degree, Closeness, and MNC) (Figure 4(c)).Figure 4
PPI network and top 10 hub genes for BYHWD anti-SCI. (a) PPI network constructed with STRING. (b) PPI network constructed with Cytoscape software (the darker the node color, the higher the number of connected proteins). (c) Top 10 hub genes for BYHWD anti-SCI were obtained by using the Degree algorithm.
(a)(b)(c)Table 1
Basic information of hub genes based on 12 topological algorithms.
UniProt IDGene symbolProtein namesAlgorithmsP31749AKT1RAC-alpha serine/threonine-protein kinase(1), (2), (3), (4), (5), (8), (10), (11), and (12)P05231IL6Interleukin-6(1), (2), (3), (4), (5), (8), (10), (11), and (12)P28482MAPK1Mitogen-activated protein kinase 1(1), (2), (3), (5), (8), (9), (10), (11), and (12)P01375TNFTumor necrosis factor(1), (2), (3), (4), (5), (8), (10), (11), and (12)P04637TP53Cellular tumor antigen p53(1), (2), (3), (4), (5), (8), (10), (11), and (12)P15692VEGFAVascular endothelial growth factor A(1), (2), (3), (4), (5), (8), (10), (11), and (12)P42574CASP3Caspase-3(1), (2), (3), (5), (8), (10), (11), and (12)P02768ALBAlbumin(1), (3), (5), (8), (10), (11), and (12)P45983MAPK8Mitogen-activated protein kinase 8(1), (2), (4), (5), (8), and (10)P05412JUNTranscription factor AP-1(1), (2), (5), (10), and (11)P04040CATCatalase(3), (4), and (12)P99999CYCSCytochrome c(3) and (12)P09238MMP10Stromelysin-2(6) and (7)P07101THTyrosine 3-monooxygenase(4) and (9)P37023ACVRL1Serine/threonine-protein kinase receptor R3(6)P01019AGTAngiotensinogen(9)P09917ALOX5Polyunsaturated fatty acid 5-lipoxygenase(6)P01185AVPVasopressin-neurophysin 2-copeptin(9)Q07812BAXApoptosis regulator BAX(6)O15392BIRC5Baculoviral IAP repeat-containing protein 5(6)P55211CASP9Caspase-9(7)P20248CCNA2Cyclin-A2(7)P14635CCNB1G2/mitotic-specific cyclin-B1(7)P13073COX4I1Cytochrome c oxidase subunit 4 isoform 1, mitochondrial(9)P02741CRPC-reactive protein(4)P02778CXCL10C-X-C motif chemokine 10(7)P00533EGFREpidermal growth factor receptor(11)P15407FOSL1Fos-related antigen 1(6)P15408FOSL2Fos-related antigen 2(6)Q13224GRIN2BGlutamate receptor ionotropic, NMDA 2B(9)Q16665HIF1AHypoxia-inducible factor 1-alpha(9)P04792HSPB1Heat shock protein beta-1(7)P05362ICAM1Intercellular adhesion molecule 1(7)P35225IL13Interleukin-13(7)P05112IL4Interleukin-4(7)P01374LTALymphotoxin-alpha(6)P01106MYCMyc proto-oncogene protein(8)P29474NOS3Nitric oxide synthase, endothelial(9)O14684PTGESProstaglandin E synthase(6)P35354PTGS2Prostaglandin G/H synthase 2(2)Q15858SCN9ASodium channel protein type 9 subunit alpha(4)P16581SELEE-selectin(7)P07204THBDThrombomodulin(6)Q8NER1TRPV1Transient receptor potential cation channel subfamily V member 1(9)P11473VDRVitamin D3 receptor(9)(1): Degree, (2): Maximal Clique Centrality (MCC), (3): Betweenness, (4): BottleNeck, (5): Closeness, (6): ClusteringCoefficient, (7): Density of Maximum Neighborhood Component (DMNC), (8): Edge Percolated Component (EPC), (9): EcCentricity, (10): Maximum Neighborhood Component (MNC), (11): Radiality, and (12): Stress


### 3.4. Results of GO and KEGG Pathway Enrichment Analysis

GO and KEGG pathway enrichment analysis are important tools to explore the biological processes and signal pathways involved in clinical drug therapy for diseases, as well as key genes. Through the ClusterProfiler package in R, we performed GO enrichment analysis on the 189 potential therapeutic targets for BYHWD anti-SCI and obtained 2842 GO items (after adjustment,P < 0.05), including 2567 BP items, 99 CC items, and 176 MF items. We combined the top 10 enrichment results for each of GO-BP, GO-CC, and GO-MF by a bar graph (Figure 5). The basic information of the top 10 GO enrichment items is shown in Table 2.Figure 5
The bar chart of top 10 GO (BP, CC, and MF) enriched items.Table 2
Top 10 items of Gene Ontology (GO) enrichment analysis.
GO itemsIDDescriptionP valueAdjustedP valueGene numberBiological processGO:0010038Response to metal ion2.68739E−451.20529E−4152Biological processGO:0031667Response to nutrient levels2.05974E−404.61896E−3754Biological processGO:0032496Response to lipopolysaccharide5.34106E−387.98488E−3545Biological processGO:0002237Response to molecule of bacterial origin3.1151E−373.4928E−3445Biological processGO:0062197Cellular response to chemical stress1.44121E−351.29277E−3244Biological processGO:0046677Response to antibiotic2.54411E−341.90172E−3142Biological processGO:0006979Response to oxidative stress3.05284E−341.956E−3147Biological processGO:0000302Response to reactive oxygen species9.13447E−345.12101E−3137Biological processGO:0007584Response to nutrient2.40193E−331.19696E−3036Biological processGO:0034599Cellular response to oxidative stress6.07803E−322.726E−2939Cellular componentGO:0045121Membrane raft3.92204E−137.60946E−1122Cellular componentGO:0098857Membrane microdomain4.18102E−137.60946E−1122Cellular componentGO:0098589Membrane region8.84519E−131.07322E−1022Cellular componentGO:0031983Vesicle lumen4.00266E−103.64242E−0819Cellular componentGO:1902911Protein kinase complex5.90811E−104.10759E−0812Cellular componentGO:0090575RNA polymerase II transcription regulator complex7.53142E−104.10759E−0814Cellular componentGO:1902554Serine/threonine-protein kinase complex7.89922E−104.10759E−0811Cellular componentGO:0034774Secretory granule lumen2.08462E−099.48503E−0818Cellular componentGO:0060205Cytoplasmic vesicle lumen2.53536E−091.02541E−0718Cellular componentGO:0005667Transcription regulator complex3.48128E−091.26718E−0720Molecular functionGO:0061629RNA polymerase II-specific DNA-binding transcription factor binding1.80793E−117.91763E−0920Molecular functionGO:0140297DNA-binding transcription factor binding2.70226E−117.91763E−0922Molecular functionGO:0046906Tetrapyrrole binding4.34881E−118.49468E−0915Molecular functionGO:0020037Heme binding1.91194E−102.801E−0814Molecular functionGO:0009055Electron transfer activity2.62858E−103.08069E−0813Molecular functionGO:0004879Nuclear receptor activity1.43492E−091.20123E−079Molecular functionGO:0098531Ligand-activated transcription factor activity1.43492E−091.20123E−079Molecular functionGO:0044389Ubiquitin-like protein ligase binding5.74415E−094.20759E−0718Molecular functionGO:0003707Steroid hormone receptor activity7.34487E−094.78232E−079Molecular functionGO:0097110Scaffold protein binding1.18514E−085.67261E−079In the same way, through the ClusterProfiler package in R, we performed KEGG pathway enrichment analysis on the 189 potential therapeutic targets for BYHWD anti-SCI and obtained 177 KEGG pathways (after adjustment,P < 0.05). The top 20 KEGG enrichment pathways are listed below, including AGE-RAGE signaling pathway in diabetic complications (hsa04933), Fluid shear stress and atherosclerosis (hsa05418), Prostate cancer (hsa05215), IL-17 signaling pathway (hsa04657), Hepatitis B (hsa05161), TNF signaling pathway (hsa04668), Kaposi sarcoma-associated herpesvirus infection (hsa05167), Human cytomegalovirus infection (hsa05163), Hepatitis C (hsa05160), Human T-cell leukemia virus 1 infection (hsa05166), Small cell lung cancer (hsa05222), Pancreatic cancer (hsa05212), Measles (hsa05162), Epstein-Barr virus infection (hsa05169), Nonalcoholic fatty liver disease (hsa04932), Influenza A (hsa05164), Endocrine resistance (hsa01522), Bladder cancer (hsa05219), Platinum drug resistance (hsa01524), and Apoptosis (hsa04210). The visualization of the abovementioned KEGG enrichment pathway is realized in Figure 6. More KEGG pathways were enriched, and another fine screening of the pathways was required. Therefore, the SCI-related literature was searched in the PubMed database, and then the retrieved potentially relevant pathways were compared with the enriched 177 pathways, and finally, 15 relatively relevant pathways for SCI were obtained. Subsequently, we used Cytoscape software to visualize the network relationship between potential therapeutic targets for BYHWD anti-SCI and pathways (Figure 7). Table 3 shows the basic information of 15 pathways that may be relevant for SCI.Figure 6
The bar chart of top 20 KEGG enriched pathways.Figure 7
Potential therapeutic targets-pathways” network. The blue diamond represents the top 10 hub genes; the cyan rectangle represents the potential therapeutic targets for BYHWD anti-SCI; the yellow ellipse represents potential SCI-related pathways supported by literature; the green ellipse represents possible pathways related to SCI.Table 3
The enriched 15 possible related pathways for SCI.
IDDescriptionP valueAdjustedP valueGene numberhsa04933AGE-RAGE signaling pathway in diabetic complications8.58479E−342.24063E−3135hsa04668TNF signaling pathway8.83261E−243.84219E−2229hsa04210Apoptosis1.13129E−151.42828E−1424hsa05022Pathways of neurodegeneration - multiple diseases1.12387E−141.0476E−1341hsa04620Toll-like receptor signaling pathway5.71758E−144.52208E−1320hsa04151PI3K-Akt signaling pathway1.43931E−131.00595E−1234hsa04066HIF-1 signaling pathway1.4646E−131.00595E−1220hsa04010MAPK signaling pathway5.84978E−123.39287E−1129hsa04064NF-kappa B signaling pathway7.32616E−124.1568E−1118hsa01521EGFR tyrosine kinase inhibitor resistance1.25843E−095.13202E−0914hsa04012ErbB signaling pathway3.21508E−081.0622E−0713hsa04630JAK-STAT signaling pathway4.40433E−071.27726E−0616hsa04310Wnt signaling pathway0.0007489370.00151529211hsa04152AMPK signaling pathway0.0012301370.0024072949hsa04150mTOR signaling pathway0.0020914930.00392719110


### 3.5. Results of Molecular Docking of Hub Genes and Key Active Compounds

Sankey diagram was used to establish one-to-one correspondence between the top 10 hub genes (AKT1, IL6, MAPK1, TNF, TP53, VEGFA, CASP3, ALB, MAPK8, and JUN) for BYHWD anti-SCI and BYHWD-related compounds, providing the basis for the next molecular docking analysis. Among them, the gene that targets the most active compounds is TNF, and the key active compounds that target the most hub genes are quercetin and luteolin (Figure8). As seen in the Sankey diagram, 6 of the 33 active compounds (quercetin, luteolin, beta-carotene, kaempferol, baicalein, and beta-sitosterol) targeted more hub genes, so we decided to perform the molecular docking between these 6 active compounds and top 10 hub genes. The results of the molecular docking score are presented as a heat map (Figure 9), with a range of −6.4 to −10.9 kcal·mol−1, representing the large molecule protein receptor (hub genes) that binds well to the small molecule ligand (key active compounds). As shown in Figure 10, the molecular docking process of each key active compound and its corresponding protein encoded by the hub gene with the best docking affinity was visualized using Pymol. Table 4 shows the basic information of active compounds targeting hub genes in BYHWD.Figure 8
Sankey diagram of “herbs-active compounds-hub genes.” The band between the two bars represents a targeting relationship.Figure 9
Heatmap of the scores of molecular docking in this study.Figure 10
Results of molecular docking between hub genes and key active compounds. (a) Quercetin to AKT1. (b) Kaempferol to TNF. (c) Beta-sitosterol to AKT1. (d) Beta-carotene to AKT1. (e) Baicalein to AKT1. (f) Luteolin to MAPK8.
(a)(b)(c)(d)(e)(f)Table 4
Basic information of active compounds targeting hub genes in BYHWD.
Molecule IDMolecule namePubChem CIDOB (%)DL2D StructureSource databaseSourceTargeted hub genesMOL000098Quercetin528034346.430.28TCMSPHonghua; HuangqiAKT1, CASP3, IL6, JUN, MAPK1, and VEGFAMOL000006Luteolin528044536.160.25TCMSPHonghuaAKT1, CASP3, IL6, JUN, MAPK1, and VEGFAMOL002773Beta-carotene528048937.180.58TCMSPHonghuaAKT1, CASP3, ALB, JUN, and VEGFAMOL000422Kaempferol528086341.880.24TCMSPHonghua; HuangqiAKT1, CASP3, JUN, and MAPK8MOL002714Baicalein528160533.520.21TCMSPChishao; HonghuaAKT1, CASP3, and VEGFAMOL000358Beta-sitosterol22228436.910.75TCMSPChishao, Danggui, Honghua, and TaorenCASP3 and JUNMOL001002Ellagic acid528185543.060.43TCMSPChishaoVEGFAMOL001924Paeoniflorin44253453.870.79TCMSPChishaoIL6MOL000392Formononetin528037869.670.21TCMSPHuangqiJUN—Anisic acid7478——BATMAN-TCMDangguiTNF and TP53—Caffeic acid dimethyl ether717531——BATMAN-TCMChuanxiongTNF—Cibarian100275——BATMAN-TCMChuanxiongALB—Linoleic acid5280450——BATMAN-TCMChuanxiongTNF—O-Cresol335——BATMAN-TCMChuanxiongTNF—O-Ethylphenol6997——BATMAN-TCMChuanxiongTNF—Oleic acid445639——BATMAN-TCMChuanxiongTNF—Thymol6989——BATMAN-TCMChuanxiongTNF—1-Methyl-2-dodecyl-4-(1h)-quinolone5319601——BATMAN-TCMDangguiJUN—4-Ethylresorcinol17927——BATMAN-TCMDangguiTNF—Angelicin10658——BATMAN-TCMDangguiAKT1—Carvacrol10364——BATMAN-TCMDangguiTNF—Dihydropinosylvin442700——BATMAN-TCMDangguiTNF—Dodecenoic acid96204——BATMAN-TCMDangguiTNF—Ethyl-P-methoxycinnamate5281783——BATMAN-TCMDangguiTNF—M-Cresol342——BATMAN-TCMDangguiTNF—M-Ethylphenol12101——BATMAN-TCMDangguiTNF—P-Cresol2879——BATMAN-TCMDangguiTNF—P-Ethylphenol31242——BATMAN-TCMDangguiTNF—Sebiferic acid5321206——BATMAN-TCMDangguiTNF—Astragaloside I13996685——BATMAN-TCMHuangqiALB—Astragaloside II13996693——BATMAN-TCMHuangqiALB—Sucrose5988——BATMAN-TCMHuangqiTNF—Cordycepin6303——BATMAN-TCMTaorenAKT1OB, oral bioavailability; DL, drug-likeness.


## 3.1. Acquisition Results of BYHWD-Related Targets and SCI-Related Targets

A total of 775 compounds (BATMAN-TCM: 219; after deduplication: 189) in BYHWD were preliminarily retrieved from the TCMSP database and BATMAN-TCM platform (score cutoff >25; adjustedP value < 0.05), of which 119 (BATMAN-TCM: 8) were from Chishao, 189 (87) were from Chuanxiong, 125 (77) were from Danggui, 87 (23) were from Huangqi, 189 (22) were from Honghua, and 66 (2) were from Taoren. Relevant compounds and targets’ information for Dilong was not retrieved in either database. After setting the ADME criteria (OB ≥ 30% and DL ≥ 0.18), a total of 103 active compounds (after removing duplicates and nontarget compounds: 60) were screened from the TCMSP database, of which 29 (After removing nontarget compounds: 14) were from Chishao, 7 (6) were from Chuanxiong, 2 (2) were from Danggui, 20 (17) were from Huangqi, 22 (17) were from Honghua, and 23 (19) were from Taoren. The BATMAN-TCM platform provides the gene names directly, while the TCMSP database first provides the target protein names of the active compounds and therefore needs to be converted to gene names via the UniProt database. Finally, after integrating the results of the two databases, we obtained 244 BYHWD-related active compounds and 744 BYHWD-related targets.After deduplication, a total of 923 SCI-related targets were obtained by mining the DisGeNET (number: 84), TTD (number: 5), CTD (number: 553), GeneCards (number: 432), and MalaCards (number: 19) databases.


## 3.2. Acquisition Results of Potential Therapeutic Targets for BYHWD Anti-SCI

The overlapping targets of BYHWD-related targets and SCI-related targets were considered as potential therapeutic targets for BYHWD anti-SCI. As shown in Figure2, we screened the overlapping targets of both by constructing a Venn diagram, and a total of 189 targets were obtained. Through multiple databases mining and construction of the Venn diagram, we obtained BYHWD-related compounds and potential therapeutic targets for BYHWD anti-SCI and finally used Cytoscape software to construct a “BYHWD-active compounds-target genes-SCI” network (Figure 3). In Figure 3, the lines between two nodes represent the existence of mutual relationships, and the larger the node is, the more relationships exist.Figure 2
Venn diagram for BYHWD-related targets and SCI-related targets.Figure 3
BYHWD-active compounds-target genes-SCI” network. The purple diamond represents the herbs contained in BYHWD; the red octagon represents SCI; the cyan rectangle represents the potential therapeutic targets for BYHWD anti-SCI; the light green ellipse represents the active compounds contained in BYHWD; the blue ellipse represents the common active compounds of the herb.


## 3.3. Acquisition Results of Hub Genes for BYHWD Anti-SCI and PPI Network Construction

In the STRING database, we first set the species to “Homo sapiens” and then entered 189 potential therapeutic targets for BYHWD anti-SCI to obtain the PPI network (Figure 4(a)), which involved 189 nodes and 3717 edges. The obtained TSV files were then imported into Cytoscape software (version 3.7.2) for further analysis and visualization (Figure 4(b)). To obtain the hub genes for BYHWD anti-SCI, based on the above the PPI network, we used the CytoHubba plug-in of Cytoscape software, which currently contains 12 topological analysis methods. We targeted the top 10 hub genes for BYHWD anti-SCI, and the results of the 12 algorithms each contained 10 genes. We integrated and analyzed the results to obtain a total of 45 different genes, which were then sorted by the number of algorithms to which they were attributed. Table 1 shows the basic information on hub genes for BYHWD anti-SCI. The top 10 hub genes (AKT1, IL6, MAPK1, TNF, TP53, VEGFA, CASP3, ALB, MAPK8, and JUN) sorted are completely consistent with some algorithms (Degree, Closeness, and MNC) (Figure 4(c)).Figure 4
PPI network and top 10 hub genes for BYHWD anti-SCI. (a) PPI network constructed with STRING. (b) PPI network constructed with Cytoscape software (the darker the node color, the higher the number of connected proteins). (c) Top 10 hub genes for BYHWD anti-SCI were obtained by using the Degree algorithm.
(a)(b)(c)Table 1
Basic information of hub genes based on 12 topological algorithms.
UniProt IDGene symbolProtein namesAlgorithmsP31749AKT1RAC-alpha serine/threonine-protein kinase(1), (2), (3), (4), (5), (8), (10), (11), and (12)P05231IL6Interleukin-6(1), (2), (3), (4), (5), (8), (10), (11), and (12)P28482MAPK1Mitogen-activated protein kinase 1(1), (2), (3), (5), (8), (9), (10), (11), and (12)P01375TNFTumor necrosis factor(1), (2), (3), (4), (5), (8), (10), (11), and (12)P04637TP53Cellular tumor antigen p53(1), (2), (3), (4), (5), (8), (10), (11), and (12)P15692VEGFAVascular endothelial growth factor A(1), (2), (3), (4), (5), (8), (10), (11), and (12)P42574CASP3Caspase-3(1), (2), (3), (5), (8), (10), (11), and (12)P02768ALBAlbumin(1), (3), (5), (8), (10), (11), and (12)P45983MAPK8Mitogen-activated protein kinase 8(1), (2), (4), (5), (8), and (10)P05412JUNTranscription factor AP-1(1), (2), (5), (10), and (11)P04040CATCatalase(3), (4), and (12)P99999CYCSCytochrome c(3) and (12)P09238MMP10Stromelysin-2(6) and (7)P07101THTyrosine 3-monooxygenase(4) and (9)P37023ACVRL1Serine/threonine-protein kinase receptor R3(6)P01019AGTAngiotensinogen(9)P09917ALOX5Polyunsaturated fatty acid 5-lipoxygenase(6)P01185AVPVasopressin-neurophysin 2-copeptin(9)Q07812BAXApoptosis regulator BAX(6)O15392BIRC5Baculoviral IAP repeat-containing protein 5(6)P55211CASP9Caspase-9(7)P20248CCNA2Cyclin-A2(7)P14635CCNB1G2/mitotic-specific cyclin-B1(7)P13073COX4I1Cytochrome c oxidase subunit 4 isoform 1, mitochondrial(9)P02741CRPC-reactive protein(4)P02778CXCL10C-X-C motif chemokine 10(7)P00533EGFREpidermal growth factor receptor(11)P15407FOSL1Fos-related antigen 1(6)P15408FOSL2Fos-related antigen 2(6)Q13224GRIN2BGlutamate receptor ionotropic, NMDA 2B(9)Q16665HIF1AHypoxia-inducible factor 1-alpha(9)P04792HSPB1Heat shock protein beta-1(7)P05362ICAM1Intercellular adhesion molecule 1(7)P35225IL13Interleukin-13(7)P05112IL4Interleukin-4(7)P01374LTALymphotoxin-alpha(6)P01106MYCMyc proto-oncogene protein(8)P29474NOS3Nitric oxide synthase, endothelial(9)O14684PTGESProstaglandin E synthase(6)P35354PTGS2Prostaglandin G/H synthase 2(2)Q15858SCN9ASodium channel protein type 9 subunit alpha(4)P16581SELEE-selectin(7)P07204THBDThrombomodulin(6)Q8NER1TRPV1Transient receptor potential cation channel subfamily V member 1(9)P11473VDRVitamin D3 receptor(9)(1): Degree, (2): Maximal Clique Centrality (MCC), (3): Betweenness, (4): BottleNeck, (5): Closeness, (6): ClusteringCoefficient, (7): Density of Maximum Neighborhood Component (DMNC), (8): Edge Percolated Component (EPC), (9): EcCentricity, (10): Maximum Neighborhood Component (MNC), (11): Radiality, and (12): Stress


## 3.4. Results of GO and KEGG Pathway Enrichment Analysis

GO and KEGG pathway enrichment analysis are important tools to explore the biological processes and signal pathways involved in clinical drug therapy for diseases, as well as key genes. Through the ClusterProfiler package in R, we performed GO enrichment analysis on the 189 potential therapeutic targets for BYHWD anti-SCI and obtained 2842 GO items (after adjustment,P < 0.05), including 2567 BP items, 99 CC items, and 176 MF items. We combined the top 10 enrichment results for each of GO-BP, GO-CC, and GO-MF by a bar graph (Figure 5). The basic information of the top 10 GO enrichment items is shown in Table 2.Figure 5
The bar chart of top 10 GO (BP, CC, and MF) enriched items.Table 2
Top 10 items of Gene Ontology (GO) enrichment analysis.
GO itemsIDDescriptionP valueAdjustedP valueGene numberBiological processGO:0010038Response to metal ion2.68739E−451.20529E−4152Biological processGO:0031667Response to nutrient levels2.05974E−404.61896E−3754Biological processGO:0032496Response to lipopolysaccharide5.34106E−387.98488E−3545Biological processGO:0002237Response to molecule of bacterial origin3.1151E−373.4928E−3445Biological processGO:0062197Cellular response to chemical stress1.44121E−351.29277E−3244Biological processGO:0046677Response to antibiotic2.54411E−341.90172E−3142Biological processGO:0006979Response to oxidative stress3.05284E−341.956E−3147Biological processGO:0000302Response to reactive oxygen species9.13447E−345.12101E−3137Biological processGO:0007584Response to nutrient2.40193E−331.19696E−3036Biological processGO:0034599Cellular response to oxidative stress6.07803E−322.726E−2939Cellular componentGO:0045121Membrane raft3.92204E−137.60946E−1122Cellular componentGO:0098857Membrane microdomain4.18102E−137.60946E−1122Cellular componentGO:0098589Membrane region8.84519E−131.07322E−1022Cellular componentGO:0031983Vesicle lumen4.00266E−103.64242E−0819Cellular componentGO:1902911Protein kinase complex5.90811E−104.10759E−0812Cellular componentGO:0090575RNA polymerase II transcription regulator complex7.53142E−104.10759E−0814Cellular componentGO:1902554Serine/threonine-protein kinase complex7.89922E−104.10759E−0811Cellular componentGO:0034774Secretory granule lumen2.08462E−099.48503E−0818Cellular componentGO:0060205Cytoplasmic vesicle lumen2.53536E−091.02541E−0718Cellular componentGO:0005667Transcription regulator complex3.48128E−091.26718E−0720Molecular functionGO:0061629RNA polymerase II-specific DNA-binding transcription factor binding1.80793E−117.91763E−0920Molecular functionGO:0140297DNA-binding transcription factor binding2.70226E−117.91763E−0922Molecular functionGO:0046906Tetrapyrrole binding4.34881E−118.49468E−0915Molecular functionGO:0020037Heme binding1.91194E−102.801E−0814Molecular functionGO:0009055Electron transfer activity2.62858E−103.08069E−0813Molecular functionGO:0004879Nuclear receptor activity1.43492E−091.20123E−079Molecular functionGO:0098531Ligand-activated transcription factor activity1.43492E−091.20123E−079Molecular functionGO:0044389Ubiquitin-like protein ligase binding5.74415E−094.20759E−0718Molecular functionGO:0003707Steroid hormone receptor activity7.34487E−094.78232E−079Molecular functionGO:0097110Scaffold protein binding1.18514E−085.67261E−079In the same way, through the ClusterProfiler package in R, we performed KEGG pathway enrichment analysis on the 189 potential therapeutic targets for BYHWD anti-SCI and obtained 177 KEGG pathways (after adjustment,P < 0.05). The top 20 KEGG enrichment pathways are listed below, including AGE-RAGE signaling pathway in diabetic complications (hsa04933), Fluid shear stress and atherosclerosis (hsa05418), Prostate cancer (hsa05215), IL-17 signaling pathway (hsa04657), Hepatitis B (hsa05161), TNF signaling pathway (hsa04668), Kaposi sarcoma-associated herpesvirus infection (hsa05167), Human cytomegalovirus infection (hsa05163), Hepatitis C (hsa05160), Human T-cell leukemia virus 1 infection (hsa05166), Small cell lung cancer (hsa05222), Pancreatic cancer (hsa05212), Measles (hsa05162), Epstein-Barr virus infection (hsa05169), Nonalcoholic fatty liver disease (hsa04932), Influenza A (hsa05164), Endocrine resistance (hsa01522), Bladder cancer (hsa05219), Platinum drug resistance (hsa01524), and Apoptosis (hsa04210). The visualization of the abovementioned KEGG enrichment pathway is realized in Figure 6. More KEGG pathways were enriched, and another fine screening of the pathways was required. Therefore, the SCI-related literature was searched in the PubMed database, and then the retrieved potentially relevant pathways were compared with the enriched 177 pathways, and finally, 15 relatively relevant pathways for SCI were obtained. Subsequently, we used Cytoscape software to visualize the network relationship between potential therapeutic targets for BYHWD anti-SCI and pathways (Figure 7). Table 3 shows the basic information of 15 pathways that may be relevant for SCI.Figure 6
The bar chart of top 20 KEGG enriched pathways.Figure 7
Potential therapeutic targets-pathways” network. The blue diamond represents the top 10 hub genes; the cyan rectangle represents the potential therapeutic targets for BYHWD anti-SCI; the yellow ellipse represents potential SCI-related pathways supported by literature; the green ellipse represents possible pathways related to SCI.Table 3
The enriched 15 possible related pathways for SCI.
IDDescriptionP valueAdjustedP valueGene numberhsa04933AGE-RAGE signaling pathway in diabetic complications8.58479E−342.24063E−3135hsa04668TNF signaling pathway8.83261E−243.84219E−2229hsa04210Apoptosis1.13129E−151.42828E−1424hsa05022Pathways of neurodegeneration - multiple diseases1.12387E−141.0476E−1341hsa04620Toll-like receptor signaling pathway5.71758E−144.52208E−1320hsa04151PI3K-Akt signaling pathway1.43931E−131.00595E−1234hsa04066HIF-1 signaling pathway1.4646E−131.00595E−1220hsa04010MAPK signaling pathway5.84978E−123.39287E−1129hsa04064NF-kappa B signaling pathway7.32616E−124.1568E−1118hsa01521EGFR tyrosine kinase inhibitor resistance1.25843E−095.13202E−0914hsa04012ErbB signaling pathway3.21508E−081.0622E−0713hsa04630JAK-STAT signaling pathway4.40433E−071.27726E−0616hsa04310Wnt signaling pathway0.0007489370.00151529211hsa04152AMPK signaling pathway0.0012301370.0024072949hsa04150mTOR signaling pathway0.0020914930.00392719110


## 3.5. Results of Molecular Docking of Hub Genes and Key Active Compounds

Sankey diagram was used to establish one-to-one correspondence between the top 10 hub genes (AKT1, IL6, MAPK1, TNF, TP53, VEGFA, CASP3, ALB, MAPK8, and JUN) for BYHWD anti-SCI and BYHWD-related compounds, providing the basis for the next molecular docking analysis. Among them, the gene that targets the most active compounds is TNF, and the key active compounds that target the most hub genes are quercetin and luteolin (Figure8). As seen in the Sankey diagram, 6 of the 33 active compounds (quercetin, luteolin, beta-carotene, kaempferol, baicalein, and beta-sitosterol) targeted more hub genes, so we decided to perform the molecular docking between these 6 active compounds and top 10 hub genes. The results of the molecular docking score are presented as a heat map (Figure 9), with a range of −6.4 to −10.9 kcal·mol−1, representing the large molecule protein receptor (hub genes) that binds well to the small molecule ligand (key active compounds). As shown in Figure 10, the molecular docking process of each key active compound and its corresponding protein encoded by the hub gene with the best docking affinity was visualized using Pymol. Table 4 shows the basic information of active compounds targeting hub genes in BYHWD.Figure 8
Sankey diagram of “herbs-active compounds-hub genes.” The band between the two bars represents a targeting relationship.Figure 9
Heatmap of the scores of molecular docking in this study.Figure 10
Results of molecular docking between hub genes and key active compounds. (a) Quercetin to AKT1. (b) Kaempferol to TNF. (c) Beta-sitosterol to AKT1. (d) Beta-carotene to AKT1. (e) Baicalein to AKT1. (f) Luteolin to MAPK8.
(a)(b)(c)(d)(e)(f)Table 4
Basic information of active compounds targeting hub genes in BYHWD.
Molecule IDMolecule namePubChem CIDOB (%)DL2D StructureSource databaseSourceTargeted hub genesMOL000098Quercetin528034346.430.28TCMSPHonghua; HuangqiAKT1, CASP3, IL6, JUN, MAPK1, and VEGFAMOL000006Luteolin528044536.160.25TCMSPHonghuaAKT1, CASP3, IL6, JUN, MAPK1, and VEGFAMOL002773Beta-carotene528048937.180.58TCMSPHonghuaAKT1, CASP3, ALB, JUN, and VEGFAMOL000422Kaempferol528086341.880.24TCMSPHonghua; HuangqiAKT1, CASP3, JUN, and MAPK8MOL002714Baicalein528160533.520.21TCMSPChishao; HonghuaAKT1, CASP3, and VEGFAMOL000358Beta-sitosterol22228436.910.75TCMSPChishao, Danggui, Honghua, and TaorenCASP3 and JUNMOL001002Ellagic acid528185543.060.43TCMSPChishaoVEGFAMOL001924Paeoniflorin44253453.870.79TCMSPChishaoIL6MOL000392Formononetin528037869.670.21TCMSPHuangqiJUN—Anisic acid7478——BATMAN-TCMDangguiTNF and TP53—Caffeic acid dimethyl ether717531——BATMAN-TCMChuanxiongTNF—Cibarian100275——BATMAN-TCMChuanxiongALB—Linoleic acid5280450——BATMAN-TCMChuanxiongTNF—O-Cresol335——BATMAN-TCMChuanxiongTNF—O-Ethylphenol6997——BATMAN-TCMChuanxiongTNF—Oleic acid445639——BATMAN-TCMChuanxiongTNF—Thymol6989——BATMAN-TCMChuanxiongTNF—1-Methyl-2-dodecyl-4-(1h)-quinolone5319601——BATMAN-TCMDangguiJUN—4-Ethylresorcinol17927——BATMAN-TCMDangguiTNF—Angelicin10658——BATMAN-TCMDangguiAKT1—Carvacrol10364——BATMAN-TCMDangguiTNF—Dihydropinosylvin442700——BATMAN-TCMDangguiTNF—Dodecenoic acid96204——BATMAN-TCMDangguiTNF—Ethyl-P-methoxycinnamate5281783——BATMAN-TCMDangguiTNF—M-Cresol342——BATMAN-TCMDangguiTNF—M-Ethylphenol12101——BATMAN-TCMDangguiTNF—P-Cresol2879——BATMAN-TCMDangguiTNF—P-Ethylphenol31242——BATMAN-TCMDangguiTNF—Sebiferic acid5321206——BATMAN-TCMDangguiTNF—Astragaloside I13996685——BATMAN-TCMHuangqiALB—Astragaloside II13996693——BATMAN-TCMHuangqiALB—Sucrose5988——BATMAN-TCMHuangqiTNF—Cordycepin6303——BATMAN-TCMTaorenAKT1OB, oral bioavailability; DL, drug-likeness.


## 4. Discussion

According to the report of the National SCI Statistics Center, there are about 17,000 new cases of SCI in the United States each year, and the annual incidence of SCI is about 54 cases per million people [1]. Over the past century, advances in understanding the mechanisms of SCI have changed the clinical management strategies for SCI, including surgical procedures, supportive measures, and rehabilitative training that have improved neurological outcomes and reduced morbidity in SCI patients [1]. However, there is still a lack of effective treatments for SCI and thus has been stimulating the desire of researchers to explore and try to find effective treatments, including drugs, surgical methods, and rehabilitation programs. As we all know, TCM has a long history of being used to treat various diseases in China and other Asian countries [28]. Among them, BYHWD, Jisuikang, single herb, and Governor Vessel electroacupuncture have all been shown in relevant studies to be used in the treatment of SCI [4, 5, 28–30].BYHWD, which comes from Qing Dynasty medical classics, consists of seven herbs: Chishao, Chuanxiong, Danggui, Dilong, Huangqi, Honghua, and Taoren [5]. Based on the theory of TCM, Huangqi is the main ingredient in this formula, which has the function of tonifying Qi and nourishing blood, eliminating blood stasis without harming the proper function; Danggui is the secondary herb, which has the effect of invigorating blood and harmonizing blood; Chishao, Chuanxiong, Dilong, Honghua, and Taoren are the adjuvant herbs in this formula; all the five herbs have the function of activating blood circulation and removing blood stasis, promoting the flow of Qi and blood circulation; the combination of the five herbs not only activates blood circulation and removes blood stasis but also helps the Qi and blood that the main herbs benefit to reach the whole body [3, 31]. Studies have shown that BYHWD in combination with neural stem cells (NSCs) attenuates the process of demyelination or improves the recovery of myelin and exerts a synergistic effect on the neurological function recovery [5]. It has been shown that BYHWD restores hindlimb motor function in rats with SCI, and the neuroprotective effect is associated with the regulation of apoptosis-related protein expression [32]. BYHWD rescues axotomized neurons and promotes functional recovery after spinal cord injury in rats. A study published in 2008 showed that BYHWD protected the injured neurons and promoted functional recovery after SCI in rats [33]. However, since the specific mechanism of BYHWD for the treatment of SCI is still unclear due to its complex composition and wide range of effects, we decided to use a network pharmacology approach to investigate in this study and to lay the foundation for the next in-depth study.In the first step, BYHWD-related compounds and related targets were obtained. We obtained relevant compounds and targets by mining two herbal-related databases, TCMSP and BATMAN-TCM. The mining methods of these two databases are different: TCMSP database (ADME principle: OB ≥ 30 and DL ≥ 0.18) first obtains BYHWD-related compounds, then obtains related target proteins, and finally needs to be converted into standard gene names by UniProt database; BATMAN-TCM platform (score cutoff >25, adjustedP value <0.05) directly obtains BYHWD-related compounds and target genes. The screening process of TCMSP database is divided into five steps: (1) the relevant compounds of each component of BYHWD are directly obtained (Chishao: 119, Chuanxiong: 189, Danggui: 125, Huangqi: 87, Honghua: 189, Taoren: 66, and Dilong: 0); (2) the relevant compounds of each component are screened by ADME criteria (Chishao: 29, Chuanxiong: 7, Danggui: 2, Huangqi: 20, Honghua: 22, Taoren: 23, and Dilong: 0); (3) the above compounds are screened again by removing the null targets (Chishao: 14, Chuanxiong: 6, Danggui: 2, Huangqi: 17, Honghua: 17, Taoren: 19, and Dilong: 0); (4) the above compounds are finally screened by removing duplicates (total number of active compounds: 60); and (5) the target proteins corresponding to the obtained active compounds are converted into gene names (total number of targets: 227). And the screening process of BATMAN-TCM platform is divided into two steps: (1) the relevant compounds and targets of each component of BYHWD are directly obtained (Chishao: 8, Chuanxiong: 87, Danggui: 77, Huangqi: 23, Honghua: 22, Taoren: 2, and Dilong: 0), and (2) the relevant compounds of each component are screened after deduplication (total number of active compounds: 189; the total number of targets: 586). The results of the two databases were integrated and deduplicated to finally obtain 244 BYHWD-related active compounds and 744 BYHWD-related targets. In the second step, a total of 923 SCI-related targets were acquired by mining 5 disease databases, including DisGeNET (number: 84), TTD (number: 5), CTD (number: 553), GeneCards (number: 432), and MalaCards (number: 19). Next, the Venn diagram of BYHWD-related targets and SCI-related targets was plotted by the website online tool to obtain 189 potential therapeutic targets of BYHWD against SCI.In total, Cytoscape software was used four times in this study: (1) the “BYHWD-active compounds-target genes-SCI” network was constructed, involving 410 nodes and 1571 interactions; (2) the PPI network obtained from the STRING database was reconstructed here, involving 189 nodes and 3717 interactions; (3) the top 10 hub genes (AKT1, IL6, MAPK1, TNF, TP53, VEGFA, CASP3, ALB, MAPK8, and JUN) of BYHWD against SCI were obtained based on 12 topological analysis methods (Degree, MCC, Stress, ClusteringCoefficient, EcCentricity, BottleNeck, Closeness, Radiality, Betweenness, EPC, DMNC, and MNC) contained in the CytoHubba plug-in; (4) the “potential therapeutic targets-pathways” network was constructed, involving 347 nodes and 2650 interactions. Studies have shown that the protein encoded by TP53 and other factors regulate the regeneration, germination, and functional recovery of axons after central nervous system injury [34]; the rapid upregulation of CASP3 mRNA seen after SCI in rats may be associated with cell death in the spinal cord [35]; inflammatory mediators such as TNF-α and IL6 mediate the recruitment of inflammatory cells to the site of injury and by targeting these cytokines may be a potential strategy to reduce secondary injury in SCI [36, 37]; by targeting MAPK1, the overexpression of miRNA-433-5p protects motor dysfunction and inflammation after SCI [38].The Sankey diagram was used once, mainly to show the correspondence between the top 10 hub genes and the corresponding active compounds contained in BYHWD (quercetin, luteolin, beta-carotene, kaempferol, baicalein, beta-sitosterol, carvacrol, O-cresol, anisic acid, 1-methyl-2-dodecyl-4-(1h)-quinolone, 4-ethylresorcinol, angelicin, astragaloside I, astragaloside II, caffeic acid dimethyl ether, cibarian, cordycepin, dihydropinosylvin, dodecenoic acid, ellagic acid, ethyl-P-methoxycinnamate, formononetin, linoleic acid, M-cresol, M-ethylphenol, O-ethylphenol, oleic acid, paeoniflorin, P-cresol, P-ethylphenol, sebiferic acid, sucrose, and thymol). Among them, quercetin targeted six hub genes, luteolin targeted six hub genes, beta-carotene targeted five hub genes, kaempferol targeted four hub genes, baicalein targeted three hub genes, beta-sitosterol targeted two hub genes, anisic acid targeted two hub genes, and the others targeted one hub gene. The experimental results showed that quercetin exerts a protective effect on the spinal cord by inhibiting the activation of the p38MAPK/iNOS signaling pathway in SCI rats and subsequently regulating secondary oxidative stress [39]; the combination of luteolin and palmitoylethanolamide reduced autophagy in SCI [40]; beta-carotene effectively reduced the course of secondary damage events after SCI by blocking NF-κB pathway activation [41]; baicalein may alleviate the harm caused by SCI by activating PI3K and inducing autophagy to reduce neuronal apoptosis [42].Based on the results of the Sankey diagram, we performed molecular docking of the top 10 genes and the top 6 compounds in terms of the number of targeted genes (quercetin, luteolin, beta-carotene, kaempferol, baicalein, and beta-sitosterol). The results of molecular docking of 60 cohorts of hub target proteins and key active compounds showed a good binding affinity for both. We used R software twice, once to present the scores of the above molecular docking as a heat map, and another time to perform GO and KEGG pathway enrichment analysis using the ClusterProfiler package.The top 10 results of the GO enrichment analysis of 189 potential therapeutic targets for BYHWD anti-SCI, based on the adjustedP value (from small to large), are as follows: response to metal ion (GO-BP:0010038), response to nutrient levels (GO-BP:0031667), response to lipopolysaccharide (GO-BP:0032496), response to molecule of bacterial origin (GO-BP:0002237), cellular response to chemical stress (GO-BP: 0062197), response to antibiotic (GO-BP:0046677), response to oxidative stress (GO-BP: 0006979), response to reactive oxygen species (GO-BP:0000302), response to nutrient (GO-BP:0007584), cellular response to oxidative stress (GO-BP: 0034599), membrane raft (GO-CC:0045121), membrane microdomain (GO-CC:0098857), membrane region (GO-CC:0098589), vesicle lumen (GO-CC: 0031983), protein kinase complex (GO-CC: 1902911), RNA polymerase II transcription regulator complex (GO-CC: 0090575), serine/threonine-protein kinase complex (GO-CC:1902554), secretory granule lumen (GO-CC: 0034774), cytoplasmic vesicle lumen (GO-CC: 0060205), transcription regulator complex (GO-CC: 0005667), RNA polymerase II-specific DNA-binding transcription factor binding (GO-MF:0061629), DNA-binding transcription factor binding (GO-MF:0140297), tetrapyrrole binding (GO-MF:0046906), heme binding (GO-MF:0020037), electron transfer activity (GO-MF: 0009055), nuclear receptor activity (GO-MF:0004879), ligand-activated transcription factor activity (GO-MF:0098531), ubiquitin-like protein ligase binding (GO-MF: 0044389), steroid hormone receptor activity (GO-MF:0003707), and scaffold protein binding (GO-MF:0097110). The first 20 results of the (adjusted P value based) KEGG pathway enrichment analysis of 189 potential therapeutic targets of BYHWD against SCI have been presented, and we screened the enriched pathways again to find possible relevant pathways for SCI. By searching the relevant literature, we screened a total of 15 pathways with strong relevance, as follows: AGE-RAGE signaling pathway in diabetic complications (hsa04933) [43], TNF signaling pathway (hsa04668), Apoptosis (hsa04210) [44], Pathways of neurodegeneration-multiple diseases (hsa05022), Toll-like receptor signaling pathway (hsa04620) [45], PI3K-Akt signaling pathway (hsa04151) [46], HIF-1 signaling pathway (hsa04066) [47], MAPK signaling pathway (hsa04010) [48], NF-kappa B signaling pathway (hsa04064) [49], EGFR tyrosine kinase inhibitor resistance (hsa01521) [50], JAK-STAT signaling pathway (hsa04630) [51], Wnt signaling pathway (hsa04310) [43], AMPK signaling pathway (hsa04152) [44], and mTOR signaling pathway (hsa04150) [46]. The above pathways and hub genes are not in a disconnected relationship; rather, many hub genes are important components of some of the pathways. Taking PI3K-Akt signaling pathway as an example, it contains AKT1, MAPK1, IL6, TP53, and VEGFA. Based on the results of the above network pharmacological analysis, we hypothesize that the key active compounds of BYHWD may exert anti-SCI effects through the above pathways and hub genes.


### 4.1. Limitations

This study still has some limitations. First, the corresponding experimental validation is lacking. Second, the compounds, targets, and pathways contained in these databases may not be exhaustive. Finally, compounds screened based on ADME principles may be missing other important compounds.


### 4.2. Future Perspectives

This study uncovered the multicompound, multitarget, and multipathway characteristics of BYHWD through a network pharmacology approach. Among them, the key compounds include quercetin, luteolin, beta-carotene, kaempferol, baicalein, and beta-sitosterol; the hub genes include AKT1, IL6, MAPK1, TNF, TP53, VEGFA, CASP3, ALB, MAPK8, and JUN; and the important pathways are the 15 related pathways mentioned above. These results clearly give a general direction for SCI-related researchers, and there is no doubt that the anti-SCI mechanism of BYHWD is much more than the above-mentioned research directions, and we will find more accurate targets and pathways while exploring them one by one.


## 4.1. Limitations

This study still has some limitations. First, the corresponding experimental validation is lacking. Second, the compounds, targets, and pathways contained in these databases may not be exhaustive. Finally, compounds screened based on ADME principles may be missing other important compounds.


## 4.2. Future Perspectives

This study uncovered the multicompound, multitarget, and multipathway characteristics of BYHWD through a network pharmacology approach. Among them, the key compounds include quercetin, luteolin, beta-carotene, kaempferol, baicalein, and beta-sitosterol; the hub genes include AKT1, IL6, MAPK1, TNF, TP53, VEGFA, CASP3, ALB, MAPK8, and JUN; and the important pathways are the 15 related pathways mentioned above. These results clearly give a general direction for SCI-related researchers, and there is no doubt that the anti-SCI mechanism of BYHWD is much more than the above-mentioned research directions, and we will find more accurate targets and pathways while exploring them one by one.


## 5. Conclusion

Using a network pharmacology approach, we explored the potential mechanism of BYHWD as a classical herbal formula for treating nerve injury in the treatment of SCI. The top 10 core genes among the potential therapeutic targets for BYHWD anti-SCI include AKT1, IL6, MAPK1, TNF, TP53, VEGFA, CASP3, ALB, MAPK8, and JUN. Key active compounds of BYHWD against SCI include quercetin, luteolin, beta-carotene, kaempferol, baicalein, and beta-sitosterol. In addition, combining the KEGG pathway enrichment results and reviewing SCI-related literature, we obtained 15 SCI potentially relevant pathways, mainly including PI3K-Akt signaling pathway, MAPK signaling pathway, NF-kappa B signaling pathway, TNF signaling pathway, Apoptosis, Toll-like receptor signaling pathway, HIF-1 signaling pathway, JAK-STAT signaling pathway, Wnt signaling pathway, AMPK signaling pathway, and mTOR signaling pathway. The above analysis results demonstrated the multitarget, multicomponent, and multipathway characteristics of BYHWD in treating SCI, which laid the foundation for us to carry out specific experimental validation in the next step, and provided new ideas for researchers dedicated to the treatment of SCI with herbal extracts.



---

*Source: 9921534-2021-04-23.xml*