Systematic Pharmacological Methodology to Explore the Pharmacological Mechanism of Siwu Decoction for Osteoporosis.

Osteoporosis is an important health problem worldwide. Siwu decoction (SWD) and its modification have a good clinical effect on osteoporosis. However, the molecular mechanism of SWD on osteoporosis has not been thoroughly explained. A systematic pharmacological methodology was utilized to predict the active compounds and potential targets of SWD, collect the genes of osteoporosis and the known targets of SWD, and analyze the osteoporosis and SWD's network. Five networks were constructed and analyzed: (1) Osteoporosis genes' protein-protein interaction (PPI) network; (2) Compound-compound target network of SWD; (3) SWD-osteoporosis PPI network; (4) Compound-known target network of SWD; and (5) SWD known target- osteoporosis PPI network. Several osteoporosis and treatment-related targets (eg.,. HSP90AB1, FGFR1, HRAS, GRB2, and PGF), clusters, biological processes, and signaling pathways (e.g., PI3K-Akt signaling pathway, insulin signaling pathway, MAPK signaling pathway and FoxO signaling pathway) were found. The therapeutic effect of SWD on osteoporosis may be achieved by interfering with the biological processes and signaling pathways related to the development of osteoporosis.

Background

Osteoporosis is an important health problem worldwide, the incidence of which is closely related to the increase of the elderly population and the decline of quality of life [1,2]. This systemic bone disease is mainly characterized by low bone mass and bone microstructural damage, which leads to increased bone fragility and fracture [1,2]. The primary goal of osteoporosis treatment is to prevent fractures, maintain or increase bone mineral density (BMD), and improve bone physiology [3,4]. Fracture prevention is a global public health priority [4], and current methods include non-drug management and drug intervention [4]. Non-drug management methods include nutritional improvement, calcium, vitamin D supplements, and exercise [4,5]. Drug interventions mainly include anti-bone resorption drugs (such as bisphosphonates) and osteogenic drugs (e.g., parathyroid hormone and other estrogen hormone therapy [raloxifene]) and other bone metabolism regulation drugs, such as teriparatide, zoledronic acid, ibandronate, and strontium ranelate [4,6]. However, long-term use of these drugs often causes adverse events; for example, anti-bone resorption drugs have liver and kidney toxicity, and it can also lead to gastrointestinal complications, musculoskeletal pain, and osteonecrosis of the jaw [7]. The adverse effects of estrogen hormone therapy are mainly endometrial hyperplasia and uterine bleeding [8]. The adverse effects of calcitonin include flushing, nausea, vomiting, dizziness, and mild facial redness, and continued application may diminish its effect; bone density will gradually decrease after stopping the drug [9].

Recent research shows that Siwu decoction (SWD) and its modification have a good clinical effect on senile osteoporosis and bone fractures [10-12]. SWD is an oral decoction, that is, a liquid dosage form prepared by soaking and boiling its 4 herbs in water and then removing the juice from the residue. Studies also show that SWD has estrogen-like effects in both in vitro and in vivo [13-15]. SWD extracts can stimulate bone formation through the PI3K/Akt/NF-κB signaling pathway in osteoblasts [16]. These facts suggest that SWD may be an alternative treatment for osteoporosis. SWD is composed of 4 herbs: Rehmanniae Radix Praeparata (Shu Di Huang), Angelicae Sinensis Radix (Dang Gui), Paeoniae Radix Alba (Bai Shao), and Chuanxiong Rhizoma (Chuan Xiong). At present, research on the effect of SWD on osteoporosis is mostly applied by traditional pharmacological methods, which is limited by the perspective of “single compound–single target–single pathway”. Therefore, this study used a systematic pharmacology approach to explore the pharmacological mechanism of SWD in treatment of osteoporosis.

Material and Methods

Data preparation

Acquisition of SWD’s compounds

The traditional Chinese Medicine (TCM) Database@Taiwan [17] (, updated in March 2014) and the Traditional Chinese Medicine Systems Pharmacology Database [18] (TCMSP™, , updated on May 31, 2014) were employed to collect the compounds of SWD.

Pharmacokinetic prediction

To predict the bioactive compounds of SWD, oral bioavailability (OB), Caco-2 permeability, and drug-likeness (DL) were applied [19-21]. The standard was OB ≥30%, Caco-2 >−0.4 and DL ≥0.18. Any compounds meeting this standard were selected for subsequent research and all others were excluded. The details of OB, Coca-4, and DL, and the methodology of acquisition were described in our previous work [19-21].

According to these indexes, we included the following compounds: beta-sitosterol, stigmasterol, sitosterol, 11alpha,12alpha-epoxy-3beta-23-dihydroxy-30-norolean-20-en-28,12beta-olide (MOL001910), paeoniflorgenone, paeoniflorin, mairin, kaempferol, (+)-catechin, (3S,5R,8R,9R,10S,14S)-3,17-dihydroxy-4,4,8,10,14-pentamethyl-2,3,5,6,7,9-hexahydro-1H-cyclopenta[a]phenanthrene-15,16-dione (Palbinone), perlolyrine, senkyunone, wallichilide, mandenol, and myricanone.

