A network pharmacology approach to explore and validate the potential targets underlying the effect of Ginsenoside on Osteoporosis

Objective This study aim to investigate the potential targets involving the effect of ginsenoside on osteoporosis using a network pharmacology approach. Methods Ginsenoside and its drug targets associated to osteoporosis (OP) were identied by using network analysis. First, the Traditional Chinese Medicine Systems Pharmacology Database and Analysis Platform (TCMSP), DrugBank database, Pharmmapper database and Cytoscape software were used to mine information relevant to Ginsenoside ingredients and Ginsenoside -related targets. Second, the Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis of ginsenoside-target gene and ginsenoside-OP target gene were performed in String, Phenopedia, DisGeNET database and Metascape software. Eventually, the protein to protein interaction (PPI) of key ginsenoside-OP targets were created by String database and Cytoscape software. The validation of the binding of ginsenoside to target proteins was plotted by utilizing SwissDock tool, UCSF Chimera and Pymol software. A total of 8 important active ingredients of ginsenosides were obtained in the TCMSP. Eighty potential targets of ginsenoside and 1304 related targets involved in OP were subjected to network analysis, and the 17 intersection targets were indicated to be linked to ginsenoside treating OP. GO and KEGG analysis showed the top 10 items of biological processes, cellular components, molecular functions and signaling pathways in the 80 targets of ginsenoside. Then, 14 key targets were determined to be the most crucial genes by protein to protein interaction (PPI) analysis. In the 14 intersection potential targets, 10 signaling pathways were dened by Metascape software. Validation plots of four target proteins IL1B, TNF, IFNG, NFKBIA binding to ginsenoside rh2 was drew, lastly. This study investigated the potential targets and signaling pathways of ginsenoside during the treatment of OP, which might be benecial to elucidate the mechanism concerned to the action of ginsenoside and might supply a better understanding of its anti-OP effects. Distribution, Metabolism, Excretion.


Introduction
Osteoporosis (OP) is determined as a chronic systemic bone disorder, characterised by low bone mass and microarchitectural degradation of bone tissue, thereby increasing the vulnerability of bone and susceptibility to fracture [1]. A recent report derived from the International Osteoporosis Foundation represents that approximately 4 million men and 16 million women are affected by OP in the six countries of the European Union [2]. Hip fractures and vertebral fractures are the most common fracture associated with bone fragility, which have led to extensive morbidity, mortality, and other adverse health consequences. It was estimated that approximately 2.7 million hip fractures happened due to OP in 2010 worldwide [1]. Pharmacological therapies, including Bisphosphonates, Denosumab, Raloxifene, Teriparatide and Abaloparatide, have been applied to reduce the fracture risk in individual osteoporotic patient [1,3]. However, these medications also brought about some adverse events, such as upper gastrointestinal adverse reactions, skin rash, breast pain and muscle cramps and so on [1]. Recently, scientists have found that traditional Chinese medicines (TCM) have great potential in preventing and treating OP on account of good curative effect and minimal side effects. Cistanche tubulosa [4], Muscone [5], Xian-Ling-Gu-Bao prescription [6], Morinda o cinalis have been explored to be bene cial for osteoporosis with different mechanisms. Although these TCMs have been examined to be potential treatment e cacy for OP, the best one is not yet found.

