Evidence-Based Complementary and Alternative Medicine
Volume 2020 (2020), Article ID 7049812, 12 pages
https://doi.org/10.1155/2020/7049812

Research Article

Deciphering the Underlying Mechanism of Eucommiae Cortex against Osteoporotic Fracture by Network Pharmacology

Yongming Shuai, Zhili Jiang, Qiuwen Yuan, Shuqiang Tu, and Fanhui Zeng

Nanchang Hongdu Hospital of Traditional Chinese Medicine, 264 Minde Road, Donghu District, Nanchang, Jiangxi 330006, China

Correspondence should be addressed to Fanhui Zeng; zengfanghui510@163.com

Received 1 May 2020; Accepted 28 August 2020; Published 7 September 2020

Academic Editor: Antonella Fioravanti

Copyright © 2020 Yongming Shuai et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

Background. Du Zhong (DZ), or Eucommiae Cortex, a traditional Chinese herbal medicine, has been used to treat osteoporosis. Although it has been reported that DZ can improve bone mass in ovariectomized rats, its pharmacological mechanisms in treating osteoporotic fractures (OPF) remain unclear. Methods. In this study, we used a network pharmacological manner to explore its potential complicated mechanism in treating OPF. We obtained DZ compounds from TCMSP and BATMAN-TCM databases and collected potential targets of these compounds through target fishing based on TCMSP and BATMAN-TCM databases. Next, we collected the OPF targets by using CTD, GeneCards, OMIM, HPO, and GenCLiP 3 databases. And then the overlapping genes between DZ potential targets and OPF targets were used to build up the protein-protein interaction (PPI) network and to analyze their interactions and find out the big hub genes in this network. Subsequently, clusterProfiler package in R language was utilized to conduct the enrichment of Gene Ontology biological process and KEGG pathways. Results. There were totally 93 active compounds and 916 related targets in DZ. After the enrichment analysis, we collected top 25 cellular biological processes and top 25 pathways based on the adjusted value and found that the DZ anti-OPF effect was mainly associated with the regulation of ROS and inflammatory response. Furthermore, 64 hub genes in PPI network, such as MAPK1 (degree = 41), SRC (degree = 39), PIK3R1 (degree = 36), VEGFA (degree = 31), TP53 (degree = 29), EGFR (degree = 29), JUN (degree = 29), AGT (degree = 29), MAPK1, SRC, PIK3R1, VEGFA, and TP53, were considered as potential therapeutic targets, implying the underlying mechanisms of DZ acting on OPF. Conclusion. We investigated the possible therapeutic mechanisms of DZ from a systemic perspective. These key targets and pathways provided promising directions for the future research to reveal the exact regulating mechanisms of DZ in treating OPF.

1. Introduction

Osteoporosis is a systemic bone disorder characterized by low bone density, poor bone quality, reduced bone strength, and an accompanying increased incidence of fractures [1]. The latest study noted that over 10 million people suffer from osteoporosis in the United States. Meanwhile, the treatment burden of osteoporosis, $22 billion in 2008, is expected to rise due to the consistently increasing aging population. Related osteoporosis treatments should be used to alleviate or reduce symptoms such as fractures, which could result in considerable cost savings [2]. Osteoporotic fractures (OPF) have serious consequences, such as declining functions and independence as well as increased morbidity and mortality [3, 4]. Hip fractures are considered to be the most expensive and debilitating of fractures because approximately 10 to 20% of patients will end up with a disability in the year following a hip fracture and half of patients will lose their independence. Osteoporotic fractures are also associated with motor function decline when occurring at other skeletal sites [5]. To curb bone loss, maintain bone mass, and decrease the risk of OPF, antiresorptive agents, anabolic agents, and bone mineral drugs have been widely used in clinical treatment [6]. However, it is noteworthy that almost all these antiosteoporosis treatments have limitations and side effects, such as increasing the risk of breast cancer, jaw necrosis, or atypical femur fracture [7, 8]. Therefore, alternative and safe intervention strategies may be promising.

Recently, traditional Chinese medicine (TCM) has attracted worldwide attention and has served as a main alternative treatment in East Asia, North America, and Europe due to its satisfactory curative effect, relatively low toxicity, and low cost [9–11]. TCM has been used in China for thousands of years, and it exerts a beneficial therapeutic effect and reduces side effects through multiple herb combinations to prevent and treat various diseases [12]. Meanwhile, TCM has increasingly been found to be effective in the treatment of osteoporosis [13, 14].

Du Zhong (DZ), or Eucommiae Cortex, is one of the most commonly used TCM herbs in the treatment of osteoporosis. DZ has been found to facilitate osteogenesis through activating osteoblasts and to inhibit osteolysis by suppressing osteoclast activity [15]. An animal experiment confirmed that Du Zhong Wan, consisting of Eucommiae Cortex and Radix Dipsaci, is capable of improving trabecular bone mineral density and bone biomechanical properties and thus plays a key role in the treatment of osteoporosis [16]. Although DZ exerted therapeutic effects against osteoporosis, the therapeutic mechanism of its exact is still unclear. Similar to TCM, DZ exerts its therapeutic efficacy by regulating multiple molecules in the human body. Therefore, it is still a major challenge to prove the mechanism of DZ through scientific trials that are used in Western medicine. Fortunately, the development of systems pharmacology offers researchers an alternative opportunity and an option to investigate the pharmacological mechanisms of TCM [17]. Recently, Wang et al. utilized a network pharmacology method to clarify the synergistic mechanism of Er-Xian decoction in treating osteoporosis [18]. Similarly, Gan et al. employed this systems pharmacology method to dissect the mechanisms of the therapeutic effect of Rhizoma Drynariae on osteoporosis [19].

In this network pharmacology work, we aim to comprehensively dissect the mechanisms of DZ in treating OPF. We collected related compounds of DZ from multiple databases and obtained the compounds’ potential targets via target fishing. Then, we matched these targets with OPF-related targets that were collected from a multisource database. Next, using overlapping targets obtained from the previous process, we built a protein-protein network to analyze their internal interactions and then screened out the hub genes. Furthermore, we used the clusterProfiler package in R to conduct biological process and KEGG pathway enrichment analyses. The protocol of our experimental procedures is shown in Figure 1.

Figure 1: The schematic map of the present study to investigate potential mechanisms of DZ in the OPF treatment.