Exceptional molecules

Since the application of biological models to predict SWD compounds has limitations [22], in order to avoid missing active compounds during the pre-screening process, we searched a large number of references and selected oral absorbable compounds with pharmacological activity.

Combined with relevant references [23-29], the following compounds were included: stachyose, senkyunolide I, rehmannioside D, rehmannioside A, paeoniflorin sulfonate, ligustilide, ferulic acid, catalpol, butylidenephthalide, albiflorin, acteoside, ligustrazine, paeonol, benzoic acid, and 5-hydroxymethylfurfural.

Compound targets and known targets for SWD

To assess the compound targets and known targets for SWD, we collected the molecular structure of each compound by searching the SciFinder (), drew them in ChemBioDraw, and saved them as “mol2” file format. We input the “mol2” files of each compound into PharmMapper () to predict their potential targets [30]. The known targets were collected from the TCMSP [18].

Protein name correction

The UniProtKB (), a database containing the accurate annotation of proteins and other substances, was used for the correction of protein’s names and the collection of official symbols with the species limited to “Homo sapiens”. The corrected target is detailed in the supplementary materials (Supplementary Tables 1, 2).

Osteoporosis genes

OMIM database and Genecards were utilized to collect the osteoporosis-related disease genes and targets [19-21]. OMIM database () is the database that catalogs all known diseases with a genetic component [31]. Genecards () is “a database about genes, their products and biomedical applications maintained by Israel’s Weizmann Institute of Science” [19-21]. Relevance scores greater than 5 were included in the network construction and analysis. The osteoporosis-related genes and their relevance scores are shown in Supplementary Table 3.

Protein–protein interaction data

The String database () and the IntAct database () were utilized to obtain data on protein–protein interaction (PPI). While searching the String database, the species was limited to “Homo sapiens” [19-21,32-33].

Network construction

Network construction method

Cytoscape 3.4.0 software was used for the network visualization and network analysis [34] (). Several networks were constructed: (1) Osteoporosis genes’ PPI network; (2) Compound-compound target network of SWD; (3) SWD-osteoporosis PPI network; (4) Compound-known target network of SWD; and (5) SWD known target- osteoporosis PPI network.

Cluster

The definition and the methodology of acquisition of clusters were described in our previous work “Exploring the Pharmacological Mechanism of Danzhi Xiaoyao Powder on ER-Positive Breast Cancer by a Network Pharmacology Approach” [19].

Enrichment analysis

The DAVID database (, ver. 6.8) was applied for Gene Ontology (GO) enrichment analysis and pathway enrichment analysis [35].

Results

Osteoporosis network analysis

Osteoporosis network

We input 357 genes into the String database to collect the PPI information so as to build this network. It contains 357 nodes and 33770 edges (Figure 1).

Figure 1

Osteoporosis network.

In this network, several genes are thought to be the key gene, as they have higher degrees: ALB (143 edges), INS (138 edges), TP53 (118 edges), IL6 (112 edges), IGF1 (108 edges), TNF (101 edges), EGF (101 edges), and JUN (101 edges). These may be the key or central genes in osteoporosis.

Clusters of osteoporosis network

After analyzing the osteoporosis network by MCODE, 7 clusters were returned (Table 1, Figure 2), and these were subjected to GO enrichment analysis. Several biological processes were obtained. Cluster 16 and 18 did not return any human biological processes. Cluster 2, 4, 5, 11, 13, 14, 15, and 17 did not return any osteoporosis-related biological processes. Taking cluster 1 as an example: cluster 1 gets GO: 0001503, GO: 0032355, GO: 0032755, GO: 0032869, GO: 0033280, GO: 0035630, GO: 0045669, GO: 0045672, GO: 0046427, GO: 0098868, and GO: 2000366. The details of each clusters and biological processes are shown in Supplementary Table 4.

Table 1

Cluster of osteoporosis network.