Page 3/24
As the primary pharmacologically active ingredients of the ginseng root, ginsenosides are a type of sterol compound including a structure with four-ring and a carbohydrate steroidal body [7]. To date, the researchers have isolated and identi ed a total of forty structurally diverging ginsenosides from the P. ginseng root (Table 1). Several studies have explored the anti-osteoporosis activity of ginsenosides in vivo, in vitro and rat [8][9][10]. However, although some progress has been achieved on the research of OP with ginsenosides, we are still lack of systematic, in-depth and comprehensive understanding on the intrinsic interaction between potential targets and signaling pathways involved this. Network pharmacology is a scienti c research method that has been extensively applied in drug study based on its strong, e cient and holistic function. The interconnectedness between "drugstargets-pathways-diseases" could be adequately illustrated according to bioinformatic algorithms, molecular biology and related database. Here, the potential molecular targets and signaling pathways on ginsenosides involving the OP treatment will be fully characterized, and the potential related mechanism will also be clearly described by establishing a network of ginsenosides and OP target interrelation.  [11] based on the keyword: renshen (ginseng). Absorption, Distribution, Metabolism, Excretion (ADME) [12] properties of each active component of medicinal herbs could be derived from the TCMSP database, which contains related information, such as Oral Bioavailability (OB), Drug-Likeness (DL), Half-Life (HL) and so on.
Firstly, we gained ginsenosides from all active substances results of ginseng by TCMSP. Then, the ingredients of ginsenosides with OB ≥ 30% and DL ≥ 0.7 were selected as main research subjects in this study. We also acquired the three-dimension (2D) chemical structures of these selected ginsenosides from the TCMSP database in mol2 format.