2. Materials and Methods

2.1. Identification of Chemical Ingredients

To collect the active compounds of DZ, we utilized the Traditional Chinese Medicine System Pharmacology Database (TCMSP™, http://lsp.nwu.edu.cn/tcmsp.php) [20], which is a frequently used platform for systems pharmacology. The BATMAN-TCM platform (http://bionet.ncpsb.org/batman-tcm/) [21], one of the largest comprehensive TCM platforms, uses the similarity of known TCM drug–target interactions to predict potential compound–target interactions. In total, ninety-three chemical compounds were collected in this part.

2.2. Compound Screening and Target Prediction

2.2.1. OB Evaluation

Oral bioavailability (OB) refers to the percentage of an orally administered drug that reaches systemic circulation, and it is one of the most important pharmacokinetic profiles for drug screening. The TCMSP platform has adopted the OBioavail1.1 system, which integrates P450, 3A4, and P-glycoprotein information to obtain the OB value [22]. To determine the active ingredients used for further steps, we set the OB threshold at 30%.

2.2.2. DL Prediction

Drug-likeness (DL) is a qualitative index that represents the degree to which the target compound is “drug-like” and is used to remove chemically unsuitable compounds. TCMSP uses the Tanimoto similarity method to calculate the DL index by comparing the target compound to all 6511 molecules in the DrugBank database [23]. In this process, the compounds that did not meet the condition that DL ≥ 0.18 were removed.

2.2.3. BATMAN-TCM

To collect candidate compounds, we set the potential TAR score cutoff ≥20 and -value cutoff <0.05.

2.2.4. Target Prediction of the Active Ingredients

The compounds have effects on the targets that induce them to exert their biological functions. Thus, we used the TCMSP and BATMAN-TCM platforms to predict the targets of active compounds.

2.3. Disease Target Identification by Multiple Databases

We used multiple databases to collect OPF-related targets, and the terms “Osteoporotic” and “Fracture” were used as the key words for the search. The databases in this step included the Comparative Toxicogenomics Database (CTD, http://ctdbase.org/) [24], GeneCards (https://www.genecards.org/) [25], OMIM (https://www.omim.org/) [26], HPO (https://hpo.jax.org/app/) [27], and GenCLiP3 (http://ci.smu.edu.cn/genclip3/analysis.php) [28]. In total, three thousand four hundred-thirteen OPF targets were found through these databases. Then, we constructed a Venn diagram to determine the overlap between OPF and the active compound targets because these overlapping targets may play a critical role in the treatment of OPF. These targets were analyzed by String (https://string-db.org/) [29], in which the minimum required interaction score was set at the highest confidence level (0.900), and then, the protein-protein interaction (PPI) data were exported.

2.4. Network Construction

2.4.1. Construction Method

Two main networks were built in this process, including the compound-hub gene network and the hub gene-pathway network. The target information was obtained from the KEGG pathway enrichment results. Cytoscape 3.6.2 (http://www.cytoscape.org/) was used for all the network visualization work. Cytoscape 3.6.2 is one of most powerful open-source software programs for constructing different networks visually [30]. However, there are several limitations faced by network pharmacology research. For example, questions including whether compounds activate or inhibit targets, how they exact their effects on targets, whether they use binding or catalysis, and so on remain unclear.

2.4.2. Topological Feature Definition

We used three parameters to describe and quantify the importance of nodes in these networks because nodes that bridge many edges with their neighbors are more likely to exert crucial mediating functions. (1) “Degree” means the number of edges shared with other nodes. Examining the node’s degree is the most straightforward method of quantifying its centrality [31]. (2) “Betweenness centrality” is used to measure node centrality based on the shortest paths. A high-value node would play a more important role in the network because more regulation information will pass through it [32]. (3) “Closeness centrality” represents the mean distance between the node and all the other nodes in the network and is the reciprocal of the sum of the length of the shortest paths between a node and all the other nodes in the network. Therefore, the central node would be more likely to be close with all other nodes [33].

2.5. Biological Process and Pathway Enrichment Analysis

We used the clusterProfiler package in R (R 3.6.2 for Windows) to perform Gene Ontology (GO) enrichment analysis and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis to identify biological processes and systemic involvement of pathways.

3. Results

3.1. Active Compounds of DZ

A total of 147 related compounds were found in DZ collected from the TCMSP database, and of these, 27 candidate compounds remained after the screening of ADME thresholds (OB ≥ 30%, DL ≥ 0.18). Using the BATMAN-TCM database, we obtained 73 active compounds that matched the filter criteria. In total, we collected 93 unique compounds.

3.2. Target Prediction and Analysis

For these unique compounds, we obtained 961 unique related targets, 104 from TCMSP and 857 from BATMAN-TCM. We integrated the OPF genes that were obtained from multisource databases, including the CTD, GeneCards, OMIM, HPO, GenCLiP 3 databases, and a total of 3,834 related genes were collected. After construction of the Venn diagram, three hundred twenty-five overlapping targets between the related targets of DZ and OPF were selected as the key targets on which DZ exerts its anti-OPF effects (Figure 2).

Figure 2: The Venn diagram for the targets of DZ and OPF. The overlap targets are the potential therapeutic genes for DZ when treating OPF.

The data obtained from the String database were used to establish the PPI network for the 325 overlapping targets. In this network, there were 249 nodes and 1006 edges in total. Then, three main parameters, namely, “degree”, “betweenness”, and “closeness”, were used as filters to select key genes and to build the large hub nodes to determine the anti-OPF effect of DZ. The first screening threshold was degree ≥ 6, closeness ≥ 0.312, and betweenness ≥ 0.000, which resulted in 122 nodes and 745 edges. Then, these 72 key nodes were further screened with a second threshold consisting of degree ≥ 11, closeness ≥ 0.340, and betweenness ≥ 0.016, and 64 nodes and 414 edges remained after this screening (Figure 3). These nodes included MAPK1, SRC, PIK3R1, VEGFA, TP53, EGFR, JUN, AGT, IL6, EGF, MAPK8, FOS, F2, FN1, EDN1, TNF, PPBP, ESR1, NCOA1, CTNNB1, NR3C1, ALB, RXRA, CBL, INS, PRKACA, AHSG, TGFB1, HSPA8, AR, FGF23, GNAS, NFKB1, GAS6, PRKCD, AVPR2, APOA1, PDGFB, HSP90AB1, IL4, HSPA1A, IGF2, SIRT1, PPARG, PPARA, CDC42, PGR, IFNG, CST3, NCOA3, TAC1, PLCB1, NOS3, CCL5, FGF2, ESR2, STAT5A, WNT5A, HTR2C, HTR2A, CCL2, CRH, APOE, and CFD (Table 1). After sorting these 64 big hub nodes in descending order, we found that MAPK1 (degree = 41), SRC (degree = 39), PIK3R1 (degree = 36), VEGFA (degree = 31), TP53 (degree = 29), EGFR (degree = 29), JUN (degree = 29), AGT (degree = 29), IL6 (degree = 28), EGF (degree = 27), MAPK8 (degree = 26), FOS (degree = 25), and F2 (degree = 25) were the most important targets in this PPI network.

Figure 3: The whole screening process for the PPI network through a topological method. In these images, the bigger size and brighter color represent higher DC value.

Table 1: Information of 64 hub targets.

UniProt ID
Gene symbol
Description
Degree

Q499G7
MAPK1
Mitogen-activated protein kinase 1
41

P12931
SRC
SRC proto-oncogene, non–receptor tyrosine kinase
39

P27986
PIK3R1
Phosphoinositide-3-kinase regulatory subunit 1
36

P15692
VEGFA
Vascular endothelial growth factor A
31

Q53GA5
TP53
Tumor protein p53
29

Q504U8
EGFR
Epidermal growth factor receptor
29

P05412
JUN
Jun proto-oncogene
29

P01019
AGT
Angiotensinogen
29

Q75MH2
IL6
Interleukin 6
28

P01133
EGF
Epidermal growth factor
27

P45983
MAPK8
Mitogen-activated protein kinase 8
26

Q6FG41
FOS
Fos proto-oncogene
25

P00734
F2
Coagulation factor II, thrombin
25

Q9UQS6
FN1
Fibronectin 1
24

Q6FH53
EDN1
Endothelin 1
24

Q5STB3
TNF
Tumor necrosis factor
24

P02775
PPBP
Pro-platelet basic protein
24

Q9UBT1
ESR1
Estrogen receptor 1
23

Q15788
NCOA1
Nuclear receptor coactivator 1
23

P35222
CTNNB1
Catenin beta 1
23

P04150
NR3C1
Nuclear receptor subfamily 3 group C member 1
23

P02768
ALB
Albumin
22

Q6P3U7
RXRA
Retinoid X receptor alpha
21

P22681
CBL
Cbl proto-oncogene
21

P01308
INS
Insulin
20

P17612
PRKACA
Protein kinase cAMP-activated catalytic subunit alpha
18

P02765
AHSG
Alpha 2-HS glycoprotein
18

P01137
TGFB1
Transforming growth factor beta 1
18

V9HW22
HSPA8
Heat shock protein family A (Hsp70) member 8
17

Q9NUA2
AR
Androgen receptor
17

Q9GZV9
FGF23
Fibroblast growth factor 23
17

Q5JWF2
GNAS
GNAS complex locus
16

P19838
NFKB1
Nuclear factor kappa B subunit 1
16

Q14393
GAS6
Growth arrest specific 6
15

Q05655
PRKCD
Protein kinase C delta
15

P30518
AVPR2
Arginine vasopressin receptor 2
15

P02647
APOA1
Apolipoprotein A1
15

P01127
PDGFB
Platelet-derived growth factor subunit B
15

Q6PK50
HSP90AB1
Heat shock protein 90 alpha family class B member 1
14

Q5FC01
IL4
Interleukin 4
14

P0DMV9
HSPA1A
Heat shock protein family A (Hsp70) member 1A
14

P01344
IGF2
Insulin-like growth factor 2
14

Q96EB6
SIRT1
Sirtuin 1
13

Q4W448
PPARG
Peroxisome proliferator activated receptor gamma
13

Q07869
PPARA
Peroxisome proliferator activated receptor alpha
13

P60953
CDC42
Cell division cycle 42
13

P06401
PGR
Progesterone receptor
13

P01579
IFNG
Interferon gamma
13

P01034
CST3
Cystatin C
13

Q9Y6Q9
NCOA3
Nuclear receptor coactivator 3
12

Q9Y494
TAC1
Tachykinin precursor 1
12

Q9NQ66
PLCB1
Phospholipase C beta 1
12

P29474
NOS3
Nitric oxide synthase 3
12

P13501
CCL5
C-C motif chemokine ligand 5
12

P09038
FGF2
Fibroblast growth factor 2
12

Q92731
ESR2
Estrogen receptor 2
11

Q59GY7
STAT5A
Signal transducer and activator of transcription 5A
11

P41221
WNT5A
Wnt family member 5A
11

P28335
HTR2C
5-Hydroxytryptamine receptor 2C
11

P28223
HTR2A
5-Hydroxytryptamine receptor 2A
11

P13500
CCL2
C-C motif chemokine ligand 2
11

P06850
CRH
Corticotropin releasing hormone
11

P02649
APOE
Apolipoprotein E
11

P00746
CFD
Complement factor D
11

Then, we further built the big hub nodes-compound network (Figure 4) based on these 64 key targets. This network included 40 compound nodes and 64 large hub target nodes. Subsequently, we reordered these compound nodes in descending order of degree and found that Civetone was relevant to 25 large hub genes, Epiquinidine to 19 genes, quercetin to 13 genes, Gamma-Sitosterol to 7 genes, and Pinosylvin to 7 genes.

Figure 4: The network for big hub genes-compounds connection. The diamond shape nodes are the hub genes and round ones represent compounds. And all nodes’ color changes according to their degree value.

3.3. GO Biological Process Enrichment Analysis

After sorting the 343 biological process (BP) terms in ascending order of adjusted value, we found that the top 25 biological processes (Figure 5) mainly concentrated on aspects of the reactive oxygen species metabolic process and its response to nutrient levels and epithelial cell proliferation. The following processes showed an enriched gene ratio: gland development (GO:0048732, 18.6%), response to nutrient levels (GO:0031667, 18.6%), response to a molecule of bacterial origin (GO:0002237, 15.5%), response to lipopolysaccharide (GO:0032496, 15.2%), response to an antibiotic (GO:0046677, 14.6%), reactive oxygen species metabolic process (GO:0072593, 13.6%), response to oxidative stress (GO:0006979, 16.4%), steroid metabolic process (GO:0008202, 14.2%), cellular response to drug (GO:0035690, 14.9%), epithelial cell proliferation (GO:0050673, 15.8%), response to a steroid hormone (GO:0048545, 14.9%), response to reactive oxygen species (GO:0000302, 11.8%), second-messenger-mediated signaling (GO:0019932, 15.2%), reproductive structure development (GO:0048608, 14.9%), calcium ion transport (GO:0006816, 14.9%), reproductive system development (GO:0061458, 14.9%), regulation of peptide secretion (GO:0002791, 15.8%), response to xenobiotic stimulus (GO:0009410, 12.4%), regulation of membrane potential (GO:0042391, 14.6%), cellular response to oxidative stress (GO:0034599, 12.4%), regulation of epithelial cell proliferation (GO:0050678, 13.6%), divalent metal ion transport (GO:0070838, 15.2%), divalent inorganic cation transport (GO:0072511, 15.2%), negative regulation of response to external stimulus (GO:0032102, 13.3%), and branching morphogenesis of an epithelial tube (GO:0048754, 9.3%). Based on these BP enrichment analyses, we found that the anti-OPF effect of DZ may result from its regulatory effect on the reactive oxygen species metabolic process, response to nutrient levels, and epithelial cell proliferation.

Figure 5: Top 25 processes of the biological process enrichment.

3.4. KEGG Pathway Enrichment Analysis

We conducted further KEGG pathway enrichment analyses of the 325 overlapping genes to determine the potential therapeutic mechanism of DZ for OPF. Then, we sorted 25 pathways based on the adjusted value (Figure 6). These pathways were the AGE-RAGE signaling pathway in diabetic complications (hsa04933, 9.6%), the estrogen signaling pathway (hsa04915, 11.1%), proteoglycans in cancer (hsa05205, 11.8%), fluid shear stress and atherosclerosis (hsa05418, 8.6%), parathyroid hormone synthesis, secretion, and action (hsa04928, 7.5%), breast cancer (hsa05224, 8.6%), Cushing’s syndrome (hsa04934, 8.6%), prostate cancer (hsa05215, 6.8%), the cAMP signaling pathway (hsa04024, 10%), morphine addiction (hsa05032, 6.4%), African trypanosomiasis (hsa05143, 4.3%), the GnRH signaling pathway (hsa04912, 6.4%), endocrine resistance (hsa01522, 6.4%), gastric cancer (hsa05226, 7.9%), neuroactive ligand-receptor interaction (hsa04080, 12.1%), the MAPK signaling pathway (hsa04010, 11.1%), the calcium signaling pathway (hsa04020, 8.6%), the oxytocin signaling pathway (hsa04921, 7.5%), malaria (hsa05144, 4.3%), Th17 cell differentiation (hsa04659, 6.1%), the HIF-1 signaling pathway (hsa04066, 6.1%), colorectal cancer (hsa05210, 5.4%), GnRH secretion (hsa04929, 4.6%), cortisol synthesis and secretion (hsa04927, 4.6%), and the TNF signaling pathway (hsa04668, 6.1%) (Table 2). We also simultaneously constructed the target-pathway network based on the DZ targets enriched in each pathway (Figure 7).

Figure 6: Top 25 pathways of the KEGG enrichment.

Table 2: Information of 25 pathways.

ID
Description
. adjust
Count

hsa04933
AGE-RAGE signaling pathway in diabetic complications
4.29750564649423e-15
27

hsa04915
Estrogen signaling pathway
4.29750564649423e-15
31

hsa05205
Proteoglycans in cancer
7.02978144058047e-12
33

hsa05418
Fluid shear stress and atherosclerosis
3.50540676797616e-09
24

hsa04928
Parathyroid hormone synthesis, secretion, and action
3.50540676797616e-09
21

hsa05224
Breast cancer
8.65129662473866e-09
24

hsa04934
Cushing syndrome
2.25061036469472e-08
24

hsa05215
Prostate cancer
2.25061036469472e-08
19

hsa04024
cAMP signaling pathway
4.25537215261066e-08
28

hsa05032
Morphine addiction
4.35290458036493e-08
18

hsa05143
African trypanosomiasis
5.25038311520481e-08
12

hsa04912
GnRH signaling pathway
5.25038311520481e-08
18

hsa01522
Endocrine resistance
1.16942061208876e-07
18

hsa05226
Gastric cancer
1.53049773606446e-07
22

hsa04080
Neuroactive ligand-receptor interaction
3.92997114777747e-07
34

hsa04010
MAPK signaling pathway
4.75554213381079e-07
31

hsa04020
Calcium signaling pathway
8.10557748077706e-07
24

hsa04921
Oxytocin signaling pathway
1.09530352609987e-06
21

hsa05144
Malaria
1.3497422679709e-06
12

hsa04659
Th17 cell differentiation
1.90497387620456e-06
17

hsa04066
HIF-1 signaling pathway
2.39525818923615e-06
17

hsa05210
Colorectal cancer
2.67545458381123e-06
15

hsa04929
GnRH secretion
2.67545458381123e-06
13

hsa04927
Cortisol synthesis and secretion
3.01425725456316e-06
13

hsa04668
TNF signaling pathway
3.01425725456316e-06
17

Figure 7: The targets-pathway network of DZ for treating OPF. The circle nodes represent big hub genes and the V-shape nodes represent the top 50 pathways. The nodes’ size and colors are dependent on DC value.

4. Discussion

Osteoporosis is becoming a serious health burden due to the aging population. Therefore, prevention of osteoporosis and of osteoporotic fracture has considerable social and economic significance [34]. Drug therapy, such as bisphosphonates and denosumab, is a practical way to treat osteoporosis by inhibiting osteoclast development, formation, and survival [35, 36]. Although these antiresorptive drugs are effective in treating osteoporosis, the long-term use or high dosage of these drugs can cause serious side effects, such as atypical long bone fractures and jaw osteonecrosis [37, 38]. In terms of fracture healing, the use of bisphosphonates remains controversial because callus remodeling is impaired in testing on ovariectomized rats [39]. Furthermore, these antiresorptive drugs cannot rebuild bone that is lost in the progression of osteoporosis, so there is a need to find a method capable of promoting bone building [40].

DZ has been found to alleviate improved osteoblast activity and reverse bone loss in ovariectomized (OVX) mice, which may result from the interaction between its active compounds and related targets. Thus, we investigated the anti-OPF mechanism of DZ through a network pharmacology method. Based on the big hub genes-compounds network (Figure 4), we found that five key compounds, including Civetone, Epiquinidine, quercetin, Gamma-Sitosterol, and Pinosylvin, play dominant roles in this network. Therefore, these active compounds may lay the foundation for the promising anti-OPF effect of DZ.

In the present study, 325 common gene targets of DZ and OPF were identified and 64 hub genes were identified that may play critical roles in this treatment. Some anti-OPF effects of these genes have been confirmed by clinical trials or animal experiments. Chow et al. found that low-magnitude high-frequency vibration stimulation is able to promote fracture healing by positively regulating the impaired innate immune response (lower expression of TNF-α and IL-6) and promoting macrophage polarization in delayed osteoporotic fracture healing [41]. Wang et al. suggested that organic gallium can improve osteoporotic fracture healing by positively modulating the OPG/RANKL ratio and inhibiting the expression of TNF-α and IL-6 [42]. Zhang et al. found that 17β-estradiol, which is used in the clinic to prevent and treat postmenopausal women with osteoporotic fracture, can significantly downregulate ROS production and the expression of TNF-α and IL-1β, as well as the phosphorylation of ERK2 (MAPK1) [43]. Wang et al. suggested that ESR1 is a susceptibility gene for osteoporotic fracture in postmenopausal Chinese women [44]. De-Ugarte et al. found that miR-320a, which is known to target CTNNB1 of osteoblasts, is overexpressed in osteoporotic fractures [45]. Swanberg et al. found that elderly women who were carriers of the IFNG variant alleles had lower BMD and a lower risk of incident fracture [46]. Based on a large population-based cohort study, Rivadeneira et al. found that ESR2 variants were related to an increased risk of fragility fracture in postmenopausal women [47]. These gene targets play key roles in osteoblast differentiation and osteoclastogenesis during fracture repair. Moreover, the results presented in Figure 4 show that one active compound can interact with multiple genes, one gene with dicompound can interact with multiple genes, and one gene can interact with different compounds, which is in accordance with the modern drug theory of “multi-ingredient, multitarget” [48].

We conducted GO biological process and KEGG enrichment analyses of the overlapping genes to identify their functions. We found 44 (13.6%) genes involved in reactive oxygen species metabolic processes. Reactive oxygen species play key roles in regulating cell proliferation, apoptosis, migration, and differentiation [49], as well as in the regulation of osteoclasts and osteoblasts [50, 51]. EGFR is the overlapping gene between reactive oxygen species (ROS) and the PI3K-Akt signaling pathway (27 genes enriched, 9.6%) based on the results of the GOBP and KEGG pathway enrichment analyses. The PI3K-Akt signaling pathway was first found in tumor cells and proved to regulate numerous cell functions, and it has been proven to be related to the promotion of osteogenesis by upregulating the proliferation and differentiation of bone marrow mesenchymal stem cells [52–54]. Furthermore, activation of the EGFR-Akt-mTOR pathway can effectively protect osteoblasts against dexamethasone, which increases the level of ROS production in osteoblastic cells [55]. MAPK1 is a common gene in reactive oxygen species, the PI3K-Akt signaling pathway, the MAPK signaling pathway (31 genes enriched, 11.1%), and the TNF signaling pathway. MAPKs play important roles in cell proliferation, differentiation, and apoptosis, as well as in regulating inflammation [56]. ROS can lead to serious periodontal tissue destruction by promoting osteoclastic bone resorption, which is closely related to the protein expression of MAPKs and NF-κB [57]. Srinivasa et al. found that the facilitation of ROS on osteoclast differentiation and resorption can be reversed by inhibiting the NF-κB and calcineurin-NFAT pathways [58]. The activation of MAPKs also plays a key role in osteoclastogenesis induced by RANKL. However, these ROS-directed upregulations of MAPK and NF-κB signaling can be attenuated by strengthening the nuclear factor-erythroid 2-related factor 2 [59]. Therefore, the anti-OPF effect of DZ may depend on the regulatory effect on ROS through these key genes and pathways.

Forty-five (13.9%) genes were involved in the regulation of the inflammatory response. In vitro and in vivo studies have shown that the proinflammatory cytokines IL-6 and TNF-α are involved in the pathogenesis of osteoporosis and influence bone metabolism [60, 61]. Several studies have also found that high expression of inflammatory markers is closely related to increased bone loss [62–64] and that these markers are also risk factors for incidents of nontraumatic fractures. Women with increased expression of inflammatory markers have a 3-fold risk of hip fractures [65]. IL-6 and TNF-α are the overlapping genes between the regulation of inflammatory response (17 genes enriched, 6.1%) and the TNF signaling pathway (45 genes enriched, 13.9%) based on the results of the BP and KEGG pathway enrichment analyses. Ali Aydin et al. found higher expression levels of TNF- α and IL-6 in OVX model rats, and these cytokines reached their apex when a fracture occurred [66]. Moreover, these cytokines also increase the production of oxygen radicals, which can enhance osteoclastic activity while decreasing osteoblastic activity [67, 68]. Therefore, the anti-OPF effect of DZ may depend on the regulation of the inflammatory response through these key genes and pathways.

5. Conclusion

By utilizing network pharmacology, we investigated the potential targets of DZ and the underlying mechanism of its anti-OPF effects, which may be based on the regulation of ROS and the inflammatory response. According to the KEGG pathway enrichment results, we found that the PI3K-Akt signaling pathway, MAPK signaling pathway, and TNF signaling pathway may be the main pathways in treating OPF. Thus, we believe that the anti-OPF effect of DZ is mainly based on its direct or indirect regulation of the abovementioned potential targets and pathways and that DZ provides promising directions for future research, which is essential to reveal its exact regulatory mechanisms.

Abbreviations

DZ: | Du Zhong
OPF: | Osteoporotic fractures
TCM: | Traditional Chinese medicine
TCMSP: | Traditional Chinese Medicine System Pharmacology Database
OB: | Oral bioavailability
DL: | Drug-likeness
CTD: | Comparative Toxicogenomics Database
OMIM: | Online Mendelian Inheritance in Man
HPO: | Human Phenotype Ontology
PPI: | Protein-protein interaction
KEGG: | Kyoto Encyclopedia of Genes and Genomes
GO: | Gene Ontology
BP: | Biological process
OVX: | Ovariectomized
ROS: | Reactive oxygen species.

Data Availability

All our main data used to support the findings of this study have been deposited in the Figshare repository (https://doi.org/10.6084/m9.figshare.12152607.v1).

Conflicts of Interest

The authors declare that they have no conflicts of interest regarding the publication of this paper.

Authors’ Contributions

Yongming Shuai and Fanhui Zeng conceived the idea of this article and supervised the research. Yongming Shuai performed the research, analyzed the data, and wrote the manuscript. Zhili Jiang and Qiuwen Yuan performed Target prediction and analysis as well as related enrichment processes. Shuqiang Tu and Fanhui Zeng participated in revising the data and improving manuscript writing. All authors reviewed the manuscript. And all authors read and approved the final version of the manuscript.

References

P. M. Camacho, S. M. Petak, N. Binkley et al., “American association of clinical endocrinologists and American college of endocrinology clinical practice guidelines for the diagnosis and treatment of postmenopausal osteoporosis–2016,” Endocrine Practice, vol. 22, no. 4, pp. 1–42, 2016.

E. M. Lewiecki, J. D. Ortendahl, J. Vanderpuye-Orgle et al., “Healthcare policy changes in osteoporosis can improve outcomes and reduce costs in the United States.,” JBMR Plus, vol. 3, Article ID e10192, 2019.

J. A. Cauley, “Public health impact of osteoporosis,” Journals of Gerontology Series A: Biomedical Sciences and Medical Sciences, vol. 68, no. 10, pp. 1243–1251, 2013.

M. A. Clynes, N. C. Harvey, E. M. Curtis, N. R. Fuggle, E. M. Dennison, and C. Cooper, “The epidemiology of osteoporosis.,” British Medical Bulletin, vol. 133, pp. 105–117, 2020.

E. M. Lewiecki, B. Chastek, K. Sundquist et al., “Osteoporotic fracture trends in a population of US managed care enrollees from 2007 to 2017,” Osteoporosis International, vol. 31, no. 7, pp. 1299–1304, 2020.

J.-G. Letarouilly, O. Broux, and A. Clabaut, “New insights into the epigenetics of osteoporosis,” Genomics, vol. 111, no. 4, pp. 793–798, 2019.

M. Kharazmi, J. Schilcher, P. Hallberg, and K. Michaëlsson, “Bisphosphonate-associated atypical fractures of the femur: an update of the current evidence,” Lakartidningen, vol. 116, 2019.

O. Tabatabaei-Malazy, P. Salari, P. Khashayar, and B. Larijani, “New horizons in treatment of osteoporosis,” Daru Journal of Pharmaceutical Sciences, vol. 25, no. 1, 2017.

L.-J. Xiao and R. Tao, “Traditional Chinese medicine (TCM) therapy,” Advances in Experimental Medicine and Biology, vol. 1010, pp. 261–280, 2017.

W. Xu, A. D. Towers, P. Li, and J.-P. Collet, “Traditional Chinese medicine in cancer care: Perspectives and experiences of patients and professionals in China,” European Journal of Cancer Care, vol. 15, no. 4, pp. 397–403, 2006.

T.-T. Sham, C.-O. Chan, Y.-H. Wang, J.-M. Yang, D. Kam-Wah Mok, and S.-W. Chan, “A review on the traditional Chinese medicinal herbs and formulae with hypolipidemic effect,” BioMed Research International, vol. 2014, Article ID 925302, pp. 1–21, 2014.

P. Li, L.-p. Yang, and Y.-w. Gong, “Application of systems biology technology in research of traditional Chinese medicine,” Journal of Traditional Chinese Medicine, vol. 29, no. 2, pp. 153–157, 2009.

L. P. Zhou, K. Y. Wong, H. T. Yeung et al., “Bone protective effects of danggui buxue tang alone and in combination with tamoxifen or raloxifene in vivo and in vitro.,” Frontiers in Pharmacology, vol. 9, p. 779, 2018.

F. Gu, J. Jiang, S. Wang et al., “An experimental research into the potential therapeutic effects of anti-osteoporosis decoction and Yougui pill on ovariectomy-induced osteoporosis,” American Journal of Translational Research, vol. 11, no. 9, pp. 6032–6039, 2019.

H. Ha, J. Ho, S. Shin et al., “Effects of Eucommiae Cortex on osteoblast-like cell proliferation and osteoclast inhibition,” Archives of Pharmacal Research, vol. 26, no. 11, pp. 929–936, 2003.

F. Li, X. Yang, J. Bi, Z. Yang, and C. Zhang, “Antiosteoporotic activity of Du-Zhong-Wan water extract in ovariectomized rats,” Pharmaceutical Biology, vol. 54, no. 9, pp. 1857–1864, 2016.

W. Zhang, Y. Huai, Z. Miao, A. Qian, and Y. Wang, “Systems pharmacology for investigation of the mechanisms of action of traditional chinese medicine in drug discovery,” Frontiers in Pharmacology, vol. 10, p. 743, 2019.

N. Wang, P. Xu, X. Wang et al., “Integrated pathological cell fishing and network pharmacology approach to investigate main active components of Er-Xian decotion for treating osteoporosis,” Journal of Ethnopharmacology, vol. 241, Article ID 111977, 2019.

D. Gan, X. Xu, D. Chen, P. Feng, and Z. Xu, “Network pharmacology-based pharmacological mechanism of the Chinese medicine rhizoma drynariae against osteoporosis,” Medical Science Monitor, vol. 25, pp. 5700–5716, 2019.

J. Ru, P. Li, J. Wang et al., “TCMSP: A database of systems pharmacology for drug discovery from herbal medicines.,” Journal of Cheminformatics, vol. 6, no. 1, p. 13, 2014.

Z. Liu, F. Guo, Y. Wang et al., “BATMAN-TCM: A bioinformatics analysis tool for molecular mechANism of traditional Chinese medicine,” Scientific Reports, vol. 6, no. 1, Article ID 21146, 2016.

X. Xu, W. Zhang, C. Huang et al., “A novel chemometric method for the prediction of human oral bioavailability,” International Journal of Molecular Sciences, vol. 13, no. 6, pp. 6964–6982, 2012.

W. Tao, X. Xu, X. Wang et al., “Network pharmacology-based prediction of the active ingredients and potential targets of Chinese herbal Radix Curcumae formula for application to cardiovascular disease,” Journal of Ethnopharmacology, vol. 145, no. 1, pp. 1–10, 2013.

C. J. Grondin, A. P. Davis, T. C. Wiegers, J. A. Wiegers, and C. J. Mattingly, “Accessing an expanded exposure science module at the comparative Toxicogenomics database.,” Environmental Health Perspectives, vol. 126, no. 1, Article ID 014501, 2018.

M. Rebhan, V. Chalifa-Caspi, J. Prilusky, and D. Lancet, “GeneCards: A novel functional genomics compendium with automated data mining and query reformulation support,” Bioinformatics, vol. 14, no. 8, pp. 656–664, 1998.

J. S. Amberger and A. Hamosh, “Searching online mendelian inheritance in man (OMIM): A knowledgebase of human genes and genetic phenotypes.,” Current Protocols in Bioinformatics, vol. 58, no. 1, 2017.

P. N. Robinson, S. Köhler, S. Bauer, D. Seelow, D. Horn, and S. Mundlos, “The Human phenotype ontology: A tool for annotating and analyzing human hereditary disease,” The American Journal of Human Genetics, vol. 83, no. 5, pp. 610–615, 2008.

J. H. Wang, L. F. Zhao, H. F. Wang et al., “GenCLiP 3: Mining human genes’ functions and regulatory networks from pubMed based on co-occurrences and natural language processing,” Bioinformatics, vol. 36, no. 6, pp. 1973–1975, 2019.

D. Szklarczyk, J. H. Morris, H. Cook et al., “The STRING database in 2017: quality-controlled protein-protein association networks, made broadly accessible,” Nucleic Acids Research, vol. 45, no. D1, pp. D362–D368, 2017.

P. Shannon, A. Markiel, and O. Ozier, “Cytoscape: A software environment for integrated models of biomolecular interaction networks,” Genome Research, vol. 13, no. 11, pp. 2498–2504, 2003.

P. V. Missiuro, K. Liu, L. Zou et al., “Information flow analysis of interactome networks,” PLOS Computational Biology, vol. 5, no. 4, Article ID e1000350, 2009.

K. Raman, N. Damaraju, and G. K. Joshi, “The organisational structure of protein networks: Revisiting the centrality-lethality hypothesis,” Systems and Synthetic Biology, vol. 8, no. 1, pp. 73–81, 2014.

Y. Tang, M. Li, J. Wang, Y. Pan, and F.-X. Wu, “CytoNCA: A cytoscape plugin for centrality analysis and evaluation of protein interaction networks,” Biosystems, vol. 127, pp. 67–72, 2015.

P. Mitchell, K. Åkesson, M. Chandran, C. Cooper, K. Ganda, and M. Schneider, “Implementation of models of care for secondary osteoporotic fracture prevention and orthogeriatric models of care for osteoporotic hip fracture,” Best Practice & Research Clinical Rheumatology, vol. 30, no. 3, pp. 536–558, 2016.

R. G. G. Russell, N. B. Watts, F. H. Ebetino, and M. J. Rogers, “Mechanisms of action of bisphosphonates: Similarities and differences and their potential influence on clinical efficacy,” Osteoporosis International, vol. 19, no. 6, pp. 733–759, 2008.

S. Fukumoto and T. Matsumoto, “Recent advances in the management of osteoporosis,” F1000Research, vol. 6, p. 625, 2017.

M. Rossini, G. Adami, S. Adami, O. Viapiana, and D. Gatti, “Safety issues and adverse reactions with osteoporosis management,” Expert Opinion on Drug Safety, vol. 15, no. 3, pp. 321–332, 2016.

M. Kharazmi, P. Hallberg, G. Warfvinge, and K. Michaëlsson, “Risk of a typical femoral fractures and osteonecrosis of the jaw associated with alendronate use compared with other oral bisphosphonates,” Rheumatology, vol. 53, no. 10, pp. 1911–1913, 2014.

T. Savaridas, R. J. Wallace, D. M. Salter, and A. H. R. W. Simpson, “Do bisphosphonates inhibit direct fracture healing?” The Bone & Joint Journal, vol. 95-B, no. 9, pp. 1263–1268, 2013.

J. Wang, K. Su, W. Sang, L. Li, and S. Ma, “Thiazide diuretics and the incidence of osteoporotic fracture: A systematic review and meta-analysis of cohort studies.,” Frontiers in Pharmacology, vol. 10, p. 1364, 2019.

S.-H. Chow, Y.-N. Chim, J. Wang et al., “Vibration treatment modulates macrophage polarisation and enhances early inflammatory response in oestrogen-deficient osteoporotic-fracture healing,” European Cells and Materials, vol. 38, pp. 228–245, 2019.

J. Wang, M. He, G. Wang, and Q. Fu, “Organic gallium treatment improves osteoporotic fracture healing through affecting the OPG/RANKL ratio and expression of serum inflammatory cytokines in ovariectomized rats,” Biological Trace Element Research, vol. 183, no. 2, pp. 270–279, 2018.

Y. Zhang, Y. He, Y. Zong et al., “17 -estradiol attenuates homocysteine-induced oxidative stress and inflammatory response as well as MAPKs cascade via activating PI3-K/Akt signal transduction pathway in Raw 264.7 cells,” Acta Biochimica et Biophysica Sinica, vol. 47, no. 2, pp. 65–72, 2015.

C. Wang, Z. Zhang, H. Zhang et al., “Susceptibility genes for osteoporotic fracture in postmenopausal Chinese women,” Journal of Bone and Mineral Research, vol. 27, no. 12, pp. 2582–2591, 2012.

L. De-Ugarte, G. Yoskovitz, S. Balcells et al., “MiRNA profiling of whole trabecular bone: Identification of osteoporosis-related changes in MiRNAs in human hip bones.,” BMC Medical Genomics, vol. 8, no. 1, p. 75, 2015.

M. Swanberg, F. E. McGuigan, K. K. Ivaska, P. Gerdhem, and K. Åkesson, “Polymorphisms in the inflammatory genes CIITA, CLEC16A and IFNG influence BMD, bone loss and fracture in elderly women,” PLoS One, vol. 7, no. 10, Article ID e47964, 2012.

F. Rivadeneira, J. B. J. van Meurs, J. Kant et al., “Estrogen receptor β (ESR2) polymorphisms in interaction with estrogen receptor α (ESR1) and insulin-like growth factor I (IGF1) variants influence the risk of fracture in postmenopausal women,” Journal of Bone and Mineral Research, vol. 21, no. 9, pp. 1443–1456, 2006.

S. Sheng, J. Wang, L. Wang et al., “Network pharmacology analyses of the antithrombotic pharmacological mechanism of Fufang Xueshuantong capsule with experimental support using disseminated intravascular coagulation rats,” Journal of Ethnopharmacology, vol. 154, no. 3, pp. 735–744, 2014.

V. Domazetovic, G. Marcucci, T. Iantomasi, M. L. Brandi, and M. T. Vincenzini, “Oxidative stress in bone remodeling: Role of antioxidants,” Clinical Cases in Mineral and Bone Metabolism, vol. 14, no. 2, pp. 209–216, 2017.

K. Schröder, “NADPH oxidases in bone homeostasis and osteoporosis,” Free Radical Biology and Medicine, vol. 132, pp. 67–72, 2019.

X. Gan, S. Huang, Q. Yu, H. Yu, and S. S. Yan, “Blockade of Drp1 rescues oxidative stress-induced osteoblast dysfunction,” Biochemical and Biophysical Research Communications, vol. 468, no. 4, pp. 719–725, 2015.

G. Luo, B. Xu, and Y. Huang, “Icariside II promotes the osteogenic differentiation of canine bone marrow mesenchymal stem cells via the PI3K/AKT/mTOR/S6K1 signaling pathways,” American Journal of Translational Research, vol. 9, pp. 2077–2087, 2017.

G. Li, J. Hu, H. Chen et al., “Enamel matrix derivative enhances the proliferation and osteogenic differentiation of human periodontal ligament stem cells on the titanium implant surface,” Organogenesis, vol. 13, no. 3, pp. 103–113, 2017.

Y. Tong, W. Feng, Y. Wu et al., “Mechano-growth factor accelerates the proliferation and osteogenic differentiation of rabbit mesenchymal stem cells through the PI3K/AKT pathway,” BMC Biochemistry, vol. 16, no. 1, p. 1, 2015.

J.-b. Fan, W. Liu, X.-h. Zhu et al., “EGFR-AKT-mTOR activation mediates epiregulin-induced pleiotropic functions in cultured osteoblasts,” Molecular and Cellular Biochemistry, vol. 398, no. 1-2, pp. 105–113, 2015.

W. V. Berghe, S. Plaisance, E. Boone et al., “p38 and extracellular signal-regulated kinase mitogen-activated protein kinase pathways are required for nuclear factor-kappaB p65 transactivation mediated by tumor necrosis factor.,” Journal of Biological Chemistry, vol. 273, no. 6, pp. 3285–3290, 1998.

Y. An, H. Zhang, C. Wang et al., “Activation of ROS/MAPK s/NF-κB/NLRP3 and inhibition of efferocytosis in osteoclast-mediated diabetic osteoporosis,” The FASEB Journal, vol. 33, no. 11, pp. 12515–12527, 2019.

L. Jiating, J. Buyun, and Z. Yinchang, “Role of metformin on osteoblast differentiation in type 2 diabetes,” Biomed Res Int, vol. 2019, Article ID 9203934, 6 pages, 2019.

D. Thummuri, V. G. M. Naidu, and P. Chaudhari, “Carnosic acid attenuates RANKL-induced oxidative stress and osteoclastogenesis via induction of Nrf2 and suppression of NF-κB and MAPK signalling,” Journal of Molecular Medicine, vol. 95, no. 10, pp. 1065–1076, 2017.

R. B. Kimble, A. B. Matayoshi, J. L. Vannice, V. T. Kung, C. Williams, and R. Pacifici, “Simultaneous block of interleukin-1 and tumor necrosis factor is required to completely prevent bone loss in the early postovariectomy period,” Endocrinology, vol. 136, no. 7, pp. 3054–3061, 1995.

M. N. Weitzmann, S. Cenci, L. Rifas, J. Haug, J. Dipersio, and R. Pacifici, “T cell activation induces human osteoclast formation via receptor activator of nuclear factor κB ligand-dependent and -independent mechanisms,” Journal of Bone and Mineral Research, vol. 16, no. 2, pp. 328–337, 2001.

L. Mosekilde, H. Beck-Nielsen, O. H. Sørensen et al., “Hormonal replacement therapy reduces forearm fracture incidence in recent postmenopausal women - results of the Danish Osteoporosis Prevention Study,” Maturitas, vol. 36, no. 3, pp. 181–193, 2000.

C. Ding, V. Parameswaran, R. Udayan, J. Burgess, and G. Jones, “Circulating levels of inflammatory markers predict change in bone mineral density and resorption in older adults: A longitudinal study.,” The Journal of Clinical Endocrinology & Metabolism, vol. 93, no. 5, pp. 1952–1958, 2008.

E. R. Gertz, N. E. Silverman, K. S. Wise et al., “Contribution of serum inflammatory markers to changes in bone mineral content and density in postmenopausal women: A 1-year investigation,” Journal of Clinical Densitometry, vol. 13, no. 3, pp. 277–282, 2010.

K. E. Barbour, L.-Y. Lui, K. E. Ensrud et al., “Inflammatory markers and risk of hip fracture in older white women: The study of osteoporotic fractures,” Journal of Bone and Mineral Research, vol. 29, no. 9, pp. 2057–2064, 2014.

A. Aydin, Z. Halici, A. Albayrak et al., “Treatment with carnitine enhances bone fracture healing under osteoporotic and/or inflammatory conditions,” Basic & Clinical Pharmacology & Toxicology, vol. 117, no. 3, pp. 173–179, Feb. 2015.

N. Mody, F. Parhami, T. A. Sarafian, and L. L. Demer, “Oxidative stress modulates osteoblastic differentiation of vascular and bone cells,” Free Radical Biology and Medicine, vol. 31, no. 4, pp. 509–519, 2001.

B. Watkins, H. E. Lippman, B. L. Le, Y. Li, and M. F. Seifert, “Bioactive fatty acids: Role in bone biology and bone cell function,” Progress in Lipid Research, vol. 40, no. 1-2, pp. 125–148, 2001.