Cluster	Score	Nodes	Edges	Genes
1	32.326	44	695	IFNG, IL2, INS, MMP1, TF, SRC, EGF, CTNNB1, CCND1, TNF, TP53, MMP2, CD40LG, ESR1, IL10, JUN, FOS, ALB, LEPR, PPARG, LEP, IGF1, CRP, IGF1R, KIT, POMC, IL1B, NR3C1, IL1A, LDLR, APOE, IL6, IGF2, CYP19A1, WT1, AR, CSF1, STAT1, TGFB1, RUNX2, MAPK3, SP1, SPP1, MMP3
2	7.714	8	27	NELF, PROK2, KAL1, PROKR2, KISS1R, TACR3, GNRHR, TAC3
3	5.7	21	57	CRH, PTH, MMP14, IL7, GGT1, VDR, FGFR1, CD79A, ESR2, TNFSF11, TERT, IL11, BMP2, PRL, BGLAP, SHBG, BMP4, TNFRSF11B, ELN, SMAD3, GNRH1
4	5.667	7	17	LEPRE1, PLOD1, COL7A1, FKBP10, CRTAP, PPIB, SERPINH1
5	5	5	10	B3GAT3, XYLT2, HSPG2, B3GALT6, B4GALT7
6	4.667	16	35	RFC2, GTF2IRD1, CYP27A1, LIMK1, CLIP2, WRAP53, TBL2, CYP27B1, DKC1, NHP2, NOP10, GC, CYP24A1, NR1I2, TINF2, BAZ1B
7	4.333	13	26	OTX2, IGFBP4, SOST, LRP5, FZD4, NOTCH3, NOG, NOTCH2, CDX2, FGF17, SP7, WNT3A, LRP6
8	4	5	8	GHR, ADIPOQ, IGFBP3, IGFBP1, GH1
9	4	4	6	SFRP1, FGF8, WNT1, DKK1
10	3.8	11	19	PTHLH, STAR, IBSP, CTSK, CYP11A1, ACP5, FSHR, NFATC1, TRAF6, COL1A1, CALCA
11	3.714	8	13	IFT122, ERCC2, IFT52, ERCC6, WDR19, WDR35, PRKACB, CDC73
12	3.091	12	17	IAPP, CYP17A1, GNAS, SOX9, NR0B1, PRKACA, NR5A1, CGA, ADCY10, CHD7, CALCR, PTH1R
13	3	3	3	CYP21A2, HSD3B2, CYP11B1
14	3	3	3	POU1F1, PROP1, HESX1
15	3	3	3	UGT2B17, AKR1D1, HSD11B1
16	3	3	3	EFEMP2, GORAB, FBLN5
17	3	3	3	BMP1, LAMA3, LAMC2
18	3	3	3	MKRN3, MAGEL2, SNRPN

Figure 2

Cluster of osteoporosis network (A, B, C, … stand for cluster 1, 2, 3, …).

Pathway of osteoporosis network

After the pathway enrichment analysis, 17 osteoporosis-related pathways were obtained (Figure 3). The PI3K-Akt signaling pathway had the largest number of genes in this analysis (30 genes); the MAPK signaling pathway includes 24 genes; Osteoclast differentiation has 22 genes; Jak-STAT signaling pathway contains 17 genes; Wnt signaling pathway has 16 genes, and so on. The details of each pathway are shown in Supplementary Table 5.

Figure 3

Pathway of osteoporosis network (blue circle stands for osteoporosis gene. Red diamond stands for pathway).

These pathways with so many osteoporosis genes may be the key pathways involved in the development of osteoporosis. Early intervention in these signaling pathways may indicate a strategy for future treatment and prevention of osteoporosis.

Compound-compound target network analysis

This network is composed of 391 nodes (361 compound target nodes and 30 compound nodes) and 3961 edges. In this network, targets near the center are regulated by more compounds than are peripheral targets. For example, GSTA1, F2, CA2, CDK2, CCNA2, PDPK1, BACE1, and GSTP1 can be controlled by all of the compounds; while NDST1, UMPS, GSS, CBS, TPI1, and CSK are regulated by only 1 compound (Figure 4).

Figure 4

Compound-compound target network of SWD. Pink hexagons stand for compound targets; Blue, orange, yellow and green circles stand for compounds of Rehmanniae Radix Praeparata, Angelicae Sinensis Radix, Chuanxiong Rhizoma, Paeoniae Radix Alba, resp. Blue triangle stands for common compound of Rehmanniae Radix Praeparata, Angelicae Sinensis Radix, Paeoniae Radix Alba. Yellow triangle stands for common compound of Chuanxiong Rhizoma, Paeoniae Radix Alba. Green triangle stands for common compound of Rehmanniae Radix Praeparata, Angelicae Sinensis Radix. Orange triangles stand for common compounds of Angelicae Sinensis Radix, Chuanxiong Rhizoma.

SWD-osteoporosis network analysis

SWD-osteoporosis network

To integrate the osteoporosis network and SWD network, the SWD-Osteoporosis network was formed. This network contains 671 nodes and 10698 edges. There are 317 nodes that are compound targets, 31 nodes are SWD-osteoporosis targets, and 323 nodes are osteoporosis genes (Figure 5).

Figure 5

SWD-osteoporosis network (pink, blue, purple circle stand for compound target, osteoporosis genes, compound-osteoporosis target, resp.).

Clusters of SWD-osteoporosis metwork

Twenty-one clusters were collected after being analyzed by MCODE. There were some osteoporosis-related genes in clusters (Blue circles), and these genes may be the key genes of SWD that are active in treating osteoporosis (Table 2, Figure 6). These clusters were subjected to GO enrichment analysis. Cluster 5, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20 did not return osteoporosis-related biological processes. As an example, in cluster 1 there were many biological processes, including: GO: 0000165, GO: 0000186, GO: 0007179, GO: 0007183, GO: 0032355, GO: 0032755, GO: 0035630, GO: 0038128, GO: 0045669, GO: 0060397, GO: 0060557, and GO: 2000366. The details of each cluster and biological process are shown in Supplementary Table 6.

Table 2

Cluster of SWD-osteoporosis network.