Prediction of potential target proteins for ginsenosides
The TCMSP database has given some potential target proteins of each active component. Some targets results are supported by DrugBank database (https://www.drugbank.ca/) [13]. We rstly stored these results as part of potential targets of selected ginsenosides. Secondly, Pharmmapper database (http://www.lilab-ecust.cn/pharmmapper/) were used to identify other potential human target proteins for the ginsenosides. Pharmmapper could identify potential target candidates for a given small molecule structure by a pharmacophore mapping approach. 3D chemical structure of each selected ginsenoside was uploaded to the Pharmmapper with default parameters of database. 300 potential target candidate results of each ginsenoside were downloaded after handling [14].
Thirdly, we chose human proteins from Pharmmapper results as target candidates. Together with results of TCMSP in the rst step, we got the potential target proteins for selected ginsenosides. Then we used Cytoscape software to visualize the relationship between selected ginsenosides and potential targets [15].
2.3 Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis for target proteins For the target proteins of selected ginsenosides we obtained above, we put to use the String database ( https://string-db.org/) to perform GO (Biological Process, Cellular Component and Molecular Function) and KEGG (Kyoto Encyclopedia of Genes and Genomes) enrichment analysis by default analysis parameters [16]. For the analysis results harvested from String, we made use of R language (version 4.0.0) to visualize them in the form of bubble charts. In order to ensure the visualization effect, we only visualized the TOP 10 items of GO and KEGG enrichment pathway. For pathways related to osteogenesis, we further applied the R package 'pathview' to visualize their pathways and related targets [17].

Screening the potential target genes related to OP
Phenopedia database (http://www.hugenavigator.net/HuGENavigator/startPagePhenoPedia.do) is an application based on web page with persistently renewed from PubMed, which is able to explore the literatures associated with human genes and diseases [18]. In this analysis, the keyword "Osteoporosis" was rstly searched in the Phenopedia database to obtain the target genes related to OP. DisGeNET database (http://www.disgenet.org/) is a standardized intellectual management platform based on multiple scienti c literatures, which covers a full range of human diseases and related genes [19]. The target genes associated with OP were also explored in this database. Then, we took the union results gained from both Phenopedia and DisGeNET databases as the potential target genes of OP. The jvenn tool [20] was used to cross-validate the ginsenoside target proteins and OP-related potential target genes to procure a possible target genes of ginsenoside treating OP.

Functional enrichment analysis for target genes of ginsenoside treating OP by Metascape
Metascape is a web-based portal characterized by functional enrichment, interactome analysis, gene annotation, and membership search to leverage over 40 independent knowledgebases. During the process of enrichment analysis, the input gene list is compared to thousands of gene sets determined by their involvement in particular biological processes, pathway membership, protein localization, enzymatic function, or other features [21]. In order to further analyze the biological function and pathway of ginsenoside target genes treating OP, we utilized Metascape database and Minimal Common Oncology Data Elements (MCODE) algorithm to perform GO(Biological Process)and KEGG enrichment analysis for key target genes [21]. The completed operation methods could be found in this published paper written by Zhou [21].

Protein-protein interaction (PPI) analysis for key target proteins
In order to study the interaction between key target proteins, we uploaded them to String database for further analysis, and the analysis parameters were set as default (combined score > 0.4, medium con dence). The result of PPI given by String was visualized by Cytoscape software. Then we made use of the MCODE plug-in of Cytoscape to analyze the key sub-network of the obtained PPI network to study the important topological sub-network structure in the entire network [22]. The analysis parameters of MCODE were set to degree cur off = 2, node score cut off = 0.2, Max depth = 100, k score = 2 and MCODE score > 5.

Prediction of the binding of ginsenoside to target proteins
In order to elucidate how the ginsenosides bind to target proteins to exert the therapeutic e cacy, we used the SwissDock tool to make docking predictions. We obtained the Protein Data Bank (PDB) identi cations (IDs) [23] of key target proteins from the PDB database. With the 3D-structure in mol2 format of selected ginsenoside, target proteins were uploaded to SwissDock website to make docking predictions [24]. The model with minimal G scores was considered as the most excellent docking model. For the optimal prediction results we got, we took advantage of UCSF Chimera and Pymol software to visualize them [25].

Pharmacologically active ingredients of ginsenosides and 2D chemical structure
We obtained 190 TCM active ingredients when we searched in TCMSP database on account of keyword 'renshen (ginseng)'.
According to lter threshold, a total of 8 important active ingredients of ginsenosides were obtained, which were listed in Table 2.
Then, the 2D chemical structures of these 8 ginsenoside active components were acquired from TCMSP. They are showed in

Prediction of target proteins in the 8 ginsenoside active components
In the TCMSP database, only ginsenoside rh2 and Ginsenoside-Rh4_qt got the predicted target proteins. Therefore, we explored the target proteins of these 8 ginsenoside active components in the Pharmmapper database. A total of 66 and 14 target proteins were secured in the Pharmmaper and TCMSP database, respectively. Finally, a network with 80 target proteins and 8 ginsenosides was constructed below (Fig. 3). The larger the red circle, the more connections. As we can see, ginsenoside rh2 contains the most connections.

GO and KEGG pathway enrichment results of target proteins of ginsenosides
The GO and KEGG pathway enrichment analysis of target proteins in these 8 ginsenosides were performed in the String database.
The results were showed as bubble chart in the Fig. 3. As we can see, the target proteins mainly involve cellular response to stimulus, organic substance and lipid in GO biological process (BP). In GO molecular function (MF), they are mainly involving in all kinds of binding, such as protein biding, ion binding, heterocyclic compound binding and so on. And, in GO cellular component (CC), intracellular organelle and membrane-bounded organelle are the most important binding part for the target proteins. Then, in the KEGG enrichment analysis, the tumor and in ammation related biological pathways are the primary pathways concerned to our target proteins, except osteoclast differentiation pathway ( Figure 3). Finally, we found ve target proteins participated in the osteoclast differentiation pathway, including BTK, IFNG, IL1B, NFKBIA and TNF. The completed pathway chart was showed in the Fig. 4.

The potential target of ginsenoside treating OP
The cross-validation was performed between OP-related genes and ginsenoside target protein. Then, we took the cross-validated intersection as the possible targets for ginsenoside in the treatment of OP. As seen in the Fig. 6, a total of 1304 OP-related genes were obtained, and a total of 17 overlapping genes were achieved after cross-validation with 80 ginsenoside target proteins.

Functional and pathway enrichment analysis in the 17 key target proteins
Metascape database was used to perform functional and pathway enrichment analysis for these 17 key target genes. The bars in the Fig. 7a are represented by discrete color and p-values of statistical signi cance. The three most signi cantly enriched ontology terms (M, GO and hsa) are combined to visualize and interpret.
As showed in the Fig. 7a, the most important pathways involved in these targeted genes was PID RXR VDR pathway, which was a pathway concerning with vitamin D receptor signaling mechanisms and bone synthesis process. The other pathway are PID lysophospholipid pathway and PID REG GR pathway. There were 5 items concerned to GO analysis in this gure, which are including positive regulation of chemokine biosynthetic process, positive regulation of cell death, gland development, organic anion transport and blood circulation. Seven pathways involved in the KEGG enrichment analysis, such as legionellpsis, pathway in cancer, Th17 cell differentiation, uid shear stress and atherosclerosis, amyotrophic lateral sclerosis, transcriptional regulation of white adipocyte differentiation and SUMO E3 ligases SUMOylate target proteins.
A mature complicated algorithm called MCODE is used by Metascape to automatically extract target complexes embedded in the large network [21]. It provides dense interactome neighborhoods and the interpretation of these biological enrichment pathways ( Figure. 7b and c).

PPI network of the 17 key targets
The 17 key target proteins were used to perform PPI analysis in the String database. Then, the PPI network obtained above was further analyzed in the plug-in MCODE of Cytoscape. There were16 nodes in this PPI network, beside HBB (Fig. 8A). Nevertheless, in the PPI sub-network, the important nodes (genes) became 14 ones (Fig. 8B), including IL1B, MMP2, TNF, IFNG, RXRG, RXRB, RARB, THRB, CASP1, CASP3, RARG, HMOX1, NFKBIA and NCOA2, beside BAX, HBB and SLC2A4. Therefore, the most important key targets were reduced to 14 ones.

Prediction of ginsenosides binding to target proteins associated with OP
The SwissDock database was utilized to predict the binding of ginsenosides and target proteins related to OP. We chose one as the representative drug from the 8 main ginsenoside active ingredient with the largest MW, the longest drug HL, and the most target proteins, which was ginsenoside rh2(MOL005344)and was currently the main ginsenoside drug on the market. Prediction analysis of four target proteins IL1B, TNF, IFNG, NFKBIA binding to ginsenoside rh2 was performed. The results with optimal docking were showed in Fig. 9 and Fig. 10 below.

Discussion
Plenty of osteoporosis-associated morbidities have been on the rise in women older than 55-yrs and men older than 65-yrs. OP has become a common disorder that has greater potential to threaten the elderly health. Many previous literatures have reported that oxidative stress response, de ciency of estrogen, bone immune disorder, chronic in ammative reaction, bone marrow microcirculation dysfunction and glucocorticoid-induced bone loss may be contributed to the pathogenesis of OP [26,27].
Network pharmacology has been broadly used to study drug and relevant disease and made full sense of drug research due to big data analysis [28]. In the clinical practice, ginsenosides have been widely utilized to treat different kinds of disease because of good effectiveness and low side effects. In this study, the network pharmacology was conducted to explore the potential key targets and mechanisms associated with OP treatment by building and analyzing drug target networks, enrichment analysis.

Page 9/24
4.1. The potential key targets Fourteen potential key targets were obtained by cross-validation of intersection and PPI network analysis, which are IL1B, MMP2,  TNF, IFNG, RXRG, RXRB, RARB, THRB, CASP1, CASP3, RARG, HMOX1, NFKBIA and NCOA2. All of these target genes or proteins were considered to be associated with ginsenosides treating OP, which have been reported in the recent literatures. IL1B was detected as an increased in ammation-relevant gene of OP in Kashin-beck disease patients, and was regulated by ginsenoside rhs in treating idiopathic pulmonary brosis [29,30]. MMP2 has been showed that might be a new possible therapeutical target in precaution and treatment of glucocorticoid-related OP [31]. TNF-α and CASP3 expression were regulated to alleviate OP by Sinomenii Caulis extract [32]. A recent Cis-eQTL analysis have elucidated that IFNG expression was positively related to OP and abdominal obesity [33]. RARB was examined to be a bone mineral density loci by genome-wide association studies in Koreans [34]. High glucose induced caspase-1 and IL-1B activation and inhibition of efferocytosis in osteoclast-mediated diabetic osteoporosis [35]. Fufang Lurong Jiangu Capsule signi cantly alleviating oxidative stress in osteoblasts induced by Nrf2/HMOX-1axis might provide new perception to the treatment of osteoporosis [36]. And, ginsenoside rb3 defends cardiomyocytes against oxidative injury via activation of PERK/Nrf2/HMOX1 axis [37]. Acteoside signi cantly in uenced osteoclastogenesis by downregulated NFKBIA expression to prevent the development of osteoporosis [38]. Most of these genes have been linked to the treatment of mechanism of OP. However, RXRG, RXRB, THRB, RARG has not been demonstrated to be associated with OP or ginsenosides in the previous studies. This implicates that these genes might be new potential targets for ginsenosides treating OP worth being explored in the future.

The Functional and pathway enrichment analysis
As shown in the Figure. 4, the enrichment analysis built by String can provide a great number of enrichment results regarding to GO and KEGG. Nevertheless, it inevitably contains some redundant information and reduces the accuracy of the outcome.
Metascape adopts a dynamic clustering procedure that could reduce the redundancy and raise the precision. Compared to Figure. 4, the Figure. 7 produced by Metascape represents much simpli ed and precise results. The pathways included in the Figure. 7 contain three most signi cantly enriched ontology terms: M,GO and hsa.
The kind of M enrichment pathway is uniquely generated by Metascape, including PID RXR VDR pathway, PID lysophospholipid pathway and PID REG GR pathway. It was reported that bone homeostasis was mediated by RXR (Retinoid X Receptor) signaling pathway and RXRs might be potential target for the treatment of OP [39]. Vitamin D3 plays a key role in human health. There is association between vitamin D3 de ciency and rising risk in OP, cancer. The presence of both vitamin D receptor (VDR) and RXR augmented the promoter activity of SLC30A10. This procedure might have potential bene t for the treatment of OP and Parkinson's disease [40]. Another report showed that 1,25(OH) D posed the potential to prevent age-relevant OP by up-regulation of Ezh2 via the VDR-regulated transcription [41]. An investigation from China evaluated the expression of both VDR and osteoprotegerin in postmenopausal women and found that these two genes were associated with increasing risk of OP [42]. All these evidences indicated that the potential targets were closely associated the PID VDR RXR pathway. Sphingosine-1-phosphate (S1P) is a biologically active lysophospholipid that is enriched in blood and involved in the bone balance [43]. A recent review has demonstrated that S1P linked to osteoclast-osteoblast coupling signals and may play a potential treatment role in OP [44]. This showed PID lysophospholipid pathway involved in the process of OP. PID_REG_GR_PATHWAY is a glucocorticoid receptor (GR) regulatory network related pathway. As is known to all, glucocorticoid is a primary factor to induce osteoporosis. Glucocorticoid receptor have been found to enhance Dexamethasone-inhibited osteogenic differentiation by silencing of let-7f-5p and agomiRlet-7f-5p signi cantly reversed OP induced by Dexamethasone. However, combining these pathways and ginsenosides, much work still needs to enrich this research eld.
There are 5 GO related pathways listed in the gure. The chemokines have been reported to play an important role in the physiological and pathological bone remodeling, which is the main pathological of OP [45]. It has been determined that CX3CL1/CX3CR1 axis is a potential target of immune treatment in OP and might be one of the most hopeful targets in OP therapy [46]. All these evidences represent positive regulation of chemokine biosynthetic process to be important in ginsenosides treating OP. Cell death is a signi cant biological process involved in OP. Many published literature have elucidated that cell death occurring from a variety of causes, such as iron-induced osteoblast cell death [47], MDM-p53-induced cell death [48], is also associated with the formation of OP.
It is so surprising that legionellpsis is listed as the rst item of KEGG pathway. In fact, we do not yet nd any legionellpsis-related information involved in OP. This may be a prediction direction about the mechanism of combining ginsenoside treating OP.
Pathway in cancer shared many common signaling pathway with OP, such as Wnt signaling pathway [49], NF-κB pathway [50], PI3K-AKT-NFATc1 pathway [51] and MAPK pathway [52], etc. Sclerostin promoted Th17 cells differentiation and reduced Treg cells differentiation, which aggravate the sclerostin induced inhibition of osteogenesis [53]. Fluid shear stress that can quickly increase release of calcium from endoplasmic reticulum to induce calcium transients, plays a crucial role in osteoblast balance [54]. The metals with neurotoxic properties, amyloid beta peptide and amyloid precursor protein found in brain and cerebrospinal uid from patients have been examined in bone tissue from patients with OP [55]. Sphingosine-1-phosphate (S1P) signaling impacts bone metabolism by markedly increasing bone formation, mass and strength and substantially reducing white adipose tissue. The medication based on S1P has been considered as a hopeful therapeutic agents for the treatment of OP [56]., The loss of SENP3 forti es interferon regulatory factor 8 (IRF8) SUMO3 modi cation at the site of K310 amino acid, which increasing osteoclastogenesis. In conclusion, IRF8 de-SUMO modi cation mediated by SENP3 inhibits osteoclast differentiation and proposes strategies to treat bone loss diseases, such as OP [57]. All of these related pathway indicated that ginsenosides treating OP involve in complicated mechanism, which needs to pay great attention to enrich and demonstrate it.

Prediction of ginsenosides binding to target proteins
The results shown in the Figure.9 and 10 represent the detail 3D optimal docking images, which indicates that ginsenosides can bind with different targets and played a role in anti-osteoporosis. What the most valuable in these outcomes is that ginsenoside not only could combine with their targets to play a role, but also present the precise docking site of the targets. This computational method has been demonstrated to play a key role in the development of drug discovery/design. Furthermore, the advantage of this molecular docking is to reduce time and save money, lastly, accelerate the identi cation of potential drug candidate [58].

The value of our study
We apply the network pharmacology approach to illustrate the molecular mechanisms of ginsenosides treating OP in this study.
The network pharmacology suggests great advantages in exploring the types of drugs and mode of action. Additionally, let me emphasize the limitations of this study. Firstly, our research is based on pharmacological network methods to study potential targets of ginsenosides in treating OP, lacking cell and animal experiments to support. Secondly, because of the di culty of obtaining data, the database might not contain all known or unknown key targets and protein-protein interactions. This could be improved in the future if more data becomes available.

Conclusion
In the present study, 14 potential targets of ginsenosides treatment for OP were determined by a network pharmacological method, which was validated by functional enrichment analysis. The mechanisms used by ginsenosides to treat OP are elucidated. According to the results of the enrichment analysis, we found that 10 signaling pathway were involved in ginsenosides treating OP. Our study laid a good theoretical foundation for further research of the ginsenosides treatment of OP. Although our research has successful outcomes, there are still some limitations in this study. Furthermore, further improvement will need to do in the future. This method only used a computational procedure to predict and analyze potential therapeutic targets and disease related pathways, and further clinical and experimental validation is required.  Flowchart of the systematic procedure to investigate the mechanisms of ginsenosides in the treatment of OP.

Figure 2
The 2D chemical structures of the 8 ginsenoside active components.
Page 17/24  The ginsenoside target proteins (red) in the osteoclast differentiation pathway.

Figure 6
The venny diagram of ginsenosides and OP intersection target genes. Note: 1287 OP non-intersection potential targets (left), 17 OP-ginsenoside intersection potential targets (middle) and 63 ginsenoside non-intersection potential targets (right).  Visualizations of functional enrichment and interactome analysis outcomes for the 17 key targets associated to ginsenosides and OP.

Figure 10
Prediction of ginsenoside rh2 binding to IFNG and NFKBIA associated with ΔG = -7.55 kcal/mol, up plots; ΔG = -8.36 kcal/mol, down plots. Note: The large structures is ginsenoside rhs, the small is target protein. The left plot is made by Chimera software, the middle and right plots are drew by Pymol software.

Figure 11
Molecular docked structure of IFNG-ginsenoside complex. Docking was performed using AutoDock software (V4.2.6) and various amino acid interactions are mentioned in the respective passage. The larger helical lamellar structure is IFNG protein and the smaller hexagon is ginsenoside.

Figure 12
Molecular docked structure of IL1B-ginsenoside complex. Docking was performed using AutoDock software (V4.2.6) and various amino acid interactions are mentioned in the respective passage. The larger helical lamellar structure is IL1B protein and the smaller hexagon is ginsenoside.