Cluster	Score	Nodes	Edges	Targets and genes
1	47.525	60	1402	MMP2, CTNNB1, IL2, CCND1, ANXA5, GGT1, GRB2, ABL1, PARP1, FOS, EGF, NOS3, INS, TP53, HRAS, MAPK8, HSP90AA1, IGF1R, JUN, MAPK14, PPARG, ACE, KDR, CALM3, CALM1, CDK2, CALM2, MAPK1, SRC, ESR1, CDC42, IGF1, MDM2, MAPK3, JAK2, EGFR, CASP3, BCL2L1, AR, ALB, CCL5, RAF1, TF, AKT1, IL10, SMAD3, IL6, F2, TGFB1, STAT1, TNF, RHOA, REN, HPGDS, IL1B, SPP1, MAP2K1, LEP, IFNG, KIT
2	11.911	46	268	ESR2, MMP9, PIK3R1, MMP3, MET, IGF2, CD40LG, PTK2, CSK, CSF1, MAP3K1, SP1, PTK2B, NOS2, HMOX1, ADIPOQ, PTPN11, LEPR, RAC2, GSK3B, PRL, GNRH1, XIAP, MMP1, CASP1, CRP, TERT, AKT2, LCK, PGR, INSR, APOE, IGFBP3, IL1A, BMP4, PTPN1, PLAU, SELE, NR3C1, LDLR, IL7, WT1, POMC, RUNX2, MAPK12, SYK
3	4.48	26	56	CD79A, TRAF6, MMP14, TNFSF11, KAT2B, CDK6, HSPA1A, LYZ, BMP2, SOD2, APAF1, ELN, PTH, EIF4E, HSPA8, IL11, CASP7, ZAP70, MMP13, TNFRSF11B, LGALS3, ERBB4, HCK, VDR, CA2, BGLAP
4	4.333	7	13	IGFBP2, IGFBP1, GH1, JAK3, SHBG, GHRL, GHR
5	4	7	12	ERCC2, NBN, CDK7, ERCC6, WRN, MCM9, RFC2
6	3.923	27	51	COL1A2, PAK7, CYP2C9, DKK1, APRT, NT5M, NOG, COL2A1, NR1I2, LAMB3, NQO1, LAMC2, GMPR, GSS, LAMA3, WNT3A, CANT1, GSTA1, CYP2C8, BMP6, PLA2G2A, AKR1B1, PDE5A, GMPR2, RAP2A, IMPDH1, GP1BA
7	3.833	13	23	PROKR2, GNRHR, TACR3, ALPL, PROK2, RXRA, FSHR, CYP11A1, SULT2A1, COL1A1, TAC3, KISS1R, IBSP
8	3.667	25	44	SHMT1, PSPH, MTHFD1, GSTZ1, THRB, ACAT1, RXRB, NHP2, AK1, THRA, NOP10, GART, BCAT2, WRAP53, HSD17B1, TINF2, NR1H2, RARB, RARG, DHODH, RPL11, OAT, CRABP2, PEPD, NDP
9	3.5	21	35	IMPDH2, PCCB, CYP24A1, BAZ1B, GTF2IRD1, DTYMK, AKR1C3, TBL2, CLIP2, CDA, ACADM, CYP27B1, HSD3B2, EPHX2, WARS, SULT2B1, PCCA, STS, AKR1D1, IVD, DCK
10	3.2	11	16	SPARC, SERPINA1, CP, BMP7, TGM2, JAG1, WNT1, NOTCH3, NOTCH2, HAMP, SP7
11	3.091	12	17	GPI, GSTA3, ANTXR2, HSD11B1, GSTM2, HADH, TNFRSF11A, GSTM1, MGP, PKLR, ADH5, LDHB
12	3	3	3	CD36, SELP, ELANE
13	3	5	6	RAB5A, RAB11A, CTSB, CTSS, HSP90AB1
14	3	5	6	EPHB4, KL, FGF23, UMPS, LIMK1
15	3	5	6	GLB1, GALK1, PYGL, HK1, AMY1A
16	3	3	3	AMY1B, AMY1C, SDS
17	3	3	3	HSD17B4, SORD, HPGD
18	3	3	3	GALNS, HEXB, PSAP
19	3	3	3	NMNAT3, BST1, NMNAT1
20	2.941	18	25	SOX9, DUT, GLI2, ATIC, PPIA, TGFB2, NME2, HPRT1, CTSK, FGF8, PTHLH, PPP5C, ACP5, TYMS, CALCA, AIP, CASR, BRAF
21	2.769	14	18	PPP1CC, ALDOA, CFD, PDPK1, AHSG, TGFBR1, TTR, HP, DPP4, G6PC, AURKA, MME, GCK, RND3

Figure 6

Cluster of SWD-Osteoporosis network (A, B, C, D … stand for cluster 1, 2, 3, 4 … pink, blue, purple circle stand for compound target, osteoporosis genes, compound-osteoporosis target, resp.)

Pathway of SWD-osteoporosis network

After pathway enrichment analysis, 17 osteoporosis-related pathways were obtained (Figure 7). The PI3K-Akt signaling pathway once again contained the largest number of genes (38 genes), while the insulin signaling pathway had 27 targets and genes, the MAPK signaling pathway had 26 targets and genes, and the FoxO signaling pathway had 24 genes. These pathways and targets may in the effect of SWD in treating osteoporosis. The details of each pathway are shown in Supplementary Table 7.

Figure 7

Pathway o of SWD-osteoporosis network (pink and blue circles stand for compound-osteoporosis target and compound target, resp.; blue, orange, yellow and green hexagon stand for compounds of Rehmanniae Radix Praeparata, Angelicae Sinensis Radix, Chuanxiong Rhizoma, Paeoniae Radix Alba, resp. Blue triangle stands for common compound of Rehmanniae Radix Praeparata, Angelicae Sinensis Radix, Paeoniae Radix Alba. Yellow triangle stands for common compound of Chuanxiong Rhizoma, Paeoniae Radix Alba. Green triangle stands for common compound of Rehmanniae Radix Praeparata, Angelicae Sinensis Radix. Orange triangles stands for common compounds of Angelicae Sinensis Radix, Chuanxiong Rhizoma. Light line stands for the relationship between compound and target, black line stand for the relationship between pathway and target).

In addition, in this network, the compounds of Rehmanniae Radix Praeparata totally regulate 103 targets (which is the most), while that of Angelicae Sinensis Radix regulate 96 targets. The compounds of Paeoniae Radix Alba regulate 92 targets, and Chuanxiong Rhizoma’s compounds regulate 84 targets. This suggests that Rehmanniae Radix Praeparata plays a major role in SWD, while other herbs assist Rehmanniae Radix Praeparata.

Compound-known target-osteoporosis network analysis

Compound-known target network

This network contains 125 nodes (103 compound targets node and 22 compound nodes) and 250 edges. This network is smaller than the compound-compound target network and has fewer targets. Some of the targets in this network have also appeared in previous networks (Figure 8).

Figure 8

Compound-known target network (blue circles stand for known targets; blue diamonds, orange diamonds, yellow diamonds and green diamonds stand for compounds of Rehmanniae Radix Praeparata, Angelicae Sinensis Radix, Chuanxiong Rhizoma, Paeoniae Radix Alba, resp. Blue hexagon stand for common compound of Rehmanniae Radix Praeparata, Angelicae Sinensis Radix, Paeoniae Radix Alba. Yellow hexagon stands for common compound of Chuanxiong Rhizoma, Paeoniae Radix Alba. Green hexagon stands for common compound of Rehmanniae Radix Praeparata, Angelicae Sinensis Radix. Orange hexagons stand for common compounds of Angelicae Sinensis Radix, Chuanxiong Rhizoma).

SWD known target-osteoporosis network

The SWD known target-osteoporosis network is composed of 440 nodes and 6262 edges (Figure 9). This network is much smaller than the SWD-osteoporosis network, but some targets that appear in Figure 5 can be found in this network, such as SWD’s predicted targets and SWD-osteoporosis targets. There are also some new targets in these areas, and these new targets may be indirectly regulated by SWD’s predicted targets.

Figure 9

SWD known target-osteoporosis network (pink, blue, purple circle stand for known target, osteoporosis genes, compound-osteoporosis target, resp.).

Cluster and pathway of SWD known target-osteoporosis network

In analyzing the network by MCODE, 11 clusters were returned. (Table 3, Figure 10). Most of them are the same as that in the clusters of the SWD-Osteoporosis network, which indirectly confirm SWD’s effects.

Table 3

Cluster of SWD known target-osteoporosis network.

Cluster	Score	Nodes	Edges	Targets and genes
1	46.621	59	1352	NR3C1, KDR, SRC, MAPK14, GSK3B, CDK2, MAPK3, IGF1, CD40LG, TF, IL10, ALB, PTGS2, IL6, TGFB1, NOS3, NOS2, IL1B, AR, SPP1, APOE, PPARG, HSP90AA1, PIK3CG, F2, LEP, IL2, IFNG, CCND1, CALM2, GGT1, CALM3, CALM1, RELA, MMP2, CTNNB1, EGF, AKT1, BCL2, BAX, TNF, JUN, CASP3, FOS, TP53, MAPK8, MMP1, STAT1, INS, HMOX1, ICAM1, SP1, SELE, VCAM1, IGF1R, LEPR, INSR, ESR1, CAT
2	9.789	20	93	SHBG, IGFBP5, ADRA1B, TAC3, AIP, KISS1R, PROKR2, FGFR1, KISS1, CASR, GNRHR, TACR3, GH1, PROK2, CHRM3, IGFBP2, CHRM1, HTR2A, GHR, ADRA1A
3	6.375	17	51	NFATC1, PTH1R, CYP11A1, CTSK, IAPP, PTHLH, FSHR, STAR, ACP5, IBSP, COL1A1, CALCA, CGA, ALPL, CALCR, ADRB1, GATA4
4	6	17	48	ESR2, GHRL, BMP4, BMP2, POMC, ADIPOQ, TNFRSF11B, TRAF6, RUNX2, BGLAP, IKBKB, CRH, IL11, PTH, GNRH1, MMP14, IGFBP1
5	5	5	10	ADRA2A, CHRM2, ADRA2C, ADRA2B, CNR2
6	4.889	19	44	TNFSF11, MAP3K1, TERT, CSF1, LDLR, KIT, IGFBP3, IL1A, IL7, SMAD3, ELN, PRL, SLC2A4, CD36, PLAU, MMP3, IGF2, CYP19A1, CRP
7	3.75	17	30	ERCC6, CYP27B1, CYP27A1, CYP1B1, ERCC2, GTF2I, GSTP1, BAZ1B, GTF2IRD1, TBL2, HSD3B7, CLIP2, CYP24A1, GSTM2, CDC73, AKR1C3, HSD3B2
8	3.6	16	27	HFE, PHGDH, CBS, NOTCH2, SP7, FZD4, LRP5, TOP2B, CP, IGFBP4, BMP6, JAG1, LRP6, WNT3A, ATP7A, DKK1
9	3	3	3	CHEK1, TOP2A, POLD1
10	3	3	3	GABRA3, GABRA1, GABRA2
11	2.667	10	12	FGF8, CDX2, DUSP6, CYP3A4, FGF17, CYP1A1, RXRA, NR5A1, NCOA1, SOX9

Figure 10

Cluster of SWD known target-osteoporosis network (A, B, C, D … stand for cluster 1, 2, 3, 4 … pink, blue, purple circle stand for known target, osteoporosis genes, known-osteoporosis target, resp.).

These clusters were subjected to GO enrichment analysis. Cluster 5, 9, 10 did not return osteoporosis-related biological processes. As an example: Cluster 1 also has many biological processes, including GO: 0000165, GO: 0038128, GO: 0042346, GO: 0043065, GO: 0043123, GO: 0043410, GO: 0046850, GO: 0048010, GO: 0071392, and GO: 0098868. The details of each cluster and biological process are shown in Supplementary Table 8.

We obtained 15 osteoporosis-related pathways after pathway enrichment analysis. They are the same as pathways in Figures 3 and 7, suggesting that these pathways may be the crucial pathways in osteoporosis development (Figure 11). The PI3K-AKT signaling pathway in Figures 7 and 11 contains the most targets and genes; these targets and genes are shown in Figure 12 (red squares) as an example. Some of the red squares in Figure 12B overlap with those in Figure 12A, which further demonstrates the reliability of the predicted network. The details are described in Supplementary Table 9.

Figure 11

Pathway o of SWD known target-osteoporosis network (pink and blue circles stand for compound-osteoporosis target and compound target, resp.; Blue, orange, yellow and green hexagon stand for compounds of Rehmanniae Radix Praeparata, Angelicae Sinensis Radix, Chuanxiong Rhizoma, Paeoniae Radix Alba, resp. Green triangle stands for common compound of Rehmanniae Radix Praeparata, Angelicae Sinensis Radix. Orange triangles stands for common compounds of Angelicae Sinensis Radix, Chuanxiong Rhizoma. Light line stands for the relationship between compound and target, black line stands for the relationship between pathway and target).

Figure 12

PI3K-AKT signaling pathway adapted from KEGG (ID: hsa04151) (A) the predicted targets and genes were signed in red; (B) the known targets and genes were signed in red.

In this compound-known target network, we found some signaling pathways and biological processes that have appeared in previous networks, such as: hsa04380, hsa04210, hsa04910, GO: 0070374, GO: 0045453, and GO: 0001957. Several biological processes are also similar in functions. This somewhat verifies the reliability of the predicted target network. However, because this method failed to directly reflect the impact of compound content on biological processes, more research is required to test the effects of SWD and compounds in relation to osteoporosis and to confirm their pharmacological and molecular mechanisms.

Discussion

Three main networks (Osteoporosis genes PPI network, SWD-osteoporosis PPI network, and SWD known target- osteoporosis PPI network) were constructed and analyzed in this research. The first network uncovered the possible mechanism of osteoporosis. The second network explored the potential mechanism of SWD treating osteoporosis. The third network confirmed the feasibility of the second network.

According to the network analysis, the compounds with high degree (with many connected predicted targets) are thought to be the core compounds of an herb; the pathways containing many targets linked by core compounds are considered to be the ones regulated by those core compounds. Based on this, the core compounds of Rehmanniae Radix Praeparata are acteoside, benzoic acid, and 5-hydroxymethylfurfural. The top 3 pathways regulated by acteoside are insulin signaling pathway, PI3K-Akt signaling pathway, and osteoclast differentiation. The top 3 pathways controlled by benzoic acid are PI3K-Akt signaling pathway, FoxO signaling pathway, and insulin signaling pathway. The top 3 pathways regulated by 5-hydroxymethylfurfural are insulin signaling pathway, PI3K-Akt signaling pathway, and osteoclast differentiation. The core compounds of Angelicae Sinensis Radix are sitosterol, ferulic acid, and senkyunolide I. The top 3 pathways regulated by sitosterol are PI3K-Akt signaling pathway, Estrogen signaling pathway, and insulin signaling pathway. The top 3 pathways regulated by ferulic acid are PI3K-Akt signaling pathway, FoxO signaling pathway, and insulin signaling pathway. The top 3 pathways regulated by senkyunolide I are PI3K-Akt signaling pathway, FoxO signaling pathway, and insulin signaling pathway. The core compounds of Chuanxiong Rhizoma are senkyunone, mandenol, and beta-sitosterol. The top 3 pathways regulated by senkyunone are PI3K-Akt signaling pathway, Insulin signaling pathway, and Estrogen signaling pathway. The top 3 pathways regulated by mandenol are PI3K-Akt signaling pathway, Insulin signaling pathway, and cAMP signaling pathway. The top 3 pathways regulated by beta-sitosterol are PI3K-Akt signaling pathway, Estrogen signaling pathway, and FoxO signaling pathway. The core compounds of Paeoniae Radix Alba are mairin, MOL001910, and paeoniflorgenone. The top 4 pathways regulated by mairin are PI3K-Akt signaling pathway, insulin signaling pathway, cAMP signaling pathway, and insulin resistance. The top 3 pathways regulated by MOL001910 are PI3K-Akt signaling pathway, Estrogen signaling pathway, and insulin signaling pathway. The top 3 pathways regulated by paeoniflorgenone are PI3K-Akt signaling pathway, FoxO signaling pathway, and insulin signaling pathway.

The development of osteoporosis is the result of a combination of multiple factors. At the molecular level, the process of bone formation and bone remodeling involves a variety of signaling pathways, and they interact with each other to play a fine regulatory role in complex regulatory network systems. Studies have shown that Wnt/β-catenin signaling pathway [36], BMP-2 signaling pathway [37], and OPG/RANKL signaling pathway [38] play important regulatory roles in bone formation and bone remodeling [39]. The BMP-2 signaling pathway process has 2 roles in the cell. One is to transfer the external growth-promoting signals to the nucleus via Smads to promote osteogenic differentiation [37]. The other is the MAPK signaling pathway, which also includes 3 signaling pathways: the extracellular signal-regulated kinase (ERK) signaling pathway, the c-Jun N-terminal kinase (JNK) signaling pathway, and the p38 signaling pathway [40-42]. These BMP-2 signaling pathways act through phosphorylation, which in turn regulates transcription of downstream target genes such as RUNX2 and Osterix, so as to promote osteogenesis [40-42]. In summary, existing studies have found that the interaction of multiple signaling pathways mediates the development of osteoporosis. Changes in the expression of these signaling pathways may reduce osteoblast formation, proliferation, and differentiation, ultimately leading to osteoporosis.

From a cytological point of view, osteoporosis is closely related to abnormal osteoclastogenesis [43]. Osteoclasts are derived from the monocyte-macrophage lineage and play an important role in normal development and in some diseases [44]. Osteoclastogenesis in vivo is controlled by various cytokines and signaling pathways [45], wherein macrophage colony-stimulating factor (M-CSF) and the receptor activator NF-κB ligand (RANKL) are the 2 most important factors for osteoclastogenesis [45]. Studies have shown that M-CSF is involved in the survival of osteoclasts and also indirectly affects the binding of RANKL to RANK, which together promote osteoclastogenesis [46]. RANKL binds to RANK to activate TNF receptor-associated factor 2/6 (TRAF 2/6) and transforming growth factor β-activated kinase 1 (TAK1), which in turn promotes IκB kinase (IKK) complex, c-Jun N-terminal kinase (JNK), extracellular signal-regulated kinase (Erk), and p38 phosphorylation, and thereby activating NF-κB and MAPK signaling to promote osteoclastogenesis [47,48]. We also found some related signaling pathways and biological processes in this disease network enrichment analysis, such as hsa04010, hsa04310, hsa04668, GO: 0070374, GO: 0000165, and GO: 0030514.

Because the OPG/RANKL/RANK signaling pathway, Wnt/β-catenin signaling pathway, and BMPs signaling pathway play an important role in regulating osteogenesis and osteoclast lineage function, as well as in regulating bone mass and bone remodeling, some protein molecules in these signaling pathways have become new targets for osteoporosis treatment. Targeted therapeutics for these targets have also entered clinical trials, but they are generally more expensive and have potential adverse effects [49-54]. This research shows that the active components of SWD can directly (purple circles in Figure 5) or indirectly (through pink circles to blue circles in Figure 5) regulate osteoporosis-related targets, thereby regulating “(hsa04010) MAPK signaling pathway”, “(hsa04668) TNF signaling pathway”, and “(hsa04064) NF-kappa B signaling pathway”, as well as “(GO: 0071773) cellular response to BMP stimulus”, “(GO: 0043410) positive regulation of MAPK cascade”, “(GO: 0060070) canonical Wnt signaling pathway”, and “(GO: 0046330) positive regulation of JNK cascade”.

Macroscopically, osteoporosis is a multi-pathogenic disease whose etiology is related to age, endocrine disorders, calcium absorption disorders, limb disuse, immunity, nutrition, genetics, and long-term use of certain drugs.

For endocrine therapy, estrogen is a bone resorption inhibitor that directly regulates bone metabolism, reduces bone resorption after menopause, relieves bone and joint pain, and reduces the incidence of fractures [49]. Therefore, selective estrogen receptor modulators play an important role as a classic drug for the treatment of osteoporosis. Our study found that SWD can regulate estrogen-related biological processes and signaling pathways to exert anti-osteoporosis effects, such as “(hsa04915) Estrogen signaling pathway”, “(GO: 0032355) response to estradiol”, and “(GO: 0043627) response to estrogen”. Liu et al. used microarray transcriptional analysis to demonstrate the discovery of phytoestrogens in SWD [15]. Wen et al. also identified SWD as Nrf2 activator and phytoestrogens [14]. Lu et al. found that SWD drug-containing serum showed an estrogen-like effect, and in vitro experiments further confirmed that SWD drug-containing serum can activate ER signaling pathways [13]. The long-term use of estrogen is likely to lead to an increased risk of estrogen-related tumors and cardiovascular adverse events. Meanwhile, because of its two-way regulation, phytoestrogens can alleviate the symptoms of low estrogen after menopause, while avoiding the risk of reproductive system tumors caused by medical treatment of estrogen replacement therapy. The estrogen-like effects of phytoestrogens affect hormone secretion, metabolic biological activity, protein synthesis, and growth factor activity. It has different functions under different physiological conditions. For example, in the case of strong estrogen physiological activity, phytoestrogens can act as antiestrogens and reduce the risk of estrogen-regulated cancers such as breast cancer. Recent studies in breast cancer have also shown a 21% reduction in mortality among patients with low intake (≥0.3 mg/d) vs. the high-intake breast cancer group (≥1.5 mg/d). This suggests that phytoestrogens are safe even in the treatment of high-risk breast cancer [55]. When estrogen has low physiological activity (such as after female menopause), phytoestrogens can enhance estrogen and reduce various symptoms caused by low estrogen (especially osteoporosis). Recent research has also shown that phytoestrogens can better prevent osteoporosis while reducing other risks associated with exogenous estrogen [56,57]. Therefore, SWD has an estrogen-like effect and may be a potential alternative treatment option [49].

In vitamin D and calcium supplementation, vitamin D is closely related to calcium and phosphorus metabolism. Vitamin D is an important hormone regulating bone metabolism in the body. It can affect bone formation through indirect and direct effects on osteoblasts [58,59]. Therefore, vitamin D plays an important role in the prevention of osteoporosis [58,59]. The active metabolite 1a,25-dihydroxyvitamin D [1,25(OH)2D] of vitamin D in the body can act on bone formation indirectly by promoting intestinal calcium absorption [60]. In addition, osteoblasts express vitamin D receptor (VDR), which is also an important target cell for 1,25(OH)2D action. In vitro osteoblast culture showed that physiological amounts of 1,25(OH)2D can stimulate osteoblast differentiation and promote bone formation [60,61]. Bone calcium stores calcium and plays a very important role in maintaining blood calcium stability. When the blood calcium concentration increases or decreases, the body will mobilize a variety of mechanisms to deal with it by adjusting the “inventory” to ensure normal blood calcium. When the blood calcium concentration is increased, calcium is sent to the bone for storage. This process is called bone formation. When the blood calcium concentration is lowered, calcium is dissolved from the bone to supplement the deficiency of blood calcium. This process is called bone resorption [62,63]. These 2 effects are regulated by hormones such as parathyroid hormone, vitamin D, and calcitonin. When people reach old age, various hormone secretions are reduced, bone metabolism is disordered, and the absorbed calcium cannot compensate for excreted calcium, and the body is in a state of negative calcium balance [64]. Hence, calcium preparations and vitamin D are often used as basic medicines for the treatment of osteoporosis. This research shows that SWD slows the progression of osteoporosis by regulating the biological processes of calcium metabolism, such as “(GO: 0035630) bone mineralization involved in bone maturation”, and “(GO: 0030278) regulation of ossification”, as well as vitamin D biological processes “(GO: 0060557) positive regulation of vitamin D biosynthetic process” and “(GO: 0042359) vitamin D metabolic process”. In terms of bone formation, Wu et al. showed that NF-κB activation contributes to SWD extract-induced bone mineralization and expression of ALP, BMP-2 and OPN in osteoblasts in vitro [14]. In addition, there are also some biological processes directly related to the development of osteoporosis in Supplementary Tables 4 and 6, such as “(GO: 0030316) osteoclast differentiation”, “(GO: 0045669) positive regulation of osteoblast differentiation”, and “(GO: 0045453) bone resorption”. This suggests that SWD may directly regulate these biological processes to achieve therapeutic effects.

Last but not least, because this method failed to directly reflect the impact of compound content on targets, pathways, and biological processes, further experimental research is required to test the effects of SWD and compounds in relation to osteoporosis and to confirm their pharmacological and molecular mechanisms. However, the present study did not produce negative results. Instead, it suggests that systemic pharmacology can save time and effort by predicting targets, signaling pathways, and biological processes to provide direction for subsequent in-depth experimental research.

Conclusions

The therapeutic effect of SWD on osteoporosis may be achieved by interfering with the above-mentioned targets (such as HSP90AB1, FGFR1, HRAS, GRB2, and PGF), biological processes, and signaling pathways (such as PI3K-Akt signaling pathway, insulin signaling pathway, MAPK signaling pathway, and FoxO signaling pathway) related to the development of osteoporosis.

Supplementary Data

Supplementary Tables 1–9 available from the corresponding author on request.